Implement a deterministic algorithm for multidimensional numerical
integration. Its algorithm uses one of several cubature rules in a globally
adaptive subdivision scheme. The subdivision algorithm is similar to
suave's.
cuhre(
f,
nComp = 1L,
lowerLimit,
upperLimit,
...,
relTol = 1e-05,
absTol = 1e-12,
minEval = 0L,
maxEval = 10^6,
flags = list(verbose = 0L, final = 1L, keep_state = 0L, level = 0L),
key = 0L,
nVec = 1L,
stateFile = NULL
)The function (integrand) to be integrated. For cuhre, it can be something as simple as a function of a single argument, say x.
The number of components of f, default 1, bears no relation to the dimension of the hypercube over which integration is performed.
The lower limit of integration, a vector for hypercubes.
The upper limit of integration, a vector for hypercubes.
All other arguments passed to the function f.
The maximum tolerance, default 1e-5.
the absolute tolerance, default 1e-12.
the minimum number of function evaluations required
The maximum number of function evaluations needed, default 10^6. Note that the actual number of function evaluations performed is only approximately guaranteed not to exceed this number.
flags governing the integration. The list here is exhaustive to keep the documentation and invocation uniform, but not all flags may be used for a particular method as noted below. List components:
verboseencodes the verbosity level, from 0 (default) to 3. Level 0 does not print any output, level 1 prints reasonable information on the progress of the integration, level 2 also echoes the input parameters, and level 3 further prints the subregion results.
finalwhen 0, all sets of samples collected on a subregion during the various iterations or phases contribute to the final result. When 1, only the last (largest) set of samples is used in the final result.
smoothApplies to Suave and Vegas only. When 0, apply additional smoothing to the importance function, this moderately improves convergence for many integrands. When 1, use the importance function without smoothing, this should be chosen if the integrand has sharp edges.
keep_statewhen nonzero, retain state file if
argument stateFile is non-null, else delete
stateFile if specified.
load_stateApplies to Vegas only. Reset the
integrator's state even if a state file is present, i.e. keep
only the grid. Together with keep_state this allows a
grid adapted by one integration to be used for another
integrand.
levelapplies only to Divonne, Suave
and Vegas. When 0, Mersenne Twister random numbers are
used. When nonzero Ranlux random numbers are used, except when
rngSeed is zero which forces use of Sobol quasi-random
numbers. Ranlux implements Marsaglia and Zaman's 24-bit RCARRY
algorithm with generation period p, i.e. for every 24 generated
numbers used, another p-24 are skipped. The luxury level for
the Ranlux generator may be encoded in level as follows:
gives very long period, passes the gap test but fails spectral test
passes all known tests, but theoretically still defective
any theoretically possible correlations have very small chance of being observed
highest possible luxury, all 24 bits chaotic
default to 3, values above 24 directly specify the period p.
level = 24.the quadrature rule key: key = 7, 9, 11, 13
selects the cubature rule of degree key. Note that the
degree-11 rule is available only in 3 dimensions, the degree-13
rule only in 2 dimensions. For other values, including the
default 0, the rule is the degree-13 rule in 2 dimensions, the
degree-11 rule in 3 dimensions, and the degree-9 rule
otherwise.
the number of vectorization points, default 1, but can be set to an integer > 1 for vectorization, for example, 1024 and the function f above needs to handle the vector of points appropriately. See vignette examples.
the name of an external file. Vegas can store its entire internal state (i.e. all the information to resume an interrupted integration) in an external file. The state file is updated after every iteration. If, on a subsequent invocation, Vegas finds a file of the specified name, it loads the internal state and continues from the point it left off. Needless to say, using an existing state file with a different integrand generally leads to wrong results. Once the integration finishes successfully, i.e. the prescribed accuracy is attained, the state file is removed. This feature is useful mainly to define ‘check-points’ in long-running integrations from which the calculation can be restarted.
A list with components:
the actual number of integrand evaluations needed
if zero, the desired accuracy was reached, if -1, dimension out of range, if 1, the accuracy goal was not met within the allowed maximum number of integrand evaluations.
vector of length nComp; the integral of
integrand over the hypercube
vector of
length nComp; the presumed absolute error of
integral
vector of length nComp; the
\(\chi^2\)-probability (not the
\(\chi^2\)-value itself!) that error is not a
reliable estimate of the true integration error.
See details in the documentation.
J. Berntsen, T. O. Espelid (1991) An adaptive algorithm for the approximate calculation of multiple integrals. ACM Transactions on Mathematical Software, 17(4), 437-451.
T. Hahn (2005) CUBA-a library for multidimensional numerical integration. Computer Physics Communications, 168, 78-95.
# NOT RUN {
integrand <- function(arg) {
x <- arg[1]
y <- arg[2]
z <- arg[3]
ff <- sin(x)*cos(y)*exp(z);
return(ff)
} # End integrand
NDIM <- 3
NCOMP <- 1
cuhre(f = integrand,
lowerLimit = rep(0, NDIM),
upperLimit = rep(1, NDIM),
relTol = 1e-3, absTol= 1e-12,
flags = list(verbose = 2, final = 0))
# }
Run the code above in your browser using DataLab