login
A118896
Number of powerful numbers <= 10^n.
6
1, 4, 14, 54, 185, 619, 2027, 6553, 21044, 67231, 214122, 680330, 2158391, 6840384, 21663503, 68575557, 217004842, 686552743, 2171766332, 6869227848, 21725636644, 68709456167, 217293374285, 687174291753, 2173105517385, 6872112993377, 21731852479862, 68722847672629, 217322225558934, 687236449779456, 2173239433013146
OFFSET
0,2
COMMENTS
These numbers agree with the asymptotic formula c*sqrt(x), with c=2.1732...(A090699). - T. D. Noe, May 09 2006
Bateman & Grosswald proved that the number of powerful numbers up to x is zeta(3/2)/zeta(3) * x^1/2 + zeta(2/3)/zeta(2) * x^1/3 + o(x^1/6). This approximates the series very closely: up to a(24), all absolute errors are less than 75. - Charles R Greathouse IV, Sep 23 2008
LINKS
Charles R Greathouse IV and Hiroaki Yamanouchi, Table of n, a(n) for n = 0..45 (terms a(0)-a(32) from Charles R Greathouse IV)
Michael Filaseta and Ognian Trifonov, The distribution of squarefull numbers in short intervals, Acta Arithmetica 67 (1994), pp. 323-333.
Paul T. Bateman and Emil Grosswald, On a theorem of Erdős and Szekeres, Illinois J. Math. 2:1 (1958), p. 88-98.
Eric Weisstein's World of Mathematics, Powerful Number
FORMULA
Pi(x) = Sum_{i=1..x^(1/3)} floor(sqrt(x/i^3)) only if i is squarefree. - Robert G. Wilson v, Aug 12 2014
MAPLE
f:= m -> nops({seq(seq(a^2*b^3, b=1..floor((m/a^2)^(1/3))), a=1..floor(sqrt(m)))}):
seq(f(10^n), n=0..10); # Robert Israel, Aug 12 2014
MATHEMATICA
f[n_] := Block[{max = 10^n}, Length@ Union@ Flatten@ Table[ a^2*b^3, {b, max^(1/3)}, {a, Sqrt[ max/b^3]}]]; Array[f, 13, 0] (* Robert G. Wilson v, Aug 11 2014 *)
powerfulNumberPi[n_] := Sum[ If[ SquareFreeQ@ i, Floor[ Sqrt[ n/i^3]], 0], {i, n^(1/3)}]; Array[ powerfulNumberPi[10^#] &, 27, 0] (* Robert G. Wilson v, Aug 12 2014 *)
PROG
(PARI) a(n)=n=10^n; sum(k=1, floor((n+.5)^(1/3)), if(issquarefree(k), sqrtint(n\k^3))) \\ Charles R Greathouse IV, Sep 23 2008
(Python)
from math import isqrt
from sympy import integer_nthroot, factorint
def A118896(n):
m = 10**n
return sum(isqrt(m//x**3) for x in range(1, integer_nthroot(m, 3)[0]+1) if max(factorint(x).values(), default=0)<=1) # Chai Wah Wu, May 13 2023
(Python)
# faster program
def A118896(n):
def squarefreepi(n): return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
m, c, l = 10**n, 0, 0
j = isqrt(m)
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)
c += squarefreepi(integer_nthroot(m, 3)[0])-l
return c # Chai Wah Wu, Sep 09 2024
CROSSREFS
Sequence in context: A308555 A000651 A192247 * A145211 A060898 A180142
KEYWORD
nonn
AUTHOR
Eric W. Weisstein, May 05 2006
EXTENSIONS
More terms from T. D. Noe, May 09 2006
a(13)-a(24) from Charles R Greathouse IV, Sep 23 2008
a(25)-a(29) from Charles R Greathouse IV, May 30 2011
a(30) from Charles R Greathouse IV, May 31 2011
STATUS
approved