Gaudi Framework, version v20r2

Generated: 18 Jul 2008

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral Class Reference

#include <NumericalIndefiniteIntegral.h>

Collaboration diagram for Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral:

Collaboration graph
[legend]
List of all members.

Detailed Description

The simple class for numerical integrations.

It allows to evaluate following indefinite integrals:

\[ {\mathcal{F}}_i \left(x_1, \dots , x_{i-1}, x_i , x_{i+1}, \dots , x_n \right) = \int\limits_{a}^{x_i} 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 , x_{i+1}, \dots , x_n \right) = \int\limits_{x_i}^{a} 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 , x_{i+1}, \dots , x_n \right) = \int\limits_{-\infty}^{x_i} 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 , x_{i+1}, \dots , x_n \right) = \int\limits_{x_i}^{+\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.Belyaev@itep.ru
Date:
2003-08-31

Definition at line 76 of file NumericalIndefiniteIntegral.h.

Public Types

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

Public Member Functions

 FUNCTION_OBJECT_DEF (NumericalIndefiniteIntegral)
 From CLHEP/GenericFunctions.
 NumericalIndefiniteIntegral (const AbsFunction &function, const size_t index, const double a, const GaudiMath::Integration::Limit limit=GaudiMath::Integration::VariableHighLimit, const GaudiMath::Integration::Type type=GaudiMath::Integration::Adaptive, const GaudiMath::Integration::KronrodRule rule=GaudiMath::Integration::Default, const double epsabs=1e-10, const double epsrel=1.e-7, const size_t size=1000)
 Standard constructor.
 NumericalIndefiniteIntegral (const AbsFunction &function, const size_t index, const double a, const Points &points, const GaudiMath::Integration::Limit limit=GaudiMath::Integration::VariableHighLimit, const double epsabs=1e-9, const double epsrel=1.e-6, const size_t size=1000)
 standard constructor
 NumericalIndefiniteIntegral (const AbsFunction &function, const size_t index, const GaudiMath::Integration::Limit limit=GaudiMath::Integration::VariableHighLimit, const double epsabs=1e-9, const double epsrel=1.e-6, const size_t size=1000)
 standard constructor the integral limt is assumed to be infinity
 NumericalIndefiniteIntegral (const NumericalIndefiniteIntegral &)
 copy constructor
virtual ~NumericalIndefiniteIntegral ()
 destructor
virtual unsigned int dimensionality () const
 dimensionality of the problem
virtual double operator() (double argument) const
 evaluate the function
virtual double operator() (const Argument &argument) const
 evaluate the function
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
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
const size_t size () const
GaudiMath::Integration::Limit limit () const
 integration limit
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

 NumericalIndefiniteIntegral ()
NumericalIndefiniteIntegral & operator= (const NumericalIndefiniteIntegral &)

Private Attributes

const AbsFunction * m_function
size_t m_DIM
size_t m_index
double m_a
GaudiMath::Integration::Limit m_limit
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

Classes

struct  _Function
struct  _Workspace

Member Typedef Documentation

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

typedef for vector of singular points

Definition at line 80 of file NumericalIndefiniteIntegral.h.

Constructor & Destructor Documentation

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

Standard constructor.

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

for value of limit = VariableHighLimit and the integral

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

for value of limit = VariableLowLimit

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 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.

Parameters:
function the base function
index the variable index
a integration limit
limit flag to distinguisch low variable limit from high variable limit
type the integration type (adaptive, non-adaptive or adaptive for singular functions
key Gauss-Kronrad integration rule
epsabs absolute precision for integration
epsrel relative precision for integration
lim bisection limit

Definition at line 71 of file NumericalIndefiniteIntegral.cpp.

References GaudiMath::Integration::Default, Exception(), GaudiMath::Integration::Fixed, m_DIM, m_index, and m_rule.

00080       : AbsFunction () 
00081       , m_function  ( function.clone()                    )
00082       , m_DIM       ( function.dimensionality()           )
00083       , m_index     ( index                               )
00084       , m_a         ( a                                   ) 
00085       , m_limit     ( limit                               )
00086       , m_type      ( type                                ) 
00087       , m_category  ( GaudiMath::Integration::Finite      ) 
00088       , m_rule      ( rule                                ) 
00089       //
00090       , m_points    (                                     ) 
00091       , m_pdata     ( 0                                   )
00092       //
00093       , m_epsabs    ( epsabs                              ) 
00094       , m_epsrel    ( epsrel                              ) 
00095       //
00096       , m_result    ( GSL_NEGINF                          ) 
00097       , m_error     ( GSL_POSINF                          ) 
00098       //
00099       , m_size      ( size                                ) 
00100       , m_ws        ( 0                                   ) 
00101       , m_argument  ( function.dimensionality()           ) 
00102     {
00103       if ( GaudiMath::Integration::Fixed == m_rule ) 
00104         { m_rule = GaudiMath::Integration::Default ; }
00105       if ( m_index >= m_DIM  ) 
00106         { Exception("::constructor: invalid variable index") ; }
00107     };

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::NumericalIndefiniteIntegral ( const AbsFunction &  function,
const size_t  index,
const double  a,
const Points &  points,
const GaudiMath::Integration::Limit  limit = GaudiMath::Integration::VariableHighLimit,
const double  epsabs = 1e-9,
const double  epsrel = 1.e-6,
const size_t  size = 1000 
)

standard constructor

The function, created with this constructor evaluates following indefinite integral:

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

for value of limit = VariableHighLimit and the integral

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

for value of limit = VariableLowLimit

The integrand is assumed to have a known discontinuities

Parameters:
function the base function
index the variable index
a integration limit
limit flag to distinguisch low variable limit from high variable limit
points list of known function singularities
epsabs absolute precision for integration
epsrel relative precision for integration

Definition at line 121 of file NumericalIndefiniteIntegral.cpp.

References std::vector< _Tp, _Alloc >::begin(), std::vector< _Tp, _Alloc >::end(), std::vector< _Tp, _Alloc >::erase(), Exception(), m_DIM, m_index, m_pdata, m_points, std::vector< _Tp, _Alloc >::push_back(), std::vector< _Tp, _Alloc >::size(), std::sort(), and std::unique().

00129       : AbsFunction () 
00130       , m_function  ( function.clone()          )
00131       , m_DIM       ( function.dimensionality() ) 
00132       , m_index     ( index                     ) 
00133       , m_a         ( a                         )
00134       , m_limit     ( limit                     )
00135       , m_type      ( GaudiMath::Integration::    Other )
00136       , m_category  ( GaudiMath::Integration:: Singular )
00137       , m_rule      ( GaudiMath::Integration::    Fixed )
00138       , m_points    ( points  ) 
00139       , m_pdata     ( 0       ) 
00140       , m_epsabs    ( epsabs  )
00141       , m_epsrel    ( epsrel  ) 
00142       //
00143       , m_result    ( GSL_NEGINF                          ) 
00144       , m_error     ( GSL_POSINF                          ) 
00145       //
00146       , m_size      ( size                                ) 
00147       , m_ws        ( 0                                   ) 
00148       , m_argument  ( function.dimensionality()           ) 
00149     {
00150       if ( m_index >= m_DIM ) 
00151         { Exception("::constructor: invalid variable index") ; } 
00152       m_pdata = new double[ 2 + m_points.size() ] ;
00153       m_points.push_back( a ) ;
00154       std::sort( m_points.begin() , m_points.end() ) ;
00155       m_points.erase ( std::unique( m_points.begin () , 
00156                                     m_points.end   () ) , m_points.end() );
00157     };

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::NumericalIndefiniteIntegral ( const AbsFunction &  function,
const size_t  index,
const GaudiMath::Integration::Limit  limit = GaudiMath::Integration::VariableHighLimit,
const double  epsabs = 1e-9,
const double  epsrel = 1.e-6,
const size_t  size = 1000 
)

standard constructor the integral limt is assumed to be infinity

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

for value of limit = VariableHighLimit and the integral

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

for value of limit = VariableLowLimit

Parameters:
function the base function
index the variable index
limit flag to distinguisch low variable limit from high variable limit
singularities list of known function singularities

Definition at line 169 of file NumericalIndefiniteIntegral.cpp.

References Exception(), m_DIM, and m_index.

00175       : AbsFunction () 
00176       , m_function  ( function.clone() ) 
00177       , m_DIM       ( function.dimensionality() )
00178       , m_index     ( index            )
00179       , m_a         ( GSL_NEGINF       ) // should not be used! 
00180       , m_limit     ( limit            )
00181       , m_type      ( GaudiMath::Integration::    Other ) 
00182       , m_category  ( GaudiMath::Integration:: Infinite )
00183       , m_rule      ( GaudiMath::Integration::    Fixed )
00184       , m_points    (            ) 
00185       , m_pdata     ( 0          )
00186       , m_epsabs    ( epsabs     ) 
00187       , m_epsrel    ( epsrel     )
00188       , m_result    ( GSL_NEGINF ) 
00189       , m_error     ( GSL_POSINF ) 
00190       , m_size      ( size       ) 
00191       , m_ws        ( 0          ) 
00192       , m_argument  ( function.dimensionality()           ) 
00193     {
00194       if ( m_index >= m_DIM ) 
00195         { Exception("::constructor: invalid variable index") ; }
00196     };

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::NumericalIndefiniteIntegral ( const NumericalIndefiniteIntegral &   ) 

copy constructor

Definition at line 204 of file NumericalIndefiniteIntegral.cpp.

References m_pdata, m_points, and std::vector< _Tp, _Alloc >::size().

00205       : AbsFunction () 
00206       , m_function  ( right.m_function->clone() ) 
00207       , m_DIM       ( right.m_DIM      ) 
00208       , m_index     ( right.m_index    ) 
00209       , m_a         ( right.m_a        ) 
00210       , m_limit     ( right.m_limit    ) 
00211       , m_type      ( right.m_type     ) 
00212       , m_category  ( right.m_category )
00213       , m_rule      ( right.m_rule     )
00214       , m_points    ( right.m_points   ) 
00215       , m_pdata     ( 0 )           // attention 
00216       , m_epsabs    ( right.m_epsabs   ) 
00217       , m_epsrel    ( right.m_epsrel   ) 
00218       , m_result    ( GSL_NEGINF       ) 
00219       , m_error     ( GSL_POSINF       )
00220       , m_size      ( right.m_size     ) 
00221       , m_ws        ( 0                )
00222       , m_argument  ( right.m_argument )
00223     {
00224       m_pdata = new double[ 2 + m_points.size() ] ; // attention! 
00225     };

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::~NumericalIndefiniteIntegral (  )  [virtual]

destructor

Definition at line 232 of file NumericalIndefiniteIntegral.cpp.

References m_function, m_pdata, m_ws, and Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Workspace::ws.

00233     {
00234       if( 0 != m_ws ) 
00235         { 
00236           gsl_integration_workspace_free ( m_ws->ws ) ;
00237           delete m_ws ;
00238           m_ws = 0 ;
00239         }
00240       if ( 0 != m_pdata    ) { delete m_pdata    ; m_pdata    = 0 ; }
00241       if ( 0 != m_function ) { delete m_function ; m_function = 0 ; }
00242     };

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::NumericalIndefiniteIntegral (  )  [private]

Member Function Documentation

Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::FUNCTION_OBJECT_DEF ( NumericalIndefiniteIntegral   ) 

From CLHEP/GenericFunctions.

virtual unsigned int Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::dimensionality (  )  const [inline, virtual]

dimensionality of the problem

Definition at line 278 of file NumericalIndefiniteIntegral.h.

References m_DIM.

00278 { return m_DIM ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::operator() ( double  argument  )  const [virtual]

evaluate the function

Definition at line 262 of file NumericalIndefiniteIntegral.cpp.

00263     {
00264       // reset the result and the error  
00265       m_result = GSL_NEGINF ;
00266       m_error  = GSL_POSINF ;
00267       
00268       // check the argument 
00269       if( 1 != m_DIM ) { Exception ( "operator(): invalid argument size " ) ; };
00270       
00271       m_argument[0] = argument ;
00272       return (*this) ( m_argument );
00273     };

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::operator() ( const Argument &  argument  )  const [virtual]

evaluate the function

Definition at line 302 of file NumericalIndefiniteIntegral.cpp.

References GaudiMath::Integration::Adaptive, GaudiMath::Integration::AdaptiveSingular, GaudiMath::Integration::Finite, Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, Genfun::GaudiMathImplementation::GSL_Adaptor(), GaudiMath::Integration::Infinite, GaudiMath::Integration::NonAdaptive, and GaudiMath::Integration::Singular.

00303     {
00304       // reset the result and the error  
00305       m_result = GSL_NEGINF ;
00306       m_error  = GSL_POSINF ;
00307       
00308       // check the argument 
00309       if( argument.dimension() != m_DIM ) 
00310         { Exception ( "operator(): invalid argument size " ) ; };
00311       
00312       // copy the argument 
00313       {for( size_t i  = 0 ; i < m_DIM ; ++i ){ m_argument[i] = argument[i];}}
00314       
00315       // create the helper object 
00316       GSL_Helper helper( *m_function , m_argument , m_index );
00317       
00318       // use GSL to evaluate the numerical derivative 
00319       gsl_function F ;
00320       F.function = &GSL_Adaptor ;
00321       F.params   = &helper                 ;
00322       _Function F1    ;
00323       F1.fn      = &F ;
00324       
00325       if        (  GaudiMath::Integration::Infinite         == category () )
00326         { return   QAGI ( &F1 ) ; }                                // RETURN
00327       else if   (  GaudiMath::Integration::Singular         == category () )
00328         { return   QAGP ( &F1 ) ; }                                // RETURN
00329       else if   (  GaudiMath::Integration::Finite           == category () )
00330         if      (  GaudiMath::Integration::NonAdaptive      == type     () ) 
00331           { return QNG  ( &F1 ) ; }                                // RETURN
00332         else if (  GaudiMath::Integration::Adaptive         == type     () ) 
00333           { return QAG  ( &F1 ) ; }                                // RETURN
00334         else if (  GaudiMath::Integration::AdaptiveSingular == type     () ) 
00335           { return QAGS ( &F1 ) ; }                                // RETURN
00336         else 
00337           { Exception ( "::operator(): invalid type "  ); }
00338       else 
00339         { Exception ( "::operator(): invalid category "  ); }
00340       
00341       return 0 ;
00342     };

virtual bool Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::hasAnalyticDerivative (  )  const [inline, virtual]

Does this function have an analytic derivative?

Definition at line 286 of file NumericalIndefiniteIntegral.h.

00286 { return true ;}

Genfun::Derivative Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::partial ( unsigned int  index  )  const [virtual]

Derivatives.

Definition at line 280 of file NumericalIndefiniteIntegral.cpp.

References Exception(), function(), m_DIM, m_index, and GaudiMath::Integration::VariableLowLimit.

00281     {
00282       if      ( idx >= m_DIM   )
00283         { Exception ( "::partial(i): invalid variable index " ) ; };
00284       if      ( idx != m_index )
00285         {
00286           const AbsFunction& aux = NumericalDerivative( *this , idx ) ;
00287           return Genfun::FunctionNoop( &aux ) ;
00288         }
00289       else if ( GaudiMath::Integration::VariableLowLimit == limit () ) 
00290         { 
00291           const AbsFunction& aux = -1 * function() ;
00292           return Genfun::FunctionNoop( &aux ) ;
00293         }
00294       const AbsFunction& aux = function() ;
00295       return Genfun::FunctionNoop( &aux ) ;
00296     };

const AbsFunction& Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::function (  )  const [inline]

accessor to the function itself

Definition at line 294 of file NumericalIndefiniteIntegral.h.

References m_function.

Referenced by partial().

00294 { return *m_function ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::a (  )  const [inline]

integration limit

Definition at line 296 of file NumericalIndefiniteIntegral.h.

References m_a.

00296 { return         m_a ; }

const Points& Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::points (  )  const [inline]

known singularities

Definition at line 298 of file NumericalIndefiniteIntegral.h.

References m_points.

Referenced by QAGP().

00298 { return    m_points ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::epsabs (  )  const [inline]

absolute precision

Definition at line 300 of file NumericalIndefiniteIntegral.h.

References m_epsabs.

00300 { return    m_epsabs ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::epsrel (  )  const [inline]

relatiove precision

Definition at line 302 of file NumericalIndefiniteIntegral.h.

References m_epsrel.

00302 { return    m_epsrel ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::result (  )  const [inline]

previous result

Definition at line 305 of file NumericalIndefiniteIntegral.h.

References m_result.

00305 { return    m_result ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::error (  )  const [inline]

evaluate of previous error

Definition at line 307 of file NumericalIndefiniteIntegral.h.

References m_error.

00307 { return     m_error ; }

const size_t Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::size (  )  const [inline]

Definition at line 310 of file NumericalIndefiniteIntegral.h.

References m_size.

00310 { return      m_size ; }

GaudiMath::Integration::Limit Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::limit (  )  const [inline]

integration limit

Definition at line 314 of file NumericalIndefiniteIntegral.h.

References m_limit.

Referenced by QAG(), QAGI(), QAGP(), QAGS(), and QNG().

00314 { return     m_limit ; }

GaudiMath::Integration::Type Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::type (  )  const [inline]

integration type

Definition at line 317 of file NumericalIndefiniteIntegral.h.

References m_type.

00317 { return      m_type ; }

GaudiMath::Integration::Category Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::category (  )  const [inline]

integration category

Definition at line 320 of file NumericalIndefiniteIntegral.h.

References m_category.

00320 { return  m_category ; }

GaudiMath::Integration::KronrodRule Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::rule (  )  const [inline]

integration rule

Definition at line 323 of file NumericalIndefiniteIntegral.h.

References m_rule.

Referenced by QAG().

00323 { return      m_rule ; }

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::QAGI ( _Function *  fun  )  const [protected]

Definition at line 364 of file NumericalIndefiniteIntegral.cpp.

References allocate(), Exception(), Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, limit(), m_argument, m_epsabs, m_epsrel, m_error, m_index, m_result, GaudiMath::Integration::VariableHighLimit, GaudiMath::Integration::VariableLowLimit, and ws().

00365     {
00366       // check the argument 
00367       if( 0 == F ) { Exception("::QAGI: invalid function"); }
00368       
00369       const double x = m_argument[m_index] ;
00370       
00371       // allocate workspace 
00372       if( 0 == ws() ) { allocate() ; }
00373       
00374       int ierror = 0 ;
00375       switch ( limit() ) 
00376         {
00377         case GaudiMath::Integration::VariableLowLimit  : 
00378           ierror = gsl_integration_qagiu ( F->fn     , x        ,  
00379                                            m_epsabs  , m_epsrel , 
00380                                            size ()   , ws()->ws , 
00381                                            &m_result , &m_error ) ; break ;
00382         case GaudiMath::Integration::VariableHighLimit :
00383           ierror = gsl_integration_qagil ( F->fn     , x        ,
00384                                            m_epsabs  , m_epsrel , 
00385                                            size ()   , ws()->ws , 
00386                                            &m_result , &m_error ) ; break ;
00387         default :
00388           Exception ( "::QAGI: invalid mode" ) ;
00389         };
00390       
00391       if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QAGI" , 
00392                                 __FILE__ , __LINE__ , ierror ) ;}
00393       
00394       return m_result ;
00395     };

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::QAGP ( _Function *  fun  )  const [protected]

Definition at line 401 of file NumericalIndefiniteIntegral.cpp.

References b, std::vector< _Tp, _Alloc >::begin(), std::copy(), std::vector< _Tp, _Alloc >::end(), Exception(), Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, std::vector< _Tp, _Alloc >::insert(), limit(), std::lower_bound(), m_a, m_argument, m_epsabs, m_epsrel, m_error, m_index, m_pdata, m_result, std::max(), std::min(), points(), QAGS(), std::vector< _Tp, _Alloc >::size(), upper(), std::upper_bound(), GaudiMath::Integration::VariableHighLimit, GaudiMath::Integration::VariableLowLimit, and ws().

00402     {
00403       if( 0 == F ) { Exception("QAGP::invalid function") ; }
00404       
00405       const double x = m_argument[m_index] ;
00406       if ( m_a == x  ) 
00407         { 
00408           m_result = 0    ;
00409           m_error  = 0    ;   // EXACT !
00410           return m_result ; 
00411         }
00412       
00413       // no known singular points ?
00414       if( points().empty() ) { return QAGS( F ) ; }
00415       
00416       // integration limits 
00417       const double a = std::min ( m_a , x ) ;
00418       const double b = std::max ( m_a , x ) ;
00419       
00420       // "active" singular points
00421       Points::const_iterator lower = 
00422         std::lower_bound ( points().begin() , points().end() , a ) ;
00423       Points::const_iterator upper = 
00424         std::upper_bound ( points().begin() , points().end() , b ) ;
00425       
00426       Points pnts ( upper - lower ) ;
00427       std::copy( lower , upper , pnts.begin() );
00428       if ( *lower != a       ) { pnts.insert( pnts.begin () , a ) ; }
00429       if ( *upper != b       ) { pnts.insert( pnts.end   () , b ) ; }
00430       std::copy( pnts.begin() , pnts.end() , m_pdata );
00431       const size_t npts = pnts.size() ;
00432       
00433       // use GSL 
00434       int ierror = 
00435         gsl_integration_qagp ( F->fn                , 
00436                                m_pdata   , npts     , 
00437                                m_epsabs  , m_epsrel ,
00438                                size ()   , ws()->ws , 
00439                                &m_result , &m_error ) ;
00440       
00441       if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QAGI " , 
00442                                 __FILE__ , __LINE__ , ierror ) ; }
00443       
00444       // sign
00445       if      ( GaudiMath::Integration::VariableHighLimit == limit() 
00446                 &&  x < m_a  ) { m_result *= -1 ; }
00447       else if ( GaudiMath::Integration::VariableLowLimit  == limit() 
00448                 &&  x > m_a  ) { m_result *= -1 ; }
00449       
00450       return m_result ;
00451     };

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::QNG ( _Function *  fun  )  const [protected]

Definition at line 457 of file NumericalIndefiniteIntegral.cpp.

References b, Exception(), Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, limit(), m_a, m_argument, m_epsabs, m_epsrel, m_error, m_index, m_result, std::max(), std::min(), GaudiMath::Integration::VariableHighLimit, and GaudiMath::Integration::VariableLowLimit.

00458     {
00459       if( 0 == F ) { Exception("QNG::invalid function") ; }
00460       
00461       const double x = m_argument[m_index] ;
00462       if ( m_a == x  ) 
00463         { 
00464           m_result = 0    ;
00465           m_error  = 0    ;   // EXACT !
00466           return m_result ; 
00467         }
00468       
00469       // integration limits 
00470       const double a = std::min ( m_a , x ) ;
00471       const double b = std::max ( m_a , x ) ;
00472       
00473       size_t neval = 0 ;
00474       int ierror = 
00475         gsl_integration_qng ( F->fn                 ,  
00476                               a         ,         b , 
00477                               m_epsabs  ,  m_epsrel , 
00478                               &m_result , &m_error  , &neval  ) ;
00479       
00480       if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QNG " , 
00481                                 __FILE__ , __LINE__ , ierror ) ; }
00482       
00483       // sign
00484       if      ( GaudiMath::Integration::VariableHighLimit == limit() 
00485                 &&  x < m_a  ) { m_result *= -1 ; }
00486       else if ( GaudiMath::Integration::VariableLowLimit  == limit() 
00487                 &&  x > m_a  ) { m_result *= -1 ; }
00488 
00489       return m_result ;
00490     };

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::QAG ( _Function *  fun  )  const [protected]

Definition at line 495 of file NumericalIndefiniteIntegral.cpp.

References allocate(), b, Exception(), Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, limit(), m_a, m_argument, m_epsabs, m_epsrel, m_error, m_index, m_result, std::max(), std::min(), rule(), GaudiMath::Integration::VariableHighLimit, GaudiMath::Integration::VariableLowLimit, and ws().

00496     {
00497       if( 0 == F ) { Exception("QAG::invalid function") ; }
00498       
00499       const double x = m_argument[m_index] ;
00500       if ( m_a == x  ) 
00501         { 
00502           m_result = 0    ;
00503           m_error  = 0    ;   // EXACT !
00504           return m_result ; 
00505         }
00506       
00507       // allocate workspace 
00508       if( 0 == ws () ) { allocate () ; }
00509       
00510       // integration limits 
00511       const double a = std::min ( m_a , x ) ;
00512       const double b = std::max ( m_a , x ) ;
00513       
00514       int ierror = 
00515         gsl_integration_qag ( F->fn                    , 
00516                               a         ,            b , 
00517                               m_epsabs  ,     m_epsrel , 
00518                               size   () , (int) rule() , ws ()->ws , 
00519                               &m_result ,     &m_error             );
00520       
00521       if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QAG " , 
00522                                 __FILE__ , __LINE__ , ierror ) ; }
00523       
00524       // sign
00525       if      ( GaudiMath::Integration::VariableHighLimit == limit() 
00526                 &&  x < m_a  ) { m_result *= -1 ; }
00527       else if ( GaudiMath::Integration::VariableLowLimit  == limit() 
00528                 &&  x > m_a  ) { m_result *= -1 ; }
00529       
00530       return m_result ;
00531     };

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::QAGS ( _Function *  fun  )  const [protected]

