Evaluating $\sum_{i=1}^{\infty} \sum_{j=1}^{\infty} \frac{1}{ij(i+j)^2}$

Question

I was playing around with double sums and encountered this problem: Evaluate $$\sum_{i=1}^{\infty} \sum_{j=1}^{\infty} \frac{1}{ij(i+j)^2}$$ It looks so simple I thought it must have been seen before, but after using ApproachZero I only found extremely similar variants of the problem:

$\sum_{j=1}^{\infty} \sum_{k=1}^{\infty} \frac{1}{j^2(j+k)^2}$ here

$\sum_{k=1}^{\infty} \sum_{n=1}^{\infty} \frac{1}{n^2k^2(n+k)^2}$ here

$\sum_{n,k \in \mathbb{N_1}} \frac{1}{n^2k^2(n+k)}$ here

Some general facts about $T(r,s,t)=\sum_{u,v=1}^{\infty} \frac{1}{u^rv^s(u+v)^t}$ here

While the last link looks promising, I don't believe that it answers my specific question(but I haven't read the paper fully, only skimmed).

I have an answer to this question(shared Q&A style), but I am curious as to whether there is another way to evaluate the sum. I have been experimenting with double sums for a few weeks now, so many techniques to attack "simple" sums such as this would be much appreciated! In terms of level, I have officially finished a Calculus 2 course but I have self-taught myself Calculus 3 and Complex Analysis, and worked through a few "integration problems" books so any further integral techniques would be fine.

Answer 1

My own answer to the above question.

We can do this by analyzing a particular $n=i+j$. Instead of summing in the row/column way(when the double sum is written out) we sum by the diagonals from the top row to the left column as shown below

where each box corresponds to a particular choice of i,j. Now, because $n=i+j$ will hit all integers from $2$ to $\infty$, we can rewrite the summation as $$S = \sum_{n=2}^{\infty} \frac{1}{n^2}f(n)$$ where $$f(n)=\sum_{i=1}^{n-1} \frac{1}{i(n-i)}$$ where we set $j$ to be $n-i$ because $i+j=n$. This is the "sum by diagonals" representation. Now, we can find an explicit formula for $f(n)$. Note that $$\frac{1}{i(n-i)}=\frac{1}{n}(\frac{1}{n-i}+\frac{1}{i})$$ so, we have that $$f(n)=\frac{1}{n} \sum_{i=1}^{n-1} \frac{1}{i}+\frac{1}{n-i}=$$ $$\frac{1}{n} (\sum_{i=1}^{n-1} \frac{1}{i}+ \sum_{i=1}^{n-1} \frac{1}{n-i})=$$ $$\frac{1}{n}(H_{n-1}+\sum_{i=1}^{n-1} \frac{1}{n-i})=$$ $$\frac{1}{n}(H_{n-1}+\sum_{i=1}^{n-1} \frac{1}{i})=$$ $$\frac{2 H_{n-1}}{n}$$ where the second-to-last equality follows from the reflection formula $\sum_{i=a}^{b} f(i)=\sum_{i=a}^{b} f(a+b-i)$. From this, we have that $$S=\sum_{n=2}^{\infty} \frac{1}{n^2} \cdot \frac{2 H_{n-1}}{n}=$$ $$2 \sum_{n=2}^{\infty} \frac{H_{n-1}}{n^3}=$$ $$2 \sum_{n=2}^{\infty} \frac{H_n-\frac{1}{n}}{n^3}=$$ $$2(\sum_{n=2}^{\infty} \frac{H_n}{n^3} - \sum_{n=2}^{\infty} \frac{1}{n^4})$$ Here, by noting that $H_n=\frac{1}{n}$ for $n=1$, we also have that $$S=2(\sum_{n=1}^{\infty} \frac{H_n}{n^3} - \sum_{n=1}^{\infty} \frac{1}{n^4})$$ as the $n=1$ terms cancel each other out. Now, the second summation is clearly $\zeta(4)=\frac{\pi^4}{90}$, so we have that $$S=2(\sum_{n=1}^{\infty} \frac{H_n}{n^3} - \frac{\pi^4}{90})=$$ $$2 \sum_{n=1}^{\infty} \frac{H_n}{n^3} - \frac{\pi^4}{45}$$ Now, all that is left is to evaluate $\sum_{n=1}^{\infty} \frac{H_n}{n^3}$. To do this, we use the method outlined in Problem 3.58 of "Limits, Series, and Fractional Part Integrals". Recall the Dilogarithm function $$\mathrm{Li_2}(z)=\sum_{n=1}^{\infty} \frac{z^n}{n^2} = - \int_{0}^{z} \frac{\ln(1-t)}{t} dt$$ For example, note that $\mathrm{Li_2}(1)=\zeta(2)=\frac{\pi^2}{6}$ and $\mathrm{Li_2}(0)=0$. It is a known result(proved in an earlier problem in the book) that $$\sum_{n=1}^{\infty} \frac{H_n}{n^3} = - \sum_{n=1}^{\infty} \frac{1}{n^2} \int_{0}^{1} x^{n-1} \ln(1-x) dx$$ and we know that $$- \sum_{n=1}^{\infty} \frac{1}{n^2} \int_{0}^{1} x^{n-1} \ln(1-x) dx=- \int_{0}^{1} \frac{\ln(1-x)}{x} (\sum_{n=1}^{\infty} \frac{x^n}{n^2}) dx=$$ $$- \int_{0}^{1} \frac{\ln(1-x)}{x} \mathrm{Li_2}(x) dx$$ Now, sub $u=\mathrm{Li_2}(x), du = \frac{-\ln(1-x)}{x} dx$, and using the special values of $\mathrm{Li_2}(x)$ above for the bounds, we see that $$\sum_{n=1}^{\infty} \frac{H_n}{n^3} = \int_{0}^{\frac{\pi^2}{6}} u du=$$ $$\frac{u^2}{2} \Big|_0^{\frac{\pi^2}{6}}=$$ $$\frac{\pi^4}{72}$$ So, we finish with $$S = 2 \cdot \frac{\pi^4}{72} - \frac{\pi^4}{45}=$$ $$\boxed{\frac{\pi^4}{180}}=$$ $$\boxed{\frac{1}{2} \zeta(4)}$$

