GammaLib  2.0.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GIntegral Class Reference

GIntegral class interface definition. More...

#include <GIntegral.hpp>

Inheritance diagram for GIntegral:
GBase

Public Member Functions

 GIntegral (void)
 Void constructor. More...
 
 GIntegral (GFunction *kernel)
 Function kernel constructor. More...
 
 GIntegral (const GIntegral &integral)
 Copy constructor. More...
 
virtual ~GIntegral (void)
 Destructor. More...
 
GIntegraloperator= (const GIntegral &integral)
 Assignment operator. More...
 
void clear (void)
 Clear integral. More...
 
GIntegralclone (void) const
 Clone integral. More...
 
std::string classname (void) const
 Return class name. More...
 
void max_iter (const int &iter)
 Set maximum number of iterations. More...
 
const int & max_iter (void) const
 Return maximum number of iterations. More...
 
void fixed_iter (const int &iter)
 Set fixed number of iterations. More...
 
const int & fixed_iter (void) const
 Return fixed number of iterations. More...
 
void eps (const double &eps)
 Set relative precision. More...
 
const double & eps (void) const
 Get relative precision. More...
 
void silent (const bool &silent)
 Set silence flag. More...
 
const bool & silent (void) const
 Get silence flag. More...
 
const int & iter (void) const
 Return number of iterations. More...
 
const int & calls (void) const
 Get number of function calls. More...
 
const bool & is_valid (void) const
 Signal if integration result is valid. More...
 
const std::string & message (void) const
 Return integration status message. More...
 
void kernel (GFunction *kernel)
 Set function kernel. More...
 
const GFunctionkernel (void) const
 Get function kernel. More...
 
double romberg (std::vector< double > bounds, const int &order=5)
 Perform Romberg integration. More...
 
double romberg (const double &a, const double &b, const int &order=5)
 Perform Romberg integration. More...
 
double trapzd (const double &a, const double &b, const int &n=1, double result=0.0)
 Perform Trapezoidal integration. More...
 
double adaptive_simpson (const double &a, const double &b) const
 Adaptive Simpson's integration. More...
 
double gauss_kronrod (const double &a, const double &b) const
 Gauss-Kronrod integration. More...
 
std::string print (const GChatter &chatter=NORMAL) const
 Print integral information. More...
 
- Public Member Functions inherited from GBase
virtual ~GBase (void)
 Destructor. More...
 

Protected Member Functions

void init_members (void)
 Initialise class members. More...
 
void copy_members (const GIntegral &integral)
 Copy class members. More...
 
void free_members (void)
 Delete class members. More...
 
double polint (double *xa, double *ya, int n, double x, double *dy)
 Perform Polynomial interpolation. More...
 
double adaptive_simpson_aux (const double &a, const double &b, const double &eps, const double &S, const double &fa, const double &fb, const double &fc, const int &bottom) const
 Auxiliary function for adaptive Simpson's method. More...
 
double rescale_error (double err, const double &result_abs, const double &result_asc) const
 Rescale errors for Gauss-Kronrod integration. More...
 

Protected Attributes

GFunctionm_kernel
 Pointer to function kernel. More...
 
double m_eps
 Requested relative integration precision. More...
 
int m_max_iter
 Maximum number of iterations. More...
 
int m_fix_iter
 Fixed number of iterations. More...
 
bool m_silent
 Suppress integration warnings in console. More...
 
int m_iter
 Number of iterations used. More...
 
int m_calls
 Number of function calls used. More...
 
bool m_isvalid
 Integration result valid (true=yes) More...
 
bool m_has_abserr
 Has absolute integration error. More...
 
bool m_has_relerr
 Has relative integration error. More...
 
double m_abserr
 Absolute integration error. More...
 
double m_relerr
 Absolute integration error. More...
 
std::string m_message
 Status message (if result is invalid) More...
 

Detailed Description

GIntegral class interface definition.

This class allows to perform integration using various methods. The integrand is implemented by a derived class of GFunction.

Definition at line 46 of file GIntegral.hpp.

Constructor & Destructor Documentation

GIntegral::GIntegral ( void  )
explicit

Void constructor.

Definition at line 231 of file GIntegral.cpp.

