NAG Library Routine Document
D01FCF
1 Purpose
D01FCF attempts to evaluate a multidimensional integral (up to $15$ dimensions), with constant and finite limits, to a specified relative accuracy, using an adaptive subdivision strategy.
2 Specification
SUBROUTINE D01FCF ( 
NDIM, A, B, MINPTS, MAXPTS, FUNCTN, EPS, ACC, LENWRK, WRKSTR, FINVAL, IFAIL) 
INTEGER 
NDIM, MINPTS, MAXPTS, LENWRK, IFAIL 
REAL (KIND=nag_wp) 
A(NDIM), B(NDIM), FUNCTN, EPS, ACC, WRKSTR(LENWRK), FINVAL 
EXTERNAL 
FUNCTN 

3 Description
D01FCF returns an estimate of a multidimensional integral over a hyperrectangle (i.e., with constant limits), and also an estimate of the relative error. You set the relative accuracy required, return values for the integrand via a routine argument
FUNCTN, and also set the minimum and maximum acceptable number of calls to
FUNCTN (in
MINPTS and
MAXPTS).
The routine operates by repeated subdivision of the hyperrectangular region into smaller hyperrectangles. In each subregion, the integral is estimated using a seventhdegree rule, and an error estimate is obtained by comparison with a fifthdegree rule which uses a subset of the same points. The fourth differences of the integrand along each coordinate axis are evaluated, and the subregion is marked for possible future subdivision in half along that coordinate axis which has the largest absolute fourth difference.
If the estimated errors, totalled over the subregions, exceed the requested relative error (or if fewer than
MINPTS calls to
FUNCTN have been made), further subdivision is necessary, and is performed on the subregion with the largest estimated error, that subregion being halved along the appropriate coordinate axis.
The routine will fail if the requested relative error level has not been attained by the time
MAXPTS calls to
FUNCTN have been made; or, if the amount
LENWRK of working storage is insufficient. A formula for the recommended value of
LENWRK is given in
Section 5. If a smaller value is used, and is exhausted in the course of execution, the routine switches to a less efficient mode of operation; only if this mode also breaks down is insufficient storage reported.
D01FCF is based on the HALF subroutine developed by
van Dooren and de Ridder (1976). It uses a different basic rule, described in
Genz and Malik (1980).
4 References
Genz A C and Malik A A (1980) An adaptive algorithm for numerical integration over an Ndimensional 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 Ndimensional cube J. Comput. Appl. Math. 2 207–217
5 Arguments
 1: $\mathrm{NDIM}$ – INTEGERInput

On entry: $n$, the number of dimensions of the integral.
Constraint:
$2\le {\mathbf{NDIM}}\le 15$.
 2: $\mathrm{A}\left({\mathbf{NDIM}}\right)$ – REAL (KIND=nag_wp) arrayInput

On entry: the lower limits of integration,
${a}_{i}$, for $\mathit{i}=1,2,\dots ,n$.
 3: $\mathrm{B}\left({\mathbf{NDIM}}\right)$ – REAL (KIND=nag_wp) arrayInput

On entry: the upper limits of integration,
${b}_{i}$, for $\mathit{i}=1,2,\dots ,n$.
 4: $\mathrm{MINPTS}$ – INTEGERInput/Output

On entry: must be set to the minimum number of integrand evaluations to be allowed.
On exit: contains the actual number of integrand evaluations used by D01FCF.
 5: $\mathrm{MAXPTS}$ – INTEGERInput

On entry: the maximum number of integrand evaluations to be allowed.
Constraints:
 ${\mathbf{MAXPTS}}\ge {\mathbf{MINPTS}}$;
 ${\mathbf{MAXPTS}}\ge \alpha $, where $\alpha ={2}^{{\mathbf{NDIM}}}+2\times {{\mathbf{NDIM}}}^{2}+2\times {\mathbf{NDIM}}+1$.
 6: $\mathrm{FUNCTN}$ – REAL (KIND=nag_wp) FUNCTION, supplied by the user.External Procedure

FUNCTN must return the value of the integrand
$f$ at a given point.
The specification of
FUNCTN is:
FUNCTION FUNCTN ( 
NDIM, Z) 
REAL (KIND=nag_wp) FUNCTN 
INTEGER 
NDIM 
REAL (KIND=nag_wp) 
Z(NDIM) 

 1: $\mathrm{NDIM}$ – INTEGERInput

On entry: $n$, the number of dimensions of the integral.
 2: $\mathrm{Z}\left({\mathbf{NDIM}}\right)$ – REAL (KIND=nag_wp) arrayInput

