1
$\begingroup$

How can I evaluate these integrals to get a smooth outcome for the ratio, preferably without any errors?

I'm trying to integrate this:

f[x_, y_, a_, 
  s_] := (Cos[a*s*y/2])^2*(a^2 + 1)^(-1/2)*(x^2 + y^2 + 1)^(-1/
     2)*((a + x)^2 + y^2 + 1)^(-1/2) (Sqrt[1 + a^2] + 
     Sqrt[x^2 + y^2 + 1] + Sqrt[(a + x)^2 + y^2 + 1])^(-2)

nintf[a_, s_] := 
 NIntegrate[
  f[x, y, a, s], {x, -Infinity, Infinity}, {y, -Infinity, Infinity} ]

ratio[a_, s_] := 2*nintf[a, s]/nintf[a, 0] - 1 

ListPlot[Table[{a, ratio[a, 1]}, {a, 5, 15, 0.5}], Joined -> True]

For ease of reading, this is f:

$f = \frac{\cos^2(asy/2)}{\sqrt{a^2+1}\sqrt{x^2+y^2+1}\sqrt{(a+x)^2+y^2+1}\left(\sqrt{a^2+1}+\sqrt{x^2+y^2+1}+\sqrt{(a+x)^2+y^2+1}\right)^2}$

This gives me errors like so:

NIntegrate::slwcon: Numerical integration converging too slowly;suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 0.005791102372788279 and 1.7598057722389672*^-7 for the integral and error estimates.

I tried using TrapezoidalRule, which gets rid of the errors, but when I use ListPlot, the graph of nintf is not as smooth. So then I just kept on using the function above because it gave smooth plots for nintf albeit with a long time and the above errors.

But hen I plot ratio[a,1] in the regime {a, 5, 15, 0.5}, Mathematica does not give a smooth graph despite giving smooth graphs of nintf[a,1] and nintf[a,0] in that regime. I just used ListPlot like above for all the graphs.

Again, I am looking for a way to evaluate ratio, preferably without any errors. Thank you!

$\endgroup$
5
  • 1
    $\begingroup$ What is nintfcommnotrap ? You haven't defined that in your code. $\endgroup$
    – codebpr
    Commented Jan 22 at 6:48
  • $\begingroup$ What are the parameter ranges of a,s? $\endgroup$ Commented Jan 22 at 9:17
  • $\begingroup$ @codebpr Sorry, edited to make it more clear. $\endgroup$
    – roamer
    Commented Jan 22 at 11:49
  • $\begingroup$ @Ulrich Neumann they are both real numbers, a is positive. $\endgroup$
    – roamer
    Commented Jan 22 at 11:51
  • $\begingroup$ Calculating nintf[3, 1] in your notations, I obtained a warning "NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small." and a warning "NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 18 recursive bisections in y near {x,y} = {-356929.,1898.59}. NIntegrate obtained 0.020228384449357567` and 1.563032080098117*^-7 for the integral and error estimates." and reasonable 0.0202284. $\endgroup$
    – user64494
    Commented Jan 23 at 9:49

1 Answer 1

2
$\begingroup$

Change the integration subroutine as follows

nintf[a_?NumericQ, s_?NumericQ] := 
NIntegrate[
Evaluate[(Cos[a*s*y/2])^2*(a^2 + 1)^(-1/2)*(x^2 + y^2 + 1)^(-1/2)*((a + x)^2 + y^2 + 1)^(-1/2) (Sqrt[1 + a^2] + Sqrt[x^2 + y^2 + 1] + Sqrt[(a + x)^2 + y^2 + 1])^(-2)]
, {x, -Infinity,Infinity}, {y, -Infinity, Infinity}, AccuracyGoal -> 6,PrecisionGoal -> 4] 

zw = With[{s = 1}, Table[{a, s, 2 nintf[a, s]/nintf[a, 0] - 1}, {a, 1, 5,.25}]]  

ListPlot[zw[[All, {1, -1}]],PlotLabel -> "s=" <> ToString[zw[[1, 2]]],  AxesLabel -> {a, ratio}]

enter image description here

Also ratio[a,s]seems to be smooth:

zw = Table[{a, s, 2 nintf[a, s]/nintf[a, 0] - 1}, {a, 1,5, .25}, {s, 0, 1.5, .1}] // Flatten[#, 1] &

ListPlot3D[zw, PlotRange -> All, AxesLabel -> {a, s, ratio[a, s]},MeshFunctions -> (#3 &) ]

enter image description here

$\endgroup$
15
  • $\begingroup$ I don't know good numeric methods for multiple improper integrals. For example, NIntegrate[ Sin[x^2 + y^2]^2/(x^2 + y^2), {x, -Infinity, Infinity}, {y, -Infinity, Infinity}, AccuracyGoal -> 6, PrecisionGoal -> 4] results in 33.0768 without any warnings/errrors, though the integral diverges. Did you pay your attention to my answer where is proven the divergence? $\endgroup$
    – user64494
    Commented Jan 22 at 21:11
  • $\begingroup$ The divergence can be establishing by switching to the polar coordinates: $\int_0^\infty \sin(r^2)^2r/r^2\,dr$ diverges and Integrate[Sin[r^2]^2/r, {r, 0, Infinity}] says "Integrate::idiv: Integral of Sin[r^2]^2/r does not converge on {0,[Infinity]}.". $\endgroup$
    – user64494
    Commented Jan 22 at 21:25
  • $\begingroup$ @user64494 I recognized your answer but couldn't follow your argumentation. You only showed that integral of the asymptotes diverges. $\endgroup$ Commented Jan 22 at 21:37
  • $\begingroup$ This implies the divergence of the integral under consideration. Look at the result of Series[(r*Cos[(a*r*s*Sin[\[Theta]])/2]^2)/(Sqrt[(1 + a^2)*(1 + r^2)]* Sqrt[1 + a^2 + r^2 + 2*a*r*Cos[\[Theta]]]), {r, Infinity, 2}] // Normal. $\endgroup$
    – user64494
    Commented Jan 22 at 21:43
  • $\begingroup$ Hi, I'm a bit confused, this all evaluates well, but I get an error for the first ListPlot: ListPlot::nonopt: Options expected (instead of Null) beyond position 1 in ListPlot[<<1>>]. An option must be a rule or a list of rules. $\endgroup$
    – roamer
    Commented Jan 23 at 6:14

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.