References init_members().

Referenced by clone().

GIntegral::GIntegral ( GFunction kernel)
explicit

Function kernel constructor.

Parameters
[in]kernelPointer to function kernel.

The function kernel constructor assigns the function kernel pointer in constructing the object.

Definition at line 249 of file GIntegral.cpp.

References init_members(), kernel(), and m_kernel.

GIntegral::GIntegral ( const GIntegral integral)

Copy constructor.

Parameters
[in]integralIntegral.

Definition at line 267 of file GIntegral.cpp.

References copy_members(), and init_members().

GIntegral::~GIntegral ( void  )
virtual

Destructor.

Definition at line 283 of file GIntegral.cpp.

References free_members().

Member Function Documentation

double GIntegral::adaptive_simpson ( const double &  a,
const double &  b 
) const

Adaptive Simpson's integration.

Parameters
[in]aLeft integration boundary.
[in]bRight integration boundary.

Integrates the function using an adaptive Simpson's rule. The initial interval [a,b] is split into two sub-intervals [a,c] and [c,b] for which the integral is computed using

\[ \frac{b-a}{6} f(a) + 4f(c) + f(b) \]

where \(c=(a+b)/2\) is the mid-point of interval [a,b]. Each sub-interval is then recursively divided into sub-interval and the process is repeated. Dividing of sub-intervals is stopped when the difference between subsequent intervals falls below the relative tolerance specified by eps(). The maximum recursion depth is set by the max_iter() method.

I almost do not dare to confess, but the code has been taken from http://en.wikipedia.org/wiki/Adaptive_Simpson%27s_method It's really pretty simple ...

< Mid-point of interval [a,b]

< Length of interval [a,b]

Definition at line 708 of file GIntegral.cpp.

References abs(), adaptive_simpson_aux(), GFunction::eval(), m_calls, m_eps, m_has_abserr, m_has_relerr, m_isvalid, m_iter, m_kernel, m_max_iter, m_message, m_silent, gammalib::str(), and gammalib::warning().

double GIntegral::adaptive_simpson_aux ( const double &  a,
const double &  b,
const double &  eps,
const double &  S,
const double &  fa,
const double &  fb,
const double &  fc,
const int &  bottom 
) const
protected

Auxiliary function for adaptive Simpson's method.

Parameters
[in]aLeft integration boundary.
[in]bRight integration boundary.
[in]epsPrecision.
[in]SIntegral of last computation.
[in]faFunction value at left integration boundary.
[in]fbFunction value at right integration boundary.
[in]fcFunction value at mid-point of interval [a,b]
[in]bottomIteration counter (stop when 0)

Implements a recursive auxiliary method for the adative_simpson() integrator.

< Mid-point of interval [a,b]

< Length of interval [a,b]

< Mid-point of interval [a,c]

< Mid-point of interval [c,b]

Definition at line 1173 of file GIntegral.cpp.

References abs(), eps(), GFunction::eval(), m_calls, m_isvalid, m_iter, and m_kernel.

Referenced by adaptive_simpson().

const int & GIntegral::calls ( void  ) const
inline

Get number of function calls.

Returns
Number of function calls.

Definition at line 233 of file GIntegral.hpp.

References m_calls.

Referenced by print(), and romberg().

std::string GIntegral::classname ( void  ) const
inlinevirtual

Return class name.

Returns
String containing the class name ("GIntegral").

Implements GBase.

Definition at line 129 of file GIntegral.hpp.

void GIntegral::clear ( void  )
virtual

Clear integral.

Implements GBase.

Definition at line 335 of file GIntegral.cpp.

References free_members(), and init_members().

GIntegral * GIntegral::clone ( void  ) const
virtual

Clone integral.

Returns
Pointer to deep copy of integral.

Implements GBase.

Definition at line 353 of file GIntegral.cpp.

References GIntegral().

void GIntegral::copy_members ( const GIntegral integral)
protected

Copy class members.

Parameters
[in]integralIntegral.

Definition at line 1044 of file GIntegral.cpp.

References m_abserr, m_calls, m_eps, m_fix_iter, m_has_abserr, m_has_relerr, m_isvalid, m_iter, m_kernel, m_max_iter, m_message, m_relerr, and m_silent.

Referenced by GIntegral(), and operator=().

const double & GIntegral::eps ( void  ) const
inline

Get relative precision.

Returns
Relative precision.

Definition at line 221 of file GIntegral.hpp.

References m_eps.

Referenced by adaptive_simpson_aux(), eps(), and print().

void GIntegral::fixed_iter ( const int &  iter)
inline

Set fixed number of iterations.

Parameters
[in]iterFixed number of iterations.

If the fixed number of iterations is set, the integration algorithm will always performed the given number of iterations, irrespectively of the precision that is reached. This feature is relevant for computing numerical derivates from numerically integrated functions.

Definition at line 183 of file GIntegral.hpp.

References iter(), and m_fix_iter.

Referenced by GCTAModelAeffBackground::aeff_integral(), cta_psf_radial_kerns_delta::eval(), GCTAModelIrfBackground::npred_roi_kern_theta::eval(), GCTAModelAeffBackground::npred_roi_kern_theta::eval(), GCTAModelSpatial::npred_roi_kern_theta::eval(), cta_irf_radial_kern_rho::eval(), GResponse::irf_radial_kern_theta::eval(), GResponse::irf_elliptical_kern_theta::eval(), cta_nroi_radial_kern_rho::eval(), cta_irf_elliptical_kern_rho::eval(), cta_nroi_elliptical_kern_rho::eval(), cta_irf_diffuse_kern_theta::eval(), cta_nroi_diffuse_kern_theta::eval(), cta_psf_radial_kern_rho::eval(), cta_psf_radial_kern_delta::eval(), cta_psf_elliptical_kern_rho::eval(), cta_psf_diffuse_kern_delta::eval(), GCTAResponseIrf::irf_diffuse(), GResponse::irf_elliptical(), GCTAResponseIrf::irf_elliptical(), GResponse::irf_radial(), GCTAResponseIrf::irf_radial(), GCTAModelIrfBackground::npred(), GCTAModelSpatial::npred(), GObservation::npred_spec(), GCTAResponseIrf::npsf(), GCTAResponseIrf::nroi(), GCTAResponseIrf::nroi_diffuse(), GCTAResponseIrf::nroi_elliptical(), GCTAResponseIrf::nroi_radial(), GCTAResponseCube::psf_diffuse(), GCTAResponseCube::psf_elliptical(), GCTAResponseCube::psf_radial(), and GCTAAeffArf::remove_thetacut().

const int & GIntegral::fixed_iter ( void  ) const
inline

Return fixed number of iterations.

Returns
Fixed number of iterations.

Definition at line 196 of file GIntegral.hpp.

References m_fix_iter.

Referenced by print().

void GIntegral::free_members ( void  )
protected

Delete class members.

Definition at line 1071 of file GIntegral.cpp.

Referenced by clear(), operator=(), and ~GIntegral().

double GIntegral::gauss_kronrod ( const double &  a,
const double &  b 
) const
void GIntegral::init_members ( void  )
protected

Initialise class members.

Definition at line 1015 of file GIntegral.cpp.

References m_abserr, m_calls, m_eps, m_fix_iter, m_has_abserr, m_has_relerr, m_isvalid, m_iter, m_kernel, m_max_iter, m_message, m_relerr, and m_silent.

Referenced by clear(), GIntegral(), and operator=().

const bool & GIntegral::is_valid ( void  ) const
inline

Signal if integration result is valid.

Returns
True is integration result is valid.

Definition at line 297 of file GIntegral.hpp.

References m_isvalid.

Referenced by print().

const int & GIntegral::iter ( void  ) const
inline

Return number of iterations.

Returns
Number of iterations.

Definition at line 141 of file GIntegral.hpp.

References m_iter.

Referenced by fixed_iter(), max_iter(), and print().

void GIntegral::kernel ( GFunction kernel)
inline

Set function kernel.

Parameters
[in]kernelFunction kernel.

Sets the function kernel for which the integral should be determined.

Definition at line 272 of file GIntegral.hpp.

References kernel(), and m_kernel.

const GFunction * GIntegral::kernel ( void  ) const
inline

Get function kernel.

