Namespace lsst::afw::math::details

namespace details

Functions

template<class T>
T norm(const T &x)

template<class T>
T real(const T &x)

template<class T>
T Epsilon()

template<class T>
T MinRep()

template<class T>
T rescale_error(T err, T const &resabs, T const &resasc)

template<class UF>
bool intGKPNA(UF const &func, IntRegion<typename UF::result_type> &reg, typename UF::result_type const epsabs, typename UF::result_type const epsrel, std::map<typename UF::result_type, typename UF::result_type> *fxmap = 0)

Non-adaptive integration of the function f over the region ‘reg’.

Note

The algorithm computes first a Gaussian quadrature value then successive Kronrod/Patterson extensions to this result. The functions terminates when the difference between successive approximations (rescaled according to rescale_error) is less than either epsabs or epsrel * I, where I is the latest estimate of the integral. The order of the Gauss/Kronron/Patterson scheme is determined by which file is included above. Currently schemes starting with order 1 and order 10 are calculated. There seems to be little practical difference in the integration times using the two schemes, so I haven’t bothered to calculate any more.

template<class UF>
void intGKP(UF const &func, IntRegion<typename UF::result_type> &reg, typename UF::result_type const epsabs, typename UF::result_type const epsrel, std::map<typename UF::result_type, typename UF::result_type> *fxmap = 0)

An adaptive integration algorithm which computes the integral of f over the region reg.

The area and estimated error are returned as reg.Area() and reg.Err() If desired, *retx and *retf return std::vectors of x,f(x) respectively They only include the evaluations in the non-adaptive pass, so they do not give an accurate estimate of the number of function evaluations.

Note

First the non-adaptive GKP algorithm is tried. If that is not accurate enough (according to the absolute and relative accuracies, epsabs and epsrel), the region is split in half, and each new region is integrated. The routine continues by successively splitting the subregion which gave the largest absolute error until the integral converges.

template<class UF>
AuxFunc1<UF> Aux1(UF uf)

Auxiliary function 1

template<class UF>
AuxFunc2<UF> Aux2(UF uf)

Auxiliary function 2

template<class BF, class Tp>
binder2_1<BF> bind21(const BF &oper, const Tp &x)

template<class TF, class Tp>
binder3_1<TF> bind31(const TF &oper, const Tp &x)

int gkp_n(int level)

template<class T>
const std::vector<T> &gkp_x(int level)

template<class T>
const std::vector<T> &gkp_wa(int level)

template<class T>
const std::vector<T> &gkp_wb(int level)

Variables

const int NGKPLEVELS = 5

template<class UF>
struct AuxFunc1 : public std::unary_function<UF::argument_type, UF::result_type>

#include <Integrate.h>

Auxiliary struct 1

template<class T>
struct ConstantReg1 : public std::unary_function<T, IntRegion<T>>

#include <Integrate.h>

Helpers for constant regions for int2d, int3d:

template<typename BinaryFunctionT>
class FunctionWrapper : public std::unary_function<BinaryFunctionT::second_argument_type, BinaryFunctionT::result_type>

#include <Integrate.h>

Wrap an integrand in a call to a 1D integrator: romberg()

When romberg2D() is called, it wraps the integrand it was given in a FunctionWrapper functor. This wrapper calls romberg() on the integrand to get a 1D (along the x-coord, for constant y) result . romberg2D() then calls romberg() with the FunctionWrapper functor as an integrand.