DiscretePlot of a function containing Nintegrate under Sum

Question

I am trying to get DiscretePlot of the following integral containing laguerre polynomial, the Mathematica code of which is given below:

A1[l0_, ρ_, n_, m_, q_, h_, M_, 
  v_] := (Sqrt[(2*(q - (Abs[m] + m)/2)!)/(l0^2*q!)])*
  Exp[-((I/(2*h)*M*v*ρ^2)/l0)]*
  Exp[-(ρ^2/(2*l0^2))]*((ρ^2/l0^2)^(Abs[m]/2))*
  LaguerreL[q - (Abs[m] + m)/2, Abs[m], (ρ^2/l0^2)]*(Sqrt[(2*(n - (Abs[m] + m)/2)!)/(l0^2*n!)])*
  Exp[-(ρ^2/(2*l0^2))]*((ρ^2/l0^2)^(Abs[m]/2))*
  LaguerreL[n - (Abs[m] + m)/2, Abs[m], (ρ^2/l0^2)]*ρ
E0[q_, m_, w0_, Ω_, 
  h_] := (q + 1/2)*h*(w0 + Ω) + (q - m + 1/2)*
   h (w0 - Ω) + (4 + 1/2)*h*10
p1[l0_, l1_, ρ_, n1_, m_, q_, h_, M_, v_, 
  w0_, Ω_, τ_] := (Sqrt[(
   2*(n1 - (Abs[m] + m)/2)!)/(l1^2*n1!)])*
  Exp[-(ρ^2/(2*l1^2))]*((ρ^2/l1^2)^(Abs[m]/2))*
  LaguerreL[n1 - (Abs[m] + m)/2, Abs[m], (ρ^2/l1^2)]*(Sqrt[(2*(q - (Abs[m] + m)/2)!)/(l1^2*q!)])*
  Exp[(I/(2*h)*M*v*ρ^2 - 
    I/h*l0*E0[q, m, w0, Ω, h]*τ)/l1]*
  Exp[-(ρ^2/(2*l1^2))]*((ρ^2/l1^2)^(Abs[m]/2))*
  LaguerreL[q - (Abs[m] + m)/2, Abs[m], (ρ^2/l1^2)]*ρ
p2[n_, q_] := 
 NIntegrate[A1[1, ρ, n, 1, q, 1, 1, 5], {ρ, 0, Infinity}]*
  NIntegrate[
   p1[1, 2, ρ, 2, 1, q, 1, 1, 5, 10, 5, 10], {ρ, 0, 
    Infinity}]
p3[n_] := Sum[p2[n, q], {q, 1, 10}]
DiscretePlot[Abs[p3[n]]^2, {n, 1, 6, 1}, PlotStyle -> Red]

But it shows the following error:

NIntegrate::errprec: Catastrophic loss of precision in the global error estimate due to insufficient WorkingPrecision or divergent integral.

I searched in the already discussed question in Mathematica stack exchange related to this kind of error and I got a particular case under heading "Catastrophic loss of precision", where it is suggested to avoid the error by using Table command, but I could not follow the same properly in my case.

Would you kindly suggest me how to amend my code such that I can get the desired plot?

Maybe you should focus on which call to NIntegrate causes the error. Alternatively, the message is implying you should try increasing WorkingPrecision. (Table seems an unlikely fix, though the separate integrals would be returned instead of their sum. That might be helpful.) — Michael E2, Sep 18 '22 at 16:42
Also, using ?NumericQ on the arguments of p2[], while it probably won't fix this problem, is good practice. — Michael E2, Sep 18 '22 at 16:48
Sir, thanks for your kind response...It sppears to me on going through your valuable comments that the order of operation matters...i.e. 2 input indexes such as n and q should be specified maintaining the order....Sir, if you dont mind, would you kindly show me the ameneded code of p2[n_,q_] inserting ?NumericQ — R. Bhattacharya, Sep 18 '22 at 18:29
p2[n_?NumericQ,q_?NumericQ] := <..rest of definition..> is how it is shown in the "Applications" section of the NumericQ documentation and the examples in the link in my previous comment. — Michael E2, Sep 18 '22 at 18:47
Sir, I have written the same thing but still it does not show anything.. p2[n_?NumericQ, q_?NumericQ] := NIntegrate[A1[1, \[Rho], n, 1, q, 1, 1, 5], {\[Rho], 0, Infinity}]* NIntegrate[ p1[1, 2, \[Rho], 2, 1, q, 1, 1, 5, 10, 5, 10], {\[Rho], 0, Infinity}]p3[n_?NumericQ] := Sum[p2[n?NumericQ, q?NumericQ], {q, 1, 10}]DiscretePlot[Abs[p3[n?NumericQ]]^2, {n, 1, 6, 1}, PlotStyle -> Red] — R. Bhattacharya, Sep 18 '22 at 19:02
(1) You seemed to have added ?NumericQ in more than the two places I indicated. (2) What does not show anything? (You said "it.") I did say I didn't think this would solve the problem with your code. To do that, you'll have to fix the NIntegrate call that gives the error. Have you determined which on that is? — Michael E2, Sep 18 '22 at 19:14

score 0 · Accepted Answer · answered Sep 19 '22 at 18:08

NIntegrate has problems, because integrand is highly osccilatroy.

But this can all be done with Integrate in 8 seconds.

A1[l0_, \[Rho]_, n_, m_, q_, h_, M_, 
  v_] := (Sqrt[(2*(q - (Abs[m] + m)/2)!)/(l0^2*q!)])*
  Exp[-((I/(2*h)*M*v*\[Rho]^2)/l0)]*
  Exp[-(\[Rho]^2/(2*l0^2))]*((\[Rho]^2/l0^2)^(Abs[m]/2))*
  LaguerreL[q - (Abs[m] + m)/2, 
   Abs[m], (\[Rho]^2/
     l0^2)]*(Sqrt[(2*(n - (Abs[m] + m)/2)!)/(l0^2*n!)])*
  Exp[-(\[Rho]^2/(2*l0^2))]*((\[Rho]^2/l0^2)^(Abs[m]/2))*
  LaguerreL[n - (Abs[m] + m)/2, Abs[m], (\[Rho]^2/l0^2)]*\[Rho]
E0[q_, m_, w0_, \[CapitalOmega]_, 
  h_] := (q + 1/2)*h*(w0 + \[CapitalOmega]) + (q - m + 1/2)*
   h (w0 - \[CapitalOmega]) + (4 + 1/2)*h*10
p1[l0_, l1_, \[Rho]_, n1_, m_, q_, h_, M_, v_, 
  w0_, \[CapitalOmega]_, \[Tau]_] := (Sqrt[(2*(n1 - (Abs[m] + m)/
           2)!)/(l1^2*n1!)])*
  Exp[-(\[Rho]^2/(2*l1^2))]*((\[Rho]^2/l1^2)^(Abs[m]/2))*
  LaguerreL[n1 - (Abs[m] + m)/2, 
   Abs[m], (\[Rho]^2/
     l1^2)]*(Sqrt[(2*(q - (Abs[m] + m)/2)!)/(l1^2*q!)])*
  Exp[(I/(2*h)*M*v*\[Rho]^2 - 
      I/h*l0*E0[q, m, w0, \[CapitalOmega], h]*\[Tau])/l1]*
  Exp[-(\[Rho]^2/(2*l1^2))]*((\[Rho]^2/l1^2)^(Abs[m]/2))*
  LaguerreL[q - (Abs[m] + m)/2, Abs[m], (\[Rho]^2/l1^2)]*\[Rho]
p2[n_, q_] := 
 Integrate[A1[1, \[Rho], n, 1, q, 1, 1, 5], {\[Rho], 0, Infinity}]*
  Integrate[
   p1[1, 2, \[Rho], 2, 1, q, 1, 1, 5, 10, 5, 10], {\[Rho], 0, 
    Infinity}]
p3[n_] := Sum[p2[n, q], {q, 1, 10}]
(tab = Table[Abs[p3[n]]^2, {n, 1, 6}]; tab // N) // AbsoluteTiming
(*   {8.2475393, {0.000440122, 0.000138649, 0.00142201, 0.00274115, 
  0.00326648, 0.0031438}}   *)

Thank you, Sir... for your kind response...Yes, I have also done the same thing using Integrate instead of NIntegrate..Actually, Sir MicaelE2 pointed out the same thing that somehow NIntegarte created the problem....Thank you once again Sir... — R. Bhattacharya, Sep 20 '22 at 05:59

DiscretePlot of a function containing Nintegrate under Sum

1 Answers1