%I #32 Oct 24 2023 12:40:11

%S 1,1,4,1,5,8,1,7,11,15,1,11,19,23,21,1,19,41,47,33,33,1,35,103,125,77,

%T 57,41,1,67,281,395,255,149,71,56,1,131,799,1373,1025,555,205,103,69,

%U 1,259,2321,5027,4503,2537,905,325,130,87,1,515,6823,18965,20657,12867,4945,1585,442,170,99

%N Square array T(n,k), n >= 1, k >= 0, read by antidiagonals downwards, where T(n,k) = Sum_{j=1..n} j^k * binomial(floor(n/j)+1,2).

%F G.f. of column k: (1/(1-x)) * Sum_{j>=1} j^k * x^j/(1 - x^j)^2.

%F T(n,k) = Sum_{j=1..n} j * sigma_{k-1}(j).

%e Square array begins:

%e 1, 1, 1, 1, 1, 1, 1, ...

%e 4, 5, 7, 11, 19, 35, 67, ...

%e 8, 11, 19, 41, 103, 281, 799, ...

%e 15, 23, 47, 125, 395, 1373, 5027, ...

%e 21, 33, 77, 255, 1025, 4503, 20657, ...

%e 33, 57, 149, 555, 2537, 12867, 68969, ...

%t T[n_, k_] := Sum[j^k * Binomial[Floor[n/j] + 1, 2], {j, 1, n}]; Table[T[k, n - k], {n, 1, 11}, {k, 1, n}] // Flatten (* _Amiram Eldar_, Jul 28 2022 *)

%o (PARI) T(n, k) = sum(j=1, n, j^k*binomial(n\j+1, 2));

%o (PARI) T(n, k) = sum(j=1, n, j*sigma(j, k-1));

%o (Python)

%o from itertools import count, islice

%o from math import isqrt

%o from sympy import bernoulli

%o def A356124_T(n,k): return ((s:=isqrt(n))*(s+1)*(bernoulli(k+1)-bernoulli(k+1,s+1))+sum(w**k*(k+1)*((q:=n//w)*(q+1))+(w*(bernoulli(k+1,q+1)-bernoulli(k+1))<<1) for w in range(1,s+1)))//(k+1)>>1

%o def A356124_gen(): # generator of terms

%o return (A356124_T(k+1,n-k-1) for n in count(1) for k in range(n))

%o A356124_list = list(islice(A356124_gen(),30)) # _Chai Wah Wu_, Oct 24 2023

%Y Column k=0..4 give A024916, A143127, A143128, A356125, A356126.

%Y T(n,n) gives A356129.

%Y T(n,n+1) gives A356128.

%Y Cf. A279394, A319649, A350106.

%K nonn,tabl

%O 1,3

%A _Seiichi Manyama_, Jul 27 2022