herwig is hosted by Hepforge, IPPP Durham
Herwig 7.3.0
ShowerAlphaQCD.h
1// -*- C++ -*-
2//
3// ShowerAlphaQCD.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_ShowerAlphaQCD_H
10#define HERWIG_ShowerAlphaQCD_H
11//
12// This is the declaration of the ShowerAlphaQCD class.
13//
14
15#include "Herwig/Shower/ShowerAlpha.h"
16
17namespace Herwig {
18
19using namespace ThePEG;
20
34
35public:
36
41 _qmin(0.630882*GeV), _asType(1), _asMaxNP(1.0),
42 _thresholds(4), _lambda(4),
43 _nloop(3),_thresopt(false),
44 _alphain(0.118),_tolerance(1e-10),
45 _maxtry(100),_alphamin(0.),_val0(1.), _optInputScale(91.1876_GeV) {}
46
47public:
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
109public:
110
118
124 void persistentInput(PersistentIStream & is, int version);
126
133 static void Init();
134
135protected:
136
143 virtual IBPtr clone() const;
144
149 virtual IBPtr fullclone() const;
151
152
153protected:
154
162 virtual void doinit();
164
165private:
166
179 Energy computeLambda(Energy match, double alpha, unsigned int nflav) const;
180
186 pair<short, Energy> getLamNfTwoLoop(Energy q) const;
188
189private:
190
196
197private:
198
202 Energy _qmin;
203
209
214 double _asMaxNP;
215
219 vector<Energy> _thresholds;
220
224 vector<Energy> _lambda;
225
229 unsigned int _nloop;
230
235
239 double _alphain;
240
245
249 unsigned int _maxtry;
250
254 double _alphamin;
255
259 double _val0;
260
266
271 vector<Energy> _quarkMasses;
272};
273
274}
275
276#endif /* HERWIG_ShowerAlphaQCD_H */
This concrete class provides the definition of the pure virtual function value() and overestimateValu...
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
double _alphain
Input value of .
double _val0
Value of at the minimum scale.
virtual double overestimateValue() const
It returns the running coupling value evaluated at the input scale multiplied by the scale factor sca...
const vector< Energy > & quarkMasses() const
Return the quark masses to be used; if not empty these masses should be considered instead of the one...
Energy _optInputScale
An optional input scale to be used for the input alphas; if zero MZ will be used out of the particle ...
double _asMaxNP
Another parameter, a possible (maximum) value of alpha in the non-perturbative region.
double _tolerance
Tolerance for discontinuities at the thresholds.
virtual void doinit()
Initialize this object after the setup phase before saving an EventGenerator to disk.
vector< Energy > _quarkMasses
The quark masses to be used; if not empty these masses should be considered instead of the ones set i...
bool _thresopt
Option for the threshold masses.
vector< Energy > _lambda
for the different number of flavours
virtual double value(const Energy2 scale) const
Methods to return the coupling.
Energy _qmin
Minimum value of the scale.
virtual void initialize()
Initialize this coupling.
Energy computeLambda(Energy match, double alpha, unsigned int nflav) const
Member functions which calculate the coupling.
unsigned int _maxtry
Maximum number of iterations for the Newton-Raphson method to converge.
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
ShowerAlphaQCD()
The default constructor.
virtual double ratio(const Energy2 scale, double factor=1.) const
Return the ratio of the coupling at the scale to the overestimated value.
string value(string)
A command to initialize the coupling and write its value at the scale given by the argument (in GeV)
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
int _asType
Parameter controlling the behaviour of in the non-perturbative region.
double _alphamin
The minimum value of the coupling.
Energy lambdaQCD(unsigned int nf)
Get the value of .
ShowerAlphaQCD & operator=(const ShowerAlphaQCD &)=delete
The assignment operator is private and must never be called.
vector< Energy > _thresholds
Thresholds for the different number of flavours.
static void Init()
The standard Init function used to initialize the interfaces.
virtual IBPtr clone() const
Make a simple clone of this object.
pair< short, Energy > getLamNfTwoLoop(Energy q) const
Return the value of and the number of flavours at the scale.
string check(string args)
Match thresholds and write alpha_s specified file; arguments are Q_low/GeV Q_high/GeV n_steps filenam...
unsigned int _nloop
Option for the number of loops.
This class is the abstract class from which all types of running couplings used in the Showering deri...
Definition: ShowerAlpha.h:49
-*- C++ -*-
ThePEG::Ptr< InterfacedBase >::pointer IBPtr