I have been trying to evaluate the following family of integrals:
$$ f:(0,\infty)^2 \rightarrow \mathbb{R} \, , \, f(\alpha,\beta) = \int \limits_0^\infty \frac{\ln (1+x^\alpha) \ln (1+x^{-\beta})}{x} \, \mathrm{d} x \, . $$
The changes of variables $\frac{1}{x} \rightarrow x$, $x^\alpha \rightarrow x$ and $x^\beta \rightarrow x$ yield the symmetry properties $$ \tag{1} f(\alpha,\beta) = f(\beta,\alpha) = \frac{1}{\alpha} f\left(1,\frac{\beta}{\alpha}\right) = \frac{1}{\alpha} f\left(\frac{\beta}{\alpha},1\right) = \frac{1}{\beta} f\left(\frac{\alpha}{\beta},1\right) = \frac{1}{\beta} f\left(1,\frac{\alpha}{\beta}\right) $$ for $\alpha,\beta > 0$ .
Using this result one readily computes $f(1,1) = 2 \zeta (3)$ . Then $(1)$ implies that $$ f(\alpha,\alpha) = \frac{2}{\alpha} \zeta (3) $$ holds for $\alpha > 0$ . Every other case can be reduced to finding $f(1,\gamma)$ for $\gamma > 1$ using $(1)$.
An approach based on xpaul's answer to this question employs Tonelli's theorem to write $$ \tag{2} f(1, \gamma) = \int \limits_0^\infty \int \limits_0^1 \int \limits_0^1 \frac{\mathrm{d}u \, \mathrm{d}v \, \mathrm{d}x}{(1+ux)(v+x^\gamma)} = \int \limits_0^1 \int \limits_0^1 \int \limits_0^\infty \frac{\mathrm{d}x \, \mathrm{d}u \, \mathrm{d}v}{(1+ux)(v+x^\gamma)} \, .$$ The special case $f(1,2) = \pi \mathrm{C} - \frac{3}{8} \zeta (3)$ is then derived via partial fraction decomposition ($\mathrm{C}$ is Catalan's constant). This technique should work at least for $\gamma \in \mathbb{N}$ (it also provides an alternative way to find $f(1,1)$), but I would imagine that the calculations become increasingly complicated for larger $\gamma$ .
Mathematica manages to evaluate $f(1,\gamma)$ in terms of $\mathrm{C}$, $\zeta(3)$ and an acceptably nice finite sum of values of the trigamma function $\psi_1$ for some small, rational values of $\gamma > 1$ (before resorting to expressions involving the Meijer G-function for larger arguments). This gives me some hope for a general formula, though I have not yet been able to recognise a pattern.
Therefore my question is:
How can we compute $f(1,\gamma)$ for general (or at least integer/rational) values of $\gamma > 1$ ?
Update 1:
Symbolic and numerical evaluations with Mathematica strongly suggest that $$ f(1, n) = \frac{1}{n (2 \pi)^{n-1}} \mathrm{G}_{n+3, n+3}^{n+3,n+1} \left(\begin{matrix} 0, 0, \frac{1}{n}, \dots, \frac{n-1}{n}, 1 , 1 \\ 0,0,0,0,\frac{1}{n}, \dots, \frac{n-1}{n} \end{matrix} \middle| \, 1 \right) $$ holds for $n \in \mathbb{N}$ . These values of the Meijer G-function admit an evaluation in terms of $\zeta(3)$ and $\psi_1 \left(\frac{1}{n}\right), \dots, \psi_1 \left(\frac{n-1}{n}\right) $ at least for small (but likely all) $n \in \mathbb{N}$ .
Interesting side note: The limit $$ \lim_{\gamma \rightarrow \infty} f(1,\gamma+1) - f(1,\gamma) = \frac{3}{4} \zeta(3) $$ follows from the definition.
Update 2:
Assume that $m, n \in \mathbb{N} $ are relatively prime (i.e. $\gcd(m,n) = 1$). Then the expression for $f(m,n)$ given in Sangchul Lee's answer can be reduced to \begin{align} f(m,n) &= \frac{2}{m^2 n^2} \operatorname{Li}_3 ((-1)^{m+n}) \\ &\phantom{=} - \frac{\pi}{4 m^2 n} \sum \limits_{j=1}^{m-1} (-1)^j \csc\left(j \frac{n}{m} \pi \right) \left[\psi_1 \left(\frac{j}{2m}\right) + (-1)^{m+n} \psi_1 \left(\frac{m + j}{2m}\right) \right] \\ &\phantom{=} - \frac{\pi}{4 n^2 m} \sum \limits_{k=1}^{n-1} (-1)^k \csc\left(k \frac{m}{n} \pi \right) \left[\psi_1 \left(\frac{k}{2n}\right) + (-1)^{n+m} \psi_1 \left(\frac{n + k}{2n}\right) \right] \\ &\equiv F(m,n) \, . \end{align} Further simplifications depend on the parity of $m$ and $n$.
This result can be used to obtain a solution for arbitrary rational arguments: For $\frac{n_1}{d_1} , \frac{n_2}{d_2} \in \mathbb{Q}^+$ equation $(1)$ yields \begin{align} f\left(\frac{n_1}{d_1},\frac{n_2}{d_2}\right) &= \frac{d_1}{n_1} f \left(1,\frac{n_2 d_1}{n_1 d_2}\right) = \frac{d_1}{n_1} f \left(1,\frac{n_2 d_1 / \gcd(n_1 d_2,n_2 d_1)}{n_1 d_2 / \gcd(n_1 d_2,n_2 d_1)}\right) \\ &= \frac{d_1 d_2}{\gcd(n_1 d_2,n_2 d_1)} f\left(\frac{n_1 d_2}{\gcd(n_1 d_2,n_2 d_1)},\frac{n_2 d_1}{\gcd(n_1 d_2,n_2 d_1)}\right) \\ &= \frac{d_1 d_2}{\gcd(n_1 d_2,n_2 d_1)} F\left(\frac{n_1 d_2}{\gcd(n_1 d_2,n_2 d_1)},\frac{n_2 d_1}{\gcd(n_1 d_2,n_2 d_1)}\right) \, . \end{align}
Therefore I consider the problem solved in the case of rational arguments. Irrational arguments can be approximated by fractions, but if anyone can come up with a general solution: you are most welcome to share it. ;)