login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

A116434
Consider the array T(r,c), the number of c-almost primes less than or equal to r^c. This is the diagonal T(r,r-1).
2
0, 1, 3, 13, 90, 726, 7089, 78369, 973404, 13377156, 201443165, 3297443264, 58304208767, 1107693755122
OFFSET
1,3
MATHEMATICA
AlmostPrimePi[k_Integer, n_] := Module[{a, i}, a[0] = 1; If[k == 1, PrimePi[n], Sum[PrimePi[n/Times @@ Prime[Array[a, k - 1]]] - a[k - 1] + 1, Evaluate[ Sequence @@ Table[{a[i], a[i - 1], PrimePi[(n/Times @@ Prime[Array[a, i - 1]])^(1/(k - i + 1))]}, {i, k - 1}]] ]]]; (* Eric W. Weisstein, Feb 07 2006 *)
Do[ Print@ AlmostPrimePi[n, (n + 1)^n], {n, 11}]
PROG
(Python)
from math import isqrt, prod
from sympy import primerange, integer_nthroot, primepi
def A116434(n):
def almostprimepi(n, k):
def g(x, a, b, c, m): yield from (((d, ) for d in enumerate(primerange(b, isqrt(x//c)+1), a)) if m==2 else (((a2, b2), )+d for a2, b2 in enumerate(primerange(b, integer_nthroot(x//c, m)[0]+1), a) for d in g(x, a2, b2, c*b2, m-1)))
return int(sum(primepi(n//prod(c[1] for c in a))-a[-1][0] for a in g(n, 0, 1, 1, k)) if k>1 else primepi(n))
return almostprimepi((n+1)**n, n) # Chai Wah Wu, Sep 02 2024
CROSSREFS
Sequence in context: A097711 A114477 A345104 * A174290 A034513 A257661
KEYWORD
hard,more,nonn
AUTHOR
EXTENSIONS
Name rephrased by R. J. Mathar, Jun 20 2021
a(13)-a(14) from Max Alekseyev, Oct 12 2023
STATUS
approved