# $\int\limits_{-\infty}^\infty \left(f(\frac{x-\mu}{1+\Psi/2})-f(\frac{x+\mu}{1-\Psi/2})\right)\frac{x\gamma}{(x-x_{0})^{2}+a}dx$

EDIT: I realized from numerical implementation that the step from \begin{align} \mathcal{I}_2=&\frac{\gamma}{2}\int\limits_{-\infty}^\infty \left(f_T(\frac{x-\mu}{1+\Psi/2})-f_T(\frac{x+\mu}{1-\Psi/2})\right)\left[\frac{1}{x-\left(x_{0}+ \mathrm{i}\frac{\gamma}{2}\right)}+\frac{1}{x-\left(x_{0}- \mathrm{i}\frac{\gamma}{2}\right)}\right]dx\nonumber \end{align} is not the same to the 'equivalent expression' $$\begin{gather} \mathcal{I}_2=\mathcal{I}^A+\mathcal{I}^R\\ \mathcal{I}^{A/R}=\frac{\gamma}{2}\int\limits_{-\infty}^\infty f_T(z)\left[\frac{1}{z-\frac{x_0-\mu\pm\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}-\frac{1}{z-\frac{x_0+\mu\pm\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2}}\right]dz \end{gather}$$ (in fact the first expression gives the black line of the plot and the second gives the red line of the plot). Can someone see the difference?

ORIGINAL QUESTION:I want to solve the following integral

\begin{align} W=\int\limits_{-\infty}^\infty \left(f_T(\frac{x-\mu}{1+\Psi/2})-f_T(\frac{x+\mu}{1-\Psi/2})\right)\frac{x\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}}dx \end{align}

with $$f_{T}(x)=\frac{1}{1+e^{\frac{x}{T}}}$$ the Fermi function and $$x_{0},\mu,\gamma, T\in\Re_{\ge0}$$ and $$\Psi\in[-2,2]$$.

I have calculated the analytical solution for the integral in terms of the digammas $$\psi(z)$$ functions (see derivation below), but when comparing the derived expression with the numerical solution using the quadrature method from the scipy routine in python for the set of parameters $$x_{0}=\mu=0,\gamma= T=1$$ and $$\Psi\in[-2,2]$$ there is an evident discrepancy between both approaches, plot at the end. I suspect that my analytical derivation is not completely correct, although when I set $$\Psi=0$$ the agreement between both approaches is perfect. Can anybody help me to understand where I'm wrong?

The required integral is derived after an auxiliary integral $$\mathcal{I}_1$$ (correctly checked). Problems may start from $$\mathcal{I}_2$$.

\begin{align} \mathcal{I}_1=\int\limits_{-\infty}^\infty f_T(x-\mu)\frac{\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}}dx \end{align}

Using the substitution $$x = \mu + T u$$ and the abbreviations $$u_0 = \frac{x_0 - \mu}{T}$$ and $$\beta = \frac{\gamma}{2T}$$, we can write the integral as

\begin{align} \mathcal{I}_{1} = 2 \beta \int \limits_{-\infty}^\infty \frac{\mathrm{d} u}{(\mathrm{e}^{u} + 1)[(u-u_0)^2 + \beta^2]} \equiv 2 \beta \int \limits_{-\infty}^\infty g_{\beta,u_0}(u) \, \mathrm{d} u \, . \end{align}

$$g_{\beta,u_0}$$ has simple poles at $$u_0 \pm \mathrm{i} \beta$$ with residues

\begin{align} \operatorname{Res}(g_{\beta,u_0},u_0 \pm \mathrm{i} \beta) = \pm \frac{1}{2 \mathrm{i} \beta (\mathrm{e}^{u_0 \pm \mathrm{i} \beta} + 1)} \end{align}

and at $$\pm (2n+1) \pi \mathrm{i}$$ (so only at the odd-integer multiples of $$\pi \mathrm{i}$$ !) with residues

