herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
GSLBisection.h
1 // -*- C++ -*-
2 //
3 // GSLBisection.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_GSLBisection_H
10 #define HERWIG_GSLBisection_H
11 //
12 // This is the declaration of the GSLBisection class.
13 //
14 
15 #include "ThePEG/Pointer/ReferenceCounted.h"
16 #include "Herwig/Utilities/GSLHelper.h"
17 #include <gsl/gsl_errno.h>
18 #include <gsl/gsl_math.h>
19 #include <gsl/gsl_roots.h>
20 
21 namespace Herwig {
22 
23 using namespace ThePEG;
24 
52 
53 public:
54 
58  struct GSLerror {};
59 
63  struct IntervalError {};
64 
70  GSLBisection() : abserr_(0), relerr_(1.E-8), maxPoints_(100) {}
71 
78  inline GSLBisection(double abserr, double relerr, int max) :
79  abserr_(abserr), relerr_(relerr), maxPoints_(max) {}
80 
82 
86  static void GSLsubstHandler(const char *, const char *,
87  int, int){
88  throw GSLerror();
89  }
90 
97  template <class T>
98  inline typename T::ArgType value(const T & function,
99  const typename T::ArgType lower,
100  const typename T::ArgType upper) const;
101 
102 private:
103 
108  GSLBisection & operator=(const GSLBisection &) = delete;
109 
110 private:
111 
115  double abserr_;
116 
120  double relerr_;
121 
126 };
127 
128 }
129 
130 #include "GSLBisection.tcc"
131 
132 #endif /* HERWIG_GSLBisection_H */
double abserr_
The parameters controlling the absolute error.
Definition: GSLBisection.h:115
int maxPoints_
The maximum number of evaluations to use.
Definition: GSLBisection.h:125
Struct that is used to throw and catch GSL errors.
Definition: GSLBisection.h:58
static void GSLsubstHandler(const char *, const char *, int, int)
Function to overwrite the default GSL error handling.
Definition: GSLBisection.h:86
GSLBisection(double abserr, double relerr, int max)
Specify all the parameters.
Definition: GSLBisection.h:78
double relerr_
The parameters controlling the relatve error.
Definition: GSLBisection.h:120
GSLBisection()
Default Constructor.
Definition: GSLBisection.h:70
This class is designed to find the root of a given function between 2 limits using bisection methods...
Definition: GSLBisection.h:51
Struct that is used to throw and catch GSL errors.
Definition: GSLBisection.h:63
-*- C++ -*-