Integrate $I=\int_e^\infty\left(\frac{\log\log y}{y(y+1)}\right)^2dy$ using residue calculus?

Question

Can this integral be done using the residue calculus?

$$I=\int_e^\infty\left(\frac{\log\log y}{y(y+1)}\right)^2dy$$ ?

My (empirical) investigative attempts have been to use a keyhole contour centred at $y=e$. My working so far indicates that the larger circle as $R\to\infty$ vanishes, as does the small circle $\varepsilon\to 0$, both centred at $y=e$. I could be wrong about this but after my manipulations numerical evidence suggest this to be the case.

That leaves the integral I'm after over $e+\varepsilon\to\infty$, and the reverse integral $\infty\to e+\varepsilon$ just under the $x$-axis. Not sure if I can combine those to get $I$...

But, If can obtain $I$ from those two and supposing all this is correct so far, I tried to work out the residues at the poles of order $2$ of $y^2(y+1)^2=0$, which gives $y=0$ and $y=-1$. I figured if the residues are not defined at these points then I can't use contour integration anyway. Can the residues at these points be computed? Or maybe I've gone completely wrong.

Answer 1

Can this integral be done using the residue calculus ?

I'm afraid this goes way beyond that... First of all, for your integral to possess a closed form, the lower integration limit should be $\color{red}1.~($No, I cannot offer any rigorous proof of this statement, but, speaking from experience, this sort of definite integrals only “make sense” when evaluated over the extent of their “natural” domain – and no, I cannot offer any rigorous definition for the two quoted terms either$).~$ Now let $y~=~e^x,~$ and expand the integrand into partial fractions. We thus have $$I~=~\int_{\color{red}0}^\infty\bigg[\frac{\ln^2x}{e^x}-\frac{\ln^2x}{e^x+1}-\frac{\ln^2x}{(e^x+1)^2}\bigg]~dx~=~I_1-I_2-I_3$$ By twice differentiating under the integral sign with regard to the parameter n the first two of the following three expressions $$\begin{align}\int_0^\infty\frac{x^n}{e^x\pm0}~dx&~=~n!\\\\\int_0^\infty\frac{x^n}{e^x+1}~dx&~=~n!\cdot\eta(n+1)\\\\\int_0^\infty\frac{x^n}{e^x-1}~dx&~=~n!\cdot\zeta(n+1)\end{align}$$ we get $I_1~=~\gamma^2+\zeta(2)~=~\gamma^2+\dfrac{\pi^2}6,~$ and $I_2~=~\bigg(\gamma^2-\zeta(2)+2\gamma_1\bigg)\cdot\ln2-\dfrac{\ln^32}3,~$ where $\gamma$ represents the Euler-Mascheroni constant, $\gamma_1$ stands for the first Stieltjes constant, while $\zeta$ and $\eta$ indicate the Riemann $\zeta$ and Dirichlet $\eta$ function. As for $I_3,~$ it can be expressed by evaluating $J'(1),~$ where $$J(a)~=~\int_0^\infty\frac{\ln^2x}{e^x+a}~dx~=~\frac{\Big(\gamma^2+\zeta(2)\Big)\cdot\ln(1+a)+2\gamma\text{Li}'_1(-a)-\text{Li}''_1(-a)}a$$ The $\text{Li}'_1(\cdot)$ and $\text{Li}''_1(\cdot)$ represent differentiation with regard to the index rather than the actual argument. Letting $a\to1,~$ we have $\text{Li}'_1(-1)~=~\dfrac{\ln^22}2-\gamma\ln2,~$ and $~\text{Li}''_1(-1)~=~2\gamma_1\ln2~+$ $+~\gamma\ln^22-\dfrac{\ln^32}3.~$ These last two constants will appear in the expression of $J'(1),~$ along with $\text{Li}'_0(-1)~=~\dfrac{\ln2-\ln\pi}2,~$ and $~\text{Li}''_0(-1)~=~\dfrac{\gamma^2}2-\dfrac{\zeta(2)}4-\ln^22+2\ln2\cdot\ln(2\pi)-\dfrac{\ln^2(2\pi)}2~+$ $+~\gamma_1.~$ Adding them all together, we get the final expression for $-I_3~=~J'(1),~$ which, when further added to the difference of the first two partial integrals, yields the final result.