Mathematica integration timing
Marc Paterno
7/23/2020
Purpose of this document
This document looks at the Genz_1 (absolute value) integral in 5D as solved by Mathematica. This integral can not be solved analytically, so we try using Mathematica’s facilities for numerical integration with controlled error bounds.
Data
The Mathematica routine NIntegrate takes an argument PrecisionGoal that specifies the number of (decimal) significant digits of precision is wanted. This is analagous to (bit not identical to) our control parameter of epsrel.
We used the routine NIntegrate on the integrand \[ | \cos(4 v +5 w + 6 x +7 y + 8 z)|\] (note the lack of a normalization constant), varying the PrecisionGoal from 1 to 9. We used the Mathematica function Timing to mreasure the CPU time taken for each evalution.
In this dataframe:
- t is the time in seconds to perform the integration
- estimate is the estimated value of the integral
- prec is the number of decimal digits of precision requested
- delta is the difference between estimate and the previous estimate
- z is the absolute value of delta times 10^prec
Do the estimates seem to be converging as the error tolerance is tightened?
If each increase in prec does in fact generate one more significant (decimal) digit of precision, we should expect z to be about 1 for each row.
| t | estimate | prec | delta | z |
|---|---|---|---|---|
| 4.690997 | 0.642876941 | 1 | NA | NA |
| 5.253748 | 0.637403806 | 2 | -0.005473135 | 0.547313503 |
| 7.310951 | 0.639144292 | 3 | 0.001740486 | 1.740486015 |
| 11.781695 | 0.636071192 | 4 | -0.003073100 | 30.731000427 |
| 45.351825 | 0.636102157 | 5 | 0.000030965 | 3.096547832 |
| 161.971142 | 0.636710653 | 6 | 0.000608496 | 608.496238309 |
| 1197.121728 | 0.636579836 | 7 | -0.000130817 | 1308.172450611 |
| 2729.919346 | 0.636597441 | 8 | 0.000017605 | 1760.527287697 |
| 24150.790993 | 0.636595455 | 9 | -0.000001986 | 1985.969537044 |
We are not seeing wild variations, once we get past the very low precision calculations:
However, we note that the we are also not seing that the value is stable to the degree it should be, if the calculation were truly done to the number of significant digits claimed. Mathematica is clearly having trouble obtaining an error estimate of high quality.
How long does the calculation take?
Mathematica is using a single core to perform the integration, and working with native double precision floating-point numbers for this calculation. We have no indication from these results that using Mathematica’s extended-precision calculation ability is needed.
For prec>4, there seems to be a (very approximate) power law relation between the time taken and the specified precision.
Concentrating only on prec > 3, we can plot the power law fit. It is not an especially good fit, but probably has at least some predictive utility.
Given the rapid increase in the time required for the calculation, it would seem like PrecisionGoal=9 is the largest feasible value to use for this calculation; 10 digits would require several days’ runtime.