SimplicialCubature (version 1.1)

adaptIntegrateSimplex: Integrate a general function over a simplex

Description

Adaptive integration of a function f(x) of n variables over an m-dimensional simplex S, 1

Usage

adaptIntegrateSimplex(f, S, fDim = 1L, maxEvals = 10000L, absError = 0, tol = 1e-05, 
     integRule = 3L, partitionInfo = FALSE, ...)
integrate.vector.fn( intervals, fDim, f, maxEvals, absError, tol, partitionInfo=FALSE  )
original.coordinates( u, S  )

Value

A list containing:integralestimated value of the integral, it is a vector if fDim > 1estAbsErrorestimated absolute errorfunctionEvaluationscount of number of function evaluationsreturnCodeinteger status: returnCode=0 is a successful return; non-zero error values are described by next variablemessagetext message explaining returnCode; "OK" for normal returnsubsimplicesif partitionInfo=TRUE, this gives an array of subsimplices, see function adsimp for more details.subsimplicesIntegralif partitionInfo=TRUE, this array gives estimated values of each component of the integral on each subsimplex, see function adsimp for more details.subsimplicesAbsErrorif partitionInfo=TRUE, this array gives estimated values of the absolute error of each component of the integral on each subsimplex, see function adsimp for more details.subsimplicesVolumeif partitionInfo=TRUE, vector of m-dim. volumes of subsimplices; this is not d-dim. volume if m < n.

Details

If m=n, then an R translation of Alan Genz's function adsimp(...) is used to evaluate the n-dimensional integral. It works by adaptively splitting the region of integration into finer partitions, always splitting the simplex with the largest estimated error. If 1 < m < n, then the integral is evaluated by mapping the m-simplex in R^n to the canonical simplex in m-dimensional space, using function adsimp on that `full' m-dimensional integral, and correcting with the Jacobian of the transformation. If m=1, the built-in R function integrate (based on QUADPACK 1-dimensional adaptive quadrature) is used to evaluate the line integral. Since it does not provide access to the final subdivision, partitionInfo=TRUE in the univariate case returns the original partition information. So, if a fine parition is desired in the m=1 case, start with a fine partition.

References

See references to Genz and Cool (2003) in SimplicialCubature-package.

Examples

n <- 4
S <- CanonicalSimplex( n )
f1 <- function( x ) { x[1]^3 } 
adaptIntegrateSimplex( f1, S )   # correct answer 0.00119047619
str( adaptIntegrateSimplex( f1, S, partitionInfo=TRUE ) ) # same result, with more info returned

# test with vector valued integrand
f2 <- function( x ) { c(x[1]^3,x[3]^4) } 
adaptIntegrateSimplex( f2, S, fDim=2 )  # correct answer 0.00119047619 0.0005952380952

# test with vector valued integrand and extra arguments
f3 <- function( x, extra.arg ) { extra.arg*c(x[1]^3,x[3]^4) } # multiple of f2 above
adaptIntegrateSimplex( f3, S, fDim=2, extra.arg=100 )  # correct answer 0.119047619 0.05952380952

# integrate over lower dimensional simplices
adaptIntegrateSimplex( f1, UnitSimplexV(4) )  # answer =  0.01666667

f4 <- function(x) { 1 }
#  2-dim integral, exact answer area of unit simplex = sqrt(3)/2 = 0.8660254...
adaptIntegrateSimplex( f4, UnitSimplexV(3) ) 

# line integral over diamond in 2-dim, exact answer=arclength=4*sqrt(2)=5.656854...
S4 <- array( c( 1,0, 0,1,  0,1, -1,0,  -1,0, 0,-1,  0,-1, 1,0) , dim=c(2,2,4)  )
adaptIntegrateSimplex( f4, S4 ) 
adaptIntegrateSimplex( f4, S4, partitionInfo=TRUE )