2
$\begingroup$

I have the following expression to integrate

$$ Overlap=\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}dt_{1}dt_{2} \left(\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{1}-t_{01})^2}}{\sqrt[4]{\pi }}\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{2}-t_{02})^2}}{\sqrt[4]{\pi }}+\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{2}-t_{01})^2}}{\sqrt[4]{\pi }}\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{1}-t_{02})^2}}{\sqrt[4]{\pi }}\right)^2 $$

This is just two Gaussian distributions and I think there is no problem for Mathematica, I have known that when $t_{01}=t_{02}$ the integral result is 4 and when $t_{02}-t_{01}\rightarrow\infty$ the integral approaches to 2. If we set $t_{01}=0$ and changes $t_{02}$ the ideal integral results would be as follows

enter image description here

However, If I directly integrate the expression using NIntegrate

f\[Tau]gaussian[\[CapitalGamma]_,t_,t0_:5]:=(E^(-(1/2) (t-t0)^2 \[CapitalGamma]^2) Sqrt[\[CapitalGamma]])/\[Pi]^(1/4);

ListPlot[Thread[{Table[t02, {t02, 0, 20, 1}], 
   Table[NIntegrate[(f\[Tau]gaussian[1, t1, 0]*
         f\[Tau]gaussian[1, t2, t02] + 
        f\[Tau]gaussian[1, t2, 0]*
         f\[Tau]gaussian[1, t1, t02])^2, {t1, -Infinity, Infinity}
        , {t2, -Infinity, Infinity}], {t02, 0, 20, 1}]}], 
 PlotRange -> Full, Joined -> True, Frame -> True,
 FrameLabel -> {"t02-t01", "Overlap"}]

enter image description here

when $t_{02}-t_{01}>9$ the results turn to be incorrect. The expression to be integrated could expand into 4 terms and If I integrate each of them independently the results would be much better

nor1 = Table[
   NIntegrate[(f\[Tau]gaussian[1, t1, 0]*
       f\[Tau]gaussian[1, t2, t02])^2, {t1, -Infinity, Infinity}
       , {t2, -Infinity, Infinity}], {t02, 0, 20, 0.5}];
nor2 = Table[
   NIntegrate[(f\[Tau]gaussian[1, t2, 0]*
       f\[Tau]gaussian[1, t1, t02])^2, {t1, -Infinity, Infinity}
       , {t2, -Infinity, Infinity}], {t02, 0, 20, 0.5}];
nor3 = Table[
   NIntegrate[(f\[Tau]gaussian[1, t1, 0]*f\[Tau]gaussian[1, t2, t02]*
      f\[Tau]gaussian[1, t2, 0]*
      f\[Tau]gaussian[1, t1, t02]), {t1, -Infinity, Infinity}
       , {t2, -Infinity, Infinity}], {t02, 0, 20, 0.5}];
nor4 = Table[
   NIntegrate[(f\[Tau]gaussian[1, t2, 0]*f\[Tau]gaussian[1, t1, t02]*
      f\[Tau]gaussian[1, t1, 0]*
      f\[Tau]gaussian[1, t2, t02]), {t1, -Infinity, Infinity}
       , {t2, -Infinity, Infinity}], {t02, 0, 20, 0.5}];
t0mat = Table[t02, {t02, 0, 20, 0.5}];
ListPlot[{Thread[{t0mat, nor1}], Thread[{t0mat, nor2}], 
  Thread[{t0mat, nor3}], Thread[{t0mat, nor4}],
  Thread[{t0mat, nor1 + nor2 + nor3 + nor4}]}, PlotRange -> Full, 
 Joined -> True, Frame -> True,
 PlotStyle -> {Dashed, Dashed, Dashed, Dashed, Thick},
 FrameLabel -> {"t02-t01", "Overlap"}, 
 PlotLegends -> {"nor1", "nor1", "nor3", "nor4", 
   "nor1+nor2+nor3+nor4"}]

enter image description here

from the integration result for each term, we can see that the integral for the expression

$$ Overlap=\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}dt_{1}dt_{2} \left(\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{1}-t_{01})^2}}{\sqrt[4]{\pi }}\frac{\sqrt{\Gamma } (e^{-\frac{1}{2} \Gamma ^2 (t_{2}-t_{02})^2}}{\sqrt[4]{\pi }}\right)^2 $$

will be incorrect for certain $t_{02}-t_{01}$ and then will affect the final result. I think this is just a very simple 2D Gaussian distribution and should be integrated easily. I don't know how to use Mathematica to do such kinds of integral and obtain relatively correct results. Can anyone help me with this?

$\endgroup$
1
  • $\begingroup$ Try adding the modes to the integration intervals: {t1, -Infinity, 0, t02, Infinity} and the same for t2. The problem looks like undersampling. The integrand probably underflows over a large portion of the domain, so you’ve got to help NIntegrate find its numeric support. It should be an otherwise easy integral. $\endgroup$
    – Michael E2
    Oct 31, 2021 at 23:36

2 Answers 2

1
$\begingroup$

Try adding the modes of the Gaussians to the integration intervals. The problem looks like undersampling. The integrand probably underflows over a large portion of the domain, so you’ve got to help NIntegrate find its numeric support. It should be an otherwise easy integral.

