login
A368107
Prime powers p^m such that p | m.
1
4, 16, 27, 64, 256, 729, 1024, 3125, 4096, 16384, 19683, 65536, 262144, 531441, 823543, 1048576, 4194304, 9765625, 14348907, 16777216, 67108864, 268435456, 387420489, 1073741824, 4294967296, 10460353203, 17179869184, 30517578125, 68719476736, 274877906944, 282429536481
OFFSET
1,1
COMMENTS
Proper subset of A072873, which in turn is a proper subset of A342090.
This sequence represents the prime power block in A072873 and A342090.
A342090 \ {a(n)} is in A126706.
A072873 \ {{1} U {a(n)}} is in A286708, in turn a proper subset of A001694.
Contains A051674.
LINKS
FORMULA
Sum_{n>=1} 1/a(n) = Sum_{n>=1} 1/A088730(n) = 0.372116188498... . - Amiram Eldar, Jan 20 2024
EXAMPLE
This sequence contains prime powers of the following form:
2^2, 2^4, i.e., 2^k such that k is even.
3^3, 3^6, 3^9, i.e., 3^k such that 3 | k.
5^5, 5^10, 5^15, i.e., 5^k such that 5 | k, etc.
MAPLE
N:= 10^13: # for terms <= N
R:= NULL:
for i from 1 do
p:= ithprime(i);
if p^p > N then break fi;
R:= R, seq(p^k, k=p..floor(log[p](N)), p);
od:
sort([R]); # Robert Israel, Jan 16 2024
MATHEMATICA
nn = 10^12; i = 1; p = 2; While[p^p <= nn, p = NextPrime[p] ];
MapIndexed[Set[S[First[#2]], #1] &, Prime@ Range@ PrimePi[p] ];
Union@ Reap[
While[j = S[i];
While[S[i]^j < nn,
Sow[S[i]^j]; j += S[i] ]; j > 2,
i++] ][[-1, 1]]
PROG
(Python)
import heapq
from itertools import islice
from sympy import nextprime
def agen(): # generator of terms
v, h, m, nextp = 4, [(4, 2)], 4, 3
while True:
v, p = heapq.heappop(h)
yield v
if v >= m:
m = nextp**nextp
heapq.heappush(h, (m, nextp))
nextp = nextprime(nextp)
heapq.heappush(h, (v*p**p, p))
print(list(islice(agen(), 31))) # Michael S. Branicky, Jan 16 2024
(Python)
from sympy import integer_nthroot, primefactors
def A368107(n):
def f(x):
c = n+x
for k in range(1, x.bit_length()):
m = integer_nthroot(x, k)[0]
c -= sum(1 for p in primefactors(k) if p<=m)
return c
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
KEYWORD
nonn,easy
AUTHOR
Michael De Vlieger, Jan 15 2024
STATUS
approved