Different ways of evaluating $\int_{0}^{\infty}\frac{x^n dx}{(e^x-1)^2}$

Question

A while ago, I asked a question: Finding a closed form for this expression

The expression in that question,

$$\frac{1}{4}\sum_{k=0}^{\infty} \frac {1}{4^k}\int_{0}^{\infty}\frac{x^ne^{-x}}{e^{kx}}(1+e^x)^{2k}dx$$

Was actually (wrongly) derived from this harmless looking integral:

$$I = \int_{0}^{\infty}\frac{x^n dx}{(e^x-1)^2}$$

Which can be tackled by:

$$I = \int_{0}^{\infty}\frac{x^ne^{-2x} dx}{1-2e^{-x}+e^{-2x}}$$ Then interpreting $2e^{-x} + e^{-2x}$ as the common ratio of a geometric progression with first term $=1$ $$\implies I = \int_{0}^{\infty}x^ne^{-2x}\sum_{r=0}^{\infty}e^{-rx}(2-e^{-x})^rdx$$

Now, after interchanging the integral and the summation sign, I applied the Binomial theorem on $(2-e^{-x})^r$

After which I interchanged the integral and summation again;

Which made $$I = \sum_{r=0}^\infty\sum_{p=0}^{r}\binom{r}{p}(-1)^p2^{r-p}\int_0^\infty x^n e^{-(p+r+2)x}dx$$

Now after a change of variable and applying the gamma function, I finally got:

$$I = n!\sum_{r=0}^\infty\sum_{p=0}^{r}\binom{r}{p}(-1)^p\frac{2^{r-p}}{(p+r+2)^{1+n}}$$

Which may not look like much, but after a close inspection, I found that $$I = n!(\zeta(n)-\zeta(1+n))$$

Keeping this result in mind, I wanted to evaluate the original integral, I, in a different way.

So I started by plugging in $x \rightarrow it$

$$\implies I = \int_{0}^{\infty}\frac{i^nt^ne^{-2it} (-i dt)}{(1-e^{-it})^2}$$ After some simplification:

$$I = \frac{-i^{1+n}}{4}\int_{0}^{\infty}\frac{t^ne^{-it} dt}{\sin^2(t/2)}$$ Writing $\sin^2(t/2) = 1-\cos^2(t/2)$ and interpreting $\cos^2t/2$ as the common ratio of a G.P with first term = $1$: $$\implies I = \frac{-i^{1+n}}{4}\int_{0}^{\infty}t^ne^{-it}\sum_{k=0}^{\infty}\cos^{2k}(t/2) dt$$ Interchanging the summation and integration sign and plugging $x = it$ back in this integral, I got:

$$I = \frac{1}{4}\sum_{k=0}^{\infty}\int_{0}^{\infty}x^ne^{-x}\cos^{2k}(ix/2) dx$$

After which I used $\cos x = (e^{ix}+e^{-ix})/2$ which led me to this: $$I = \frac{1}{4}\sum_{k=0}^{\infty} \frac {1}{4^k}\int_{0}^{\infty}\frac{x^ne^{-x}}{e^{kx}}(1+e^x)^{2k}dx$$ Where I got stuck.

My questions are:

1. How do I approach this integral using the second method I tried? (As my approach is definitely wrong.)
2. What actually happens when I substitute $x \rightarrow it$?

I just assumed that as $x \rightarrow \infty, t \rightarrow \infty$

But am I thinking about this the wrong way?

Edit: Different ways to approach this problem are appreciated, but please think about ways to solve this after substituting $x \rightarrow it$

Answer 1

Consider the generalized integral: $$J(a)=\int_0^{\infty} \frac{1}{e^{ax}-1} \; dx$$ Then, $$J'(a)=\int_0^{\infty} \frac{-xe^{ax}}{{\left(e^{ax}-1\right)}^2} \; dx$$ The integral you're looking for is quite similar to $J'(a)$ above. We want to construct a new integral, $I(a)$, so that $I'(1)$ is the integral you're looking for in your question. $$I(a)=-\int_0^{\infty} \frac{x^{n-1}e^{-x}}{{\left(e^{ax}-1\right)}} \; dx$$ $$I'(a)=\int_0^{\infty} \frac{x^{n}e^{-x}e^{ax}}{{\left(e^{ax}-1\right)}^2} \; dx$$ This works because $I'(1)$ equals the integral in your question. Okay, now begins the fun. We will begin to integrate $I(a)$: $$I(a)=-\int_0^{\infty} \frac{x^{n-1}e^{-x}e^{-ax}}{{\left(1-e^{-ax}\right)}} \; dx$$ $$I(a)=-\int_0^{\infty} x^{n-1}e^{-x(1+a)} \sum_{k=0}^{\infty} e^{-axk} \; dx$$ Here, we can interchange the summation and integral sign: $$I(a)=-\sum_{k=0}^{\infty} \int_0^{\infty} x^{n-1}e^{-x(1+a)} e^{-axk} \; dx$$ $$I(a)=-\sum_{k=0}^{\infty} \int_0^{\infty} x^{n-1}e^{-x(1+a+ak)} \; dx$$ Let $u=x(1+a+ak)$: $$I(a)=-\sum_{k=0}^{\infty} \frac{1}{{\left(1+a+ak\right)}^n} \int_0^{\infty} u^{n-1}e^{-u} \; du$$ Recognize this integral as the gamma function: $$I(a)=-(n-1)! \sum_{k=0}^{\infty} \frac{1}{{\left(1+a+ak\right)}^n} $$ Now, differentiate both sides with respect to $a$ to obtain $I'(a)$: $$I'(a)=-(n-1)! \sum_{k=0}^{\infty} \frac{-n(k+1)}{{\left(1+a+ak\right)}^{n+1}} $$ $$I'(a)=n! \sum_{k=0}^{\infty} \frac{k+1}{{\left(1+a+ak\right)}^{n+1}} $$ Plug in $a=1$ to find what the integral in your question evaluates to: $$I'(1)=n! \sum_{k=0}^{\infty} \frac{k+1}{{\left(k+2\right)}^{n+1}} $$ $$I'(1)=n! \sum_{k=0}^{\infty} \frac{k+2-1}{{\left(k+2\right)}^{n+1}} $$ $$I'(1)=n! \left(\sum_{k=0}^{\infty} \frac{1}{{\left(k+2\right)}^{n}} -\sum_{k=0}^{\infty} \frac{1}{{\left(k+2\right)}^{n+1}}\right)$$ $$I'(1)=n! \left(\sum_{k=1}^{\infty} \frac{1}{{k}^{n}} -1 -\sum_{k=1}^{\infty} \frac{1}{{k}^{n+1}}+1\right)$$ $$\boxed{I'(1)=\int_0^{\infty} \frac{x^n}{{\left(e^x-1\right)}^2} \; dx = n! \left(\zeta(n) -\zeta(n+1)\right)}$$

According to desmos, it appears that the integral converges for all $n >1$ and blows up at around $n=8$.

Answer 2

$$\int_{0}^{+\infty}\frac{x^n}{(e^x-1)^2}\,dx = \int_{0}^{+\infty} x^n e^{-x}\frac{e^x}{(e^x-1)^2}\,dx\stackrel{\text{IBP}}{=}\int_{0}^{+\infty}\frac{nx^{n-1}-x^n}{e^x(e^x-1)}\,dx $$ and since $\int_{0}^{+\infty}(n x^{n-1}-x^n)e^{-x}\,dx = n(n-1)!-n!=0$ and $\frac{1}{e^x(e^x-1)}=\frac{1}{e^x-1}-\frac{1}{e^x}$ we have $$\int_{0}^{+\infty}\frac{x^n}{(e^x-1)^2}\,dx = n\int_{0}^{+\infty}\frac{x^{n-1}}{e^x-1}\,dx - \int_{0}^{+\infty}\frac{x^n}{e^x-1}\,dx=n(n-1)!\zeta(n)-n!\zeta(n+1) $$ i.e. $\int_{0}^{+\infty}\frac{x^n}{(e^x-1)^2}\,dx = n!(\zeta(n)-\zeta(n+1))$ as was to be shown.

Answer 3

The first step of your second approach is invalid. You should have $\mathrm dx=i~\mathrm dt$ and the upper bound of your integral should be $-i\infty$. Notice that along real $t$, since $e^{it}=1$ at every multiple of $2\pi$, you have infinitely many poles. Furthermore, the rest of the integrand just grows unboundedly in magnitude.

Instead of expanding the denominator as you have, it is easier to use the derivative of the geometric series:

\begin{align}I(n)&=\int_0^\infty\frac{x^n}{(e^x-1)^2}~\mathrm dx\\&=\int_0^\infty\frac{x^ne^{-2x}}{(1-e^{-x})^2}~\mathrm dx\\&=\int_0^\infty x^ne^{-x}\left(\frac{\mathrm d}{\mathrm dx}\frac{-1}{1-e^{-x}}\right)~\mathrm dx\\&=\int_0^\infty x^ne^{-x}\left(\frac{\mathrm d}{\mathrm dx}\sum_{k=0}^\infty-e^{-kx}\right)~\mathrm dx\\&=\sum_{k=0}^\infty\int_0^\infty x^ne^{-x}\left(\frac{\mathrm d}{\mathrm dx}-e^{-kx}\right)~\mathrm dx\\&=\sum_{k=0}^\infty k\int_0^\infty x^ne^{-(k+1)x}~\mathrm dx\\&=\sum_{k=0}^\infty\frac k{(k+1)^{n+1}}\int_0^\infty x^ne^{-x}~\mathrm dx\tag{$x\mapsto x/(k+1)$}\\&=n!\sum_{k=0}^\infty\frac{(k+1)-1}{(k+1)^{n+1}}\\&=n![\zeta(n)-\zeta(n+1)]\end{align}

Answer 4

Consider a similar looking integral, with $a$ being a positive real:

\begin{equation} L(a)=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}+1} \,dx \hspace{.3cm}\text{if}\hspace{.3cm} n\in\mathbb{Z^{+}} \end{equation}

It would be nice to get rid of the $e^{ax}+1$ below, and we could do this with the geometric series, but it would not converge as it is, but $e^{-ax}+1$ would, then:

\begin{equation} L(a)=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}}\frac{1}{1+e^{-ax}} \,dx=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}}\frac{1}{1-(-e^{-ax})} \,dx \end{equation}

Plugging the geometric series leaves us with:

\begin{equation} L(a)=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}}\sum_{k=0}^{+\infty}(-e^{-ax})^{k} \,dx=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}}\sum_{k=0}^{+\infty}-e^{-axk} \,dx \end{equation}

After interchanging the sum and the integral, we have the following:

\begin{equation} L(a)=-\sum_{k=0}^{+\infty}\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}}e^{-axk} \,dx \end{equation}

\begin{equation} L(a)=-\sum_{k=0}^{+\infty}\int\limits_{0}^{+\infty} x^{n}e^{-ax(k+1)}\,dx \end{equation}

The integral almost has the form of the gamma function. Let $ax(k+1)=t$, which implies that: $dx=dt/a(k+1)$. The limits remain the same, thus:

\begin{equation} L(a)=-\sum_{k=0}^{+\infty}\int\limits_{0}^{+\infty} \left(\frac{t}{a(k+1)}\right)^{n}e^{-t}\,\frac{dt}{a(k+1)} \end{equation}

\begin{equation} L(a)=-\frac{1}{a^{n+1}}\sum_{k=0}^{+\infty}\frac{1}{(k+1)^{n+1}}\int\limits_{0}^{+\infty} t^{n}e^{-t}\,dt \end{equation}

The last integral is $\Gamma(n+1)$, which is equal to $n!$, then:

\begin{equation} L(a)=-\frac{n!}{a^{n+1}}\sum_{k=0}^{+\infty}\frac{1}{(k+1)^{n+1}} \end{equation}

Let $k+1=s$, which implies that $s(k=0)=1$. We now have that:

\begin{equation} L(a)=-\frac{n!}{a^{n+1}}\sum_{s=1}^{+\infty}\frac{1}{s^{n+1}} \end{equation}

From the definition of the Riemann zeta function, the last summation is equal to $\zeta(n+1)$, then:

\begin{equation} \boxed{L(a)=\int\limits_{0}^{+\infty} \frac{x^{n}}{e^{ax}+1} \,dx =-\frac{n!}{a^{n+1}}\zeta(n+1)} \end{equation}

If we differentiate once with respect to $a$, we obtain:

\begin{equation} L'(a)=\int\limits_{0}^{+\infty} \frac{x^{n+1}e^{ax}}{(e^{ax}+1)^{2}} \,dx =-(n+1)\frac{n!}{a^{n+2}}\zeta(n+1) \end{equation}

If we add and subtract $1$ in the $e^{ax}$ term, things get simpler:

\begin{equation} L'(a)=\int\limits_{0}^{+\infty} \frac{x^{n+1}(e^{ax}+1-1)}{(e^{ax}+1)^{2}} \,dx =-(n+1)\frac{n!}{a^{n+2}}\zeta(n+1) \end{equation}

\begin{equation} L'(a)=\int\limits_{0}^{+\infty} \frac{x^{n+1}(e^{ax}+1-1)}{(e^{ax}+1)^{2}} \,dx =\int\limits_{0}^{+\infty} \frac{x^{n+1}}{e^{ax}+1} \,dx-\underbrace{\int\limits_{0}^{+\infty} \frac{x^{n+1}}{(e^{ax}+1)^{2}} \,dx}_{I} \end{equation}

We have already derived the result for the first integral and the second integral is the one we want to compute, denoted as $I$, then:

\begin{equation} -(n+1)\frac{n!}{a^{n+2}}\zeta(n+1)=-\frac{(n+1)!}{a^{n+2}}\zeta(n+2)-I \end{equation}

\begin{equation} (n+1)\frac{n!}{a^{n+2}}\zeta(n+1)=\frac{(n+1)!}{a^{n+2}}\zeta(n+2)+I \end{equation}

Let $a=1$, then:

\begin{equation} I=\int\limits_{0}^{+\infty} \frac{x^{n+1}}{(e^{x}+1)^{2}} \,dx=\underbrace{(n+1)n!}_{(n+1)!}\,\zeta(n+1)-(n+1)!\,\zeta(n+2) \end{equation}

\begin{equation} I=\int\limits_{0}^{+\infty} \frac{x^{n+1}}{(e^{x}+1)^{2}} \,dx=(n+1)!\bigl[\zeta(n+1)-\zeta(n+2)\bigr] \end{equation}

If we set $n+1=s$, then we arrive at the desired result:

\begin{equation} \boxed{I=\int\limits_{0}^{+\infty} \frac{x^{s}}{(e^{x}+1)^{2}} \,dx=s!\bigl[\zeta(s)-\zeta(s+1)\bigr]} \end{equation}

with $s$ being a positive integer.