mcintegral {cooltools} | R Documentation |

## Monte Carlo and Quasi-Monte Carlo integration in any dimension

### Description

Numerical integration using a Monte Carlo (MC) or Quasi-Monte Carlo (QMC) algorithm, based on a Halton sequence. These algorithms are of low order (1/sqrt(n) for MC, log(n)/n for QMC in one dimension) compared to the typical orders of 1D deterministic integrators, such as those available in the `integrate`

function. The MC and QMC integrators are suitable to compute D-dimensional integrals with D>>1, since the order of most deterministic methods deteriorates exponentially with D, whereas the order of MC remains 1/sqrt(n), irrespective of D, and the order of QMC only deteriorates slowly with D as log(n)^D/n.

### Usage

```
mcintegral(f, a, b, n = 1e+05, qmc = FALSE, seed = NULL, warn = TRUE)
```

### Arguments

`f` |
scalar function of a D-vector to be integrated numerically; for fast performance, this function should be vectorized, such that it returns an N-element vector if it is given an N-by-D matrix as argument. An automatic warning is produced if the function is not vectorized in this manner. |

`a` |
D-vector with lower limit(s) of the integration |

`b` |
D-vector with upper limit(s) of the integration |

`n` |
approximate number of random evaluations. (The exact number is max(1,round(sqrt(n)))^2.) |

`qmc` |
logical flag. If false (default), pseudo-random numbers are used; if true, quasi-random numbers from a D-dimensional Halton sequence are used. |

`seed` |
optional seed for random number generator. Only used if |

`warn` |
logical flag. If true (default), a warning is produced if the function f is not vectorized. |

### Value

Returns a list of items:

`value` |
the best estimate of the integral. |

`error` |
an estimate of the statistical 1-sigma uncertainty. |

`iterations` |
exact number of evaluations (close to n). |

### Author(s)

Danail Obreschkow

### Examples

```
## Numerically integrate sin(x)
f = function(x) sin(x)
m = mcintegral(f,0,pi)
cat(sprintf('Integral = %.3f\u00B1%.3f (true value = 2)\n',m$value,m$sigma))
## Numerically compute the volume of a unit sphere
sphere = function(x) as.numeric(rowSums(x^2)<=1) # this is vectorized
vmc = mcintegral(sphere,rep(-1,3),rep(1,3),seed=1)
vqmc = mcintegral(sphere,rep(-1,3),rep(1,3),qmc=TRUE)
cat(sprintf('Volume of unit sphere = %.3f\u00B1%.3f (MC)\n',vmc$value,vmc$error))
cat(sprintf('Volume of unit sphere = %.3f\u00B1%.3f (QMC)\n',vqmc$value,vqmc$error))
cat(sprintf('Volume of unit sphere = %.3f (exact)\n',4*pi/3))
```

*cooltools*version 2.4 Index]