On entry: the coordinates of the point at which the integrand $f$ must be evaluated.
FUNCTN must either be a module subprogram USEd by, or declared as EXTERNAL in, the (sub)program from which D01FCF is called. Arguments denoted as
Input must
not be changed by this procedure.
 7: $\mathrm{EPS}$ – REAL (KIND=nag_wp)Input

On entry: the relative error acceptable to you. When the solution is zero or very small relative accuracy may not be achievable but you may still set
EPS to a reasonable value and check for the error exit
${\mathbf{IFAIL}}={\mathbf{2}}$.
Constraint:
${\mathbf{EPS}}>0.0$.
 8: $\mathrm{ACC}$ – REAL (KIND=nag_wp)Output

On exit: the estimated relative error in
FINVAL.
 9: $\mathrm{LENWRK}$ – INTEGERInput

On entry: the dimension of the array
WRKSTR as declared in the (sub)program from which D01FCF is called.
Suggested value:
for maximum efficiency,
${\mathbf{LENWRK}}\ge \left({\mathbf{NDIM}}+2\right)\times \left(1+{\mathbf{MAXPTS}}/\alpha \right)$ (see argument
MAXPTS for
$\alpha $).
If
LENWRK is less than this, D01FCF will usually run less efficiently and may fail.
Constraint:
${\mathbf{LENWRK}}\ge 2\times {\mathbf{NDIM}}+4$.
 10: $\mathrm{WRKSTR}\left({\mathbf{LENWRK}}\right)$ – REAL (KIND=nag_wp) arrayWorkspace

 11: $\mathrm{FINVAL}$ – REAL (KIND=nag_wp)Output

On exit: the best estimate obtained for the integral.
 12: $\mathrm{IFAIL}$ – INTEGERInput/Output

On entry:
IFAIL must be set to
$0$,
$1\text{ or}1$. If you are unfamiliar with this argument you should refer to
Section 3.4 in How to Use the NAG Library and its Documentation for details.
For environments where it might be inappropriate to halt program execution when an error is detected, the value
$1\text{ or}1$ is recommended. If the output of error messages is undesirable, then the value
$1$ is recommended. Otherwise, because for this routine the values of the output arguments may be useful even if
${\mathbf{IFAIL}}\ne {\mathbf{0}}$ on exit, the recommended value is
$1$.
When the value $\mathbf{1}\text{ or}1$ is used it is essential to test the value of IFAIL on exit.
On exit:
${\mathbf{IFAIL}}={\mathbf{0}}$ unless the routine detects an error or a warning has been flagged (see
Section 6).
6 Error Indicators and Warnings
If on entry
${\mathbf{IFAIL}}=0$ or
$1$, explanatory error messages are output on the current error message unit (as defined by
X04AAF).
Note: D01FCF may return useful information for one or more of the following detected errors or warnings.
Errors or warnings detected by the routine:
 ${\mathbf{IFAIL}}=1$

On entry,  ${\mathbf{NDIM}}<2$, 
or  ${\mathbf{NDIM}}>15$, 
or  MAXPTS is too small, 
or  ${\mathbf{LENWRK}}<2\times {\mathbf{NDIM}}+4$, 
or  ${\mathbf{EPS}}\le 0.0$. 
 ${\mathbf{IFAIL}}=2$

MAXPTS was too small to obtain the required relative accuracy
EPS. On soft failure,
FINVAL and
ACC contain estimates of the integral and the relative error, but
ACC will be greater than
EPS.
 ${\mathbf{IFAIL}}=3$

LENWRK was too small. On soft failure,
FINVAL and
ACC contain estimates of the integral and the relative error, but
ACC will be greater than
EPS.
 ${\mathbf{IFAIL}}=99$
An unexpected error has been triggered by this routine. Please
contact
NAG.
See
Section 3.9 in How to Use the NAG Library and its Documentation for further information.
 ${\mathbf{IFAIL}}=399$
Your licence key may have expired or may not have been installed correctly.
See
Section 3.8 in How to Use the NAG Library and its Documentation for further information.
 ${\mathbf{IFAIL}}=999$
Dynamic memory allocation failed.
See
Section 3.7 in How to Use the NAG Library and its Documentation for further information.
7 Accuracy
A relative error estimate is output through the argument
ACC.
8 Parallelism and Performance
D01FCF is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
D01FCF makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the
X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this routine. Please also consult the
Users' Note for your implementation for any additional implementationspecific information.
Execution time will usually be dominated by the time taken to evaluate
FUNCTN, and hence the maximum time that could be taken will be proportional to
MAXPTS.
10 Example
This example estimates the integral
The accuracy requested is one part in $10000$.
10.1 Program Text
Program Text (d01fcfe.f90)
10.2 Program Data
None.
10.3 Program Results
Program Results (d01fcfe.r)