Integrals of the form $\int_0^1\text{arctanh}^a(y) y^{2b} dy$

Question

Recently I have been trying to calculate integrals of the form: $$ I(a,b)=\int_0^1\text{arctanh}^a(y) y^{2b} dy$$ for some positive integer-valued $a$ and $b$. The values $a=0$ or 1 have quite trivial solutions. For example, with $a=1$, I can simply expand the arctanh function to get: $$ I(1,b) = \sum_{k=0}^\infty \int_0^1 \frac{y^{2(k+b)+1}}{2k+1} =\sum_{k=0}^\infty\frac{1}{(2b+2k+2)(2k+1)}=\frac{H_{b}-2\ln2}{2(2b+1)}, $$ where $H_b$ is the Harmonic number. For $a=2$ or $a=3$, however, the situation is quite different. The first a few values provided by Mathematica are: $$I(2,0)=\frac{\pi^2}{12},\quad I(2,1) = \frac{1}{3}+\frac{\pi^2}{36},\quad I(2,2) = \frac{1}{3}+\frac{\pi^2}{60},\\ \quad I(2,3) = \frac{14}{45}+\frac{\pi^2}{84},\quad I(2,4) = \frac{818}{2835}+\frac{\pi^2}{108} $$ $$I(3,0) =\frac{9\zeta(3)}{8},\quad I(3,1) =\ln2+\frac{3\zeta(3)}{8},\quad I(3,2) =\frac{1}{20}+\ln2+\frac{9\zeta(3)}{40},\quad I(3,4)=\frac{17}{210}+ \frac{14}{15}\ln2 +\frac{9}{56}\zeta(3) $$ However, I had no idea how to prove even the simple case $I(2,0)$, which I imagine should be related to the zeta function since $\zeta(2) = \pi^2/6$. I also don't know whether there's any general relation between those rational numbers and $b$. Any thoughts and ideas are highly appreciated.

---

Here are some of my ideas to prove $I(2,b)$, which may or may not be helpful. We have the Maclaurin series: $$\text{arctanh}(y)^2=\sum_{k=0}^\infty \left(\sum_{n=0}^k\frac{1}{2n+1} \right)\frac{y^{2(k+1)}}{k+1} $$ Therefore, integrating the above series multiplied by $y^{2b}$ gives: $$I(2,b)= \sum_{k=0}^\infty \left(\sum_{n=0}^k\frac{1}{2n+1} \right)\frac{1}{(k+1)(2b+2k+3)}, $$ which should somehow simplify to get the above analytical values. However, I could not find a way to simplify the double sum. I would imagine this becomes a triple sum when $a=3$ where everything quickly turns intractable...

Answer 1

Here is a procedure to evaluate the integral systematically. Denote

\begin{align} J(m,n) &= \int_0^1 y^n \text{arctanh}^m y \ dy\\ K(m,n) &= \int_0^1 \frac{y^n \text{arctanh}^m y }{1+y}\ dy\\ \end{align} and then establish the following recursive formulas \begin{align} K(m,n)&=J(m,n-1)-K(m,n-1)\\ J(m,n) &= \frac m{n+1}K(m-1,n) +\frac n{n+1}J(m,n-1) \end{align} along with the starting values \begin{align} & K(0,n) =\int_0^1 \frac{y^n}{1+y}dy = (-1)^n \ln2+\sum_{k=0}^{n-1}\frac{(-1)^k}{n-k}\\ \end{align} For $m=1$ \begin{align} &J(1,0)=K(0,0)=\ln2\\ & J(1,1)=\frac12 K(0,1)+\frac12J(1,0)=\frac12\\ & J(1,2)=\frac13 K(0,2) +\frac23 J(1,1)=\frac16+\frac13\ln2 \\ & J(1,3)=\frac14 K(0,3) + \frac34 J(1,2)=\frac13 \\ & J(1,4)=\frac15 K(0,4)+\frac45 J(1,3)=\frac3{20}+\frac15\ln2\\ \\ &K(1,0)=\int_0^1 \frac{\text{arctanh} \ y}{1+y}dy =\frac{\pi^2}{24}\\ &K(1,1)=J(1,0)-K(1,0)=\ln2-\frac{\pi^2}{24} \\ &K(1,2)=J(1,1)-K(1,1)=\frac12-\ln2+ \frac{\pi^2}{24} \\ &K(1,3)=J(1,2)-K(1,2)= -\frac13+\frac43\ln2-\frac{\pi^2}{24} \\ &K(1,4)=J(1,3)-K(1,3)=\frac23-\frac43\ln2+ \frac{\pi^2}{24}\\ \end{align} For $m=2$ \begin{align} &J(2,0)=2K(1,0)= \frac{\pi^2}{12}\\ &J(2,1)=K(1,1)+\frac12 J(2,0)=\ln2\\ &J(2,2)=\frac23 K(1,2) +\frac23 J(2,1)=\frac13+ \frac{\pi^2}{36}\\ &J(2,3)=\frac12 K(1,3) + \frac34 J(2,2)=\frac1{12}+\frac23\ln2\\ &J(2,4)=\frac25 K(1,4)+\frac45 J(2,3)=\frac13+\frac{\pi^2}{60}\\ \\ &K(2,0)=\int_0^1 \frac{\text{arctanh}^2 y}{1+y}dy =\frac{3}{8}\zeta(3) \\ &K(2,1)=J(2,0)-K(2,0)=\frac{\pi^2}{12} -\frac{3}{8}\zeta(3)\\ &K(2,2)=J(2,1)-K(2,1)=\ln2- \frac{\pi^2}{12}+ \frac{3}{8}\zeta(3) \\ &K(2,3)=J(2,2)-K(2,2)= \frac13-\ln2+\frac{\pi^2}{9} -\frac{3}{8}\zeta(3)\\ &K(2,4)=J(2,3)-K(2,3)=-\frac14+\frac53 \ln2-\frac{\pi^2}{9}+ \frac{3}{8}\zeta(3)\\ \end{align}

Answer 2

$\newcommand{\d}{\,\mathrm{d}}\newcommand{\artanh}{\operatorname{artanh}}$A practical (without closed form) algorithm for inductively determining the solution for even $a$. It's just mechanical enough that I could code a python script to do this. The method is very similar to this except the residues at infinity are harder to get in general (I think).

By using the logarithm with $0\le\arg<2\pi$ we have an analytic branch of $\artanh:\Bbb C\setminus[-1,1]\to\Bbb C,z\mapsto\frac{1}{2}\log\left(\frac{1+z}{1-z}\right)$, where for $-1<x<1$ we have the branch jumps: $$\lim_{t\to0^+}\artanh(x+it)=\artanh(x)\\\lim_{t\to0^-}\artanh(x+it)=\artanh(x)+\pi i$$

Now consider an integrand of type $f(z)=\artanh(z)^a\cdot z^{2b}$ for positive integer $a,b$. This thing is going to be analytic on $\Bbb C\setminus[-1,1]$ but it will have a nontrivial residue at infinity. We have: $$\artanh(z)=\frac{1}{2}\log\left(1-\frac{2}{1-\frac{1}{z}}\right)=\frac{1}{2}\log(-1-2z^{-1}-z^{-2}-\cdots)\\=\frac{\pi i}{2}+\sum_{n=1}^\infty\frac{(-1)^{n-1}\cdot 2^{n-1}}{n}\left(\sum_{m=1}^\infty\frac{1}{mz^m}\right)^n,\quad\quad\text{ for sufficiently large $z$}$$And we are faced with the dauting task of extracting the $z^{-(2b+1)}$th coefficient of this raised to the power $a$. However, in practice this is a finite computation that we can always do. Assume we can do this; call this coefficient $\varsigma$.

Then it is true by integrating $\artanh(z)^{a+1}z^{2b}$ around a counterclockwise oriented dogbone about $[-1,1]$ that, bearing in mind a branch jump $(?)^{a+1}-(?+\pi i)^{a+1}$ gives us something binomially: $$2\pi i\cdot\varsigma_{a+1,b}=-2\pi i\operatorname{Res}(\infty)=2\sum_{\,\,k=0\\k\text{ even}}^a\binom{a+1}{k}(\pi i)^{a+1-k}I_{k,b}$$Thus: $$I_{2a',b}=\frac{1}{2a'+1}\left(\varsigma_{2a'+1,b}-\sum_{m=0}^{a'-1}\binom{2a'+1}{2m}(-1)^{a'-m}\pi^{2(a'-m)}I_{2m,b}\right)$$For even $a=2a'$.

For a really straightforward example: let's determine $I_{2,0}$ with this method. Cubing the asymptotic series for $\artanh(z)$ I find: $$\artanh(z)^3=-\frac{\pi^3i}{8}-\frac{3\pi^2}{4z}+O(z^{-2})$$Thus: $$I_{2,0}=\frac{1}{3}\left(-\frac{3\pi^2}{4}-(-1)\pi^2\cdot1\right)=\frac{\pi^2}{12}$$