A300909 - OEIS

%I #25 Jan 24 2024 10:15:46

%S 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,17,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,17,1,

%T 1,1,1,1,1,1,1,1,1,1,1,1,1,1,17,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,17,1,1,

%U 1,1,1,1,1,1,1,1,1,1,1,1,1,17,82,1,1,1,1,1,1,1,1,1,1,1,1,1,1,17,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,17

%N Sum of 4th powers dividing n.

%C Multiplicative with a(p^e) = (p^(4*(1+floor(e/4)))-1)/(p^4-1). - _Robert Israel_, Mar 15 2018

%H Antti Karttunen, <a href="/A300909/b300909.txt">Table of n, a(n) for n = 1..65537</a>

%H A. Dixit, B. Maji, A. Vatwani, <a href="https://arxiv.org/abs/2303.09937">Voronoi summation formula for the generalized divisor function sigma_z^k(n)</a>, arXiv:2303.09937 [math.NT], 2023, sigma(z=4,k=4,n).

%F G.f.: Sum_{k>=1} k^4*x^(k^4)/(1 - x^(k^4)).

%F L.g.f.: -log(Product_{k>=1} (1 - x^(k^4))) = Sum_{n>=1} a(n)*x^n/n.

%F D.g.f.: zeta(s)*zeta(4s-4). - _Robert Israel_, Mar 15 2018

%F Sum_{k=1..n} a(k) ~ zeta(5/4)*n^(5/4)/5 - n/2. - _Vaclav Kotesovec_, Dec 01 2020

%e a(16) = 17 because 16 has 5 divisors {1, 2, 4, 8, 16} among which 2 divisors {1, 16} are 4th powers and 1 + 16 = 17.

%e L.g.f.: L(x) = x + x^2/2 + x^3/3 + x^4/4 + x^5/5 + x^6/6 + x^7/7 + x^8/8 + x^9/9 + x^10/10 + x^11/11 + x^12/12 + x^13/13 + x^14/14 + x^15/15 + 17*x^16/16 + x^17/17 + ...

%e exp(L(x)) = 1 + x + x^2 + x^3 + x^4 + x^5 + x^6 + x^7 + x^8 + x^9 + x^10 + x^11 + x^12 + x^13 + x^14 + x^15 + 2*x^16 + 2*x^17 + ... + A046042(n)*x^n + ...

%p N:= 1000: # for a(1)..a(N)

%p V:= Vector(N,1):

%p for m from 2 to floor(N^(1/4)) do

%p   R:= [seq(i,i=m^4 .. N, m^4)];

%p   V[R]:= map(`+`,V[R],m^4)

%p od:

%p convert(V,list); # _Robert Israel_, Mar 15 2018

%t Table[DivisorSum[n, # &, IntegerQ[#^(1/4)] &], {n, 112}]

%t nmax = 112; Rest[CoefficientList[Series[Sum[k^4 x^k^4/(1 - x^k^4), {k, 1, 10}], {x, 0, nmax}], x]]

%t nmax = 112; Rest[CoefficientList[Series[-Log[Product[(1 - x^k^4), {k, 1, 10}]], {x, 0, nmax}], x] Range[0, nmax]]

%t f[p_, e_] := (p^(4*(1 + Floor[e/4])) - 1)/(p^4 - 1); a[1] = 1; a[n_] := Times @@ f @@@ FactorInteger[n]; Array[a, 100] (* _Amiram Eldar_, May 01 2020 *)

%o (PARI) a(n) = sumdiv(n, d, d*ispower(d, 4)); \\ _Michel Marcus_, Mar 15 2018

%Y Cf. A000583, A001159, A035316, A046042, A046100 (positions of ones), A063775, A113061.

%K nonn,mult

%O 1,16

%A _Ilya Gutkovskiy_, Mar 15 2018