Sums of powers, read off a staircase

The odd numbers stack into squares, two triangles make a rectangle, six stepped pyramids make a box, and the cubes add up to the square of a triangular number. Each of those identities gives a sum of powers — $1 + 2 + \cdots + n$, or $1^2 + 2^2 + \cdots + n^2$, or the cubes — as a closed formula, and each has its own picture. By the fourth powers the pictures run out: a sum of fourth powers would need a five-dimensional stack, and nobody draws those.

The formulas do not run out. For every power $m$ the sum $1^m + 2^m + \cdots + n^m$ is a polynomial in $n$ of degree $m + 1$, and there is a single picture that explains all of them at once. It is not a stack. It is a staircase standing on a curve.

The staircase and the curve #

Draw a bar of width one and height $k^m$ over each interval from $k - 1$ to $k$. The bars’ total area is the sum. Draw the curve $y = x^m$ across the same interval from $0$ to $n$. The curve passes through the top right-hand corner of every bar, and between those corners it bends under the bars’ flat tops.

For the first powers the picture is exact and familiar. The region under $y = x$ from $0$ to $n$ is a triangle of area $n^2/2$. Each bar pokes above the line by a small triangle of area $\tfrac12$, and there are $n$ of them, so the slivers hold $n/2$. The sum is $n^2/2 + n/2$, which is the triangular number $n(n+1)/2$ — the same formula the doubled staircase gave, now read as area under the line plus the slivers above it.

For squares the curve bends and the slivers are no longer triangles.

The area under $y = x^2$ from $0$ to $6$ is $72$, which is the dominant part of 91. The slivers hold 19. Half the last bar is 18, which accounts for almost all of them. And the leftover 1 is a sixth of 6. Put together:

$$1^2 + 2^2 + \cdots + n^2 = \frac{n^3}{3} + \frac{n^2}{2} + \frac{n}{6},$$

which is $n(n+1)(2n+1)/6$, the pyramid formula, multiplied out. The three terms are three different things: the area under the curve, a correction of half a step, and a correction to the correction.

Why the second term is always one half #

The slivers above the curve are nearly triangles, and their size can be estimated without knowing the power.

Over the interval from $k - 1$ to $k$ the bar’s top is at height $k^m$ and the curve rises from $(k-1)^m$ to $k^m$. If the curve were straight across that interval, the sliver would be a triangle of base one and height $k^m - (k-1)^m$, with area half of that. Adding those half-differences from $k = 1$ to $n$ telescopes: everything cancels except $\tfrac12 n^m$. So the slivers hold about half the last bar, for any power, and the second term of every formula is $n^m/2$.

The approximation is that the curve is straight across each unit interval, and it is not quite. A curve that bends upward — as $x^m$ does for $m \ge 2$ — sags below its chord, so each sliver is slightly larger than the triangle, by an amount that depends on how sharply the curve bends. That excess is the third term, and it is where the power starts to matter.

For cubes the numbers are 441 in the bars, 324 under the curve, 117 in the slivers and 108 for the half-step, leaving 9. The polynomial is

$$1^3 + 2^3 + \cdots + n^3 = \frac{n^4}{4} + \frac{n^3}{2} + \frac{n^2}{4},$$

and $6^2/4 = 9$. The picture of the cubes as the square of a triangular number is still true — $(n(n+1)/2)^2$ multiplies out to exactly this — but now it has a second reading: an integral, a half-step, and a curvature correction.

Every coefficient, exactly #

The pattern continues for every power, and the coefficients can be computed exactly as fractions.

Three diagonals of the table are forced. The leading coefficient is $1/(m+1)$, the area under $x^m$. The next is $1/2$, the half-step. The third is $m/12$ — for squares $2/12 = 1/6$, for cubes $3/12 = 1/4$, for fourth powers $4/12 = 1/3$. After that the coefficients alternate between nought and a fraction, and the fractions look unrelated: $-1/30$, $-1/12$, $1/42$, $-1/6$.

They are related. Every coefficient in the table is a binomial coefficient times one of a single sequence of numbers:

$$1^m + \cdots + n^m = \frac{1}{m+1}\sum_{j=0}^{m} \binom{m+1}{j} B_j\, n^{m+1-j},$$

where $B_0 = 1$, $B_1 = \tfrac12$, $B_2 = \tfrac16$, $B_3 = 0$, $B_4 = -\tfrac1{30}$, $B_5 = 0$, $B_6 = \tfrac1{42}$, and so on. Those are the Bernoulli numbers. Johann Faulhaber had worked out the polynomials up to the seventeenth power by 1631 without recognising the pattern in their coefficients; Jacob Bernoulli found it and published it in 1713, remarking that he had computed the sum of the tenth powers of the first thousand numbers in under half a quarter of an hour.

The number he got is 91,409,924,241,424,243,424,241,924,242,500, and its strange near-repetitions are the polynomial showing through. With $n = 1000$ each power of $n$ is a one followed by threes of zeros, so the terms of the degree-eleven polynomial land in overlapping blocks of digits, and the answer is the decimal expansions of the coefficients $1/11$, $1/2$, $5/6$, $-1$, $1$, $-1/2$ and $5/66$ laid over one another, each shifted three digits for every power of $n$ it drops. A sum that would take thousands of multiplications becomes the coefficients of one row of the table written out side by side. Seki Takakazu found the same numbers independently in Japan, in a book published in 1712, a year before Bernoulli’s.

The third term is curvature #

The picture says what each term is, and it says it most clearly for the third.

The sliver over the interval from $k-1$ to $k$ exceeds the chord’s triangle by the area between the chord and the curve. For a smooth curve that area is about $\tfrac1{12}$ times the curve’s second derivative there — the familiar error of the trapezium rule. The second derivative of $x^m$ is $m(m-1)x^{m-2}$, and adding $\tfrac1{12}m(m-1)k^{m-2}$ from $k = 1$ to $n$ gives, to leading order, $\tfrac1{12} m\, n^{m-1}$. The third coefficient is the trapezium rule’s error summed along the staircase, and $m/12$ is what it comes to.

The general statement is the Euler–Maclaurin formula. For any smooth function $f$, the sum $f(1) + \cdots + f(n)$ equals the integral of $f$, plus half the difference of its end values, plus $B_2/2!$ times the difference of its first derivative at the ends, plus $B_4/4!$ times the difference of its third derivative, and so on — one term for each odd derivative, weighted by a Bernoulli number, plus a remainder. For $f(x) = x^m$ the derivatives run out after $m$ steps, the remainder vanishes, and the formula is Faulhaber’s polynomial exactly. For any other function the derivatives never run out and the formula is an infinite series of corrections.

The odd-indexed Bernoulli numbers past $B_1$ are all nought for a reason the picture also supplies. The function $x/(e^x - 1) + x/2$, whose Taylor coefficients generate the Bernoulli numbers, is even — and an even function has no odd powers. In the staircase that symmetry is the observation that a correction for curvature cannot depend on which way the staircase is read, so only even derivatives of the error can appear.

The corrections that eventually grow #

The corrections for a polynomial stop. For a function like $\log x$ or $1/x$ they do not, and whether an infinite list of corrections converges depends on how big the Bernoulli numbers get.

They get very big. After a few small values the sizes grow like $2\,(2k)!/(2\pi)^{2k}$, which outruns every exponential because of the factorial. So in the Euler–Maclaurin series for a function like $1/x$, the terms first shrink — each correction smaller than the last — reach a smallest term, and then grow without limit. The series diverges for every function that is not a polynomial.

That does not make it useless; it makes it an asymptotic series. Stopping at the smallest term gives an approximation whose error is about the size of that term, and for large $n$ that error is fantastically small. Euler used exactly this to compute the constant $\gamma = 0.5772\ldots$ — the amount by which the harmonic series exceeds $\log n$ — to sixteen decimal places, adding ten terms by hand and correcting with a handful of Bernoulli numbers, at a time when adding up enough terms directly would have taken longer than the age of the universe.

The twelve in Stirling’s formula #

The same corrections, applied to a function whose sum everybody needs, produce a constant that is usually presented as a mystery.

The logarithm of $n!$ is $\log 1 + \log 2 + \cdots + \log n$, a staircase over the curve $y = \log x$. The integral of $\log x$ from 1 to $n$ is $n\log n - n + 1$. The half-step is $\tfrac12\log n$. The curvature correction involves the first derivative of $\log x$, which is $1/x$, weighted by $B_2/2! = 1/12$. Putting the pieces together, and collecting everything that does not depend on $n$ into one constant,

$$\log n! = n\log n - n + \tfrac12\log n + C + \frac{1}{12n} - \frac{1}{360n^3} + \cdots.$$

The constant turns out to be $\tfrac12\log 2\pi$, which is where the $\sqrt{2\pi n}$ in Stirling’s formula comes from, and the $\frac{1}{12n}$ is the first correction that makes the formula accurate enough to use for small $n$. That twelve is the same twelve as in the $m/12$ of the table — the trapezium rule’s curvature error — and the $\frac{1}{360}$ after it is $B_4/(4 \cdot 3)$. Stirling’s series is the Euler–Maclaurin series for the logarithm, and like every such series for a non-polynomial it diverges: its terms shrink for a while at any fixed $n$ and then grow, with the turning point further out the larger $n$ is.

That series is what sits behind the $\sqrt{\pi d/2}$ in the average room at which a shared birthday first appears, where factorials in the all-distinct product had to be estimated. A problem about people and days reaches its $\pi$ through a staircase over a logarithm.

Who needed these sums, and for what #

The formulas were not wanted for their own sake. Every one of them was first needed to find an area or a volume, and the history runs exactly parallel to the picture.

Archimedes needed the sum of the squares to find the area inside a spiral and the area of a parabolic segment, and he proved it by the same kind of bounding from above and below that he used to pin the area of a circle between two sequences. Around the year 1000 Ibn al-Haytham needed the sum of the fourth powers to find the volume of a paraboloid, and derived it — an argument that, written in modern notation, is the step from one power to the next that Faulhaber’s table encodes. In the seventeenth century Cavalieri, Fermat and Pascal each used sums of powers to find the area under $y = x^m$, and it was the recognition that the answer is $n^{m+1}/(m+1)$ for every $m$ that became one of the first general theorems of integration, and the slicing method Cavalieri built on it.

