2026-06-23
This manuscript provides an algorithm for deriving closed formulas for multifold sums of powers of integers by combining variations of Newton’s interpolation formula with hockey-stick family identities for binomial coefficients.
The story begins in ancient Greece, where discrete approximations were employed to obtain continuous areas through the method of exhaustion. In particular, Archimedes (287–212 BC) evaluated sums of squares while determining the area enclosed by his spiral (Pengelley 2002), leading to calculations analogous to the modern identity \[\begin{align*} 1^2 + 2^2 + 3^2 + \cdots + n^2 = \frac{n(n+1)(2n+1)}{6}. \end{align*}\]
Moving to early modern Europe, Pierre de Fermat (1601–1665) described the problem of finding formulas for sums of powers as “perhaps the most beautiful problem of all arithmetic” in his correspondence with Blaise Pascal (1623–1662). Fermat’s contribution to the subject was substantial: he revived the study of figurate numbers and investigated their arithmetic properties, including what is now known as Fermat’s polygonal number theorem. Applications of polygonal numbers to sums of powers are discussed extensively in (Marko and Litvinov 2020, (2020)).
Meanwhile, Johann Faulhaber (1580–1635) investigated sums of powers in his celebrated work Academia Algebrae (Faulhaber 1631, (1631)). His achievements were extraordinary, as he obtained formulas up to \(\Sigma^{1}\,{n}^{25}=1^{25}+2^{25}+\cdots+n^{25}\). Faulhaber explicitly tabulated formulas through the seventeenth power and provided sufficient information for determining higher cases, including the twenty-fifth power; the latter reconstruction was later carried out by Knuth. More importantly, Faulhaber discovered deep structural patterns, particularly the representation of odd-power sums as polynomials in triangular numbers. What remained missing was a completely general symbolic formula.
A major breakthrough was achieved by Jacob Bernoulli (1655–1705) in his posthumously published work Ars Conjectandi (Bernoulli 1713, (1713)). There, Bernoulli introduced the numbers that now bear his name, obtaining a general formula for arbitrary power sums \[\begin{align*} \sum_{k=1}^{n} k^p = \frac{1}{p+1} \sum_{r=0}^{p} \binom{p+1}{r} B_r n^{p+1-r}. \end{align*}\] Here it is assumed that \(B_1=+\frac{1}{2}\). The choice of convention for \(B_1\) has a long historical background; see (Luschny, Peter HN 2021, 2020) for a detailed discussion. Bernoulli’s formula supplied the missing symbolical framework for Faulhaber’s observations, thus the theory were complete. Although the expression above is commonly known as Faulhaber’s formula, its modern Bernoulli-number formulation originates from Bernoulli’s work.
During the eighteenth and nineteenth centuries, the theory was further developed by Euler, Stirling, Jacobi, and many others through connections with finite differences, interpolation, generating functions, and special number sequences.
In 1993, Donald Knuth revisited the theory of sums of powers and opened a new direction of development. Building upon ideas of Carl Gustav Jacobi (1804–1851) and the work of Gessel and Viennot (Gessel and Viennot 1989, (1989)), Knuth reexamined Faulhaber’s results through the lens of modern enumerative combinatorics. Remarkably, Bernoulli numbers are not central to Knuth’s approach. Instead, he derived formulas for power sums using special number sequences, including Stirling and central factorial numbers, together with binomial coefficients. Knuth also introduced the convenient notation \(\Sigma^{r}\,{n}^{m}\) for \(r\)-fold sums of powers, which will be used throughout this paper.
In this manuscript, we show the role of Newton’s interpolation formula in the derivation of identities for sums of powers. Our main result is an algorithm for obtaining closed-form expressions for power sums by combining variants of Newton’s interpolation formula with binomial identities from the hockey-stick family. The results demonstrate that infinitely many distinct closed-form representations of sums of powers can be constructed in this manner.
Allow us to define the mathematical objects we utilize throughout the manuscript. For multifold sums of powers, we use the notation provided by Donald Knuth in (Knuth, Donald E. 1993), because it is simple and beautiful.
Definition 1 (Multifold sums of powers). For non-negative integers \(r,n,m\) \[\begin{align*} \Sigma^{0}\,{n}^{m} &= n^m, \\ \Sigma^{1}\,{n}^{m} &= \Sigma^{0}\,{1}^{m} + \Sigma^{0}\,{2}^{m} + \cdots + \Sigma^{0}\,{n}^{m}, \\ \Sigma^{r+1}\,{n}^{m} &= \Sigma^{r}\,{1}^{m} + \Sigma^{r}\,{2}^{m} + \cdots + \Sigma^{r}\,{n}^{m}. \end{align*}\]
We utilize the following notation for falling factorials
Definition 2 (Falling factorials). For integers \(x,n\) \[\begin{align*} \left(n\right)_{k} = \begin{cases} 1, & \mathrm{if \; } k = 0, \\ n(n-1)(n-2)\cdots(n-k+1) = \prod_{j=0}^{k-1} (n-j), & \mathrm{if \; } k > 0, \\ 0, & \mathrm{else}. \end{cases} \end{align*}\]
We utilize the following notation for central factorials
Definition 3 (Central factorials). For integers \(n,k\) \[\begin{align*} n^{[k]} = \begin{cases} 0, \quad & \mathrm{if \; } k < 0, \\ 1, \quad & \mathrm{if \; } k = 0, \\ n \left( n + \frac{k}{2} -1 \right)\left( n + \frac{k}{2} -2 \right) \cdots \left( n - \frac{k}{2} +1 \right) = n \left(n+\frac{k}{2}-1\right)_{k-1}, \quad & \mathrm{if \; } k>0, \end{cases} \end{align*}\] where \(\left(n+\frac{k}{2}-1\right)_{k-1} = \prod_{j=0}^{k-2} \left( n+\frac{k}{2}-1 -j \right)\) are falling factorials.
We encourage the audience to read J. F. Steffensen’s work on the definition of generalized factorials (Steffensen 1933).
Lemma 4 (Central Newton’s formula). For a function \(f\) defined on integers and arbitrary integers \(x,t\) \[\begin{align*} f(x) = \sum_{k=0}^{\infty} \frac{(x-t)^{[k]}}{k!} \delta^{k} f(t), \end{align*}\] where \(\delta^{k} f(t)\) denotes the \(k\)-th central finite difference \[\begin{align*} \delta^{k} f(t) = \mathop{\textstyle\sum}_{j=0}^{k} (-1)^j \tbinom{k}{j} f\!\left(t+\tfrac{k}{2} - j\right), \end{align*}\] and \((x-t)^{[k]}\) are central factorials defined by [def:central-factorials].
Proof. See (Butzer et al. 1989, 462), and (Steffensen, Johan Frederik 1927, sec. 2, eq. (19)), and (Riordan 1968, sec. 6, eq. (21)). ◻
Lemma 5 (Newton’s formula for powers). For integers \(n,m \geq 0\), and an arbitrary integer \(t\), \[\begin{align*} n^m = \sum_{k=0}^{m} \frac{(n-t)^{[k]}}{k!} \delta^{k} t^m. \end{align*}\]
Proof. Special case of Lemma [lem:central-newtons-formula] for \(f(n) = n^m\). ◻
Lemma 6 (Central factorials binomial form). For integers \(n\) and \(k\neq0\), \[\begin{align*} \frac{n^{[k]}}{k!} = \frac{n}{k} \binom{n+\frac{k}{2}-1}{k-1}. \end{align*}\]
Proof. We have \[\begin{align*} \frac{n^{[k]}}{k!} &= \frac{n}{k!}\left(n+\frac{k}{2}-1\right)_{k-1} = \frac{n}{k (k-1)!} \left(n+\frac{k}{2}-1\right)_{k-1} = \frac{n}{k} \binom{n+\frac{k}{2}-1}{k-1}, \end{align*}\] because of the identity in falling factorials \(\frac{\left(x\right)_{n}}{n!} = \binom{x}{n}\). ◻
Lemma 7 (Central factorials binomial form decomposition). For non-negative integers \(n,k\), \[\begin{align*} \frac{n^{[k]}}{k!} = \frac{n}{k} \binom{n+\frac{k}{2}-1}{k-1} = \frac{1}{2} \left[ \binom{n+\frac{k}{2}}{k} + \binom{n+\frac{k}{2}-1}{k} \right]. \end{align*}\]
Proof. We have \(\binom{a}{k} = \frac{a}{k} \binom{a-1}{k-1}\). Let \(a=n + \frac{k}{2}\), then \[\begin{align*} \tbinom{n + \frac{k}{2}}{k} = \tfrac{n + \tfrac{k}{2}}{k} \tbinom{n+\frac{k}{2}-1}{k-1} = \tfrac{n}{k} \tbinom{n+\frac{k}{2}-1}{k-1} + \tfrac{1}{2} \tbinom{n+\frac{k}{2}-1}{k-1}. \end{align*}\] Thus, \[\begin{align*} \tfrac{n}{k} \tbinom{n+\frac{k}{2}-1}{k-1} = \tbinom{n + \frac{k}{2}}{k} - \tfrac{1}{2} \tbinom{n+\frac{k}{2}-1}{k-1}. \end{align*}\] By moving under common denominator, and applying binomial recurrence yields \[\begin{align*} \tfrac{n}{k} \tbinom{n+\frac{k}{2}-1}{k-1} = \tfrac{1}{2} \left[ 2 \tbinom{n + \frac{k}{2}}{k} - \tbinom{n+\frac{k}{2}-1}{k-1} \right] = \tfrac{1}{2} \left[ \tbinom{n + \frac{k}{2}}{k} + \tbinom{n + \frac{k}{2}-1}{k} + \tbinom{n + \frac{k}{2}-1}{k-1} - \tbinom{n+\frac{k}{2}-1}{k-1} \right]. \end{align*}\] Thus, the statement follows. ◻
Lemma 8 (Newton’s formula for powers decomposition). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} n^m = \frac{1}{2} \sum_{k=0}^{m} \left[ \binom{n-t+\frac{k}{2}}{k} + \binom{n-t+\frac{k}{2}-1}{k} \right] \delta^{k} t^m. \end{align*}\]
Proof. By Lemma [lem:central-newtons-formula-for-powers], and by Lemma [lem:central-factorials-binomial-form-decomposition]. ◻
Lemma 9 (Generalized hockey-stick identity). For arbitrary integers \(a,b\) and \(j\), \[\begin{align*} \sum_{k=a}^{b} \binom{k}{j} = \binom{b+1}{j+1} - \binom{a}{j+1}. \end{align*}\]
Proof. We have, \(\sum_{k=a}^{b} \binom{k}{j} = \binom{a}{j} + \binom{a+1}{j} + \cdots + \binom{b}{j} = \left( \sum_{k=0}^{b} \binom{k}{j} \right) - \left( \sum_{k=0}^{a-1} \binom{k}{j} \right)\). By hockey-stick identity \(\sum_{k=0}^{n} \binom{k}{j} = \binom{n+1}{j+1}\) yields, \[\begin{align*} \mathop{\textstyle\sum}_{k=a}^{b} \tbinom{k}{j} = \left( \mathop{\textstyle\sum}_{k=0}^{b} \tbinom{k}{j} \right) - \left( \mathop{\textstyle\sum}_{k=0}^{a-1} \tbinom{k}{j} \right) = \tbinom{b+1}{j+1} - \tbinom{a}{j+1}. \end{align*}\] This completes the proof. ◻
Lemma 10 (Centered hockey-stick identity). For integer \(n\geq1\), and arbitrary integers \(t,k,r\) \[\begin{align*} \sum_{j=1}^{n} \binom{j -t + \frac{k}{2} +r}{k} = \binom{n -t + \frac{k}{2} + r +1}{k+1} - \binom{1 -t + \tfrac{k}{2} +r}{k+1}. \end{align*}\]
Proof. By setting \(a=1 -t + \tfrac{k}{2} +r\), and \(b=n-t+ \tfrac{k}{2} +r\) into Lemma [lem:generalized-hockey-stick-identity] yields \[\begin{align*} \mathop{\textstyle\sum}_{j=1}^{n} \tbinom{j -t + \tfrac{k}{2} +r}{k} = \mathop{\textstyle\sum}_{j=1 -t + \tfrac{k}{2} +r}^{n -t + \tfrac{k}{2} +r} \tbinom{j}{k} = \tbinom{n -t + \tfrac{k}{2} + r +1}{k+1} - \tbinom{1 -t + \tfrac{k}{2} +r}{k+1}. \end{align*}\] This completes the proof. ◻
Lemma 11 (Multifold sums of zero powers). For non-negative integers \(r,n\) \[\begin{align*} \Sigma^{r}\,{n}^{0} = \binom{r+n-1}{r}. \end{align*}\]
Proof.
Let \(r=0\), then we have \(\Sigma^{0}\,{n}^{0} = 1\), by definition [def:multifold-sums-of-powers-recurrence].
Let \(r=1\), then we have \(\Sigma^{1}\,{n}^{0} = \sum_{k=1}^{n} 1 = \binom{n}{1}\).
Let \(r=2\), then we have \(\Sigma^{2}\,{n}^{0} = \sum_{k=1}^{n} \binom{k}{1} = \binom{n+1}{2}\).
Let \(r=3\), then we have \(\Sigma^{3}\,{n}^{0} = \sum_{k=1}^{n} \binom{k+1}{2} = \binom{n+2}{3}\).
Similarly, for arbitrary integer \(r\), by hockey-stick identity \(\sum_{k=1}^{n} \binom{k}{r} = \binom{n+1}{r+1}\), yields \(\Sigma^{r}\,{n}^{0} = \sum_{k=1}^{n} \Sigma^{r-1}\,{k}^{0} = \sum_{k=1}^{n} \binom{k+r-2}{r-1} = \binom{r+n-1}{r}\).
This completes the proof. ◻
Convention 12. For all \(x\) \[\begin{align*} x^0 = 1. \end{align*}\]
Donald Knuth give extensive discussion on the convention \(x^0=1\) for all \(x\), including zero, in Concrete Mathematics (Graham, Ronald L. and Knuth, Donald E. and Patashnik, Oren 1994, 162), and in (Knuth 1992).
By decomposition of Newton’s formula for powers [lem:central-newtons-formula-for-powers-decomposition], for non-negative integers \(n,m\), and an arbitrary integer \(t\), we have \[\begin{align*} \Sigma^{1}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \left[ \sum_{j=1}^{n} \binom{j-t+\frac{k}{2}}{k} + \sum_{j=1}^{n} \binom{j-t+\frac{k}{2}-1}{k} \right] \delta^{k} t^m. \end{align*}\] Now, the closed form of ordinary sums of powers becomes straightforward. By setting \(r=0\), and \(r=-1\) into Lemma [lem:centered-hockey-stick-identity], we get closed forms of the sums \(\sum_{j=1}^{n} \binom{j -t + \frac{k}{2}}{k}\), and \(\sum_{j=1}^{n} \binom{j -t + \frac{k}{2}-1}{k}\), respectively. Therefore,
Proposition 13 (Ordinary sums of powers decomposition). For integers \(n,m \geq 0\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{1}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \left[ \binom{n -t + \frac{k}{2} +1}{k+1} - \binom{1 -t + \tfrac{k}{2}}{k+1} + \binom{n -t + \frac{k}{2}}{k+1} - \binom{-t + \tfrac{k}{2}}{k+1} \right] \delta^{k} t^m. \end{align*}\]
Note that we do not bother ourselves about fractions in upper index of binomial coefficients \(\binom{n -t + \frac{k}{2} +1}{k+1}\), because binomial coefficient \(\binom{n}{k}\) is a polynomial in \(n\). The well-known identity in unsigned Stirling numbers of the first kind \(\genfrac{[}{]}{0pt}{}{n}{k}\) shows it explicitly \[\begin{align*} \binom{n}{k} = \sum_{j=0}^{k} \frac{(-1)^{k-j}}{k!} \genfrac{[}{]}{0pt}{}{k}{j} n^j. \end{align*}\] For example, for \(\Sigma^{1}\,{n}^{2}\) we have
Example 14 (Sums of squares in central differences). \[\begin{align*} t=0:\quad \Sigma^{1}\,{n}^{2} &= \frac{t^2}{2}(2n) + \frac{\delta t^2}{2}\, n(1+n) + \frac{\delta^2 t^2}{2}\,\frac{n(1+3n+2n^2)}{6}, \\ % t=1:\quad \Sigma^{1}\,{n}^{2} &= \frac{t^2}{2}(2n) + \frac{\delta t^2}{2}\, n(n-1) + \frac{\delta^2 t^2}{2}\,\frac{n(1-3n+2n^2)}{6}, \\ % t=2:\quad \Sigma^{1}\,{n}^{2} &= \frac{t^2}{2}(2n) + \frac{\delta t^2}{2}\, n(n-3) + \frac{\delta^2 t^2}{2}\,\frac{n(13-9n+2n^2)}{6}, \\ % t=3:\quad \Sigma^{1}\,{n}^{2} &= \frac{t^2}{2}(2n) + \frac{\delta t^2}{2}\, n(n-5) + \frac{\delta^2 t^2}{2}\,\frac{n(37-15n+2n^2)}{6}, \\ % t=4:\quad \Sigma^{1}\,{n}^{2} &= \frac{t^2}{2}(2n) + \frac{\delta t^2}{2}\, n(n-7) + \frac{\delta^2 t^2}{2}\,\frac{n(73-21n+2n^2)}{6}. \end{align*}\]
For example, for \(\Sigma^{1}\,{n}^{3}\) we have
Example 15 (Sums of cubes in central differences). \[\begin{align*} t=0:\quad \Sigma^{1}\,{n}^{3} &= \tfrac{t^3}{2}(2n) + \tfrac{\delta t^3}{2}\, n(1+n) + \tfrac{\delta^2 t^3}{2}\,\tfrac{n(1+3n+2n^2)}{6} + \tfrac{\delta^3 t^3}{2}\,\tfrac{n(-1+n+4n^2+2n^3)}{24}, \\ % t=1:\quad \Sigma^{1}\,{n}^{3} &= \tfrac{t^3}{2}(2n) + \tfrac{\delta t^3}{2}\, n(n-1) + \tfrac{\delta^2 t^3}{2}\,\tfrac{n(1-3n+2n^2)}{6} + \tfrac{\delta^3 t^3}{2}\,\tfrac{n(1+n-4n^2+2n^3)}{24}, \\ % t=2:\quad \Sigma^{1}\,{n}^{3} &= \tfrac{t^3}{2}(2n) + \tfrac{\delta t^3}{2}\, n(n-3) + \tfrac{\delta^2 t^3}{2}\,\tfrac{n(13-9n+2n^2)}{6} + \tfrac{\delta^3 t^3}{2}\,\tfrac{n(-21+25n-12n^2+2n^3)}{24}, \\ % t=3:\quad \Sigma^{1}\,{n}^{3} &= \tfrac{t^3}{2}(2n) + \tfrac{\delta t^3}{2}\, n(n-5) + \tfrac{\delta^2 t^3}{2}\,\tfrac{n(37-15n+2n^2)}{6} + \tfrac{\delta^3 t^3}{2}\,\tfrac{n(-115+73n-20n^2+2n^3)}{24}, \\ % t=4:\quad \Sigma^{1}\,{n}^{3} &= \tfrac{t^3}{2}(2n) + \tfrac{\delta t^3}{2}\, n(n-7) + \tfrac{\delta^2 t^3}{2}\,\tfrac{n(73-21n+2n^2)}{6} + \tfrac{\delta^3 t^3}{2}\,\tfrac{n(-329+145n-28n^2+2n^3)}{24}. \end{align*}\]
It is indeed interesting to notice the connection between Proposition [prop:ordinary-sums-of-powers-decomposition], and central factorial numbers of the second kind \(T(n,k)\)
Proposition 16 (Central factorial numbers sums). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{1}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \left[ \binom{n + \frac{k}{2} +1}{k+1} - \binom{1 + \tfrac{k}{2}}{k+1} + \binom{n + \frac{k}{2}}{k+1} - \binom{\tfrac{k}{2}}{k+1} \right] k! T(m,k), \end{align*}\] where \(T(m,k)\) are central factorial numbers of the second kind (in Riordan’s notation). See (Riordan 1968, 213), and (Carlitz and Riordan 1963).
These central factorial numbers of the second kind \(T(n,k)\) are defined by central Newton’s formula for powers [lem:central-newtons-formula-for-powers] evaluated at zero, such that \[\begin{align*} n^m = \sum_{k=0}^{\infty} T(m,k) n^{[k]}. \end{align*}\] Which implies \[\begin{align*} T(n,k) = \frac{\delta^{k} 0^m}{k!}, \end{align*}\] by Newton’s formula [lem:central-newtons-formula-for-powers].
In (Knuth, Donald E. 1993) Donald Knuth provides a remarkable formula for sums of odd powers, in terms of central factorial numbers of the second kind \(T(n,k)\)
Proposition 17 (Knuth’s odd powers sum). For integers \(n \geq 0\), and \(m\geq 1\) \[\begin{align*} \Sigma^{1}\,{n}^{2m-1} = \sum_{k=1}^{m} (2k-1)! \binom{n+k}{2k} T(2m,2k), \end{align*}\] where \(T(2m,2k)\) are central factorial numbers of the second kind.
It is quite interesting to notice that Proposition [prop:knuth-odd-powers-sum] originates from centered Newton’s formula for powers [lem:central-newtons-formula-for-powers] evaluated at zero. By Newton’s formula, for \(m\geq 1\), we have \[\begin{align*} n^{2m} = \sum_{k=1}^{2m} \frac{n^{[k]}}{k!} \cdot \delta^{k} 0^{2m} \end{align*}\] The iteration for \(k=0\) can be omitted because \(\delta^{0} 0^{2m} = 0\). Now, we observe that central difference \(\delta^k 0^m\) satisfies the parity of arguments \(m\) and \(k\), such that \[\begin{align*} \begin{cases} \delta^{k} 0^m \neq 0, \quad &\text{when} \quad m \equiv k \pmod{2}, \\ \delta^{k} 0^m = 0, \quad &\text{when} \quad m \not\equiv k \pmod{2}, \end{cases} \end{align*}\] meaning that \(\delta^{k} 0^m\) is zero whether \(k,m\) are odd and even, and vise versa. Note that the parity property holds only for central difference of powers evaluated at zero, by contrast \(\delta^{3} 1^{4} = 24 \neq 0\). Thus, the parity of \(\delta^{k} 0^{2m}\) allows us to omit odd iterations of \(k\), because \(\delta^{k} 0^{2m}\) is zero for odd \(k\) \[\begin{align*} n^{2m} = \sum_{k=1}^{m} \frac{n^{[2k]}}{(2k)!} \cdot \delta^{2k} 0^{2m}. \end{align*}\] By definition of central factorial numbers of the second kind, yields \[\begin{align*} n^{2m} = \sum_{k=1}^{m} n^{[2k]} T(2m,2k). \end{align*}\] Thus, \[\begin{align*} n^{2m} = \sum_{k=1}^{m} (2k)! \frac{n^{[2k]}}{(2k)!} T(2m,2k) \end{align*}\] By simplifying \(\frac{n^{[2k]}}{(2k)!}\) using Lemma [lem:central-factorials-binomial-form], we get \[\begin{align*} n^{2m} = \sum_{k=1}^{m} (2k)! \frac{n}{2k} \binom{n+k-1}{2k-1} T(2m,2k). \end{align*}\] Therefore, \[\begin{align*} n^{2m} = \sum_{k=1}^{m} (2k-1)! n \binom{n+k-1}{2k-1} T(2m,2k). \end{align*}\] By factoring out \(n\), we obtain \[\begin{align*} n^{2m-1} &= \sum_{k=1}^{m} (2k-1)! \binom{n+k-1}{2k-1} T(2m,2k). \end{align*}\] By using hockey-stick identity \(\sum_{k=0}^{n} \binom{k}{j} = \binom{n+1}{j+1}\) yields \[\begin{align*} \Sigma^{1}\,{n}^{2m-1} = \sum_{k=1}^{m} (2k-1)! \binom{n+k}{2k} T(2m,2k). \end{align*}\] The formula above matches the Proposition [prop:knuth-odd-powers-sum] precisely. Clearly, the formula for sums of odd powers can be generalized to multifold case, by applying hockey-stick identity \(\sum_{k=0}^{n} \binom{k}{j} = \binom{n+1}{j+1}\) repeatedly
Proposition 18 (Knuth’s odd powers sum multifold). For non-negative integers \(r,n\), and \(m \geq 1\) \[\begin{align*} \Sigma^{r+1}\,{n}^{2m-1} = \sum_{k=1}^{m} (2k-1)! \binom{n+k+r}{2k+r} T(2m,2k), \end{align*}\] where \(T(2m,2k)\) are central factorial numbers of the second kind.
For example, for \(r=0\) we have
\[\begin{align*}
\Sigma^{1}\,{n}^{1} &= \tbinom{n+1}{2} 1, \\
\Sigma^{1}\,{n}^{3} &= \tbinom{n+2}{4} 6 + \tbinom{n+1}{2} 1,
\\
\Sigma^{1}\,{n}^{5} &= \tbinom{n+3}{6} 120 + \tbinom{n+2}{4}
30 + \tbinom{n+1}{2} 1, \\
\Sigma^{1}\,{n}^{7} &= \tbinom{n+4}{8} 5040 + \tbinom{n+3}{6} 1680
+ \tbinom{n+2}{4} 126 + \tbinom{n+1}{2} 1.
\end{align*}\] For example, for arbitrary integer \(r\) we have \[\begin{align*}
\Sigma^{r}\,{n}^{1} &= \tbinom{n+0+r}{1+r} 1, \\
\Sigma^{r}\,{n}^{3} &= \tbinom{n+1+r}{3+r} 6 +
\tbinom{n+0+r}{1+r} 1, \\
\Sigma^{r}\,{n}^{5} &= \tbinom{n+2+r}{5+r} 120 +
\tbinom{n+1+r}{3+r} 30 + \tbinom{n+0+r}{1+r} 1, \\
\Sigma^{r}\,{n}^{7} &= \tbinom{n+3+r}{7+r} 5040 +
\tbinom{n+3+r}{5+r} 1680 + \tbinom{n+1+r}{3+r} 126 + \tbinom{n+0+r}{1+r}
1.
\end{align*}\] The coefficients \(1, 6,
1, 120, 30, 1, 5040, 1680,\ldots\) is the sequence A303675 in the
OEIS (Sloane, Neil
J.A. and others 2003). However, Proposition [prop:knuth-odd-powers-sum-multifold]
serves as the source of inspiration for one more remarkable identity. By
Newton’s formula [lem:central-newtons-formula-for-powers],
and by Lemma [lem:central-factorials-binomial-form],
we have \[\begin{align*}
n^m = \sum_{k=1}^{m} \frac{n}{k} \binom{n+\frac{k}{2}-1}{k-1}
\delta^{k} 0^m
\end{align*}\] By factoring out \(n\) \[\begin{align*}
n^{m-1} = \sum_{k=1}^{m} \frac{1}{k} \binom{n+\frac{k}{2}-1}{k-1}
\delta^{k} 0^m
\end{align*}\] By shifting \[\begin{align*}
n^{m} = \sum_{k=1}^{m+1} \frac{1}{k} \binom{n+\frac{k}{2}-1}{k-1}
\delta^{k} 0^{m+1}
\end{align*}\] By re-indexing \[\begin{align}
\label{eq:power-identity-at-zero}
n^{m} = \sum_{k=0}^{m} \frac{1}{k+1}
\binom{n+\frac{k}{2}-\frac{1}{2}}{k} \delta^{k+1} 0^{m+1}
\end{align}\] Thus, by setting \(r=0\) and \(t=\frac{1}{2}\) into Lemma [lem:centered-hockey-stick-identity]
yields
Proposition 19 (Ordinary sums of powers at zero). For non-negative integers \(n,m\) \[\begin{align*} \Sigma^{1}\,{n}^{m} = \sum_{k=0}^{m} \frac{1}{k+1} \left[ \binom{n+\frac{k}{2}+\frac{1}{2}}{k+1} - \binom{\frac{k}{2}+\frac{1}{2}}{k+1}\right] \delta^{k+1} 0^{m+1}. \end{align*}\]
Which can also be rewritten in terms of central factorial numbers of the second kind \[\begin{align*} \Sigma^{1}\,{n}^{m} = \sum_{k=0}^{m} k! \left[ \binom{n+\frac{k}{2}+\frac{1}{2}}{k+1} - \binom{\frac{k}{2}+\frac{1}{2}}{k+1}\right] T(k+1, m+1). \end{align*}\] By setting \(r=1\) and \(t=\frac{1}{2}\) into Lemma [lem:centered-hockey-stick-identity], we get formula for double sums of powers
Proposition 20 (Double sums of powers at zero). For non-negative integers \(n,m\) \[\begin{align*} \Sigma^{2}\,{n}^{m} = \sum_{k=0}^{m} \frac{1}{k+1} \left[ \binom{n+\frac{k}{2}+\frac{3}{2}}{k+2} - \binom{\frac{k}{2}+\frac{1}{2}}{k+1} \Sigma^{1}\,{n}^{0} - \binom{\frac{k}{2}+\frac{3}{2}}{k+2} \Sigma^{0}\,{n}^{0} \right] \delta^{k+1} 0^{m+1}. \end{align*}\]
In general, we have
Theorem 21 (Multifold sums of powers at zero). For non-negative integers \(r,n,m\) \[\begin{align*} \Sigma^{r}\,{n}^{m} = \sum_{k=0}^{m} \left[ \binom{n+\frac{k}{2}+\frac{2r-1}{2}}{k+r} - \sum_{s=0}^{r-1} \binom{\frac{k}{2} + \frac{2(r-s)-1}{2}}{k+r-s} \Sigma^{s}\,{n}^{0} \right] \frac{\delta^{k+1} 0^{m+1}}{k+1}. \end{align*}\]
Proof. By repeatedly applied Lemma [lem:centered-hockey-stick-identity], and Equation [eq:power-identity-at-zero]. ◻
By Lemma [lem:multifold-sum-of-zero-powers], from Theorem [thm:multifold-sums-of-powers-at-zero] follows
Proposition 22 (Multifold sums of powers at zero binomial form). For non-negative integers \(r,n,m\) \[\begin{align*} \Sigma^{r}\,{n}^{m} = \sum_{k=0}^{m} \left[ \binom{n+\frac{k}{2}+\frac{2r-1}{2}}{k+r} - \sum_{s=0}^{r-1} \binom{\frac{k}{2} + \frac{2(r-s)-1}{2}}{k+r-s} \binom{s+n-1}{s} \right] \frac{\delta^{k+1} 0^{m+1}}{k+1}. \end{align*}\]
Proof. By repeatedly applied Lemma [lem:centered-hockey-stick-identity], and Equation [eq:power-identity-at-zero]. ◻
Which, as well, can be expressed in terms of central factorial numbers of the second kind \[\begin{align*} \Sigma^{r}\,{n}^{m} = \sum_{k=0}^{m} \left[ \binom{n+\frac{k}{2}+\frac{2r-1}{2}}{k+r} - \sum_{s=0}^{r-1} \binom{\frac{k}{2} + \frac{2(r-s)-1}{2}}{k+r-s} \binom{s+n-1}{s} \right] k! T(k+1, m+1). \end{align*}\] Thus, we have discussed quite remarkable results, that follow from centered Newton’s formula evaluated at zero.
Definition 23 (Generalized central factorial numbers of the second kind). For integers \(0 \leq k \leq m\), and an arbitrary integer \(t\), the generalized central factorial numbers of the second kind are defined by the relation \[\begin{align*} n^m = \sum_{k=0}^{\infty} T_t(m,k) \cdot (n-t)^{[k]}. \end{align*}\]
It is also straightforward that \[\begin{align*} T_t (m,k) = \frac{\delta^{k} t^m}{k!}, \end{align*}\] by Newton’s formula [lem:central-newtons-formula-for-powers]. For example, for \(T_0(n,k)\) we have
| \(n/k\) | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | ||||||||
| 1 | 0 | 1 | |||||||
| 2 | 0 | 0 | 1 | ||||||
| 3 | 0 | \(\tfrac{1}{4}\) | 0 | 1 | |||||
| 4 | 0 | 0 | 1 | 0 | 1 | ||||
| 5 | 0 | \(\tfrac{1}{16}\) | 0 | \(\tfrac{5}{2}\) | 0 | 1 | |||
| 6 | 0 | 0 | 1 | 0 | 5 | 0 | 1 | ||
| 7 | 0 | \(\tfrac{1}{64}\) | 0 | \(\tfrac{91}{16}\) | 0 | \(\tfrac{35}{4}\) | 0 | 1 | |
| 8 | 0 | 0 | 1 | 0 | 21 | 0 | 14 | 0 | 1 |
For example, for \(T_1(n,k)\) we have
| \(n/k\) | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | ||||||||
| 1 | 1 | 1 | |||||||
| 2 | 1 | 2 | 1 | ||||||
| 3 | 1 | \(\tfrac{13}{4}\) | 3 | 1 | |||||
| 4 | 1 | 5 | 7 | 4 | 1 | ||||
| 5 | 1 | \(\tfrac{121}{16}\) | 15 | \(\tfrac{25}{2}\) | 5 | 1 | |||
| 6 | 1 | \(\tfrac{91}{8}\) | 31 | 35 | 20 | 6 | 1 | ||
| 7 | 1 | \(\tfrac{1093}{64}\) | 63 | \(\tfrac{1491}{16}\) | 70 | \(\tfrac{119}{4}\) | 7 | 1 | |
| 8 | 1 | \(\tfrac{205}{8}\) | 127 | \(\tfrac{483}{2}\) | 231 | 126 | 42 | 8 | 1 |
For example, for \(T_2(n,k)\) we have
| \(n/k\) | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | ||||||||
| 1 | 2 | 1 | |||||||
| 2 | 4 | 4 | 1 | ||||||
| 3 | 8 | \(\tfrac{49}{4}\) | 6 | 1 | |||||
| 4 | 16 | 34 | 25 | 8 | 1 | ||||
| 5 | 32 | \(\tfrac{1441}{16}\) | 90 | \(\tfrac{85}{2}\) | 10 | 1 | |||
| 6 | 64 | \(\tfrac{931}{4}\) | 301 | 190 | 65 | 12 | 1 | ||
| 7 | 128 | \(\tfrac{37969}{64}\) | 966 | \(\tfrac{12411}{16}\) | 350 | \(\tfrac{371}{4}\) | 14 | 1 | |
| 8 | 256 | \(\tfrac{6001}{4}\) | 3025 | 3003 | 1701 | 588 | 126 | 16 | 1 |
For example, for \(T_3(n,k)\) we have
| \(n/k\) | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 1 | ||||||||
| 1 | 3 | 1 | |||||||
| 2 | 9 | 6 | 1 | ||||||
| 3 | 27 | \(\tfrac{109}{4}\) | 9 | 1 | |||||
| 4 | 81 | 111 | 55 | 12 | 1 | ||||
| 5 | 243 | \(\tfrac{6841}{16}\) | 285 | \(\tfrac{185}{2}\) | 15 | 1 | |||
| 6 | 729 | \(\tfrac{12753}{8}\) | 1351 | 585 | 140 | 18 | 1 | ||
| 7 | 2187 | \(\tfrac{372709}{64}\) | 6069 | \(\tfrac{53011}{16}\) | 1050 | \(\tfrac{791}{4}\) | 21 | 1 | |
| 8 | 6561 | \(\tfrac{167943}{8}\) | 26335 | \(\tfrac{35049}{2}\) | 6951 | 1722 | 266 | 24 | 1 |
For example, for \(T_t (2n, 2k)\), such that \(0 \leq t \leq 3\), we have
| \(t\) | \(T_t(2m,2k)\) | OEIS ID |
|---|---|---|
| \(0\) | \(1, 0, 1, 0, 1, 1, 0, 1, 5, 1, 0, 1, 21, 14, 1, 0, 1, 85, 147, 30, 1, 0, 1, \cdots\) | A269945 |
| \(1\) | \(1, 1, 1, 1, 7, 1, 1, 31, 20, 1, 1, 127, 231, 42, 1, 1, 511, 2290, 987,\cdots\) | A394692 |
| \(2\) | \(1, 4, 1, 16, 25, 1, 64, 301, 65, 1, 256, 3025, 1701, 126, 1, 1024, \cdots\) | A395456 |
| \(3\) | \(1, 9, 1, 81, 55, 1, 729, 1351, 140, 1, 6561, 26335, 6951, 266, 1, \cdots\) | A395457 |
Thus, yields the formula for sums of powers in terms of generalized central factorial numbers of the second kind
Proposition 24 (Generalized Central factorial sums). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{1}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \left[ \binom{n -t + \frac{k}{2} +1}{k+1} - \binom{1 -t + \tfrac{k}{2}}{k+1} + \binom{n -t + \frac{k}{2}}{k+1} - \binom{-t + \tfrac{k}{2}}{k+1} \right] k! T_t (m,k). \end{align*}\]
By following similarly, it is easy to compute closed forms of double sums of powers, by using Lemma [lem:centered-hockey-stick-identity] along with Proposition [prop:ordinary-sums-of-powers-decomposition]. Thus, we have to compute closed forms for sums \(\sum_{j=1}^{n} \binom{j -t + \frac{k}{2} +1}{k+1}\), and \(\sum_{j=1}^{n} \binom{j -t + \frac{k}{2}}{k+1}\), hence
Proposition 25 (Double sums of powers decomposition). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{2}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \Bigg[ &+ \binom{n -t + \frac{k}{2} +2}{k+2} - \binom{2 -t + \tfrac{k}{2}}{k+2} \Sigma^{0}\,{n}^{0} - \binom{1 -t + \tfrac{k}{2}}{k+1} \Sigma^{1}\,{n}^{0} \\ &+ \binom{n -t + \frac{k}{2}+1}{k+2} - \binom{1-t + \tfrac{k}{2}}{k+2} \Sigma^{0}\,{n}^{0} - \binom{0-t + \tfrac{k}{2}}{k+1} \Sigma^{1}\,{n}^{0} \Bigg] \delta^{k} t^m. \end{align*}\]
For example, for \(\Sigma^{2}\,{n}^{2}\) we have
Example 26 (Double sum of squares in central differences). For integer \(0 \leq t \leq 4\) \[\begin{align*} t=0:\quad \Sigma^{2}\,{n}^{2} &= \frac{t^2}{2}\, n(1+n) + \frac{\delta t^2}{2}\,\frac{n(2+3n+n^2)}{3} + \frac{\delta^2 t^2}{2}\,\frac{n(1+n)^2(2+n)}{12}, \\ % t=1:\quad \Sigma^{2}\,{n}^{2} &= \frac{t^2}{2}\, n(1+n) + \frac{\delta t^2}{2}\,\frac{n(-1+n^2)}{3} + \frac{\delta^2 t^2}{2}\,\frac{n^2(-1+n^2)}{12}, \\ % t=2:\quad \Sigma^{2}\,{n}^{2} &= \frac{t^2}{2}\, n(1+n) + \frac{\delta t^2}{2}\,\frac{n(-4-3n+n^2)}{3} + \frac{\delta^2 t^2}{2}\,\frac{n(10+5n-4n^2+n^3)}{12}, \\ % t=3:\quad \Sigma^{2}\,{n}^{2} &= \frac{t^2}{2}\, n(1+n) + \frac{\delta t^2}{2}\,\frac{n(-7-6n+n^2)}{3} + \frac{\delta^2 t^2}{2}\,\frac{n(32+23n-8n^2+n^3)}{12}, \\ % t=4:\quad \Sigma^{2}\,{n}^{2} &= \frac{t^2}{2}\, n(1+n) + \frac{\delta t^2}{2}\,\frac{n(-10-9n+n^2)}{3} + \frac{\delta^2 t^2}{2}\,\frac{n(66+53n-12n^2+n^3)}{12}. \end{align*}\]
Repeating the procedure above, we easily obtain the formula for triple sums of powers
Proposition 27 (Triple sums of powers decomposition). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} &\Sigma^{3}\,{n}^{m} = \frac{1}{2} \sum_{k=0}^{m} \Bigg[ \\ &+ \binom{n -t + \frac{k}{2} +3}{k+3} - \binom{3 -t + \tfrac{k}{2}}{k+3} \Sigma^{0}\,{n}^{0} - \binom{2 -t + \tfrac{k}{2}}{k+2} \Sigma^{1}\,{n}^{0} - \binom{1 -t + \tfrac{k}{2}}{k+1} \Sigma^{2}\,{n}^{0} \\ &+ \binom{n -t + \frac{k}{2}+2}{k+3} - \binom{2-t + \tfrac{k}{2}}{k+3} \Sigma^{0}\,{n}^{0} - \binom{1-t + \tfrac{k}{2}}{k+2} \Sigma^{1}\,{n}^{0} - \binom{0-t + \tfrac{k}{2}}{k+1} \Sigma^{2}\,{n}^{0} \Bigg] \delta^{k} t^m. \end{align*}\]
For example, for \(\Sigma^{3}\,{n}^{2}\) we have
Example 28 (Triple sum of squares in central differences). For integer \(0 \leq t \leq 4\) \[\begin{align*} t=0:\quad \Sigma^{3}\,{n}^{2} &= \tfrac{t^2}{2}\,\tfrac{n(2+3n+n^2)}{3} + \tfrac{\delta t^2}{2}\,\tfrac{n(6+11n+6n^2+n^3)}{12} + \tfrac{\delta^2 t^2}{2}\,\tfrac{n(18+45n+40n^2+15n^3+2n^4)}{120}, \\ % t=1:\quad \Sigma^{3}\,{n}^{2} &= \tfrac{t^2}{2}\,\tfrac{n(2+3n+n^2)}{3} + \tfrac{\delta t^2}{2}\,\tfrac{n(-2-n+2n^2+n^3)}{12} + \tfrac{\delta^2 t^2}{2}\,\tfrac{n(-2-5n+5n^3+2n^4)}{120}, \\ % t=2:\quad \Sigma^{3}\,{n}^{2} &= \tfrac{t^2}{2}\,\tfrac{n(2+3n+n^2)}{3} + \tfrac{\delta t^2}{2}\,\tfrac{n(-10-13n-2n^2+n^3)}{12} + \tfrac{\delta^2 t^2}{2}\,\tfrac{n(58+65n-5n^3+2n^4)}{120}, \\ % t=3:\quad \Sigma^{3}\,{n}^{2} &= \tfrac{t^2}{2}\,\tfrac{n(2+3n+n^2)}{3} + \tfrac{\delta t^2}{2}\,\tfrac{n(-18-25n-6n^2+n^3)}{12} + \tfrac{\delta^2 t^2}{2}\,\tfrac{n(198+255n+40n^2-15n^3+2n^4)}{120}, \\ % t=4:\quad \Sigma^{3}\,{n}^{2} &= \tfrac{t^2}{2}\,\tfrac{n(2+3n+n^2)}{3} + \tfrac{\delta t^2}{2}\,\tfrac{n(-26-37n-10n^2+n^3)}{12} + \tfrac{\delta^2 t^2}{2}\,\tfrac{n(418+565n+120n^2-25n^3+2n^4)}{120}. \end{align*}\]
By applying Lemma [lem:centered-hockey-stick-identity] on formula for ordinary sums of powers [prop:ordinary-sums-of-powers-decomposition] repeatedly, the formula for multifold sums of powers follows
Theorem 29 (Multifold sums of powers). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \frac{1}{2} \sum_{k=0}^{m} \Bigg \{ \binom{n-t+\frac{k}{2}+r}{k+r} + \binom{n-t+\frac{k}{2}+r-1}{k+r} \\ &- \sum_{s=0}^{r-1} \left[ \binom{r-s-t+\frac{k}{2}}{k+r-s} + \binom{r-1-s-t+\frac{k}{2}}{k+r-s} \right] \Sigma^{s}\,{n}^{0} \Bigg \} \delta^{k} t^m. \end{align*}\]
For example,
Example 30 (Examples of multifold sums). For integer \(0 \leq t \leq 2\) \[\begin{align*} t=0:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n+\tfrac{k}{2}+r}{k+r} + \tbinom{n+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s+\tfrac{k}{2}}{k+r-s} + \tbinom{r-1-s+\tfrac{k}{2}}{k+r-s} \Big] \Sigma^{s}\,{n}^{0} \Bigg\} \,\delta^k 0^m , \\ % t=1:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n-1+\tfrac{k}{2}+r}{k+r} + \tbinom{n-1+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s-1+\tfrac{k}{2}}{k+r-s} + \tbinom{r-2-s+\tfrac{k}{2}}{k+r-s} \Big] \Sigma^{s}\,{n}^{0} \Bigg\} \,\delta^k 1^m , \\ % t=2:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n-2+\tfrac{k}{2}+r}{k+r} + \tbinom{n-2+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s-2+\tfrac{k}{2}}{k+r-s} + \tbinom{r-3-s+\tfrac{k}{2}}{k+r-s} \Big] \Sigma^{s}\,{n}^{0} \Bigg\} \,\delta^k 2^m . \end{align*}\]
By definition [def:generalized-central-factorial-numbers], from Theorem [thm:multifold-sums-of-powers] follows
Proposition 31 (Multifold central factorial sums of powers). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \frac{1}{2} \sum_{k=0}^{m} \Bigg \{ \binom{n-t+\frac{k}{2}+r}{k+r} + \binom{n-t+\frac{k}{2}+r-1}{k+r} \\ &- \sum_{s=0}^{r-1} \left[ \binom{r-s-t+\frac{k}{2}}{k+r-s} + \binom{r-1-s-t+\frac{k}{2}}{k+r-s} \right] \Sigma^{s}\,{n}^{0} \Bigg \} k! T_t(m,k). \end{align*}\]
For example, for \(\Sigma^{1}\,{n}^{2}\) we have
Example 32 (Sums of squares in central factorial numbers). For integer \(0 \leq t \leq 5\) \[\begin{align*} t=0:\quad \Sigma^{1}\,{n}^{2} &= n\cdot0 + \tfrac12 n(1+n)\cdot0 + \tfrac16 n(1+3n+2n^2)\cdot1, \\ % t=1:\quad \Sigma^{1}\,{n}^{2} &= n\cdot1 + \tfrac12 n(n-1)\cdot2 + \tfrac16 n(1-3n+2n^2)\cdot1, \\ % t=2:\quad \Sigma^{1}\,{n}^{2} &= n\cdot4 + \tfrac12 n(n-3)\cdot4 + \tfrac16 n(13-9n+2n^2)\cdot1, \\ % t=3:\quad \Sigma^{1}\,{n}^{2} &= n\cdot9 + \tfrac12 n(n-5)\cdot6 + \tfrac16 n(37-15n+2n^2)\cdot1, \\ % t=4:\quad \Sigma^{1}\,{n}^{2} &= n\cdot16 + \tfrac12 n(n-7)\cdot8 + \tfrac16 n(73-21n+2n^2)\cdot1, \\ % t=5:\quad \Sigma^{1}\,{n}^{2} &= n\cdot25 + \tfrac12 n(n-9)\cdot10 + \tfrac16 n(121-27n+2n^2)\cdot1. \end{align*}\]
From the example above, we may observe that
For \(t=0\) the coefficients \(0,0,1\) correspond to the second row of the triangle [tab:central-factorial-numbers-t0].
For \(t=1\) the coefficients \(1,2,1\) correspond to the second row of the triangle [tab:central-factorial-numbers-t1].
For \(t=2\) the coefficients \(4,4,1\) correspond to the second row of the triangle [tab:central-factorial-numbers-t2].
For \(t=3\) the coefficients \(9,6,1\) correspond to the second row of the triangle [tab:central-factorial-numbers-t3].
For example, for \(\Sigma^{1}\,{n}^{4}\) we have
Example 33 (Sums of fourth powers in central factorial numbers). For integer \(0 \leq t \leq 3\) \[\begin{align*} t=0:\quad \Sigma^{1}\,{n}^{4} &= n\cdot0 + \tfrac12 n(1+n)\cdot0 + \tfrac16 n(1+3n+2n^2)\cdot1 \\ &+ \tfrac18 n(-1+n+4n^2+2n^3)\cdot0 + \tfrac1{10} n(-2-5n+5n^3+2n^4)\cdot1, \\ % t=1:\quad \Sigma^{1}\,{n}^{4} &= n\cdot1 + \tfrac12 n(n-1)\cdot5 + \tfrac16 n(1-3n+2n^2)\cdot7 \\ &+ \tfrac18 n(1+n-4n^2+2n^3)\cdot4 + \tfrac1{10} n(-2+5n-5n^3+2n^4)\cdot1, \\ % t=2:\quad \Sigma^{1}\,{n}^{4} &= n\cdot16 + \tfrac12 n(n-3)\cdot34 + \tfrac16 n(13-9n+2n^2)\cdot25 \\ &+ \tfrac18 n(-21+25n-12n^2+2n^3)\cdot8 \\ &+ \tfrac1{10} n(18-45n+40n^2-15n^3+2n^4)\cdot1, \\ % t=3:\quad \Sigma^{1}\,{n}^{4} &= n\cdot81 + \tfrac12 n(n-5)\cdot111 + \tfrac16 n(37-15n+2n^2)\cdot55 \\ &+ \tfrac18 n(-115+73n-20n^2+2n^3)\cdot12 \\ &+ \tfrac1{10} n(298-275n+120n^2-25n^3+2n^4)\cdot1. \end{align*}\]
From the example above, we may observe that
For \(t=0\) the coefficients \(0,0,1,0,1\) correspond to the fourth row of the triangle [tab:central-factorial-numbers-t0].
For \(t=1\) the coefficients \(1,5,7,4,1\) correspond to the fourth row of the triangle [tab:central-factorial-numbers-t1].
For \(t=2\) the coefficients \(16,34,25,8,1\) correspond to the fourth row of the triangle [tab:central-factorial-numbers-t2].
For \(t=3\) the coefficients \(81,111,55,12,1\) correspond to the fourth row of the triangle [tab:central-factorial-numbers-t3].
By Lemma [lem:multifold-sum-of-zero-powers], and by Theorem [thm:multifold-sums-of-powers] yields
Proposition 34 (Multifold binomial sums of powers). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \frac{1}{2} \sum_{k=0}^{m} \Bigg \{ \binom{n-t+\frac{k}{2}+r}{k+r} + \binom{n-t+\frac{k}{2}+r-1}{k+r} \\ &- \sum_{s=0}^{r-1} \left[ \binom{r-s-t+\frac{k}{2}}{k+r-s} + \binom{r-1-s-t+\frac{k}{2}}{k+r-s} \right] \binom{s+n-1}{s} \Bigg \} \delta^{k} t^m. \end{align*}\]
For example,
Example 35 (Examples of multifold binomial sums). For integer \(0 \leq t \leq 2\) \[\begin{align*} t=0:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n+\tfrac{k}{2}+r}{k+r} + \tbinom{n+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s+\tfrac{k}{2}}{k+r-s} + \tbinom{r-1-s+\tfrac{k}{2}}{k+r-s} \Big] \tbinom{s+n-1}{s} \Bigg\} \,\delta^k 0^m , \\ % t=1:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n-1+\tfrac{k}{2}+r}{k+r} + \tbinom{n-1+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s-1+\tfrac{k}{2}}{k+r-s} + \tbinom{r-2-s+\tfrac{k}{2}}{k+r-s} \Big] \tbinom{s+n-1}{s} \Bigg\} \,\delta^k 1^m , \\ % t=2:\quad \Sigma^{r}\,{n}^{m} &= \tfrac12 \mathop{\textstyle\sum}_{k=0}^{m} \Bigg\{ \tbinom{n-2+\tfrac{k}{2}+r}{k+r} + \tbinom{n-2+\tfrac{k}{2}+r-1}{k+r} \\ &\qquad- \mathop{\textstyle\sum}_{s=0}^{r-1} \Big[ \tbinom{r-s-2+\tfrac{k}{2}}{k+r-s} + \tbinom{r-3-s+\tfrac{k}{2}}{k+r-s} \Big] \tbinom{s+n-1}{s} \Bigg\} \,\delta^k 2^m . \end{align*}\]
By setting \(t \rightarrow -t\) we get negated version of Proposition [prop:multifold-binomial-sums-of-powers]
Proposition 36 (Multifold binomial sums of powers negated). For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \frac{1}{2} \sum_{k=0}^{m} (-1)^{m+k} \Bigg \{ \binom{n+t+\frac{k}{2}+r}{k+r} + \binom{n+t+\frac{k}{2}+r-1}{k+r} \\ &- \sum_{s=0}^{r-1} \left[ \binom{r-s+t+\frac{k}{2}}{k+r-s} + \binom{r-1-s+t+\frac{k}{2}}{k+r-s} \right] \binom{s+n-1}{s} \Bigg \} \delta^{k} t^m. \end{align*}\]
Proof. By setting \(t \rightarrow -t\) into Proposition [prop:multifold-binomial-sums-of-powers], and by \(\delta^k (-t)^m = (-1)^{m+k} \delta^{k} t^{m}\). ◻
There are indeed many (infinitely) variations of such formulas as above, at least for every integer \(t\). However, this entropy can be even wilder, by recalling the fundamental properties of binomial coefficients, such as \[\begin{align*} \binom{-x}{k} &= (-1)^{k} \binom{x+k-1}{k}, \\ \binom{x}{k} &= \binom{x}{x-k}. \end{align*}\]
In this manuscript, we derive formulas for sums of powers by combining centered Newton’s interpolation formula for powers [lem:central-newtons-formula-for-powers] with hockey-stick identities [lem:generalized-hockey-stick-identity], and [lem:centered-hockey-stick-identity]. This approach can be generalized further. In particular, while using Newton’s formulas \[\begin{align} \label{eq:newtons-formulas-casess} \begin{cases} f(x) &= \sum_{k=0}^{\infty} \frac{(x-a)^{[k]}}{k!} \delta^{k} f(a) = \sum_{k=0}^{\infty} \frac{x-a}{k} \binom{x-a+\frac{k}{2}-1}{k-1} \delta^{k} f(a), \\ f(x) &= \sum_{k=0}^{\infty} \frac{\left(x-a\right)_{k}}{k!} \Delta^k f(a) = \sum_{k=0}^{\infty} \binom{x-a}{k} \Delta^k f(a), \\ f(x) &= \sum_{k=0}^{\infty} \frac{(x-a)^{\left(k\right)}}{k!} \nabla^k f(a) = \sum_{k=0}^{\infty} \binom{x-a+k-1}{k} \nabla^k f(a), \\ f(x) &= \sum_{k=0}^{\infty} \frac{(x-a)^k}{k!} \frac{\mathrm{d}^{k} f(a)}{\mathrm{d} x^{k}}. \end{cases} \end{align}\] one may replace linear difference operators \(\delta, \Delta, \nabla\) with their non-linear dynamic analogs, allowing more generic approach for sums of powers. It is indeed so remarkable that Taylor’s series \(f(x)= \sum_{k=0}^{\infty} \frac{(x-a)^k}{k!} \frac{\mathrm{d}^{k} f(a)}{\mathrm{d} x^{k}}\) is a continuous analogue of Newton’s formula! Back to the topic, such non-linear difference operators arise naturally in the theory of dynamic equations on time scales, developed by Bohner and Peterson in (Bohner, Martin and Peterson, Allan 2001). Thus, the discrete change \(x+h\) of a function \(f(x+h)\) is replaced by function \(g(h)\) which implies that dynamic difference operator is \[\begin{align} D_g f = \frac{f(g(x)) - f(x)}{g(x) - x}. \label{eq:dynamic-difference-operator} \end{align}\] A well-known example is Jackson’s \(q\)-difference (Jackson, Frederick H. 1909) \[\begin{align*} D_q f(x) = \frac{f(qx) - f(x)}{qx - x}, \end{align*}\] where \(g(x) = qx\). In particular, for Newton’s formula in linear differences \(\delta, \Delta, \nabla\), their respective falling, rising, and central factorials are related to binomial coefficients via \[\begin{align} \label{eq:factorial-identities} \begin{cases} \frac{\left(x\right)_{n}}{n!} &= \frac{1}{n!} x(x-1)(x-2)\cdots(x-n+1) =\binom{x}{n}, \\ \frac{x^{\left(n\right)}}{n!} &= \frac{1}{n!} x(x+1)(x+2)\cdots(x+n-1) =\binom{x+n-1}{n}, \\ \frac{x^{[n]}}{n!} &= \frac{1}{n!} \left( n + \frac{k}{2} -1 \right)\left( n + \frac{k}{2} -2 \right) \cdots \left( n - \frac{k}{2} +1 \right) =\frac{x}{n} \binom{x+\frac{n}{2}-1}{n-1}. \end{cases} \end{align}\] Motivated by this viewpoint, we can see that dynamic version of Newton’s formula is \[\begin{align} \label{eq:generalized-newtons-formula} f(x) = \sum_{k=0}^{\infty} \frac{\left(x\mid \alpha, \beta\right)_{k}}{k!} \cdot D_g f(a), \end{align}\] where \(a\) is an arbitrary integer, \(D_g f(a)\) is dynamic difference operator, and \(\left(x\mid \alpha, \beta\right)_{n}\) is generalized factorial which uses the combination of Nörlund’s notation (Nörlund 1924), and definition given by Steffensen (Steffensen 1933)
Definition 37 (Generalized factorial). For non-negative integers \(x,n,\alpha,\beta\) \[\begin{align*} \left(x\mid \alpha, \beta\right)_{n} = x(x-n\alpha - \beta)(x-n\alpha - 2\beta)\cdots(x-n\alpha - n - I\beta). \end{align*}\]
For instance, \[\begin{align*} \begin{cases} \left(x\right)_{n} &= \left(x\mid \alpha, \beta\right)_{n}, \quad \alpha=0, \quad \beta =1, \quad \text{(forward Newton's formula)} \\ x^{\left(n\right)} &= \left(x\mid \alpha, \beta\right)_{n}, \quad \alpha=0, \quad \beta =-1, \quad \text{(backward Newton's formula)} \\ x^{[n]} &= \left(x\mid \alpha, \beta\right)_{n}, \quad \alpha=\frac{1}{2}, \quad \beta =1, \quad \text{(central Newton's formula)}\\ x^n &= \left(x\mid \alpha, \beta\right)_{n}, \quad \alpha=0, \quad \beta =0, \quad \text{(Taylor's formula)} \end{cases} \end{align*}\] More precisely, the algorithm for closed formulas for sums of powers is
Algorithm 38 (Algorithm for closed formulas).
Define the function \(f(x) = x^m\).
Pick generalized difference operator \(D_g f(x)\), see [eq:dynamic-difference-operator]. In this paper, we use central difference \(\delta^k f(x)\).
Derive generalized Newton’s formula [eq:generalized-newtons-formula] in terms of \(D_g f(x)\), for example one of [eq:newtons-formulas-casess].
Find an identity that connects respective factorials \(\left(x\mid \alpha, \beta\right)_{n}\) with binomial coefficients of the form \(\binom{x}{k}\), so that upper index is expression in \(x\). In this manuscript, we use Lemma [lem:central-factorials-binomial-form], and Lemma [lem:central-factorials-binomial-form-decomposition], however they may vary, based on selected difference operator, which is central difference \(\delta^k f(x)\) for this work.
Apply hockey-stick identities, to derive closed form for sums of powers \(\Sigma^{1}\,{n}^{m}\). The precise hockey-stick identities may vary, depending on factorials \(\left(x\mid \alpha, \beta\right)_{n}\), in this work we use Lemma [lem:generalized-hockey-stick-identity], and Lemma [lem:centered-hockey-stick-identity].
Apply hockey-stick identities repeatedly, to get formula for multifold sums of powers \(\Sigma^{r}\,{n}^{m}\), for example, Theorem [thm:multifold-sums-of-powers-at-zero], and Theorem [thm:multifold-sums-of-powers].
The algorithm above has already been utilized in practice. For instance, the work (Kolosov 2026a) provides remarkable formulas for sums of powers, based on forward difference Newton’s formula,
Proposition 39. For non-negative integers \(r,n,m\) and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \sum_{j=0}^{m} \left[ \binom{n-t+r}{j+r} + \sum_{s=1}^{r} (-1)^{j+s-1} \binom{j+t-1}{j+s} \binom{r-s+n-1}{r-s} \right] \Delta^{j} t^{m}. \end{align*}\]
Also, in terms of Eulerian numbers \(\genfrac{\langle}{\rangle}{0pt}{}{n}{k}\)
Proposition 40. For non-negative integers \(r,n,m\) and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} = \sum_{j=0}^{m} \sum_{k=0}^{m} \left[ \binom{n-t+r}{j+r} - \sum_{s=1}^{r} \binom{s-t}{j+s} \binom{r-s+n-1}{r-s} \right] \binom{t+k}{m-j} \genfrac{\langle}{\rangle}{0pt}{}{m}{k}. \end{align*}\]
The work entitled Sums of powers of integers and generalized Stirling numbers of the second kind (Cereceda 2022) discusses closed formulas for sums of powers in terms of generalized Stirling numbers of the second kind \(\genfrac{\{}{\}}{0pt}{}{k}{j}_{r}\) which originates from forward Newton’s formula for powers, evaluated at arbitrary integer \(t\). It is quite interesting that both (Cereceda 2022), and (Kolosov 2026a) reached the same result independently. In (Cereceda 2022), we read
Proposition 41 (Cereceda (2022)). For non-negative integers \(n,k\), and an arbitrary integer \(r\) \[\begin{align} \label{eq:cereceda-sums-1} \Sigma^{1}\,{n}^{k} = \sum_{j=0}^{k} j! \left[ \binom{n+1-r}{j+1} + (-1)^j \binom{r+j-1}{j+1} \right] \genfrac{\{}{\}}{0pt}{}{k}{j}_{r}, \end{align}\] \[\begin{align} \label{eq:cereceda-sums-2} \Sigma^{1}\,{n}^{k} = \sum_{j=0}^{k} j! \left[ \binom{n+1+r}{j+1} + \binom{r+1}{j+1} \right] \genfrac{\{}{\}}{0pt}{}{k}{j}_{-r}. \end{align}\]
Formula [eq:cereceda-sums-1] is a special case of proposition (Kolosov 2026a, Prop. (1.9))
Proposition 42. For non-negative integers \(n,m\), and an arbitrary integer \(t\) \[\begin{align*} \Sigma^{1}\,{n}^{m} = \sum_{j=0}^{m} \left[ \binom{n-t+1}{j+1} + (-1)^j \binom{j+t-1}{j+1} \right] \Delta^{j} t^m. \end{align*}\]
in terms of generalized Stirling numbers \(\genfrac{\{}{\}}{0pt}{}{k}{j}_{r}\). This implies that finite differences of powers can be expressed in terms of generalized Stirling numbers of the second kind, \(\Delta^{j} t^m = j! \genfrac{\{}{\}}{0pt}{}{m}{j}_{t}\). Formula [eq:cereceda-sums-2] is a special case of proposition (Kolosov 2026a, Prop. (1.8)) with setting \(t=-r\)
Proposition 43. For non-negative integers \(n,m\), and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{1}\,{n}^{m} = \sum_{j=0}^{m} \left[ \binom{n+t+1}{j+1} - \binom{1+t}{j+1} \right] \Delta^{j} (-t)^m. \end{align*}\]
Moving discussion to backward Newton’s formula, the work (Kolosov 2026b) utilizes algorithm [alg:sums-of-powers] in terms of backward differences, providing another remarkable formula
Proposition 44. For non-negative integers \(r,n,m\) and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r+1}\,{n}^{m} = \sum_{j=0}^{m} \left[ \binom{j-t+n+r}{j+r+1} + \sum_{s=0}^{r} (-1)^{j+s} \binom{t}{j+s+1} \binom{r-s+n-1}{r-s} \right] \nabla^{j} t^{m}. \end{align*}\]
Finally, the work (Kolosov 2026c) explores Newton’s formula in central differences, allowing it to be applied with algorithm [alg:sums-of-powers]. Essentially, it achieves similar formulas for sums of powers as current paper, but using different binomial identities to reach the closed form. For instance,
Proposition 45. For non-negative integers \(r,n,m\) and an arbitrary integer \(t\), \[\begin{align*} \Sigma^{r}\,{n}^{m} &= \binom{r+n-1}{r} t^m + \sum_{k=1}^{m} \frac{\delta^{k} t^{m}}{2} \Bigg \{ \left[ \binom{n-t+\frac{k}{2}+r}{k+r} + \binom{n-t+\frac{k}{2}+r-1}{k+r} \right] \\ &- \sum_{s=0}^{r-1} \left[ \binom{-t+\frac{k}{2} + r -s}{k+r-s} + \binom{-t+\frac{k}{2} + r -s-1}{k+r-s} \right] \binom{s+n-1}{s} \Bigg \}. \end{align*}\]
Which matches Theorem [thm:multifold-sums-of-powers], derived however using the following binomial identity, for non-negative integers \(r,n,m\) \[\begin{align*} n \binom{n+r}{m} = (m+1) \binom{n+r}{m+1} - (r-m) \binom{n+r}{m}, \end{align*}\] instead of hockey-stick approach we utilize in this manuscript.
The main results establish an algorithm [alg:sums-of-powers] for deriving closed formulas for multifold sums of powers of integers by combining variations of Newton’s interpolation formula with hockey-stick identities for binomial coefficients.
Using this algorithm, we obtain several notable results, including a closed formula for sums of powers in Proposition [prop:ordinary-sums-of-powers-decomposition], a formula for \(r\)-fold sums of powers in Theorem [thm:multifold-sums-of-powers], and a pure binomial representation of \(r\)-fold sums of powers, in Proposition [prop:multifold-binomial-sums-of-powers].
Furthermore, we show that Knuth’s famous formula for sums of odd powers [prop:knuth-odd-powers-sum] originates from Newton’s interpolation formula [lem:central-newtons-formula-for-powers], thereby obtaining its \(r\)-fold generalization in Proposition [prop:knuth-odd-powers-sum-multifold].
Motivated by Knuth’s formula for odd powers, we also employ the centered Newton interpolation formula evaluated at zero to derive another closed formula for sums of powers in Proposition [prop:ordinary-sums-of-powers-at-zero]. This approach is extended to the \(r\)-fold case in Theorem [thm:multifold-sums-of-powers-at-zero], and its pure binomial counterpart is established in Proposition [prop:multifold-sums-of-powers-at-zero-binomial-form].
Use this GitHub
repository to validate the results using Mathematica programs.
| Mathematica Function | Validates / Prints |
|---|---|
MultifoldSumOfPowersRecurrence[r, n, m] |
Computes \(\Sigma^{r}\,{n}^{m}\), [def:multifold-sums-of-powers-recurrence] |
CentralFactorial[n, k] |
Computes \(n^{[k]}\), [def:central-factorials] |
ValidateCentralFactorialsBinomialDecomposition.txt |
Validates Lem. [lem:central-factorials-binomial-form-decomposition] |
ValidateNewtonsFormulaForPowers.txt |
Validates Lem. [lem:central-newtons-formula-for-powers-decomposition] |
ValidateCenteredHockeyStickIdentity.txt |
Validates Lem. [lem:centered-hockey-stick-identity] |
ValidateOrdinarySumsOfPowersDecomposition.txt |
Validates Prop. [prop:ordinary-sums-of-powers-decomposition] |
ValidateMultifoldSumsKnuthFormula.txt |
Validates Prop. [prop:knuth-odd-powers-sum-multifold] |
ValidateOrdinarySumsOfPowersAtZero.txt |
Validates Prop. [prop:ordinary-sums-of-powers-at-zero] |
ValidateDoubleSumsOfPowersAtZero.txt |
Validates Prop. [prop:double-sums-of-powers-at-zero] |
ValidateMultifoldSumsOfPowersAtZero.txt |
Validates Thm. [thm:multifold-sums-of-powers-at-zero] |
ValidateDoubleSumsOfPowersDecomposition.txt |
Validates Prop. [prop:double-sums-of-powers-decomposition] |
ValidateTripleSumsOfPowersDecomposition.txt |
Validates Prop. [prop:triple-sums-of-powers-decomposition] |
ValidateMultifoldSumsOfPowersTheorem.txt |
Validates Thm. [thm:multifold-sums-of-powers] |
ValidateMultifoldSumsOfZeroPowersLemma.txt |
Validates Lem. [lem:multifold-sum-of-zero-powers] |
ValidateMultifoldBinomialSumsOfPowersProposition.txt |
Validates Prop. [prop:multifold-binomial-sums-of-powers] |
ValidateMultifoldBinomialSumsOfPowersNegated.txt |
Validates Prop. [prop:multifold-binomial-sums-of-powers-negated] |
Metadata
Initial release date: June 04, 2026.
Current release date: 2026-06-23.
Version: 1.1.0+main.c9139c4
.
MSC2010: 05A10, 11B68, 11B73, 11B83 .
Keywords: Sums of powers, Newton’s interpolation formula, Finite differences, Binomial coefficients, Faulhaber’s formula, Bernoulli numbers, Bernoulli polynomials, Interpolation, Central factorial numbers, Stirling numbers, Eulerian numbers, OEIS .
License: This work is licensed under a CC BY 4.0 License.
Web Version: kolosovpetro.github.io/math/sums-of-powers-framework-for-closed-forms
GitHub: github.com/kolosovpetro/SumsOfPowersACompleteFrameworkForClosedForms
ORCID: 0000-0002-6544-8880
Email: kolosovp94@gmail.com