%I #12 Sep 08 2022 08:46:19
%S 2,3,7,19,37,47,83,89,107,137,181,251,257,349,379,569,631,653,677,691,
%T 797,823,839,863,883,919,1009,1021,1223,1229,1361,1423,1571,1609,1831,
%U 1873,1907,1993,2053,2113,2207,2239,2293,2309,2579,2833,3137,3319,3593,3673
%N Primes p such that p^3*q^3 + p^3 + q^3 is prime, where q is the next prime after p.
%H Charles R Greathouse IV, <a href="/A291339/b291339.txt">Table of n, a(n) for n = 1..10000</a>
%e a(2) = 3 is prime; 5 is the next prime: 3^3*5^3 + 3^3 + 5^3 = 27*125 + 27 + 125 = 3527 that is a prime.
%e a(3) = 7 is prime; 11 is the next prime: 7^3*11^3 + 7^3 + 11^3 = 343*1331 + 343 + 1331 = 458207 that is a prime.
%p select(p -> andmap(isprime,[p,(p^3*nextprime(p)^3+p^3+nextprime(p)^3)]), [seq(p, p=1..10^4)]);
%t Prime@Select[Range[1000], PrimeQ[Prime[#]^3*Prime[# + 1]^3 + Prime[#]^3 + Prime[# + 1]^3] &]
%o (Magma) [p: p in PrimesUpTo (5000) | IsPrime(p^3*q^3+p^3+q^3)];
%o (PARI) is(n) = my(q=nextprime(n+1)); ispseudoprime(n^3*q^3+n^3+q^3)
%o forprime(p=1, 3700, if(is(p), print1(p, ", "))) \\ _Felix Fröhlich_, Aug 22 2017
%o (PARI) list(lim)=my(v=List(),p=2,p3=8,q3); forprime(q=3,nextprime(lim\1+1), q3=q^3; if(isprime(p3*q3+p3+q3), listput(v,p)); p=q; p3=q3); Vec(v) \\ _Charles R Greathouse IV_, Aug 23 2017
%Y Cf. A000040, A001043, A006094, A030078, A096342, A120398, A126148, A152241.
%K nonn
%O 1,1
%A _K. D. Bajpai_, Aug 22 2017