It is easy to see that the integral is equivalent to
$$
\begin{align*}
\int_0^\infty \frac{1}{x\sqrt{2}+\sqrt{2x^2+1}}\frac{\log x}{\sqrt{1+x^2}}dx &= \sqrt{2}\int_0^\infty \frac{\sqrt{x^2+\frac{1}{2}}-x}{\sqrt{1+x^2}}\log x\; dx\tag{1}
\end{align*}
$$
This integral is a special case of the following generalised equation:
$$\begin{align*}\mathcal{I}(k) :&= \int_0^\infty \frac{\sqrt{x^2+k^2}-x}{\sqrt{1+x^2}}\log x\; dx \\ &= E'(k)-\left(\frac{1+k^2}{2} \right)K'(k)+\left(k^2 K'(k)-E'(k) \right)\frac{\log k}{2}+\log 2-1 \tag{2}\end{align*}$$
where $K'(k)$ and $E'(k)$ are complementary elliptic integrals of the first and second kind respectively.
Putting $k=\frac{1}{\sqrt{2}}$ in equation $(2)$,
$$
\begin{align*}
\mathcal{I}\left(\frac{1}{\sqrt{2}}\right)&=E'\left(\frac{1}{\sqrt{2} }\right)-\frac{3}{4}K'\left(\frac{1}{\sqrt{2}} \right)-\left\{\frac{1}{2} K'\left(\frac{1}{\sqrt{2}} \right)-E'\left(\frac{1}{\sqrt{2}} \right)\right\}\frac{\log 2}{4}+\log 2-1
\end{align*}
$$
Using the special values,
$$
\begin{align*}
E'\left(\frac{1}{\sqrt2} \right) &= \frac{\Gamma\left(\frac{3}{4} \right)^2}{2\sqrt\pi}+\frac{\sqrt{\pi^3}}{4\Gamma\left(\frac{3}{4} \right)^2}\\
K'\left(\frac{1}{\sqrt2} \right) &= \frac{\sqrt{\pi^3}}{2\Gamma\left(\frac{3}{4} \right)^2}
\end{align*}
$$
we get
$$
\mathcal{I}\left(\frac{1}{\sqrt{2}}\right)=\frac{1+\log\sqrt[4]2}{2\sqrt{\,\pi}}\Gamma\left(\frac34\right)^2-\frac{\sqrt{\,\pi^3}}8\Gamma\left(\frac34\right)^{-2}+(\log 2-1)\, \tag{3}
$$
Putting this in equation $(1)$, we get the answer that Cleo posted.
How to prove Equation $(2)$?
We begin with Proposition 7.1 of "The integrals in Gradshteyn and Ryzhik:
Part 16" by Boettner and Moll.
$$\int_0^\infty \frac{\log x}{\sqrt{(1+x^2)(m^2+x^2)}}dx = \frac{1}{2}K'(m)\log m$$
Multiplying both sides by $m$ and integrating from $0$ to $k$:
$$
\begin{align*}
\int_0^\infty \frac{\sqrt{x^2+k^2}-x}{\sqrt{1+x^2}}\log x\; dx &= \frac{1}{2}\int_0^k m K'(m)\log(m)\; dm
\end{align*}
$$
The result follows since
$$\begin{align*} \int m K'(m)\log(m)\; dm &= 2E'(m)-\left(1+m^2 \right)K'(m)+\left(m^2 K'(m)-E'(m) \right)\log m\\ &\quad +\text{constant} \tag{4}
\end{align*}$$
One can verify equation $(4)$ easily by differentiating both sides with respect to $m$ and using the identities
$$
\begin{align*}
\frac{dE'(k)}{dk}&= \frac{k}{k^{'2}}(K'(k)-E'(k))\\
\frac{dK'(k)}{dk}&= \frac{k^2 K'(k)-E^{'}(k)}{kk^{'2}}
\end{align*}
$$
Start by splitting the integral into two convergent parts:
$$I=\int_0^\infty \frac{\tanh x-x+x-x\operatorname{sech} x}{x^3}dx=\int_0^\infty \frac{\tanh x-x}{x^3}dx+\int_0^\infty \frac{1-\operatorname{sech} x}{x^2}dx$$
$$I_1=\int_0^\infty \frac{\tanh x-x}{x^3}dx\overset{IBP}=\frac12 \int_0^\infty \frac{\operatorname{sech}^2 x -1}{x^2}dx$$
$$\overset{IBP}=-\int_0^\infty \frac{\tanh x\operatorname{sech}^2 x}{x}dx=-\frac{7\zeta(3)}{\pi^2}$$
$$I_2=\int_0^\infty \frac{1-\operatorname{sech} x}{x^2}dx\overset{IBP}=\int_0^\infty \frac{\tanh x\operatorname{sech} x}{x}dx=\frac{4G}{\pi}$$
$$I_2=\int_0^\infty \frac{\tanh x \operatorname{sech} x}{x}dx\overset{x=\ln t}=2\int_1^\infty \frac{(t^2-1)}{(t^2+1)\ln t}dt\overset{t=\frac{1}{x}}=\int_0^\infty \frac{x^2-1}{(x^2+1)^2\ln x}dx$$
Above the two integrals were averaged after the reciprocal subtitution was done.
Now we will use Feynman's trick alongside beta function: $$I(a)=\int_0^\infty \frac{x^a-1}{(x^2+1)^2 \ln x}dx\Rightarrow I'(a)=\int_0^\infty \frac{x^a}{(x^2+1)^2}dx=\frac12 \left(\frac{1-a}{2}\right)\frac{\pi}{\sin\left(\frac{\pi(a+1)}{2}\right)}$$
$$I(0)=0\Rightarrow I_2=\frac{\pi}{4}\int_0^2 \frac{1-a}{\sin\frac{\pi(a+1)}{2}}da =\frac{2}{\pi} \int_{0}^\frac{\pi}{2}\frac{t}{\sin t}dt=\frac{4 G}{\pi}$$
$I_1$ can be found here.
Best Answer
\begin{align}\mathcal{I}&=\frac{1}{2}\int_0^\infty \frac{x^n}{\cosh(x)}\mathrm dx\\&=\int_0^\infty \frac{x^n}{e^x+e^{-x}}\mathrm dx\\&=\int_0^\infty \frac{x^ne^{-x}}{1+e^{-2x}}\mathrm dx.\end{align}
Proposition: $$\int_0^\infty \frac{x^{s-1}e^{-x}}{1+e^{-2x}}\mathrm dx=\beta(s)\Gamma(s),$$ where $\beta(s)$ is the Dirichlet beta function defined as $\beta(s)=\sum_{n=0}^\infty \frac{(-1)^n}{(2n+1)^s}$.
Proof:\begin{align}\int_0^\infty x^{s-1}\left(\frac{e^{-x}}{1-(-e^{-2x})}\right)\mathrm dx&=\int_0^\infty x^{s-1}\left(\sum_{k=0}^\infty (-1)^ke^{-(2k+1)x}\right)\mathrm dx \text{, using geometric series}\\&=\sum_{k=0}^\infty (-1)^k\int_0^\infty x^{s-1} e^{-(2k+1)x}\mathrm dx\\&=\sum_{k=0}^\infty \frac{(-1)^k}{(2k+1)^s}\displaystyle\int_0^\infty x^{s-1}e^{-x}\mathrm dx\text{, substituting $(2k+1)x\mapsto x$}\\&=\beta(s)\Gamma(s). \end{align}
Therefore, $$\mathcal{I}=\beta(n+1)\Gamma(n+1).$$
We can evaluate for different values of $n\ge 0$. For $n=1$, $\frac{1}{2}\int_0^\infty x\operatorname{sech}(x)=\beta(2)\Gamma(2)=G$, where $G$ is Catalan's constant. For $n=2$, $\int_0^\infty x^2\operatorname{sech}(x)=2\beta(3)\Gamma(3)=\frac{\pi^3}{8}$