# Speed up this NIntegrate

Is there any trick to speed up this numerical integral:

Clear[Cef2];
Cef2[b1_, b2_, x_, y_, r_] = E^(-r - E^-r Cosh[r (x^2 + y^2)]) ((E^(x y E^(-2 r)) Sin[(y b1 + x b2)] + E^(x y ) Cos[(y b1 - x b2)]));
NIntegrate[Cef2[0.1243, 1.1321, x, y, 2.51], {x,0, Infinity}, {y, 0, Infinity}]


EDIT: for the final result I do not care to have a precision better than $$10^{-5}$$

• Work it out by hand and then implement it using weighting functions applied to the scalar portions of your functions you will integrate! So you'll just be assembling the final summation using slots and algebraic operations instead of Integrating to Infinity. Alternately & additionally, provide some assumptions to your inputs, so the system does not have to take time by assuming they are anything but numerical values. – CA Trevillian May 4 '19 at 17:08
• Related. – corey979 May 4 '19 at 17:27

NIntegrate can spend a distressing amount of time trying to simplify the integrand symbolically if you allow it. In some cases, like this one, it is significantly faster to simply force it to treat it as a black-box numeric function by defining the function to only take numeric arguments. For direct comparison, here is the original definition of the function with AbsoluteTiming stuck on the end:

Clear[Cef2];
Cef2[b1_, b2_, x_, y_, r_] =
E^(-r - E^-r Cosh[
r (x^2 + y^2)]) ((E^(x y E^(-2 r)) Sin[(y b1 + x b2)] +
E^(x y) Cos[(y b1 - x b2)]));
NIntegrate[
Cef2[0.1243, 1.1321, x, y, 2.51], {x, 0, Infinity}, {y, 0,
Infinity}] // AbsoluteTiming


{121.665, 0.103711}

This takes 122 seconds and finds the result to be 0.103711 (though it does throw a few warnings out that the integral is converging slowly).

Redefining the function to use ?NumericQ (and := instead of =, as that can cause some unexpected issues):

Clear[Cef2];
Cef2[b1_?NumericQ, b2_?NumericQ, x_?NumericQ, y_?NumericQ,
r_?NumericQ] :=
E^(-r - E^-r Cosh[
r (x^2 + y^2)]) ((E^(x y E^(-2 r)) Sin[(y b1 + x b2)] +
E^(x y) Cos[(y b1 - x b2)]));
NIntegrate[
Cef2[0.1243, 1.1321, x, y, 2.51], {x, 0, Infinity}, {y, 0,
Infinity}] // AbsoluteTiming


{0.0916949, 0.103711}

This takes less than a tenth of a second to achieve the same result, though it still complains about the convergence rate being slow.

For the final result I care to have a precision of at most 10^-5, can this be useful to make it even faster?

Yes, you can use this piece of information to make this integral even faster. Just add PrecisionGoal -> 4 to the NIntegrate. This tells NIntegrate that you only need 4 digits of precision (covering $$10^{-1}$$ through $$10^{-4}$$ and then the $$10^{-5}$$ digit is a matter of luck).

This provides another roughly 50% speedup (0.0532 seconds compared to the above 0.092), though the evaluation times are small enough already that this isn't going to be a reliable or consistent measure of that.

• Wow! Question: Did you clear the cache before you ran the second result? If so, I LOVE the speed-up we get when we tell the computer just a liiiiiittle bit more about our problem, instead of letting it think everything about it on its own! And Nice catch on the set-delayed, I did not even notice!!! @m137, if you're running this a lot, for many different functional inputs, then accessing them again multiple times, memoization will valuable, as for the convergence rate (barring some group misunderstanding on our parts) just silence the message from being outputted. – CA Trevillian May 4 '19 at 17:19
• I didn't explicitly clear the cache, but it's about the same speed on a completely fresh kernel. While the convergence rate error can be silenced, I would highly recommend not doing so, as that's putting a band-aid over the limitations of the numerical methods. I'd much rather the reminder that the result may not be perfectly accurate than be lulled into a false sense of security by Quiet. – eyorble May 4 '19 at 17:24
• Ah yes, I agree with this view definitely. Mt silencing suggestion comes from my normal use of internal implementations of functions I (usually) know enough about to either ignore the message with silence, or because it is an internal function, there is a slowdown associated with the message being outputted during evaluation, and then it just doesn't look very pretty when its being printed out in the middle of my outputs. – CA Trevillian May 4 '19 at 17:30
• @eyorble thanks, your answer is super useful. The speed up is amazing. P.S. For the final result I care to have a precision of at most 10^-5, can this be useful to make it even faster? – m137 May 4 '19 at 20:36
• @m137 Yes, you can use PrecisionGoal to lower NIntegrate's tolerances if you don't need the precision. If PrecisionGoal -> 4 isn't quite what you meant, it's still a speed up for PrecisionGoal -> 5, just less dramatic. They turn out to get the same result for this problem. – eyorble May 4 '19 at 23:09

Use the option setting "SymbolicProcessing"->0. (Given/prompted at Wolfram Support's "How do I accelerate NIntegrate evaluations?".)

Clear[Cef2];
Cef2[b1_, b2_, x_, y_, r_] :=
E^(-r - E^-r Cosh[
r (x^2 + y^2)]) ((E^(x y E^(-2 r)) Sin[(y b1 + x b2)] +
E^(x y) Cos[(y b1 - x b2)]));

AbsoluteTiming[
NIntegrate[
Cef2[0.1243, 1.1321, x, y, 2.51], {x, 0, Infinity}, {y, 0, Infinity},
Method -> {Automatic, "SymbolicProcessing" -> 0}]
]

(* During evaluation of In[6]:= 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. *)

(* {0.025366, 0.103711} *)


## Faster results (when using higher working precision and precision goal)

P.S. For the final result I care to have a precision of at most 10^-5, can this be useful to make it even faster?

It seems it is a good idea to compute this integral with a Cartesian product rule. (See Method->"GaussKronrodRule" used below.)

This makes sense -- we get "NIntegrate::slwcon" with the standard multi-dimensional rule, hence using the Cartesian rule (which fills more densely the integration domain with sampling points) we might get results faster.

AbsoluteTiming[
NIntegrate[
Cef2[1243/10000, 11321/10000, x, y, 251/100], {x, 0, Infinity}, {y, 0, Infinity},
Method -> {"GlobalAdaptive", Method -> "GaussKronrodRule",
"SymbolicProcessing" -> 0}, WorkingPrecision -> 30,
PrecisionGoal -> 10]
]

(* {0.632539, 0.103711103761559644223925578653} *)


Note that I rationalized the numerical arguments given to Cef2.

Also, in the method specification above replacing "GlobalAdaptive" with "LocalAdaptive" is not going to produce results for at least 5 minutes.

• I was about to post the same thing. The option "LocalAdaptive" is somewhat faster though. – AccidentalFourierTransform May 4 '19 at 20:06
• @AntonAntonov thank you for your answer, also this is super useful. The speed up is amazing. P.S. For the final result I care to have a precision of at most 10^-5, can this be useful to make it even faster? – m137 May 4 '19 at 20:37
• @m137 See my update. – Anton Antonov May 4 '19 at 21:03