Double harmonic sum $\sum_{n\geq 1}\frac{H^{(p)}_nH_n}{n^q}$ Are there any general formula for the following series 
$$\tag{1}\sum_{n\geq 1}\frac{H^{(p)}_nH_n}{n^q}$$
Where we define 
$$H^{(p)}_n= \sum_{k=1}^n \frac{1}{k^p}\,\,\,\,\,H^{(1)}_n\equiv H_n =\sum_{k=1}^n\frac{1}{k} $$
For the special case $p=q=2$ in (1) I found the following paper
Stating that 
$$\sum_{n\geq 1}\frac{H^{(2)}_nH_n}{n^2}=\zeta(5)+\zeta(2)\zeta(3)$$
See equation (3a) .
Is there any other paper in the literature discussing (1) or any special cases ?
 A: A couple of places to start for further looking:
The paper Further summation formulae related to generalized harmonic numbers by Zheng mentions the $p = q = 2$ case in example 2.3, and has several  related results, but not explicitly of the form you're looking for.
Appendix B of the paper On some log-cosine integrals related to $\zeta(3), \zeta(4)$ and $\zeta(6)$ by Mark Coffey has some similar things, and maybe looking through its references and papers that cite it would yield more.
A: Here we provide an answer for $p=1$ and arbitrary $q \ge 2$. Let us denote the quantity in question as follows:
\begin{equation}
{\mathfrak S}^{(p,1)}_q := \sum\limits_{m=1}^\infty \frac{H^{(p)}_m H_m}{m^q}
\end{equation}
In order to calculate the above we start from the generating function. We have:
\begin{eqnarray}
\sum\limits_{m=1}^\infty H^{(p)}_m H_m \frac{t^m}{m} &= & - \int\limits_0^t \frac{\log(1-\xi/t)}{\xi} \cdot \frac{Li_p(\xi)}{1-\xi} d\xi \\
&=& \left\{
\begin{array}{rr}
-\frac{1}{3} [\log(1-t)]^3 -\log(1-t) Li_2(t) + Li_3(t) & \mbox{for $p=1$} \\
\cdots
\end{array}
\right.
\end{eqnarray}
Now we use the common trick,i.e. we divide by $t$ and multiply by an appropriate power of $\log(1/t)$ and integrate appropriately. We have:
\begin{eqnarray}
{\mathfrak S}^{(1,1)}_q &=&
\underbrace{-\frac{1}{3} \int\limits_0^1 \frac{[\log(1/\xi)]^{q-2}}{(q-2)!} \cdot \frac{[\log(1-\xi)]^3}{\xi} d\xi}_{I_1} +\\
&&\underbrace{\int\limits_0^1 \frac{[\log(1/\xi)]^{q-2}}{(q-2)!} \cdot \frac{Li_1(\xi) Li_2(\xi)}{\xi} d\xi}_{I_2}+\\
&&\underbrace{\int\limits_0^1 \frac{[\log(1/\xi)]^{q-2}}{(q-2)!} \cdot \frac{ Li_3(\xi)}{\xi} d\xi}_{I_3}
\end{eqnarray}
Clearly all the inthree integrals in here have already been dealt with in MSE. The last integral is trivial . We have:
\begin{equation}
I_3= Li_{q+2}(+1)
\end{equation}
The second integral has been evaluated in An integral involving product of poly-logarithms and a power of a logarithm.. We have:
\begin{equation}
I_2= \frac{1}{2}\left[
\zeta(2) \zeta(q) + \sum\limits_{j=1}^{q-3} j \zeta(q-j-1) \zeta(j+3)-\sum\limits_{j=1}^{q-2} j {\bf H}^{(q-j-1)}_{j+3}(+1)
\right]
\end{equation}
where the quantities highlighted in bold font have been calculated in Calculating alternating Euler sums of odd powers .They all reduce to single zeta values for $q+2\le 7$ and otherwise involve a few additional two-dimensional zeta values. 
Finally the first integral has been evaluated in Compute an integral containing a product of powers of logarithms. .It always reduces to single zeta values. We have:
\begin{eqnarray}
&&I_1 = \\
&&\frac{(-1)^{q-1}}{(q-1)!} \left[
-\frac{1}{3} \Psi^{(q+1)}(1) + \frac{1}{2} \sum\limits_{j=1}^{q-2} \binom{q-1}{j} \left\{\Psi^{(j+1)}(1)\Psi^{(q-1-j)}(1) + \Psi^{(j+0)}(1)\Psi^{(q-j)}(1)\right\}+\right.\\
&&
\left.
-\frac{1}{3} \sum\limits_{1 \le j < j_1 \le q-2} \frac{(q-1)!}{j! (j_1-j)!(q-1-j_1)!} \Psi^{(j)}(1) \Psi^{(j_1-j)}(1) \Psi^{(q-1-j_1)}(1)
\right]
\end{eqnarray}
where $\Psi^{(j)}(1)= (-1)^{j+1} j! \zeta(j+1)$ for $j=1,2,3,\cdots$.
A: You can refer to the following papers 
$1$. Evaluations of some quadratic Euler sums
$2$. Euler sums and integrals of polylogarithm functions
$3$. Multiple zeta values and Euler sums
$4$. Tornheim type series and nonlinear Euler sums
A: In here we provide a generating function of the quantities in question. Let us define:
\begin{equation}
{\bf H}^{(p,r)}_q(t) := \sum\limits_{m=1}^\infty H_m^{(p)} H_m^{(r)} \frac{t^m}{m^q}
\end{equation}
In here we take $q\ge1$. We have:
\begin{eqnarray}
&&{\bf H}^{(p,1)}_q(t) = Li_p(1) \cdot \frac{1}{2} [\log(1-t)]^2 \cdot 1_{q=1}+\\
&&\frac{(-1)^{q}}{2} \sum\limits_{l=(q-2)}^{p+q-3}  \left(\binom{l}{q-2} 1_{l < p+q-3} + ({\mathcal A}^{(p)}_{q-2}) 1_{l=p+q-3}\right) \cdot \underbrace{\int\limits_0^1 \frac{[Li_1(t \xi)]^2}{\xi} Li_{l+1}(\xi)  \frac{[\log(1/\xi)]^{p+q-3-l}}{(p+q-l-3)!}d\xi}_{I_1}+\\
&& \frac{(-1)^{q-1}}{2} \sum\limits_{j=0}^{q-3} \left({\mathcal A}^{(p)}_{q-2-j}\right) \cdot \zeta(p+q-2-j) \underbrace{\int\limits_0^1 \frac{[Li_1(t \xi)]^2}{\xi} \frac{[\log(\xi)]^j}{j!} d\xi}_{I_2}+\\
&& \sum\limits_{l=1}^p  \underbrace{\int\limits_0^1 \frac{Li_q(t \xi)}{\xi} Li_l(\xi) \frac{[\log(1/\xi)]^{p-l}}{(p-l)!} d\xi}_{I_3}
\end{eqnarray}
Here $t\in (-1,1)$ and $p=1,2,\cdots$ and
\begin{equation}
{\mathcal A}^{(p)}_{q} := p+\sum\limits_{j=2}^{q} \binom{p+j-2}{j}
= p \cdot 1_{p=1} + \frac{p+q-1}{p-1} \binom{p+q-2}{q}\cdot 1_{p > 1}
\end{equation}
Note 1: The quantities in the right hand side all contain products of poly-logarithms and a power of logarithm. Those quantities, in principle, have been already dealt with in An integral involving product of poly-logarithms and a power of a logarithm. for example.
Note 2: Now that we have the generating functions we will find recurrence relations for the sums in question  and hopefully provide some closed form expressions .
Now we have:
\begin{eqnarray}
&&I_1 =\\
&&\sum\limits_{l_1=2}^{l+1} \binom{p+q-2-l_1}{l+1-l_1} (-1)^{l+1-l_1} \zeta(l_1) \left({\bf H}^{(1)}_{p+q-l_1}(t) - Li_{p+q+1-l_1}(t) \right)+\\
&&\sum\limits_{l_1=2}^{p+q-2-l} \binom{p+q-2-l_1}{l} (-1)^{l-1} \zeta(l_1) \left({\bf H}^{(1)}_{p+q-l_1}(t) - Li_{p+q+1-l_1}(t) \right)+\\
&&\sum\limits_{l_1=1}^{p+q-2-l} \binom{p+q-2-l_1}{l} (-1)^{l-0} \left( {\bf H}^{(l_1,1)}_{p+q-l_1}(t) - {\bf H}^{(l_1)}_{p+q+1-l_1}(t) \right)
\end{eqnarray}
and
\begin{eqnarray}
&&I_2=2 (-1)^j \left[ {\bf H}^{(1)}_{j+2}(t) - Li_{j+3}(t)\right]
\end{eqnarray}
and
\begin{eqnarray}
&&I_3=\\
&&\sum\limits_{l_1=2}^l \binom{p-l_1}{p-l}(-1)^{l-l_1} \zeta(l_1) Li_{p+q+1-l_1}(t) +\\
&& \sum\limits_{l_1=2}^{p-l+1} \binom{p-l_1}{l-1}(-1)^{l} \zeta(l_1) Li_{p+q+1-l_1}(t)-
\sum\limits_{l_1=1}^{p-l+1} \binom{p-l_1}{l-1} (-1)^l {\bf H}^{(l_1)}_{p+q+1-l_1}(t)
\end{eqnarray}
