Logo ROOT   6.10/09
Reference Guide
AdaptiveIntegratorMultiDim.h
Go to the documentation of this file.
1 // @(#)root/mathcore:$Id$
2 // Author: M. Slawinska 08/2007
3 
4 /**********************************************************************
5  * *
6  * Copyright (c) 2007 , LCG ROOT MathLib Team *
7  * *
8  * *
9  **********************************************************************/
10 
11 // Header source file for class AdaptiveIntegratorMultiDim
12 
13 
14 #ifndef ROOT_Math_AdaptiveIntegratorMultiDim
15 #define ROOT_Math_AdaptiveIntegratorMultiDim
16 
17 #include "Math/IFunctionfwd.h"
18 
19 #include "Math/VirtualIntegrator.h"
20 
21 namespace ROOT {
22 namespace Math {
23 
24 /** \class AdaptiveIntegratorMultiDim
25  \ingroup Integration
26 
27 Class for adaptive quadrature integration in multi-dimensions using rectangular regions.
28 Algorithm from A.C. Genz, A.A. Malik, An adaptive algorithm for numerical integration over
29 an N-dimensional rectangular region, J. Comput. Appl. Math. 6 (1980) 295-302.
30 
31 Converted/adapted by R.Brun to C++ from Fortran CERNLIB routine RADMUL (D120)
32 The new code features many changes compared to the Fortran version.
33 
34 Control parameters are:
35 
36  - \f$ minpts \f$: Minimum number of function evaluations requested. Must not exceed maxpts.
37  if minpts < 1 minpts is set to \f$ 2^n +2n(n+1) +1 \f$ where n is the function dimension
38  - \f$ maxpts \f$: Maximum number of function evaluations to be allowed.
39  \f$ maxpts >= 2^n +2n(n+1) +1 \f$
40  if \f$ maxpts<minpts \f$, \f$ maxpts \f$ is set to \f$ 10minpts \f$
41  - \f$ epstol \f$, \f$ epsrel \f$ : Specified relative and absolute accuracy.
42 
43 The integral will stop if the relative error is less than relative tolerance OR the
44 absolute error is less than the absolute tolerance
45 
46 The class computes in addition to the integral of the function in the desired interval:
47 
48  - an estimation of the relative accuracy of the result.
49  - number of function evaluations performed.
50  - status code:
51  0. Normal exit. . At least minpts and at most maxpts calls to the function were performed.
52  1. maxpts is too small for the specified accuracy eps.
53  The result and relerr contain the values obtainable for the
54  specified value of maxpts.
55  2. size is too small for the specified number MAXPTS of function evaluations.
56  3. n<2 or n>15
57 
58 ### Method:
59 
60 An integration rule of degree seven is used together with a certain
61 strategy of subdivision.
62 For a more detailed description of the method see References.
63 
64 ### Notes:
65 
66  1..Multi-dimensional integration is time-consuming. For each rectangular
67  subregion, the routine requires function evaluations.
68  Careful programming of the integrand might result in substantial saving
69  of time.
70  2..Numerical integration usually works best for smooth functions.
71  Some analysis or suitable transformations of the integral prior to
72  numerical work may contribute to numerical efficiency.
73 
74 ### References:
75 
76  1. A.C. Genz and A.A. Malik, Remarks on algorithm 006:
77  An adaptive algorithm for numerical integration over
78  an N-dimensional rectangular region, J. Comput. Appl. Math. 6 (1980) 295-302.
79  2. A. van Doren and L. de Ridder, An adaptive algorithm for numerical
80  integration over an n-dimensional cube, J.Comput. Appl. Math. 2 (1976) 207-217.
81 
82 */
83 
85 
86 public:
87 
88  /**
89  Construct given optionally tolerance (absolute and relative), maximum number of function evaluation (maxpts) and
90  size of the working array.
91  The integration will stop when the absolute error is less than the absolute tolerance OR when the relative error is less
92  than the relative tolerance. The absolute tolerance by default is not used (it is equal to zero).
93  The size of working array represents the number of sub-division used for calculating the integral.
94  Higher the dimension, larger sizes are required for getting the same accuracy.
95  The size must be larger than
96  \f$ (2n + 3) (1 + maxpts/(2^n + 2n(n + 1) + 1))/2) \f$.
97  For smaller value passed, the minimum allowed will be used
98  */
99  explicit
100  AdaptiveIntegratorMultiDim(double absTol = 0.0, double relTol = 1E-9, unsigned int maxpts = 100000, unsigned int size = 0);
101 
102  /**
103  Construct with a reference to the integrand function and given optionally
104  tolerance (absolute and relative), maximum number of function evaluation (maxpts) and
105  size of the working array.
106  */
107  explicit
108  AdaptiveIntegratorMultiDim(const IMultiGenFunction &f, double absTol = 0.0, double relTol = 1E-9, unsigned int maxcall = 100000, unsigned int size = 0);
109 
110  /**
111  destructor (no operations)
112  */
114 
115 
116  /**
117  evaluate the integral with the previously given function between xmin[] and xmax[]
118  */
119  double Integral(const double* xmin, const double * xmax) {
120  return DoIntegral(xmin,xmax, false);
121  }
122 
123 
124  /// evaluate the integral passing a new function
125  double Integral(const IMultiGenFunction &f, const double* xmin, const double * xmax);
126 
127  /// set the integration function (must implement multi-dim function interface: IBaseFunctionMultiDim)
128  void SetFunction(const IMultiGenFunction &f);
129 
130  /// return result of integration
131  double Result() const { return fResult; }
132 
133  /// return integration error
134  double Error() const { return fError; }
135 
136  /// return relative error
137  double RelError() const { return fRelError; }
138 
139  /// return status of integration
140  /// - status = 0 successful integration
141  /// - status = 1
142  /// MAXPTS is too small for the specified accuracy EPS.
143  /// The result contain the values
144  /// obtainable for the specified value of MAXPTS.
145  /// - status = 2
146  /// size is too small for the specified number MAXPTS of function evaluations.
147  /// - status = 3
148  /// wrong dimension , N<2 or N > 15. Returned result and error are zero
149  int Status() const { return fStatus; }
150 
151  /// return number of function evaluations in calculating the integral
152  int NEval() const { return fNEval; }
153 
154  /// set relative tolerance
155  void SetRelTolerance(double relTol);
156 
157  /// set absolute tolerance
158  void SetAbsTolerance(double absTol);
159 
160  ///set workspace size
161  void SetSize(unsigned int size) { fSize = size; }
162 
163  ///set min points
164  void SetMinPts(unsigned int n) { fMinPts = n; }
165 
166  ///set max points
167  void SetMaxPts(unsigned int n) { fMaxPts = n; }
168 
169  /// set the options
171 
172  /// get the option used for the integration
174 
175 protected:
176 
177  // internal function to compute the integral (if absVal is true compute abs value of function integral
178  double DoIntegral(const double* xmin, const double * xmax, bool absVal = false);
179 
180  private:
181 
182  unsigned int fDim; // dimensionality of integrand
183  unsigned int fMinPts; // minimum number of function evaluation requested
184  unsigned int fMaxPts; // maximum number of function evaluation requested
185  unsigned int fSize; // max size of working array (explode with dimension)
186  double fAbsTol; // absolute tolerance
187  double fRelTol; // relative tolerance
188 
189  double fResult; // last integration result
190  double fError; // integration error
191  double fRelError; // Relative error
192  int fNEval; // number of function evaluation
193  int fStatus; // status of algorithm (error if not zero)
194 
195  const IMultiGenFunction* fFun; // pointer to integrand function
196 
197 };
198 
199 }//namespace Math
200 }//namespace ROOT
201 
202 #endif /* ROOT_Math_AdaptiveIntegratorMultiDim */
double DoIntegral(const double *xmin, const double *xmax, bool absVal=false)
float xmin
Definition: THbookFile.cxx:93
const double absTol
void SetFunction(const IMultiGenFunction &f)
set the integration function (must implement multi-dim function interface: IBaseFunctionMultiDim) ...
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
ROOT::Math::IntegratorMultiDimOptions Options() const
get the option used for the integration
int Status() const
return status of integration
void SetMaxPts(unsigned int n)
set max points
void SetAbsTolerance(double absTol)
set absolute tolerance
virtual ~AdaptiveIntegratorMultiDim()
destructor (no operations)
double Result() const
return result of integration
void SetRelTolerance(double relTol)
set relative tolerance
int NEval() const
return number of function evaluations in calculating the integral
double Integral(const double *xmin, const double *xmax)
evaluate the integral with the previously given function between xmin[] and xmax[] ...
Documentation for the abstract class IBaseFunctionMultiDim.
Definition: IFunction.h:62
Numerical multi dimensional integration options.
Interface (abstract) class for multi numerical integration It must be implemented by the concrete Int...
double Error() const
return integration error
float xmax
Definition: THbookFile.cxx:93
constexpr Double_t E()
Definition: TMath.h:74
double f(double x)
Namespace for new Math classes and functions.
void SetSize(unsigned int size)
set workspace size
double RelError() const
return relative error
void SetOptions(const ROOT::Math::IntegratorMultiDimOptions &opt)
set the options
const Int_t n
Definition: legend1.C:16
AdaptiveIntegratorMultiDim(double absTol=0.0, double relTol=1E-9, unsigned int maxpts=100000, unsigned int size=0)
Construct given optionally tolerance (absolute and relative), maximum number of function evaluation (...
Class for adaptive quadrature integration in multi-dimensions using rectangular regions.
void SetMinPts(unsigned int n)
set min points