\begin{align} \operatorname{Res}(g_{\beta,u_0},\pm (2n+1) \pi \mathrm{i}) = - \frac{1}{[(2n+1) \pi \mathrm{i} \mp u_0]^2 + \beta^2} \end{align}

for $$n \in \mathbb{N}_0$$ . The integrals of $$g_{\beta,u_0}$$ along semi-circles (avoiding the poles on the imaginary axis) vanish in the limit of large radii, so we can use the residue theorem to evaluate the integral.

Closing the contour in the upper half-plane (the lower half-plane works just as well) yields

\begin{align} \mathcal{I}_{1} =& 2 \beta \, 2 \pi \mathrm{i} \left[\operatorname{Res}(g_{\beta,u_0},u_0 + \mathrm{i} \beta) + \sum \limits_{n=0}^\infty \operatorname{Res}(g_{\beta,u_0},(2n+1) \pi \mathrm{i}) \right] \\ =&4 \pi \beta\mathrm{i} \left[ \frac{1}{2 \mathrm{i} \beta (\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1)}- \sum \limits_{n=0}^\infty \frac{1}{[(2n+1) \pi \mathrm{i} - u_0]^2 + \beta^2} \right]\\ =&4 \pi \beta\mathrm{i} \left[ \frac{1}{2 \mathrm{i} \beta (\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1)}+ \sum \limits_{n=0}^\infty\frac{\frac{1}{4\pi\beta}}{\left(n+\left[\frac{1}{2}+\frac{\mathrm{i}u_0+\beta}{2\pi}\right]\right)}-\sum \limits_{n=0}^\infty\frac{\frac{1}{4\pi\beta}}{\left(n+\left[\frac{1}{2}+\frac{\mathrm{i}u_0-\beta}{2\pi}\right]\right)} \right] \end{align}

Now we use the series formula for the digamma function $$\psi$$

\begin{align} \psi(z)=\sum \limits_{n=0}^\infty\left(\frac{1}{n+1}-\frac{1}{n+z}\right)-\gamma \end{align}

to find

\begin{align} \mathcal{I} = 2 \pi \left[\frac{1}{\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1} - \frac{1}{2 \pi \mathrm{i}} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) - \psi \left(\frac{1}{2} + \frac{-\beta + \mathrm{i} u_0}{2 \pi}\right) \right)\right]\end{align}

Finally, we apply the reflection formula

\begin{align} \psi(1-z)=\psi(z)+\pi\cot{\left(\pi z\right)} \end{align}

to the second digamma function together with the properties

$$\begin{gather} \operatorname{Im} \left[\psi(x+ \mathrm{i}y)\right]=\frac{ \mathrm{i}}{2}\left(\psi(x- \mathrm{i}y)-\psi(x+ \mathrm{i}y)\right)\\ \cot{\left(\frac{\pi}{2}+z\right)}=-\tan(z)\\ \frac{1}{e^{z}+1}+\frac{\mathrm{i}}{2}\tan\left(\frac{-\mathrm{i}z}{2}\right)=\frac{1}{2} \end{gather}$$

and simplify the result: \begin{align} \mathcal{I}&= 2 \pi \left[\frac{1}{\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1} - \frac{1}{2 \mathrm{i}} \tan\left(\frac{\beta - \mathrm{i} u_0}{2}\right) -\frac{1}{2 \pi \mathrm{i}} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) - \psi \left(\frac{1}{2} + \frac{\beta - \mathrm{i} u_0}{2 \pi}\right) \right)\right] \\ &= 2 \pi \left[\frac{1}{2} - \frac{1}{\pi} \operatorname{Im} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) \right) \right] = \pi - 2 \operatorname{Im} \left[\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) \right] \, . \end{align}

Returning to the original parameters, we end up with

