In Mathematica 9, I want to evaluate the following integral:
$\qquad \int_{0}^{\infty}x^{\frac{-9}{8}}*e^{-16*x^{\frac{-1}{8}}-\frac{2*x}{0.5^{9}}}dx$
The result should be positive. However, Mathematica gives me -0.0824979.
This is the code from Mathematica in InputForm:
Integrate[E^(-16/x^8^(-1) - (2*x)/0.5^9)/x^(9/8), {x, 0, Infinity}]
I don't think I typed a wrong formula.
Could someone help me answering why did Mathematica give me a wrong answer?
This is a precision issue.
$Version
(* "12.0.0 for Mac OS X x86 (64-bit) (April 7, 2019)" *)
Since you used a machine precision number (i.e., 0.5) in the integrand, the integration is done with only machine precision.
Integrate[E^(-16/x^8^(-1) - (2*x)/(0.5)^9)/x^(9/8), {x, 0, Infinity}]
(* -0.0824979 *)
Using exact numbers
int = Integrate[E^(-16/x^8^(-1) - (2*x)/(1/2)^9)/x^(9/8), {x, 0, Infinity}]
(* (1/(2^(3/4) \[Pi]^(
7/2)))MeijerG[{{}, {}}, {{-(1/8), 0, 1/8, 1/4, 3/8, 1/2, 5/8, 3/4, 7/
8}, {}}, 262144] *)
If this complicated result is converted to a number using machine precision
int // N
(* 0.0270485 *)
However, using arbitrary-precision
int // N[#, 20] &
(* 5.6160550672798812373*10^-16 *)
Comparing with numerical integration
int = NIntegrate[E^(-16/x^8^(-1) - (2*x)/(1/2)^9)/x^(9/8), {x, 0, Infinity}]
(* 5.61606*10^-16 *)
This value is consistent with the earlier arbitrary-precision result and does not change with increased precision.
int = NIntegrate[E^(-16/x^8^(-1) - (2*x)/(1/2)^9)/x^(9/8), {x, 0, Infinity},
WorkingPrecision -> 20]
(* 5.6160550672798812757*10^-16 *)
Looking at the integrand
max = NMaximize[
{E^(-16/x^8^(-1) - (2*x)/(1/2)^9)/x^(9/8),
x > 0}, x, WorkingPrecision -> 20]
(* {1.3826054145707193900*10^-13, {x -> 0.0029477253331967126781}} *)
LogPlot[E^(-16/x^8^(-1) - (2*x)/(1/2)^9)/x^(9/8), {x, 1*^-4, 2*^-2},
WorkingPrecision -> 20,
Epilog -> {Red, AbsolutePointSize[4],
Point[{x /. max[[2]], Log@max[[1]]}]}]
Consequently, a small value for the integration is expected.
Another case of using non-exact numbers with exact function.
See the difference:
integrand = x^(-9/8) Exp[-16 x^(-1/8) - 2 x/(1/2)^9];
Integrate[integrand, {x, 0, Infinity}] // N
integrand = x^(-9/8) Exp[-16 x^(-1/8) - 2 x/(0.5)^9];
Integrate[integrand, {x, 0, Infinity}]
Rule of thumb: use exact numbers when calling exact functions like DSolve, Integrate, etc.. as sometimes when using non-exact numbers, Mathematica result can go wrong.
1/2 is exact and 0.5 is approximate. When Mathematica sees even one approximate number in an expression, it falls back to numerical approximation techniques. If all terms are exact, then exact symbolic techniques are used (where possible). In the case at hand, the presence of 0.5 term causes a machine precision approximation to be calculated which differs substantially from the exact symbolic result.
$\endgroup$