%I #19 Jun 04 2022 21:16:56
%S 0,2,3,6,5,11,7,12,12,17,11,27,13,23,23,28,17,38,19,41,31,35,23,57,30,
%T 41,36,55,29,71,31,60,47,53,47,90,37,59,55,87,41,95,43,83,77,71,47,
%U 121,56,92,71,97,53,116,71,117,79,89,59,167,61,95,103,120,83,143,67,125,95
%N Sum of the divisors of n whose cube does not divide n.
%H Harvey P. Dale, <a href="/A345321/b345321.txt">Table of n, a(n) for n = 1..1000</a>
%F a(n) = Sum_{k=1..n} k * (ceiling(n/k^3) - floor(n/k^3)) * (1 - ceiling(n/k) + floor(n/k)).
%F a(n) = A000203(n) - A333843(n). - _Rémy Sigrist_, Jun 14 2021
%e a(16) = 28; The divisors of 16 whose cube does not divide 16 are: 4, 8 and 16. The sum of these divisors is then 4 + 8 + 16 = 28.
%t Table[Sum[k (Ceiling[n/k^3] - Floor[n/k^3]) (1 - Ceiling[n/k] + Floor[n/k]), {k, n}], {n, 80}]
%t Table[Total[Select[Divisors[n],Mod[n,#^3]!=0&]],{n,100}] (* _Harvey P. Dale_, May 01 2022 *)
%o (PARI) a(n) = sumdiv(n, d, if (n % d^3, d)); \\ _Michel Marcus_, Jun 13 2021
%o (Python 3.8+)
%o from math import prod
%o from sympy import factorint
%o def A345321(n):
%o f = factorint(n).items()
%o return prod((p**(q+1)-1)//(p-1) for p, q in f) - prod((p**(q//3+1)-1)//(p-1) for p, q in f) # _Chai Wah Wu_, Jun 14 2021
%Y Cf. A000203, A333843.
%K nonn
%O 1,2
%A _Wesley Ivan Hurt_, Jun 13 2021
|