fτgaussian[Γ_, t_, 
   t0_ : 
    5] := (E^(-(1/2) (t - 
           t0)^2 Γ^2) Sqrt[Γ])/π^(1/4);

ListPlot[Thread[{Table[t02, {t02, 0, 20, 1}], Table[
    NIntegrate[(fτgaussian[1, t1, 0]*
         fτgaussian[1, t2, t02] + 
        fτgaussian[1, t2, 0]*fτgaussian[1, t1, t02])^2,
     {t1, -Infinity, 0, t02, Infinity},
     {t2, -Infinity, 0, t02, Infinity},
     Method -> {"MultidimensionalRule", "Generators" -> 9}],
    {t02, 0, 20, 1}]}],
 PlotRange -> Full, Joined -> True, Frame -> True, 
 FrameLabel -> {"t02-t01", "Overlap"}]

enter image description here

$\endgroup$
5
  • $\begingroup$ The documentation to "MultidimensionalRule" says only about the hypercube $ \left[-\frac{1}{2},\frac{1}{2}\right]^d \text{, } d\in \mathbb{N},d>1$. Is your usage to an improper integral documented? $\endgroup$
    – user64494
    Nov 1, 2021 at 5:58
  • $\begingroup$ Every method has its limitations: ListPlot[Thread[{Table[t02, {t02, 0, 160, 1}], Table[NIntegrate[(f\[Tau]gaussian[1, t1, 0]* f\[Tau]gaussian[1, t2, t02] + f\[Tau]gaussian[1, t2, 0]* f\[Tau]gaussian[1, t1, t02])^2, {t1, -Infinity, 0, t02, Infinity}, {t2, -Infinity, 0, t02, Infinity}, Method -> {"MultidimensionalRule", "Generators" -> 9}], {t02, 0, 160, 1}]}], PlotRange -> Full, Joined -> True, Frame -> True, FrameLabel -> {"t02-t01", "Overlap"}] ( 160 instead of 20) produces a bad plot $\endgroup$
    – user64494
    Nov 1, 2021 at 7:05
  • $\begingroup$ @user64494 To follow the advice to "find its numeric support" when $t_{01}-t_{02}$ is too large for the standard recursive subdivision, usually one subdivides further at $\mu\pm a\sigma$, with $a$ around $5$ or $6$ for a precision goal of around 6. As you are a user who has complained at others' interest in ridiculously large or small numbers, I wonder why 160 is at all useful. f\[Tau]gaussian[1, 0, 160] is on the order of $10^{-5660}$. If you just want to stymie Mma, use 10^8. --- For the (centered) unit cube and infinite intervals, read the rest of the NIntegrate documentation. $\endgroup$
    – Michael E2
    Nov 1, 2021 at 12:58
  • $\begingroup$ Thank you for your reply. You wrote " For the (centered) unit cube and infinite intervals, read the rest of the NIntegrate documentation". I repeat I don't find any example of the usage of "MultidimensionalRule" to improper integrals in the documentation to NIntegrate. Can you give a concrete reference? TIA. $\endgroup$
    – user64494
    Nov 1, 2021 at 18:13
  • $\begingroup$ @user64494 This is too short to completely document the workings of NIntegrate, but if you read it carefully, you will discover how NIntegrate handles infinite intervals. You may have to experiment to fill in missing details, but such is the nature of the documentation that we have. The option IntegrationMonitor, described on this site, is a useful tool for investigating the operation of NIntegrate. Good luck with your researches! $\endgroup$
    – Michael E2
    Nov 2, 2021 at 18:55
1
$\begingroup$

I don't know good numeric methods for improper multiple integrals. My best is

f\[Tau]gaussian[\[CapitalGamma]_, t_,t0_] := (E^(-(1/2) (t - 
       t0)^2 \[CapitalGamma]^2) Sqrt[\[CapitalGamma]])/\[Pi]^(1/
  4);

ListPlot[Thread[{Table[t02, {t02, 0, 20, 1}], Table[NIntegrate[(f\[Tau]gaussian[1, t1, 0]*
     f\[Tau]gaussian[1, t2, t02] + 
    f\[Tau]gaussian[1, t2, 0]*
     f\[Tau]gaussian[1, t1, t02])^2, {t1, -Infinity, 
  Infinity}, {t2, -Infinity, Infinity}, 
 Method -> {"AdaptiveQuasiMonteCarlo", "BisectionDithering" -> 0},
  PrecisionGoal -> 2, AccuracyGoal -> 3, MaxRecursion -> 100, 
 WorkingPrecision -> 20], {t02, 0, 20, 1}]}], PlotRange -> Full, 
 Joined -> True, Frame -> True, FrameLabel -> {"t02-t01", "Overlap"}]

enter image description here

, not so bad when comparing with the exact result

ListPlot[Thread[{Table[t02, {t02, 0, 20, 1}], 
   Table[Integrate[(f\[Tau]gaussian[1, t1, 0]*f\[Tau]gaussian[1, t2, t02] +   f\[Tau]gaussian[1, t2, 0]*f\[Tau]gaussian[1, t1, t02])^2, {t1, -Infinity, Infinity}, {t2, -Infinity, 
  Infinity}], {t02, 0, 20, 1}]}], PlotRange -> Full, 
 Joined -> True, Frame -> True, FrameLabel -> {"t02-t01", "Overlap"}]

enter image description here

$\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.