Definition at line 537 of file NumericalIndefiniteIntegral.cpp.

References allocate(), b, Exception(), Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Function::fn, limit(), m_a, m_argument, m_epsabs, m_epsrel, m_error, m_index, m_result, std::max(), std::min(), GaudiMath::Integration::VariableHighLimit, GaudiMath::Integration::VariableLowLimit, and ws().

Referenced by QAGP().

00538     {
00539       if( 0 == F ) { Exception("QAG::invalid function") ; }
00540       
00541       const double x = m_argument[m_index] ;
00542       if ( m_a == x  ) 
00543         { 
00544           m_result = 0    ;
00545           m_error  = 0    ;   // EXACT !
00546           return m_result ; 
00547         }
00548       
00549       // allocate workspace 
00550       if( 0 == ws () ) { allocate () ; }
00551       
00552       // integration limits 
00553       const double a = std::min ( m_a , x ) ;
00554       const double b = std::max ( m_a , x ) ;
00555       
00556       int ierror = 
00557         gsl_integration_qags ( F->fn                , 
00558                                a         , b        , 
00559                                m_epsabs  , m_epsrel , 
00560                                size   () , ws()->ws , 
00561                                &m_result , &m_error );
00562       
00563       if( ierror ) { gsl_error( "NumericalIndefiniteIntegral::QAGS " , 
00564                                 __FILE__ , __LINE__ , ierror ) ; }
00565       
00566       // sign
00567       if      ( GaudiMath::Integration::VariableHighLimit == limit() 
00568                 &&  x < m_a  ) { m_result *= -1 ; }
00569       else if ( GaudiMath::Integration::VariableLowLimit  == limit() 
00570                 &&  x > m_a  ) { m_result *= -1 ; }
00571       
00572       return m_result ;
00573     };

