OFFSET
1,2
MAPLE
select(n -> GaussInt:-GIprime(n*(1+I)^n+1), [$1..50000]); # Robert Israel, May 08 2023
MATHEMATICA
Do[ If[ PrimeQ[ n * (1 + I)^n + 1, GaussianIntegers -> True], Print[n] ], {n, 1, 4000} ]
PROG
(Python)
from itertools import count, islice
from sympy import I
from sympy.ntheory import is_gaussian_prime
def A058770_gen(startvalue=1): # generator of terms
x = (1+I)**(m:=max(startvalue, 1))
for k in count(m):
if is_gaussian_prime(k*x+1):
yield k
x *= (1+I)
CROSSREFS
KEYWORD
nonn,hard,more
AUTHOR
Robert G. Wilson v, Jan 02 2001
EXTENSIONS
Corrected by Robert Israel, May 08 2023
STATUS
approved