A062762 - OEIS

%I #27 Sep 16 2024 12:47:30

%S 1,1,2,3,5,8,11,18,26,38,55,80,116,166,240,345,497,710,1016,1453,2073,

%T 2955,4211,5992,8523,12111,17202,24423,34648,49152,69694,98795,140009,

%U 198378,281016,398002,563612,797999,1129737,1599166,2263457,3203381

%N Number of powerful numbers not exceeding 2^n.

%C Number of terms x from A001694 for which x <= 2^n.

%H Chai Wah Wu, <a href="/A062762/b062762.txt">Table of n, a(n) for n = 0..127</a> (terms 0..90 from Daniel Suteu)

%F a(n) = Sum_{k=0..n} A062761(k). - _Daniel Suteu_, Feb 18 2020

%e Below 128, the 18 powerful numbers {1,4,8,9,16,25,...,100,108,121,125,128} can be found, so a(7)=18.

%t nn = 41; s = Union@ Flatten@ Table[a^2*b^3, {b, (2^nn)^(1/3)}, {a, Sqrt[(2^nn)/b^3]}]; Table[FirstPosition[s, 2^k][[1]], {k, 2, nn}] (* _Michael De Vlieger_, Oct 29 2023 *)

%o (PARI) a(n) = my(s=0,N=2^n); forsquarefree(k=1, sqrtnint(N, 3), s += sqrtint(N\k[1]^3)); s; \\ _Daniel Suteu_, Feb 18 2020

%o (Python)

%o from math import isqrt

%o from sympy import mobius, integer_nthroot

%o def A062762(n):

%o     def squarefreepi(n): return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))

%o     m = 1<<n

%o     c, l, j = squarefreepi(integer_nthroot(m,3)[0]), 0, isqrt(m)

%o     while j>1:

%o         k2 = integer_nthroot(m//j**2,3)[0]+1

%o         w = squarefreepi(k2-1)

%o         c += j*(w-l)

%o         l, j = w, isqrt(m//k2**3)

%o     return c-l # _Chai Wah Wu_, Sep 13 2024

%Y Cf. A001694, A029837, A036380, A036386, A062761.

%K nonn

%O 0,3

%A _Labos Elemer_, Jul 16 2001

%E a(19)-a(41) from _Donovan Johnson_, Oct 01 2009