\begin{align} \boxed{\mathcal{I}_{1}(\gamma,x_{0}, \mu,T) =\int \limits_{-\infty}^\infty \frac{1}{\mathrm{e}^{(x - \mu)/T} + 1} \frac{\gamma}{(x-x_0)^2 + \gamma^2/4} \, \mathrm{d} x = \pi - 2 \operatorname{Im} \left[\psi \left(\frac{1}{2} + \frac{\frac{\gamma}{2} + \mathrm{i} (x_0 - \mu)}{2 \pi T}\right) \right] }\, . \end{align}

Note that the dependence on the parameters is explicit and in the special case $$x_0 = \mu$$ the result is simply $$\pi$$. This can also be shown using elementary methods, so the implicit assumption $$u_0 + \mathrm{i} \beta \not\in (2 \mathbb{Z} + 1) \pi \mathrm{i}$$ used in the computation of the residues is justified.

The integral we are interested in is

\begin{align} \mathcal{I}_2=\int\limits_{-\infty}^\infty \left(f_T(\frac{x-\mu}{1+\Psi/2})-f_T(\frac{x+\mu}{1-\Psi/2})\right)\frac{x\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}}dx \label{eq:2} \end{align} We can rewrite the integral $$I_{2}$$ decomposing the second factor as the sum of two complex parts. This is related to the advanced and retarded Green's functions

\begin{align} \frac{x\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}} = &\frac{\gamma}{2}\left[\frac{1}{x-\left(x_{0}+ \mathrm{i}\frac{\gamma}{2}\right)}+\frac{1}{x-\left(x_{0}- \mathrm{i}\frac{\gamma}{2}\right)}\right] - \mathrm{i}x_{0}\left[\frac{1}{x-\left(x_{0}+ \mathrm{i}\frac{\gamma}{2}\right)}-\frac{1}{x-\left(x_{0}- \mathrm{i}\frac{\gamma}{2}\right)}\right]\\ =&\frac{\gamma}{2}\left[G^{A}+G^{R}\right]- \mathrm{i}x_{0}\left[G^{A}-G^{R}\right] \end{align}

On the other hand, we can rewrite the Cauchy distribution used in the previous integral as function of the Green's functions too and using that the second contribution is related to the previously calculated integral $$\mathcal{I}_{1}$$

\begin{align} \frac{\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}} =- \mathrm{i}\left[\frac{1}{x-\left(x_{0}+ \mathrm{i}\frac{\gamma}{2}\right)}-\frac{1}{x-\left(x_{0}- \mathrm{i}\frac{\gamma}{2}\right)}\right] =- \mathrm{i}\left[G^{A}-G^{R}\right]. \end{align}

So

