Bernoulli corrections for polynomial sums (source code)

= Bernoulli corrections for polynomial sums
{c}
{title2=$\sum_{k=1}^{N}k^4=\frac{N^5}{5}+\frac{N^4}{2}+\frac{N^3}{3}-\frac{N}{30}$}

For a concrete comparison, take $f(x)=x^2$ on $[0,6]$. The <integral> is $I=\int_0^6x^2\,dx=72$, whereas the unit-width left and right <Riemann sums> are $L=55$ and $R=91$. The exact <integer power sum> formula gives
$$
L=I-\frac{6^2}{2}+\frac66=55,\qquad
R=I+\frac{6^2}{2}+\frac66=91.
$$
The endpoint choice explains the opposite signs of the $18$ correction. The remaining correction is one in both cases.

\Image[/bernoulli-square-sums-versus-integral.png]
{title=Left and right <Riemann sums> for $x^2$ differ from its <integral> by exact <Bernoulli number> corrections}
{description=Six rectangles of width one sample $x^2$ at $0,\ldots,5$ on the left and $1,\ldots,6$ on the right. Their areas are $55$ and $91$; the smooth area is $72$. The endpoint terms are $-18$ and $+18$, while the $B_2$ correction adds one in each case.}

The same method is useful far beyond <squares of numbers>. Substituting $B_0$, $B_1$, $B_2$ and $B_4$ into the <integer power sum> formula gives the exact fourth-power result
$$
\sum_{k=1}^{N}k^4=\frac{N^5}{5}+\frac{N^4}{2}+\frac{N^3}{3}-\frac{N}{30}.
$$
At $N=10$, this is
$$
\underbrace{20000}_{\text{integral}}+\underbrace{5000}_{\text{endpoint}}+
\underbrace{\frac{1000}{3}}_{B_2\text{ term}}-
\underbrace{\frac13}_{B_4\text{ term}}=25333.
$$
Thus a finite <polynomial> calculation replaces ten additions, or a million additions if $N$ is a million. Although the individual corrections can be <rational numbers>, the result is the exact <integer> sum.

\Image[/bernoulli-fourth-power-exact-corrections.png]
{title=Exact <Bernoulli number> corrections turn an <integral> into the fourth-power sum $1^4+\cdots+10^4$}
{description=The first panel builds the sum from integral $20000$, endpoint term $5000$, and $B_2$ term $1000/3$. The enlarged second panel resolves the final $B_4$ correction $-1/3$, reaching $25333$, checked independently by direct addition.}

There is also a useful <numerical integration> interpretation. The average of the left and right <Riemann sums> is the composite <trapezoidal rule>. For $x^2$, the endpoint corrections cancel but its value is $73$, one above the true <integral>. Subtracting the $B_2$ correction makes this particular <polynomial> example exact. For $f(x)=x^4$, the two interior corrections above can be written $B_2(f'(N)-f'(0))/2!$ and $B_4(f'''(N)-f'''(0))/4!$. These are exact identities for these <polynomials>; extending such corrections to general <functions> requires a separate error analysis.