10
$\begingroup$

Bug introduced in 10.0 and fixed in 10.0.2


Trying to do the integral

Integrate[
E^(-((201 x1^2)/101))
x1^2 (15 - 20 x1^2 + 4 x1^4) (5913508078417951503 + 
1124782662003060300000 x1^2 - 2983951574394000000000 x1^4 + 
2591462040000000000000 x1^6 - 797900000000000000000 x1^8 + 
80000000000000000000 x1^10), {x1, -Infinity, Infinity}]

Mathmatica tells me that it does not converge. However it certainly does due to the Gaussian term upfront. It seems to be an issue with the large oscillations due to the polynomials and I was wondering what options I need to supply so that Mathematica can evaluate the integral. (I will have to do this with hundreds of similar polynomials, so I am looking for a general solution to this problem)

Thanks so much

EDIT: Based on the answers and comments I contacted Wolfram and they acknowledged the issue. Here is their response:

Thank you for your email.

I agree that Mathematica should return a result in this case and have filed a report with our developers on this issue. This may allow them to fix the problem in a future version of Mathematica. I believe that you have received several workarounds to this problem in the StackExchange thread that you cite, so I will not repeat them here.

Please let me know if you have any further questions.

$\endgroup$
  • $\begingroup$ The minimal example: Integrate[Exp[-x^2] x^8 (1 + x^2 + x^4 + x^8), {x, -∞, ∞}]. I think it's a bug. Could anyone confirm it? $\endgroup$ – ybeltukov Oct 11 '14 at 17:26
  • 1
    $\begingroup$ @rhermans: Gaussians shrink faster than any polynomial, so the integral converges. $\endgroup$ – DumpsterDoofus Oct 11 '14 at 17:27
  • $\begingroup$ @ybeltukov I am using the 64 bit Linux Version of Mathematica 10. I will double check the minor version tomorrow or on Monday. You think this should be reported to wolfram? $\endgroup$ – ftiaronsem Oct 11 '14 at 19:27
  • $\begingroup$ Yes, I think it is good idea. $\endgroup$ – ybeltukov Oct 11 '14 at 19:30
  • 4
    $\begingroup$ Definitely a bug. $\endgroup$ – Daniel Lichtblau Oct 12 '14 at 19:46
5
$\begingroup$

Expectation, interestingly, does not have the same problem as Integrate, and it is much faster. If you "have to do this with hundreds of similar polynomials," then I recommend using Expectation.

First the exponential factor of the integrand corresponds to the pdf of a normal distribution, scaled by a constant:

1 / Sqrt[(201/(101 π))] * PDF[NormalDistribution[0, Sqrt[101/201/2]]][x1]

(* E^(-((201 x1^2)/101)) *)

The Expectation returns quite quickly, about 160 times faster than either Integrate method in the present answers, suggesting that the expectation of polynomials over normal distributions is a special case:

1/Sqrt[(201/(101 π))] * Expectation[ 
   x1^2 (15 - 20 x1^2 + 4 x1^4) (5913508078417951503 + 
      1124782662003060300000 x1^2 - 2983951574394000000000 x1^4 + 
      2591462040000000000000 x1^6 - 797900000000000000000 x1^8 + 
      80000000000000000000 x1^10), 
   x1 \[Distributed] NormalDistribution[0, Sqrt[101/201/2]]]

(* (368559503694222488015947879582225000 Sqrt[(101 π)/201])/32891481882087921 *)

One can check that this is equal to the result obtained by Fred Simons.

Update -- Further explanation

Examining the Trace,

integrand = 
  E^(-((201 x1^2)/101)) x1^2 (15 - 20 x1^2 + 
     4 x1^4) (5913508078417951503 + 1124782662003060300000 x1^2 - 
     2983951574394000000000 x1^4 + 2591462040000000000000 x1^6 - 
     797900000000000000000 x1^8 + 80000000000000000000 x1^10);
trace = With[{f = integrand /. E^_ -> 1},
   Trace[Expectation[f, 
     x1 \[Distributed] NormalDistribution[0, Sqrt[101/201/2]]], 
    TraceInternal -> True]
   ];

one observes that the integrand is detected to be a polynomial and that for each power x1^n, the coefficient of x1^n is multiplied by Moment[NormalDistribution[0, Sqrt[101/402]], n].

Timing

Needs["GeneralUtilities`"]

Integrate[integrand, {x1, -Infinity, Infinity}, PrincipalValue -> True] // AccurateTiming

(* 0.368640 *)

Integrate[#, {x1, -Infinity, Infinity}] & /@ Expand[integrand] // AccurateTiming

(* 0.306265 *)

Integrate[integrand, {x1, -Infinity, Infinity}, GenerateConditions -> False] // AccurateTiming

(* 0.216459 *)

1/Sqrt[(201/(101 \[Pi]))] *
  Expectation[integrand /. E^_ -> 1, 
   x1 \[Distributed] NormalDistribution[0, Sqrt[101/201/2]]] // AccurateTiming

(* 0.00135676 *)
$\endgroup$
  • $\begingroup$ Expectation might be calling Integrate, but with the setting GenerateConditions->False. Maybe. $\endgroup$ – Daniel Lichtblau Oct 14 '14 at 21:14
  • $\begingroup$ Thanks, Expectation is indeed very useful and fast. I didn't know about it before. $\endgroup$ – ftiaronsem Oct 20 '14 at 2:58
  • $\begingroup$ @ftiaronsem You're welcome. :) $\endgroup$ – Michael E2 Oct 20 '14 at 3:01
11
$\begingroup$

I have no idea why Mathematica complains about the convergence. When you expand the integrand, then integrate each of the terms and add, you will find the exact result:

E^(-((201 x1^2)/101)) x1^2 (15 - 20 x1^2 + 
4 x1^4) (5913508078417951503 + 1124782662003060300000 x1^2 - 
2983951574394000000000 x1^4 + 2591462040000000000000 x1^6 - 
797900000000000000000 x1^8 + 80000000000000000000 x1^10) // Expand;
Integrate[#, {x1, -\[Infinity], \[Infinity]}] & /@ %


(*-((209432594936834177932282080417913954555 Sqrt[(101 \[Pi])/201])/65782963764175842)+(9561258132402170537415665 Sqrt[(303 \[Pi])/67])/8978*)

You could also use NIntegrate, with the option WorkingPrecision->20 to force Mathematica to use arbitrary precision numbers, and the obtain numerically the same result as above.

$\endgroup$
  • $\begingroup$ Nice workaround! +1 $\endgroup$ – ybeltukov Oct 11 '14 at 17:24
4
$\begingroup$
Integrate[E^(-((201 x1^2)/101))
x1^2 (15 - 20 x1^2 + 4 x1^4) (5913508078417951503 + 
1124782662003060300000 x1^2 - 2983951574394000000000 x1^4 + 
2591462040000000000000 x1^6 - 797900000000000000000 x1^8 + 
80000000000000000000 x1^10), {x1, -Infinity, Infinity}, PrincipalValue -> True]

(*(368559503694222488015947879582225000 Sqrt[(
 101 \[Pi])/201])/32891481882087921*)
$\endgroup$
  • $\begingroup$ Most interesting. This is not even a principal value integral. There are no poles in the complex plane. But it works. Marvelous $\endgroup$ – ftiaronsem Oct 11 '14 at 21:04
  • $\begingroup$ @ftiaronsem, in my experience PrincipalValue can force an integral "convergent" to Mathematica (because, like you pointed out, sometimes MMA thinks a convergent integral diverges), which is also used internally in FourierTransform I guess. However the result you get may not always be correct, which is sometimes annoying... $\endgroup$ – Leo Fang Oct 14 '14 at 1:31
2
$\begingroup$

Mathmatica tells me that it does not converge

This has been fixed in V 10.0.2. On windows 7, 64 bit.

Mathematica graphics

Integrate[E^(-((201 x1^2)/101)) x1^2 (15 - 20 x1^2 + 4 x1^4) 
     (5913508078417951503 + 1124782662003060300000 x1^2 - 
    2983951574394000000000 x1^4 + 2591462040000000000000 x1^6 - 
    797900000000000000000 x1^8 + 
    80000000000000000000 x1^10), {x1, -Infinity, Infinity}]

Mathematica graphics

 Integrate[Exp[-x^2] x^8 (1 + x^2 + x^4 + x^8), {x, -\[Infinity], \[Infinity]}]

Mathematica graphics

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.