\begin{align} \mathcal{I}_2=&\frac{\gamma}{2}\int\limits_{-\infty}^\infty \left(f_T(\frac{x-\mu}{1+\Psi/2})-f_T(\frac{x+\mu}{1-\Psi/2})\right)\left[\frac{1}{x-\left(x_{0}+ \mathrm{i}\frac{\gamma}{2}\right)}+\frac{1}{x-\left(x_{0}- \mathrm{i}\frac{\gamma}{2}\right)}\right]dx\nonumber\\&+x_{0}\left(I_{1}\left(\gamma,x_{0}, \mu,T(1+\Psi/2)\right)-I_{1}(\gamma,x_{0}, -\mu,T(1-\Psi/2))\right)\\ =&\frac{\gamma}{2}\int\limits_{-\infty}^\infty f_T(z)\left[\left(\frac{1}{z-\frac{x_0-\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}+\frac{1}{z-\frac{x_0-\mu-\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}\right)-\left(\frac{1}{z-\frac{x_0+\mu+\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2}}+\frac{1}{z-\frac{x_0+\mu-\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2}}\right)\right]dz\nonumber\\&+x_{0}\left(\mathcal{I}_{1}\left(\gamma,x_{0}, \mu,T_L\right)-\mathcal{I}_{1}(\gamma,x_{0}, -\mu,T_R)\right) \end{align}

Note that the integral $$\mathcal{I}_2$$ is convergent because the difference of the Fermi functions contributes asymptotically as $$x^{-1}$$ and the Green's function contributes another $$x^{-1}$$. The integral of a single Green's function with a single Fermi function would be logarithmically divergent. In the second equality the Fermi functions arguments have been shifted to the Green's functions terms we can split the integral in two parts such that

$$\begin{gather} \mathcal{I}_2=\mathcal{I}^A+\mathcal{I}^R+x_{0}\left(\mathcal{I}_{1}\left(\gamma,x_{0}, \mu,T_L\right)-\mathcal{I}_{1}(\gamma,x_{0}, -\mu,T_R)\right)\\ \mathcal{I}^{A/R}=\frac{\gamma}{2}\int\limits_{-\infty}^\infty f_T(z)\left[\frac{1}{z-\frac{x_0-\mu\pm\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}-\frac{1}{z-\frac{x_0+\mu\pm\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2}}\right]dz=\frac{\gamma}{2}\int\limits_{-\infty}^\infty g^{A/R}(z)dz \end{gather}$$

closing the integration contour in the upper half and using now the series relation for the Fermi function

\begin{align} f_{T}(z)=\frac{1}{1+e^{\frac{z}{T}}} =\frac{1}{2}- \frac{\mathrm{i}}{2\pi}\sum \limits_{n=-\infty}^\infty \frac{1}{n+\frac{1}{2} +\mathrm{i}\frac{z}{2\pi T}} \end{align}

we have

$$\begin{gather} \mathcal{I}^{A}=\frac{\gamma}{2}2\pi\mathrm{i}\left( \operatorname{Res}(g^{A},\frac{x_0-\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2})-\operatorname{Res}(g^{A},\frac{x_0+\mu+\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2})+\sum\limits_{n=0}^\infty \operatorname{Res}(g^{A},(2n+1)\pi\mathrm{i}T)\right)\\ =\gamma\pi\mathrm{i}\left(f_{T}(\frac{x_0-\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2})-f_{T}(\frac{x_0+\mu+\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2})-T\sum\limits_{n=0}^\infty\frac{1}{(2n+1)i\pi T-\frac{x_{0}-\mu+i\gamma/2}{1+\Psi/2}}+T\sum\limits_{n=0}^\infty\frac{1}{(2n+1)i\pi T-\frac{x_{0}+\mu+i\gamma/2}{1-\Psi/2}}\right)\\ =\gamma\pi\mathrm{i}\left(f_{T}(\frac{x_0-\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2})-f_{T}(\frac{x_0+\mu+\mathrm{i}\frac{\gamma}{2}}{1-\Psi/2})+\frac{i}{2\pi}\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}-\mu)}{2\pi T(1+\Psi/2)}}-\frac{i}{2\pi}\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}+\mu)}{2\pi T(1-\Psi/2)}}\right)\\ =\gamma\pi\mathrm{i}\frac{i}{2\pi}\left(-\sum \limits_{n=-\infty}^\infty \frac{1}{n+\frac{1}{2} +\mathrm{i}\frac{\frac{x_0-\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}{2\pi T}}+\sum \limits_{n=-\infty}^\infty \frac{1}{n+\frac{1}{2} +\mathrm{i}\frac{\frac{x_0+\mu+\mathrm{i}\frac{\gamma}{2}}{1+\Psi/2}}{2\pi T}}+\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}-\mu)}{2\pi T(1+\Psi/2)}}-\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}+\mu)}{2\pi T(1-\Psi/2)}}\right)\\ =\frac{\gamma}{2\pi}\left(\sum \limits_{n=-\infty}^\infty \frac{1}{n+\frac{1}{2} +\frac{-\gamma/2+i(x_{0}-\mu)}{2\pi T_L}}-\sum \limits_{n=-\infty}^\infty \frac{1}{n+\frac{1}{2} +\frac{-\gamma/2+i(x_{0}+\mu)}{2\pi T_R}}-\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}-\mu)}{2\pi T_L}}+\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{-\gamma/2+i(x_{0}+\mu)}{2\pi T_R}}\right)\\ =\frac{\gamma}{2}\left(\sum \limits_{n=-\infty}^{-1} \frac{1}{n+\frac{1}{2} +\frac{-\gamma/2+i(x_{0}-\mu)}{2\pi T_L}}-\sum \limits_{n=-\infty}^{-1} \frac{1}{n+\frac{1}{2} +\frac{-\gamma/2+i(x_{0}+\mu)}{2\pi T_R}}\right)\\ =\frac{\gamma}{2}\left(\sum \limits_{n=1}^{\infty} \frac{-1}{n-\frac{1}{2} +\frac{\gamma/2-i(x_{0}-\mu)}{2\pi T_L}}+\sum \limits_{n=1}^{\infty} \frac{1}{n-\frac{1}{2} +\frac{\gamma/2-i(x_{0}+\mu)}{2\pi T_R}}\right)\\ =\frac{\gamma}{2}\left(\psi(\frac{1}{2} +\frac{\gamma/2-i(x_{0}-\mu)}{2\pi T_L})-\psi(\frac{1}{2} +\frac{\gamma/2-i(x_{0}+\mu)}{2\pi T_R})\right) \end{gather}$$

