login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

A349711
a(n) = Sum_{d|n} sopfr(d) * sopfr(n/d).
2
0, 0, 0, 4, 0, 12, 0, 16, 9, 20, 0, 44, 0, 28, 30, 40, 0, 54, 0, 68, 42, 44, 0, 104, 25, 52, 36, 92, 0, 124, 0, 80, 66, 68, 70, 147, 0, 76, 78, 152, 0, 164, 0, 140, 108, 92, 0, 200, 49, 110, 102, 164, 0, 144, 110, 200, 114, 116, 0, 298, 0, 124, 144, 140, 130, 244, 0, 212, 138, 236, 0, 300, 0, 148, 140
OFFSET
1,4
COMMENTS
Dirichlet convolution of A001414 with itself.
LINKS
FORMULA
Dirichlet g.f.: ( zeta(s) * Sum_{p prime} p/(p^s-1) )^2.
a(p^k) = (k^3-k)*p^2/6 = A000292(k-1)*p^2 for p prime. - Chai Wah Wu, Nov 28 2021
MAPLE
b:= proc(n) option remember; add(i[1]*i[2], i=ifactors(n)[2]) end:
a:= n-> add(b(d)*b(n/d), d=numtheory[divisors](n)):
seq(a(n), n=1..75); # Alois P. Heinz, Nov 26 2021
MATHEMATICA
sopfr[1] = 0; sopfr[n_] := Plus @@ Times @@@ FactorInteger@n; a[n_] := Sum[sopfr[d] sopfr[n/d], {d, Divisors[n]}]; Table[a[n], {n, 1, 75}]
PROG
(PARI) sopfr(n) = (n=factor(n))[, 1]~*n[, 2]; \\ A001414
a(n) = sumdiv(n, d, sopfr(d)*sopfr(n/d)); \\ Michel Marcus, Nov 26 2021
(Python)
from itertools import product
from sympy import factorint
def A349711(n):
f = factorint(n)
plist, m = list(f.keys()), sum(f[p]*p for p in f)
return sum((lambda x: x*(m-x))(sum(d[i]*p for i, p in enumerate(plist))) for d in product(*(list(range(f[p]+1)) for p in plist))) # Chai Wah Wu, Nov 27 2021
CROSSREFS
KEYWORD
nonn,look
AUTHOR
Ilya Gutkovskiy, Nov 26 2021
STATUS
approved