Logo ROOT   6.12/07
Reference Guide
RichardsonDerivator.h
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Authors: David Gonzalez Maline 01/2008
3 
4 /**********************************************************************
5  * *
6  * Copyright (c) 2006 , LCG ROOT MathLib Team *
7  * *
8  * *
9  **********************************************************************/
10 
11 // Header file for RichardsonDerivator
12 //
13 // Created by: David Gonzalez Maline : Mon Feb 4 2008
14 //
15 
16 #ifndef ROOT_Math_RichardsonDerivator
17 #define ROOT_Math_RichardsonDerivator
18 
19 #include <Math/IFunction.h>
20 
21 /**
22  @defgroup Deriv Numerical Differentiation
23  Classes for numerical differentiation
24  @ingroup NumAlgo
25 */
26 
27 
28 namespace ROOT {
29 namespace Math {
30 
31 //___________________________________________________________________________________________
32 /**
33  User class for calculating the derivatives of a function. It can calculate first (method Derivative1),
34  second (method Derivative2) and third (method Derivative3) of a function.
35 
36  It uses the Richardson extrapolation method for function derivation in a given interval.
37  The method use 2 derivative estimates (one computed with step h and one computed with step h/2)
38  to compute a third, more accurate estimation. It is equivalent to the
39  <a href = http://en.wikipedia.org/wiki/Five-point_stencil>5-point method</a>,
40  which can be obtained with a Taylor expansion.
41  A step size should be given, depending on x and f(x).
42  An optimal step size value minimizes the truncation error of the expansion and the rounding
43  error in evaluating x+h and f(x+h). A too small h will yield a too large rounding error while a too large
44  h will give a large truncation error in the derivative approximation.
45  A good discussion can be found in discussed in
46  <a href=http://www.nrbook.com/a/bookcpdf/c5-7.pdf>Chapter 5.7</a> of Numerical Recipes in C.
47  By default a value of 0.001 is uses, acceptable in many cases.
48 
49  This class is implemented using code previously in TF1::Derivate{,2,3}(). Now TF1 uses this class.
50 
51  @ingroup Deriv
52 
53  */
54 
56 public:
57 
58  /** Destructor: Removes function if needed. */
60 
61  /** Default Constructor.
62  Give optionally the step size for derivation. By default is 0.001, which is fine for x ~ 1
63  Increase if x is in average larger or decrease if x is smaller
64  */
65  RichardsonDerivator(double h = 0.001);
66 
67  /** Construct from function and step size
68  */
69  RichardsonDerivator(const ROOT::Math::IGenFunction & f, double h = 0.001, bool copyFunc = false);
70 
71  /**
72  Copy constructor
73  */
75 
76  /**
77  Assignment operator
78  */
80 
81 
82  /** Returns the estimate of the absolute Error of the last derivative calculation. */
83  double Error () const { return fLastError; }
84 
85 
86  /**
87  Returns the first derivative of the function at point x,
88  computed by Richardson's extrapolation method (use 2 derivative estimates
89  to compute a third, more accurate estimation)
90  first, derivatives with steps h and h/2 are computed by central difference formulas
91  \f[
92  D(h) = \frac{f(x+h) - f(x-h)}{2h}
93  \f]
94  the final estimate
95  \f[
96  D = \frac{4D(h/2) - D(h)}{3}
97  \f]
98  "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
99 
100  the argument eps may be specified to control the step size (precision).
101  the step size is taken as eps*(xmax-xmin).
102  the default value (0.001) should be good enough for the vast majority
103  of functions. Give a smaller value if your function has many changes
104  of the second derivative in the function range.
105 
106  Getting the error via TF1::DerivativeError:
107  (total error = roundoff error + interpolation error)
108  the estimate of the roundoff error is taken as follows:
109  \f[
110  err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
111  \f]
112  where k is the double precision, ai are coefficients used in
113  central difference formulas
114  interpolation error is decreased by making the step size h smaller.
115  */
116  double Derivative1 (double x) { return Derivative1(*fFunction,x,fStepSize); }
117  double operator() (double x) { return Derivative1(*fFunction,x,fStepSize); }
118 
119  /**
120  First Derivative calculation passing function object and step-size
121  */
122  double Derivative1(const IGenFunction & f, double x, double h);
123 
124  /// Computation of the first derivative using a forward formula
125  double DerivativeForward(double x) {
127  }
128 
129  /// Computation of the first derivative using a forward formula
130  double DerivativeForward(const IGenFunction &f, double x, double h);
131 
132  /// Computation of the first derivative using a backward formula
133  double DerivativeBackward(double x) {
135  }
136 
137  /// Computation of the first derivative using a forward formula
138  double DerivativeBackward(const IGenFunction &f, double x, double h) {
139  return DerivativeForward(f, x, -h);
140  }
141 
142  /**
143  Returns the second derivative of the function at point x,
144  computed by Richardson's extrapolation method (use 2 derivative estimates
145  to compute a third, more accurate estimation)
146  first, derivatives with steps h and h/2 are computed by central difference formulas
147  \f[
148  D(h) = \frac{f(x+h) - 2f(x) + f(x-h)}{h^{2}}
149  \f]
150  the final estimate
151  \f[
152  D = \frac{4D(h/2) - D(h)}{3}
153  \f]
154  "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
155 
156  the argument eps may be specified to control the step size (precision).
157  the step size is taken as eps*(xmax-xmin).
158  the default value (0.001) should be good enough for the vast majority
159  of functions. Give a smaller value if your function has many changes
160  of the second derivative in the function range.
161 
162  Getting the error via TF1::DerivativeError:
163  (total error = roundoff error + interpolation error)
164  the estimate of the roundoff error is taken as follows:
165  \f[
166  err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
167  \f]
168  where k is the double precision, ai are coefficients used in
169  central difference formulas
170  interpolation error is decreased by making the step size h smaller.
171  */
172  double Derivative2 (double x) {
173  return Derivative2( *fFunction, x, fStepSize);
174  }
175 
176  /**
177  Second Derivative calculation passing function and step-size
178  */
179  double Derivative2(const IGenFunction & f, double x, double h);
180 
181  /**
182  Returns the third derivative of the function at point x,
183  computed by Richardson's extrapolation method (use 2 derivative estimates
184  to compute a third, more accurate estimation)
185  first, derivatives with steps h and h/2 are computed by central difference formulas
186  \f[
187  D(h) = \frac{f(x+2h) - 2f(x+h) + 2f(x-h) - f(x-2h)}{2h^{3}}
188  \f]
189  the final estimate
190  \f[
191  D = \frac{4D(h/2) - D(h)}{3}
192  \f]
193  "Numerical Methods for Scientists and Engineers", H.M.Antia, 2nd edition"
194 
195  the argument eps may be specified to control the step size (precision).
196  the step size is taken as eps*(xmax-xmin).
197  the default value (0.001) should be good enough for the vast majority
198  of functions. Give a smaller value if your function has many changes
199  of the second derivative in the function range.
200 
201  Getting the error via TF1::DerivativeError:
202  (total error = roundoff error + interpolation error)
203  the estimate of the roundoff error is taken as follows:
204  \f[
205  err = k\sqrt{f(x)^{2} + x^{2}deriv^{2}}\sqrt{\sum ai^{2}},
206  \f]
207  where k is the double precision, ai are coefficients used in
208  central difference formulas
209  interpolation error is decreased by making the step size h smaller.
210  */
211  double Derivative3 (double x) {
212  return Derivative3(*fFunction, x, fStepSize);
213  }
214 
215  /**
216  Third Derivative calculation passing function and step-size
217  */
218  double Derivative3(const IGenFunction & f, double x, double h);
219 
220  /** Set function for derivative calculation (copy the function if option has been enabled in the constructor)
221 
222  \@param f Function to be differentiated
223  */
224  void SetFunction (const IGenFunction & f);
225 
226  /** Set step size for derivative calculation
227 
228  \@param h step size for calculation
229  */
230  void SetStepSize (double h) { fStepSize = h; }
231 
232 protected:
233 
234  bool fFunctionCopied; // flag to control if function is copied in the class
235  double fStepSize; // step size used for derivative calculation
236  double fLastError; // error estimate of last derivative calculation
237  const IGenFunction* fFunction; // pointer to function
238 
239 };
240 
241 } // end namespace Math
242 
243 } // end namespace ROOT
244 
245 #endif /* ROOT_Math_RichardsonDerivator */
246 
double Derivative3(double x)
Returns the third derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
Interface (abstract class) for generic functions objects of one-dimension Provides a method to evalua...
Definition: IFunction.h:135
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
double DerivativeForward(double x)
Computation of the first derivative using a forward formula.
TH1 * h
Definition: legend2.C:5
double Derivative2(double x)
Returns the second derivative of the function at point x, computed by Richardson&#39;s extrapolation meth...
Double_t x[n]
Definition: legend1.C:17
void SetFunction(const IGenFunction &f)
Set function for derivative calculation (copy the function if option has been enabled in the construc...
double Derivative1(double x)
Returns the first derivative of the function at point x, computed by Richardson&#39;s extrapolation metho...
double DerivativeBackward(double x)
Computation of the first derivative using a backward formula.
~RichardsonDerivator()
Destructor: Removes function if needed.
double DerivativeBackward(const IGenFunction &f, double x, double h)
Computation of the first derivative using a forward formula.
RichardsonDerivator(double h=0.001)
Default Constructor.
Namespace for new Math classes and functions.
RichardsonDerivator & operator=(const RichardsonDerivator &rhs)
Assignment operator.
void SetStepSize(double h)
Set step size for derivative calculation.
double Error() const
Returns the estimate of the absolute Error of the last derivative calculation.
User class for calculating the derivatives of a function.