While tinkering with Mathematica I found that it could evaluate sums of the form $$\sum\limits_{k = 1}^{\infty} \frac{H_k H_{k+m}}{k^2}$$ for some nonnegative integer $m$. As it turns out (or the evidence is overwhelmingly in favor of it), this sum is always a linear combination of a constant term and the first three $\zeta$-values, i.e. $\zeta(2)$ through $\zeta(4)$. The coefficients of the $\zeta$-values were not too hard to figure out with the help of the OEIS and a bit of intuition. We have $$\sum\limits_{k = 1}^{\infty} \frac{H_k H_{k+m}}{k^2} = -C -\left(\sum\limits_{k = 1}^{m}\frac{1}{k^2}\right)\zeta(2) + 2\left(\sum\limits_{k = 1}^{m}\frac{1}{k}\right)\zeta(3) + \frac{17}{4}\zeta(4)$$ with just the non-zeta term $C$ being undetermined. For $m = 1,\ldots,9$ this constant takes on the values $$C = 0,\frac{1}{4},\frac{4}{9},\frac{341}{576},\frac{679}{960},\frac{25921}{32400},\frac{19879}{22680},\frac{95594629}{101606400},\frac{182134073}{182891520}$$ which... isn't all that insightful, to say the least. Visually, this looks like

This looks like logarithmic growth which is typical for harmonic sums. But it does not seem to be a simple combination in this case.

Thanks to the comment by ho boon suan (and verification by Dr. Wolfgang Hintze) I was able to verify that indeed $$C = \frac{1}{2} \sum\limits_{j = 1}^{m} \frac{H_{j-1}^2+H_{j-1}^{(2)}}{j^2}$$

and thus

$$\sum\limits_{k = 1}^{\infty} \frac{H_k H_{k+m}}{k^2} = -\sum\limits_{j = 1}^{m} \underbrace{\frac{H_{j-1}^2+H_{j-1}^{(2)}}{2j^2}}_{=\frac{1}{j^2(j-1)!} \begin{bmatrix}j \\ 3\end{bmatrix}} -H_m^{(2)}\zeta(2) + 2H_m\zeta(3) + \frac{17}{4}\zeta(4)$$

As a side product of an auxiliary calculation, we also have the identity

$$\sum\limits_{k = 1}^{\infty} \frac{H_k}{k^2(k+j)} = \frac{H_{j-1}^2+H_{j-1}^{(2)}}{j^2} - \frac{1}{j^2} \zeta(2) + \frac{2}{j}\zeta(3)$$

for arbitrary $j \in \mathbb{N}$.

  • $\begingroup$ Whole thing is constant only. $\endgroup$
    – Z Ahmed
    Feb 9 at 13:48
  • $\begingroup$ Yes, let me rephrase that part a bit. $\endgroup$
    – TheOutZ
    Feb 9 at 13:48
  • 1
    $\begingroup$ Try use: $H_{k+m}=H_k+\sum _{j=1}^m \frac{1}{j+k}$ ? $\endgroup$ Feb 9 at 15:12
  • $\begingroup$ Well, this reduces the problem to evaluating $\sum_{k = 1}^{\infty} \frac{H_k}{k^3+jk^2}$, so it certainly helps but these new coefficients seem to be just as obscure as in the original problem. $\endgroup$
    – TheOutZ
    Feb 9 at 15:31
  • 1
    $\begingroup$ I believe this question gives a general method for evaluating sums of the form $\sum_{k\ge1}H_k/p(k)$, where $p(k)$ is an arbitrary polynomial, in terms of the digamma function $\psi$. $\endgroup$ Feb 9 at 16:25

1 Answer 1


EDIT 13.02.23

I am grateful to Ali Shadhar for a very helpful conversation which allowed me to improve part of the discussion section.

EDIT 11.02.23

I am grateful to Mariusz Iwaniuk who pointed out a sign mistake in $(3)$, which I have now corrected. The original author of the OP, TheOutZ, already had the term correct.

Original post

The sum

$$s(m)=\sum _{k=1}^{\infty } \frac{H_k H_{k+m}}{k^2}\tag{1}$$

can be simplified in terms of finite sums writing

$$H_{k+m}=\sum _{i=1}^m \frac{1}{i+k}+H_k\tag{2}$$

then expanding the product and doing the sums with the result

$$s(m) =-\zeta (2) H_m^{(2)}+2 \zeta (3) H_m+\frac{17 \zeta (4)}{4} -\frac{1}{2} \sum _{s=1}^m \frac{1}{s^2}( H_{s-1}{}^2+H_{s-1}^{(2)})\tag{3}$$

so that your "constant" is the rational number given by

$$C_m = \frac{1}{2} \sum _{s=1}^m \frac{1}{s^2}( H_{s-1}{}^2+H_{s-1}^{(2)})\tag{4}$$

Mathematica finds $C_\infty = \frac{7\pi^4}{360} \simeq 1.89407$


Some possibly useful relations are.

The combinations of $H$ in the summands in $C$ can be written as

$$c(n)=\left(H_{n}\right){}^2+H_{n}^{(2)} = 2 \sum_{i=1}^{n}\frac{1}{i^2} + 2\sum_{1\le i\lt j \le n}\frac{1}{i j}$$

its generating function

$$\sum_{n=1}^ {\infty}(\left(H_{n}\right){}^2+H_{n}^{(2)}) x^n =\frac{2 \text{Li}_2(x)}{1-x}+\frac{\log ^2(1-x)}{1-x} \tag{5}$$

Here we have combined

$$\sum_{n=1}^ {\infty}(H_{n}^{(2)}) x^n =\frac{ \text{Li}_2(x)}{1-x} $$


$$\sum_{n=1}^ {\infty}(H_n^2 - H_{n}^{(2)}) x^n=\frac{\log ^2(1-x)}{1-x} $$

(see Ali's book, 2.1.5)

Generalization to non-integer $m$

As a first test I tried to calculate $s(\frac {1}{2}) \simeq 5.01804$.

With the help of Ali and drawing from his book we could find the closed expression

$$s(\frac{1}{2}) = 8 \text{Li}_4\left(\frac{1}{2}\right)+4 \zeta (3)+\zeta (3) \log (8)-\frac{\pi ^4}{90}\\ +\frac{\log ^4(2)}{3}-\frac{1}{3} \pi ^2 \log ^2(2)-8 \log ^2(2)\tag{6}$$

Indeed, writing e.g. (https://math.stackexchange.com/a/3035646/198592) $$H_{n+\frac{1}{2}}=2 H_{2 n+1}-H_n-\log (4)\tag{7}$$

we encounter the in $s(\frac{1}{2})$ the sums over $ H_n H_{2n}/n^2, H_n^2/n^2, H_n/(n^2(2n+1)), H_n/n^2$ which are all solved in Ali's book.

For a general small $m=\alpha$ we can form a power series expansion of $s(\alpha)$ based on

$$\frac{\partial H_{n+\alpha }}{\partial \alpha }\text{/.}\, \alpha \to 0 = \frac{\pi ^2}{6}-H_n^{(2)}\tag{8}$$

leads to the sum

$$s_2 = \sum _{n=1}^{\infty } \frac{H_n H_n^{(2)}}{n^2}=\zeta(5)+\zeta(3)\zeta(2)\tag{9}$$

which was calculated in Ali's book, formula $(4.107)$.

Summarizing, the case of non-integer $m$ promises to lead to interesting problems in terms of sums (or the equivalent integrals).

The case of rational $\alpha$ can perhaps be solved.

  • $\begingroup$ God damn! Just noticed my error in the calculation and got the same result. Nice job! $\endgroup$
    – TheOutZ
    Feb 9 at 18:32
  • $\begingroup$ Also very interesting that the sum on the right hand side can be expressed as a weighted sum of the Stirling numbers of the first kind... $\endgroup$
    – TheOutZ
    Feb 9 at 18:50
  • $\begingroup$ I tried for m=2 aka. s(2) give me: $$3 \zeta (3)-\frac{5 \pi ^2}{24}+\frac{17 \pi ^4}{360}$$,but correct answer is:$$3 \zeta (3)-\frac{1}{4}-\frac{5 \pi ^2}{24}+\frac{17 \pi ^4}{360}$$ ? Mathematica code: -Zeta[2]*HarmonicNumber[m, 2] + 2 Zeta[3] HarmonicNumber[m] + 17/4*Zeta[4] - 1/2*Sum[1/ s^2*(HarmonicNumber[s - 1]^2 - HarmonicNumber[s - 1, 2]), {s, 1, m}] /. m -> 2 $\endgroup$ Feb 10 at 16:22
  • 1
    $\begingroup$ It should be: $c(m)=\frac{1}{2} \sum _{s=1}^m \frac{\left(H_{s-1}\right){}^2+H_{s-1}^{(2)}}{s^2}$ :) . $\endgroup$ Feb 10 at 16:39
  • 2
    $\begingroup$ @Dr.WolfgangHintze Using Mathematica with MZIntegrate package see here:math.stackexchange.com/questions/4495798/… $$\int_0^1 \frac{\text{Li}_3(1-x) \log (1-x)}{x} \, dx=\frac{\pi ^2 \zeta (3)}{3}-\frac{9 \zeta (5)}{2}$$ and $$\int_0^1 \frac{\text{Li}_3(x) \log (1-x)}{x} \, dx=\frac{\pi ^2 \zeta (3)}{6}-3 \zeta (5)$$ $\endgroup$ Feb 13 at 10:37

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .