cp's OEIS Frontend

A076407 Sum of perfect powers <= n.

%I A076407 #22 Feb 16 2025 08:32:47
%S A076407 1,1,1,5,5,5,5,13,22,22,22,22,22,22,22,38,38,38,38,38,38,38,38,38,63,
%T A076407 63,90,90,90,90,90,122,122,122,122,158,158,158,158,158,158,158,158,
%U A076407 158,158,158,158,158,207,207,207,207,207,207,207,207,207,207,207,207,207
%N A076407 Sum of perfect powers <= n.
%H A076407 Robert Israel, <a href="/A076407/b076407.txt">Table of n, a(n) for n = 1..10000</a>
%H A076407 Eric Weisstein's World of Mathematics, <a href="https://mathworld.wolfram.com/PerfectPower.html">Perfect Powers</a>.
%F A076407 a(n) = 1 - Sum_{k=2..floor(log_2(n))} mu(k) * (F(k, floor(n^(1/k))) - 1), where F(k, n) = Sum_{j=1..n} j^k = (Bernoulli(k+1, n+1) - Bernoulli(k+1, 1))/(k+1). - _Daniel Suteu_, Aug 19 2023
%e A076407 Sum of the 8 perfect powers <= 32: a(32) = 1+4+8+9+16+25+27+32 = 122.
%p A076407 N:= 100: # for a(1)..a(N)
%p A076407 V:= Vector(N,1):
%p A076407 pps:= {seq(seq(x^k,k=2..floor(log[x](N))),x=2..floor(sqrt(N)))}:
%p A076407 for y in pps do
%p A076407    V[y..N]:= V[y..N] +~ y
%p A076407 od:
%p A076407 convert(V,list); # _Robert Israel_, Oct 19 2023
%o A076407 (PARI)
%o A076407 F(k,n) = (subst(bernpol(k+1), x, n+1) - subst(bernpol(k+1), x, 1)) / (k+1);
%o A076407 a(n) = 1 - sum(k=2, logint(n,2), moebius(k) * (F(k, sqrtnint(n,k)) - 1)); \\ _Daniel Suteu_, Aug 19 2023
%Y A076407 Cf. A001597, A076408, A069623.
%K A076407 nonn
%O A076407 1,4
%A A076407 _Reinhard Zumkeller_, Oct 09 2002