A344992 a(n) = Sum_{1 <= i <= j <= k <= m <= n} gcd(i,j,k,m).

1, 6, 18, 44, 83, 159, 249, 401, 592, 867, 1163, 1655, 2122, 2796, 3594, 4594, 5579, 7046, 8394, 10328, 12339, 14699, 17021, 20441, 23526, 27317, 31379, 36323, 40846, 47300, 52786, 59954, 67191, 75380, 83720, 94662, 103837, 115137, 126851, 141059, 153440

Offset

1,2

Comments

In general, if g.f.: 1/(1-x) * Sum_{j>=1} phi(j) * x^j/(1 - x^j)^k, where k > 2 and phi is the Euler totient function (A000010), then a(n) ~ zeta(k-1) * n^k / (k! * zeta(k)).

Links

Vaclav Kotesovec, Table of n, a(n) for n = 1..10000

Formula

G.f.: 1/(1-x) * Sum_{k>=1} phi(k) * x^k/(1 - x^k)^4, where phi is the Euler totient function (A000010).

a(n) = Sum_{k=1..n} Sum_{d|k} phi(k/d) * d*(d+1)*(d+2)/6.

a(n) ~ 15 * zeta(3) * n^4 / (4*Pi^4).

Mathematica

Table[Sum[Sum[Sum[Sum[GCD[i, j, k, m], {i, 1, j}], {j, 1, k}], {k, 1, m}], {m, 1, n}], {n, 1, 100}]
nmax = 100; Rest[CoefficientList[Series[1/(1-x) * Sum[EulerPhi[k]*x^k/(1 - x^k)^4, {k, 1, nmax}], {x, 0, nmax}], x]]
Accumulate[Table[Sum[EulerPhi[n/d] * d*(d+1)*(d+2)/6, {d, Divisors[n]}], {n, 1, 100}]] (* faster *)

Prog

(PARI) a(n) = sum(i=1, n, sum(j=i, n, sum(k=j, n, sum(m=k, n, gcd([i, j, k, m]))))); \\ Michel Marcus, Jun 06 2021

Crossrefs

Cf. A000010, A059358, A272718, A309322, A309323, A344521.

Sequence in context: A191101 A184630 A009957 * A011929 A070735 A136028

Adjacent sequences: A344989 A344990 A344991 * A344993 A344994 A344995

Keyword

nonn

Author

Vaclav Kotesovec, Jun 05 2021

Status

approved