A096165 Prime powers with exponents that are themselves prime powers.
2, 3, 4, 5, 7, 8, 9, 11, 13, 16, 17, 19, 23, 25, 27, 29, 31, 32, 37, 41, 43, 47, 49, 53, 59, 61, 67, 71, 73, 79, 81, 83, 89, 97, 101, 103, 107, 109, 113, 121, 125, 127, 128, 131, 137, 139, 149, 151, 157, 163, 167, 169, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227
Offset: 1
Keywords
Examples
512=2^9=2^(3^2), A000961(118)=A000040(1)^A000961(118), therefore 512 is a term; 64=2^6, but 6 is not a prime power, therefore 64 is not a term.
Links
- Reinhard Zumkeller, Table of n, a(n) for n = 1..10000
Programs
-
Haskell
a096165 n = a096165_list !! (n-1) a096165_list = filter ((== 1) . a010055 . a001222) $ tail a000961_list -- Reinhard Zumkeller, Nov 17 2011
-
Maple
F:= proc(t) local P; P:= ifactors(t)[2]; nops(P) = 1 and (P[1][2]=1 or nops(numtheory:-factorset(P[1][2]))=1) end proc: select(F, [$2..1000]); # Robert Israel, Jul 20 2015
-
Mathematica
Select[Range@ 240, Or[PrimeQ@ #, PrimePowerQ@ # && PrimePowerQ@ FactorInteger[#][[1, 2]]] &] (* Michael De Vlieger, Jul 20 2015 *)
-
PARI
is(n)=while(1,if(!(n=isprimepower(n)),return(0),if(n==1,return(1)))) \\ Anders Hellström, Jul 19 2015
-
PARI
ispp(n)=n==1 || isprimepower(n) is(n)=ispp(isprimepower(n)) \\ Charles R Greathouse IV, Oct 19 2015
-
Python
from sympy import primepi, integer_nthroot, factorint def A096165(n): def f(x): return int(n+x-sum(primepi(integer_nthroot(x,k)[0]) for k in range(1,x.bit_length()) if len(factorint(k))<=1)) def bisection(f,kmin=0,kmax=1): while f(kmax) > kmax: kmax <<= 1 while kmax-kmin > 1: kmid = kmax+kmin>>1 if f(kmid) <= kmid: kmax = kmid else: kmin = kmid return kmax return bisection(f,n,n) # Chai Wah Wu, Sep 12 2024
Formula
a(n) ~ n log n. - Charles R Greathouse IV, Oct 19 2015
Comments