Returns
Function kernel.

Definition at line 285 of file GIntegral.hpp.

References m_kernel.

Referenced by GIntegral(), and kernel().

void GIntegral::max_iter ( const int &  iter)
inline

Set maximum number of iterations.

Parameters
[in]iterMaximum number of iterations.

Definition at line 153 of file GIntegral.hpp.

References iter(), and m_max_iter.

const int & GIntegral::max_iter ( void  ) const
inline

Return maximum number of iterations.

Returns
Maximum number of iterations.

Definition at line 166 of file GIntegral.hpp.

References m_max_iter.

Referenced by print(), and romberg().

const std::string & GIntegral::message ( void  ) const
inline

Return integration status message.

Returns
Integration status message.

Definition at line 309 of file GIntegral.hpp.

References m_message.

Referenced by print().

GIntegral & GIntegral::operator= ( const GIntegral integral)

Assignment operator.

Parameters
[in]integralIntegral.
Returns
Integral.

Definition at line 305 of file GIntegral.cpp.

References copy_members(), free_members(), and init_members().

double GIntegral::polint ( double *  xa,
double *  ya,
int  n,
double  x,
double *  dy 
)
protected

Perform Polynomial interpolation.

Parameters
[in]xaPointer to array of X values.
[in]yaPointer to array of Y values.
[in]nNumber of elements in arrays.
[in]xX value at which interpolations should be performed.
[out]dyError estimate for interpolated values.

Given arrays xa[1,..,n] and ya[1,..,n], and given a value x, this method returns a value y, and an error estimate dy. If P(x) is the polynomial of degree n-1, then the returned value y=P(x).

Todo:

Implement exceptions instead of screen dump.

Use std::vector for xa and ya and start at 0

Definition at line 1094 of file GIntegral.cpp.

References abs(), G_POLINT, m_isvalid, m_message, m_silent, gammalib::str(), and gammalib::warning().

Referenced by romberg().

std::string GIntegral::print ( const GChatter chatter = NORMAL) const
virtual

Print integral information.

Parameters
[in]chatterChattiness.
Returns
String containing integral information.

Implements GBase.

Definition at line 946 of file GIntegral.cpp.

References calls(), eps(), fixed_iter(), is_valid(), iter(), m_abserr, m_fix_iter, m_has_abserr, m_has_relerr, m_relerr, max_iter(), message(), gammalib::parformat(), SILENT, silent(), and gammalib::str().

Referenced by cta_irf_diffuse_kern_theta::eval(), and cta_nroi_diffuse_kern_theta::eval().

double GIntegral::rescale_error ( double  err,
const double &  result_abs,
const double &  result_asc 
) const
protected

Rescale errors for Gauss-Kronrod integration.

Parameters
[in]errError estimate.
[in]result_abs???.
[in]result_asc???.
Returns
Rescaled error estimate.

Definition at line 1236 of file GIntegral.cpp.

References abs(), and pow().

Referenced by gauss_kronrod().

double GIntegral::romberg ( std::vector< double >  bounds,
const int &  order = 5 
)

Perform Romberg integration.

Parameters
[in]boundsIntegration boundaries.
[in]orderIntegration order (default: 5)
Exceptions
GException::invalid_argumentIntegration order incompatible with number of iterations.

Returns the integral of the integrand, computed over a number of intervals [a0,a1], [a1,a2], ... that are given as an unordered vector by the bounds argument.

Integration is performed by Romberg's method of order 2*order, where

order=1 is equivalent to the trapezoidal rule,
order=2 is equivalent to Simpson's rule, and
order=3 is equivalent to Boole's rule.

The number of iterations is limited by m_max_iter. m_eps specifies the requested fractional accuracy. By default it is set to 1e-6.

Definition at line 381 of file GIntegral.cpp.

References calls(), and m_calls.

Referenced by GCTAModelAeffBackground::aeff_integral(), GModelSpectralMultiplicative::eflux(), GModelSpectralExponential::eflux(), GModelSpectralLogParabola::eflux(), GModelSpectralSmoothBrokenPlaw::eflux(), GModelSky::eflux(), cta_psf_radial_kerns_delta::eval(), GCTAModelIrfBackground::npred_roi_kern_theta::eval(), GCTAModelAeffBackground::npred_roi_kern_theta::eval(), GCTAModelSpatial::npred_roi_kern_theta::eval(), GModelSpatial::circle_int_kern_rho::eval(), cta_irf_radial_kern_rho::eval(), GResponse::irf_radial_kern_theta::eval(), GResponse::irf_elliptical_kern_theta::eval(), cta_nroi_radial_kern_rho::eval(), cta_irf_elliptical_kern_rho::eval(), cta_nroi_elliptical_kern_rho::eval(), cta_irf_diffuse_kern_theta::eval(), cta_nroi_diffuse_kern_theta::eval(), cta_psf_radial_kern_rho::eval(), cta_psf_radial_kern_delta::eval(), cta_psf_elliptical_kern_rho::eval(), cta_psf_diffuse_kern_delta::eval(), GModelSpectralMultiplicative::flux(), GModelSpectralExponential::flux(), GModelSpectralLogParabola::flux(), GModelSpatial::flux(), GModelSpectralSmoothBrokenPlaw::flux(), GModelSky::flux(), GLATPsfV3::integrate_psf(), GCTAResponseIrf::irf_diffuse(), GResponse::irf_elliptical(), GCTAResponseIrf::irf_elliptical(), GResponse::irf_radial(), GCTAResponseIrf::irf_radial(), GCOMIaq::location_smearing(), GLATPsfV3::normalize_psf(), GCTAEdisp2D::normalize_table(), GCTAModelIrfBackground::npred(), GCTAModelRadialAcceptance::npred(), GCTAModelSpatial::npred(), GObservation::npred(), GObservation::npred_spec(), GCTAResponseIrf::npsf(), GCTAResponseIrf::nroi(), GCTAResponseIrf::nroi_diffuse(), GCTAResponseIrf::nroi_elliptical(), GCTAResponseIrf::nroi_radial(), GCTAModelRadialGauss::omega(), GCTAModelRadialPolynom::omega(), GCTAModelRadialProfile::omega(), GLATPsfV1::psf(), GCTAResponseCube::psf_diffuse(), GCTAResponseCube::psf_elliptical(), GCTAResponseCube::psf_radial(), GCTAAeffArf::remove_thetacut(), GCTABackgroundPerfTable::solidangle(), and GCOMD2Response::update_response_vector().

double GIntegral::romberg ( const double &  a,
const double &  b,
const int &  order = 5 
)

Perform Romberg integration.

Parameters
[in]aLeft integration boundary.
[in]bRight integration boundary.
[in]orderIntegration order (default: 5)
Exceptions
GException::invalid_argumentIntegration order incompatible with number of iterations.

Returns the integral of the integrand from a to b. Integration is performed by Romberg's method of order 2*order, where

order=1 is equivalent to the trapezoidal rule,
order=2 is equivalent to Simpson's rule, and
order=3 is equivalent to Boole's rule.

The number of iterations is limited by m_max_iter. m_eps specifies the requested fractional accuracy. By default it is set to 1e-6.

Definition at line 426 of file GIntegral.cpp.

References abs(), G_ROMBERG, gammalib::is_infinite(), gammalib::is_notanumber(), m_abserr, m_calls, m_eps, m_fix_iter, m_has_abserr, m_has_relerr, m_isvalid, m_iter, m_max_iter, m_message, m_relerr, m_silent, max_iter(), polint(), gammalib::str(), trapzd(), and gammalib::warning().

void GIntegral::silent ( const bool &  silent)
inline

Set silence flag.

Parameters
[in]silentSilence flag.

Definition at line 245 of file GIntegral.hpp.

References m_silent, and silent().

Referenced by GModelSpatial::flux(), GCOMIaq::location_smearing(), and GCOMD2Response::update_response_vector().

const bool & GIntegral::silent ( void  ) const
inline

Get silence flag.

Returns
True is class is silent, false otherwise.

Definition at line 258 of file GIntegral.hpp.

References m_silent.

Referenced by print(), and silent().

double GIntegral::trapzd ( const double &  a,
const double &  b,
const int &  n = 1,
double  result = 0.0 
)