Answer 2

$$\sum_{j=1}^{\infty} \frac{1}{i\,j\,(i+j)^2}=\frac{\gamma+\psi ^{(0)}(i+1)-i \psi ^{(1)}(i+1) }{i^3}$$

$$S_p=\sum_{i=1}^{p}\frac{\gamma}{i^3}=\gamma\, H_p^{(3)}$$

$$\sum_{i=1}^{p}\frac{\psi ^{(1)}(i+1)}{i^2}=\frac{1}{120} \left(\pi^4-60 \psi ^{(1)}(p+1)^2+10 \psi^{(3)}(p+1)\right)$$

$$\sum_{i=1}^{\infty}\frac{\psi ^{(0)}(i+1)}{i^3}=\frac{\pi ^4}{72}-\gamma \zeta (3)$$

Expanding as series

$$S_p=\frac{\pi ^4}{180}+\frac{1-\gamma }{2 p^2}+\frac{3 \gamma -4}{6 p^3}+\frac{13-6 \gamma }{24 p^4} +O\left(\frac{1}{p^5}\right)$$ which is in a relative error smaller than $1.0$% as soon as $p>17$ and smaller than $0.1$% as soon as $p>65$.

Edit

What is interesting is that $\gamma$ appears only in the even powers after $p^{-4}$. Isolating it, we have $$S_p=\frac{\pi ^4}{180}-\frac\gamma{2p^2} \Bigg(1-\frac{1}{p}+\frac{1}{2 p^2}-\frac{1}{6 p^4}+\frac{1}{6 p^6}+O\left(\frac{1}{p^{8}}\right) \Bigg)+$$ $$\frac1{2p^2}\Bigg(1-\frac{4}{3 p}+\frac{13}{12 p^2}-\frac{1}{2 p^3}-\frac{7}{180 p^4}+\frac{1}{5 p^5}+\frac{23}{630 p^6}-\frac{31}{126 p^7}+O\left(\frac{1}{p^{8}}\right) \Bigg)$$

