OFFSET
1,3
LINKS
Michael De Vlieger, Table of n, a(n) for n = 1..10000
FORMULA
EXAMPLE
a(16) is 2 since A000961(16) is 27 which is 3^3 = (p_2)^3, i.e., the prime index of 3 is 2.
MAPLE
b:= proc(n) option remember; local k; for k from
1+b(n-1) while nops(ifactors(k)[2])>1 do od; k
end: b(1):=1:
a:= n-> `if`(n=1, 0, numtheory[pi](ifactors(b(n))[2, 1$2])):
seq(a(n), n=1..100); # Alois P. Heinz, Feb 20 2020
MATHEMATICA
mxn = 500; Join[{0}, Transpose[ Sort@ Flatten[ Table[ {Prime@n^ex, n}, {n, PrimePi@ mxn}, {ex, Log[Prime@n, mxn]}], 1]][[2]]]
PROG
(PARI) lista(nn) = {print1(0); for(n=2, nn, if(isprimepower(n, &p), print1(", ", primepi(p)))); } \\ Jinyuan Wang, Feb 19 2020
(Python)
from sympy import primepi, integer_nthroot, primefactors
def A330669(n):
if n == 1: return 0
def f(x): return int(n-2+x-sum(primepi(integer_nthroot(x, k)[0]) for k in range(1, x.bit_length())))
kmin, kmax = 1, 2
while f(kmax) >= kmax:
kmax <<= 1
while True:
kmid = kmax+kmin>>1
if f(kmid) < kmid:
kmax = kmid
else:
kmin = kmid
if kmax-kmin <= 1:
break
return int(primepi(primefactors(kmax)[0])) # Chai Wah Wu, Aug 20 2024
CROSSREFS
KEYWORD
easy,nonn
AUTHOR
Grant E. Martin and Robert G. Wilson v, Dec 23 2019
STATUS
approved