WorkingPrecision and AccuracyGoal for using with NIntegrate

Question

I'm modeling a quantum tunneling barrier. Because I want to express my quantities in mks units, that means that some of the variables I'm using are pretty small (10^-9).

I use NIntegrate for part of the program. If I try and use it naively:

Vapp = 1;
tunder = TRaw[energy, 0, bar - V/2, -V, s] /. V -> Vapp;
fermilevel = 2;
integral = 
 NIntegrate[Vapp*tunder, {energy, 0, Max[fermilevel - Vapp, 0]}]

Where tunder is this horrendous function of energy (sorry for the pastebin, it was taking up too much space here): http://pastebin.com/aJzea7Nb

I get these two errors:

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

NIntegrate::ncvb: "NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in energy near {energy} = {0.912094}. NIntegrate obtained 3.3633806566925494`*^-18 and 2.6134212824171803`*^-17 for the integral and error estimates.

It stops giving me those errors if I specify some values of WorkingPrecision and AccuracyGoal, but I don't really know why or how to choose proper values:

accgoal = 12;
WP = 12;
integral = 
 NIntegrate[Vapp*tunder, {energy, 0, Max[fermilevel - Vapp, 0]}, 
  AccuracyGoal -> accgoal, WorkingPrecision -> WP]

It also gives a fairly different answer for the integral in either case.

I apologize if this is a simple question, but I can't figure out what do to. I thought that if I increase the AccuracyGoal, it would take longer to calculate, but be more accurate and not give me errors, so I tried AccuracyGoal->40, but that actually causes more errors.

What could be causing this? thank you!

edit: Okay, I read a few more threads, but I'm still confused. Here's what I tried. I manually set the argument of NIntegrate to have a precision (30) higher than the workingprecision and accuracygoal (25), and that gets rid of errors. I still don't really know what's an appropriate value for them, so I plotted the result of my integral for many values of the precision, where I set the precision, working precision, and accuracy goal to the same:

inttable = 
  Table[{WP, 
    NIntegrate[
     SetPrecision[Vapp*tunder, WP], {energy, 0, 
      Max[fermilevel - Vapp, 0]}, AccuracyGoal -> WP, 
     WorkingPrecision -> WP]}, {WP, 10, 30}];
ListPlot@inttable

Which gives me:

Where the horizontal axis is the workingprecision/accuracygoal/setprecision value and the vertical axis is the integral value.

My naive interpretation of this is that, because it levels off at around 22, it needs to be about that, but not much higher. Is that right?

However, this is for just one integral. It's possible that for other values of the parameters in it, I'd need more precision.

Answer 1

Try

f0 = Rationalize[tunder];
NIntegrate[f0, {energy, 0, Max[fermilevel - Vapp, 0]}, WorkingPrecision -> 20]
(*  1.0319078410138864350*10^-17  *)

---

FWIW...

NIntegrate obtained 3.3633806566925494`*^-18 and 2.6134212824171803`*^-17 for the integral and error estimates.

Since

-Log10[2.6134212824171803`*^-17]
(* 16.5828  *)

this means that an AccuracyGoal higher than about 16.5 will result in a warning related to the goal. (The AccuracyGoal sets a limit below which errors will be treated as "zero" or "acceptable.")

Since

-Log10[3.3633806566925494`*^-18/2.6134212824171803`*^-17]
(*  0.890433  *)

this means that a PrecisionGoal higher than about 0.9 will result is a warning related to the goal.

NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in energy near {energy} = {0.912094}.

Since the default recursion setting is MaxRecursion -> 9, you might consider raising the setting to achieve the desired accuracy indicated by AccuracyGoal.

...suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

One might check the first three of the four possibilities by plotting the integrand (when low dimensional).

Generally, one can tell that WorkingPrecision is high enough, when the computed value converges without error warnings. In the OP's code, this happens at about a setting of 21. For me, I get convergence at a setting of 20:

vals = Table[
   NIntegrate[f0, {energy, 0, Max[fermilevel - Vapp, 0]}, WorkingPrecision -> wp],
   {wp, 15, 30}];

Differences[vals]/Rest[vals]  (* precision *)
(*
  {Indeterminate, Indeterminate, 1.00000000000000000, 
   0.431977933335144191, 0.019191433540064037, 0.*10^-20, 0.*10^-21, 
   0.*10^-22, 0.*10^-23, 0.*10^-24, 0.*10^-25, 0.*10^-26, 0.*10^-27, 
   0.*10^-28, 0.*10^-29}
*)

In practice, I usually don't check successive values. That can be slow. Doubling PrecisionGoal or increasing it by 50% usually what I do. Note that increasing WorkingPrecision increases PrecisionGoal, so sometimes one needs to fix it at an appropriate setting to avoid convergence warnings.