Solve the differential equation in alphaS using an implementation of RK4. More...
#include <AlphaS.h>
Public Member Functions | |
std::string | type () const |
Implementation type of this solver. | |
double | alphasQ2 (double q2) const |
Calculate alphaS(Q2) | |
void | setMZ (double mz) |
Set MZ, and also the caching flag. | |
void | setAlphaSMZ (double alphas) |
Set alpha_s(MZ), and also the caching flag. | |
void | setMassReference (double mref) |
Set reference mass, and also the caching flag. | |
void | setAlphaSReference (double alphas) |
Set alpha_s(MZ), and also the caching flag. | |
void | setQValues (const std::vector< double > &qs) |
Set the array of Q values for interpolation, and also the caching flag. | |
void | setQ2Values (std::vector< double > q2s) |
Set the array of Q2 values for interpolation, and also the caching flag. More... | |
![]() | |
AlphaS () | |
Base class constructor for default param setup. | |
virtual | ~AlphaS () |
Destructor. | |
double | alphasQ (double q) const |
Calculate alphaS(Q) | |
int | numFlavorsQ (double q) const |
Calculate the number of active flavours at energy scale Q. | |
virtual int | numFlavorsQ2 (double q2) const |
Calculate the number of active flavours at energy scale Q2. | |
double | quarkMass (int id) const |
Get a quark mass by PDG code. | |
void | setQuarkMass (int id, double value) |
Set quark masses by PDG code. More... | |
double | quarkThreshold (int id) const |
Get a flavor scale threshold by PDG code. More... | |
void | setQuarkThreshold (int id, double value) |
Set a flavor threshold by PDG code (= quark masses by default) More... | |
int | orderQCD () |
void | setOrderQCD (int order) |
Set the order of QCD (expressed as number of loops) More... | |
void | setMZ (double mz) |
Set the Z mass used in this alpha_s. More... | |
void | setAlphaSMZ (double alphas) |
Set the alpha_s(MZ) used in this alpha_s. More... | |
void | setMassReference (double mref) |
Set the Z mass used in this alpha_s. More... | |
void | setAlphaSReference (double alphas) |
Set the alpha_s(MZ) used in this alpha_s. More... | |
virtual void | setLambda (unsigned int, double) |
Set the {i}th Lambda constant for i active flavors. More... | |
void | setFlavorScheme (FlavorScheme scheme, int nf=-1) |
Set flavor scheme of alpha_s solver. | |
FlavorScheme | flavorScheme () const |
Get flavor scheme. | |
Private Member Functions | |
double | _derivative (double t, double y, const std::vector< double > &beta) const |
Calculate the derivative at Q2 = t, alpha_S = y. | |
double | _decouple (double y, double t, unsigned int ni, unsigned int nf) const |
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. | |
void | _solve (double q2, double &t, double &y, const double &allowed_relative, double h, double accuracy) const |
Solve alpha_s for q2 using RK4. | |
void | _interpolate () const |
Create interpolation grid. | |
Private Attributes | |
std::vector< double > | _q2s |
Vector of Q2s in case specific anchor points are used. | |
bool | _calculated |
Whether or not the ODE has been solved yet. | |
AlphaS_Ipol | _ipol |
The interpolation used to get Alpha_s after the ODE has been solved. | |
Additional Inherited Members | |
![]() | |
enum | FlavorScheme { FIXED, VARIABLE } |
Enum of flavor schemes. | |
![]() | |
double | _beta (int i, int nf) const |
std::vector< double > | _betas (int nf) const |
![]() | |
int | _qcdorder |
Order of QCD evolution (expressed as number of loops) | |
double | _mz |
Mass of the Z boson in GeV. | |
double | _alphas_mz |
Value of alpha_s(MZ) | |
double | _mreference |
Reference mass in GeV. | |
double | _alphas_reference |
Value of alpha_s(reference mass) | |
bool | _customref |
Decides whether to use custom reference values or fall back on MZ/AlphaS_MZ. | |
std::map< int, double > | _quarkmasses |
std::map< int, double > | _flavorthresholds |
FlavorScheme | _flavorscheme |
The flavor scheme in use. | |
int | _fixflav |
The allowed numbers of flavours in a fixed scheme. | |
Solve the differential equation in alphaS using an implementation of RK4.
|
private |
Calculate the decoupling relation when going from num. flav. = ni -> nf abs(ni - nf) must be = 1
|
inline |
Set the array of Q2 values for interpolation, and also the caching flag.
Writes to the same internal array as setQValues, appropriately transformed.