# Wrong sign with NIntegrate

I'm using Mathematica to integrate this function:

w[x_, z_] :=
E^x/(E^x + 1)^2 Log[(E^(z^2/(4 x)) + E^-x)/(E^(z^2/(4 x)) - 1)]
W[z_] := NIntegrate[w[x, z], {x, 0, \[Infinity]}]


where z > 0 is a positive parameter. Analytically, I know that the integrated is non-negative, so the integral itself should be non-negative. However, I'm getting some incorrect results due to a lack of accuracy:

W[100]


-6.05023*10^-20


My question: How can I get a reliable result? Thanks!

w[x_, z_] := E^x/(E^x + 1)^2 Log[(E^(z^2/(4 x)) + E^-x)/(E^(z^2/(4 x)) - 1)]


For larger z the integrand is very small and high precision is required in the integration.

Plot3D[w[x, z], {x, 0, 40}, {z, 0, 30}, WorkingPrecision -> 15,
PlotPoints -> 50,
MaxRecursion -> 5,
ClippingStyle -> None]


W2[z_?NumericQ] := NIntegrate[w[x, z], {x, 0, ∞},
WorkingPrecision -> 40,

W2[100]

(* 4.679853458969239635780655689865016458810*10^-43 *)

• Thank you so much! Nov 11, 2021 at 22:33

Simplifying and rearranging the integrand, along with high precision, solves the problem in this case:

With[{z = 100},
NIntegrate[1/4 Sech[x/2]^2 Log[(1 + E^-x E^(-(z^2/(4 x))))/(1 - E^(-(z^2/(4 x))))],
{x, 0, ∞}, WorkingPrecision -> 25]]
4.679853735636909286562544*10^-43


In particular, the important part here was to try to express the integrand in terms of functions that decay at infinity, like $$\exp(-x)$$ or $$\operatorname{sech} x$$.

• Point taken, thanks a lot. May 21, 2022 at 22:11