OFFSET
1,1
COMMENTS
a(49) = 1247136427. For n > 45, a(n) > 1.5*10^9 except n = 49. - David A. Corneth, Feb 15 2023
LINKS
David A. Corneth, Table of n, a(n) for n = 1..45
PROG
(PARI) a(n) = {my(p=prime(n), s=1); while(p!=lift(znprimroot(prime(s))), s++); s; } \\ Modified by Jinyuan Wang, Apr 03 2020
(PARI) upto(u, {maxn = 100}) = { my(t = 1, m = Map(), res = []); forprime(p = 2, oo, mapput(m, p, t); t++; if(t > maxn, break ) ); t = 1; u = prime(u); forprime(p = 2, u, c = lift(znprimroot(p)); if(mapisdefined(m, c), ind = mapget(m, c); if(ind > #res, res = concat(res, vector(ind - #res)) ); if(res[ind] == 0, res[ind] = t; ) ); t++ ); res } \\ David A. Corneth, Feb 15 2023
(Python)
from sympy import nextprime, primitive_root
def a(n):
k, pk, pn = 1, 2, prime(n)
while primitive_root(pk) != pn: k += 1; pk = nextprime(pk)
return k
print([a(n) for n in range(1, 19)]) # Michael S. Branicky, Feb 13 2023
(Python) # faster version for segments of sequence
from itertools import count, islice
from sympy import isprime, nextprime, prime, primepi, primitive_root
def agen(startk=1, startn=1): # generator of terms
p, vdict, adict, n = prime(startk), dict(), dict(), startn
for k in count(startk):
v = primitive_root(p)
if v not in vdict and isprime(v):
vdict[v] = k
adict[primepi(v)] = k
while n in adict: yield adict[n]; n += 1
p = nextprime(p)
print(list(islice(agen(), 18))) # Michael S. Branicky, Feb 14 2023
CROSSREFS
KEYWORD
nonn
AUTHOR
Benoit Cloitre, Feb 02 2003
EXTENSIONS
a(17)-a(18) from Jinyuan Wang, Apr 03 2020
a(19)-a(36) from Michael S. Branicky, Feb 14 2023
STATUS
approved