herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
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
19namespace Herwig {
20
21using namespace ThePEG;
22
34
35public:
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
79private:
80
86
90 void Init();
91
92private:
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
118
119};
120
121}
122
123#include "GaussianIntegrator.tcc"
124
125#endif /* HERWIG_GaussianIntegrator_H */
GaussianIntegrator & operator=(const GaussianIntegrator &)=delete
The assignment operator is private and must never be called.
GaussianIntegrator()
Default Constructor.
std::vector< std::vector< double > > _weights
The weights for the gaussian quadrature.
std::vector< std::vector< double > > _abscissae
The abscissae.
ValT< T > value(const T &, const typename T::ArgType lower, const typename T::ArgType upper) const
The value of the integral.
double _abserr
The parameters controlling the error.
void Init()
Initialise the weights and abscissae.
int _maxeval
Maximum number of function evaluations.
GaussianIntegrator(double abserr, double relerr, double binwidth, int maxeval)
Specify all the parameters.
decltype(std::declval< typename T::ValType >() *std::declval< typename T::ArgType >()) ValT
helper type for the integration result
double _binwidth
The minimum width of a bin as a fraction of the integration region.
-*- C++ -*-