with $$T_{L/R}=T(1\pm\Psi/2)$$. In the last step we have used

\begin{align} \psi(z+1)=-\gamma+\sum\limits_{n=1}^\infty\left(\frac{1}{n}-\frac{1}{n+z}\right) \end{align}

for the second contribution

$$\begin{gather} \mathcal{I}^{R}=\frac{\gamma}{2}2\pi\mathrm{i}\left(\sum\limits_{n=0}^\infty \operatorname{Res}(g^{R},(2n+1)\pi\mathrm{i}T)\right)\\ =\gamma\pi\mathrm{i}\left(\frac{i}{2\pi}\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T(1+\Psi/2)}}-\frac{i}{2\pi}\sum\limits_{n=0}^\infty\frac{1}{n+\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T)}}\right)\\ =\frac{\gamma}{2}\left(\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T_L})-\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T_R})\right) \end{gather}$$

So

$$\begin{gather} \mathcal{I}^{A}+\mathcal{I}^{R}=\frac{\gamma}{2}\left(\psi(\frac{1}{2}+\frac{\gamma/2-i(x_{0}-\mu)}{2\pi T_L})-\psi(\frac{1}{2}+\frac{\gamma/2-i(x_{0}+\mu)}{2\pi T_R})+\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T_L})-\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T_R})\right)\\ =\gamma\operatorname{Re} \left(\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T_L})\right)-\gamma\operatorname{Re} \left( \psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T_R}\right) \end{gather}$$

Finally

$$\begin{gather} \mathcal{I}_2=\gamma\operatorname{Re} \left(\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T_L})\right)-\gamma\operatorname{Re} \left( \psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T_{R}})\right)+x_{0}\left(\mathcal{I}_{1}\left(\gamma,x_{0}, \mu,T_L\right)-\mathcal{I}_{1}(\gamma,x_{0}, -\mu,T_R)\right) \end{gather}$$ So

$$\begin{gather} \boxed{W=\gamma\operatorname{Re} \left(\psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}-\mu)}{2\pi T_L})\right)-\gamma\operatorname{Re} \left( \psi(\frac{1}{2}+\frac{\gamma/2+i(x_{0}+\mu)}{2\pi T_{R}})\right)+x_{0}\left(\mathcal{I}_{1}\left(\gamma,x_{0}, \mu,T_L\right)-\mathcal{I}_{1}(\gamma,x_{0}, -\mu,T_R)\right)} \end{gather}$$

(The motivation is related to a quantum problem, but people downvoted previous question arguing that this is related to physics, so if somebody wants me to explain the physical interpretation of the parameters please do not hesitate.)

• To debug this, how about evaluating some of the intermediate steps numerically as well - that might help to pin down where the problem is. May 6, 2021 at 2:49
• A tricky point, not least in terms of numerics, is how you let the integration regions go to infinity. Although, as you say, your combination of terms converges, this comes about through a cancellation between divergent terms, and so the convergence is conditional. If you swap things around so that the approach to infinite integration range is changed between the terms relative to each other, you may change the result by a finite amount. You'd have to carefully say how exactly the divergent terms are combined to form a convergent result, and take care to stick with that in your manipulations. May 6, 2021 at 18:14
• So, concretely, I suspect the problem happens when you rescale the integration variable by $1+\Psi /2$ and $1-\Psi /2$ respectively in two divergent terms that need to cancel. I'd focus on that step. May 6, 2021 at 18:34
• To boil it down to essentials, the following type of thing seems to be happening. Consider $f(x)=\{ 1$ for $x<1$ and $1/x$ for $x>1 \}$, and the reasoning (substituting $x\rightarrow 2x$ in one of the integrals) $0 = \int_{0}^{\infty } dx f(x) - \int_{0}^{\infty } dx f(x) = \int_{0}^{\infty } dx f(x) - \int_{0}^{\infty } dx 2f(2x) = \int_{0}^{1/2} dx (-1) + \int_{1/2}^{1} dx (1-1/x) = -\ln 2 \neq 0$ May 6, 2021 at 20:40
• I could imagine the difference between your red and black plots is something rather simple, like a linear combination of $\ln (1+\Psi /2)$ and $\ln (1-\Psi /2)$. May 6, 2021 at 21:20

Solution by @ComplexYetTrivial in math.stackexchange: The step in the edit is indeed the only mistake, the rest of your calculation seems to be correct. In order to see what the problem is we write $$\mathcal{J}^{A/R} = \frac{\gamma}{2} \int \limits_{-\infty}^{\infty} \left[f_T \left(\frac{x - \mu}{1 + \Psi/2}\right) - f_T \left(\frac{x + \mu}{1 - \Psi/2}\right)\right] \, \frac{\mathrm{d} x}{x - x_0 \mp \mathrm{i} \frac{\gamma}{2}} \, .$$ You correctly noted that these integrals are convergent due to the difference of the two Fermi functions. If there was only one Fermi function, they would diverge logarithmically at the lower limit.
However, you then proceed to perform two different changes of variables ($$x = \pm \mu + (1 \pm \Psi/2) z$$) for the first and the second term. This amounts to splitting the convergent integral into two divergent integrals, doing the substitutions and then combining them again, which is not allowed and changes the result!
I will show the second possibility here, since it illustrates nicely how the error arises. We have \begin{align} \mathcal{J}^{A/R} &= \lim_{r \to \infty} \frac{\gamma}{2} \int \limits_{-r}^{\infty} \left[f_T \left(\frac{x - \mu}{1 + \Psi/2}\right) - f_T \left(\frac{x + \mu}{1 - \Psi/2}\right)\right] \, \frac{\mathrm{d} x}{x - x_0 \mp \mathrm{i} \frac{\gamma}{2}} \\ &= \lim_{r \to \infty} \frac{\gamma}{2} \left[~\int \limits_{-r}^{\infty} \, \frac{f_T \left(\frac{x - \mu}{1 + \Psi/2}\right)}{x - x_0 \mp \mathrm{i} \frac{\gamma}{2}} \, \mathrm{d} x - \int \limits_{-r}^{\infty} \, \frac{f_T \left(\frac{x + \mu}{1 - \Psi/2}\right)}{x - x_0 \mp \mathrm{i} \frac{\gamma}{2}} \, \mathrm{d} x \right] \\ &= \lim_{r \to \infty} \frac{\gamma}{2} \left[~\int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{\infty} \, \frac{f_T \left(z\right)}{z - \frac{x_0 - \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 + \Psi/2}} \, \mathrm{d} z - \int \limits_{- \frac{r - \mu}{1 - \Psi/2}}^{\infty} \, \frac{f_T \left(z\right)}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} \, \mathrm{d} z\right] \\ &= \lim_{r \to \infty} \frac{\gamma}{2} \left[~\int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{\infty} f_T (z) \left(\frac{1}{z - \frac{x_0 - \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 + \Psi/2}} - \frac{1}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}}\right) \, \mathrm{d} z + \int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{- \frac{r - \mu}{1 - \Psi/2}} \frac{f_T \left(z\right)}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} \, \mathrm{d} z \right]. \end{align} Letting $$r \to \infty$$ in the first integral we simply obtain your result ($$\mathcal{I}^{A/R}$$). The second term is the correction we are looking for: \begin{align} \mathcal{J}^{A/R} - \mathcal{I}^{A/R} &= \lim_{r \to \infty} \frac{\gamma}{2} \int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{- \frac{r - \mu}{1 - \Psi/2}} \frac{f_T \left(z\right)}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} \, \mathrm{d} z \\ &= \lim_{r \to \infty} \frac{\gamma}{2} \int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{- \frac{r - \mu}{1 - \Psi/2}} \frac{\mathrm{d} z}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} - \lim_{r \to \infty} \frac{\gamma}{2} \int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{- \frac{r - \mu}{1 - \Psi/2}} \frac{1 - f_T \left(z\right)}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} \, \mathrm{d} z \, . \end{align} Since $$1 - f_T(z) \leq \mathrm{e}^{z/T}$$, the second limit vanishes and we are left with \begin{align} \mathcal{J}^{A/R} - \mathcal{I}^{A/R} &= \lim_{r \to \infty} \frac{\gamma}{2} \int \limits_{- \frac{r + \mu}{1 + \Psi/2}}^{- \frac{r - \mu}{1 - \Psi/2}} \frac{\mathrm{d} z}{z - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}} = \lim_{r \to \infty} \frac{\gamma}{2} \log \left(\frac{- \frac{r - \mu}{1 - \Psi/2} - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}}{- \frac{r + \mu}{1 + \Psi/2} - \frac{x_0 + \mu \pm \mathrm{i} \frac{\gamma}{2}}{1 - \Psi/2}}\right) \\ &= \frac{\gamma}{2} \log \left(\frac{1 + \Psi/2}{1-\Psi/2}\right) = \gamma \operatorname{artanh} \left(\frac{\Psi}{2}\right) \, . \end{align}
With this addition the correct final result is \begin{align} W &= \int\limits_{-\infty}^\infty \left[f_T\left(\frac{x-\mu}{1+\Psi/2}\right)-f_T\left(\frac{x+\mu}{1-\Psi/2}\right)\right] \frac{\gamma x}{(x-x_0)^2 + \gamma^2/4} \, \mathrm{d} x \\ &=2 \gamma \operatorname{artanh} \left(\frac{\Psi}{2}\right) \!+ \!(\gamma \operatorname{Re} {} \! -\! 2 x_0 \operatorname{Im}) \left[\operatorname{\psi} \left(\frac{1}{2} \! + \! \frac{\gamma/2 \! + \! \mathrm{i} (x_0 \! - \! \mu)}{2 \pi T (1 \! + \! \Psi/2)}\right) \! - \! \operatorname{\psi} \left(\frac{1}{2} \! + \! \frac{\gamma/2 \! + \! \mathrm{i} (x_0 \! + \! \mu)}{2 \pi T (1 \! - \! \Psi/2)} \right) \right] . \end{align}