Perform Trapezoidal integration.

Parameters
[in]aLeft integration boundary.
[in]bRight integration boundary.
[in]nNumber of steps.
[in]resultResult from a previous trapezoidal integration step.
Returns
Integration results.
Exceptions
GException::invalid_valueFunction kernel not set.

The method performs a trapezoidal integration of the function kernel for the integration interval [a,b].

If n = 1 the integral is computed using

\[ \int_a^b f(x) \, dx = \frac{1}{2} (b-a) (f(a) + f(b)) \]

For n > 1 the integral is computed using

\[ \int_a^b f(x) \, dx = \frac{1}{2} \left[{\tt result} + \frac{b-a}{2^{n-1}} \sum_{i=0}^{2^{n-1}-1} f \left( a + (0.5+i) \frac{b-a}{2^{n-1}} \right) \right] \]

where \({\tt result}\) is the integration result from a previous call to the method with n = n - 1.

Definition at line 585 of file GIntegral.cpp.

References GFunction::eval(), G_TRAPZD, m_calls, m_isvalid, m_kernel, m_message, m_silent, gammalib::str(), sum(), and gammalib::warning().

Referenced by GCTAResponseIrf::npsf(), and romberg().

Member Data Documentation

double GIntegral::m_abserr
mutableprotected

Absolute integration error.

Definition at line 117 of file GIntegral.hpp.

Referenced by copy_members(), gauss_kronrod(), init_members(), print(), and romberg().

int GIntegral::m_calls
mutableprotected

Number of function calls used.

Definition at line 113 of file GIntegral.hpp.

Referenced by adaptive_simpson(), adaptive_simpson_aux(), calls(), copy_members(), gauss_kronrod(), init_members(), romberg(), and trapzd().

double GIntegral::m_eps
protected

Requested relative integration precision.

Definition at line 106 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), eps(), gauss_kronrod(), init_members(), and romberg().

int GIntegral::m_fix_iter
protected

Fixed number of iterations.

Definition at line 108 of file GIntegral.hpp.

Referenced by copy_members(), fixed_iter(), init_members(), print(), and romberg().

bool GIntegral::m_has_abserr
mutableprotected

Has absolute integration error.

Definition at line 115 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), gauss_kronrod(), init_members(), print(), and romberg().

bool GIntegral::m_has_relerr
mutableprotected

Has relative integration error.

Definition at line 116 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), gauss_kronrod(), init_members(), print(), and romberg().

bool GIntegral::m_isvalid
mutableprotected

Integration result valid (true=yes)

Definition at line 114 of file GIntegral.hpp.

Referenced by adaptive_simpson(), adaptive_simpson_aux(), copy_members(), gauss_kronrod(), init_members(), is_valid(), polint(), romberg(), and trapzd().

int GIntegral::m_iter
mutableprotected

Number of iterations used.

Definition at line 112 of file GIntegral.hpp.

Referenced by adaptive_simpson(), adaptive_simpson_aux(), copy_members(), gauss_kronrod(), init_members(), iter(), and romberg().

GFunction* GIntegral::m_kernel
protected

Pointer to function kernel.

Definition at line 105 of file GIntegral.hpp.

Referenced by adaptive_simpson(), adaptive_simpson_aux(), copy_members(), gauss_kronrod(), GIntegral(), init_members(), kernel(), and trapzd().

int GIntegral::m_max_iter
protected

Maximum number of iterations.

Definition at line 107 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), init_members(), max_iter(), and romberg().

std::string GIntegral::m_message
mutableprotected

Status message (if result is invalid)

Definition at line 119 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), gauss_kronrod(), init_members(), message(), polint(), romberg(), and trapzd().

double GIntegral::m_relerr
mutableprotected

Absolute integration error.

Definition at line 118 of file GIntegral.hpp.

Referenced by copy_members(), gauss_kronrod(), init_members(), print(), and romberg().

bool GIntegral::m_silent
protected

Suppress integration warnings in console.

Definition at line 109 of file GIntegral.hpp.

Referenced by adaptive_simpson(), copy_members(), gauss_kronrod(), init_members(), polint(), romberg(), silent(), and trapzd().


The documentation for this class was generated from the following files: