<s>[INST] Evaluate $\int_0^1 \frac{\operatorname{arctanh}^3(x)}{x}dx$ I'm trying to evaluate the following integral:
$$\int_0^1 \frac{\operatorname{arctanh}^3(x)}{x}dx$$
I was playing around trying to numerically approximate the answer with known constants and found that the integral is almost exactly $\frac{\pi^4}{64}$. The integral seems to break down after about $11$ decimal places.
I have a suspicion that this stems from the integral:
$$\int_0^1 \frac{\operatorname{arctanh}(x)}{x}dx$$
since this is exactly equal to $\frac{\pi^2}{8}$. 
Also for the integral:
$$\int_0^1 \frac{\operatorname{arctanh}^5(x)}{x}dx$$
This is suspicously close to $\frac{\pi^6}{128}$ but not exactly.
For some reason the above integrals diverges slightly from the some from of $\frac{\pi^n}{2^m}$ for some $n$ and $m$.

My question is: Why does this happen? And what are the true values of those integrals?

 [/INST] As MrTaurho pointed out in the comments, we can rewrite $\,\displaystyle{\operatorname{arctanh}x=\frac12 \ln\left(\frac{1+x}{1-x}\right)}$, this gives:
$$I=\int_{0}^1 \frac{\operatorname{arctanh}(x)^3}{x}dx=-\frac18\int_0^1 \frac{\ln^3\left(\frac{1-x}{1+x}\right)}{x}dx$$
And we will substitute $\displaystyle{\frac{1-x}{1+x}=t\Rightarrow dx=-\frac{2}{(1+t)^2}dt}$ in order to get:
$$I=-\frac18\int_0^1 \frac{\ln^3 t}{\frac{1-t}{1+t}}\frac{2}{(1+t)^2}dt=-\frac14 \int_0^1 \frac{\ln^3 t}{1-t^2}dt=-\frac14 \sum_{n=0}^\infty \int_0^1 t^{2n}\ln^3t dt$$
Now we consider the following integral: $$\int_0^1 t^a dt=\frac{1}{a+1}\Rightarrow \int_0^1 t^a \ln^3 tdt=\frac{d^3}{da^3} \left(\frac{1}{a+1}\right)=-\frac{6}{(a+1)^4}$$
$$\Rightarrow I=\frac{6}{4}\sum_{n=0}^\infty \frac{1}{(2n+1)^4}=\frac32\cdot\frac{\pi^4}{96}=\frac{\pi^4}{64}$$

Of course this can be generalized for any power, so let's do that. Consider:
$$I(k)=\int_0^1 \frac{\text{arctanh}^kx}{x}dx=\frac{(-1)^k}{2^k}\int_0^1 \frac{\ln^k\left(\frac{1-x}{1+x}\right)}{x}dx\overset{\large\frac{1-x}{1+x}=t}=\frac{2(-1)^k}{2^k}\int_0^1 \frac{\ln^k t}{1-t^2}dt$$
$$=\frac{(-1)^k}{2^{k-1}}\sum_{n=0}^\infty \int_0^1 x^{2n}\ln^k xdx=\frac{(-1)^k}{2^{k-1}}\sum_{n=0}^\infty \frac{(-1)^k k!}{(2n+1)^{k+1}}=\frac{k!}{2^{k-1}} \sum_{n=0}^\infty \frac{1}{(2n+1)^{k+1}}$$
Where above we used the following result: $$\int_0^1 t^a dt=\frac{1}{a+1}\Rightarrow \int_0^1 t^a \ln^k tdt=\frac{d^k}{da^k} \left(\frac{1}{a+1}\right)=\frac{(-1)^k k!}{(a+1)^{k+1}}$$
And finally it reduces to: $$I(k)=\frac{k!}{2^{k-1}} \left(1-\frac{1}{2^{k+1}}\right)\zeta(k+1)=\boxed{\frac{k!\left(2^{k+1}-1\right)}{4^k}\zeta(k+1)}$$

One can verify the result by comparing to the one announced by Maple: $$I(5)=\int_0^1 \frac{\text{arctanh}^5 x}{x}dx=\frac{5!(2^6-1)}{4^5}\zeta(6)=\frac{945}{128}\cdot\frac{\pi^6}{945}=\frac{\pi^6}{128}$$
</s>