|
||||
File indexing completed on 2025-01-18 10:10:21
0001 // @(#)root/mathcore:$Id$ 0002 // Authors: David Gonzalez Maline 01/2008 0003 0004 /********************************************************************** 0005 * * 0006 * Copyright (c) 2006 , LCG ROOT MathLib Team * 0007 * * 0008 * * 0009 **********************************************************************/ 0010 0011 // Header file for RichardsonDerivator 0012 // 0013 // Created by: David Gonzalez Maline : Mon Feb 4 2008 0014 // 0015 0016 #ifndef ROOT_Math_RichardsonDerivator 0017 #define ROOT_Math_RichardsonDerivator 0018 0019 #include <Math/IFunction.h> 0020 0021 /** 0022 @defgroup Deriv Numerical Differentiation 0023 Classes for numerical differentiation 0024 @ingroup NumAlgo 0025 */ 0026 0027 0028 namespace ROOT { 0029 namespace Math { 0030 0031 //___________________________________________________________________________________________ 0032 /** 0033 User class for calculating the derivatives of a function. It can calculate first (method Derivative1), 0034 second (method Derivative2) and third (method Derivative3) of a function. 0035 0036 It uses the Richardson extrapolation method for function derivation in a given interval. 0037 The method use 2 derivative estimates (one computed with step h and one computed with step h/2) 0038 to compute a third, more accurate estimation. It is equivalent to the 0039 <a href = http://en.wikipedia.org/wiki/Five-point_stencil>5-point method</a>, 0040 which can be obtained with a Taylor expansion. 0041 A step size should be given, depending on x and f(x). 0042 An optimal step size value minimizes the truncation error of the expansion and the rounding 0043 error in evaluating x+h and f(x+h). A too small h will yield a too large rounding error while a too large 0044 h will give a large truncation error in the derivative approximation. 0045 A good discussion can be found in discussed in 0046 <a href=http://www.nrbook.com/a/bookcpdf/c5-7.pdf>Chapter 5.7</a> of Numerical Recipes in C. 0047 By default a value of 0.001 is uses, acceptable in many cases. 0048 0049 This class is implemented using code previously in TF1::Derivate{,2,3}(). Now TF1 uses this class. 0050 0051 @ingroup Deriv 0052 0053 */ 0054 0055 class RichardsonDerivator { 0056 public: 0057 0058 /** Destructor: Removes function if needed. */ 0059 ~RichardsonDerivator(); 0060 0061 /** Default Constructor. 0062 Give optionally the step size for derivation. By default is 0.001, which is fine for x ~ 1 0063 Increase if x is in average larger or decrease if x is smaller 0064 */ 0065 RichardsonDerivator(double h = 0.001); 0066 0067 /** Construct from function and step size 0068 */ 0069 RichardsonDerivator(const ROOT::Math::IGenFunction & f, double h = 0.001, bool copyFunc = false); 0070 0071 /** 0072 Copy constructor 0073 */ 0074 RichardsonDerivator(const RichardsonDerivator & rhs); 0075 0076 /** 0077 Assignment operator 0078 */ 0079 RichardsonDerivator & operator= ( const RichardsonDerivator & rhs); 0080 0081 0082 /** Returns the estimate of the absolute Error of the last derivative calculation. */ 0083 double Error () const { return fLastError; } 0084 0085 0086 /** 0087 Returns the first derivative of the function at point x, 0088 computed by Richardson's extrapolation method (use 2 derivative estimates 0089 to compute a third, more accurate estimation) 0090 first, derivatives with steps h and h/2 are computed by central difference formulas 0091 \f[ 0092 D(h) = \frac{f(x+h) - f(x-h)}{2h} 0093 \f] 0094 the final estimate 0095 \f[ 0096 D = \frac{4D(h/2) - D(h)}{3} 0097 \f] 0098 "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition" 0099 0100 the argument eps may be specified to control the step size (precision). 0101 the step size is taken as eps*(xmax-xmin). 0102 the default value (0.001) should be good enough for the vast majority 0103 of functions. Give a smaller value if your function has many changes 0104 of the second derivative in the function range. 0105 0106 Getting the error via TF1::DerivativeError: 0107 (total error = roundoff error + interpolation error) 0108 the estimate of the roundoff error is taken as follows: 0109 \f[ 0110 err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}}, 0111 \f] 0112 where k is the double precision, ai are coefficients used in 0113 central difference formulas 0114 interpolation error is decreased by making the step size h smaller. 0115 */ 0116 double Derivative1 (double x) { return Derivative1(*fFunction,x,fStepSize); } 0117 double operator() (double x) { return Derivative1(*fFunction,x,fStepSize); } 0118 0119 /** 0120 First Derivative calculation passing function object and step-size 0121 */ 0122 double Derivative1(const IGenFunction & f, double x, double h); 0123 0124 /// Computation of the first derivative using a forward formula 0125 double DerivativeForward(double x) { 0126 return DerivativeForward(*fFunction, x, fStepSize); 0127 } 0128 0129 /// Computation of the first derivative using a forward formula 0130 double DerivativeForward(const IGenFunction &f, double x, double h); 0131 0132 /// Computation of the first derivative using a backward formula 0133 double DerivativeBackward(double x) { 0134 return DerivativeBackward(*fFunction, x, fStepSize); 0135 } 0136 0137 /// Computation of the first derivative using a forward formula 0138 double DerivativeBackward(const IGenFunction &f, double x, double h) { 0139 return DerivativeForward(f, x, -h); 0140 } 0141 0142 /** 0143 Returns the second derivative of the function at point x, 0144 computed by Richardson's extrapolation method (use 2 derivative estimates 0145 to compute a third, more accurate estimation) 0146 first, derivatives with steps h and h/2 are computed by central difference formulas 0147 \f[ 0148 D(h) = \frac{f(x+h) - 2f(x) + f(x-h)}{h^{2}} 0149 \f] 0150 the final estimate 0151 \f[ 0152 D = \frac{4D(h/2) - D(h)}{3} 0153 \f] 0154 "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition" 0155 0156 the argument eps may be specified to control the step size (precision). 0157 the step size is taken as eps*(xmax-xmin). 0158 the default value (0.001) should be good enough for the vast majority 0159 of functions. Give a smaller value if your function has many changes 0160 of the second derivative in the function range. 0161 0162 Getting the error via TF1::DerivativeError: 0163 (total error = roundoff error + interpolation error) 0164 the estimate of the roundoff error is taken as follows: 0165 \f[ 0166 err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}}, 0167 \f] 0168 where k is the double precision, ai are coefficients used in 0169 central difference formulas 0170 interpolation error is decreased by making the step size h smaller. 0171 */ 0172 double Derivative2 (double x) { 0173 return Derivative2( *fFunction, x, fStepSize); 0174 } 0175 0176 /** 0177 Second Derivative calculation passing function and step-size 0178 */ 0179 double Derivative2(const IGenFunction & f, double x, double h); 0180 0181 /** 0182 Returns the third derivative of the function at point x, 0183 computed by Richardson's extrapolation method (use 2 derivative estimates 0184 to compute a third, more accurate estimation) 0185 first, derivatives with steps h and h/2 are computed by central difference formulas 0186 \f[ 0187 D(h) = \frac{f(x+2h) - 2f(x+h) + 2f(x-h) - f(x-2h)}{2h^{3}} 0188 \f] 0189 the final estimate 0190 \f[ 0191 D = \frac{4D(h/2) - D(h)}{3} 0192 \f] 0193 "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition" 0194 0195 the argument eps may be specified to control the step size (precision). 0196 the step size is taken as eps*(xmax-xmin). 0197 the default value (0.001) should be good enough for the vast majority 0198 of functions. Give a smaller value if your function has many changes 0199 of the second derivative in the function range. 0200 0201 Getting the error via TF1::DerivativeError: 0202 (total error = roundoff error + interpolation error) 0203 the estimate of the roundoff error is taken as follows: 0204 \f[ 0205 err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}}, 0206 \f] 0207 where k is the double precision, ai are coefficients used in 0208 central difference formulas 0209 interpolation error is decreased by making the step size h smaller. 0210 */ 0211 double Derivative3 (double x) { 0212 return Derivative3(*fFunction, x, fStepSize); 0213 } 0214 0215 /** 0216 Third Derivative calculation passing function and step-size 0217 */ 0218 double Derivative3(const IGenFunction & f, double x, double h); 0219 0220 /** Set function for derivative calculation (copy the function if option has been enabled in the constructor) 0221 0222 \@param f Function to be differentiated 0223 */ 0224 void SetFunction (const IGenFunction & f); 0225 0226 /** Set step size for derivative calculation 0227 0228 \@param h step size for calculation 0229 */ 0230 void SetStepSize (double h) { fStepSize = h; } 0231 0232 protected: 0233 0234 bool fFunctionCopied; ///< flag to control if function is copied in the class 0235 double fStepSize; ///< step size used for derivative calculation 0236 double fLastError; ///< error estimate of last derivative calculation 0237 const IGenFunction* fFunction; ///< pointer to function 0238 0239 }; 0240 0241 } // end namespace Math 0242 0243 } // end namespace ROOT 0244 0245 #endif /* ROOT_Math_RichardsonDerivator */ 0246
[ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |