Logo ROOT   6.10/09
Reference Guide
GaussIntegrator.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 GaussIntegrator
12 //
13 // Created by: David Gonzalez Maline : Wed Jan 16 2008
14 //
15 
16 #ifndef ROOT_Math_GaussIntegrator
17 #define ROOT_Math_GaussIntegrator
18 
19 #include "Math/IFunction.h"
20 
21 #include "Math/VirtualIntegrator.h"
22 
23 
24 namespace ROOT {
25 namespace Math {
26 
27 
28 //___________________________________________________________________________________________
29 /**
30  User class for performing function integration.
31 
32  It will use the Gauss Method for function integration in a given interval.
33  This class is implemented from TF1::Integral().
34 
35  @ingroup Integration
36 
37  */
38 
40 
41 
42 public:
43 
44  /** Destructor */
45  virtual ~GaussIntegrator();
46 
47  /** Default Constructor.
48  If the tolerance are not given, use default values specified in ROOT::Math::IntegratorOneDimOptions
49  */
50  GaussIntegrator(double absTol = -1, double relTol = -1);
51 
52 
53  /** Static function: set the fgAbsValue flag.
54  By default TF1::Integral uses the original function value to compute the integral
55  However, TF1::Moment, CentralMoment require to compute the integral
56  using the absolute value of the function.
57  */
58  void AbsValue(bool flag);
59 
60 
61  // Implementing VirtualIntegrator Interface
62 
63  /** Set the desired relative Error. */
64  virtual void SetRelTolerance (double eps) { fEpsRel = eps; }
65 
66  /** This method is not implemented. */
67  virtual void SetAbsTolerance (double eps) { fEpsAbs = eps; }
68 
69  /** Returns the result of the last Integral calculation. */
70  double Result () const;
71 
72  /** Return the estimate of the absolute Error of the last Integral calculation. */
73  double Error () const;
74 
75  /** return the status of the last integration - 0 in case of success */
76  int Status () const;
77 
78  // Implementing VirtualIntegratorOneDim Interface
79 
80  /**
81  Returns Integral of function between a and b.
82  Based on original CERNLIB routine DGAUSS by Sigfried Kolbig
83  converted to C++ by Rene Brun
84 
85  This function computes, to an attempted specified accuracy, the value
86  of the integral.
87 
88  Method:
89  For any interval [a,b] we define g8(a,b) and g16(a,b) to be the 8-point
90  and 16-point Gaussian quadrature approximations to
91  \f[
92  I = \int^{b}_{a} f(x)dx
93  \f]
94  and define
95  \f[
96  r(a,b) = \frac{\left|g_{16}(a,b)-g_{8}(a,b)\right|}{1+\left|g_{16}(a,b)\right|}
97  \f]
98  Then,
99  \f[
100  G = \sum_{i=1}^{k}g_{16}(x_{i-1},x_{i})
101  \f]
102  where, starting with \f$x_{0} = A\f$ and finishing with \f$x_{k} = B\f$,
103  the subdivision points \f$x_{i}(i=1,2,...)\f$ are given by
104  \f[
105  x_{i} = x_{i-1} + \lambda(B-x_{i-1})
106  \f]
107  \f$\lambda\f$ is equal to the first member of the
108  sequence 1,1/2,1/4,... for which \f$r(x_{i-1}, x_{i}) < EPS\f$.
109  If, at any stage in the process of subdivision, the ratio
110  \f[
111  q = \left|\frac{x_{i}-x_{i-1}}{B-A}\right|
112  \f]
113  is so small that 1+0.005q is indistinguishable from 1 to
114  machine accuracy, an error exit occurs with the function value
115  set equal to zero.
116 
117  Accuracy:
118  The user provides absolute and relative error bounds (epsrel and epsabs) and the
119  algorithm will stop when the estimated error is less than the epsabs OR is less
120  than |I| * epsrel.
121  Unless there is severe cancellation of positive and negative values of
122  f(x) over the interval [A,B], the relative error may be considered as
123  specifying a bound on the <I>relative</I> error of I in the case
124  |I|&gt;1, and a bound on the absolute error in the case |I|&lt;1. More
125  precisely, if k is the number of sub-intervals contributing to the
126  approximation (see Method), and if
127  \f[
128  I_{abs} = \int^{B}_{A} \left|f(x)\right|dx
129  \f]
130  then the relation
131  \f[
132  \frac{\left|G-I\right|}{I_{abs}+k} < EPS
133  \f]
134  will nearly always be true, provided the routine terminates without
135  printing an error message. For functions f having no singularities in
136  the closed interval [A,B] the accuracy will usually be much higher than
137  this.
138 
139  Error handling:
140  The requested accuracy cannot be obtained (see Method).
141  The function value is set equal to zero.
142 
143  Note 1:
144  Values of the function f(x) at the interval end-points A and B are not
145  required. The subprogram may therefore be used when these values are
146  undefined
147  */
148  double Integral (double a, double b);
149 
150  /** Returns Integral of function on an infinite interval.
151  This function computes, to an attempted specified accuracy, the value of the integral:
152  \f[
153  I = \int^{\infty}_{-\infty} f(x)dx
154  \f]
155  Usage:
156  In any arithmetic expression, this function has the approximate value
157  of the integral I.
158 
159  The integral is mapped onto [0,1] using a transformation then integral computation is surrogated to DoIntegral.
160  */
161  double Integral ();
162 
163  /** Returns Integral of function on an upper semi-infinite interval.
164  This function computes, to an attempted specified accuracy, the value of the integral:
165  \f[
166  I = \int^{\infty}_{A} f(x)dx
167  \f]
168  Usage:
169  In any arithmetic expression, this function has the approximate value
170  of the integral I.
171  - A: lower end-point of integration interval.
172 
173  The integral is mapped onto [0,1] using a transformation then integral computation is surrogated to DoIntegral.
174  */
175  double IntegralUp (double a);
176 
177  /** Returns Integral of function on a lower semi-infinite interval.
178  This function computes, to an attempted specified accuracy, the value of the integral:
179  \f[
180  I = \int^{B}_{-\infty} f(x)dx
181  \f]
182  Usage:
183  In any arithmetic expression, this function has the approximate value
184  of the integral I.
185  - B: upper end-point of integration interval.
186 
187  The integral is mapped onto [0,1] using a transformation then integral computation is surrogated to DoIntegral.
188  */
189  double IntegralLow (double b);
190 
191 
192  /** Set integration function (flag control if function must be copied inside).
193  \@param f Function to be used in the calculations.
194  */
195  void SetFunction (const IGenFunction &);
196 
197  /** This method is not implemented. */
198  double Integral (const std::vector< double > &pts);
199 
200  /** This method is not implemented. */
201  double IntegralCauchy (double a, double b, double c);
202 
203  /// get the option used for the integration
205 
206  // set the options
207  virtual void SetOptions(const ROOT::Math::IntegratorOneDimOptions & opt);
208 
209 private:
210 
211  /**
212  Integration surrogate method. Return integral of passed function in interval [a,b]
213  Derived class (like GaussLegendreIntegrator) can re-implement this method to modify to use
214  an improved algorithm
215  */
216  virtual double DoIntegral (double a, double b, const IGenFunction* func);
217 
218 protected:
219 
220  static bool fgAbsValue; // AbsValue used for the calculation of the integral
221  double fEpsRel; // Relative error.
222  double fEpsAbs; // Absolute error.
223  bool fUsedOnce; // Bool value to check if the function was at least called once.
224  double fLastResult; // Result from the last estimation.
225  double fLastError; // Error from the last estimation.
226  const IGenFunction* fFunction; // Pointer to function used.
227 
228 };
229 
230 /**
231  Auxiliary inner class for mapping infinite and semi-infinite integrals
232 */
234 public:
235  enum ESemiInfinitySign {kMinus = -1, kPlus = +1};
236  IntegrandTransform(const IGenFunction* integrand);
237  IntegrandTransform(const double boundary, ESemiInfinitySign sign, const IGenFunction* integrand);
238  double operator()(double x) const;
239  double DoEval(double x) const;
240  IGenFunction* Clone() const;
241 private:
244  double fBoundary;
246  double DoEval(double x, double boundary, int sign) const;
247 };
248 
249 
250 
251 } // end namespace Math
252 
253 } // end namespace ROOT
254 
255 #endif /* ROOT_Math_GaussIntegrator */
Interface (abstract) class for 1D numerical integration It must be implemented by the concrate Integr...
virtual void SetAbsTolerance(double eps)
This method is not implemented.
double IntegralUp(double a)
Returns Integral of function on an upper semi-infinite interval.
double Error() const
Return the estimate of the absolute Error of the last Integral calculation.
Interface (abstract class) for generic functions objects of one-dimension Provides a method to evalua...
Definition: IFunction.h:134
const double absTol
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
double Integral()
Returns Integral of function on an infinite interval.
GaussIntegrator(double absTol=-1, double relTol=-1)
Default Constructor.
TArc * a
Definition: textangle.C:12
virtual double DoIntegral(double a, double b, const IGenFunction *func)
Integration surrogate method.
TRObject operator()(const T1 &t1) const
double IntegralLow(double b)
Returns Integral of function on a lower semi-infinite interval.
virtual void SetRelTolerance(double eps)
Set the desired relative Error.
Double_t x[n]
Definition: legend1.C:17
double IntegralCauchy(double a, double b, double c)
This method is not implemented.
int Status() const
return the status of the last integration - 0 in case of success
void AbsValue(bool flag)
Static function: set the fgAbsValue flag.
const IGenFunction * fIntegrand
virtual ~GaussIntegrator()
Destructor.
Numerical one dimensional integration options.
virtual void SetOptions(const ROOT::Math::IntegratorOneDimOptions &opt)
set the options (should be re-implemented by derived classes -if more options than tolerance exist ...
User class for performing function integration.
double Result() const
Returns the result of the last Integral calculation.
void SetFunction(const IGenFunction &)
Set integration function (flag control if function must be copied inside).
virtual ROOT::Math::IntegratorOneDimOptions Options() const
get the option used for the integration
double func(double *x, double *p)
Definition: stressTF1.cxx:213
Auxiliary inner class for mapping infinite and semi-infinite integrals.
Namespace for new Math classes and functions.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
const IGenFunction * fFunction