LHAPDF  6.2.1
AlphaS.h
1 // -*- C++ -*-
2 //
3 // This file is part of LHAPDF
4 // Copyright (C) 2012-2016 The LHAPDF collaboration (see AUTHORS for details)
5 //
6 #pragma once
7 #ifndef LHAPDF_AlphaS_H
8 #define LHAPDF_AlphaS_H
9 
10 #include "LHAPDF/Utils.h"
11 #include "LHAPDF/Exceptions.h"
12 #include "LHAPDF/KnotArray.h"
13 
14 namespace LHAPDF {
15 
23  class AlphaS {
24  public:
25 
27  AlphaS();
28 
30  virtual ~AlphaS() {};
31 
33 
34 
36  double alphasQ(double q) const { return alphasQ2(q*q); }
37 
40  virtual double alphasQ2(double q2) const = 0;
41 
43 
44 
46 
47 
49  int numFlavorsQ(double q) const { return numFlavorsQ2(q*q); }
50 
52  virtual int numFlavorsQ2(double q2) const;
53 
55  double quarkMass(int id) const;
56 
60  void setQuarkMass(int id, double value);
61 
65  double quarkThreshold(int id) const;
66 
70  void setQuarkThreshold(int id, double value);
71 
75  int orderQCD() { return _qcdorder; }
76 
80  void setOrderQCD(int order) { _qcdorder = order; }
81 
85  void setMZ(double mz) { _mz = mz; }
86 
90  void setAlphaSMZ(double alphas) { _alphas_mz = alphas; }
91 
95  void setMassReference(double mref) { _mreference = mref; _customref = true; }
96 
100  void setAlphaSReference(double alphas) { _alphas_reference = alphas; _customref = true; }
101 
105  virtual void setLambda(unsigned int, double) {};
106 
108 
110  enum FlavorScheme { FIXED, VARIABLE };
111 
113  virtual std::string type() const = 0;
114 
116  void setFlavorScheme(FlavorScheme scheme, int nf = -1);
117 
120 
121 
122  protected:
123 
125 
126 
130  double _beta(int i, int nf) const;
131 
134  std::vector<double> _betas(int nf) const;
135 
137 
138 
139  protected:
140 
143 
145  double _mz;
146 
148  double _alphas_mz;
149 
151  double _mreference;
152 
155 
158 
162  std::map<int, double> _quarkmasses, _flavorthresholds;
163 
166 
168  int _fixflav;
169 
170  };
171 
172 
173 
175  class AlphaS_Analytic : public AlphaS {
176  public:
177 
179  std::string type() const { return "analytic"; }
180 
182  double alphasQ2(double q2) const;
183 
185  int numFlavorsQ2(double q2) const;
186 
188  void setLambda(unsigned int i, double lambda);
189 
190 
191  private:
192 
194  double _lambdaQCD(int nf) const;
195 
197  void _setFlavors();
198 
199 
201  std::map<int, double> _lambdas;
202 
204  int _nfmax;
206  int _nfmin;
207 
208  };
209 
210 
211 
215  class AlphaS_Ipol : public AlphaS {
216  public:
217 
219  std::string type() const { return "ipol"; }
220 
222  double alphasQ2(double q2) const;
223 
227  void setQValues(const std::vector<double>& qs);
228 
234  void setQ2Values(const std::vector<double>& q2s) { _q2s = q2s; }
235 
241  void setAlphaSValues(const std::vector<double>& as) { _as = as; }
242 
243 
244  private:
245 
247  double _interpolateCubic(double T, double VL, double VDL, double VH, double VDH) const;
249  double _ddq_central( size_t i ) const;
251  double _ddq_forward( size_t i ) const;
253  double _ddq_backward( size_t i ) const;
254 
257  void _setup_grids() const;
258 
259 
262  mutable std::map<double, AlphaSArray> _knotarrays;
263 
265  std::vector<double> _q2s;
267  std::vector<double> _as;
268 
269  };
270 
271 
272 
274  class AlphaS_ODE : public AlphaS {
275  public:
276 
278  std::string type() const { return "ode"; }
279 
281  double alphasQ2( double q2 ) const;
282 
284  void setMZ( double mz ) { _mz = mz; _calculated = false; }
285 
287  void setAlphaSMZ( double alphas ) { _alphas_mz = alphas; _calculated = false; }
288 
290  void setMassReference( double mref ) { _mreference = mref; _calculated = false; _customref = true; }
291 
293  void setAlphaSReference( double alphas ) { _alphas_reference = alphas; _calculated = false; _customref = true; }
294 
296  void setQValues(const std::vector<double>& qs);
297 
301  void setQ2Values( std::vector<double> q2s ) { _q2s = q2s; _calculated = false; }
302 
303 
304  private:
305 
307  double _derivative(double t, double y, const std::vector<double>& beta) const;
308 
311  double _decouple(double y, double t, unsigned int ni, unsigned int nf) const;
312 
314  void _rk4(double& t, double& y, double h, const double allowed_change, const vector<double>& bs) const;
315 
317  void _solve(double q2, double& t, double& y, const double& allowed_relative, double h, double accuracy) const;
318 
320  void _interpolate() const;
321 
322 
324  mutable std::vector<double> _q2s;
325 
327  mutable bool _calculated;
328 
331 
332  };
333 
334 
335 }
336 #endif
LHAPDF::AlphaS_Ipol::_interpolateCubic
double _interpolateCubic(double T, double VL, double VDL, double VH, double VDH) const
Standard cubic interpolation formula.
LHAPDF::AlphaS_Analytic::_nfmin
int _nfmin
Min number of flavors.
Definition: AlphaS.h:206
LHAPDF::AlphaS_Ipol::setAlphaSValues
void setAlphaSValues(const std::vector< double > &as)
Definition: AlphaS.h:241
LHAPDF::AlphaS::setMZ
void setMZ(double mz)
Set the Z mass used in this alpha_s.
Definition: AlphaS.h:85
LHAPDF::AlphaS::_betas
std::vector< double > _betas(int nf) const
LHAPDF::AlphaS_Ipol::_ddq_backward
double _ddq_backward(size_t i) const
Get the gradient for a patch at the high end of the grid.
LHAPDF::AlphaS::_mz
double _mz
Mass of the Z boson in GeV.
Definition: AlphaS.h:145
LHAPDF::AlphaS_ODE::_decouple
double _decouple(double y, double t, unsigned int ni, unsigned int nf) const
LHAPDF::AlphaS_ODE::alphasQ2
double alphasQ2(double q2) const
Calculate alphaS(Q2)
LHAPDF::AlphaS::_alphas_reference
double _alphas_reference
Value of alpha_s(reference mass)
Definition: AlphaS.h:154
LHAPDF::AlphaS_Analytic::_lambdaQCD
double _lambdaQCD(int nf) const
Get lambdaQCD for nf.
LHAPDF::AlphaS_Analytic::setLambda
void setLambda(unsigned int i, double lambda)
Set lambda_i (for i = flavour number)
LHAPDF::AlphaS_ODE::setQValues
void setQValues(const std::vector< double > &qs)
Set the array of Q values for interpolation, and also the caching flag.
LHAPDF::AlphaS_ODE::_calculated
bool _calculated
Whether or not the ODE has been solved yet.
Definition: AlphaS.h:327
LHAPDF::AlphaS_ODE::_ipol
AlphaS_Ipol _ipol
The interpolation used to get Alpha_s after the ODE has been solved.
Definition: AlphaS.h:330
LHAPDF::AlphaS::_quarkmasses
std::map< int, double > _quarkmasses
Definition: AlphaS.h:162
LHAPDF::AlphaS_Ipol::_as
std::vector< double > _as
Array of alpha_s values for the Q2 knots.
Definition: AlphaS.h:267
LHAPDF::AlphaS_ODE::setMassReference
void setMassReference(double mref)
Set reference mass, and also the caching flag.
Definition: AlphaS.h:290
LHAPDF::AlphaS::_beta
double _beta(int i, int nf) const
LHAPDF::AlphaS_ODE
Solve the differential equation in alphaS using an implementation of RK4.
Definition: AlphaS.h:274
LHAPDF::AlphaS::setQuarkThreshold
void setQuarkThreshold(int id, double value)
Set a flavor threshold by PDG code (= quark masses by default)
LHAPDF::AlphaS_Analytic
Calculate alpha_s(Q2) by an analytic approximation.
Definition: AlphaS.h:175
LHAPDF::AlphaS::setOrderQCD
void setOrderQCD(int order)
Set the order of QCD (expressed as number of loops)
Definition: AlphaS.h:80
LHAPDF::AlphaS_ODE::_interpolate
void _interpolate() const
Create interpolation grid.
LHAPDF::AlphaS_Ipol::alphasQ2
double alphasQ2(double q2) const
Calculate alphaS(Q2)
LHAPDF::AlphaS::AlphaS
AlphaS()
Base class constructor for default param setup.
LHAPDF::AlphaS_Ipol::setQValues
void setQValues(const std::vector< double > &qs)
LHAPDF::AlphaS::setFlavorScheme
void setFlavorScheme(FlavorScheme scheme, int nf=-1)
Set flavor scheme of alpha_s solver.
LHAPDF::AlphaS::setMassReference
void setMassReference(double mref)
Set the Z mass used in this alpha_s.
Definition: AlphaS.h:95
LHAPDF::AlphaS_Ipol::type
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:219
LHAPDF::AlphaS_Ipol::_ddq_central
double _ddq_central(size_t i) const
Get the gradient for a patch in the middle of the grid.
LHAPDF::AlphaS::flavorScheme
FlavorScheme flavorScheme() const
Get flavor scheme.
Definition: AlphaS.h:119
LHAPDF::AlphaS_Analytic::type
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:179
LHAPDF::AlphaS_ODE::_solve
void _solve(double q2, double &t, double &y, const double &allowed_relative, double h, double accuracy) const
Solve alpha_s for q2 using RK4.
LHAPDF::AlphaS_ODE::_derivative
double _derivative(double t, double y, const std::vector< double > &beta) const
Calculate the derivative at Q2 = t, alpha_S = y.
LHAPDF::AlphaS::setQuarkMass
void setQuarkMass(int id, double value)
Set quark masses by PDG code.
LHAPDF::AlphaS_Analytic::alphasQ2
double alphasQ2(double q2) const
Calculate alphaS(Q2)
LHAPDF::AlphaS::_customref
bool _customref
Decides whether to use custom reference values or fall back on MZ/AlphaS_MZ.
Definition: AlphaS.h:157
LHAPDF::AlphaS::numFlavorsQ2
virtual int numFlavorsQ2(double q2) const
Calculate the number of active flavours at energy scale Q2.
LHAPDF::AlphaS_Ipol
Definition: AlphaS.h:215
LHAPDF::AlphaS_Analytic::_lambdas
std::map< int, double > _lambdas
LambdaQCD values.
Definition: AlphaS.h:201
LHAPDF::AlphaS_Analytic::_setFlavors
void _setFlavors()
Recalculate min/max flavors in case lambdas have changed.
LHAPDF::AlphaS::FlavorScheme
FlavorScheme
enum of flavor schemes
Definition: AlphaS.h:110
LHAPDF::AlphaS_ODE::type
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:278
LHAPDF::AlphaS::setAlphaSMZ
void setAlphaSMZ(double alphas)
Set the alpha_s(MZ) used in this alpha_s.
Definition: AlphaS.h:90
LHAPDF::AlphaS::setLambda
virtual void setLambda(unsigned int, double)
Set the {i}th Lambda constant for i active flavors.
Definition: AlphaS.h:105
LHAPDF::AlphaS_ODE::setAlphaSReference
void setAlphaSReference(double alphas)
Set alpha_s(MZ), and also the caching flag.
Definition: AlphaS.h:293
LHAPDF::AlphaS::quarkMass
double quarkMass(int id) const
Get a quark mass by PDG code.
LHAPDF::AlphaS_Ipol::_knotarrays
std::map< double, AlphaSArray > _knotarrays
Definition: AlphaS.h:262
LHAPDF::AlphaS::type
virtual std::string type() const =0
Get the implementation type of this AlphaS.
LHAPDF::AlphaS_Ipol::_setup_grids
void _setup_grids() const
LHAPDF::AlphaS_Analytic::_nfmax
int _nfmax
Max number of flavors.
Definition: AlphaS.h:204
LHAPDF::AlphaS::_alphas_mz
double _alphas_mz
Value of alpha_s(MZ)
Definition: AlphaS.h:148
LHAPDF::AlphaS::_flavorscheme
FlavorScheme _flavorscheme
The flavor scheme in use.
Definition: AlphaS.h:165
LHAPDF::AlphaS::setAlphaSReference
void setAlphaSReference(double alphas)
Set the alpha_s(MZ) used in this alpha_s.
Definition: AlphaS.h:100
LHAPDF::AlphaS::_fixflav
int _fixflav
The allowed numbers of flavours in a fixed scheme.
Definition: AlphaS.h:168
LHAPDF::AlphaS_ODE::setMZ
void setMZ(double mz)
Set MZ, and also the caching flag.
Definition: AlphaS.h:284
LHAPDF::AlphaS_Ipol::_q2s
std::vector< double > _q2s
Array of ipol knots in Q2.
Definition: AlphaS.h:265
LHAPDF::AlphaS::alphasQ
double alphasQ(double q) const
Calculate alphaS(Q)
Definition: AlphaS.h:36
LHAPDF::AlphaS_ODE::_rk4
void _rk4(double &t, double &y, double h, const double allowed_change, const vector< double > &bs) const
Calculate the next step using RK4 with adaptive step size.
LHAPDF::AlphaS::quarkThreshold
double quarkThreshold(int id) const
Get a flavor scale threshold by PDG code.
LHAPDF::AlphaS_Ipol::_ddq_forward
double _ddq_forward(size_t i) const
Get the gradient for a patch at the low end of the grid.
LHAPDF
Namespace for all LHAPDF functions and classes.
Definition: AlphaS.h:14
LHAPDF::AlphaS::alphasQ2
virtual double alphasQ2(double q2) const =0
LHAPDF::AlphaS_ODE::setAlphaSMZ
void setAlphaSMZ(double alphas)
Set alpha_s(MZ), and also the caching flag.
Definition: AlphaS.h:287
LHAPDF::AlphaS_ODE::_q2s
std::vector< double > _q2s
Vector of Q2s in case specific anchor points are used.
Definition: AlphaS.h:324
LHAPDF::AlphaS::_qcdorder
int _qcdorder
Order of QCD evolution (expressed as number of loops)
Definition: AlphaS.h:142
LHAPDF::AlphaS_Ipol::setQ2Values
void setQ2Values(const std::vector< double > &q2s)
Definition: AlphaS.h:234
LHAPDF::AlphaS::orderQCD
int orderQCD()
Definition: AlphaS.h:75
LHAPDF::AlphaS
Calculator interface for computing alpha_s(Q2) in various ways.
Definition: AlphaS.h:23
LHAPDF::AlphaS::~AlphaS
virtual ~AlphaS()
Destructor.
Definition: AlphaS.h:30
LHAPDF::AlphaS::numFlavorsQ
int numFlavorsQ(double q) const
Calculate the number of active flavours at energy scale Q.
Definition: AlphaS.h:49
LHAPDF::AlphaS::_mreference
double _mreference
Reference mass in GeV.
Definition: AlphaS.h:151
LHAPDF::AlphaS_Analytic::numFlavorsQ2
int numFlavorsQ2(double q2) const
Analytic has its own numFlavorsQ2 which respects the min/max nf set by the Lambdas.
LHAPDF::AlphaS_ODE::setQ2Values
void setQ2Values(std::vector< double > q2s)
Set the array of Q2 values for interpolation, and also the caching flag.
Definition: AlphaS.h:301