herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
GaussianIntegrator.h
1 // -*- C++ -*-
2 //
3 // GaussianIntegrator.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2019 The Herwig Collaboration
5 //
6 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef HERWIG_GaussianIntegrator_H
10 #define HERWIG_GaussianIntegrator_H
11 //
12 // This is the declaration of the GaussianIntegrator class.
13 //
14 
15 #include "ThePEG/Pointer/ReferenceCounted.h"
16 #include "ThePEG/Repository/CurrentGenerator.h"
17 #include <vector>
18 
19 namespace Herwig {
20 
21 using namespace ThePEG;
22 
34 
35 public:
36 
43  : _abserr(1.E-35), _relerr(5.E-5), _binwidth(1.E-5), _maxeval(100000) {
44  // setup the weights and abscissae
45  Init();
46  }
47 
55  GaussianIntegrator(double abserr, double relerr, double binwidth,
56  int maxeval)
57  : _abserr(abserr), _relerr(relerr),
58  _binwidth(binwidth),
59  _maxeval(maxeval) {
60  // setup the weights and abscissae
61  Init();
62  }
63 
65  template <class T>
66  using ValT = decltype(std::declval<typename T::ValType>()
67  * std::declval<typename T::ArgType>());
68 
74  template <class T>
75  inline ValT<T> value(const T &,
76  const typename T::ArgType lower,
77  const typename T::ArgType upper) const;
78 
79 private:
80 
85  GaussianIntegrator & operator=(const GaussianIntegrator &) = delete;
86 
90  void Init();
91 
92 private:
93 
97  std::vector< std::vector<double> > _weights;
98 
102  std::vector< std::vector <double> > _abscissae;
103 
107  double _abserr,_relerr;
108 
112  double _binwidth;
113 
117  int _maxeval;
118 
119 };
120 
121 }
122 
123 #include "GaussianIntegrator.tcc"
124 
125 #endif /* HERWIG_GaussianIntegrator_H */
double _binwidth
The minimum width of a bin as a fraction of the integration region.
GaussianIntegrator()
Default Constructor.
std::vector< std::vector< double > > _abscissae
The abscissae.
std::vector< std::vector< double > > _weights
The weights for the gaussian quadrature.
decltype(std::declval< typename T::ValType >() *std::declval< typename T::ArgType >()) ValT
helper type for the integration result
double _abserr
The parameters controlling the error.
int _maxeval
Maximum number of function evaluations.
-*- C++ -*-
GaussianIntegrator(double abserr, double relerr, double binwidth, int maxeval)
Specify all the parameters.