Integral $\int_0^1 \frac{\ln(1-x)\ln(1+x^2)}{x}dx$

I am trying to solve by a different approach the fourth sum from here, namely: $$S= \sum_{n=1}^{\infty} \sum_{m=1}^{\infty} \frac{1}{nm(4n+m)} =\int_0^1 \frac{\ln(1-x)\ln(1-x^4)}{x}dx= \frac{67}{32} \zeta(3) -\frac{\pi}{2}G$$

One way to solve it is similarly to my answer from there: $$S=\int_0^1 \frac{\ln(1-x)\ln(1-x^2)}{x}dx+\int_0^1 \frac{\ln(1-x)\ln(1+x^2)}{x}dx$$ From here we know that: $$\small \int_0^1 \frac{[m\ln(1+x)+n\ln(1-x)][q\ln(1+x)+p\ln(1-x)]}{x}dx=\left(\frac{mq}{4}-\frac{5}{8}(mp+nq)+2np\right)\zeta(3)$$ Thus by setting $m=0,n,p,q=1$ in the first integral we get that: $$S=\frac{11}{8}\zeta(3)+\int_0^1 \frac{\ln(1-x)\ln(1+x^2)}{x}dx=\frac{11}{8}\zeta(3)+I$$ $$I=\sum_{n=1}^\infty \frac{(-1)^{n-1}}{n} \int_0^1 x^{2n-1} \ln(1-x)dx=\frac12\sum_{n=1}^\infty \frac{(-1)^n H_{2n}}{n^2}=\frac{23}{32}\zeta(3)-\frac{\pi}{2}G$$ And the result for $S$ follows. The last sum appears to be known, see $(659)$ from here, or alternatively since $I=2\Re\left( S(i)\right)$ just use the following identity: $$S(x)=\sum_{n=1}^\infty \frac{x^n}{n^2}H_n=\operatorname{Li}_3(x)-\operatorname{Li}_3(1-x)+\operatorname{Li}_2(1-x)\ln(1-x)+\frac{1}{2}\ln x \ln^2(1-x)+\zeta(3)$$ However I am trying to find a different method since the result is quite nice and I believe there's a nicer way to solve the integral without using such sums.

Thus I would appreciate to get some help with the following problem:

Prove without using Euler's sum or polylogs that $$\int_0^1 \frac{\ln(1-x)\ln(1+x^2)}{x}dx=\frac{23}{32}\zeta(3)-\frac{\pi}{2}G$$

I also tried to consider the following integral: $$J=\int_0^1 \frac{\ln(1+x)\ln(1+x^2)}{x}dx$$ $$\Rightarrow I+J=\int_0^1 \frac{\ln(1-x^2)\ln(1+x^2)}{x}dx\overset{x^2=t}=\frac12 \int_0^1\frac{\ln(1-t)\ln(1+t)}{t}dt=-\frac{5}{16}\zeta(3)$$

So now I am after the following integral: $$I-J=\int_0^1 \frac{\ln\left(\frac{1-x}{1+x}\right)\ln(1+x^2)}{x}dx=\frac74 \zeta(3)-\pi G$$

\begin{align}I&=\int_0^1 \frac{\ln(1-x)\ln(1+x^2)}{x}dx\\ &=\Big[\ln x\ln(1-x)\ln(1+x^2)\Big]_0^1+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx-\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx\\ &=\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx-\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx \end{align} Let $R$ the function defined for $[0;1]$ by, \begin{align} R(x)&=\int_0^x \frac{2t\ln t}{1+t^2}\,dt\\ &=\int_0^1 \frac{2tx^2\ln(tx)}{1+t^2x^2}\,dt \end{align}For $0<A<1$, \begin{align}\int_0^A \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=\Big[R(x)\ln(1-x)\Big]_0^A+\int_0^A \frac{R(x)}{1-x}\,dx\\ &=R(A)\ln(1-A)+\int_0^A \left(\int_0^1\frac{2tx^2\ln(tx)}{(1-x)(1+t^2x^2)}\,dt\right)\,dx\\ &=R(A)\ln(1-A)+\int_0^1 \left(\int_0^A\frac{2tx^2\ln t}{(1-x)(1+t^2x^2)}\,dx\right)\,dt+\\ &\int_0^A \left(\int_0^1\frac{2tx^2\ln x}{(1-x)(1+t^2x^2)}\,dt\right)\,dx\\ &=R(A)\ln(1-A)-\int_0^1 \frac{\ln t\ln(1+A^2t^2)}{(1+t^2)t}\,dt-2\int_0^1\frac{\ln t\arctan t }{1+t^2}\,dt-\\ &2\ln(1-A)\int_0^1 \frac{t\ln t}{1+t^2}\,dt+\int_0^A \frac{\ln x\ln(1+x^2)}{1-x)}\,dx \end{align} Take the limit at $A=1$, \begin{align}\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=-\int_0^1 \frac{\ln t\ln(1+t^2)}{(1+t^2)t}\,dt-2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx\end{align}In the first integral perform the change of variable $y=x^2$, \begin{align}\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=-\frac{1}{4}\int_0^1 \frac{\ln t\ln(1+t)}{(1+t)t}\,dt-2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx\\ &=\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{1+x}\,dx-\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{x}\,dx-2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\\ &\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx \end{align}In the second integral perform integration by parts, \begin{align}\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{1+x}\,dx+\frac{1}{8}\int_0^1\frac{\ln^2 x}{1+x}\,dx-2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\\ &\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx\\ &=\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{1+x}\,dx-\frac{1}{8}\int_0^1\frac{2x\ln^2 x}{1-x^2}\,dx+\frac{1}{8}\int_0^1\frac{\ln^2 x}{1-x}\,dx-\\ &2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx \end{align} In the second integral perform the change of variable $y=x^2$, \begin{align}\int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{1+x}\,dx+\frac{3}{32}\int_0^1\frac{\ln^2 x}{1-x}\,dx-\\ &2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx\\ &=\frac{1}{4}\int_0^1\frac{\ln x\ln(1+x)}{1+x}\,dx+\frac{3}{16}\zeta(3)-2\int_0^1\frac{\ln t\arctan t}{1+t^2}\,dt+\\ &\int_0^1 \frac{\ln x\ln(1+x^2)}{1-x}\,dx \\J&=\int_0^1 \frac{\ln(1+x)\ln x}{1+x}\\ A&=\int_0^1 \frac{\ln^2 x}{1-x^2}\,dx\\ &=\int_0^1 \frac{\ln^2 x}{1-x}\,dx-\int_0^1 \frac{x\ln^2 x}{1-x^2}\,dx \end{align}In the latter integral perform the change of variable $y=x^2$: \begin{align}A&=\int_0^1 \frac{\ln^2 x}{1-x}\,dx-\frac{1}{4}\int_0^1 \frac{\ln^2 x}{1-x}\,dx\\ &=\frac{7}{8}\int_0^1 \frac{\ln^2 x}{1-x}\,dx\\ &=\frac{7}{4}\zeta(3) \end{align}On the other hand, perform the change of variable $y=\dfrac{1-x}{1+x}$, \begin{align}A&=\frac{1}{2}\int_0^1 \frac{\ln^2\left(\frac{1-x}{1+x}\right) }{x}\,dx\\ B&=\frac{1}{2}\int_0^1 \frac{\ln^2\left(1-x^2\right) }{x}\,dx \end{align}In the latter integral perform the change of variable $y=1-x^2$,\begin{align}B&=\frac{1}{4}\int_0^1 \frac{\ln^2 x}{1-x}\,dx\\ B&=\frac{1}{2}\zeta(3)\\ A+B&=\int_0^1 \frac{\ln^2\left(1-x\right) }{x}\,dx+\int_0^1 \frac{\ln^2\left(1+x\right) }{x}\,dx\\ &=\int_0^1 \frac{\ln^2\left(1-x\right) }{x}\,dx+\Big[\ln x\ln(1+x)^2\Big]_0^1-2\int_0^1 \frac{\ln(1+x)\ln x}{1+x}\,dx \end{align}In the first integral perform the change of variable $y=1-x$,\begin{align}A+B&=\int_0^1 \frac{\ln^2 x}{1-x}\,dx-2J\end{align}But,\begin{align}A+B&=\frac{9}{4}\zeta(3)\end{align}Therefore,\begin{align}J&=\boxed{-\dfrac{1}{8}\zeta(3)}\\ K&=\int_0^1 \frac{\ln x\arctan x}{1+x^2}\,dx\\ 2K&=\int_0^1 \frac{\ln x\arctan x}{1+x^2}\,dx-\int_1^\infty \frac{\ln x\arctan\left(\frac{1}{x}\right)}{1+x^2}\,dx\\ &=\int_0^\infty \frac{\ln x\arctan x}{1+x^2}\,dx+\frac{\pi}{2}\int_0^1 \frac{\ln x}{1+x^2}\,dx\\ &=\int_0^\infty \frac{\ln x\arctan x}{1+x^2}\,dx-\frac{1}{2}\text{G}\pi \end{align} Let $S$ the function defined on $[0;\infty]$ by, \begin{align} S(x)&=\int_0^x\frac{\ln t}{1+t^2}\,dt\\ &=\int_0^1\frac{x\ln(tx)}{1+t^2x^2}\,dt \end{align}Observe that, \begin{align}S(0)&=0,\lim_{x\rightarrow \infty} S(x)=0\\ \int_0^\infty \frac{\ln x\arctan x}{1+x^2}\,dx&=\Big[S(x)\arctan x\Big]_0^\infty-\int_0^\infty \frac{S(x)}{1+x^2}\,dx\\ &=-\int_0^\infty\left(\int_0^1 \frac{x\ln(tx)}{(1+x^2)(1+t^2x^2)}\,dt\right)\,dx\\ &=-\int_0^1\left(\int_0^\infty \frac{x\ln t}{(1+x^2)(1+t^2x^2)}dx\right)dt-\int_0^\infty\left(\int_0^1 \frac{x\ln x}{(1+x^2)(1+t^2x^2)}dt\right)dx\\ &=A-\int_0^\infty \frac{\ln x\arctan x}{1+x^2}\,dx \end{align} Therefore,\begin{align} \int_0^\infty \frac{\ln x\arctan x}{1+x^2}\,dx&=\frac{7}{8}\zeta(3)\\ K&=\boxed{\frac{7}{16}\zeta(3)-\frac{1}{4}\text{G}\pi}\\ \int_0^1 \frac{2x\ln x\ln(1-x)}{1+x^2}\,dx&=\frac{1}{2}\text{G}\pi-\frac{23}{32}\zeta(3)+\int_0^1\frac{\ln x\ln(1+x^2)}{1-x}\,dx\\ I&=\boxed{\frac{23}{32}\zeta(3)-\frac{1}{2}\text{G}\pi} \end{align} NB: I assume, \begin{align}\int_0^1 \frac{\ln^2 x}{1-x}\,dx=2\zeta(3)\end{align} I have computed $\displaystyle \int_0^1 \frac{\ln x\ln(1+x)}{1+x}\,dx$ using only univariate changes of variable and performing integration by parts.

PS: $\text{I}$ is linked to Evaluate $\int_{0}^{\frac{\pi}{2}}\frac{x^2}{ \sin x}dx$ see: https://math.stackexchange.com/a/2716753/186817

Here is my attempt to compute $I-J$ using complex analysis method (I hope you don't mind.) Let $f(z) $ be an analytic function defined by $\displaystyle f(z) = \frac{2\log(1+z)\log(1+iz)}z$ on the unit disk. By Cauchy's integral theorem, we have that \begin{align*} \int_{[0,i]} f(z)dz - \int_{[0,1]} f(z) dz =& \int_{[1,i]} f(z)dz. \end{align*} Then the LHS is ($[0,i]$ is parametrized by $z = ix, x\in [0,1]$) \begin{align*} \int_{[0,i]} f(z)dz - \int_{[0,1]} f(z) dz =&\int_0^1 \frac{2\log(1+ix)\log(1-x)}{x} dx -\int_0^1 \frac{2\log(1+x)\log(1+ix)}x dx \\ =& \color{red}{\int_0^1 \frac{2\log\left(\frac{1-x}{1+x}\right)\log(1+ix)}x dx}. \end{align*} On the other hand, the RHS is ($[1,i]$ is parametrized by $z = e^{i\theta}, \theta \in [0,\frac\pi 2]$) \begin{align*} \int_{[1,i]} f(z)dz =&2i \int_0^{\frac\pi 2} \log(1+e^{i\theta})\log(1+ie^{i\theta}) d\theta \\ =&\color{blue}{2i\int_0^{\frac \pi 2}\Big[\log\left(2\cos(\theta/2)\right) + i\theta/2\Big]\Big[\log(2\cos\left(\theta/2 +\pi /4\right)+i(\theta/2+\pi/4)\Big]d\theta } \end{align*} where we have used $\log(1+e^{i\theta}) = \log(2\cos (\theta/2)) + i\theta/2$ for $|\theta|<\pi$.

Note that for all real $x$, it holds that $2\Re[\log(1+ix)] = \ln(1+x^2)$. So by equating the real parts of $\color{red} {\text{red}}$ and $\color{blue} {\text{blue}}$ integrals, we get \begin{align*} I-J =& \Re\left[\int_0^1 \frac{2\log\left(\frac{1-x}{1+x}\right)\log(1+ix)}x dx\right]\\ =&-\int_0^{\frac \pi 2} \left(\theta+ \frac \pi 2\right)\log(2\cos (\theta/2)) -\int_0^{\frac \pi 2} \theta \log(2\cos(\theta/2 + \pi /4)) d\theta\\ =& -\int_0^{\frac \pi 2} \left(\theta+ \frac \pi 2\right)\log(2\cos (\theta/2)) -\int_0^{\frac \pi 2} \left(\frac \pi 2 -\theta\right) \log(2\sin(\theta/2)) d\theta \\ =& \int_0^{\frac \pi 2} \theta \log (\tan (\theta/2))d\theta -\frac \pi 2\left(\int_0^{\frac \pi 2} \log(2\cos(\theta/2)) d\theta+\int_0^{\frac \pi 2} \log(2\sin(\theta/2)) d\theta\right)\\ =& \int_0^{\frac \pi 2} \theta \log (\tan (\theta/2))d\theta \end{align*} because \begin{align*} \int_0^{\frac \pi 2} \log(2\cos(\theta/2)) d\theta+\int_0^{\frac \pi 2} \log(2\sin(\theta/2)) d\theta =& \int_0^{\frac \pi 2} \log(2\cos(\theta/2)) d\theta+\int_{\frac \pi 2}^\pi \log(2\cos(\theta/2)) d\theta\\ =& \int_0^{\pi } \log(2\cos(\theta/2)) d\theta \\=& 0. \end{align*} Finally, using the Fourier series of $\displaystyle \log\left(\tan\left(\theta/2\right)\right) = \sum_{k=1}^\infty \frac{(-1)^k-1}{k}\cos(k\theta)$, we get \begin{align*} I - J =&\sum_{k=1}^\infty \frac{(-1)^k-1}{k}\int_0^{\frac \pi 2}\theta\cos(k\theta)d\theta\\ =&\sum_{k=1}^\infty \frac{(-1)^k-1}{k}\left(\frac{\pi\sin(k\pi /2)}{2k}+\frac{\cos(k\pi /2) - 1}{k^2}\right)\\ \overset{k=2j+1}=&\sum_{j=0}^\infty \left[\frac{\pi (-1)^{j+1}}{(2j+1)^2} +\frac 2{(2j+1)^3}\right]\\ =& -\pi \text{G} + \frac 7 4\zeta(3). \end{align*}