%I #22 Feb 16 2021 02:10:05
%S 2,3,7,43,79,547,3319,6163,36979,42667,258847,1553119,1573207,1834639,
%T 1854763,11131927,20224159,20451679,124027567,141569107,141588763,
%U 467477683,1840398379,3278780359,5276533183,6089163523,6155955079,11168428363,11185512199,31655671459
%N Primes p such that p-1 is squarefree and all prime divisors of p-1 other than 13 are also in the sequence.
%C Is this sequence infinite?
%H Andrew Howroyd, <a href="/A267505/b267505.txt">Table of n, a(n) for n = 1..500</a>
%t fa = FactorInteger; is[2, p_] = True; is[2, p_];
%t is[n_, p_] := PrimeQ[n] && MoebiusMu[n - 1] ≠ 0 && Union@Table[is[fa[n - 1][[i, 1]], p] || fa[n - 1][[ i, 1]] == p , {i, Length[fa[n - 1]]}] == {True}; Select[Prime[Range[100000]], is[#, 13] &]
%o (PARI)
%o leastdiv(v, pred, inf)={ \\ finds least divisor d satisfying pred(d) && d>=inf
%o my(recurse(k,d,lim)= if(d >= lim, lim, if(d>=inf && pred(d), d, k++; if(k<=#v, lim=self()(k, d*v[k], lim); self()(k, d, lim), lim))));
%o my(stop=vecprod(v), lim=inf, m=4);
%o while(lim<=stop, lim*=m; my(d=recurse(0,1,lim)); if(d<lim, return(d), m*=4*sqrtint(m))); oo;
%o }
%o lista(n, S=[13])={my(t=2); print1(t, ", "); for(i=2, n, S=concat(S, [t]); t=leastdiv(S, d->isprime(d+1), S[#S]); if(t==oo, break); t++; print1(t, ", "))} \\ _Andrew Howroyd_, Nov 13 2018
%Y Cf. A267503, A267504, A267506, A267507, A005117, A227455, A227007, A227006.
%K nonn
%O 1,1
%A _José María Grau Ribas_, Jan 16 2016
%E Terms a(16) and beyond from _Andrew Howroyd_, Nov 13 2018
|