herwig is hosted by Hepforge, IPPP Durham
Herwig  7.1.5
ShowerAlphaQCD.h
1 // -*- C++ -*-
2 //
3 // ShowerAlphaQCD.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2017 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_ShowerAlphaQCD_H
10 #define HERWIG_ShowerAlphaQCD_H
11 //
12 // This is the declaration of the ShowerAlphaQCD class.
13 //
14 
15 #include "Herwig/Shower/Core/Couplings/ShowerAlpha.h"
16 
17 namespace Herwig {
18 
19 using namespace ThePEG;
20 
33 class ShowerAlphaQCD: public ShowerAlpha {
34 
35 public:
36 
41  _qmin(0.630882*GeV), _asType(1), _asMaxNP(1.0),
42  _thresholds(4), _lambda(4),
43  _nloop(3),_lambdaopt(false),_thresopt(false),
44  _lambdain(0.208364*GeV),_alphain(0.118),_inopt(true),_tolerance(1e-10),
45  _maxtry(100),_alphamin(0.),_val0(1.), _optInputScale(ZERO) {}
46 
47 public:
48 
59  virtual double value(const Energy2 scale) const;
60 
65  virtual double overestimateValue() const;
66 
70  virtual double ratio(const Energy2 scale,double factor =1.) const;
71 
75  virtual void initialize() { doinit(); }
76 
81  string value(string);
82 
88  string check(string args);
89 
91 
96  Energy lambdaQCD(unsigned int nf) {
97  if (nf <= 3) return _lambda[0];
98  else if (nf==4 || nf==5) return _lambda[nf-3];
99  else return _lambda[3];
100  }
101 
107  const vector<Energy>& quarkMasses() const { return _quarkMasses; }
108 
109 public:
110 
117  void persistentOutput(PersistentOStream & os) const;
118 
124  void persistentInput(PersistentIStream & is, int version);
126 
133  static void Init();
134 
135 protected:
136 
143  virtual IBPtr clone() const;
144 
149  virtual IBPtr fullclone() const;
151 
152 
153 protected:
154 
162  virtual void doinit();
164 
165 private:
166 
177  double alphaS(Energy q, Energy lam, int nf) const;
178 
185  double derivativealphaS(Energy q, Energy lam, int nf) const;
186 
195  Energy computeLambda(Energy match, double alpha, unsigned int nflav) const;
196 
202  pair<short, Energy> getLamNfTwoLoop(Energy q) const;
204 
205 private:
206 
211  ShowerAlphaQCD & operator=(const ShowerAlphaQCD &) = delete;
212 
213 private:
214 
218  Energy _qmin;
219 
224  int _asType;
225 
230  double _asMaxNP;
231 
235  vector<Energy> _thresholds;
236 
240  vector<Energy> _lambda;
241 
245  unsigned int _nloop;
246 
252 
256  bool _thresopt;
257 
261  Energy _lambdain;
262 
266  double _alphain;
267 
271  bool _inopt;
272 
276  double _tolerance;
277 
281  unsigned int _maxtry;
282 
286  double _alphamin;
287 
291  double _val0;
292 
298 
303  vector<Energy> _quarkMasses;
304 };
305 
306 }
307 
308 #endif /* HERWIG_ShowerAlphaQCD_H */
Energy _qmin
Minimum value of the scale.
bool _lambdaopt
Option for the translation between and .
double _asMaxNP
Another parameter, a possible (maximum) value of alpha in the non-perturbative region.
unsigned int _nloop
Option for the number of loops.
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
ShowerAlphaQCD()
The default constructor.
double _alphamin
The minimum value of the coupling.
vector< Energy > _quarkMasses
The quark masses to be used; if not empty these masses should be considered instead of the ones set i...
int _asType
Parameter controlling the behaviour of in the non-perturbative region.
Energy lambdaQCD(unsigned int nf)
Get the value of .
double _alphain
Input value of .
virtual void initialize()
Initialize this coupling.
double _val0
Value of at the minimum scale.
double _tolerance
Tolerance for discontinuities at the thresholds.
Energy _optInputScale
An optional input scale to be used for the input alphas; if zero MZ will be used out of the particle ...
const vector< Energy > & quarkMasses() const
Return the quark masses to be used; if not empty these masses should be considered instead of the one...
bool _inopt
Option for the calculation of Lambda from input parameters.
This concrete class provides the definition of the pure virtual function value() and overestimateValu...
This class is the abstract class from which all types of running couplings used in the Showering deri...
Definition: ShowerAlpha.h:49
bool _thresopt
Option for the threshold masses.
-*- C++ -*-
vector< Energy > _thresholds
Thresholds for the different number of flavours.
matcher< Density > match(const Generator< Density > &gen)
Indicate that the argument density should be matched to the previous one in a piecewise definition...
Energy _lambdain
Input value of Lambda.
vector< Energy > _lambda
for the different number of flavours
constexpr ZeroUnit ZERO
unsigned int _maxtry
Maximum number of iterations for the Newton-Raphson method to converge.