Infinite sum of harmonic number over polynomial of 2nd degreee

Question

Inspired by the recent problem How to evaluate the sum $\sum_{n=1}^{\infty}\frac{1}{(2n+1)(2n+2)}\left(1+\frac{1}{2}+...+\frac{1}{n}\right)$ I asked for a generalization, and I found it with Mathematica.

Can you prove (or disprove) this heuristic finding?

Let

$$s(p,q) =\sum _{k=1}^{\infty } \frac{H_k}{(k+p) (k+q)}\tag{1} $$

Then

$$s(p,q)=\frac{1}{2p-2q}\left(H_{p-1}^2-H_{q-1}^2-\psi ^{(1)}(p)+\psi ^{(1)}(q)\right)\tag{2}$$

for $p\neq q$ and

$$s(p,q)=H_{p-1} \psi ^{(1)}(p)-\frac{1}{2}\psi ^{(2)}(p)\tag{3}$$ for $q =p $

Answer 1

Assuming $q>p>2$, from $\sum_{n\geq 1}H_n x^n = \frac{-\log(1-x)}{1-x}$ we have

$$\begin{eqnarray*} s(p,q)=\sum_{n\geq 1}\frac{H_n}{(n+p)(n+q)} &=& \frac{1}{q-p}\int_{0}^{1}(x^{p-1}-x^{q-1})\frac{-\log(1-x)}{1-x}\,dx\\&\stackrel{\text{IBP}}{=}&\frac{1}{2(q-p)}\left[(q-1)A_{q-2}-(p-1)A_{p-2}\right] \end{eqnarray*}$$ where $$\begin{eqnarray*} A_m = \int_{0}^{1}x^m\log^2(1-x)\,dx&=&\left.\frac{d^2}{d\beta^2}\int_{0}^{1}x^m(1-x)^{\beta}\,dx\right|_{\beta=0}\\&=&\left.\frac{d^2}{d\beta^2}\left(\frac{\Gamma(m+1)\Gamma(\beta+1)}{\Gamma(m+\beta+2)}\right)\right|_{\beta=0}\\&=&\frac{\pi^2+6H_{m+1}^2-6\psi'(m+2)}{6(m+1)} \end{eqnarray*}$$ leads to $(p-1)A_{p-2}=\zeta(2)+H_{p-1}^2-\psi'(p)$ and to OP's $(2)$.
Then $\lim_{p\to q}\frac{\psi'(q)-\psi'(q)}{q-p}=\psi''(q)$ and $H_{p-1}=\psi(p)+\gamma$ easily prove $(3)$, too.

Answer 2

Self answer posted as a new question

As the self answer below is not really an answer to the OP but a non-trivial generalization of the OP I have posted the general formula as a new question Infinite sum of harmonic number over general polynomial

Original self answer (21.10.17)

In the meantime Jack d'Aurizio provided the proof which used the hints given by himself and reuns. Hence the problem of the OP is completely solved.

But I was not satisfied with my general formula for only a 2nd degree polynomial. So I attempted to extend it to a polynomial of arbitrary degree in $k$ using the same hints.

The attempt was successful and here is the result:

Let

$$s(p_1,\dotsc,p_n)=\sum _{k=1}^{\infty } \frac{H_k}{\prod _{i=1}^n (p_i +k)}\tag{A}$$

and

$$f(p_1,\dotsc,p_n)=\frac{1}{2} (-1)^n \sum _{i=1}^n \frac{h(p_i)}{\prod _{j=1,j\neq i}^n (p_i-p_j)}\tag{B}$$

where

$$h(a)=H_{a-1}^2-\psi ^{(1)}(a)+\zeta (2)\tag{C}$$

Then

$$s(p_1,\dotsc,p_n)=f(p_1,\dotsc,p_n)\tag{D}$$

In the case where two or more parameters $p_i$ are meant to be equal appropriate limits have to be taken of $f(p)$.

As a non trivial check and an application I rederived the Borwein-Borwein formulas (19) - (21) given in http://mathworld.wolfram.com/HarmonicNumber.html which consider sums of the type

$$\sum_{k=1}^\infty \frac{H_k}{k^n}$$

Many specific formulas can be written down using the relation established here.

Here is one more example for now

Setting $a=i$, $b=-i$ we obtain

$$\sum _{k=1}^{\infty } \frac{H_k}{k^2+1}=\frac{1}{4} \left(2 \gamma (1+\pi \coth (\pi ))+i \left(\psi ^{(0)}(-i)^2-\psi ^{(0)}(i)^2-\psi ^{(1)}(1-i)+\psi ^{(1)}(1+i)\right)\right)$$