The picture reverses the direction of that history. Those mathematicians found the sums of powers in order to compute areas under curves. Once the areas were known the relationship could be run the other way, and the area became the leading term of the sum, with the Bernoulli numbers keeping the accounts for the difference.

Where the dashed curve comes from #

The dashed curve in the figure is not a fit, and the reason it matches the Bernoulli numbers to within two per cent from $B_6$ on is one of the more surprising connections in elementary mathematics.

The exact relation is

$$B_{2k} = (-1)^{k+1}\,\frac{2\,(2k)!}{(2\pi)^{2k}}\;\zeta(2k), \qquad \zeta(s) = 1 + \frac1{2^s} + \frac1{3^s} + \cdots.$$

The sum $\zeta(2k)$ is between 1 and $\pi^2/6$, and it tends to 1 quickly as $k$ grows, because the terms after the first become negligible. So the Bernoulli numbers are the factorial curve times a factor that barely differs from one. Read the other way, the relation gives the values of the zeta function at even integers: $\zeta(2) = \pi^2/6$, which Fourier’s coefficients deliver by a completely different route, $\zeta(4) = \pi^4/90$, and every $\zeta(2k)$ as a rational multiple of $\pi^{2k}$ with a Bernoulli number in it.

So the coefficients that correct a staircase for its corners are the same numbers that sum the reciprocals of the even powers. The $2\pi$ in the formula is there because the generating function $x/(e^x - 1)$ has poles at $\pm 2\pi i$, which set the radius of convergence of its Taylor series and therefore the growth rate of its coefficients — a fact about complex numbers controlling the size of rational numbers that were first met adding up whole-number powers.

What nothing drawn here establishes #

The staircases are drawn for powers one, two and three at six to eight steps, and the table checks each polynomial against direct sums for $n$ up to 25. That is evidence that the formulas are right and a proof of nothing: the polynomial identity for a given power is proved by induction on $n$, or by the telescoping of $(k+1)^{m+1} - k^{m+1}$, and no finite table does that work.

The explanation of the second and third terms as a half-step and a curvature correction is an argument about areas, and it is exact only in the limit. At six steps the slivers are visibly not triangles, and the decomposition of 19 into 18 and 1 is a statement about the polynomial rather than something the drawn shapes can be measured to show.

And the divergence of the Euler–Maclaurin series cannot be seen in any picture of a staircase, because it concerns the corrections for a non-polynomial function continued to very high order. The Bernoulli numbers’ growth is the evidence, and it is shown only up to $B_{30}$.

A curve that explains the triangular numbers again #

The staircase reading also connects back to the figurate numbers where this started. The triangular number $n(n+1)/2$ is the first staircase, and the eight-triangles picture that turns every odd square into eight triangles and a dot is a statement about how staircases of that kind fit inside squares. The pyramid numbers are the second staircase, and the six pyramids in a box are the polynomial $n^3/3 + n^2/2 + n/6$ multiplied by six and read as a solid. Each dimension has its own trick, and the tricks stop at three.

What the staircase over a curve replaces them with is a single mechanism: an integral for the bulk, a half-step for the edge, and an infinite family of corrections for the corners, each one a finer measurement of how the flat top of a bar differs from the curve it sits on. It works for every power and extends to every smooth function, and the same comparison of bars with a curve that defines the integral turns out, carried to its corrections, to be the thing that computes the sums the integral was defined to approximate.

Still open: whether the odd zeta values are like the even ones #

The even zeta values are rational multiples of powers of $\pi$, with a Bernoulli number fixing the rational, and all of them are transcendental because $\pi$ is. Nothing like that is known for the odd values. $\zeta(3) = 1.2020569\ldots$ was proved irrational by Roger Apéry in 1978, in an argument so unexpected that it was checked line by line before anyone believed it. Since then it has been shown that infinitely many of $\zeta(5), \zeta(7), \zeta(9), \dots$ are irrational, and that at least one of $\zeta(5), \zeta(7), \zeta(9), \zeta(11)$ is.

Nobody knows whether $\zeta(5)$ is irrational, and nobody knows whether $\zeta(3)$ is a rational multiple of $\pi^3$ — the analogue of the even case — though it is universally expected not to be. The Bernoulli numbers answer the even half of the question completely, and there is no known sequence of rational numbers that plays their part for the odd half.

The polynomial in every staircase #

A sum of whole-number powers looks like the most discrete thing imaginable, a list of integers added up. Its formula turns out to be a continuous object — the area under a curve — with a sequence of corrections, and the corrections are organised by one list of rational numbers that shows up again in the zeta function, in the trapezium rule, and in the Taylor series of $x/(e^x - 1)$.

The staircase is the reason. A bar has a flat top and a curve does not, and every term after the first is a more precise account of that one difference: half a step, then the bend, then the bend of the bend — a list that stops for polynomials and, for everything else, first shrinks and then grows past every bound.