NumericalIndefiniteIntegral::_Workspace * Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::allocate (  )  const [protected]

allocate the integration workspace

Definition at line 349 of file NumericalIndefiniteIntegral.cpp.

References Exception(), m_ws, and Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::_Workspace::ws.

Referenced by QAG(), QAGI(), and QAGS().

00350     {
00351       if ( 0 != m_ws ) { return m_ws; }
00352       gsl_integration_workspace* aux = 
00353         gsl_integration_workspace_alloc( size () );
00354       if ( 0 == aux ) { Exception ( "allocate()::invalid workspace" ) ; };
00355       m_ws = new _Workspace() ;
00356       m_ws->ws = aux ;
00357       return m_ws ;
00358     };

_Workspace* Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::ws (  )  const [inline, protected]

Definition at line 341 of file NumericalIndefiniteIntegral.h.

References m_ws.

Referenced by QAG(), QAGI(), QAGP(), and QAGS().

00342       { return m_ws ; };

StatusCode Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::Exception ( const std::string &  message,
const StatusCode &  sc = StatusCode::FAILURE 
) const [protected]

Definition at line 249 of file NumericalIndefiniteIntegral.cpp.

Referenced by allocate(), NumericalIndefiniteIntegral(), partial(), QAG(), QAGI(), QAGP(), QAGS(), and QNG().

00251     {
00252       throw GaudiException( "NumericalIndefiniteIntegral::" + message , 
00253                             "*GaudiMath*" , sc ) ;
00254       return sc ;
00255     };

NumericalIndefiniteIntegral& Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::operator= ( const NumericalIndefiniteIntegral &   )  [private]

Member Data Documentation

const AbsFunction* Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_function [private]

Definition at line 359 of file NumericalIndefiniteIntegral.h.

Referenced by function(), and ~NumericalIndefiniteIntegral().

size_t Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_DIM [private]

Definition at line 360 of file NumericalIndefiniteIntegral.h.

Referenced by dimensionality(), NumericalIndefiniteIntegral(), and partial().

size_t Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_index [private]

Definition at line 361 of file NumericalIndefiniteIntegral.h.

Referenced by NumericalIndefiniteIntegral(), partial(), QAG(), QAGI(), QAGP(), QAGS(), and QNG().

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_a [private]

Definition at line 363 of file NumericalIndefiniteIntegral.h.

Referenced by a(), QAG(), QAGP(), QAGS(), and QNG().

GaudiMath::Integration::Limit Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_limit [private]

Definition at line 365 of file NumericalIndefiniteIntegral.h.

Referenced by limit().

GaudiMath::Integration::Type Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_type [private]

Definition at line 366 of file NumericalIndefiniteIntegral.h.

Referenced by type().

GaudiMath::Integration::Category Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_category [private]

Definition at line 367 of file NumericalIndefiniteIntegral.h.

Referenced by category().

GaudiMath::Integration::KronrodRule Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_rule [private]

Definition at line 368 of file NumericalIndefiniteIntegral.h.

Referenced by NumericalIndefiniteIntegral(), and rule().

Points Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_points [private]

Definition at line 370 of file NumericalIndefiniteIntegral.h.

Referenced by NumericalIndefiniteIntegral(), and points().

double* Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_pdata [private]

Definition at line 371 of file NumericalIndefiniteIntegral.h.

Referenced by NumericalIndefiniteIntegral(), QAGP(), and ~NumericalIndefiniteIntegral().

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_epsabs [private]

Definition at line 373 of file NumericalIndefiniteIntegral.h.

Referenced by epsabs(), QAG(), QAGI(), QAGP(), QAGS(), and QNG().

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_epsrel [private]

Definition at line 374 of file NumericalIndefiniteIntegral.h.

Referenced by epsrel(), QAG(), QAGI(), QAGP(), QAGS(), and QNG().

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_result [mutable, private]

Definition at line 376 of file NumericalIndefiniteIntegral.h.

Referenced by QAG(), QAGI(), QAGP(), QAGS(), QNG(), and result().

double Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_error [mutable, private]

Definition at line 377 of file NumericalIndefiniteIntegral.h.

Referenced by error(), QAG(), QAGI(), QAGP(), QAGS(), and QNG().

size_t Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_size [private]

Definition at line 379 of file NumericalIndefiniteIntegral.h.

Referenced by size().

_Workspace* Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_ws [mutable, private]

Definition at line 380 of file NumericalIndefiniteIntegral.h.

Referenced by allocate(), ws(), and ~NumericalIndefiniteIntegral().

Argument Genfun::GaudiMathImplementation::NumericalIndefiniteIntegral::m_argument [mutable, private]

Definition at line 382 of file NumericalIndefiniteIntegral.h.

Referenced by QAG(), QAGI(), QAGP(), QAGS(), and QNG().

The documentation for this class was generated from the following files:
Generated at Fri Jul 18 12:10:40 2008 for Gaudi Framework, version v20r2 by Doxygen version 1.5.1 written by Dimitri van Heesch, © 1997-2004