Integrals

Utilities for exploiting active subspaces when estimating integrals.

active_subspaces.integrals.av_integrate(avfun, avmap, N)

Approximate the integral of a function of active variables.

Parameters:
  • avfun (function) – a function of the active variables
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of points in the quadrature rule
Returns:

mu – an estimate of the integral

Return type:

float

Notes

This function is usually used when one has already constructed a response surface on the active variables and wants to estimate its integral.

active_subspaces.integrals.av_quadrature_rule(avmap, N)

Get a quadrature rule on the space of active variables.

Parameters:
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of quadrature nodes in the active variables
Returns:

  • Yp (ndarray) – quadrature nodes on the active variables
  • Yw (ndarray) – quadrature weights on the active variables

See also

integrals.quadrature_rule()

active_subspaces.integrals.integrate(fun, avmap, N, NMC=10)

Approximate the integral of a function of m variables.

Parameters:
  • fun (function) – an interface to the simulation that returns the quantity of interest given inputs as an 1-by-m ndarray
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of points in the quadrature rule
  • NMC (int, optional) – the number of points in the Monte Carlo estimates of the conditional expectation and conditional variance (default 10)
Returns:

  • mu (float) – an estimate of the integral of the function computed against the weight function on the simulation inputs
  • lb (float) – a central-limit-theorem 95% lower confidence from the Monte Carlo part of the integration
  • ub (float) – a central-limit-theorem 95% upper confidence from the Monte Carlo part of the integration

See also

integrals.quadrature_rule()

Notes

The CLT-based bounds lb and ub are likely poor estimators of the error. They only account for the variance from the Monte Carlo portion. They do not include any error from the integration rule on the active variables.

active_subspaces.integrals.interval_quadrature_rule(avmap, N, NX=10000)

Quadrature rule on a one-dimensional interval.

Quadrature when the dimension of the active subspace is 1 and the simulation parameter space is bounded.

Parameters:
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of quadrature nodes in the active variables
  • NX (int, optional) – the number of samples to use to estimate the quadrature weights (default 10000)
Returns:

  • Yp (ndarray) – quadrature nodes on the active variables
  • Yw (ndarray) – quadrature weights on the active variables

See also

integrals.quadrature_rule()

active_subspaces.integrals.quadrature_rule(avmap, N, NMC=10)

Get a quadrature rule on the space of simulation inputs.

Parameters:
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of quadrature nodes in the active variables
  • NMC (int, optional) – the number of samples in the simple Monte Carlo over the inactive variables (default 10)
Returns:

  • Xp (ndarray) – (N*NMC)-by-m matrix containing the quadrature nodes on the simulation input space
  • Xw (ndarray) – (N*NMC)-by-1 matrix containing the quadrature weights on the simulation input space
  • ind (ndarray) – array of indices identifies which rows of Xp correspond to the same fixed value of the active variables

See also

integrals.av_quadrature_rule()

Notes

This quadrature rule uses an integration rule on the active variables and simple Monte Carlo on the inactive variables.

If the simulation inputs are bounded, then the quadrature nodes on the active variables is constructed with a Delaunay triangulation of a maximin design. The weights are computed by sampling the original variables, mapping them to the active variables, and determining which triangle the active variables fall in. These samples are used to estimate quadrature weights. Note that when the dimension of the active subspace is one-dimensional, this reduces to operations on an interval.

If the simulation inputs are unbounded, the quadrature rule on the active variables is given by a tensor product Gauss-Hermite quadrature rule.

active_subspaces.integrals.zonotope_quadrature_rule(avmap, N, NX=10000)

Quadrature rule on a zonotope.

Quadrature when the dimension of the active subspace is greater than 1 and the simulation parameter space is bounded.

Parameters:
  • avmap (ActiveVariableMap) – a domains.ActiveVariableMap
  • N (int) – the number of quadrature nodes in the active variables
  • NX (int, optional) – the number of samples to use to estimate the quadrature weights (default 10000)
Returns:

  • Yp (ndarray) – quadrature nodes on the active variables
  • Yw (ndarray) – quadrature weights on the active variables

See also

integrals.quadrature_rule()