7
$\begingroup$

I have integrals of the form $$\int_0^\infty\mathrm d x\; e^{-x^2}\;\;(\textrm{Polynomial of degree n})\times(\textrm{Polynomial of degree m}),$$ where $n,m$ can be as large as 300.

When I use Integrate[f[x], {x, 0, Infinity}], the calculation is too slow for large polynomial degrees.

Is there a way to speed up this calculation? (I tried to use NIntegrate, but unfortunately, since the integrand is highly oscillatory at large $n,m$, I'm unable to get reliable results from numerical methods).

$\endgroup$
10
$\begingroup$

Perhaps this?:

Table[1/2 Gamma[(1 + n)/2], {n, 0, Length@# - 1}] . # &@
 CoefficientList[poly1*poly2, x]

Low-order check:

poly1 = LaguerreL[4, x];
poly2 = LaguerreL[5, x];

Table[1/2 Gamma[(1 + n)/2], {n, 0, Length@# - 1}] . # &@
 CoefficientList[poly1*poly2, x]

(*  -(7981/240) + (38645 Sqrt[π])/2048  *)

Integrate[Exp[-x^2] poly1*poly2, {x, 0, Infinity}]

(*  -(7981/240) + (38645 Sqrt[π])/2048  *)

High-order timing:

poly1 = LaguerreL[300, x];
poly2 = LaguerreL[301, x];

Table[1/2 Gamma[(1 + n)/2], {n, 0, Length@# - 1}] . # &@
   CoefficientList[poly1*poly2, x]; // AbsoluteTiming

(*  {0.957942, Null}  *)
$\endgroup$
5
  • $\begingroup$ Thank you for your quick answer. Can you add some description of what the code actually does? I'm not very experienced in Mathematica, so I don't actually understand your code :D $\endgroup$ – temporal99x Mar 4 at 3:20
  • 2
    $\begingroup$ @temporal99x $\int_0^\infty e^{-x^2} x^n \, dx = {1\over2} \Gamma\left[{n+1 \over 2}\right]$, so $\int_0^\infty e^{-x^2} \sum_n a_n x^n \, dx = \sum_n a_n\,{1\over2} \Gamma\left[{n+1 \over 2}\right]$. CoefficientList gives us the $a_n$ and the . makes the Dot product with the list/table of gamma function values. -- It could be broken down into three steps: an = CoefficientList[poly1*poly2, x]; gammafactors = Table[1/2 Gamma[(1 + n)/2], {n, 0, Length@an - 1}]; an . gammafactors $\endgroup$ – Michael E2 Mar 4 at 3:23
  • $\begingroup$ @temporal99x A couple of other ways to speed up Integrate, but they're not as fast as the method above. Both assume the integral is convergent: (1) Integrate[Exp[-x^2] poly1*poly2, {x, 0, Infinity}, GenerateConditions -> False] (2) Integrate[Exp[-x^2] poly1*poly2, x] /. {{x -> 0}, {Exp[-x^2] -> 0, Erf[x] -> 1}} // Differences // First $\endgroup$ – Michael E2 Mar 4 at 3:45
  • $\begingroup$ I now fully understand your solution. It is just beautiful, and fast! You may not know how helpful this is to me. $\endgroup$ – temporal99x Mar 4 at 3:57
  • 1
    $\begingroup$ @temporal99x I'm glad to be of help. Here is a slightly faster way to calculate the coefficients of the product of two polynomials: an = ListConvolve[CoefficientList[poly1, x], CoefficientList[poly2, x], {1, -1}, 0]; $\endgroup$ – Michael E2 Mar 4 at 5:17

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.