Calculating Harmonic Sums with residues. Evaluate: 
$$\sum_{n=1}^{\infty} \frac{H_n}{(n+1)^2}$$
A user stated: "most of the time sum up the residues of $(\gamma+\psi(z))^2\cdot r(z)$. To determine the residues, just expand the digamma function as a Laurent series about the positive integers." where $r(z)$ is the rational function.
Firstly, what is the Laurent series of $\psi(z)$??
And how would I find the residues? Any help is appreciated.
 A: This is a partial answer to your question and I do not know how to get the Laurent series of $\psi(z)$. Note that 
\begin{eqnarray} \sum_{n=1}^\infty\frac{H_n}{(n+1)^2}&=&\sum_{n=1}^\infty\frac{1}{(n+1)^2}\left(H_n+\frac{1}{n+1}\right)-\sum_{n=1}^\infty\frac{1}{(n+1)^3}\\
&=&\sum_{n=1}^\infty\frac{H_{n+1}}{(n+1)^2}-\sum_{n=1}^\infty\frac{1}{(n+1)^3}\\
&=&\sum_{n=1}^\infty\frac{H_{n}}{n^2}-\sum_{n=1}^\infty\frac{1}{n^3}\\
&=&\sum_{n=1}^\infty\frac{H_{n}}{n^2}-\zeta(3).\\
\end{eqnarray}
The latter series can be evaluated using this generating function
$$\sum_{n=1}^\infty \frac{H_n\, x^n}{n^2}=\zeta(3)+\frac{\ln^2(1-x)\ln(x)}{2}+\ln(1-x)\operatorname{Li}_2(1-x)+\operatorname{Li}_3(x)-\operatorname{Li}_3(1-x)$$
Taking its limit as $x\to1$, we get
\begin{eqnarray} \sum_{n=1}^\infty\frac{H_n}{n^2}=2\zeta(3)\end{eqnarray}
and hence
$$ \sum_{n=1}^\infty\frac{H_n}{(n+1)^2}=\zeta(3).$$
A: Start from the singular expansion of the digamma function:
$$\psi(z+1)+\gamma = \sum_{n\ge 1} \frac{z}{n(n+z)}
= \sum_{n\ge 1} \left(\frac{1}{n}-\frac{1}{z+n}\right).$$
This immediately  produces the Laurent expansion  at negative integers
$-q$ where $q\ge 1$:
$$\psi(z+1)+\gamma = -\frac{1}{z+q}
+ H_{q-1} + \cdots$$
This is because the sum minus the singular term is
$$\frac{1}{q} +
\sum_{n=1, n\ne q} \left(\frac{1}{n}-\frac{1}{(z+q)+n-q}\right)
\\ = \frac{1}{q} +
\sum_{n=1, n\ne q} 
\left(\frac{1}{n}-\frac{1}{n-q}\frac{1}{(z+q)/(n-q)+1}\right)$$
The constant term here is
$$\frac{1}{q} +
\sum_{n=1, n\ne q} 
\left(\frac{1}{n}-\frac{1}{n-q}\right)$$
The outside  term and the first  term inside the sum  iterate over the
inverses  of all natural  $n.$ The  second term  inside the  sum first
collects  $$\sum_{n=1}^{q-1}  \left(-\frac{1}{n-q}\right)$$  which  is
$H_{q-1}$  and thereafter iterates  over the  inverses of  all natural
$n,$  canceling  the  first   contribution,  for  a  final  result  of
$H_{q-1}.$
This immediately implies that
$$\mathrm{Res}_{z=-(q-1)} (\psi(z)+\gamma)^2
= -2 H_{q-1}.$$
The sum being evaluated is
$$\sum_{n=1}^\infty \frac{H_n}{(n+1)^2}
= \sum_{n=2}^\infty \frac{H_{n-1}}{n^2}.$$
Since the residue formula also holds for $q=1$ we can re-write this as
$$\sum_{n=1}^\infty \frac{H_{n-1}}{n^2}.$$
Now integrating  $$f(z) = \frac{1}{z^2}(\psi(z)+\gamma)^2$$  along a
large  circle of  radius $R$  with $R$ going to infinity and passing
between the poles we  pick up exactly  $-2$ times the  value of  the
sum and the residue at zero.
But at zero we have the expansion
$$\psi(z+1)+\gamma  = -\sum_{k\ge 1} \zeta(k+1) (-1)^k z^k$$
which yields
$$\psi(z)+\gamma = -\frac{1}{z} 
- \sum_{k\ge 1} \zeta(k+1) (-1)^k z^k.$$
This immediately implies that
$$\mathrm{Res}_{z=0} \frac{1}{z^2}(\psi(z)+\gamma)^2
= 2\zeta(3).$$
Because the integral on the circular contour vanishes in the limit (an
exercise that is left to  the reader that follows from the observation
that $2\pi R \times (\log R)^2/R^2$ vanishes in the limit) we get
$$2\zeta(3) - 2 \sum_{n=1}^\infty \frac{H_{n-1}}{n^2} = 0$$
or
$$\sum_{n=1}^\infty \frac{H_{n-1}}{n^2} = \zeta(3).$$
The digamma formulae are from the
Wikipedia entry.
