D01EAF computes approximations to the integrals of a vector of similar functions, each defined over the same multidimensional hyper-rectangular region. The routine uses an adaptive subdivision strategy, and also computes absolute error estimates.
Upon entry, unless MINCLS has been set to a value less than or equal to , D01EAF divides the integration region into a number of subregions with randomly selected volumes. Inside each subregion the integrals and their errors are estimated. The initial number of subregions is chosen to be as large as possible without using more than MINCLS calls to FUNSUB. The results are stored in a partially ordered list (a heap). The routine then proceeds in stages. At each stage the subregion with the largest error (measured using the maximum norm) is halved along the coordinate axis where the integrands have largest absolute fourth differences. The basic rule is applied to each half of this subregion and the results are stored in the list. The results from the two halves are used to update the global integral and error estimates (FINEST and ABSEST) and the routine continues unless where the norm is the maximum norm, or further subdivision would use more than MAXCLS calls to FUNSUB. If at some stage there is insufficient working storage to keep the results for the next subdivision, the routine switches to a less efficient mode; only if this mode of operation breaks down is insufficient storage reported.
Genz A C and Malik A A (1980) An adaptive algorithm for numerical integration over an N-dimensional rectangular region J. Comput. Appl. Math.6 295–302
van Dooren P and de Ridder L (1976) An adaptive algorithm for numerical integration over an N-dimensional cube J. Comput. Appl. Math.2 207–217
1: NDIM – INTEGERInput
On entry: , the number of dimensions of the integrals.
On entry: the upper limits of integration,
, for .
4: MINCLS – INTEGERInput/Output
On entry: must be set either to the minimum number of FUNSUB calls to be allowed, in which case or to a negative value. In this case, the routine continues the calculation started in a previous call with the same integrands and integration limits: no parameters other than MINCLS, MAXCLS, ABSREQ, RELREQ or IFAIL must be changed between the calls.
On exit: gives the number of FUNSUB calls actually used by D01EAF. For the continuation case ( on entry) this is the number of new FUNSUB calls on the current call to D01EAF.
5: MAXCLS – INTEGERInput
On entry: the maximum number of FUNSUB calls to be allowed. In the continuation case this is the number of new FUNSUB calls to be allowed.
6: NFUN – INTEGERInput
On entry: , the number of integrands.
7: FUNSUB – SUBROUTINE, supplied by the user.External Procedure
FUNSUB must evaluate the integrands at a given point.
On exit: specifies the estimated absolute accuracy of , for .
14: IFAIL – INTEGERInput/Output
On entry: IFAIL must be set to , . If you are unfamiliar with this parameter you should refer to Section 3.3 in the Essential Introduction for details.
For environments where it might be inappropriate to halt program execution when an error is detected, the value is recommended. If the output of error messages is undesirable, then the value is recommended. Otherwise, because for this routine the values of the output parameters may be useful even if on exit, the recommended value is . When the value is used it is essential to test the value of IFAIL on exit.
On exit: unless the routine detects an error or a warning has been flagged (see Section 6).
6 Error Indicators and Warnings
If on entry or , explanatory error messages are output on the current error message unit (as defined by X04AAF).
Errors or warnings detected by the routine:
MAXCLS was too small for D01EAF to obtain the required accuracy. The arrays FINEST and ABSEST respectively contain current estimates for the integrals and errors.
LENWRK is too small for the routine to continue. The arrays FINEST and ABSEST respectively contain current estimates for the integrals and errors.
On a continuation call, MAXCLS was set too small to make any progress. Increase MAXCLS before calling D01EAF again.
An absolute error estimate for each integrand is output in the array ABSEST. The routine exits with if
8 Further Comments
Usually the running time for D01EAF will be dominated by the time in FUNSUB, so the maximum time that could be used by D01EAF will be proportional to MAXCLS multiplied by the cost of a call to FUNSUB.
On a normal call, you should set on entry.
For some integrands, particularly those that are poorly behaved in a small part of the integration region, D01EAF may terminate prematurely with values of ABSEST that are significantly smaller than the actual absolute errors. This behaviour should be suspected if the returned value of MINCLS is small relative to the expected difficulty of the integrals. When this occurs D01EAF should be called again, but with an entry value of , (see specification of MAXCLS) and the results compared with those from the previous call.
If the routine is called with , the exact values of FINEST and ABSEST on return will depend (within statistical limits) on the sequence of random numbers generated internally within D01EAF by calls to G05SAF. Separate runs will produce identical answers unless the part of the program executed prior to calling D01EAF also calls (directly or indirectly) routines from Chapter G05, and, in addition, the series of such calls differs between runs.
Because of moderate instability in the application of the basic integration rule, approximately the last decimal digits may be inaccurate when using D01EAF for large values of .
This example computes
where , . The program is intended to show how to exploit the continuation facility provided with D01EAF: the routine exits with (printing an explanatory error message) and is re-entered with MAXCLS reset to a larger value. The program can be used with any values of NDIM and NFUN, except that the expression for must be changed if (see specification of MAXCLS).