QAND
Integrates a function on a hyper-rectangle.
Required Arguments
F — User-supplied FUNCTION to be integrated. The form is
F(N, X), where
N – The dimension of the hyper-rectangle. (Input)
X – The independent variable of dimension N. (Input)
F – The value of the integrand at X. (Output)
F must be declared EXTERNAL in the calling program.
N — The dimension of the hyper-rectangle. (Input)
N must be less than or equal to 20.
A — Vector of length N. (Input)
Lower limits of integration.
B — Vector of length N. (Input)
Upper limits of integration.
RESULT — Estimate of the integral from A to B of F. (Output)
The integral of F is approximated over the N-dimensional hyper-rectangle A.LE.X.LE.B.
Optional Arguments
ERRABS — Absolute accuracy desired. (Input)
Default: ERRABS = 1.e-3 for single precision and 1.d-8 for double precision.
ERRREL — Relative accuracy desired. (Input)
Default: ERRREL = 1.e-3 for single precision and 1.d-8 for double precision.
MAXFCN — Approximate maximum number of function evaluations to be permitted. (Input)
MAXFCN cannot be greater than 256N or IMACH(5) if N is greater than 3.
Default: MAXFCN = 32**N.
ERREST — Estimate of the absolute value of the error. (Output)
FORTRAN 90 Interface
Generic: CALL QAND (F, N, A, B, RESULT [, …])
Specific: The specific interface names are S_QAND and D_QAND.
FORTRAN 77 Interface
Single: CALL QAND (F, N, A, B, ERRABS, ERRREL, MAXFCN, RESULT, ERREST)
Double: The double precision name is DQAND.
Description
The routine QAND approximates the n-dimensional iterated integral
with the approximation returned in RESULT. An estimate of the error is returned in ERREST. The approximation is achieved by iterated applications of product Gauss formulas. The integral is first estimated by a two-point tensor product formula in each direction. Then for i = 1, …, n the routine calculates a new estimate by doubling the number of points in the i-th direction, but halving the number immediately afterwards if the new estimate does not change appreciably. This process is repeated until either one complete sweep results in no increase in the number of sample points in any dimension, or the number of Gauss points in one direction exceeds 256, or the number of function evaluations needed to complete a sweep would exceed MAXFCN.
Comments
1. Informational errors
Type | Code | Description |
|---|
3 | 1 | MAXFCN was set greater than 256N. |
4 | 2 | The maximum number of function evaluations has been reached, and convergence has not been attained. |
2. If EXACT is the exact value, QAND attempts to find RESULT such that ABS(EXACT ‑ RESULT).LE.MAX(ERRABS, ERRREL * ABS(EXACT)). To specify only a relative error, set ERRABS to zero. Similarly, to specify only an absolute error, set ERRREL to zero.
Example
In this example, we approximate the integral of
on an expanding cube. The values of the error estimates are machine dependent. The exact integral over
USE QAND_INT
USE UMACH_INT
IMPLICIT NONE
INTEGER I, J, MAXFCN, N, NOUT
REAL A(3), B(3), CNST, ERRABS, ERREST, ERRREL, F, RESULT
EXTERNAL F
! Get output unit number
CALL UMACH (2, NOUT)
!
N = 3
MAXFCN = 100000
! Set error tolerances
ERRABS = 0.0001
ERRREL = 0.001
!
DO 20 I=1, 6
CNST = I/2.0
! Set limits of integration
! As CNST approaches infinity, the
! answer approaches PI**1.5
DO 10 J=1, 3
A(J) = -CNST
B(J) = CNST
10 CONTINUE
CALL QAND (F, N, A, B, RESULT, ERRABS, ERRREL, MAXFCN, ERREST)
WRITE (NOUT,99999) CNST, RESULT, ERREST
20 CONTINUE
99999 FORMAT (1X, 'For CNST = ', F4.1, ', result = ', F7.3, ' with ', &
'error estimate ', 1PE10.3)
END
!
REAL FUNCTION F (N, X)
INTEGER N
REAL X(N)
REAL EXP
INTRINSIC EXP
F = EXP(-(X(1)*X(1)+X(2)*X(2)+X(3)*X(3)))
RETURN
END
Output
For CNST = 0.5, result = 0.785 with error estimate 3.934E-06
For CNST = 1.0, result = 3.332 with error estimate 2.100E-03
For CNST = 1.5, result = 5.021 with error estimate 1.192E-05
For CNST = 2.0, result = 5.491 with error estimate 2.413E-04
For CNST = 2.5, result = 5.561 with error estimate 4.232E-03
For CNST = 3.0, result = 5.568 with error estimate 2.580E-04
