TransWikia.com

NIntegrate of a convergent integral working with large integration limits, but not with infinite integration limits

Mathematica Asked by Henry Shackleton on March 31, 2021

Code for reproducing is below:

Integrand[x_] := ((x + 1) Abs[x] Log[(x + 1)^2/Abs[x]^3])/((x + 1)^2 - Abs[x]^3);
    
NIntegrate[Integrand[x], {x, -10000000, 10000000}]
NIntegrate[Integrand[x], {x, -Infinity, Infinity}]

The first instance of NIntegrate gives $2.48398$, and one can check by plotting the result of NIntegrate as a function of the integration limits that the integral appears to converge to this answer for very large integration limits. However, the second instance of NIntegrate gives a completely different answer of $1.75434 times 10^{8}$. What’s going on here, and how can I make the two integrals agree? I think something that NIntegrate might be having trouble with is that the convergence of the integral requires the fact that the integrand is essentially odd for very large $|x|$ and these contributions tend to cancel out. Ideally I’d like to find some way to get Mathematica to deal with this properly with the infinite integration range, so that I don’t have to keep remembering to exchange the infinite integration limits with large finite values every time I need to do an integral like this.

2 Answers

Update

Assuming that the principal value is wanted, the method below is probably the best way, which computes the PV at infinity by reflecting the integration over the negative axis onto the positive axis.

NIntegrate[Integrand[x] + Integrand[-x], {x, 0, Infinity}]

(*  2.8130684254425558`  *)

@user64494 pointed out that it's a divergent integral (at infinity) and that Method -> "PrincipalValue" computes the principal value only at finite points and not at infinity. Add the singularity at 1 to the interval to get higher precision:

    NIntegrate[Integrand[x] + Integrand[-x], {x, 0, 1, Infinity},...]

(See older part below for discussion of the singularities.)

The OP's NIntegrate[Integrand[x], {x, -10000000, 10000000}] in essence approximates the PV.


I can get them to agree in this way:

sing = SortBy[
   x /. Solve[FunctionSingularities[Integrand[x], x], x, Reals], N];
N@sing

(*  {-1., -0.56984, 0., 2.1479}  *)

NIntegrate[integrand[x], 
 Evaluate@{x, -10000000, Sequence @@ sing, 10000000}]
NIntegrate[integrand[x],
 Evaluate@{x, -Infinity, Sequence @@ sing, Infinity}, 
 WorkingPrecision -> 32, PrecisionGoal -> 6, AccuracyGoal -> 20, 
 Method -> "PrincipalValue"]

(*
  2.813065409383123       (no warnings)
  2.81306842613267712727  (warnings omitted)
*)

At -1 and 0, there seem to be infinite derivatives, so-called "weak" singularities. At the other two "singularities" there might be numerical difficulties from subtractive cancellation (e.g. they're where the denominator (1 + x)^2 - Abs[x]^3 is zero). Weak singularities tend to make convergence harder.

Correct answer by Michael E2 on March 31, 2021

Numerical integration can be without any problems by splitting the integration into two parts:

f1 = FullSimplify[Integrand[x], Assumptions -> x > 0];
f2 = FullSimplify[Integrand[-x], Assumptions -> x > 0];
NIntegrate[f1 + f2, {x, 0, ∞}]
(*2.81307*)

Higher accuracy can be reached (again without any warnings as follows)

NIntegrate[f1 + f2, {x, 0, 1, ∞}, WorkingPrecision -> 100]
 (*2.81306842544297991357319331946335010348754523474644575568568457957304
   5154609418329813658075049081527*)

Answered by yarchik on March 31, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP