cp's OEIS Frontend

A333753 Sum of prime power divisors of n that are <= sqrt(n).

%I A333753 #9 Feb 15 2023 13:43:34
%S A333753 0,0,0,2,0,2,0,2,3,2,0,5,0,2,3,6,0,5,0,6,3,2,0,9,5,2,3,6,0,10,0,6,3,2,
%T A333753 5,9,0,2,3,11,0,5,0,6,8,2,0,9,7,7,3,6,0,5,5,13,3,2,0,14,0,2,10,14,5,5,
%U A333753 0,6,3,14,0,17,0,2,8,6,7,5,0,19,12,2,0,16,5,2,3,14,0,19
%N A333753 Sum of prime power divisors of n that are <= sqrt(n).
%H A333753 Robert Israel, <a href="/A333753/b333753.txt">Table of n, a(n) for n = 1..10000</a>
%F A333753 G.f.: Sum_{p prime, k>=1} p^k * x^(p^(2*k)) / (1 - x^(p^k)).
%p A333753 f:= proc(n) local F,i,j,t;
%p A333753   F:= ifactors(n)[2];
%p A333753   t:= 0;
%p A333753   for i from 1 to nops(F) do
%p A333753     j:= min(F[i,2],ilog[F[i,1]^2](n));
%p A333753     t:= t + (F[i,1]^j-1)*F[i,1]/(F[i,1]-1)
%p A333753   od;
%p A333753   t
%p A333753 end proc:
%p A333753 map(f, [$1..100]); # _Robert Israel_, Feb 15 2023
%t A333753 Table[DivisorSum[n, # &, # <= Sqrt[n] && PrimePowerQ[#] &], {n, 1, 90}]
%t A333753 nmax = 90; CoefficientList[Series[Sum[Boole[PrimePowerQ[k]] k x^(k^2)/(1 - x^k), {k, 1, nmax}], {x, 0, nmax}], x] // Rest
%o A333753 (PARI) a(n) = sumdiv(n, d, if ((d^2<=n) && isprimepower(d), d)); \\ _Michel Marcus_, Apr 03 2020
%Y A333753 Cf. A023889, A066839, A069289, A069293, A097974, A246655, A333750, A333751, A333752.
%K A333753 nonn
%O A333753 1,4
%A A333753 _Ilya Gutkovskiy_, Apr 03 2020