%I #10 Jul 11 2024 13:15:52
%S 1,5,10,17,26,50,50,85,82,130,122,170,170,250,260,257,290,410,362,442,
%T 500,610,530,850,626,850,820,850,842,1300,962,1285,1220,1450,1300,
%U 1394,1370,1810,1700,2210,1682,2500,1850,2074,2132,2650,2210,2570,2402,3130,2900
%N The sum of the squares of the infinitary divisors of n.
%C Also the sum of the infinitary divisors of n^2.
%H Amiram Eldar, <a href="/A374539/b374539.txt">Table of n, a(n) for n = 1..10000</a>
%F a(n) = A049417(n^2).
%F a(n) <= A001157(n), with equality if and only if n is in A036537.
%F Multiplicative with a(p^e) = Product{k>=1, e_k=1} (p^(2^(k+1)) + 1), where e = Sum_{k} e_k * 2^k is the binary representation of e, i.e., e_k is bit k of e.
%F Sum_{k=1..n} a(k) ~ c * n^3 / 3, where c = Product_{P} (1 + 1/(P^2*(P+1))) = 1.14142906130350119631..., and P are numbers of the form p^(2^k) where p is prime and k >= 0 (A050376).
%t f[p_, e_] := p^(2^Position[Reverse@IntegerDigits[e, 2], _?(# == 1 &)]); a[1] = 1; a[n_] := Times @@ (Flatten@(f @@@ FactorInteger[n]) + 1); Array[a, 100]
%o (PARI) a(n) = {my(f = factor(n), b); prod(i = 1, #f~, b = binary(2 * f[i, 2]); prod(k=1, #b, if(b[k], 1+f[i, 1]^(2^(#b-k)), 1)));}
%o (Python)
%o from math import prod
%o from sympy import factorint
%o def A374539(n): return prod(p**(1<<i)+1 for p, e in factorint(n).items() for i, j in enumerate(bin(e)[-1:1:-1],1) if j=='1') # _Chai Wah Wu_, Jul 11 2024
%Y Cf. A049417, A050376, A077609.
%Y Similar sequences: A001157, A034676, A050999, A067558, A076577, A351265, A351307, A351647.
%K nonn,easy,mult
%O 1,2
%A _Amiram Eldar_, Jul 11 2024