3
$\begingroup$

I am trying to solve a 6D integral over this integrated

integrand=ConditionalExpression[(9*(1 + E^(R*\[Eta]1))^2*(1 +E^(R*\[Eta]2))^2*R^6*(\[Eta]1 - \[Xi]1)*(\[Eta]1 + \[Xi]1)*(\[Eta]2 - \[Xi]2)*(\[Eta]2 + \[Xi]2))/(E^(R*(-2 + \[Eta]1 + \[Eta]2 + \[Xi]1 + \[Xi]2))*(128*Pi^2*(3 + 3*E^R + R*(3 + R))^2*Sqrt[R^2*(-2 + \[Eta]1^2 + \[Eta]2^2 + \[Xi]1^2 - 2*\[Eta]1*\[Eta]2*\[Xi]1*\[Xi]2 + \[Xi]2^2 - 2*Sqrt[-((-1 + \[Eta]1^2)*(-1 + \[Xi]1^2))]*Sqrt[-((-1 + \[Eta]2^2)*(-1 + \[Xi]2^2))]*Cos[\[Phi]1 - \[Phi]2])])), Re[R] > 0]

Vee = Integrate[integrand, {\[Xi]1, 1, \[Infinity]}, {\[Eta]1, -1, 1}, {\[Phi]1, 0, 2*\[Pi]}, {\[Xi]2, 1, \[Infinity]}, {\[Eta]2, -1, 1}, {\[Phi]2, 0, 2*\[Pi]}]

with the limits $\xi_1,\xi_2$ from 1 to $\infty$, $\eta_1$,$\eta_2$ from -1 to 1 and $\phi_1,\phi_2$ from 0 to $2\pi$. By the way, this is calculating the electrostatic energy of the H$_2$ molecule in confocal, elliptic coordinates. The problem is, this integral is running for a long time and does not finish. I am not sure if it is working properly or not, since there is no response, it just keeps running. Any idea how to solve or speed up this integral?

$\endgroup$
7
  • $\begingroup$ People here generally like users to post code as copyable Mathematica code instead of just images or TeX, so they can copy-paste it. It makes it convenient for them and more likely you will get someone to help you. You may find the meta Q&A, How to copy code from Mathematica so it looks good on this site, helpful $\endgroup$
    – Michael E2
    Commented May 22, 2023 at 0:47
  • $\begingroup$ I tried to follow the procedure. Is it better now? $\endgroup$
    – Guiste
    Commented May 22, 2023 at 0:54
  • $\begingroup$ yes, thanks. Not sure whether it will be easy for someone to help or not. If the conditional expression is meant to add an assumption that Re[R] > 0, the typical way to do that is adding the option Assumptions -> Re[R] > 0 to Integrate. $\endgroup$
    – Michael E2
    Commented May 22, 2023 at 1:01
  • $\begingroup$ I will try that, definitely R is a real number and larger than 0 $\endgroup$
    – Guiste
    Commented May 22, 2023 at 1:02
  • $\begingroup$ "R is a real number and larger than 0" -- then use Assumptions -> R > 0. Variables that appear in inequalities are assumed to be real. Re[R] > 0 signals that R is complex, with a positive real part. (Integrate[] still seems to run forever, though.) $\endgroup$
    – Michael E2
    Commented May 22, 2023 at 1:05

1 Answer 1

2
$\begingroup$

We can reduce dimension from 6 to 4 using

 Integrate[
 1/Sqrt[a Cos[p1 - p2] + b], {p1, 0, 2 Pi}, {p2, 0, 2 Pi}, 
 Assumptions -> {a > 0, b > 0}]

(*Out[1]= (8 \[Pi] EllipticK[(2 a)/(a + b)])/Sqrt[a + b] *)
 % /. {a -> -2*Sqrt[-((\[Eta]1^2 - 1)*(\[Xi]1^2 - 1))]*Sqrt[-((\[Eta]2^2 - 1)*(\[Xi]2^2 - 1))], b -> \[Eta]1^2 - 2*\[Eta]1*\[Eta]2*\[Xi]1*\[Xi]2 + \[Eta]2^2 + \[Xi]1^2 + \[Xi]2^2 - 2}

(*Out[2]= (8 \[Pi] EllipticK[-((
   4 Sqrt[-((-1 + \[Eta]1^2) (-1 + \[Xi]1^2))]
     Sqrt[-((-1 + \[Eta]2^2) (-1 + \[Xi]2^2))])/(-2 + \[Eta]1^2 + \
\[Eta]2^2 + \[Xi]1^2 - 2 \[Eta]1 \[Eta]2 \[Xi]1 \[Xi]2 + \[Xi]2^2 - 
    2 Sqrt[-((-1 + \[Eta]1^2) (-1 + \[Xi]1^2))]
      Sqrt[-((-1 + \[Eta]2^2) (-1 + \[Xi]2^2))]))])/Sqrt[-2 + \
\[Eta]1^2 + \[Eta]2^2 + \[Xi]1^2 - 
 2 \[Eta]1 \[Eta]2 \[Xi]1 \[Xi]2 + \[Xi]2^2 - 
 2 Sqrt[-((-1 + \[Eta]1^2) (-1 + \[Xi]1^2))]
   Sqrt[-((-1 + \[Eta]2^2) (-1 + \[Xi]2^2))]]*)

Therefore integral can be evaluated in 4 dimensional space as follows

f[R_?NumericQ] := (9*R^5)/(128*
     Pi^2*(R*(R + 3) + 3*E^R + 
        3)^2) NIntegrate[(8*(1 + E^(R*\[Eta]1))^2*(1 + 
         E^(R*\[Eta]2))^2*
      Pi*(\[Eta]1 - \[Xi]1)*(\[Eta]1 + \[Xi]1)*(\[Eta]2 - \[Xi]2)*(\
\[Eta]2 + \[Xi]2)*
      EllipticK[-((4*Sqrt[-((-1 + \[Eta]1^2)*(-1 + \[Xi]1^2))]*
            Sqrt[-((-1 + \[Eta]2^2)*(-1 + \[Xi]2^2))])/(-2 + \
\[Eta]1^2 + \[Eta]2^2 + \[Xi]1^2 - 
            2*\[Eta]1*\[Eta]2*\[Xi]1*\[Xi]2 + \[Xi]2^2 - 
            2*Sqrt[-((-1 + \[Eta]1^2)*(-1 + \[Xi]1^2))]*
             Sqrt[-((-1 + \[Eta]2^2)*(-1 + \[Xi]2^2))]))])/(E^(R*(-2 \
+ \[Eta]1 + \[Eta]2 + \[Xi]1 + \[Xi]2))*
      Sqrt[-2 + \[Eta]1^2 + \[Eta]2^2 + \[Xi]1^2 - 
        2*\[Eta]1*\[Eta]2*\[Xi]1*\[Xi]2 + \[Xi]2^2 - 
        2*Sqrt[-((-1 + \[Eta]1^2)*(-1 + \[Xi]1^2))]*
         Sqrt[-((-1 + \[Eta]2^2)*(-1 + \[Xi]2^2))]]), {\[Xi]1, 1, 
    Infinity}, {\[Eta]1, -1, 1}, {\[Xi]2, 1, Infinity}, {\[Eta]2, -1, 
    1}, Method -> "AdaptiveQuasiMonteCarlo"]; 

Example of usage

f[1] // AbsoluteTiming

Out[]= {0.412131, 0.588934}

Also this integral can be computed in 6D as follows

integrand = (9*R^6*(\[Eta]1 - \[Xi]1)*(\[Eta]1 + \[Xi]1)*(\[Eta]2 - \[Xi]2)*(\[Eta]2 + \[Xi]2)*(E^(\[Eta]1*R) + 1)^2*
     (E^(\[Eta]2*R) + 1)^2)/(E^(R*(\[Eta]1 + \[Eta]2 + \[Xi]1 + \[Xi]2 - 2))*(128*Pi^2*(R*(R + 3) + 3*E^R + 3)^2*
      Sqrt[R^2*(-(2*Sqrt[-((\[Eta]1^2 - 1)*(\[Xi]1^2 - 1))]*Sqrt[-((\[Eta]2^2 - 1)*(\[Xi]2^2 - 1))]*
           Cos[\[Phi]1 - \[Phi]2]) + \[Eta]1^2 - 2*\[Eta]1*\[Eta]2*\[Xi]1*\[Xi]2 + \[Eta]2^2 + \[Xi]1^2 + \[Xi]2^2 - 2)])); 
lst = Table[{x, 
   NIntegrate[
    integrand /. R -> x, {\[Xi]1, 1, Infinity}, {\[Eta]1, -1, 
     1}, {\[Phi]1, 0, 2*Pi}, {\[Xi]2, 1, Infinity}, {\[Eta]2, -1, 
     1}, {\[Phi]2, 0, 2*Pi}, 
    Method -> "AdaptiveQuasiMonteCarlo"]}, {x, 1/10, 5, 1/10}]

Visualization

ListPlot[{lst, {{1, f[1]}}}, PlotStyle -> {Blue, Red}, Frame -> True, 
 FrameLabel -> {"R", "Vee"}, PlotLegends -> {"6D", "4D"}]

Figure 1

$\endgroup$
2
  • $\begingroup$ Yeah I guess numerical integration is better in this case. I was hoping for an analytic solution, but Mathematica seems not to be able to give it to me. $\endgroup$
    – Guiste
    Commented May 31, 2023 at 0:34
  • $\begingroup$ @Guiste Yes, you are right. Analytical expression in this case maybe is not so useful. Even reduced from 6D in 4D looks very cumbersome. $\endgroup$ Commented May 31, 2023 at 4:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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