Help with the integral $\int_{0}^{\infty}\frac{y^{2}e^{y}}{e^{sy}+e^{-sy}-2}dy$ I want to do the integral :
$$I(s)=\int_{0}^{\infty}\frac{y^{2}e^{y}}{e^{sy}+e^{-sy}-2} \, \mathrm{d}y$$
$s$ being a complex parameter. I tried expanding the dominator of the integrand, but this way we lose the symmetry $I(s)=I(-s)$. I tried converting the integral into a contour integral, but I don't know how to close the contour. 
 A: According to Mathematica:

$$I(s) = \frac{2 s \, \psi_{1}\left(-\frac{1}{s}\right)-\psi_{2}\left(-\frac{1}{s}\right)}{s^4}, \quad \mathrm{Re}(s) > 1,  \quad \psi(z) = \Gamma'(z) / \Gamma(z), \quad \psi_n = \psi^{(n)}$$

A: For the integral
\begin{align}
I(s) = \int_{0}^{\infty} \frac{y^{2} \, e^{y} \, dy}{e^{sy} - 2 + e^{-sy}}
\end{align}
it is evident that when $s \to -s$ the result yields $I(-s) = I(s)$. For the evaluation of the integral consider the following.
By making use of
\begin{align}
\frac{1}{(1-x)^{2}} = \sum_{n=0}^{\infty} (n+1) \, x^{n}
\end{align}
then
\begin{align}
I(s) &= \int_{0}^{\infty} \frac{y^{2} \, e^{y} \, dy}{(e^{sy/2} - e^{-sy/2})^{2}} \\
&= \int_{0}^{\infty} \frac{y^{2} \, e^{-(s-1)y} \, dy}{(1 - e^{-sy})^{2}} \\
&= \sum_{n=0}^{\infty} (n+1) \, \int_{0}^{\infty} e^{-(sn + s-1)y} \, y^{2} \, dy \\
&= \sum_{n=0}^{\infty} \frac{n+1}{(sn + s -1)^{3}} \, \int_{0}^{\infty} e^{-u} \, u^{2} \, du \\
&= \frac{2}{s^{3}} \, \sum_{n=0}^{\infty} \frac{n+1}{(n+p)^{3}}  \hspace{10mm} p = 1 - \frac{1}{s} \\
&= \frac{2}{s^{3}} \left[ \sum_{n=0}^{\infty} \frac{1}{(n+p)^{2}} + \frac{1}{s} \sum_{n=0}^{\infty} \frac{1}{(n+p)^{3}} \right]  \\
&= \frac{1}{s^{4}} \left[ 2s \, \psi^{(1)}\left(1 - \frac{1}{s}\right) - \psi^{(2)}\left(1 - \frac{1}{s}\right) \right]
\end{align}
where $\psi^{(m)}(x)$ are $m^{th}$ derivative of the digamma function. 
