7 #ifndef LHAPDF_AlphaS_H 8 #define LHAPDF_AlphaS_H 10 #include "LHAPDF/Utils.h" 11 #include "LHAPDF/Exceptions.h" 12 #include "LHAPDF/KnotArray.h" 40 virtual double alphasQ2(
double q2)
const = 0;
113 virtual std::string
type()
const = 0;
130 double _beta(
int i,
int nf)
const;
134 std::vector<double>
_betas(
int nf)
const;
179 std::string
type()
const {
return "analytic"; }
188 void setLambda(
unsigned int i,
double lambda);
194 double _lambdaQCD(
int nf)
const;
219 std::string
type()
const {
return "ipol"; }
227 void setQValues(
const std::vector<double>& qs);
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;
257 void _setup_grids()
const;
278 std::string
type()
const {
return "ode"; }
284 void setMZ(
double mz ) {
_mz = mz; _calculated =
false; }
296 void setQValues(
const std::vector<double>& qs);
301 void setQ2Values( std::vector<double> q2s ) { _q2s = q2s; _calculated =
false; }
307 double _derivative(
double t,
double y,
const std::vector<double>& beta)
const;
311 double _decouple(
double y,
double t,
unsigned int ni,
unsigned int nf)
const;
314 void _rk4(
double& t,
double& y,
double h,
const double allowed_change,
const vector<double>& bs)
const;
317 void _solve(
double q2,
double& t,
double& y,
const double& allowed_relative,
double h,
double accuracy)
const;
320 void _interpolate()
const;
324 mutable std::vector<double>
_q2s;
virtual void setLambda(unsigned int, double)
Set the {i}th Lambda constant for i active flavors.
Definition: AlphaS.h:105
void setMassReference(double mref)
Set the Z mass used in this alpha_s.
Definition: AlphaS.h:95
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:219
double quarkMass(int id) const
Get a quark mass by PDG code.
void setMZ(double mz)
Set MZ, and also the caching flag.
Definition: AlphaS.h:284
AlphaS_Ipol _ipol
The interpolation used to get Alpha_s after the ODE has been solved.
Definition: AlphaS.h:330
void setAlphaSMZ(double alphas)
Set alpha_s(MZ), and also the caching flag.
Definition: AlphaS.h:287
FlavorScheme flavorScheme() const
Get flavor scheme.
Definition: AlphaS.h:119
void setAlphaSMZ(double alphas)
Set the alpha_s(MZ) used in this alpha_s.
Definition: AlphaS.h:90
std::map< double, AlphaSArray > _knotarrays
Definition: AlphaS.h:262
FlavorScheme _flavorscheme
The flavor scheme in use.
Definition: AlphaS.h:165
bool _customref
Decides whether to use custom reference values or fall back on MZ/AlphaS_MZ.
Definition: AlphaS.h:157
virtual int numFlavorsQ2(double q2) const
Calculate the number of active flavours at energy scale Q2.
FlavorScheme
enum of flavor schemes
Definition: AlphaS.h:110
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:278
void setAlphaSReference(double alphas)
Set the alpha_s(MZ) used in this alpha_s.
Definition: AlphaS.h:100
int _nfmax
Max number of flavors.
Definition: AlphaS.h:204
virtual std::string type() const =0
Get the implementation type of this AlphaS.
int orderQCD()
Definition: AlphaS.h:75
std::map< int, double > _lambdas
LambdaQCD values.
Definition: AlphaS.h:201
void setAlphaSReference(double alphas)
Set alpha_s(MZ), and also the caching flag.
Definition: AlphaS.h:293
double _alphas_mz
Value of alpha_s(MZ)
Definition: AlphaS.h:148
std::vector< double > _q2s
Vector of Q2s in case specific anchor points are used.
Definition: AlphaS.h:324
int numFlavorsQ(double q) const
Calculate the number of active flavours at energy scale Q.
Definition: AlphaS.h:49
double alphasQ(double q) const
Calculate alphaS(Q)
Definition: AlphaS.h:36
double quarkThreshold(int id) const
Get a flavor scale threshold by PDG code.
int _fixflav
The allowed numbers of flavours in a fixed scheme.
Definition: AlphaS.h:168
int _qcdorder
Order of QCD evolution (expressed as number of loops)
Definition: AlphaS.h:142
void setAlphaSValues(const std::vector< double > &as)
Definition: AlphaS.h:241
std::vector< double > _q2s
Array of ipol knots in Q2.
Definition: AlphaS.h:265
std::vector< double > _betas(int nf) const
double _mreference
Reference mass in GeV.
Definition: AlphaS.h:151
virtual double alphasQ2(double q2) const =0
std::vector< double > _as
Array of alpha_s values for the Q2 knots.
Definition: AlphaS.h:267
double _alphas_reference
Value of alpha_s(reference mass)
Definition: AlphaS.h:154
void setQ2Values(const std::vector< double > &q2s)
Definition: AlphaS.h:234
bool _calculated
Whether or not the ODE has been solved yet.
Definition: AlphaS.h:327
virtual ~AlphaS()
Destructor.
Definition: AlphaS.h:30
std::string type() const
Implementation type of this solver.
Definition: AlphaS.h:179
Solve the differential equation in alphaS using an implementation of RK4.
Definition: AlphaS.h:274
void setMZ(double mz)
Set the Z mass used in this alpha_s.
Definition: AlphaS.h:85
Namespace for all LHAPDF functions and classes.
Definition: AlphaS.h:14
std::map< int, double > _quarkmasses
Definition: AlphaS.h:162
AlphaS()
Base class constructor for default param setup.
void setQ2Values(std::vector< double > q2s)
Set the array of Q2 values for interpolation, and also the caching flag.
Definition: AlphaS.h:301
void setFlavorScheme(FlavorScheme scheme, int nf=-1)
Set flavor scheme of alpha_s solver.
Calculator interface for computing alpha_s(Q2) in various ways.
Definition: AlphaS.h:23
double _mz
Mass of the Z boson in GeV.
Definition: AlphaS.h:145
void setMassReference(double mref)
Set reference mass, and also the caching flag.
Definition: AlphaS.h:290
void setQuarkThreshold(int id, double value)
Set a flavor threshold by PDG code (= quark masses by default)
double _beta(int i, int nf) const
Calculate alpha_s(Q2) by an analytic approximation.
Definition: AlphaS.h:175
void setOrderQCD(int order)
Set the order of QCD (expressed as number of loops)
Definition: AlphaS.h:80
int _nfmin
Min number of flavors.
Definition: AlphaS.h:206
void setQuarkMass(int id, double value)
Set quark masses by PDG code.