A difficult improper integral 
What is the value of the $$\int_{1}^{\infty}\frac{dx}{ x^2(e^x+1)}?$$

I got  $2$ as a result.  Wolfram has given $0.11111$ as a result. I am confused.
Is it correct? Or incorrect ? Please tell me...
Thanks in advance
 A: The result can't be $2$ because
$$\int_{1}^{\infty}\frac{dx}{ x^2(e^x+1)}\leq \frac{1}{e+1}\int_{1}^{\infty}\frac{dx}{x^2}=\frac{1}{e+1}<\frac{1}{3}<2.$$
Note that 
$$\int_{1}^{\infty}\frac{dx}{ x^2(e^x+1)}=\int_{1}^{\infty}\frac{e^{-x}}{ x^2(1+e^{-x})}dx=\sum_{k=1}^{\infty}(-1)^{k-1}\int_{1}^{\infty}\frac{e^{-kx}}{x^2}\,dx$$
which can be written by using the exponential integral. 
WA gives an approximate value of 0.119161046435885.I don't think that it has a closed form.
A: As an alternative expression, we have
$$ \int_{1}^{\infty} \frac{dx}{x^2(e^x + 1)}
= C + \int_{0}^{1} \left( \frac{1}{e^x + 1} \right)'' \log x \, dx $$
where $C$ is given by
\begin{align*}
C
&= \frac{1}{4} + \frac{1}{1+e} + \frac{\gamma}{4} + \frac{\log 2}{3} - 3\log A \\
&\approx 0.14803096668067398473\cdots.
\end{align*}
and $\gamma$ is the Euler-Mascheroni constant and $A$ is the Glaisher-Kinkelin constant. So if we plug the Taylor expansion of $\frac{1}{e^x+1} = \frac{1}{2}(1-\tanh(x/2))$, we obtain
$$ \int_{1}^{\infty} \frac{dx}{x^2(e^x + 1)}
= C - \sum_{n=2}^{\infty} (-1)^n \frac{2n-1}{n-1} \cdot \frac{2^{2n}-1}{(2\pi)^{2n}} \zeta(2n). $$
Finally, expanding $\zeta(2n) = \sum_{k=1}^{\infty} \frac{1}{k^{2n}} $ and interchanging the order of summation, we find that
\begin{align*}
\sum_{n=2}^{\infty} (-1)^n \frac{2n-1}{n-1} \cdot \frac{2^{2n}-1}{(2\pi)^{2n}} \zeta(2n)
&= \frac{3-e}{4(1+e)} + \sum_{k=1}^{\infty} \frac{1}{(\pi k)^2} \log\left(1 + \frac{1}{(\pi k)^2} \right) \\
&\hspace{6em} - \sum_{k=1}^{\infty} \frac{1}{(2\pi k)^2} \log\left(1 + \frac{1}{(2\pi k)^2} \right)
\end{align*}
which might be better suited for studying this quantity.

Addendum. Here is a numerical computation:

A: We are dealing with
$$ \int_{1}^{+\infty}\frac{1-\tanh(x/2)}{2x^2}\,dx=\frac{1}{2}-\frac{1}{4}\int_{1/2}^{+\infty}\frac{\tanh(z)}{z^2}\,dz=\frac{1}{2}-\int_{1/2}^{+\infty}\sum_{n\geq 0}\frac{2\,dz}{z(4z^2+(2n+1)^2\pi^2)}$$
or
$$ \frac{1}{2}-\sum_{n\geq 0}\frac{\log(1+(2n+1)^2 \pi^2)}{(2n+1)^2 \pi^2}=\frac{1}{2}+\frac{\gamma}{4}+\frac{\log 2}{3}-3\log A-\sum_{n\geq 0}\frac{\log\left(1+\frac{1}{(2n+1)^2\pi^2}\right)}{(2n+1)^2 \pi^2} $$
where the last series can be numerically approximated by exploiting the Padé relation $\log(1+z^2)\approx\frac{z^2}{1+\frac{z^2}{2}}$ in a neighborhood of the origin. This leads to
$$ \int_{1}^{+\infty}\frac{dx}{x^2(e^x+1)}\,dx\approx \frac{1}{12} \left(3+3\gamma+4\log 2-36\log A+6 \sqrt{2}\tanh\frac{1}{2\sqrt{2}}\right)=\color{green}{0.11916}9\ldots$$
