QMC
Integrates a function over a hyper-rectangle using a quasi-Monte Carlo method.
Required Arguments
FCN — User-supplied FUNCTION to be integrated. The form is FCN(X), where
X − The independent variable. (Input)
FCN − The value of the integrand at X. (Output)
FCN must be declared EXTERNAL in the calling program.
A — Vector containing lower limits of integration. (Input)
B — Vector containing upper limits of integration. (Input)
RESULT — The value of
is returned, where n is the dimension of X. If no value can be computed, then NaN is returned. (Output)
Optional Arguments
ERRABS — Absolute accuracy desired. (Input)
Default: 1.0e-2.
ERRREL — Relative accuracy desired. (Input)
Default: 1.0e-2.
ERREST — Estimate of the absolute value of the error. (Output)
MAXEVALS — Number of evaluations allowed. (Input)
Default: No limit.
BASE — The base of the Faure sequence. (Input)
Default: The smallest prime number greater than or equal to the number of dimensions (length of a and b).
SKIP — The number of points to be skipped at the beginning of the Faure sequence. (Input)
Default: 
, where 
and B is the largest representable integer.
FORTRAN 90 Interface
Generic: CALL QMC (FCN, A, B, RESULT [, …])
Specific: The specific interface names are S_QMC and D_QMC.
Description
Integration of functions over hyper rectangle by direct methods, such as QAND, is practical only for fairly low dimensional hypercubes. This is because the amount of work required increases exponentially as the dimension increases.
An alternative to direct methods is QMC, in which the integral is evaluated as the value of the function averaged over a sequence of randomly chosen points. Under mild assumptions on the function, this method will converge like
where k is the number of points at which the function is evaluated.
It is possible to improve on the performance of QMC by carefully choosing the points at which the function is to be evaluated. Randomly distributed points tend to be non-uniformly distributed. The alternative to a sequence of random points is a low-discrepancy sequence. A low-discrepancy sequence is one that is highly uniform.
This function is based on the low-discrepancy Faure sequence as computed by FAURE_NEXT, see the Fortran Stat Library, Chapter 18, “Random Number Generation”.
Example
This example evaluates the n-dimensional integral
with n=10.
use qmc_int
implicit none
integer, parameter :: ndim=10
real(kind(1d0)) :: a(ndim)
real(kind(1d0)) :: b(ndim)
real(kind(1d0)) :: result
integer :: I
external fcn
a = 0.d0
b = 1.d0
call qmc(fcn, a, b, result)
write (*,*) 'result = ', result
end
real(kind(1d0)) function fcn(x)
implicit none
real(kind(1d0)), dimension(:) :: x
integer :: i, j
real(kind(1d0)) :: prod, sum, sign
sign = -1.d0
sum = 0.d0
do i=1, size(x)
prod = 1.d0
prod = product(x(1:i))
sum = sum + (sign * prod)
sign = -sign
end do
fcn = sum
end function fcn
result = -0.3334789
