2
$\begingroup$

I have a very complicated 2D integral that I need to calculate repeatedly, and I'm trying to speed it up a bit, since at the moment it's taking a couple of days to complete.

One thing I've noticed is that when NIntegrate samples the integrand, it sometimes dives into very high precision numbers, even though WorkingPrecision is set to $MachinePrecision.

Here's a fairly minimal example to demonstrate:

f[x_?NumericQ] := (result = Exp[Sin[Sqrt[x]]]; Sow[{x, result}]; result)

{result, {samples}} = Reap[NIntegrate[(f[x]), {x, 0, 1}, 
    WorkingPrecision -> MachinePrecision]]

Tally[Precision /@ samples]

The last result is then

{{MachinePrecision, 122}, {65.9546, 7}}

I suspect that the high-precision is necessary to achieve the requested PrecisionGoal or AccuracyGoal, but in the case of my real integrand, the underlying dataset simply doesn't have data at that level of precision (it's mostly linear interpolations of experimental data, with some complicated transformations). Since MachinePrecision is much faster than arbitrary precision, I'd rather stick to it if possible.

Interestingly, there doesn't seem to be any particular features at the points that are being focussed on. In the example above, it samples multiple points very close to 0.06, and I can't see any special significance to that part of the function. I have a similar problem with my real integrand.

I've tried setting $MaxExtraPrecision, but I'm not sure that it has any impact on NIntegrate.

Is there some other way to tell Mathematica that it should stick to $MachinePrecision?

$\endgroup$
  • $\begingroup$ Use N on x inside f? $\endgroup$ – Michael E2 Mar 23 '15 at 13:05
  • $\begingroup$ One can see from ListPlot[samples[[2 ;;, 1]]] that after 4 applications of the "GlobalAdaptive" strategy with the "GaussKronrodRule", a different rule or approach is applied to the interval {0, 1/16}. I don't know why NIntegrate evaluates the function for x near and at the upper end point 1/16 with high precision. It could be estimating the derivative in order to estimate the error or check for a singularity/oscillatory behavior, maybe. $\endgroup$ – Michael E2 Mar 23 '15 at 13:25
  • $\begingroup$ There are strategies for integrating interpolating functions. Your use-case might be slow for reasons other than a few high-precision evaluations. $\endgroup$ – Michael E2 Mar 23 '15 at 13:31
  • $\begingroup$ Using N on x inside f changes the precision of the output, but the integration routine continues to request absurdly high-precision x. Changing the integration methods seem to have little effect - even the trapezoidal rule gives this behavior. Even the LocalAdaptive strategy exhibits the same problem. Unfortunately, my real integrand is a lot more than interpolating functions - it contains a mixture of theory and experiment, so there are a few Piecewises and both 1 and 2-D interpolations with Guassian fits for extrapolation - in short, it would be very difficult to handle it in a specific way. $\endgroup$ – Widjet Mar 24 '15 at 6:11
  • $\begingroup$ Why not f2[x_?NumericQ, y_?NumericQ] := f[N[x], N[y]]? This bumps down the precision before the function is evaluated and it should not take a lot of time. BTW, the absurdly high precision is $MachinePrecision + $MaxExtraPrecision. $\endgroup$ – Michael E2 Mar 24 '15 at 11:15

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.