%I #24 May 10 2022 13:48:21

%S 0,1,5,9,16,20,36,40,50,57,73,77,105,109,125,141,154,158,186,190,218,

%T 234,250,254,294,301,317,327,355,359,423,427,443,459,475,491,540,544,

%U 560,576,616,620,684,688,716,744,760,764,816,823,851,867,895,899,939,955,995

%N a(n) = Sum_{k=1..n} tau(k^3), where tau is the number of divisors function A000005.

%H Karl-Heinz Hofmann, <a href="/A353551/b353551.txt">Table of n, a(n) for n = 0..10000</a>

%F a(n) = Sum_{k=1..n} tau(k^3).

%F a(n) = a(n-1) + A048785(n) for n >= 1, a(0) = 0.

%e A048785(0) = 0

%e + A048785(1) = 1

%e + A048785(2) = 4

%e + A048785(3) = 4

%e ------------------

%e = A353551(3) = 9

%p a:= proc(n) option remember; `if`(n=0, 0, a(n-1)+numtheory[tau](n^3)) end:

%p seq(a(n), n=0..100); # _Alois P. Heinz_, May 08 2022

%t Accumulate[Join[{0}, Table[DivisorSigma[0, k^3], {k, 1, 50}]]] (* _Amiram Eldar_, May 08 2022 *)

%o (Python) from sympy import divisor_count

%o def A048785(n): return divisor_count(n**3)

%o def A353551(n): return sum(A048785(n) for n in range(1, n))

%o print([A353551(n) for n in range(1, 58)])

%o (PARI) a(n) = sum(k=1, n, numdiv(k^3)); \\ _Michel Marcus_, May 08 2022

%o (Python)

%o from math import prod

%o from sympy import factorint

%o def A353551(n): return sum(prod(3*e+1 for e in factorint(k).values()) for k in range(1,n+1)) # _Chai Wah Wu_, May 10 2022

%Y Partial sums of A048785.

%Y Cf. A000005, A006218, A061503 (squares).

%K nonn

%O 0,3

%A _Karl-Heinz Hofmann_, May 07 2022