A381496 Number of powerful numbers that are not prime powers that do not exceed 10^n.

0, 0, 3, 28, 133, 510, 1790, 5997, 19639, 63541, 204037, 652173, 2078320, 6609816, 20993381, 66612867, 211222374, 669428537, 2120835892, 6717184256, 21270247404, 67341572823, 213173925948, 674739560651, 2135491756895, 6758117426102, 21385762133815, 67670426242420

Offset: 0

Author

Michael De Vlieger, Feb 25 2025

Keywords

• nonn

Comments

Number of k such that omega(k) > 1 and rad(k)^2 | k (i.e., in A286708) that do not exceed 10^n, where omega = A001221 and rad = A007947.

Examples

			Let S = A286708 = A001694 \ A246547 = A126706 \ A001694.
a(0) = a(1) = 0 since 36 is the smallest term in S.
a(2) = 3 since S(1..3) = {36, 72, 100}.
a(3) = 28 since S(4..28) = {108, 144, ..., 972, 1000}.
a(4) = 133 since S(29..133) = {1089, 1125, ..., 9801, 10000}, etc.
		

Crossrefs

Cf. A001694, A118896, A126706, A246547, A267574, A286708, A380430, A380431, A381391.

Programs

Mathematica

Table[Sum[Boole[SquareFreeQ[k]]*Floor[Sqrt[10^n/k^3]], {k, 10^(n/3)}] - Sum[PrimePi[10^(n/k)], {k, 2, n*Log2[10]}] - 1, {n, 0, 12}]

Python

from math import isqrt
from sympy import primepi, integer_nthroot, mobius
def A381496(n):
    def squarefreepi(n): return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
    m, l = 10**n, 0
    j, c = isqrt(m), -1-sum(primepi(integer_nthroot(m,k)[0]) for k in range(2,m.bit_length())),
    while j>1:
        k2 = integer_nthroot(m//j**2,3)[0]+1
        w = squarefreepi(k2-1)
        c += j*(w-l)
        l, j = w, isqrt(m//k2**3)
    return c+squarefreepi(integer_nthroot(m,3)[0])-l # Chai Wah Wu, Feb 25 2025

Formula

a(n) = -1 + Sum_{k=1..10^(n/3)} [rad(k)=k]*floor(sqrt(10^n/k^3)) - Sum_{k=2..n*log_2(10)} pi(10^(n/k)).

a(n) = -1 + A118896(n) - A267574(n).

a(n) < A381391(n) for n > 0 since A286708 is a proper subset of A126706.