Suppose that we have the integral $$ I(a) = \int_0^\infty dx \frac{e^{x^2-a}(x^2-a)^2}{(e^{x^2-a}+1)^2} $$ where $a$ is a real value. It can be positive or negative.
I can deal with the negative case. In that case, the integral may be $$ I(a)=\frac{\sqrt{\pi}}{2} \left[\frac34f_{3/2}(\lambda)-(\ln\lambda) f_{1/2}(\lambda) + (\ln\lambda)^2 f_{-1/2}(\lambda)\right] $$ where $\lambda=e^a$, and $f_n(\lambda)=\sum_{k=1}^\infty (-1)^{k+1}\frac{\lambda^k}{k^n}$. Here, $f_n(\lambda)$ could be regarded as a special function since it relates to the Dirichlet eta function.
My question is, how to solve the positive $a$ case?
For $\,a>0\,$ the result is the $\underline{\text{same}}$!
The calculation starts with:
$\displaystyle I(a) = \int_0^\sqrt{a} \frac{e^{x^2-a}(x^2-a)^2}{(e^{x^2-a}+1)^2} + \int_\sqrt{a} ^\infty \frac{e^{x^2-a}(x^2-a)^2}{(e^{x^2-a}+1)^2}$
As a simplification I've used $\enspace\displaystyle \frac{d}{dz}\frac{x^2-a}{e^{z(x^2-a)}+1}=-\frac{e^{z(x^2-a)}(x^2-a)^2}{(e^{z(x^2-a)}+1)^2}\enspace$ but that's not really necessary.
The result is:
You only have to use the analytical continuation of $\,f_n(e^a)\equiv -Li_n(-e^a)\,$ , that’s all.
See polylogarithm e.g. here (Integral representations, Series representations).
E.g.: $\enspace\displaystyle f_n(e^a) = e^a \left(\frac{1}{2} + \int\limits_0^\infty \frac{\sin(n \arctan t - at)}{(1+t^2)^{n/2} \sinh(\pi t)} dt\right)\enspace$ for $\,a\in\mathbb{R}\,,\,n\in\mathbb{C}$
An example.
$\displaystyle I(0.3)= \int_0^\infty \frac{e^{x^2-0.3}(x^2-0.3)^2}{(e^{x^2-0.3}+1)^2} \approx 0.481255$
$\displaystyle \frac{\sqrt{\pi}}{2}\left( -\frac{3}{4}Li_{3/2}(-e^{0.3})+0.3 Li_{1/2}(-e^{0.3})- 0.3^2 Li_{-1/2}(-e^{0.3}) \right) \approx 0.481255$