Derive a formula $\sum_{j=1}^{n-1}j^a \sim \zeta(-a)+\frac{n^{a+1}}{a+1}\sum_{s=0}^\infty \binom{a+1}{s}\frac{B_s}{n^s} $

Question

Detailed derivations on this formula:

There is an exercise problem (Ex. 3.2) in Olver's book, on page 292. (or you can find it on this website: Eq.2.10.7).

$$\sum_{j=1}^{n-1}j^a \sim \zeta(-a)+\frac{n^{a+1}}{a+1}\sum_{s=0}^\infty \binom{a+1}{s}\frac{B_s}{n^s} \tag{*}$$

I try to derive it. I begin with:

$$\sum_{j=n_0}^n f(j)=\int_{n_0}^n f(x) dx+\frac{f(n_0)+f(n)}{2}+\sum_{s=1}^{m-1}\frac{B_{2s}}{(2s)!}\left( f^{(2s-1)}(n)-f^{(2s-1)}(n_0) \right)+R_m(n) $$

let $n_0=1$ and $f(x)=x^a$, for left-hand-side:

$$\sum_{j=n_0}^n f(j)=\sum_{j=1}^n j^a$$

for right-hand-side:

$$\begin{align} \int_{n_0}^n f(x) dx&=\int_1^n x^a dx=\frac{1}{a+1}x^{a+1}|_1^n=\frac{n^{a+1}}{a+1}-\frac{1}{a+1}\tag{1}\\ \\ \frac{f(n_0)+f(n)}{2}&=\frac{f(1)+f(n)}{2}=\frac{1+n^a}{2}\tag{2}\\ \\ f^{(2s-1)}(x)&=a(a-1)...(a-2s+2)x^{a-2s+1}=\frac{a!}{(a-2s+1)!}x^{a-2s+1}\\ \\ f^{(2s-1)}(n)&=\frac{a!}{(a-2s+1)!}n^{a-2s+1}\\ \\ f^{(2s-1)}(n_0)&=f^{(2s-1)}(1)=\frac{a!}{(a-2s+1)!} \end{align}$$

$$\begin{align} \sum_{s=1}^{m-1}\frac{B_{2s}}{(2s)!}\left( f^{(2s-1)}(n)-f^{(2s-1)}(n_0) \right)&=\sum_{s=1}^{m-1}\frac{B_{2s}}{(2s)!}\frac{a!}{(a-2s+1)!}\left(n^{a-2s+1}-1\right)\\ \\ &=\frac{1}{a+1}\sum_{s=1}^{m-1}\frac{B_{2s}}{(2s)!}\frac{(a+1)!}{(a-2s+1)!}\left(n^{a-2s+1}-1\right)\\ \\ &=\frac{1}{a+1} \sum_{s=1}^{m-1}B_{2s}\binom{a+1}{2s }\left(n^{a-2s+1}-1\right) \end{align}$$

Next, substitute: $s'=2s$ and use the fact $B_{2k+1}=0$ for $k=1,2,3,...$

$$\begin{align} \sum_{s=1}^{m-1}\frac{B_{2s}}{(2s)!}\left( f^{(2s-1)}(n)-f^{(2s-1)}(n_0) \right)&=\frac{1}{a+1} \sum_{s'=2}^{2m-2}B_{s'}\binom{a+1}{s'}\left(n^{a-s'+1}-1\right)\\ \\ &=\frac{n^{a+1}}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}\frac{B_{s'}}{n^{s'}}-\frac{1}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}B_{s'}\tag{3} \end{align}$$

Combine $(1)(2)(3)$,

$$\begin{align} \sum_{j=1}^n j^a&=\frac{n^{a+1}}{a+1}-\frac{1}{a+1}+ \frac{1+n^a}{2}+\frac{n^{a+1}}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}\frac{B_{s'}}{n^{s'}}-\frac{1}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}B_{s'}+R_m(n)\\ \\ \sum_{j=1}^n j^a&=n^a+\frac{n^{a+1}}{a+1}-\frac{n^a}{2}+\frac{n^{a+1}}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}\frac{B_{s'}}{n^{s'}}\\ &~~~~~~~~~-\frac{1}{a+1}+\frac{1}{2}-\frac{1}{a+1} \sum_{s'=2}^{2m-2}\binom{a+1}{s'}B_{s'}+R_m(n)\\ \end{align}$$

Use the fact $B_0=1,~B_1=-\frac{1}{2}$:

$$\begin{align} \sum_{j=1}^n j^a=n^a+\frac{n^{a+1}}{a+1}\sum_{s'=0}^{2m-2}\binom{a+1}{s'}\frac{B_{s'}}{n^{s'}}-\frac{1}{a+1}\sum_{s'=0}^{2m-2}\binom{a+1}{s'}B_{s'}+R_m(n)\tag{4} \end{align}$$

The remainder becomes:

$$\begin{align} R_m(n)&=\int_{n_0}^n \frac{B_{2m}-B_{2m}(x-\lfloor x\rfloor)}{(2m)!}f^{(2m)}(x)dx\\ \\ &=\color{red}{\int_{n_0}^\infty \frac{B_{2m}-B_{2m}(x-\lfloor x\rfloor)}{(2m)!}f^{(2m)}(x)dx}\color{blue}{-\int_{n}^\infty \frac{B_{2m}-B_{2m}(x-\lfloor x\rfloor)}{(2m)!}f^{(2m)}(x)dx}\tag{5}\end{align}$$

where $$f^{(2m)}(x)=a(a-1)...(a-2m+1)x^{a-2m}=\frac{a!}{(a-2m)!}x^{a-2m}$$

Since $|B_{2m}-B_{2m}(x-\lfloor x\rfloor)|\le 2|B_{2m}|$, the $\color{blue}{\text{blue-colored term}}$ is bounded by

$$\begin{align}\left|-\int_{n}^\infty \frac{B_{2m}-B_{2m}(x-\lfloor x\rfloor)}{(2m)!}f^{(2m)}(x)dx\right|&\le 2|B_{2m}|\int_{n}^\infty \frac{1}{(2m)!}|f^{(2m)}(x)|dx\\ \\ &=\left|\frac{2B_{2m}}{(2m)!}\cdot\frac{a!}{(a-2m)!}\cdot \frac{-1}{a-2m+1}\cdot n^{a-2m+1}\right|\\ \\ &=\left|\frac{2B_{2m}}{a+1}\cdot\binom{a+1}{2m}\right|\cdot \frac{n^{a+1}}{n^{2m}}\\ \\ &=\mathcal{O}\left( \frac{n^{a+1}}{n^{2m}}\right)\tag{6}\end{align}$$

The $\color{red}{\text{red-colored term}}$ becomes

$$ \begin{align} \int_{n_0}^\infty \frac{B_{2m}-B_{2m}(x-\lfloor x\rfloor)}{(2m)!}f^{(2m)}(x)dx&=\frac{B_{2m}}{(2m)!}f^{(2m-1)}(x)\bigg|_{n_0}^\infty-\frac1{(2m)!}\sum_{k=n_0}^\infty\int_{k}^{k+1} B_{2m}(x-k)f^{(2m)}(x)dx\\ \\ &=-\frac{B_{2m}}{(2m)!}f^{(2m-1)}(n_0)-\frac{1}{(2m)!}\sum_{k=n_0}^\infty I_{k,2m}\tag{7} \end{align}$$

where

$$I_{k,p}=\int_{k}^{k+1} B_{p}(x-k)f^{(p)}(x)dx$$

perform the integration by part, and we get

$$I_{k,p}=B_p(1)f^{(p-1)}(k+1)-B_p(0)f^{(p-1)}(k)-\int_{k}^{k+1} B'_{p}(x-k)f^{(p-1)}(x)dx$$

Use the property: $B'_s(t)=s\cdot B_{s-1}(t),~~s=1,2,3,\dots$, we get recursive equation

$$I_{k,p}=B_p(1)f^{(p-1)}(k+1)-B_p(0)f^{(p-1)}(k)-p\cdot I_{k,p-1}\tag{8}$$

Note the following properties between Bernoulli polynomial and Bernoulli number:

$$B_p(1)=B_{p}(0)=B_p,~~p=2,3,\dots~~~~~\text{But}~~~~B_1(1)=-B_1=\frac12,~~ B_1(0)=B_1=-\frac12$$

hence, we need to deal with the $p=1$ cases separately. Then, eq.(8) becomes

$$I_{k,p}=B_p\cdot\left[f^{(p-1)}(k+1)-f^{(p-1)}(k)\right]-p\cdot I_{k,p-1},~~~~p=2,3,\dots\tag{9}$$

Evaluate eq.(9) recursively and we get

$$I_{k,p}=\sum_{i=2}^p (-1)^{p-i}\cdot \frac{p!}{i!}\cdot B_i\cdot\left[f^{(i-1)}(k+1)-f^{(i-1)}(k)\right]+(-1)^{p-1}\cdot p!\cdot I_{k,1},~~p\ge2\tag{10}$$

Not done yet, I will update it later.

Not done yet, I will update it later.

Not done yet, I will update it later.

Answer 1

Exact formula

We can show a simpler version with $a$ a non-negative integer, which gives a formula that is exact:

$$\sum_{j=1}^{n-1} j^a = \zeta(-a) + \frac{n^{a+1}}{a+1} \sum_{m\ge 0} {a+1\choose m} \frac{B_m}{n^m}.$$

We start with the sum term

$$\frac{n^{a+1}}{a+1} \sum_{m\ge 0} {a+1\choose m} \frac{B_m}{n^m} = \frac{n^{a+1}}{a+1} (a+1)! [z^{a+1}] \frac{z/n}{\exp(z/n)-1} \exp(z) \\ = a! [z^{a+1}] \frac{z \exp(nz)}{\exp(z)-1} \\ = a! [z^{a+1}] \frac{z}{\exp(z)-1} \sum_{q=0}^n {n\choose q} (\exp(z)-1)^q \\ = \frac{B_a}{a+1} + a! [z^a] \sum_{q=1}^n {n\choose q} (\exp(z)-1)^{q-1}.$$

The first term cancels with the Zeta function value, leaving

$$a! [z^a] \sum_{q=1}^n {n\choose q} \sum_{j=0}^{q-1} {q-1\choose j} \exp(jz) (-1)^{q-1-j} \\ = a! [z^a] \sum_{q=0}^{n-1} {n\choose q+1} \sum_{j=0}^q {q\choose j} \exp(jz) (-1)^{q-j} \\ = \sum_{j=1}^{n-1} j^a (-1)^j \sum_{q=j}^{n-1} (-1)^q {n\choose q+1} {q\choose j} \\ = \sum_{j=1}^{n-1} j^a (-1)^j [z^{n-1}] (1+z)^n \sum_{q\ge j} {q\choose j} (-1)^q z^q \\ = \sum_{j=1}^{n-1} j^a [z^{n-1}] (1+z)^n z^j \sum_{q\ge 0} {q+j\choose j} (-1)^q z^q \\ = \sum_{j=1}^{n-1} j^a [z^{n-1}] (1+z)^n z^j \frac{1}{(1+z)^{j+1}} = \sum_{j=1}^{n-1} j^a [z^{n-1-j}] (1+z)^{n-1-j} \\ = \sum_{j=1}^{n-1} j^a.$$

This is the claim.

Asymptotic formula

There is a standard technique that produces the complete asymptotic expansion for this sum and many others like it, which is to use harmonic sums and Mellin transforms. We will deploy this here to prove the formula for $a \lt -1$, with $a$ a real number.

Put $a=-\alpha$ so that $\alpha \gt 1$ and introduce the telescoping sum

$$S(x) = \sum_{k\ge 1} \left(\frac{1}{k^\alpha}-\frac{1}{(x+k)^\alpha}\right).$$ This sum has the property that $$S(n-1) = S(n)-\frac{1}{n^\alpha} = \sum_{q=1}^{n-1} \frac{1}{q^\alpha},$$ so that $S(n)$ is the value we are looking for.

Re-write the sum as follows: $$S(x) = \sum_{k\ge 1} \frac{1}{k^\alpha} \left(1-\frac{1}{(x/k+1)^\alpha}\right).$$

The sum term is harmonic and may be evaluated by inverting its Mellin transform.

Recall the harmonic sum identity $$\mathfrak{M}\left(\sum_{k\ge 1} \lambda_k g(\mu_k x);s\right) = \left(\sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} \right) g^*(s)$$ where $g^*(s)$ is the Mellin transform of $g(x).$

In the present case we have $$\lambda_k = \frac{1}{k^\alpha}, \quad \mu_k = \frac{1}{k} \quad \text{and} \quad g(x) = 1 - \frac{1}{(1+x)^\alpha}.$$

It follows that $$\sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} = \sum_{k\ge 1} \frac{1}{k^\alpha} \times k^s = \zeta(\alpha-s)$$ which has half-plane of convergence $\alpha-s > 1$ or $s < \alpha-1.$

We need the Mellin transform $g^*(s)$ of $g(x)$ which is $$\int_0^\infty \left(1 - \frac{1}{(1+x)^\alpha}\right) x^{s-1} dx$$

which is immediately seen to be a beta function integral with value $$g^*(s) = - \frac{1}{\Gamma(\alpha)} \Gamma(s)\Gamma(\alpha-s)$$ and fundamental strip $\langle -1, 0 \rangle,$ which is covered by the half-plane of convergence of the Zeta function term.

It follows that the Mellin transform $Q(s)$ of the harmonic sum $S(x)$ is given by $$Q(s) = - \frac{1}{\Gamma(\alpha)} \Gamma(s)\Gamma(\alpha-s) \zeta(\alpha-s).$$ The Mellin inversion integral here is $$\frac{1}{2\pi i} \int_{-1/2-i\infty}^{-1/2+i\infty} Q(s)/x^s ds$$ which we evaluate by shifting it to the right for an expansion at infinity.

First treat the pole from the zeta function term at $s=\alpha-1$, which has

$$\mathrm{Res}(Q(s)/x^s; s=\alpha-1) = -\frac{1}{\Gamma(\alpha)} \Gamma(\alpha-1) \Gamma(1)\times -1 \times x^{1-\alpha} \\ = -\frac{1}{1-\alpha} x^{1-\alpha}.$$

For the pole at $s=-1$ we get no contribution as it is to the left of the inversion integral.

For the pole at $s=0$ from the simple gamma function term we obtain $$\mathrm{Res}(Q(s)/x^s; s=0) = -\frac{1}{\Gamma(\alpha)} \Gamma(\alpha) \zeta(\alpha) = -\zeta(\alpha).$$

For the pole at $s=\alpha$ from the compound gamma function term we obtain $$\mathrm{Res}(Q(s)/x^s; s=\alpha) = -\frac{1}{\Gamma(\alpha)} \Gamma(\alpha) \times -1 \times\zeta(0) x^{-\alpha} = -\frac{1}{2} x^{-\alpha}.$$

The remaining poles are at $s = q+\alpha$ where $q\ge 1$ and contribute $$\mathrm{Res}(Q(s)/x^s; s=q+\alpha) = -\frac{1}{\Gamma(\alpha)} \Gamma(q+\alpha) \frac{(-1)^{q+1}}{q!} \zeta(-q) \frac{1}{x^{q+\alpha}} \\ = - \prod_{p=0}^{q-1} (p+\alpha) \times \frac{(-1)^{q+1}}{q!} (-1)^q \frac{B_{q+1}}{q+1} \frac{1}{x^{q+\alpha}} = B_{q+1} \frac{\prod_{p=0}^{q-1} (p+\alpha)}{(q+1)!} \frac{1}{x^{q+\alpha}} \\ = (-1)^q B_{q+1} \frac{\prod_{p=0}^{q-1} (-\alpha-p)}{(q+1)!} \frac{1}{x^{q+\alpha}} = \frac{(-1)^q}{1-\alpha} {-\alpha+1\choose q+1} \frac{B_{q+1}}{x^{q+\alpha}}.$$

The zero values of the Bernoulli numbers correctly represent cancelation of the gamma function poles by the trivial zeros of the zeta function.

Setting $x=n$ and observing that the shift to the right produces a minus sign we obtain the following asymptotic expansion: $$S(n) = \sum_{k=1}^n \frac{1}{k^\alpha} \sim \frac{1}{1-\alpha} n^{1-\alpha} + \frac{1}{2} n^{-\alpha} + \zeta(\alpha) + \frac{1}{1-\alpha} \sum_{q\ge 2} (-1)^{q} {1-\alpha\choose q} \frac{B_q}{n^{q+\alpha-1}}.$$

Note however that the cited formula is $S(n)-\frac{1}{n^\alpha}$ so we get

$$S(n)-\frac{1}{n^\alpha} = \sum_{k=1}^{n-1} \frac{1}{k^\alpha} \sim \frac{1}{1-\alpha} n^{1-\alpha} - \frac{1}{2} n^{-\alpha} + \zeta(\alpha) + \frac{n^{1-\alpha}}{1-\alpha} \sum_{q\ge 2} {1-\alpha\choose q} \frac{B_q}{n^q}.$$

Here we have made use of the fact that with $q\ge 2$ the odd-index Bernoulli numbers are zero so we may drop the $(-1)^q$ term.

Now observe that for $q=1$ the sum term will produce

$$\frac{n^{1-\alpha}}{1-\alpha} \times (1-\alpha) \times -\frac{1}{2} \frac{1}{n} = - \frac{1}{2} n^{-\alpha}$$

Furthermore for $q=0$ we find

$$\frac{n^{1-\alpha}}{1-\alpha} \times 1 \times 1 = \frac{n^{1-\alpha}}{1-\alpha}.$$

This means we may merge the two leading terms into the sum and we find with $a\lt -1$, $\alpha = -a$,

$$\bbox[5px,border:2px solid #00A000]{ \sum_{j=1}^{n-1} j^a \sim \zeta(-a) + \frac{n^{a+1}}{a+1} \sum_{q\ge 0} {a+1\choose q} \frac{B_q}{n^q}.}$$

as claimed.

Answer 2

The problem is that you did not deal properly with the terms coming from the lower bound of the summation and you did not investigate the remainder. Let $\alpha \neq -1$ fixed and choose a positive integer $m$ so that $2m - 1 > \alpha$. Using the Euler–Maclaurin formula, we have \begin{align*} & \!\!\!\!\!\!\!\sum\limits_{j = 1}^{n - 1} {j^\alpha } = - n^\alpha \!+ \!\sum\limits_{j = 1}^n {j^\alpha } = - n^\alpha \! +\! \int_1^n\! {t^\alpha dt} + \frac{{n^\alpha + 1}}{2} + \!\sum\limits_{s = 1}^{m - 1} {\frac{{B_{2s} }}{{(2s)!}}\frac{{\alpha !}}{{(\alpha - 2s + 1)!}}\left[ {n^{\alpha - 2s + 1} \! -\! 1} \right]} \\ & + \int_1^n {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt} \\ = \; & \frac{{n^{\alpha + 1} }}{{\alpha + 1}} - \frac{{n^\alpha }}{2} + n^{\alpha + 1} \sum\limits_{s = 1}^{m - 1} {\frac{{B_{2s} }}{{(2s)!}}\frac{{\alpha !}}{{(\alpha - 2s + 1)!}}n^{ - 2s} } \\ & - \int_n^{ + \infty } {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt} \\ & - \frac{1}{{\alpha + 1}} + \frac{1}{2} - \!\sum\limits_{s = 1}^{m - 1} {\frac{{B_{2s} }}{{(2s)!}}\frac{{\alpha !}}{{(\alpha - 2s + 1)!}}} +\! \int_1^{ + \infty } {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt} \\ = \; & \frac{{n^{\alpha + 1} }}{{\alpha + 1}}\sum\limits_{s = 0}^{2m - 2} {\frac{{B_s }}{{s!}}\frac{{(\alpha + 1)!}}{{(\alpha - s + 1)!}}n^{ - s} } -\! \int_n^{ + \infty } {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt} \\ & - \frac{1}{{\alpha + 1}} + \frac{1}{2} - \!\sum\limits_{s = 1}^{m - 1} {\frac{{B_{2s} }}{{(2s)!}}\frac{{\alpha !}}{{(\alpha - 2s + 1)!}}} + \!\int_1^{ + \infty } {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt}. \end{align*} It is easy to show that $$ \int_n^{ + \infty } {\frac{{B_{2m} - B_{2m} (\left\{ t \right\})}}{{(2m)!}}\frac{{\alpha !}}{{(\alpha - 2m)!}}t^{\alpha - 2m} dt} = n^{\alpha + 1} \mathcal{O}(n^{ - 2m} ) $$ as $n\to +\infty$. The last line we got from the Euler–Maclaurin formula does not depend on $n$ but seemingly depends on $m$. You have to show that it is independent of $m$ and its value is $\zeta(-\alpha)$. If $\alpha<-1$, this follows from the Euler–Maclaurin formula applied to the infinite version of the sum. To extend it to $\alpha<2m-1$, you can appeal to analytic continuation.