This class allows the numerical evaluation of the following functions: More...

#include <NumericalDefiniteIntegral.h>

Inheritance diagram for Genfun::GaudiMathImplementation::NumericalDefiniteIntegral:

Collaboration diagram for Genfun::GaudiMathImplementation::NumericalDefiniteIntegral:

Classes
struct	_Function

struct	_Workspace

Public Types
typedef std::vector< double >	Points
	typedef for vector of singular points

Public Member Functions
	NumericalDefiniteIntegral (const AbsFunction &function, const size_t index, const double a, const double b, const GaudiMath::Integration::Type type=GaudiMath::Integration::Adaptive, const GaudiMath::Integration::KronrodRule rule=GaudiMath::Integration::Default, const double epsabs=1.e-10, const double epsrel=1.e-7, const size_t size=1000)
	From CLHEP/GenericFunctions.

	NumericalDefiniteIntegral (const AbsFunction &function, const size_t index, const double a, const double b, const Points &points, const double epsabs=1e-9, const double epsrel=1.e-6, const size_t size=1000)
	Standard constructor The function created with this constructor compute the following integral:

	NumericalDefiniteIntegral (const AbsFunction &function, const size_t index, const double a, const GaudiMath::Integration::Inf b=GaudiMath::Integration::Infinity, const double epsabs=1e-9, const double epsrel=1.e-6, const size_t size=1000)
	Standard constructor The function created with this constructor compute the following integral:

	NumericalDefiniteIntegral (const AbsFunction &function, const size_t index, const GaudiMath::Integration::Inf a, const double b, const double epsabs=1e-9, const double epsrel=1.e-6, const size_t size=1000)
	Standard constructor The function created with this constructor compute the following integral:

	NumericalDefiniteIntegral (const AbsFunction &function, const size_t index, const float epsabs=1e-9, const float epsrel=1.e-6, const size_t size=1000)
	Standard constructor The function created with this constructor compute the following integral:

	NumericalDefiniteIntegral (const NumericalDefiniteIntegral &)
	copy constructor

virtual	~NumericalDefiniteIntegral ()
	destructor

virtual unsigned int	dimensionality () const
	dimensionality of the problem

virtual double	operator() (double argument) const
	Function value.

virtual double	operator() (const Argument &argument) const
	Function value.

virtual bool	hasAnalyticDerivative () const
	Does this function have an analytic derivative?

virtual Genfun::Derivative	partial (unsigned int index) const
	Derivatives.

const AbsFunction &	function () const
	accessor to the function itself

double	a () const
	integration limit

double	b () const

const Points &	points () const
	known singularities

double	epsabs () const
	absolute precision

double	epsrel () const
	relatiove precision

double	result () const
	previous result

double	error () const
	evaluate of previous error

size_t	size () const

GaudiMath::Integration::Type	type () const
	integration type

GaudiMath::Integration::Category	category () const
	integration category

GaudiMath::Integration::KronrodRule	rule () const
	integration rule

Protected Member Functions
double	QAGI (_Function *fun) const

double	QAGP (_Function *fun) const

double	QNG (_Function *fun) const

double	QAG (_Function *fun) const

double	QAGS (_Function *fun) const

_Workspace *	allocate () const
	allocate the integration workspace

_Workspace *	ws () const

StatusCode	Exception (const std::string &message, const StatusCode &sc=StatusCode::FAILURE) const

Private Member Functions
	NumericalDefiniteIntegral ()

NumericalDefiniteIntegral &	operator= (const NumericalDefiniteIntegral &)

Private Attributes
const AbsFunction *	m_function

size_t	m_DIM

size_t	m_index

double	m_a

double	m_b

bool	m_ia

bool	m_ib

GaudiMath::Integration::Type	m_type

GaudiMath::Integration::Category	m_category

GaudiMath::Integration::KronrodRule	m_rule

Points	m_points

double *	m_pdata

double	m_epsabs

double	m_epsrel

double	m_result

double	m_error

size_t	m_size

_Workspace *	m_ws

Argument	m_argument

Argument	m_argF

Detailed Description

This class allows the numerical evaluation of the following functions:

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{b} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{+\infty} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{-\infty}^{b} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{-\infty}^{+\infty} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

Author: Vanya BELYAEV Ivan..nosp@m.Bely.nosp@m.aev@i.nosp@m.tep..nosp@m.ru

Date: 2003-08-31

Definition at line 72 of file NumericalDefiniteIntegral.h.

Member Typedef Documentation

typedef std::vector<double> Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::Points

typedef for vector of singular points

Definition at line 76 of file NumericalDefiniteIntegral.h.

Constructor & Destructor Documentation

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral	(	const AbsFunction &	function,
		const size_t	index,
		const double	a,
		const double	b,
		const GaudiMath::Integration::Type	type = `GaudiMath::Integration::Adaptive`,
		const GaudiMath::Integration::KronrodRule	rule = `GaudiMath::Integration::Default`,
		const double	epsabs = `1.e-10`,
		const double	epsrel = `1.e-7`,
		const size_t	size = `1000`
	)

From CLHEP/GenericFunctions.

from CLHEP/GenericFunctions

Standard constructor The function created with this constructor compute the following integral:

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{b} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

If function contains singularities, the type = Type::AdaptiveSingular need to be used

For faster integration of smooth function non-adaptive integration can be used: type = Type::NonAdaptive need to be used

For adaptive integration type = Type>>Adaptive one can specify the order of Gauss-Kronrad integration rule rule = KronrodRule::Gauss15 The higher-order rules give better accuracy for smooth functions, while lower-order rules save the time when the function contains local difficulties, such as discontinuites.

The GSL routine gsl_integration_qng is used for type = Type:NonAdaptive :

The GSL routine gsl_integration_qag is used for type = Type:Adaptive :

The GSL routine gsl_integration_qags is used for type = Type:AdaptiveSingular :

Parameters

function	funtion to be integrated
index	variable index
a	lower intgeration limit
b	high integration limit
type	integration type
rule	integration rule (for adaptive integration)
epsabs	required absolute precision
epsrel	required relative precision
size	maximal number of bisections for adaptive integration

Definition at line 110 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction ()
      , m_function  ( function.clone          ()    )
      , m_DIM       ( 0                             )
      , m_index     ( index                         )
      , m_a         ( a      )
      , m_b         ( b      )
      , m_ia        ( false  )
      , m_ib        ( false  )
      , m_type      ( type   )
      , m_category  ( GaudiMath::Integration::Finite )
      , m_rule      ( rule       )
      , m_points    (            )
      , m_pdata     ( 0          )
      , m_epsabs    ( epsabs     )
      , m_epsrel    ( epsrel     )
      , m_result    ( GSL_NEGINF )
      , m_error     ( GSL_POSINF )
      , m_size      ( size       )
      , m_ws        ()
      , m_argument  ()
      , m_argF      ()
    {
      if ( GaudiMath::Integration::Fixed == m_rule )
        { m_rule = GaudiMath::Integration::Default ; }
      if ( function.dimensionality() < 2          )
        { Exception("::constructor: invalid dimensionality ") ; }
      if ( m_index >= function.dimensionality()   )
        { Exception("::constructor: invalid variable index") ; }
      m_DIM = function.dimensionality() - 1 ;
      m_argument = Argument( m_DIM      ) ;
      m_argF     = Argument( m_DIM + 1  ) ;
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral	(	const AbsFunction &	function,
		const size_t	index,
		const double	a,
		const double	b,
		const Points &	points,
		const double	epsabs = `1e-9`,
		const double	epsrel = `1.e-6`,
		const size_t	size = `1000`
	)

Standard constructor The function created with this constructor compute the following integral:

standard constructor

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{b} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

where function have known singulatities

The GSL routine gsl_integration_qagp is used for integration

Parameters

function	funtion to be integrated
index	variable index
a	lower intgeration limit
b	high integration limit
pointns	vector of know singular points
epsabs	required absolute precision
epsrel	required relative precision
size	maximal number of bisections for adaptive integration
function	the base function
index	the variable index
a	integration limit
b	integration limit
points	list of known function singularities
epsabs	absolute precision for integration
epsrel	relative precision for integration

Definition at line 165 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction ()
      , m_function  ( function.clone()             )
      , m_DIM       ( 0                            )
      , m_index     ( index                        )
      , m_a         ( a                            )
      , m_b         ( b                            )
      , m_ia        ( false  )
      , m_ib        ( false  )
      , m_type      ( GaudiMath::Integration::    Other )
      , m_category  ( GaudiMath::Integration:: Singular )
      , m_rule      ( GaudiMath::Integration::    Fixed )
      , m_points    ( points  )
      , m_pdata     ( 0       )
      , m_epsabs    ( epsabs  )
      , m_epsrel    ( epsrel  )
      //
      , m_result    ( GSL_NEGINF                          )
      , m_error     ( GSL_POSINF                          )
      //
      , m_size      ( size                                )
      , m_ws        ( 0                                   )
      , m_argument  ()
      , m_argF      ()
    {
      if ( function.dimensionality() < 2          )
        { Exception("::constructor: invalid dimensionality ") ; }
      if ( m_index >= function.dimensionality()   )
        { Exception("::constructor: invalid variable index") ; }
      m_DIM = function.dimensionality() - 1 ;
      m_argument = Argument( m_DIM      ) ;
      m_argF     = Argument( m_DIM + 1  ) ;
      const double l1 = std::min ( a , b ) ;
      const double l2 = std::max ( a , b ) ;
      m_points.push_back ( l1 ) ;
      m_points.push_back ( l2 ) ;
      std::sort ( m_points.begin() , m_points.end() ) ;
      m_points.erase( std::unique( m_points.begin () ,
                                   m_points.end   () ) ,
                      m_points.end() );
      Points::iterator lower =
        std::lower_bound ( m_points.begin () , m_points.end () , l1 ) ;
      m_points.erase     ( m_points.begin () , lower                ) ;
      Points::iterator upper =
        std::upper_bound ( m_points.begin () , m_points.end () , l2 ) ;
      m_points.erase     ( upper             , m_points.end ()      ) ;
      m_pdata = new double[ m_points.size() ] ;
      std::copy( m_points.begin() , m_points.end() , m_pdata );
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral	(	const AbsFunction &	function,
		const size_t	index,
		const double	a,
		const GaudiMath::Integration::Inf	b = `GaudiMath::Integration::Infinity`,
		const double	epsabs = `1e-9`,
		const double	epsrel = `1.e-6`,
		const size_t	size = `1000`
	)

Standard constructor The function created with this constructor compute the following integral:

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{+\infty} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

The GSL routine gsl_integration_qagiu is used for integration

Parameters

function	funtion to be integrated
index	variable index
a	lower intgeration limit
b	indication for infinity
epsabs	required absolute precision
epsrel	required relative precision
size	maximal number of bisections for adaptive integration

Definition at line 252 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction()
      , m_function  ( function.clone()             )
      , m_DIM       ( 0                            )
      , m_index     ( index                        )
      , m_a         ( a                            )
      , m_b         ( GSL_POSINF                   )
      , m_ia        ( false  )
      , m_ib        ( true   )
      , m_type      ( GaudiMath::Integration::    Other )
      , m_category  ( GaudiMath::Integration:: Infinite )
      , m_rule      ( GaudiMath::Integration::    Fixed )
      , m_points    (         )
      , m_pdata     ( 0       )
      , m_epsabs    ( epsabs  )
      , m_epsrel    ( epsrel  )
      //
      , m_result    ( GSL_NEGINF                          )
      , m_error     ( GSL_POSINF                          )
      //
      , m_size      ( size                                )
      , m_ws        ( 0                                   )
      , m_argument  ()
      , m_argF      ()
    {
      if ( function.dimensionality() < 2          )
        { Exception("::constructor: invalid dimensionality ") ; }
      if ( m_index >= function.dimensionality()   )
        { Exception("::constructor: invalid variable index") ; }
      m_DIM = function.dimensionality() - 1 ;
      m_argument = Argument( m_DIM      ) ;
      m_argF     = Argument( m_DIM + 1  ) ;
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral	(	const AbsFunction &	function,
		const size_t	index,
		const GaudiMath::Integration::Inf	a,
		const double	b,
		const double	epsabs = `1e-9`,
		const double	epsrel = `1.e-6`,
		const size_t	size = `1000`
	)

Standard constructor The function created with this constructor compute the following integral:

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{-\infty}^{b} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

The GSL routine gsl_integration_qagiu is used for integration

Parameters

function	funtion to be integrated
index	variable index
a	lower intgeration limit
b	indication for infinity
epsabs	required absolute precision
epsrel	required relative precision
size	maximal number of bisections for adaptive integration

Definition at line 321 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction()
      , m_function  ( function.clone()             )
      , m_DIM       ( 0                            )
      , m_index     ( index                        )
      , m_a         ( GSL_NEGINF                   )
      , m_b         ( b                            )
      , m_ia        ( true   )
      , m_ib        ( false  )
      , m_type      ( GaudiMath::Integration::    Other )
      , m_category  ( GaudiMath::Integration:: Infinite )
      , m_rule      ( GaudiMath::Integration::    Fixed )
      , m_points    (         )
      , m_pdata     ( 0       )
      , m_epsabs    ( epsabs  )
      , m_epsrel    ( epsrel  )
      //
      , m_result    ( GSL_NEGINF                          )
      , m_error     ( GSL_POSINF                          )
      //
      , m_size      ( size                                )
      , m_ws        ( 0                                   )
      , m_argument  ()
      , m_argF      ()
    {
      if ( function.dimensionality() < 2          )
        { Exception("::constructor: invalid dimensionality ") ; }
      if ( m_index >= function.dimensionality()   )
        { Exception("::constructor: invalid variable index") ; }
      m_DIM = function.dimensionality() - 1 ;
      m_argument = Argument( m_DIM      ) ;
      m_argF     = Argument( m_DIM + 1  ) ;
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral	(	const AbsFunction &	function,
		const size_t	index,
		const float	epsabs = `1e-9`,
		const float	epsrel = `1.e-6`,
		const size_t	size = `1000`
	)

Standard constructor The function created with this constructor compute the following integral:

${\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_{i+1}, \dots , x_n \right) = \int\limits_{-\infty}^{+\infty} f \left(x_1, \dots , x_{i-1}, \hat{x}_i , x_{i+1}, \dots , x_n \right) \, {\mathrm{d}} \hat{x}_i$

The GSL routine gsl_integration_qagi is used for integration

Parameters

function	funtion to be integrated
index	variable index
epsabs	required absolute precision
epsrel	required relative precision
size	maximal number of bisections for adaptive integration

Definition at line 386 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction()
      , m_function  ( function.clone()             )
      , m_DIM       ( 0      )
      , m_index     ( index                        )
      , m_a         ( GSL_NEGINF                   )
      , m_b         ( GSL_POSINF                   )
      , m_ia        ( true   )
      , m_ib        ( true   )
      , m_type      ( GaudiMath::Integration::    Other )
      , m_category  ( GaudiMath::Integration:: Infinite )
      , m_rule      ( GaudiMath::Integration::    Fixed )
      , m_points    (         )
      , m_pdata     ( 0       )
      , m_epsabs    ( epsabs  )
      , m_epsrel    ( epsrel  )
      //
      , m_result    ( GSL_NEGINF                          )
      , m_error     ( GSL_POSINF                          )
      //
      , m_size      ( size                                )
      , m_ws        ( 0                                   )
      , m_argument  ()
      , m_argF      ()
    {
      if ( function.dimensionality() < 2          )
        { Exception("::constructor: invalid dimensionality ") ; }
      if ( m_index >= function.dimensionality()   )
        { Exception("::constructor: invalid variable index") ; }
      m_DIM = function.dimensionality() - 1 ;
      m_argument = Argument( m_DIM      ) ;
      m_argF     = Argument( m_DIM + 1  ) ;
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral ( const NumericalDefiniteIntegral & right )

copy constructor

Definition at line 428 of file NumericalDefiniteIntegral.cpp.

      : AbsFunction ()
      , m_function  ( right.m_function->clone() )
      , m_DIM       ( right.m_DIM      )
      , m_index     ( right.m_index    )
      , m_a         ( right.m_a        )
      , m_b         ( right.m_b        )
      , m_ia        ( right.m_ia       )
      , m_ib        ( right.m_ib       )
      , m_type      ( right.m_type     )
      , m_category  ( right.m_category )
      , m_rule      ( right.m_rule     )
      , m_points    ( right.m_points   )
      , m_pdata     ( 0                )
      , m_epsabs    ( right.m_epsabs   )
      , m_epsrel    ( right.m_epsrel   )
      , m_result    ( GSL_NEGINF       )
      , m_error     ( GSL_POSINF       )
      , m_size      ( right.m_size     )
      , m_ws        ( 0                )
      , m_argument  ( right.m_argument )
      , m_argF      ( right.m_argF     )
    {
      if( 0 != right.m_pdata )
        {
          m_pdata = new double[m_points.size()] ;
          std::copy( m_points.begin() , m_points.end() , m_pdata );
        }
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::~NumericalDefiniteIntegral ( )

virtual

destructor

Definition at line 458 of file NumericalDefiniteIntegral.cpp.

    {
      if( 0 != m_function ) { delete   m_function ; m_function = 0 ; }
      if( 0 != m_pdata    ) { delete[] m_pdata    ; m_pdata    = 0 ; }
    }

Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::NumericalDefiniteIntegral ( )

private

Member Function Documentation

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::a ( ) const

inline

integration limit

Definition at line 322 of file NumericalDefiniteIntegral.h.

{ return m_a ; }

NumericalDefiniteIntegral::_Workspace * Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::allocate ( ) const

protected

allocate the integration workspace

Definition at line 574 of file NumericalDefiniteIntegral.cpp.

    {
      if ( 0 != m_ws ) { return m_ws; }
      gsl_integration_workspace* aux =
        gsl_integration_workspace_alloc( size () );
      if ( 0 == aux ) { Exception ( "allocate()::invalid workspace" ) ; };
      m_ws = new _Workspace() ;
      m_ws->ws = aux ;
      return m_ws ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::b ( ) const

inline

Definition at line 323 of file NumericalDefiniteIntegral.h.

{ return m_b ; }

GaudiMath::Integration::Category Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::category ( ) const

inline

integration category

Definition at line 344 of file NumericalDefiniteIntegral.h.

{ return m_category ; }

virtual unsigned int Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::dimensionality ( ) const

inlinevirtual

dimensionality of the problem

Definition at line 304 of file NumericalDefiniteIntegral.h.

{ return m_DIM ; }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::epsabs ( ) const

inline

absolute precision

Definition at line 327 of file NumericalDefiniteIntegral.h.

{ return m_epsabs ; }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::epsrel ( ) const

inline

relatiove precision

Definition at line 329 of file NumericalDefiniteIntegral.h.

{ return m_epsrel ; }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::error ( ) const

inline

evaluate of previous error

Definition at line 334 of file NumericalDefiniteIntegral.h.

{ return m_error ; }

StatusCode Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::Exception	(	const std::string &	message,
		const StatusCode &	sc = `StatusCode::FAILURE`
	)		const

protected

Definition at line 468 of file NumericalDefiniteIntegral.cpp.

    {
      throw GaudiException( "NumericalDefiniteIntegral::" + message ,
                            "*GaudiMath*" , sc ) ;
      return sc ;
    }

const AbsFunction& Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::function ( ) const

inline

accessor to the function itself

Definition at line 320 of file NumericalDefiniteIntegral.h.

{ return *m_function ; }

virtual bool Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::hasAnalyticDerivative ( ) const

inlinevirtual

Does this function have an analytic derivative?

Definition at line 312 of file NumericalDefiniteIntegral.h.

{ return true ;}

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::operator() ( double argument ) const

virtual

Function value.

evaluate the function

Definition at line 481 of file NumericalDefiniteIntegral.cpp.

    {
      // reset the result and the error
      m_result = GSL_NEGINF ;
      m_error  = GSL_POSINF ;
      // check the argument
      if( 1 != m_DIM ) { Exception ( "operator(): invalid argument size " ) ; };
      m_argument[0] = argument ;
      return (*this) ( m_argument );
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::operator() ( const Argument & argument ) const

virtual

Function value.

evaluate the function

Definition at line 514 of file NumericalDefiniteIntegral.cpp.

    {
      // reset the result and the error
      m_result = GSL_NEGINF ;
      m_error  = GSL_POSINF ;
      // check the argument
      if( argument.dimension() != m_DIM )
        { Exception ( "operator(): invalid argument size " ) ; };
      // copy the argument
      for( size_t i  = 0 ; i < m_DIM ; ++i )
        {
          m_argument [i] = argument [i] ;
          const size_t j =  i < m_index ? i : i + 1 ;
          m_argF     [j] = argument [i] ;
        };
      // create the helper object
      GSL_Helper helper( *m_function , m_argF , m_index );
      // use GSL to evaluate the numerical derivative
      gsl_function F ;
      F.function = &GSL_Adaptor ;
      F.params   = &helper                 ;
      _Function F1    ;
      F1.fn      = &F ;
      if        (  GaudiMath::Integration::Infinite         == category () )
        { return   QAGI ( &F1 ) ; }                                // RETURN
      if ( m_a == m_b )
        {
          m_result = 0    ; m_error  = 0    ;                      // EXACT
          return m_result ;                                       // RETURN
        }
      if        (  GaudiMath::Integration::Singular         == category () )
         { return   QAGP ( &F1 ) ; }                                // RETURN
      else if   (  GaudiMath::Integration::Finite           == category () )
        if      (  GaudiMath::Integration::NonAdaptive      == type     () )
          { return QNG  ( &F1 ) ; }                                // RETURN
        else if (  GaudiMath::Integration::Adaptive         == type     () )
          { return QAG  ( &F1 ) ; }                                // RETURN
        else if (  GaudiMath::Integration::AdaptiveSingular == type     () )
          { return QAGS ( &F1 ) ; }                                // RETURN
        else
          { Exception ( "::operator(): invalid type "  ); }
      else
        { Exception ( "::operator(): invalid category "  ); }
      return 0 ;
    }

NumericalDefiniteIntegral& Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::operator= ( const NumericalDefiniteIntegral & )

private

Genfun::Derivative Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::partial ( unsigned int index ) const

virtual

Derivatives.

Definition at line 499 of file NumericalDefiniteIntegral.cpp.

    {
      if      ( idx >= m_DIM   )
        { Exception ( "::partial(i): invalid variable index " ) ; };
      //
      const AbsFunction& aux = NumericalDerivative( *this , idx ) ;
      return Genfun::FunctionNoop( &aux ) ;
    }

const Points& Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::points ( ) const

inline

known singularities

Definition at line 325 of file NumericalDefiniteIntegral.h.

{ return m_points ; }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::QAG ( _Function * fun ) const

protected

Definition at line 692 of file NumericalDefiniteIntegral.cpp.

    {
      if( 0 == F ) { Exception("QAG::invalid function") ; }
      // allocate workspace
      if( 0 == ws () ) { allocate () ; }
      // integration limits
      const double low  = std::min ( m_a , m_b ) ;
      const double high = std::max ( m_a , m_b ) ;
      int ierror =
        gsl_integration_qag ( F->fn                    ,
                              low       ,         high ,
                              m_epsabs  ,     m_epsrel ,
                              size   () , (int) rule() , ws ()->ws ,
                              &m_result ,     &m_error             );
      if( ierror ) { gsl_error( "NumericalDefiniteIntegral::QAG " ,
                                __FILE__ , __LINE__ , ierror ) ; }
      // sign
      if ( m_a > m_b ) { m_result *= -1 ; }
      return m_result ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::QAGI ( _Function * fun ) const

protected

Definition at line 589 of file NumericalDefiniteIntegral.cpp.

    {
      // check the argument
      if ( 0 == F    ) { Exception("::QAGI: invalid function"); }
      // allocate workspace
      if ( 0 == ws() ) { allocate() ; }
      int ierror = 0 ;
      if ( m_ia && m_ib )
        {
          ierror = gsl_integration_qagi  ( F->fn                ,
                                           m_epsabs  , m_epsrel ,
                                           size ()   , ws()->ws ,
                                           &m_result , &m_error ) ;
        }
      else if ( m_ia )
        {
          ierror = gsl_integration_qagil ( F->fn     , m_b      ,
                                           m_epsabs  , m_epsrel ,
                                           size ()   , ws()->ws ,
                                           &m_result , &m_error ) ;
        }
      else if ( m_ib )
        {
          ierror = gsl_integration_qagiu ( F->fn     , m_a      ,
                                           m_epsabs  , m_epsrel ,
                                           size ()   , ws()->ws ,
                                           &m_result , &m_error ) ;
        }
      else
        { Exception ( "::QAGI: invalid mode" ) ; };
      if( ierror ) { gsl_error( "NumericalDefiniteIntegral::QAGI" ,
                                __FILE__ , __LINE__ , ierror ) ;}
      return m_result ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::QAGP ( _Function * fun ) const

protected

Definition at line 633 of file NumericalDefiniteIntegral.cpp.

    {
      if( 0 == F ) { Exception("QAGP::invalid function") ; }
      // no known singular points ?
      if( points().empty() || 0 == m_pdata ) { return QAGS( F ) ; }
      const size_t npts = points().size();
      // use GSL
      int ierror =
        gsl_integration_qagp ( F->fn                ,
                               m_pdata   , npts     ,
                               m_epsabs  , m_epsrel ,
                               size ()   , ws()->ws ,
                               &m_result , &m_error ) ;
      if( ierror ) { gsl_error( "NumericalDefiniteIntegral::QAGI " ,
                                __FILE__ , __LINE__ , ierror ) ; }
      // the sign
      if ( m_a > m_b ) { m_result *= -1 ; }
      return m_result ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::QAGS ( _Function * fun ) const

protected

Definition at line 723 of file NumericalDefiniteIntegral.cpp.

    {
      if( 0 == F ) { Exception("QAG::invalid function") ; }
      if ( m_a == m_b )
        {
          m_result = 0    ;
          m_error  = 0    ;   // EXACT !
          return m_result ;
        }
      // allocate workspace
      if( 0 == ws () ) { allocate () ; }
      // integration limits
      const double low  = std::min ( m_a , m_b ) ;
      const double high = std::max ( m_a , m_b ) ;
      int ierror =
        gsl_integration_qags ( F->fn                ,
                               low       , high     ,
                               m_epsabs  , m_epsrel ,
                               size   () , ws()->ws ,
                               &m_result , &m_error );
      if( ierror ) { gsl_error( "NumericalDefiniteIntegral::QAGS " ,
                                __FILE__ , __LINE__ , ierror ) ; }
      // sign
      if ( m_a > m_b ) { m_result *= -1 ; }
      return m_result ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::QNG ( _Function * fun ) const

protected

Definition at line 663 of file NumericalDefiniteIntegral.cpp.

    {
      if( 0 == F ) { Exception("QNG::invalid function") ; }
      // integration limits
      const double low  = std::min ( m_a , m_b ) ;
      const double high = std::max ( m_a , m_b ) ;
      size_t neval = 0 ;
      int ierror =
        gsl_integration_qng ( F->fn                 ,
                              low       ,      high ,
                              m_epsabs  ,  m_epsrel ,
                              &m_result , &m_error  , &neval  ) ;
      if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QNG " ,
                                __FILE__ , __LINE__ , ierror ) ; }
      // sign
      if ( m_a > m_b ) { m_result *= -1 ; }
      return m_result ;
    }

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::result ( ) const

inline

previous result

Definition at line 332 of file NumericalDefiniteIntegral.h.

{ return m_result ; }

GaudiMath::Integration::KronrodRule Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::rule ( ) const

inline

integration rule

Definition at line 347 of file NumericalDefiniteIntegral.h.

{ return m_rule ; }

size_t Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::size ( ) const

inline

Definition at line 337 of file NumericalDefiniteIntegral.h.

{ return m_size ; }

GaudiMath::Integration::Type Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::type ( ) const

inline

integration type

Definition at line 341 of file NumericalDefiniteIntegral.h.

{ return m_type ; }

_Workspace* Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::ws ( ) const

inlineprotected

Definition at line 365 of file NumericalDefiniteIntegral.h.

{ return m_ws ; };

Member Data Documentation

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_a

private

Definition at line 386 of file NumericalDefiniteIntegral.h.

Argument Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_argF

mutableprivate

Definition at line 408 of file NumericalDefiniteIntegral.h.

Argument Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_argument

mutableprivate

Definition at line 407 of file NumericalDefiniteIntegral.h.

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_b

private

Definition at line 387 of file NumericalDefiniteIntegral.h.

GaudiMath::Integration::Category Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_category

private

Definition at line 392 of file NumericalDefiniteIntegral.h.

size_t Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_DIM

private

Definition at line 383 of file NumericalDefiniteIntegral.h.

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_epsabs

private

Definition at line 398 of file NumericalDefiniteIntegral.h.

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_epsrel

private

Definition at line 399 of file NumericalDefiniteIntegral.h.

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_error

mutableprivate

Definition at line 402 of file NumericalDefiniteIntegral.h.

const AbsFunction* Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_function

private

Definition at line 382 of file NumericalDefiniteIntegral.h.

bool Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_ia

private

Definition at line 388 of file NumericalDefiniteIntegral.h.

bool Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_ib

private

Definition at line 389 of file NumericalDefiniteIntegral.h.

size_t Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_index

private

Definition at line 384 of file NumericalDefiniteIntegral.h.

double* Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_pdata

private

Definition at line 396 of file NumericalDefiniteIntegral.h.

Points Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_points

private

Definition at line 395 of file NumericalDefiniteIntegral.h.

double Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_result

mutableprivate

Definition at line 401 of file NumericalDefiniteIntegral.h.

GaudiMath::Integration::KronrodRule Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_rule

private

Definition at line 393 of file NumericalDefiniteIntegral.h.

size_t Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_size

private

Definition at line 404 of file NumericalDefiniteIntegral.h.

GaudiMath::Integration::Type Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_type

private

Definition at line 391 of file NumericalDefiniteIntegral.h.

_Workspace* Genfun::GaudiMathImplementation::NumericalDefiniteIntegral::m_ws

mutableprivate

Definition at line 405 of file NumericalDefiniteIntegral.h.

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

/afs/.cern.ch/sw/Gaudi/releases/GAUDI/GAUDI_v23r9/GaudiGSL/GaudiMath/NumericalDefiniteIntegral.h
/afs/.cern.ch/sw/Gaudi/releases/GAUDI/GAUDI_v23r9/GaudiGSL/src/Lib/NumericalDefiniteIntegral.cpp

Gaudi Framework, version v23r9

Classes

Public Types

Public Member Functions

Protected Member Functions

Private Member Functions

Private Attributes

Detailed Description

Member Typedef Documentation

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation