herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
GSLIntegrator.h
1// -*- C++ -*-
2//
3// GSLIntegrator.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_GSLIntegrator_H
10#define HERWIG_GSLIntegrator_H
11//
12// This is the declaration of the GSLIntegrator class.
13//
14
15#include "ThePEG/Pointer/ReferenceCounted.h"
16#include "ThePEG/Repository/CurrentGenerator.h"
17#include "gsl/gsl_integration.h"
18#include "gsl/gsl_errno.h"
19
20namespace Herwig {
21
22using namespace ThePEG;
23
41
42public:
43
49 GSLIntegrator() : _abserr(1.0E-35), _relerr(1.0E-3), _nbins(1000) {}
50
57 GSLIntegrator(double abserr, double relerr, int nbins) :
58 _abserr(abserr), _relerr(relerr), _nbins(nbins) {}
60
62 template <class T>
63 using ValT = decltype(std::declval<typename T::ValType>()
64 * std::declval<typename T::ArgType>());
65
72 template <class T>
73 inline ValT<T>
74 value(const T & function,
75 const typename T::ArgType lower,
76 const typename T::ArgType upper) const;
77
85 template <class T>
86 inline ValT<T>
87 value(const T & function,
88 const typename T::ArgType lower,
89 const typename T::ArgType upper,
90 ValT<T> & error) const;
91
92private:
93
99
100private:
101
105 double _abserr;
106
110 double _relerr;
111
116};
117
118}
119
120#include "GSLIntegrator.tcc"
121
122#endif /* HERWIG_GSLIntegrator_H */
This class is designed to integrate a given function between 2 limits using the gsl QAGS integration ...
Definition: GSLIntegrator.h:40
decltype(std::declval< typename T::ValType >() *std::declval< typename T::ArgType >()) ValT
helper type for the integration result
Definition: GSLIntegrator.h:64
int _nbins
The maximum number of intervals to use.
ValT< T > value(const T &function, const typename T::ArgType lower, const typename T::ArgType upper) const
The value of the integral.
double _relerr
The parameters controlling the relative error.
GSLIntegrator()
Default Constructor uses values in GSL manual as parameters.
Definition: GSLIntegrator.h:49
double _abserr
The parameters controlling the absolute error.
ValT< T > value(const T &function, const typename T::ArgType lower, const typename T::ArgType upper, ValT< T > &error) const
The value of the integral.
GSLIntegrator(double abserr, double relerr, int nbins)
Specify all the parameters.
Definition: GSLIntegrator.h:57
GSLIntegrator & operator=(const GSLIntegrator &)=delete
The assignment operator is private and must never be called.
-*- C++ -*-