Adaptive Simpson Integration#

template<typename MemorySpace = Kokkos::HostSpace>
class AdaptiveSimpson : public mpart::RecursiveQuadratureBase<Kokkos::HostSpace>#

Adaptive Simpson-rule integration based on applying a simple Simpson 1/3 rule recursively on subintervals.

Public Functions

inline AdaptiveSimpson(unsigned int maxSub, unsigned int fdim, double absTol, double relTol, QuadError::Type errorMetric, unsigned int minSub = 0)#
inline AdaptiveSimpson(unsigned int maxSub, unsigned int fdim, double *workspace, double absTol, double relTol, QuadError::Type errorMetric, unsigned int minSub = 0)#
inline void SetDim(unsigned int fdim)#
template<class FunctionType>
inline void Integrate(FunctionType const &f, double lb, double ub, double *res) const#

Approximates the integral \(\int_{x_L}^{x_U} f(x) dx\).

Parameters:
  • f[in] The integrand. Can return any type that overloads multiplication with a double and the += operator. doubles and Eigen::VectorXd are examples.

  • lb[in] The lower bound \(x_L\) in the integration.

  • lb[in] The upper bound \(x_U\) in the integration.

Template Parameters:

FunctionType – The type of the integrand. Must have an operator()(double x) function.

Returns:

An approximation of \(\int_{x_L}^{x_U} f(x) dx\). The return type will be the same as the type returned by the integrand function.

template<class FunctionType>
inline void Integrate(double *workspace, FunctionType const &f, double lb, double ub, double *res) const#

Public Static Functions

static inline unsigned int GetWorkspaceSize(unsigned int maxSub, unsigned int fdim)#