Now, the absolute relative error is $0.020$ percent for $p=2$.

Answer 3

$$S=\sum_{n=1}^{\infty} \sum_{k=1}^{\infty} \frac{1}{nk(n+k)^2}$$

Using the standard result from Laplace transform,

$$\color{red}{\int_0^\infty x e^{-sx}\,dx=\frac{1}{s^2}}$$

$$\implies \int_0^\infty x e^{-(n+k)x}\,dx=\frac{1}{(n+k)^2}$$

$$S=\sum_{n=1}^{\infty} \sum_{k=1}^{\infty} \frac{1}{nk}\int_0^\infty x e^{-(n+k)x}\,dx$$

$$S=\sum_{n=1}^{\infty} \sum_{k=1}^{\infty} \frac{1}{nk}\int_0^\infty x e^{-nx}e^{-kx}\,dx$$

$$S=\int_0^\infty x\sum_{n=1}^{\infty} \frac{e^{-nx}}{n}\sum_{k=1}^{\infty} \frac{e^{-kx}}{k} \,dx$$

Since both summations are the same,

$$S=\int_0^\infty x\left(\sum_{n=1}^{\infty} \frac{e^{-nx}}{n}\right)^2 \,dx$$

Using the below result,

$$\sum_{n=1}^{\infty} \frac{z^n}{n} = -\ln(1 - z)\implies\sum_{n=1}^{\infty} \frac{e^{-n x}}{n} = -\ln(1 - e^{-x})$$

Squaring both sides,

$$\color{red}{\left(\sum_{n=1}^{\infty} \frac{e^{-nx}}{n}\right)^2 = \log^2(1-e^{-x})}$$

$$S=\int_0^\infty x\log^2(1-e^{-x})\,dx\underset{e^{-x}\to x}\implies -\int_{0}^{1}\frac{\log^2(1-x)\log(x)}{x}\,dx$$

$$S=-\int_{0}^{1}\frac{\log^2(1-x)\log(x)}{x}\,dx\underset{1-x\to x}\implies -\int_{0}^{1}\frac{\log^2(x)\log(1-x)}{1-x}\,dx$$

Using the below formula,

$$\color{red}{\sum_{n= 1}^\infty H_n x^n=-\frac{\log(1-x)}{1-x}}$$

$$S=\sum_{n= 1}^\infty H_n\int_0^1 x^n \log^2(x)\,dx$$

Using the below result, derived by integrating by parts,

$$\color{red}{\int_0^1 x^n \log^2(x)\,dx=\frac{2}{(n+1)^3}}$$

$$S=2\sum_{n= 1}^\infty \frac{H_n}{\left(n + 1\right)^{3}}$$

Using the below property,

$$\color{red}{\sum_{n\ge 1} \frac{H_n}{(n+1)^k} z^{n+1}=\sum_{n\ge 1} \frac{H_n}{n^k} z^{n}-\sum_{n\ge 1} \frac{1}{n^{k+1}} z^{n}}$$

Set $z=1,k=3$,

$$\sum_{n\ge 1} \frac{H_n}{(n+1)^3} =\sum_{n\ge 1} \frac{H_n}{n^3} -\sum_{n\ge 1} \frac{1}{n^4}$$

The first sum can be quickly looked up here and second sum from here

$$\color{red}{ \sum_{n=1}^{\infty}\frac{H_n}{n^q} = \frac12(q+2)\zeta(q+1)- \frac12\sum_{m=1}^{q-2}\zeta(m+1)\zeta(q-m)}$$

$$\sum_{n\ge 1} \frac{H_n}{(n+1)^3}=\frac{\pi^4}{72}-\frac{\pi^4}{90}=\frac{\pi^4}{360}$$

$$S=2\sum_{n\ge 1} \frac{H_n}{(n+1)^3}=\frac{\pi^4}{180}$$

$$\therefore \sum_{n=1}^{\infty} \sum_{k=1}^{\infty} \frac{1}{nk(n+k)^2}=\frac{\pi^4}{180}=\frac{\zeta(4)}{2}$$