| The Gaudi Framework v28r2p1 (f1a77ff4)

- Main Page
- Related Pages
- Modules
- Namespaces
- Classes
- Files

- File List
- File Members

NumericalDerivative.cpp
Go to the documentation of this file.
1 // ============================================================================
2 // Include files
3 // ============================================================================
4 // GaudiKernel
5 // ============================================================================
6 #include "GaudiKernel/GaudiException.h"
7 // ============================================================================
8 // local
9 // ============================================================================
10 #include "GaudiMath/NumericalDerivative.h"
11 #include "GaudiMath/Constant.h"
12 // ============================================================================
13 // GSL
14 // ============================================================================
15 #include "gsl/gsl_errno.h"
16 #include "gsl/gsl_diff.h"
17 #include "gsl/gsl_math.h"
18 // ============================================================================
19 // local
20 // ============================================================================
21 #include "Helpers.h"
22 // ============================================================================
23
24 // ============================================================================
31 // ============================================================================
32
33 namespace Genfun
34 {
35 namespace GaudiMathImplementation
36 {
37 // ========================================================================
38
39 // ========================================================================
46 // ========================================================================
47 NumericalDerivative::NumericalDerivative
48 ( const AbsFunction& function ,
49 const size_t index ,
50 const NumericalDerivative::Type& type )
51 : AbsFunction ()
52 , m_function ( function.clone() )
53 , m_index ( index )
54 , m_DIM ( function.dimensionality() )
55 , m_type ( type )
56 , m_argument ( function.dimensionality() )
57 , m_result ( GSL_NEGINF )
58 , m_error ( GSL_POSINF )
59 {
60 if( m_index >= m_DIM )
61 { Exception ( "::constructor invalid variable index " ) ; };
62 }
63 // ========================================================================
64
65 // ========================================================================
67 // ========================================================================
68 NumericalDerivative::NumericalDerivative
69 ( const NumericalDerivative& right )
70 : AbsFunction ()
71 , m_function ( right.m_function->clone() )
72 , m_index ( right.m_index )
73 , m_DIM ( right.m_DIM )
74 , m_type ( right.m_type )
75 , m_argument ( right.m_DIM )
76 , m_result ( GSL_NEGINF )
77 , m_error ( GSL_POSINF )
78 {}
79 // ========================================================================
80
81 // ========================================================================
83 // ========================================================================
84 const NumericalDerivative::Type&
85 NumericalDerivative::setType ( const NumericalDerivative::Type& value )
86 { m_type = value ; return type() ; }
87 // ========================================================================
88
89 // ========================================================================
91 // ========================================================================
92 Genfun::Derivative NumericalDerivative::partial( unsigned int idx ) const
93 {
94 if( idx >= m_DIM )
95 { Exception ( "::partial(i): invalid variable index" ) ; }
96 const AbsFunction& aux =
97 NumericalDerivative( *this , idx , type() ) ;
98 return FunctionNoop( &aux ) ;
99 }
100 // ========================================================================
101
102 // ========================================================================
104 // ========================================================================
105 double NumericalDerivative::operator()
106 ( const Argument& argument ) const
107 {
108 // reset the result and the error
109 m_result = GSL_NEGINF ;
110 m_error = GSL_POSINF ;
111
112 // check the argument
113 if( argument.dimension() != m_DIM )
114 { Exception ( "::operator():invalid argument size" ) ; };
115
117 {for( size_t i = 0 ; i < m_DIM ; ++i ){ m_argument[i] = argument[i];}}
118
119 // create the helper object
120 GSL_Helper helper( *m_function , m_argument , m_index );
121
123 gsl_function F ;
124 F.function = &GSL_Adaptor ;
125 F.params = &helper ;
126
127 double x = argument[m_index];
128 int ierrno = 0 ;
129 switch ( type() )
130 {
131 case Central :
132 ierrno = gsl_diff_central ( &F , x , &m_result , &m_error ) ; break ;
133 case Forward :
134 ierrno = gsl_diff_forward ( &F , x , &m_result , &m_error ) ; break ;
135 case Backward :
136 ierrno = gsl_diff_backward ( &F , x , &m_result , &m_error ) ; break ;
137 default:
138 Exception ( "::operator(): invalid diffrentiation type " ) ;
139 }
140
141 if( ierrno )
142 { gsl_error ( " NumericalDerivative:: the error from gsl_diff_XXXX" ,
143 __FILE__ , __LINE__ , ierrno ) ; }
144
145 return result() ;
146 }
147 // ========================================================================
148
149 // ========================================================================
151 // ========================================================================
152 double NumericalDerivative::operator() ( const double argument ) const
153 {
154 // reset the result and the error
155 m_result = GSL_NEGINF ;
156 m_error = GSL_POSINF ;
157
158 // check the argument
159 if( 1 != m_DIM ) { Exception ( "operator(): invalid argument size " ); }
160
161 Argument arg(1) ;
162 arg[0] = argument ;
163 return (*this)( arg ) ;
164 }
165 // ========================================================================
166
167 // ========================================================================
169 // ========================================================================
170 StatusCode NumericalDerivative::Exception
171 ( const std::string& message ,
172 const StatusCode& sc ) const
173 {
174 throw GaudiException( "NumericalDerivative" + message ,
175 "*GaudiMath*" , sc ) ;
176 return sc ;
177 }
178 // ========================================================================
179
180
181 } // end of namespace GaudiMathImplementation
182
183 } // end of namespace Genfun (due to GSL)
184
185 // ============================================================================
186 // The END
187 // ============================================================================
GaudiException
Define general base for Gaudi exception.
Definition: GaudiException.h:23
Genfun::GaudiMathImplementation::NumericalDerivative::type
const Type & type() const
type
Definition: NumericalDerivative.h:124
Genfun::GaudiMathImplementation::NumericalDerivative::Type
Type
the type of numerical differentiation
Definition: NumericalDerivative.h:46
Genfun::GaudiMathImplementation::NumericalDerivative::m_argument
Argument m_argument
Definition: NumericalDerivative.h:153
Genfun::GaudiMathImplementation::NumericalDerivative::m_function
std::unique_ptr< const AbsFunction > m_function
Definition: NumericalDerivative.h:148
Genfun::GaudiMathImplementation::NumericalDerivative::NumericalDerivative
NumericalDerivative()
GaudiMath::Derivative
Genfun::GaudiMathImplementation::NumericalDerivative Derivative
Definition: GaudiMath.h:29
Genfun::GaudiMathImplementation::NumericalDerivative::operator()
double operator()(double argument) const override
Function value.
Definition: NumericalDerivative.cpp:152
Genfun::GaudiMathImplementation::NumericalDerivative::partial
Derivative partial(unsigned int index) const override
Derivatives.
Definition: NumericalDerivative.cpp:92
std::string
STL class.
Genfun::GaudiMathImplementation::NumericalDerivative::m_error
double m_error
Definition: NumericalDerivative.h:155
NumericalDerivative.h
Genfun::GaudiMathImplementation::NumericalDerivative::setType
const Type & setType(const Type &value)
change the type of the adaptive differentiation
Definition: NumericalDerivative.cpp:85
Genfun::GaudiMathImplementation::NumericalDerivative::Forward
Definition: NumericalDerivative.h:48
StatusCode
This class is used for returning status codes from appropriate routines.
Definition: StatusCode.h:26
GaudiException.h
gaudirun.type
type
Definition: gaudirun.py:125
Genfun::GaudiMathImplementation::NumericalDerivative::m_index
size_t m_index
Definition: NumericalDerivative.h:149
Genfun::GaudiMathImplementation::GSL_Adaptor
double GSL_Adaptor(double x, void *params)
Definition: Helpers.cpp:37
Genfun::GaudiMathImplementation::NumericalDerivative::m_DIM
size_t m_DIM
Definition: NumericalDerivative.h:150
Genfun::GaudiMathImplementation::NumericalDerivative
Numerical derivative (using GSL adaptive numerical differentiation)
Definition: NumericalDerivative.h:42
Genfun::GaudiMathImplementation::NumericalDerivative::m_type
Type m_type
Definition: NumericalDerivative.h:151
Genfun::GaudiMathImplementation::NumericalDerivative::Central
Definition: NumericalDerivative.h:47
Genfun::GaudiMathImplementation::NumericalDerivative::m_result
double m_result
Definition: NumericalDerivative.h:154
Helpers.h
Constant.h
Genfun::GaudiMathImplementation::NumericalDerivative::Backward
Definition: NumericalDerivative.h:49
Genfun::GaudiMathImplementation::NumericalDerivative::Exception
StatusCode Exception(const std::string &message, const StatusCode &sc=StatusCode::FAILURE) const
Exception.
Definition: NumericalDerivative.cpp:171
Genfun
CLHEP.
Definition: IEqSolver.h:13
Genfun::GaudiMathImplementation::GSL_Helper
the simple structure to be used for adaption interface Genfun::AbsFunction to gsl_function structure ...
Definition: Helpers.h:26
Genfun::GaudiMathImplementation::NumericalDerivative::result
double result() const
the result of the last call
Definition: NumericalDerivative.h:130