Compute $\sum\limits_{n=1}^\infty\frac{1}{(n(n+1))^p}$ where $p\geq 1$

I was recently told to compute some integral, and the result turned out to be a scalar multiple of the series $$\sum\limits_{n=1}^\infty\frac{1}{(n(n+1))^p},$$ where $p\geq 1$. I know it converges by comparison for $$\dfrac{1}{(n(n+1))^p}\leq\dfrac{1}{n(n+1)}<\dfrac{1}{n^2},$$ and we know thanks to Euler that $$\sum\limits_{n=1}^\infty\frac1{n^2}=\frac{\pi^2}6.$$ I managed to work out the cases where $p=1$ and $p=2$. With $p=1$ being a telescoping sum, and my solution for $p=2$ being $$\frac13\pi^2-3,$$ which I obtained based on Euler's solution to the Basel Problem. I see no way to generalize the results to values to arbitrary values of $p$ however. Any advice on where to start would be much appreciated.

Also, in absence of another formula, is the series itself a valid answer? Given that it converges of course.

Solution 1:

I was doing this exercise not a long ago with a student of mine, nice thing I have a chance to discuss it here. Assuming $p\in\mathbb{N}^+$, we just have to perform a partial fraction decomposition of $\frac{1}{x^p(x+1)^p}$. By stars and bars we have $\frac{1}{(1+x)^p}=\sum_{n\geq 0}\binom{p+n-1}{p-1}(-1)^n x^{n}$, hence the singularity at the origin gives a contribution equal to $$ \sum_{n=0}^{p-1}\binom{p+n-1}{p-1}\frac{(-1)^n}{x^{p-n}}=(-1)^p\sum_{n=1}^{p}\binom{2p-n-1}{p-1}\frac{(-1)^{n}}{x^{n}} $$ to the partial fraction decomposition of $\frac{1}{x^p(x+1)^p}$. By symmetry, the contribution provided by the singularity at $x=-1$ equals $$ (-1)^p\sum_{n=1}^{p}\binom{2p-n-1}{p-1}\frac{1}{(x+1)^n}$$ hence $\sum_{x\geq 1}\frac{1}{x^p(x+1)^p}$ splits as the following linear combination of values of the $\zeta$ function at even integers and telescopic series:

$$(-1)^p\sum_{k=1}^{\lfloor p/2\rfloor}\binom{2p-2k-1}{p-1}(2\zeta(2k)-1)-(-1)^p\sum_{k=0}^{\lfloor p/2\rfloor}\binom{2p-2k-2}{p-1} $$ which can be simplified as $$\boxed{\sum_{n\geq 1}\frac{1}{n^p(n+1)^p}= 2(-1)^p\sum_{k=1}^{\lfloor p/2 \rfloor}\binom{2p-2k-1}{p-1}\zeta(2k)-\frac{(-1)^p}{2}\binom{2p}{p}.}$$

Solution 2:

$\newcommand{\msc}[2]{\left(\binom{#1}{#2}\right)}$ Assume that $p$ is a positive integer, and write $\msc{n}{k}$ for the multiset coefficient. Then by the method of Heaviside we have the partial fraction decomposition $$\frac{1}{n^p(n+1)^p}=(-1)^p\sum_{k=1}^{p}\msc{p}{p-k}\left(\frac{(-1)^k}{n^k}+\frac{1}{(1+n)^k}\right)$$ So that $$\begin{split}\sum_{n=1}^{\infty}\frac{1}{n^p(n+1)^p}&=(-1)^{p-1}\msc{p}{p-1} + (-1)^p\sum_{k=2}^p\msc{p}{p-k}(-1+(1+(-1)^k)\zeta(k))\\ &=(-1)^{p-1}\sum_{k=1}^p\msc{p}{p-k} + 2(-1)^p\sum_{k=1}^{\lfloor p/2\rfloor}\msc{p}{p-2k}\zeta(2k) \end{split}$$ $$\boxed{\sum_{n=1}^{\infty}\frac{1}{n^p(n+1)^p}=(-1)^{p-1}\sum_{k=1}^p\msc{p}{p-k} + 2(-1)^p\sum_{k=1}^{\lfloor p/2\rfloor}\msc{p}{p-2k}\frac{(2\pi)^{2k}\lvert B_{2k}\rvert }{2(2k)!}}$$ where $\zeta(s)$ denotes the Riemann zeta function and $B_n$ denotes the Bernoulli numbers.