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