To show your identity holds, make the one-dimensional change of variables $x=ty$ to get
\begin{align*}
\int_{0}^\infty \int_0^\infty \frac{y \ln(y) \ln (x)}{(y^2+1)(x^2+y^2)} \ dy \ dx &= \int_{0}^\infty \int_0^\infty \frac{y^2 \ln(y) \ln (ty)}{(y^2+1)(t^2+1) y^2} \ dy \ dt \\
&=\int_{0}^\infty \int_0^\infty \frac{ \ln(y) \ln (t)}{(y^2+1)(t^2+1)} \ dy \ dt \\
&\qquad + \int_{0}^\infty \int_0^\infty \frac{ \ln^2(y)}{(y^2+1)(t^2+1)} \ dt \ dy \\
&= \left(\int_0^\infty \frac{ \ln (t)}{t^2+1} \ dt \right)^2 + \frac{\pi}{2} \int_0^\infty \frac{\ln^2(y)}{y^2+1} \ dy.
\end{align*}
Note upon a change of variables $t=1/u,$ we have
\begin{align*}
\int_0^1\frac{ \ln (t)}{t^2+1} \ dt &= -\int_1^\infty \frac{ \ln (t)}{t^2+1} \ dt,
\end{align*}
which implies $\int_0^\infty\frac{ \ln (t)}{t^2+1} \ dt =0.$ Hence,
\begin{align*}
\int_{0}^\infty \int_0^\infty \frac{y \ln(y) \ln (x)}{(y^2+1)(x^2+y^2)} \ dy \ dx &= \frac{\pi}{2} \int_0^\infty \frac{\ln^2(y)}{y^2+1} \ dy.
\end{align*}
In order to evaluate the right hand side, consider the triple integral
$$I=\int_0^\infty \int_0^\infty \int_0^\infty \frac{xz}{(x^2+1)(z^2+1)(x^2y^2z^2+1)} \ dy \ dx \ dz.$$ By carrying out the integration directly, we can find $I=\frac{\pi^3}{8}.$ On the other hand, using Fubini's Theorem to integrate with respect to $z$ first, we get by partial fractions (or Mathematica)
\begin{align*}
I &=\int_0^\infty \int_0^\infty \int_0^\infty \frac{xz}{(x^2+1)(z^2+1)(x^2y^2z^2+1)} \ dz \ dx \ dy \\
&= \int_0^\infty \int_0^\infty \frac{x \ln \left(xy\right)}{\left(x^2+1\right) \left( x^2 y^2-1\right)} \ dx \ dy \\
&= \int_0^\infty \int_0^\infty \frac{x \ln \left(x\right)}{\left(x^2+1\right) \left( x^2 y^2-1\right)} \ dy \ dx + \int_0^\infty \int_0^\infty \frac{x \ln \left(y\right)}{\left(x^2+1\right) \left( x^2 y^2-1\right)} \ dx \ dy \\
&= \int_0^\infty \int_0^\infty \frac{x \ln \left(y\right)}{\left(x^2+1\right) \left( x^2 y^2-1\right)} \ dx \ dy \\
&= \int_0^\infty \frac{\ln^2(y)}{y^2+1} \ dy,
\end{align*}
in which the first double integral in the third to last equality has Cauchy Principal Value $0$ and thus vanishes, and we carried out the double integral with respect to $x$ in the second to last equality using partial fractions once again.
Hence, $$I=\int_0^\infty \frac{\ln^2(y)}{y^2+1} \ dy=\frac{\pi^3}{8},$$ and so the desired integral is equal to $$\frac{\pi}{2} I=\frac{\pi^4}{16}.$$
Remark: These type of integral calculations are prominent in my joint AMS paper with Daniele Ritelli. See here: https://www.ams.org/journals/qam/2018-76-03/S0033-569X-2018-01499-3/