QCD Class Referenceabstract

A class for parameters related to QCD, hadrons and quarks. More...

#include <QCD.h>

Inheritance diagram for QCD:
[legend]
Collaboration diagram for QCD:
[legend]

Detailed Description

A class for parameters related to QCD, hadrons and quarks.

Author
HEPfit Collaboration

This class is a Model class that assigns and updates parameters related to and derived from QCD. A complete list of parameters in the QCD class can be found below. This class includes, but is not limited to, the running of the strong coupling constant (Full NNLO), running of the quark masses and conversions between pole mass and \(\overline{\mathrm{MS}}\) mass. All hadronization parameters like the bag parameters for the mesons and their decay constants are assigned and updated by this class.

Initialization

The constructor QCD() sets the charge and isospin of the quarks. It also sets the mass scale of the light quarks UP, DOWN and STRANGE to 2 GeV. The cache is initialized too along with the computation of \(\zeta(2)\) and \(\zeta(3)\).

The initializations and updates of the model parameters and flags are explained below.

Model parameters

The model parameters of QCD are summarized below:

Label LaTeX symbol Description
AlsM \(\alpha_s(M_{\alpha_s})\) The strong coupling constant at the scale \(M_{\alpha_s}\).
MAls \(M_{\alpha_s}\) The mass scale in GeV at which the strong coupling constant measurement is provided. Must be in the 5 flavour energy range.
mup \(m_{u}(2\,\mathrm{GeV})\) The \(\overline{\mathrm{MS}}\) mass of the up quark at 2 GeV, \(m_u(2\,\mathrm{GeV})\), in GeV.
mdown \(m_{d}(2\,\mathrm{GeV})\) The \(\overline{\mathrm{MS}}\) mass of the down quark at 2 GeV, \(m_d(2\,\mathrm{GeV})\), in GeV.
mcharm \(m_{c}(m_c)\) The \(\overline{\mathrm{MS}}\) scale-invariant mass of the charm quark, \(m_c(m_c)\), in GeV.
mstrange \(m_{s}(2\,\mathrm{GeV})\) The \(\overline{\mathrm{MS}}\) mass of the strange quark at 2 GeV , \(m_s(2\,\mathrm{GeV})\), in GeV.
mtop \(m_{t}\) The pole mass of the top quark in GeV.
mbottom \(m_{b}(m_b)\) The \(\overline{\mathrm{MS}}\) scale-invariant mass of the bottom quark, \(m_b(m_b)\), in GeV.
muc \(\mu_c\) The threshold between three- and four-flavour theory in GeV.
mub \(\mu_b\) The threshold between four- and five-flavour theory in GeV.
mut \(\mu_t\) The threshold between five- and six-flavour theory in GeV.
Attention
The parameters AlsM and MAls are not used in StandardModel and the model classes inherited from it.

The parameters below are associated with flavour observables

MK0 \(M_{K^0}\) The mass of the \( K^0 \) meson in GeV.
MKp \(M_{K^\pm}\) The mass of the \( K^\pm \) meson in GeV.
MKstar \(M_{K^*}\) The mass of the \( K^* \) meson in GeV.
Mphi \(M_{\phi}\) The mass of the \( \phi \) meson in GeV.
MD \(M_{D^0}\) The mass of the \( D^0 \) meson in GeV.
MBd \(M_{B_d}\) The mass of the \( B_d \) meson in GeV.
MBp \(M_{B^\pm}\) The mass of the \( B^\pm \) meson in GeV.
MBs \(M_{B_s}\) The mass of the \( B_s \) meson in GeV.
tKl \(\tau_{K_L}\) The lifetime of the \( K_L \) meson in \(\mathrm{ps}^{-1} \).
tKp \(\tau_{K^\pm}\) The lifetime of the \( K^\pm \) meson in \(\mathrm{ps}^{-1} \).
tKstar \(\tau_{K^*}\) The lifetime of the \( K^* \) meson in \(\mathrm{ps}^{-1} \).
tphi \(\tau_{\phi}\) The lifetime of the \( \phi \) meson in \(\mathrm{ps}^{-1} \).
tBd \(\tau_{B_d}\) The lifetime of the \( B_d \) meson in \(\mathrm{ps}^{-1} \).
tBp \(\tau_{B^\pm}\) The lifetime of the \( B^\pm \) meson in \(\mathrm{ps}^{-1} \).
tBs \(\tau_{B_s}\) The lifetime of the \( B_s \) meson in \(\mathrm{ps}^{-1} \).
DGs_Gs
FK \(F_{K^-}\) The decay constant of the \( K^- \) meson in GeV.
FKstar \(F_{K^*}\) The decay constant of the \( K^* \) meson in GeV.
FKstarp \(F_{K^*}^{\perp}\) The decay constant of a transversely polarized \( K^* \) meson in GeV.
Fphi \(F_{\phi}\) The decay constant of the \( \phi \) meson in GeV.
FD \(F_{D^0}\) The decay constant of the \( D^0 \) meson in GeV.
FBs \(F_{B_s}\) The decay constant of the \( B_s \) meson in GeV.
FBsoFBd \(F_{B_d}/F_{B_d}\) The ratio \( F_{B_s}/F_{B_d} \) necessary to compute \( F_{B_s} \).
BK1 - BK5 \(B^1_{K} - B^5_{K}\) The bag parameter for \( O_1 - O_5\) in \( \Delta s = 2 \) processes in \( K^0 \).
BKscale \(\mu_K\) The scale at which the bag parameters are specified for the \( K^0 \) system.
BKscheme \(\) The scheme in which the bag parameters are specified for the \( K^0 \) system.
BD1 - BD5 \(B^1_{D} - B^5_{D}\) The bag parameter for \( O_1 - O_5\) in \( \Delta c = 2 \) processes in \( D^0 \).
BDscale \(\mu_D\) The scale at which the bag parameters are specified for the \( D_0 \) system.
BDscheme \(\) The scheme in which the bag parameters are specified for the \( D_0 \) system.
BBsoBBd \(B_{B_s}/B_{B_d}\) The ratio \( B_{B_s}/B_{B_d} \) necessary to compute \( B_{B_s} \).
BBs1 - BBs5 \(B^1_{B_s} - B^5_{B_s}\) The bag parameter for \( O_1 - O_5 \) in \( \Delta b = 2 \) processes in \( B_s \).
BBsscale \(\mu_{B_{s}}\) The scale at which the bag parameters are specified for the \( B_s \) system.
BBsscheme \(\) The scheme in which the bag parameters are specified for the \( B_s \) system.
BK(1/2)1 - BK(1/2)10 \(\)
BK(3/2)1 - BK(3/2)10
BKd_scale \(\)
BKd_scheme \(\)
ReA0_Kd \({\cal Re}(A_0(K\to\pi\pi))\) The experimental value of the real part of the amplitude for \(K^0\to\pi\pi\) with \(\Delta I=0\).
ReA2_Kd \({\cal Re}(A_2(K\to\pi\pi))\) The experimental value of the real part of the amplitude for \(K^0\to\pi\pi\) with \(\Delta I=2\).
Omega_eta_etap \(\Omega_{\eta/\eta'}\) The isospin breaking contribution in \(K^0\to\pi\pi\).
Br_Kp_P0enu \(\mathrm{BR}(K^+\to\pi^0e^+\nu)\) The experimental value for the branching ratio of \(K^+\to\pi^0e^+\nu\).
Br_Kp_munu \(\mathrm{BR}(K^+\to\mu^+\nu)\) The experimental value for the branching ratio of \(K^+\to\mu^+\nu\).
Br_B_Xcenu \(\mathrm{BR}(B\to X_ce\nu)\) The experimental value for the branching ratio of \(B\to X_c e\nu\).
DeltaP_cu \(\) The long-distance correction to the charm contribution of \(K^+\to\pi^+\nu\bar{\nu}\).
IB_Kl \(\) The isospin breaking corrections between \(K_L\to\pi^0\nu\bar{\nu}\) and \(K^+\to\pi^0 e^+\nu\).
IB_Kp \(\) The isospin breaking corrections between \(K^+\to\pi^+ \nu\bar{\nu}\) and \(K^+\to\pi^0 e^+\nu\).
a_0V, a_1V, a_2V, MRV \(a_0^V, a_1^V, a_2^V, \Delta m^V\) The fit parameters for the form factor \(V\) of the \(B\to K^*\).
a_0A0, a_1A0, a_2A0, MRA0 \(a_0^{A_0}, a_1^{A_0}, a_2^{A_0}, \Delta m^{A_0}\) The fit parameters for the form factor \(A_0\) of the \(B\to K^*\).
a_0A1, a_1A1, a_2A1, MRA1 \(a_0^{A_1}, a_1^{A_1}, a_2^{A_1}, \Delta m^{A_1}\) The fit parameters for the form factor \(A_1\) of the \(B\to K^*\).
a_0A12, a_1A12, a_2A12, MRA12 \(a_0^{A_{12}}, a_1^{A_{12}}, a_2^{A_{12}}, \Delta m^{A_{12}}\) The fit parameters for the form factor \(A_{12}\) of the \(B\to K^*\).
a_0T1, a_1T1, a_2T1, MRA0 \(a_0^{T_1}, a_1^{T_1}, a_2^{T_1}, \Delta m^{T_1}\) The fit parameters for the form factor \(T_1\) of the \(B\to K^*\).
a_0T2, a_1T2, a_2T2, MRA1 \(a_0^{T_2}, a_1^{T_2}, a_2^{T_2}, \Delta m^{T_2}\) The fit parameters for the form factor \(T_2\) of the \(B\to K^*\).
a_0T23, a_1T23, a_2T23, MRA1 \(a_0^{T_{23}}, a_1^{T_{23}}, a_2^{T_{23}}, \Delta m^{T_{23}}\) The fit parameters for the form factor \(T_{23}\) of the \(B\to K^*\).
a_0Vphi, a_1Vphi, a_2Vphi, MRVphi \(a_0^V, a_1^V, a_2^V, \Delta m^V\) The fit parameters for the form factor \(V\) of the \(B\to\phi\).
a_0A0phi, a_1A0phi, a_2A0phi, MRA0phi \(a_0^{A_0}, a_1^{A_0}, a_2^{A_0}, \Delta m^{A_0}\) The fit parameters for the form factor \(A_0\) of the \(B\to\phi\).
a_0A1phi, a_1A1phi, a_2A1phi, MRA1phi \(a_0^{A_1}, a_1^{A_1}, a_2^{A_1}, \Delta m^{A_1}\) The fit parameters for the form factor \(A_1\) of the \(B\to\phi\).
a_0A1phi, a_1A1phi, a_2A1phi, MRA1phi \(a_0^{A_{12}}, a_1^{A_{12}}, a_2^{A_{12}}, \Delta m^{A_{12}}\) The fit parameters for the form factor \(A_{12}\) of the \(B\to\phi\).
a_0T1phi, a_1T1phi, a_2T1phi, MRA0phi \(a_0^{T_1}, a_1^{T_1}, a_2^{T_1}, \Delta m^{T_1}\) The fit parameters for the form factor \(T_1\) of the \(B\to\phi\).
a_0T2phi, a_1T2phi, a_2T2phi, MRA1phi \(a_0^{T_2}, a_1^{T_2}, a_2^{T_2}, \Delta m^{T_2}\) The fit parameters for the form factor \(T_2\) of the \(B\to\phi\).
a_0T23phi, a_1T23phi, a_2T23phi, MRA1phi \(a_0^{T_{23}}, a_1^{T_{23}}, a_2^{T_{23}}, \Delta m^{T_{23}}\) The fit parameters for the form factor \(T_{23}\) of the \(B\to\phi\).
absh_0, absh_0_1, absh_0_2 \(\vert h_0^{(0)}\vert, \vert h_0^{(1)}\vert, \vert h_0^{(2)}\vert\) The constant, linear and quadratic terms of the absolute value of the hadronic parameter \(h_0\) of the \(B\to K^*\).
argh_0, argh_0_1, argh_0_2 \(\mathrm{Arg}\,h_0^{(0)}, \mathrm{Arg}\,h_0^{(1)}, \mathrm{Arg}\,h_0^{(2)}\) The constant, linear and quadratic terms of the argument of the hadronic parameter \(h_0\) of the \(B\to K^*\).
absh_1, absh_1_1, absh_1_2 \(\vert h_+^{(0)}\vert, \vert h_+^{(1)}\vert, \vert h_+^{(2)}\vert\) The constant, linear and quadratic terms of the absolute value of the hadronic parameter \(h_+\) of the \(B\to K^*\).
argh_1, argh_1_1, argh_1_2 \(\mathrm{Arg}\,h_+^{(0)}, \mathrm{Arg}\,h_+^{(1)}, \mathrm{Arg}\,h_+^{(2)}\) The constant, linear and quadratic terms of the argument of the hadronic parameter \(h_+\) of the \(B\to K^*\).
absh_2, absh_2_1, absh_2_2 \(\vert h_-^{(0)}\vert, \vert h_-^{(1)}\vert, \vert h_-^{(2)}\vert\) The constant, linear and quadratic terms of the absoulute value of the hadronic parameter \(h_-\) of the \(B\to K^*\).
argh_2, argh_2_1, argh_2_2 \(\mathrm{Arg}\,h_-^{(0)}, \mathrm{Arg}\,h_-^{(1)}, \mathrm{Arg}\,h_-^{(2)}\) The constant, linear and quadratic terms of the argument of the hadronic parameter \(h_-\) of the \(B\to K^*\).
alpha1kst, alpha2kst \(\alpha_1(\bar{K}^*), \alpha_2(\bar{K}^*)\) The Gegenbauer coefficients for the \(\bar{K}^*\) meson.
lambdaB \(\Lambda_{B,+}\) The integrated leading twist light-cone distribution amplitudes of the B meson divided by the integral variable.
r_1_fplus, r_2_fplus, m_fit2_fplus \(r_1^{f_+}, r_2^{f_+}, m_{fit}^{2,f_+}\) The fit parameters for the LCSR form factor \(f_+\) of the \(B\to K\).
r_1_fT, r_2_fT, m_fit2_fT \(r_1^{f_T}, r_2^{f_T}, m_{fit}^{2,f_T}\) The fit parameters for the LCSR form factor \(f_T\) of the \(B\to K\).
r_2_f0, m_fit2_f0 \(r_2^{f_0}, m_{fit}^{2,f_0}\) The fit parameters for the LCSR form factor \(f_0\) of the \(B\to K\).
absh_0_MP, absh_0_1_MP \(\vert h_0^{(0)}\vert, \vert h_0^{(1)}\vert\) The constant and linear terms of the absolute value of the hadronic parameter \(h_0\) of the \(B\to K\).
argh_0_MP, argh_0_1_MP \(\mathrm{Arg}\,h_0^{(0)}, \mathrm{Arg}\,h_0^{(1)}\) The constant and linear terms of the argument of the hadronic parameter \(h_0\) of the \(B\to K\).
bsgamma_E0 \(E_0^{b\to s\gamma}\) The experimental energy cutoff \(E_0\) of the \(b\to s\gamma\).
BLNPcorr \(N_{b\to s\gamma}\) The non perturbative uncertainty associated to the \(b\to s\gamma\) BR.
Gambino_mukin \(\mu^{\rm kin}\) The mass scale for the b quark in the kinetic scheme, employed in the \(b\to s\gamma\) BR.
Gambino_BRsem \(\mathrm{BR}(B\to X_ce\nu)^{\rm Gambino}\) The branching ratio of \(B\to X_c e\nu\) fitted by Paolo Gambino, employed in the \(b\to s\gamma\) BR.
Gambino_Mbkin \(m_b^{\rm kin,\,Gambino}\) The b quark mass in the kinetic scheme fitted by Paolo Gambino, employed in the \(b\to s\gamma\) BR.
Gambino_Mcatmuc \(m_c(\mu_c)^{\rm Gambino}\) The c quark at \(\mu_c\) fitted by Paolo Gambino, employed in the \(b\to s\gamma\) BR.
Gambino_mupi2, Gambino_rhoD3, Gambino_muG2, Gambino_rhoLS3 \(\mu_{\pi}^{2,\,\rm Gambino}, \rho_{D}^{3,\,\rm Gambino}, \mu_{G}^{2,\,\rm Gambino}, \rho_{LS}^{3,\,\rm Gambino}\) The B meson expecation values for the relevant dim. 5 and 6 operators fitted by Paolo Gambino, employed in the \(b\to s\gamma\) BR.

The set of the model parameters are initialized and updated with the methods Init() and Update(), respectively, where the former calls the latter. In Update(), the methods PreUpdate() and PostUpdate() are called to run all the procedures that are need to be executed before and after the model parameters are updated. The \(\overline{\mathrm{MS}}\) mass for the top quark is computed and the scale set in PostUpdate() with the updated parameters. Inside the Update() method, the individual model parameter is assigned with the protected member function setParameter().

Computation of the strong coupling constant \(\alpha_s\):

The strong coupling constant \(\alpha_s\) at an arbitrary scale can be computed with the member functions:

  • AlsWithInit(const double mu, const double alsi, const double mu_i, const orders order),
  • AlsWithLambda(const double mu, const orders order),

where another function

  • Als(const double mu, const orders order = FULLNLO)

calls AlsWithInit() for order=LO/FULLNLO(NLO), and AlsWithLambda() for order=FULLNNLO(NNLO).

The function AlsWithInit() computes \(\alpha_s(\mu)\) with a given initial value \(\alpha_s(\mu_i)\):

\[ \alpha_s(\mu)=\frac{\alpha_s(\mu_i)}{v(\mu)},\qquad \alpha_s(\mu)=\frac{\alpha_s(\mu_i)}{v(\mu)} \left[1-\frac{\beta_1}{\beta_0}\frac{\alpha_s(\mu_i)}{4\pi} \frac{\ln v(\mu)}{v(\mu)}\right], \]

at LO and FULLNLO (NLO) respectively, where

\[ v(\mu) = 1- \beta_0\frac{\alpha_s(\mu_i)}{2\pi}\ln\frac{\mu_i}{\mu}, \]

and the one-loop and two-loop beta functions are given by

\[ \beta_0 = \frac{11N_c-2N_f}{3},\qquad \beta_1 = \frac{34}{3}N_c^2-\frac{10}{3}N_cN_f-2C_FN_f. \]

The function AlsWithLambda() computes \(\alpha_s(\mu)\) with the use of \(\Lambda_{\rm QCD}\), where the value of \(\Lambda_{\rm QCD}\) for \(N_f=5\) is derived from \(\alpha_s(M_{\alpha_s})\) by solving the equation (see e.g., [56], which follows the convention in [45] for \(\ln(\mu^2/\Lambda^2)\)):

\[ \frac{\alpha_s(\mu)}{4\pi}=\frac{1}{\beta_0L}-\frac{\beta_1\ln L}{\beta_0^3L^2}+ \frac{1}{\left(\beta_0L\right)^3}\left[\frac{\beta_1^2}{\beta_0^2}\left(\ln^2L-\ln L-1\right)+ \frac{\beta_2}{\beta_0}\right]+O\left(\frac{\ln^3L}{L^4}\right), \]

where \(L\equiv\ln(\mu^2/\Lambda^2)\), and the three-loop beta function is given by

\[ \beta_2 = \frac{2857}{54}N_c^3+C_F^2N_f-\frac{205}{18}C_FN_cN_f -\frac{1415}{54}N_c^2N_f+\frac{11}{9}C_FN_f^2 + \frac{79}{54}N_cN_f^2. \]

For \(N_f < 5\), \(\Lambda_{QCD}\) can be obtained by solving the following matching condition at \(\mu\) [56] [52] :

\begin{eqnarray} \beta_0^{'}\ln\frac{\Lambda^{'2}}{\Lambda^2} &=& (\beta_0^{'}-\beta_0)L+\left(\frac{\beta_1^{'}}{\beta_0^{'}}- \frac{\beta_1}{\beta_0}\right)\ln L-\frac{\beta_1^{'}}{\beta_0^{'}}\ln\frac{\beta_0^{'}}{\beta_0} - C_1 \\ &&+ \frac{1}{\beta_0L}\left[\frac{\beta_1}{\beta_0}\left(\frac{\beta_1^{'}}{\beta_0^{'}}-\frac{\beta_1}{\beta_0} \right)\ln L + \frac{\beta_1^{'2}}{\beta_0^{'2}}-\frac{\beta_1^2}{\beta_0^2}-\frac{\beta_2^{'2}}{\beta_0^{'2}} +\frac{\beta_2^2}{\beta_0^2}+\frac{\beta_1^{'}}{\beta_0^{'}}C_1-C_1^2-C_2\right]+ O\left(\frac{\ln^2L}{L^2}\right), \end{eqnarray}

where the primed (unprimed) quantities refer to those pertaining to \(N_f-1\) \((N_f)\). The terms \(C_1\) and \(C_2\) are given by

\[ C_1 = \frac{2}{3}\ln\frac{\mu^2}{\mu_f^2},\qquad C_2 = -16\left(\frac{1}{36}\ln^2\frac{\mu^2}{\mu_f^2}-\frac{19}{24}\ln\frac{\mu^2}{\mu_f^2}+\frac{11}{72}\right). \]

where \(\mu_f=m_f(\mu_f)\) denotes the \(\overline{\mathrm{MS}}\) invariant mass of the \(N_f\)-th flavour. Moreover, the matching condition at a flavour threshold is given by

\[ \alpha_s^{(N_f-1)}(\mu) = (\zeta_g^f)^2\alpha_s^{(N_f)}(\mu), \]

where

\[ (\zeta_g^f)^2 = 1+\frac{\alpha_s^{(N_f)}(\mu)}{\pi}\left(-\frac{1}{6}\ln\frac{\mu^2}{\mu_f^2}\right)+ \left(\frac{\alpha_s^{(N_f)}(\mu)}{\pi}\right)^2\left(\frac{1}{36}\ln^2\frac{\mu^2}{\mu_f^2} - \frac{19}{24}\ln\frac{\mu^2}{\mu_f^2} + \frac{11}{72}\right) + O\left(\left(\frac{\alpha_s}{\pi} \right)^3\right). \]

For the top quark mass, the pole mass \(m_t\) is used as an input instead of the \(\overline{\mathrm{MS}}\) invariant mass. Then the \(\overline{\mathrm{MS}}\) invariant mass is computed from the pole mass with the computed value of \(\alpha_s^{(6)}(m_t)\), which in turn is computed from \(\alpha_s^{(5)}(M_Z)\). Hence, the matching condition from \(N_f = 6\) to \(N_f = 5\) has to be free from the \(\overline{\mathrm{MS}}\) invariant mass of the top quark. In this case we use the following matching condition

\[ \alpha_s^{(6)}(\mu_t) = \frac{1}{(\zeta_g^{OS,t})^2}\alpha_s^{(5)}(\mu_t), \]

where

\[ \frac{1}{(\zeta_g^{OS,t})^2} = 1+\frac{\alpha_s^{(5)}(\mu)}{\pi}\left(\frac{1}{6}\ln\frac{\mu^2}{m_t^2}\right)+ \left(\frac{\alpha_s^{(5)}(\mu)}{\pi}\right)^2\left(\frac{1}{36}\ln^2\frac{\mu^2}{m_t^2} + \frac{19}{24}\ln\frac{\mu^2}{m_t^2} + \frac{7}{24}\right) + O\left(\left(\frac{\alpha_s}{\pi}\right)^3\right). \]

Besides, \(\Lambda_{QCD}\) for \(N_f=6\) is derived from that for \(N_f=5\) with the relation:

\begin{eqnarray} \beta_0\ln\frac{\Lambda^{2}}{\Lambda^{'2}} &=& (\beta_0-\beta_0^{'})L^{'}+\left(\frac{\beta_1}{\beta_0}- \frac{\beta_1^{'}}{\beta_0^{'}}\right)\ln L^{'}-\frac{\beta_1}{\beta_0}\ln\frac{\beta_0}{\beta_0^{'}}-C_1^{'} \\ && + \frac{1}{\beta_0^{'}L^{'}}\left[\frac{\beta_1^{'}}{\beta_0^{'}}\left(\frac{\beta_1}{\beta_0}-\frac{\beta_1^{'}} {\beta_0^{'}}\right)\ln L^{'}+\frac{\beta_1^{2}}{\beta_0^{2}}-\frac{\beta_1^{'2}}{\beta_0^{'2}}- \frac{\beta_2^{2}}{\beta_0^{2}}+\frac{\beta_2^{'2}}{\beta_0^{'2}}+\frac{\beta_1}{\beta_0}C_1^{'}-C_1^{'2}- C_2^{'}\right]+O\left(\frac{\ln^2L}{L^2}\right), \end{eqnarray}

where

\[ C_1^{'} = -\frac{2}{3}\ln\frac{\mu^2}{m_t^2},\qquad C_2^{'} = -16\left(\frac{1}{36}\ln^2\frac{\mu^2}{m_t^2}+\frac{19}{24}\ln\frac{\mu^2}{m_t^2}+\frac{7}{24}\right). \]

Computation of the running of the quark masses:

In the \(\overline{\mathrm{MS}}\) scheme the quark mass at a scale \(\mu\) is given by (see e.g., [57])

\[ m_q(\mu) = m_q(\mu_0)\left[\frac{\alpha_s(\mu)}{\alpha_s(\mu_0)}\right]^{\frac{\gamma^{(0)}_m}{2\beta_0}} \left\{1+A_1\frac{\alpha_s(\mu)-\alpha_s(\mu_0)}{4\pi}+\frac{A_1^2}{2}\left(\frac{\alpha_s(\mu)-\alpha_s(\mu_0)}{4\pi} \right)^2+\frac{A_2}{2}\left[\left(\frac{\alpha_s(\mu)}{4\pi}\right)^2-\left(\frac{\alpha_s(\mu_0)}{4\pi} \right)^2\right]\right\}, \]

where

\[ A_1 = \frac{\gamma^{(1)}_m}{2\beta_0} - \frac{\beta_1\gamma^{(0)}_m}{2\beta_0^2},\qquad A_2 = \frac{\beta_1^2\gamma^{(0)}_m}{2\beta_0^3}-\frac{\beta_2\gamma^{(0)}_m}{2\beta_0^2}- \frac{\beta_1\gamma^{(1)}_m}{2\beta_0^2}+\frac{\gamma^{(2)}_m}{2\beta_0}, \]

and

\[ \gamma_m^{(0)} = 6C_F,\qquad \gamma_m^{(1)} = C_F\left(3C_F+\frac{97}{3}N_c-\frac{10}{3}N_f\right),\\ \gamma_m^{(2)} = 129C_F^3-\frac{129}{2}C_F^2N_c+\frac{11413}{54}C_FN_c^2+C_F^2N_f(-46+48\zeta(3))+ C_FN_CN_f\left(-\frac{556}{27}-48\zeta(3)\right)-\frac{70}{27}C_FN_f^2, \]

in the \(\overline{\mathrm{MS}}\) scheme. The threshold conditions are given by [52]

\[ m_q^{(N_f-1)}(\mu)=\zeta_m^f\, m_q^{(N_f)}(\mu),\qquad m_q^{(N_f)}(\mu)=\frac{1}{\zeta_m^f}\, m_q^{(N_f-1)}(\mu), \]

where

\[ \zeta^f_m=1+\left(\frac{\alpha_s^{(N_f)}(\mu)}{\pi}\right)^2\left(\frac{1}{12}\ln^2\frac{\mu^2}{\mu_f^2} -\frac{5}{36}\ln\frac{\mu^2}{\mu_f^2}+\frac{89}{432}\right),\\ \frac{1}{\zeta^f_m}=1+\left(\frac{\alpha_s^{(N_f-1)}(\mu)}{\pi}\right)^2\left(-\frac{1}{12}\ln^2\frac{\mu^2} {\mu_f^2}+\frac{5}{36}\ln\frac{\mu^2}{\mu_f^2}-\frac{89}{432}\right), \]

with \(\mu_f=m_f(\mu_f)\).

Pole mass vs. \( \overline{\mathrm{MS}} \) scale invariant mass:

The pole mass \(M_q\) of a heavy quark is related to the \(\overline{\mathrm{MS}}\) scale invariant mass \(m_q(m_q)\) as [52]

\[ m_q(m_q)=M_q\left\{1-\frac{4}{3}\frac{\alpha_s(M_q)}{\pi}+\left[-\frac{2251}{268}-2\zeta(2) -\frac{2}{3}\zeta(2)\ln2+\frac{\zeta(3)}{6}+\frac{n_l}{3}\left(\zeta(2)+\frac{71}{48}\right) -\frac{4}{3}\sum_{1\le i \le n_l}\Delta\left(\frac{M_i}{M_q}\right)\right]\left(\frac{\alpha_s(M_q)}{\pi} \right)^2\right\}, \]

where \(n_l\) is the number of the light quarks, and \(\Delta(x)\) is given by

\[ \Delta(x)=\frac{\pi^2}{8}x-0.597x^2+0.230x^3+O\left(x^4\right). \]

It should be noted that the above formula requires the light quark pole masses which are not well defined. The pole mass in terms of the \(\overline{\mathrm{MS}}\) mass is given by

\[ M_q=m_q(m_q)\left\{1+\frac{4}{3}\frac{\alpha_s(m_q)}{\pi}+\left[\frac{307}{32}+2\zeta(2) +\frac{2}{3}\zeta(2)\ln2-\frac{\zeta(3)}{6}-\frac{n_l}{3}\left(\zeta(2)+\frac{71}{48}\right) +\frac{4}{3}\sum_{1\le i \le n_l}\Delta\left(\frac{m_i}{m_q}\right)\right]\left(\frac{\alpha_s(m_q)}{\pi} \right)^2\right\} \]

To get the pole mass of the light quarks we solve the above equation numerically and derive it from the corresponding \(\overline{\mathrm{MS}}\) mass.

Definition at line 707 of file QCD.h.

Public Types

enum  meson { P_0, P_P, K_0, K_P, D_0, B_D, B_P, B_S, PHI, K_star, MESON_END }
 An enum type for mesons. More...
 
enum  quark { UP, DOWN, CHARM, STRANGE, TOP, BOTTOM }
 An enum type for quarks. More...
 

Public Member Functions

double AboveTh (const double mu) const
 The active flavour threshold above the scale \(\mu\) as defined in QCD::Thresholds(). More...
 
virtual double alphaMz () const =0
 
double Als (const double mu, const orders order=FULLNLO) const
 Computes the running strong coupling \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme. In the cases of LO, NLO and FULLNNLO, the coupling is computed with AlsWithInit(). On the other hand, in the cases of NNLO and FULLNNLO, the coupling is computed with AlsWithLambda(). More...
 
double Als4 (const double mu) const
 The value of \(\alpha_s^{\mathrm{FULLNLO}}\) at any scale \(\mu\) with the number of flavours \(n_f = 4\). More...
 
double AlsWithInit (const double mu, const double alsi, const double mu_i, const orders order) const
 Computes the running strong coupling \(\alpha_s(\mu)\) from \(\alpha_s(\mu_i)\) in the \(\overline{\mathrm{MS}}\) scheme, where it is forbidden to across a flavour threshold in the RG running from \(\mu_i\) to \(\mu\). More...
 
double AlsWithLambda (const double mu, const orders order) const
 Computes the running strong coupling \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme with the use of \(\Lambda_{\rm QCD}\). More...
 
double BelowTh (const double mu) const
 The active flavour threshold below the scale \(\mu\) as defined in QCD::Thresholds(). More...
 
double Beta0 (const double nf) const
 The \(\beta_0(n_f)\) coefficient for a certain number of flavours \(n_f\). More...
 
double Beta1 (const double nf) const
 The \(\beta_1(n_f)\) coefficient for a certain number of flavours \(n_f\). More...
 
double Beta2 (const double nf) const
 The \(\beta_2(n_f)\) coefficient for a certain number of flavours \(n_f\). More...
 
virtual bool CheckFlags () const
 A method to check the sanity of the set of model flags. More...
 
virtual bool CheckParameters (const std::map< std::string, double > &DPars)
 A method to check if all the mandatory parameters for StandardModel have been provided in model initialization. More...
 
double geta_0A0 () const
 
double geta_0A0phi () const
 
double geta_0A1 () const
 
double geta_0A12 () const
 
double geta_0A12phi () const
 
double geta_0A1phi () const
 
double geta_0T1 () const
 
double geta_0T1phi () const
 
double geta_0T2 () const
 
double geta_0T23 () const
 
double geta_0T23phi () const
 
double geta_0T2phi () const
 
double geta_0V () const
 
double geta_0Vphi () const
 
double geta_1A0 () const
 
double geta_1A0phi () const
 
double geta_1A1 () const
 
double geta_1A12 () const
 
double geta_1A12phi () const
 
double geta_1A1phi () const
 
double geta_1T1 () const
 
double geta_1T1phi () const
 
double geta_1T2 () const
 
double geta_1T23 () const
 
double geta_1T23phi () const
 
double geta_1T2phi () const
 
double geta_1V () const
 
double geta_1Vphi () const
 
double geta_2A0 () const
 
double geta_2A0phi () const
 
double geta_2A1 () const
 
double geta_2A12 () const
 
double geta_2A12phi () const
 
double geta_2A1phi () const
 
double geta_2T1 () const
 
double geta_2T1phi () const
 
double geta_2T2 () const
 
double geta_2T23 () const
 
double geta_2T23phi () const
 
double geta_2T2phi () const
 
double geta_2V () const
 
double geta_2Vphi () const
 
double getAlsM () const
 A get method to access the value of \(\alpha_s(M_{\alpha_s})\). More...
 
BParameter getBBd () const
 For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta b = 2\) process in the \(B_d\) meson system. More...
 
BParameter getBBs () const
 For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta b = 2\) process in the \(B_s\) meson system. More...
 
BParameter getBD () const
 For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta c = 2\) process in the \(D^0\) meson system. More...
 
BParameter getBK () const
 For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta s = 2\) process in the \(K^0\) meson system. More...
 
BParameter getBKd1 () const
 
BParameter getBKd3 () const
 
double getBLNPcorr () const
 
double getBr_B_Xcenu () const
 
double getBr_Kp_munu () const
 
double getBr_Kp_P0enu () const
 
double getbsgamma_E0 () const
 
double getCF () const
 A get method to access the Casimir factor of QCD. More...
 
double getDeltaP_cu () const
 
double getFKstarp () const
 
double getGambino_BRsem () const
 
double getGambino_Mbkin () const
 
double getGambino_Mcatmuc () const
 
double getGambino_muG2 () const
 
double getGambino_mukin () const
 
double getGambino_mupi2 () const
 
double getGambino_rhoD3 () const
 
double getGambino_rhoLS3 () const
 
gslpp::complex geth_0 () const
 
gslpp::complex geth_0_1 () const
 
gslpp::complex geth_0_1_MP () const
 
gslpp::complex geth_0_2 () const
 
gslpp::complex geth_0_MP () const
 
gslpp::complex geth_m () const
 
gslpp::complex geth_m_1 () const
 
gslpp::complex geth_m_2 () const
 
gslpp::complex geth_p () const
 
gslpp::complex geth_p_1 () const
 
gslpp::complex geth_p_2 () const
 
double getIB_Kl () const
 
double getIB_Kp () const
 
double getm_fit2_f0 () const
 
double getm_fit2_fplus () const
 
double getm_fit2_fT () const
 
double getMAls () const
 A get method to access the mass scale \(M_{\alpha_s}\) at which the strong coupling constant measurement is provided. More...
 
Meson getMesons (const meson m) const
 A get method to access a meson as an object of the type Meson. More...
 
double getMRA0 () const
 
double getMRA0phi () const
 
double getMRA1 () const
 
double getMRA12 () const
 
double getMRA12phi () const
 
double getMRA1phi () const
 
double getMRT1 () const
 
double getMRT1phi () const
 
double getMRT2 () const
 
double getMRT23 () const
 
double getMRT23phi () const
 
double getMRT2phi () const
 
double getMRV () const
 
double getMRVphi () const
 
double getMtpole () const
 A get method to access the pole mass of the top quark. More...
 
double getMub () const
 A get method to access the threshold between five- and four-flavour theory in GeV. More...
 
double getMuc () const
 A get method to access the threshold between four- and three-flavour theory in GeV. More...
 
double getMut () const
 A get method to access the threshold between six- and five-flavour theory in GeV. More...
 
double getNc () const
 A get method to access the number of colours \(N_c\). More...
 
double getOmega_eta_etap () const
 
Particle getQuarks (const quark q) const
 A get method to access a quark as an object of the type Particle. More...
 
double getr_1_fplus () const
 
double getr_1_fT () const
 
double getr_2_f0 () const
 
double getr_2_fplus () const
 
double getr_2_fT () const
 
double getReA0_Kd () const
 
double getReA2_Kd () const
 
virtual bool Init (const std::map< std::string, double > &DPars)
 Initializes the QCD parameters found in the argument. More...
 
double logLambda (const double nf, orders order) const
 Computes \(\ln\Lambda_\mathrm{QCD}\) with nf flavours in GeV. More...
 
double Nf (const double mu) const
 The number of active flavour at scale \(\mu\). More...
 
std::string orderToString (const orders order) const
 Converts an object of the enum type "orders" to the corresponding string. More...
 
virtual bool PostUpdate ()
 The post-update method for QCD. More...
 
virtual bool PreUpdate ()
 The pre-update method for QCD. More...
 
 QCD ()
 Constructor. More...
 
virtual bool setFlag (const std::string name, const bool value)
 A method to set a flag of QCD. More...
 
virtual bool setFlagStr (const std::string name, const std::string value)
 A method to set a flag of QCD. More...
 
void setNc (double Nc)
 A set method to change the number of colours \(N_c\). More...
 
double Thresholds (const int i) const
 For accessing the active flavour threshold scales. More...
 
virtual bool Update (const std::map< std::string, double > &DPars)
 The update method for QCD. More...
 
- Public Member Functions inherited from Model
const double & getModelParam (std::string name) const
 
bool IsModelInitialized () const
 A method to check if the model is initialized. More...
 
bool isModelParam (std::string name) const
 
bool isModelSUSY () const
 
bool isModelTHDM () const
 
bool IsUpdateError () const
 A method to check if there was any error in the model update process. More...
 
 Model ()
 The default constructor. More...
 
std::string ModelName () const
 A method to fetch the name of the model. More...
 
void setModelInitialized (bool ModelInitialized)
 A set method to fix the failure or success of the initialization of the model. More...
 
void setModelName (const std::string name)
 A method to set the name of the model. More...
 
void setModelSUSY ()
 
void setModelTHDM ()
 
void setUpdateError (bool UpdateError)
 A set method to fix the update status as success or failure. More...
 
virtual ~Model ()
 The default destructor. More...
 

Static Public Attributes

static const int NQCDvars = 186
 The number of model parameters in QCD. More...
 
static const std::string QCDvars [NQCDvars]
 An array containing the labels under which all QCD parameters are stored in a vector of ModelParameter via InputParser::ReadParameters(). More...
 

Additional Inherited Members

- Protected Member Functions inherited from Model
virtual void setParameter (const std::string name, const double &value)=0
 A method to set the value of a parameter of the model. More...
 
- Protected Attributes inherited from Model
std::map< std::string, boost::reference_wrapper< const double > > ModelParamMap
 
bool UpdateError
 A boolean set to false if update is successful. More...
 

Member Enumeration Documentation

enum QCD::meson

An enum type for mesons.

Enumerator
P_0 

\(\pi^0\) meson

P_P 

\(\pi^\pm\) meson

K_0 

\(K^0\) meson

K_P 

\(K^\pm\) meson

D_0 

\(D^0\) meson

B_D 

\(B_d\) meson

B_P 

\(B^\pm\) meson

B_S 

\(B_s\) meson

PHI 

\(\phi\) meson

K_star 

\(K^*\) meson

MESON_END 

The size of this enum.

Definition at line 713 of file QCD.h.

713  {
714  P_0,
715  P_P,
716  K_0,
717  K_P,
718  D_0,
719  B_D,
720  B_P,
721  B_S,
722  PHI,
723  K_star,
724  MESON_END
725  };
Definition: QCD.h:717
Definition: QCD.h:716
Definition: QCD.h:718
Definition: QCD.h:715
Definition: QCD.h:720
Definition: QCD.h:721
Definition: QCD.h:722
Definition: QCD.h:719
Definition: QCD.h:714
enum QCD::quark

An enum type for quarks.

Enumerator
UP 

Up quark

DOWN 

Down quark

CHARM 

Charm quark

STRANGE 

Strange quark

TOP 

Top quark

BOTTOM 

Bottom quark

Definition at line 730 of file QCD.h.

730  {
731  UP,
732  DOWN,
733  CHARM,
734  STRANGE,
735  TOP,
736  BOTTOM
737  };
Definition: QCD.h:731
Definition: QCD.h:735
Definition: QCD.h:732

Constructor & Destructor Documentation

QCD::QCD ( )

Constructor.

Definition at line 56 of file QCD.cpp.

57 : BBs(5), BBd(5), BD(5), BK(5), BKd1(10), BKd3(10)
58 {
59  Nc = 3.;
60  CF = Nc / 2. - 1. / (2. * Nc);
61  // Particle(std::string name, double mass, double mass_scale = 0., double width = 0., double charge = 0.,double isospin = 0.);
62  quarks[UP] = Particle("UP", 0., 2., 0., 2. / 3., .5);
63  quarks[CHARM] = Particle("CHARM", 0., 0., 0., 2. / 3., .5);
64  quarks[TOP] = Particle("TOP", 0., 0., 0., 2. / 3., .5);
65  quarks[DOWN] = Particle("DOWN", 0., 2., 0., -1. / 3., -.5);
66  quarks[STRANGE] = Particle("STRANGE", 0., 2., 0., -1. / 3., -.5);
67  quarks[BOTTOM] = Particle("BOTTOM", 0., 0., 0., -1. / 3., -.5);
68  zeta2 = gsl_sf_zeta_int(2);
69  zeta3 = gsl_sf_zeta_int(3);
70  for (int i = 0; i < CacheSize; i++) {
71  for (int j = 0; j < 8; j++)
72  als_cache[j][i] = 0.;
73  for (int j = 0; j < 4; j++)
74  logLambda5_cache[j][i] = 0.;
75  for (int j = 0; j < 10; j++)
76  mrun_cache[j][i] = 0.;
77  for (int j = 0; j < 5; j++)
78  mp2mbar_cache[j][i] = 0.;
79  }
80 
81  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("AlsM", boost::cref(AlsM)));
82  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MAls", boost::cref(MAls)));
83  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mup", boost::cref(quarks[UP].getMass())));
84  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mdown", boost::cref(quarks[DOWN].getMass())));
85  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mcharm", boost::cref(quarks[CHARM].getMass())));
86  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mstrange", boost::cref(quarks[STRANGE].getMass())));
87  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mtop", boost::cref(mtpole)));
88  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mbottom", boost::cref(quarks[BOTTOM].getMass())));
89  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("muc", boost::cref(muc)));
90  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mub", boost::cref(mub)));
91  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("mut", boost::cref(mut)));
92  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MK0", boost::cref(mesons[K_0].getMass())));
93  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MKp", boost::cref(mesons[K_P].getMass())));
94  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MD", boost::cref(mesons[D_0].getMass())));
95  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MBd", boost::cref(mesons[B_D].getMass())));
96  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MBp", boost::cref(mesons[B_P].getMass())));
97  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MBs", boost::cref(mesons[B_S].getMass())));
98  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MKstar", boost::cref(mesons[K_star].getMass())));
99  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Mphi", boost::cref(mesons[PHI].getMass())));
100  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tKl", boost::cref(mesons[K_0].getWidth())));
101  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tKp", boost::cref(mesons[K_P].getWidth())));
102  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tBd", boost::cref(mesons[B_D].getWidth())));
103  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tBp", boost::cref(mesons[B_P].getWidth())));
104  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tBs", boost::cref(mesons[B_S].getWidth())));
105  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("DGs_Gs", boost::cref(mesons[B_S].getDgamma_gamma())));
106  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tKstar", boost::cref(mesons[K_star].getWidth())));
107  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("tphi", boost::cref(mesons[PHI].getWidth())));
108  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FK", boost::cref(mesons[K_0].getDecayconst())));
109  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FD", boost::cref(mesons[D_0].getDecayconst())));
110  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FBs", boost::cref(mesons[B_S].getDecayconst())));
111  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FKstar", boost::cref(mesons[K_star].getDecayconst())));
112  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FKstarp", boost::cref(FKstarp)));
113  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Fphi", boost::cref(mesons[PHI].getDecayconst())));
114  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("FBsoFBd", boost::cref(FBsoFBd)));
115  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK1", boost::cref(BK.getBpars()(0))));
116  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK2", boost::cref(BK.getBpars()(1))));
117  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK3", boost::cref(BK.getBpars()(2))));
118  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK4", boost::cref(BK.getBpars()(3))));
119  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK5", boost::cref(BK.getBpars()(4))));
120  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BKscale", boost::cref(BK.getMu())));
121  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BD1", boost::cref(BD.getBpars()(0))));
122  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BD2", boost::cref(BD.getBpars()(1))));
123  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BD3", boost::cref(BD.getBpars()(2))));
124  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BD4", boost::cref(BD.getBpars()(3))));
125  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BD5", boost::cref(BD.getBpars()(4))));
126  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BDscale", boost::cref(BD.getMu())));
127  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBsoBBd", boost::cref(BBsoBBd)));
128  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBs1", boost::cref(BBs.getBpars()(0))));
129  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBs2", boost::cref(BBs.getBpars()(1))));
130  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBs3", boost::cref(BBs.getBpars()(2))));
131  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBs4", boost::cref(BBs.getBpars()(3))));
132  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBs5", boost::cref(BBs.getBpars()(4))));
133  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BBsscale", boost::cref(BBs.getMu())));
134  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)1", boost::cref(BKd1.getBpars()(0))));
135  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)2", boost::cref(BKd1.getBpars()(1))));
136  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)3", boost::cref(BKd1.getBpars()(2))));
137  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)4", boost::cref(BKd1.getBpars()(3))));
138  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)5", boost::cref(BKd1.getBpars()(4))));
139  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)6", boost::cref(BKd1.getBpars()(5))));
140  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)7", boost::cref(BKd1.getBpars()(6))));
141  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)8", boost::cref(BKd1.getBpars()(7))));
142  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)9", boost::cref(BKd1.getBpars()(8))));
143  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(1/2)10", boost::cref(BKd1.getBpars()(9))));
144  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)1", boost::cref(BKd3.getBpars()(0))));
145  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)2", boost::cref(BKd3.getBpars()(1))));
146  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)3", boost::cref(BKd3.getBpars()(2))));
147  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)4", boost::cref(BKd3.getBpars()(3))));
148  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)5", boost::cref(BKd3.getBpars()(4))));
149  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)6", boost::cref(BKd3.getBpars()(5))));
150  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)7", boost::cref(BKd3.getBpars()(6))));
151  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)8", boost::cref(BKd3.getBpars()(7))));
152  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)9", boost::cref(BKd3.getBpars()(8))));
153  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BK(3/2)10", boost::cref(BKd3.getBpars()(9))));
154  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BKd_scale", boost::cref(BKd1.getMu())));
155  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("ReA2_Kd", boost::cref(ReA2_Kd)));
156  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("ReA0_Kd", boost::cref(ReA0_Kd)));
157  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Omega_eta_etap", boost::cref(Omega_eta_etap)));
158  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Br_Kp_P0enu", boost::cref(Br_Kp_P0enu)));
159  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Br_Kp_munu", boost::cref(Br_Kp_munu)));
160  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Br_B_Xcenu", boost::cref(Br_B_Xcenu)));
161  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("DeltaP_cu", boost::cref(DeltaP_cu)));
162  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("IB_Kl", boost::cref(IB_Kl)));
163  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("IB_Kp", boost::cref(IB_Kp)));
164  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_0", boost::cref(absh_0)));
165  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_p", boost::cref(absh_p)));
166  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_m", boost::cref(absh_m)));
167  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_0", boost::cref(argh_0)));
168  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_p", boost::cref(argh_p)));
169  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_m", boost::cref(argh_m)));
170  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_0_1", boost::cref(absh_0_1)));
171  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_p_1", boost::cref(absh_p_1)));
172  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_m_1", boost::cref(absh_m_1)));
173  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_0_1", boost::cref(argh_0_1)));
174  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_p_1", boost::cref(argh_p_1)));
175  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_m_1", boost::cref(argh_m_1)));
176  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_0_2", boost::cref(absh_0_2)));
177  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_p_2", boost::cref(absh_p_2)));
178  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_m_2", boost::cref(absh_m_2)));
179  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_0_2", boost::cref(argh_0_2)));
180  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_p_2", boost::cref(argh_p_2)));
181  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_m_2", boost::cref(argh_m_2)));
182  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_0_MP", boost::cref(absh_0_MP)));
183  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_0_MP", boost::cref(argh_0_MP)));
184  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("absh_0_1_MP", boost::cref(absh_0_1_MP)));
185  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("argh_0_1_MP", boost::cref(argh_0_1_MP)));
186  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0V", boost::cref(a_0V)));
187  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1V", boost::cref(a_1V)));
188  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2V", boost::cref(a_2V)));
189  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRV", boost::cref(MRV)));
190  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A0", boost::cref(a_0A0)));
191  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A0", boost::cref(a_1A0)));
192  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A0", boost::cref(a_2A0)));
193  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA0", boost::cref(MRA0)));
194  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A1", boost::cref(a_0A1)));
195  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A1", boost::cref(a_1A1)));
196  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A1", boost::cref(a_2A1)));
197  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA1", boost::cref(MRA1)));
198  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A12", boost::cref(a_0A12)));
199  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A12", boost::cref(a_1A12)));
200  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A12", boost::cref(a_2A12)));
201  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA12", boost::cref(MRA12)));
202  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T1", boost::cref(a_0T1)));
203  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T1", boost::cref(a_1T1)));
204  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T1", boost::cref(a_2T1)));
205  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT1", boost::cref(MRT1)));
206  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T2", boost::cref(a_0T2)));
207  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T2", boost::cref(a_1T2)));
208  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T2", boost::cref(a_2T2)));
209  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT2", boost::cref(MRT2)));
210  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T23", boost::cref(a_0T23)));
211  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T23", boost::cref(a_1T23)));
212  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T23", boost::cref(a_2T23)));
213  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT23", boost::cref(MRT23)));
214  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0Vphi", boost::cref(a_0Vphi)));
215  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1Vphi", boost::cref(a_1Vphi)));
216  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2Vphi", boost::cref(a_2Vphi)));
217  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRVphi", boost::cref(MRVphi)));
218  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A0phi", boost::cref(a_0A0phi)));
219  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A0phi", boost::cref(a_1A0phi)));
220  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A0phi", boost::cref(a_2A0phi)));
221  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA0phi", boost::cref(MRA0phi)));
222  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A1phi", boost::cref(a_0A1phi)));
223  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A1phi", boost::cref(a_1A1phi)));
224  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A1phi", boost::cref(a_2A1phi)));
225  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA1phi", boost::cref(MRA1phi)));
226  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0A12phi", boost::cref(a_0A12phi)));
227  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1A12phi", boost::cref(a_1A12phi)));
228  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2A12phi", boost::cref(a_2A12phi)));
229  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRA12phi", boost::cref(MRA12phi)));
230  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T1phi", boost::cref(a_0T1phi)));
231  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T1phi", boost::cref(a_1T1phi)));
232  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T1phi", boost::cref(a_2T1phi)));
233  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT1phi", boost::cref(MRT1phi)));
234  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T2phi", boost::cref(a_0T2phi)));
235  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T2phi", boost::cref(a_1T2phi)));
236  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T2phi", boost::cref(a_2T2phi)));
237  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT2phi", boost::cref(MRT2phi)));
238  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_0T23phi", boost::cref(a_0T23phi)));
239  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_1T23phi", boost::cref(a_1T23phi)));
240  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("a_2T23phi", boost::cref(a_2T23phi)));
241  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("MRT23phi", boost::cref(MRT23phi)));
242  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("r_1_fplus", boost::cref(r_1_fplus)));
243  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("r_2_fplus", boost::cref(r_2_fplus)));
244  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("m_fit2_fplus", boost::cref(m_fit2_fplus)));
245  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("r_1_fT", boost::cref(r_1_fT)));
246  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("r_2_fT", boost::cref(r_2_fT)));
247  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("m_fit2_fT", boost::cref(m_fit2_fT)));
248  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("r_2_f0", boost::cref(r_2_f0)));
249  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("m_fit2_f0", boost::cref(m_fit2_f0)));
250  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("bsgamma_E0", boost::cref(bsgamma_E0)));
251  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("BLNPcorr", boost::cref(BLNPcorr)));
252  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_mukin", boost::cref(Gambino_mukin)));
253  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_BRsem", boost::cref(Gambino_BRsem)));
254  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_Mbkin", boost::cref(Gambino_Mbkin)));
255  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_Mcatmuc", boost::cref(Gambino_Mcatmuc)));
256  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_mupi2", boost::cref(Gambino_mupi2)));
257  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_rhoD3", boost::cref(Gambino_rhoD3)));
258  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_muG2", boost::cref(Gambino_muG2)));
259  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("Gambino_rhoLS3", boost::cref(Gambino_rhoLS3)));
260  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("lambdaB", boost::cref(mesons[B_D].getLambdaM())));
261  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("alpha1kst", boost::cref(mesons[K_star].getGegenalpha(0))));
262  ModelParamMap.insert(std::pair<std::string, boost::reference_wrapper<const double> >("alpha2kst", boost::cref(mesons[K_star].getGegenalpha(1))));
263 
264  unknownParameterWarning = true;
265 }
Definition: QCD.h:717
A class for particles.
Definition: Particle.h:26
Definition: QCD.h:716
Definition: QCD.h:718
Definition: QCD.h:731
Definition: QCD.h:735
Definition: QCD.h:720
Definition: QCD.h:732
Definition: QCD.h:721
Definition: QCD.h:722
Definition: QCD.h:719
std::map< std::string, boost::reference_wrapper< const double > > ModelParamMap
Definition: Model.h:200

Member Function Documentation

double QCD::AboveTh ( const double  mu) const

The active flavour threshold above the scale \(\mu\) as defined in QCD::Thresholds().

Parameters
[in]mua scale \(\mu\) in GeV
Returns
the higher active flavour threshold

Definition at line 849 of file QCD.cpp.

850 {
851  int i;
852  for (i = 4; i >= 0; i--)
853  if (mu < Thresholds(i)) return (Thresholds(i));
854 
855  throw std::runtime_error("Error in QCD::AboveTh()");
856 }
double Thresholds(const int i) const
For accessing the active flavour threshold scales.
Definition: QCD.cpp:835
virtual double QCD::alphaMz ( ) const
pure virtual

Implemented in StandardModel.

double QCD::Als ( const double  mu,
const orders  order = FULLNLO 
) const

Computes the running strong coupling \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme. In the cases of LO, NLO and FULLNNLO, the coupling is computed with AlsWithInit(). On the other hand, in the cases of NNLO and FULLNNLO, the coupling is computed with AlsWithLambda().

Parameters
[in]mua scale \(\mu\) in GeV
[in]orderLO, NLO, FULLNLO, NNLO or FULLNNLO in the \(\alpha_s\) expansion defined in OrderScheme
Returns
the strong coupling constant \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme

Definition at line 1004 of file QCD.cpp.

1005 {
1006  int i;
1007  for (i = 0; i < CacheSize; ++i)
1008  if ((mu == als_cache[0][i]) && ((double) order == als_cache[1][i]) &&
1009  (AlsM == als_cache[2][i]) && (MAls == als_cache[3][i]) &&
1010  (mut == als_cache[4][i]) && (mub == als_cache[5][i]) &&
1011  (muc == als_cache[6][i]))
1012  return als_cache[7][i];
1013 
1014  double nfmu = Nf(mu), nfz = Nf(MAls), mu_thre1, mu_thre2, Als_tmp, mf;
1015  double als;
1016 
1017  switch (order) {
1018  case LO:
1019  case FULLNLO:
1020  case NLO:
1021  if (nfmu == nfz)
1022  als = AlsWithInit(mu, AlsM, MAls, order);
1023  else if (nfmu > nfz) {
1024  if (order == NLO)
1025  throw std::runtime_error("NLO is not implemented in QCD::Als(mu,order).");
1026  if (nfmu == nfz + 1.) {
1027  mu_thre1 = AboveTh(MAls); // mut
1028  Als_tmp = AlsWithInit(mu_thre1 - MEPS, AlsM, MAls, order);
1029  if (order == FULLNLO) {
1030  mf = mtpole;
1031  Als_tmp = (1. + Als_tmp / M_PI * log(mu_thre1 / mf) / 3.) * Als_tmp;
1032  }
1033  als = AlsWithInit(mu, Als_tmp, mu_thre1 + MEPS, order);
1034  } else
1035  throw std::runtime_error("Error in QCD::Als(mu,order)");
1036  } else {
1037  if (order == NLO)
1038  throw std::runtime_error("NLO is not implemented in QCD::Als(mu,order).");
1039  if (nfmu == nfz - 1.) {
1040  mu_thre1 = BelowTh(MAls); // mub
1041  Als_tmp = AlsWithInit(mu_thre1 + MEPS, AlsM, MAls, order);
1042  if (order == FULLNLO) {
1043  mf = getQuarks(BOTTOM).getMass();
1044  Als_tmp = (1. - Als_tmp / M_PI * log(mu_thre1 / mf) / 3.) * Als_tmp;
1045  }
1046  als = AlsWithInit(mu, Als_tmp, mu_thre1 - MEPS, order);
1047  } else if (nfmu == nfz - 2.) {
1048  mu_thre1 = BelowTh(MAls); // mub
1049  mu_thre2 = AboveTh(mu); // muc
1050  Als_tmp = Als(mu_thre1 + MEPS, order);
1051  if (order == FULLNLO) {
1052  mf = getQuarks(BOTTOM).getMass();
1053  Als_tmp = (1. - Als_tmp / M_PI * log(mu_thre1 / mf) / 3.) * Als_tmp;
1054  }
1055  Als_tmp = AlsWithInit(mu_thre2, Als_tmp, mu_thre1 - MEPS, order);
1056  if (order == FULLNLO) {
1057  mf = getQuarks(CHARM).getMass();
1058  Als_tmp = (1. - Als_tmp / M_PI * log(mu_thre2 / mf) / 3.) * Als_tmp;
1059  }
1060  als = AlsWithInit(mu, Als_tmp, mu_thre2 - MEPS, order);
1061  } else
1062  throw std::runtime_error("Error in QCD::Als(mu,order)");
1063  }
1064  break;
1065  case FULLNNLO:
1066  case NNLO:
1067  /* alpha_s(mu) computed with Lambda_QCD for Nf=nfmu */
1068  als = AlsWithLambda(mu, order);
1069  break;
1070  default:
1071  throw std::runtime_error(orderToString(order) + " is not implemented in QCD::Als(mu,order).");
1072  }
1073 
1074  CacheShift(als_cache, 8);
1075  als_cache[0][0] = mu;
1076  als_cache[1][0] = (double) order;
1077  als_cache[2][0] = AlsM;
1078  als_cache[3][0] = MAls;
1079  als_cache[4][0] = mut;
1080  als_cache[5][0] = mub;
1081  als_cache[6][0] = muc;
1082  als_cache[7][0] = als;
1083 
1084  return als;
1085 }
double AlsWithLambda(const double mu, const orders order) const
Computes the running strong coupling in the scheme with the use of .
Definition: QCD.cpp:970
double Als(const double mu, const orders order=FULLNLO) const
Computes the running strong coupling in the scheme. In the cases of LO, NLO and FULLNNLO...
Definition: QCD.cpp:1004
double BelowTh(const double mu) const
The active flavour threshold below the scale as defined in QCD::Thresholds().
Definition: QCD.cpp:858
double AboveTh(const double mu) const
The active flavour threshold above the scale as defined in QCD::Thresholds().
Definition: QCD.cpp:849
Definition: OrderScheme.h:33
Particle getQuarks(const quark q) const
A get method to access a quark as an object of the type Particle.
Definition: QCD.h:869
#define MEPS
Definition: QCD.h:14
double AlsWithInit(const double mu, const double alsi, const double mu_i, const orders order) const
Computes the running strong coupling from in the scheme, where it is forbidden to across a flavour...
Definition: QCD.cpp:903
const double & getMass() const
A get method to access the particle mass.
Definition: Particle.h:61
complex log(const complex &z)
double Nf(const double mu) const
The number of active flavour at scale .
Definition: QCD.cpp:867
std::string orderToString(const orders order) const
Converts an object of the enum type "orders" to the corresponding string.
Definition: QCD.cpp:267
double QCD::Als4 ( const double  mu) const

The value of \(\alpha_s^{\mathrm{FULLNLO}}\) at any scale \(\mu\) with the number of flavours \(n_f = 4\).

Parameters
[in]muthe scale at which \(\alpha_s\) has to be computed
Returns
\(\alpha_s^{\mathrm{FULLNLO}}(\mu)\) with \(n_f = 4\)
Attention
Temporary function waiting for the implementation of NNLO etact.

Definition at line 923 of file QCD.cpp.

924 {
925  double v = 1. - Beta0(4.) * AlsM / 2. / M_PI * log(MAls / mu);
926  return (AlsM / v * (1. - Beta1(4.) / Beta0(4.) * AlsM / 4. / M_PI * log(v) / v));
927 }
double Beta1(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:892
complex log(const complex &z)
double Beta0(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:887
double QCD::AlsWithInit ( const double  mu,
const double  alsi,
const double  mu_i,
const orders  order 
) const

Computes the running strong coupling \(\alpha_s(\mu)\) from \(\alpha_s(\mu_i)\) in the \(\overline{\mathrm{MS}}\) scheme, where it is forbidden to across a flavour threshold in the RG running from \(\mu_i\) to \(\mu\).

Parameters
[in]mua scale \(\mu\) in GeV
[in]alsithe initial value for the coupling at the scale given below
[in]mu_ithe initial scale \(\mu_i\) in GeV
[in]orderLO, NLO or FULLNLO in the \(\alpha_s\) expansion defined in OrderScheme
Returns
the strong coupling constant \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme

Definition at line 903 of file QCD.cpp.

905 {
906  double nf = Nf(mu);
907  if (nf != Nf(mu_i))
908  throw std::runtime_error("Error in QCD::AlsWithInit().");
909 
910  double v = 1. - Beta0(nf) * alsi / 2. / M_PI * log(mu_i / mu);
911  switch (order) {
912  case LO:
913  return (alsi / v);
914  case FULLNLO:
915  return (alsi / v * (1. - Beta1(nf) / Beta0(nf) * alsi / 4. / M_PI * log(v) / v));
916  case NLO:
917  return (alsi / v * (-Beta1(nf) / Beta0(nf) * alsi / 4. / M_PI * log(v) / v));
918  default:
919  throw std::runtime_error(orderToString(order) + " is not implemented in QCD::Als(mu,alsi,mi,order).");
920  }
921 }
double Beta1(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:892
Definition: OrderScheme.h:33
complex log(const complex &z)
double Nf(const double mu) const
The number of active flavour at scale .
Definition: QCD.cpp:867
std::string orderToString(const orders order) const
Converts an object of the enum type "orders" to the corresponding string.
Definition: QCD.cpp:267
double Beta0(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:887
double QCD::AlsWithLambda ( const double  mu,
const orders  order 
) const

Computes the running strong coupling \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme with the use of \(\Lambda_{\rm QCD}\).

Parameters
[in]muA scale \(\mu\) in GeV
[in]orderLO, NLO, FULLNLO, NNLO or FULLNNLO in the \(\alpha_s\) expansion defined in OrderScheme
Returns
the strong coupling constant \(\alpha_s(\mu)\) in the \(\overline{\mathrm{MS}}\) scheme

Definition at line 970 of file QCD.cpp.

972 {
973  double nf = Nf(mu);
974  double L = 2. * (log(mu) - logLambda);
975 
976  // LO contribution
977  double b0 = Beta0(nf);
978  double b0L = b0*L;
979  double alsLO = 4. * M_PI / b0L;
980  if (order == LO) return alsLO;
981 
982  // NLO contribution
983  double b1 = Beta1(nf);
984  double log_L = log(L);
985  double alsNLO = 4. * M_PI / b0L * (-b1 * log_L / b0 / b0L);
986  if (order == NLO) return alsNLO;
987  if (order == FULLNLO) return (alsLO + alsNLO);
988 
989  // NNLO contribution
990  double b2 = Beta2(nf);
991  double alsNNLO = 4. * M_PI / b0L * (1. / b0L / b0L
992  * (b1 * b1 / b0 / b0 * (log_L * log_L - log_L - 1.) + b2 / b0));
993  if (order == NNLO) return alsNNLO;
994  if (order == FULLNNLO) return (alsLO + alsNLO + alsNNLO);
995 
996  throw std::runtime_error(orderToString(order) + " is not implemented in QCD::AlsWithLambda().");
997 }
double logLambda(const double nf, orders order) const
Computes with nf flavours in GeV.
Definition: QCD.cpp:1195
double Beta1(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:892
Definition: OrderScheme.h:33
complex log(const complex &z)
double Nf(const double mu) const
The number of active flavour at scale .
Definition: QCD.cpp:867
std::string orderToString(const orders order) const
Converts an object of the enum type "orders" to the corresponding string.
Definition: QCD.cpp:267
double Beta0(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:887
double Beta2(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:897
double QCD::BelowTh ( const double  mu) const

The active flavour threshold below the scale \(\mu\) as defined in QCD::Thresholds().

Parameters
[in]mua scale \(\mu\) in GeV
Returns
the lower active flavour threshold

Definition at line 858 of file QCD.cpp.

859 {
860  int i;
861  for (i = 0; i < 5; i++)
862  if (mu >= Thresholds(i)) return (Thresholds(i));
863 
864  throw std::runtime_error("Error in QCD::BelowTh()");
865 }
double Thresholds(const int i) const
For accessing the active flavour threshold scales.
Definition: QCD.cpp:835
double QCD::Beta0 ( const double  nf) const

The \(\beta_0(n_f)\) coefficient for a certain number of flavours \(n_f\).

Parameters
[in]nfthe number of active flavours \(n_f\)
Returns
\(\beta_0(n_f)\)

Definition at line 887 of file QCD.cpp.

888 {
889  return ( (11. * Nc - 2. * nf) / 3.);
890 }
double QCD::Beta1 ( const double  nf) const

The \(\beta_1(n_f)\) coefficient for a certain number of flavours \(n_f\).

Parameters
[in]nfthe number of active flavours \(n_f\)
Returns
\(\beta_1(n_f)\)

Definition at line 892 of file QCD.cpp.

893 {
894  return ( 34. / 3. * Nc * Nc - 10. / 3. * Nc * nf - 2. * CF * nf);
895 }
double QCD::Beta2 ( const double  nf) const

The \(\beta_2(n_f)\) coefficient for a certain number of flavours \(n_f\).

Parameters
[in]nfthe number of active flavours \(n_f\)
Returns
\(\beta_2(n_f)\)

Definition at line 897 of file QCD.cpp.

898 {
899  return ( 2857. / 54. * Nc * Nc * Nc + CF * CF * nf - 205. / 18. * CF * Nc * nf
900  - 1415. / 54. * Nc * Nc * nf + 11. / 9. * CF * nf * nf + 79. / 54. * Nc * nf * nf);
901 }
bool QCD::CheckFlags ( ) const
virtual

A method to check the sanity of the set of model flags.

Returns
a boolean that is true if the set of model flags is sane

Implements Model.

Reimplemented in StandardModel.

Definition at line 828 of file QCD.cpp.

829 {
830  return (true);
831 }
bool QCD::CheckParameters ( const std::map< std::string, double > &  DPars)
virtual

A method to check if all the mandatory parameters for StandardModel have been provided in model initialization.

Parameters
[in]DParsa map of the parameters that are being updated in the Monte Carlo run (including parameters that are varied and those that are held constant)
Returns
a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, NPEffectiveGIMR, NPEffectiveBS, HiggsKvKfgen, NPZbbbarLinearized, NPZbbbar, NPHiggs, HiggsKvgenKf, NPEpsilons, HiggsKvKf, NPSTUZbbbarLR, THDM, NPEpsilons_pureNP, NPSTU, SUSY, SUSYMassInsertion, FlavourWilsonCoefficient, myModel, and GeneralSUSY.

Definition at line 804 of file QCD.cpp.

805 {
806  for (int i = 0; i < NQCDvars; i++)
807  if (DPars.find(QCDvars[i]) == DPars.end()) {
808  std::cout << "missing mandatory QCD parameter " << QCDvars[i] << std::endl;
809  return false;
810  }
811  return true;
812 }
static const int NQCDvars
The number of model parameters in QCD.
Definition: QCD.h:739
static const std::string QCDvars[NQCDvars]
An array containing the labels under which all QCD parameters are stored in a vector of ModelParamete...
Definition: QCD.h:745
double QCD::geta_0A0 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A0}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1200 of file QCD.h.

1201  {
1202  return a_0A0;
1203  }
double QCD::geta_0A0phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A0}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1424 of file QCD.h.

1425  {
1426  return a_0A0phi;
1427  }
double QCD::geta_0A1 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1232 of file QCD.h.

1233  {
1234  return a_0A1;
1235  }
double QCD::geta_0A12 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A12}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1264 of file QCD.h.

1265  {
1266  return a_0A12;
1267  }
double QCD::geta_0A12phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A12}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1488 of file QCD.h.

1489  {
1490  return a_0A12phi;
1491  }
double QCD::geta_0A1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{A1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1456 of file QCD.h.

1457  {
1458  return a_0A1phi;
1459  }
double QCD::geta_0T1 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1296 of file QCD.h.

1297  {
1298  return a_0T1;
1299  }
double QCD::geta_0T1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1520 of file QCD.h.

1521  {
1522  return a_0T1phi;
1523  }
double QCD::geta_0T2 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T2}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1328 of file QCD.h.

1329  {
1330  return a_0T2;
1331  }
double QCD::geta_0T23 ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T23}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1360 of file QCD.h.

1361  {
1362  return a_0T23;
1363  }
double QCD::geta_0T23phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T23}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1584 of file QCD.h.

1585  {
1586  return a_0T23phi;
1587  }
double QCD::geta_0T2phi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^{T2}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1552 of file QCD.h.

1553  {
1554  return a_0T2phi;
1555  }
double QCD::geta_0V ( ) const
inline
Returns
the LCSR fit parameter \(a_0^V\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1168 of file QCD.h.

1169  {
1170  return a_0V;
1171  }
double QCD::geta_0Vphi ( ) const
inline
Returns
the LCSR fit parameter \(a_0^V\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1392 of file QCD.h.

1393  {
1394  return a_0Vphi;
1395  }
double QCD::geta_1A0 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A0}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1208 of file QCD.h.

1209  {
1210  return a_1A0;
1211  }
double QCD::geta_1A0phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A0}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1432 of file QCD.h.

1433  {
1434  return a_1A0phi;
1435  }
double QCD::geta_1A1 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1240 of file QCD.h.

1241  {
1242  return a_1A1;
1243  }
double QCD::geta_1A12 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A12}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1272 of file QCD.h.

1273  {
1274  return a_1A12;
1275  }
double QCD::geta_1A12phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A12}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1496 of file QCD.h.

1497  {
1498  return a_1A12phi;
1499  }
double QCD::geta_1A1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{A1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1464 of file QCD.h.

1465  {
1466  return a_1A1phi;
1467  }
double QCD::geta_1T1 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1304 of file QCD.h.

1305  {
1306  return a_1T1;
1307  }
double QCD::geta_1T1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1528 of file QCD.h.

1529  {
1530  return a_1T1phi;
1531  }
double QCD::geta_1T2 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T2}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1336 of file QCD.h.

1337  {
1338  return a_1T2;
1339  }
double QCD::geta_1T23 ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T23}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1368 of file QCD.h.

1369  {
1370  return a_1T23;
1371  }
double QCD::geta_1T23phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T23}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1592 of file QCD.h.

1593  {
1594  return a_1T23phi;
1595  }
double QCD::geta_1T2phi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^{T2}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1560 of file QCD.h.

1561  {
1562  return a_1T2phi;
1563  }
double QCD::geta_1V ( ) const
inline
Returns
the LCSR fit parameter \(a_1^V\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1176 of file QCD.h.

1177  {
1178  return a_1V;
1179  }
double QCD::geta_1Vphi ( ) const
inline
Returns
the LCSR fit parameter \(a_1^V\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1400 of file QCD.h.

1401  {
1402  return a_1Vphi;
1403  }
double QCD::geta_2A0 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A0}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1216 of file QCD.h.

1217  {
1218  return a_2A0;
1219  }
double QCD::geta_2A0phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A0}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1440 of file QCD.h.

1441  {
1442  return a_2A0phi;
1443  }
double QCD::geta_2A1 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1248 of file QCD.h.

1249  {
1250  return a_2A1;
1251  }
double QCD::geta_2A12 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A12}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1280 of file QCD.h.

1281  {
1282  return a_2A12;
1283  }
double QCD::geta_2A12phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A12}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1504 of file QCD.h.

1505  {
1506  return a_2A12phi;
1507  }
double QCD::geta_2A1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{A1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1472 of file QCD.h.

1473  {
1474  return a_2A1phi;
1475  }
double QCD::geta_2T1 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1312 of file QCD.h.

1313  {
1314  return a_2T1;
1315  }
double QCD::geta_2T1phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1536 of file QCD.h.

1537  {
1538  return a_2T1phi;
1539  }
double QCD::geta_2T2 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T2}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1344 of file QCD.h.

1345  {
1346  return a_2T2;
1347  }
double QCD::geta_2T23 ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T23}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1376 of file QCD.h.

1377  {
1378  return a_2T23;
1379  }
double QCD::geta_2T23phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T23}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1600 of file QCD.h.

1601  {
1602  return a_2T23phi;
1603  }
double QCD::geta_2T2phi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^{T2}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1568 of file QCD.h.

1569  {
1570  return a_2T2phi;
1571  }
double QCD::geta_2V ( ) const
inline
Returns
the LCSR fit parameter \(a_2^V\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1184 of file QCD.h.

1185  {
1186  return a_2V;
1187  }
double QCD::geta_2Vphi ( ) const
inline
Returns
the LCSR fit parameter \(a_2^V\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1408 of file QCD.h.

1409  {
1410  return a_2Vphi;
1411  }
double QCD::getAlsM ( ) const
inline

A get method to access the value of \(\alpha_s(M_{\alpha_s})\).

Returns
the strong coupling constant at \(M_{\alpha_s}\), \(\alpha_s(M_{\alpha_s})\)

Definition at line 878 of file QCD.h.

879  {
880  return AlsM;
881  }
BParameter QCD::getBBd ( ) const
inline

For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta b = 2\) process in the \(B_d\) meson system.

Returns
the vector of bag parameters

Definition at line 943 of file QCD.h.

944  {
945  return BBd;
946  }
BParameter QCD::getBBs ( ) const
inline

For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta b = 2\) process in the \(B_s\) meson system.

Returns
the vector of bag parameters

Definition at line 954 of file QCD.h.

955  {
956  return BBs;
957  }
BParameter QCD::getBD ( ) const
inline

For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta c = 2\) process in the \(D^0\) meson system.

Returns
the vector of bag parameters

Definition at line 965 of file QCD.h.

966  {
967  return BD;
968  }
BParameter QCD::getBK ( ) const
inline

For getting the bag parameters corresponding to the operator basis \(O_1 -O_5\) in \(\Delta s = 2\) process in the \(K^0\) meson system.

Returns
the vector of bag parameters

Definition at line 976 of file QCD.h.

977  {
978  return BK;
979  }
BParameter QCD::getBKd1 ( ) const
inline
Returns

Definition at line 984 of file QCD.h.

985  {
986  return BKd1;
987  }
BParameter QCD::getBKd3 ( ) const
inline
Returns

Definition at line 992 of file QCD.h.

993  {
994  return BKd3;
995  }
double QCD::getBLNPcorr ( ) const
inline
Returns
the Benzke, Lee, Neubert, Paz non perturbative correction for \(b\to s \gamma\) from arXiv:1003.5012.

Definition at line 1688 of file QCD.h.

1689  {
1690  return BLNPcorr;
1691  }
double QCD::getBr_B_Xcenu ( ) const
inline
Returns
the experimental value for the branching ratio of \(B\to X_c e\nu\)

Definition at line 1042 of file QCD.h.

1043  {
1044  return Br_B_Xcenu;
1045  }
double QCD::getBr_Kp_munu ( ) const
inline
Returns
the experimental value for the branching ratio of \(K^+\to\mu^+\nu\)

Definition at line 1034 of file QCD.h.

1035  {
1036  return Br_Kp_munu;
1037  }
double QCD::getBr_Kp_P0enu ( ) const
inline
Returns
the experimental value for the branching ratio of \(K^+\to\pi^0e^+\nu\)

Definition at line 1026 of file QCD.h.

1027  {
1028  return Br_Kp_P0enu;
1029  }
double QCD::getbsgamma_E0 ( ) const
inline
Returns
the experimental energy cutoff for \(b\to s \gamma\)

Definition at line 1680 of file QCD.h.

1681  {
1682  return bsgamma_E0;
1683  }
double QCD::getCF ( ) const
inline

A get method to access the Casimir factor of QCD.

Returns
the Casimir factor

Definition at line 932 of file QCD.h.

933  {
934  return CF;
935  }
double QCD::getDeltaP_cu ( ) const
inline
Returns
the long-distance correction to the charm contribution of \(K^+\to\pi^+\nu\bar{\nu}\)

References: [Isidori et al.(2005)], [Buras et al.(2006)]

Definition at line 1054 of file QCD.h.

1055  {
1056  return DeltaP_cu;
1057  }
double QCD::getFKstarp ( ) const
inline
Returns
the decay constant of a transversely polarized \(K^*\) meson at 1 GeV

Definition at line 1760 of file QCD.h.

1761  {
1762  return FKstarp;
1763  }
double QCD::getGambino_BRsem ( ) const
inline
Returns
the fit value for the branching ratio of \(B\to X_c e\nu\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1704 of file QCD.h.

1705  {
1706  return Gambino_BRsem;
1707  }
double QCD::getGambino_Mbkin ( ) const
inline
Returns
the fit value for the kinetic b mass \(M_b^{\rm kin}(\mu^{\rm kin})\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1712 of file QCD.h.

1713  {
1714  return Gambino_Mbkin;
1715  }
double QCD::getGambino_Mcatmuc ( ) const
inline
Returns
the fit value for the MSbar mass \(M_c(\mu_c)\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1720 of file QCD.h.

1721  {
1722  return Gambino_Mcatmuc;
1723  }
double QCD::getGambino_muG2 ( ) const
inline
Returns
the fit value for \(\mu_G^2\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1744 of file QCD.h.

1745  {
1746  return Gambino_muG2;
1747  }
double QCD::getGambino_mukin ( ) const
inline
Returns
the kinetic scale used in arXiv:1411.6560.

Definition at line 1696 of file QCD.h.

1697  {
1698  return Gambino_mukin;
1699  }
double QCD::getGambino_mupi2 ( ) const
inline
Returns
the fit value for \(\mu_{\pi}^2\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1728 of file QCD.h.

1729  {
1730  return Gambino_mupi2;
1731  }
double QCD::getGambino_rhoD3 ( ) const
inline
Returns
the fit value for \(\rho_D^3\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1736 of file QCD.h.

1737  {
1738  return Gambino_rhoD3;
1739  }
double QCD::getGambino_rhoLS3 ( ) const
inline
Returns
the fit value for \(\rho_{LS}^3\) computed as in arXiv:1411.6560, but with \(\mu_c=2GeV\).

Definition at line 1752 of file QCD.h.

1753  {
1754  return Gambino_rhoLS3;
1755  }
gslpp::complex QCD::geth_0 ( ) const
inline
Returns
the constant term of the charm loop long distance contribution \(h_0\) to \(M\to V\)

Definition at line 1080 of file QCD.h.

1081  {
1082  return myh_0;
1083  }
gslpp::complex QCD::geth_0_1 ( ) const
inline
Returns
the linear term of the charm loop long distance contribution \(h_0\) to \(M\to V\)

Definition at line 1104 of file QCD.h.

1105  {
1106  return myh_0_1;
1107  }
gslpp::complex QCD::geth_0_1_MP ( ) const
inline
Returns
the linear term of the charm loop long distance contribution \(h_0\) to \(M\to P\)

Definition at line 1160 of file QCD.h.

1161  {
1162  return gslpp::complex(absh_0_1_MP,argh_0_1_MP,true);
1163  }
A class for defining operations on and functions of complex numbers.
Definition: gslpp_complex.h:35
gslpp::complex QCD::geth_0_2 ( ) const
inline
Returns
the quadratic term of the charm loop long distance contribution \(h_0\) to \(M\to V\)

Definition at line 1128 of file QCD.h.

1129  {
1130  return myh_0_2;
1131  }
gslpp::complex QCD::geth_0_MP ( ) const
inline
Returns
the constant term of the charm loop long distance contribution \(h_0\) to \(M\to P\)

Definition at line 1152 of file QCD.h.

1153  {
1154  return gslpp::complex(absh_0_MP,argh_0_MP,true);
1155  }
A class for defining operations on and functions of complex numbers.
Definition: gslpp_complex.h:35
gslpp::complex QCD::geth_m ( ) const
inline
Returns
the constant term of the charm loop long distance contribution \(h_-\) to \(M\to V\)

Definition at line 1096 of file QCD.h.

1097  {
1098  return myh_m;
1099  }
gslpp::complex QCD::geth_m_1 ( ) const
inline
Returns
the linear term of the charm loop long distance contribution \(h_-\) to \(M\to V\)

Definition at line 1120 of file QCD.h.

1121  {
1122  return myh_m_1;
1123  }
gslpp::complex QCD::geth_m_2 ( ) const
inline
Returns
the quadratic term of the charm loop long distance contribution \(h_-\) to \(M\to V\)

Definition at line 1144 of file QCD.h.

1145  {
1146  return myh_m_2;
1147  }
gslpp::complex QCD::geth_p ( ) const
inline
Returns
the constant term of the charm loop long distance contribution \(h_+\) to \(M\to V\)

Definition at line 1088 of file QCD.h.

1089  {
1090  return myh_p;
1091  }
gslpp::complex QCD::geth_p_1 ( ) const
inline
Returns
the linear term of the charm loop long distance contribution \(h_+\) to \(M\to V\)

Definition at line 1112 of file QCD.h.

1113  {
1114  return myh_p_1;
1115  }
gslpp::complex QCD::geth_p_2 ( ) const
inline
Returns
the quadratic term of the charm loop long distance contribution \(h_+\) to \(M\to V\)

Definition at line 1136 of file QCD.h.

1137  {
1138  return myh_p_2;
1139  }
double QCD::getIB_Kl ( ) const
inline
Returns
the isospin breaking corrections between \(K_L\to\pi^0\nu\bar{\nu}\) and \(K^+\to\pi^0 e^+\nu\)

Definition at line 1063 of file QCD.h.

1064  {
1065  return IB_Kl;
1066  }
double QCD::getIB_Kp ( ) const
inline
Returns
the isospin breaking corrections between \(K^+\to\pi^+ \nu\bar{\nu}\) and \(K^+\to\pi^0 e^+\nu\)

Definition at line 1072 of file QCD.h.

1073  {
1074  return IB_Kp;
1075  }
double QCD::getm_fit2_f0 ( ) const
inline
Returns
the LCSR fit parameter \(m_{fit}^{2,f_0}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1672 of file QCD.h.

1673  {
1674  return m_fit2_f0;
1675  }
double QCD::getm_fit2_fplus ( ) const
inline
Returns
the LCSR fit parameter \(m_{fit}^{2,f_+}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1632 of file QCD.h.

1633  {
1634  return m_fit2_fplus;
1635  }
double QCD::getm_fit2_fT ( ) const
inline
Returns
the LCSR fit parameter \(m_{fit}^{2,f_T}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1656 of file QCD.h.

1657  {
1658  return m_fit2_fT;
1659  }
double QCD::getMAls ( ) const
inline

A get method to access the mass scale \(M_{\alpha_s}\) at which the strong coupling constant measurement is provided.

Returns
the mass scale in GeV \(M_{\alpha_s}\) at which the strong coupling constant measurement is provided

Definition at line 887 of file QCD.h.

888  {
889  return MAls;
890  }
Meson QCD::getMesons ( const meson  m) const
inline

A get method to access a meson as an object of the type Meson.

Parameters
[in]mthe name of a meson
Returns
the object of the meson specified in the argument

Definition at line 859 of file QCD.h.

860  {
861  return mesons[m];
862  }
double QCD::getMRA0 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A0}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1224 of file QCD.h.

1225  {
1226  return MRA0;
1227  }
double QCD::getMRA0phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A0}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1448 of file QCD.h.

1449  {
1450  return MRA0phi;
1451  }
double QCD::getMRA1 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1256 of file QCD.h.

1257  {
1258  return MRA1;
1259  }
double QCD::getMRA12 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A12}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1288 of file QCD.h.

1289  {
1290  return MRA12;
1291  }
double QCD::getMRA12phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A12}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1512 of file QCD.h.

1513  {
1514  return MRA12phi;
1515  }
double QCD::getMRA1phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{A1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1480 of file QCD.h.

1481  {
1482  return MRA1phi;
1483  }
double QCD::getMRT1 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T1}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1320 of file QCD.h.

1321  {
1322  return MRT1;
1323  }
double QCD::getMRT1phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T1}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1544 of file QCD.h.

1545  {
1546  return MRT1phi;
1547  }
double QCD::getMRT2 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T2}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1352 of file QCD.h.

1353  {
1354  return MRT2;
1355  }
double QCD::getMRT23 ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T23}\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1384 of file QCD.h.

1385  {
1386  return MRT23;
1387  }
double QCD::getMRT23phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T23}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1608 of file QCD.h.

1609  {
1610  return MRT23phi;
1611  }
double QCD::getMRT2phi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^{T2}\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1576 of file QCD.h.

1577  {
1578  return MRT2phi;
1579  }
double QCD::getMRV ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^V\) for \(B\to K^*\) from arXiv:1503.05534v1

Definition at line 1192 of file QCD.h.

1193  {
1194  return MRV;
1195  }
double QCD::getMRVphi ( ) const
inline
Returns
the LCSR fit parameter \(\Delta m^V\) for \(B\to\phi\) from arXiv:1503.05534v1

Definition at line 1416 of file QCD.h.

1417  {
1418  return MRVphi;
1419  }
double QCD::getMtpole ( ) const
inline

A get method to access the pole mass of the top quark.

Returns
the pole mass of the top quark \(m_t^{pole}\)

Definition at line 923 of file QCD.h.

924  {
925  return mtpole;
926  }
double QCD::getMub ( ) const
inline

A get method to access the threshold between five- and four-flavour theory in GeV.

Returns
the threshold \(\mu_b\)

Definition at line 905 of file QCD.h.

906  {
907  return mub;
908  }
double QCD::getMuc ( ) const
inline

A get method to access the threshold between four- and three-flavour theory in GeV.

Returns
the threshold \(\mu_c\)

Definition at line 914 of file QCD.h.

915  {
916  return muc;
917  }
double QCD::getMut ( ) const
inline

A get method to access the threshold between six- and five-flavour theory in GeV.

Returns
the threshold \(\mu_t\)

Definition at line 896 of file QCD.h.

897  {
898  return mut;
899  }
double QCD::getNc ( ) const
inline

A get method to access the number of colours \(N_c\).

Returns
the number of colours

Definition at line 840 of file QCD.h.

841  {
842  return Nc;
843  }
double QCD::getOmega_eta_etap ( ) const
inline
Returns
the isospin breaking contribution in \(K^0\to\pi\pi\)

Definition at line 1018 of file QCD.h.

1019  {
1020  return Omega_eta_etap;
1021  }
Particle QCD::getQuarks ( const quark  q) const
inline

A get method to access a quark as an object of the type Particle.

Parameters
[in]qthe name of a quark
Returns
the object of the quark found in the argument

Definition at line 869 of file QCD.h.

870  {
871  return quarks[q];
872  }
double QCD::getr_1_fplus ( ) const
inline
Returns
the LCSR fit parameter \(r_1^{f_+}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1616 of file QCD.h.

1617  {
1618  return r_1_fplus;
1619  }
double QCD::getr_1_fT ( ) const
inline
Returns
the LCSR fit parameter \(r_1^{f_T}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1640 of file QCD.h.

1641  {
1642  return r_1_fT;
1643  }
double QCD::getr_2_f0 ( ) const
inline
Returns
the LCSR fit parameter \(r_2^{f_0}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1664 of file QCD.h.

1665  {
1666  return r_2_f0;
1667  }
double QCD::getr_2_fplus ( ) const
inline
Returns
the LCSR fit parameter \(r_2^{f_+}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1624 of file QCD.h.

1625  {
1626  return r_2_fplus;
1627  }
double QCD::getr_2_fT ( ) const
inline
Returns
the LCSR fit parameter \(r_2^{f_T}\) for \(B\to K\) from arXiv:hep-ph/0406232v1

Definition at line 1648 of file QCD.h.

1649  {
1650  return r_2_fT;
1651  }
double QCD::getReA0_Kd ( ) const
inline
Returns
the experimental value of the real part of the amplitude for \(K^0\to\pi\pi\) with \(\Delta I=0\)

Definition at line 1001 of file QCD.h.

1002  {
1003  return ReA0_Kd;
1004  }
double QCD::getReA2_Kd ( ) const
inline
Returns
the experimental value of the real part of the amplitude for \(K^0\to\pi\pi\) with \(\Delta I=2\)

Definition at line 1010 of file QCD.h.

1011  {
1012  return ReA2_Kd;
1013  }
bool QCD::Init ( const std::map< std::string, double > &  DPars)
virtual

Initializes the QCD parameters found in the argument.

Parameters
[in]DParsa map containing the parameters (all as double) to be used in Monte Carlo

Implements Model.

Reimplemented in StandardModel, THDM, SUSY, SUSYMassInsertion, FlavourWilsonCoefficient, myModel, and GeneralSUSY.

Definition at line 287 of file QCD.cpp.

288 {
289  bool check = CheckParameters(DPars);
290  if (!check) return (check);
291  check *= Update(DPars);
292  unknownParameterWarning = false;
293  return (check);
294 }
virtual bool Update(const std::map< std::string, double > &DPars)
The update method for QCD.
Definition: QCD.cpp:307
virtual bool CheckParameters(const std::map< std::string, double > &DPars)
A method to check if all the mandatory parameters for StandardModel have been provided in model initi...
Definition: QCD.cpp:804
double QCD::logLambda ( const double  nf,
orders  order 
) const

Computes \(\ln\Lambda_\mathrm{QCD}\) with nf flavours in GeV.

Parameters
[in]nfthe number of active flavours \(n_f\)
[in]orderLO, NLO, FULLNLO, NNLO or FULLNNLO in the \(\alpha_s\) expansion defined in OrderScheme
Returns
\(\ln\Lambda_\mathrm{QCD}\) with nf flavours in GeV

Definition at line 1195 of file QCD.cpp.

1198 {
1199  if (fabs(nfNEW - nfORG) != 1.)
1200  throw std::runtime_error("Error in QCD::logLambda()");
1201  if (order == NLO) order = FULLNLO;
1202  if (order == NNLO) order = FULLNNLO;
1203 
1204  /* We do not use the codes below for FULLNLO, since threshold corrections
1205  * can be regarded as an NNLO effect as long as setting the matching scale
1206  * to be close to the mass scale of the decoupling quark. In order to use
1207  * the relation als^{nf+1} = als^{nf} exactly, we use logLambdaNLO method.
1208  */
1209  if (order == FULLNLO)
1210  return logLambdaNLO(nfNEW, nfORG, logLambdaORG);
1211 
1212  double logMuMatching = log(muMatching);
1213  double L = 2. * (logMuMatching - logLambdaORG);
1214  double rNEW = 0.0, rORG = 0.0, log_mu2_mf2 = 0.0, log_L = 0.0;
1215  double C1 = 0.0, C2 = 0.0; // threshold corrections
1216  double logLambdaNEW;
1217 
1218  // LO contribution
1219  logLambdaNEW = 1. / 2. / Beta0(nfNEW)
1220  *(Beta0(nfNEW) - Beta0(nfORG)) * L + logLambdaORG;
1221 
1222  // NLO contribution
1223  if (order == FULLNLO || order == FULLNNLO) {
1224  rNEW = Beta1(nfNEW) / Beta0(nfNEW);
1225  rORG = Beta1(nfORG) / Beta0(nfORG);
1226  log_mu2_mf2 = 2. * (logMuMatching - log(mf));
1227  log_L = log(L);
1228  if (nfNEW < nfORG)
1229  C1 = 2. / 3. * log_mu2_mf2;
1230  else
1231  C1 = -2. / 3. * log_mu2_mf2;
1232  logLambdaNEW += 1. / 2. / Beta0(nfNEW)
1233  *((rNEW - rORG) * log_L
1234  - rNEW * log(Beta0(nfNEW) / Beta0(nfORG)) - C1);
1235  }
1236 
1237  // NNLO contribution
1238  if (order == FULLNNLO) {
1239  if (nfNEW == 5. && nfORG == 6.)
1240  C2 = -16. * (log_mu2_mf2 * log_mu2_mf2 / 36. - 19. / 24. * log_mu2_mf2 - 7. / 24.);
1241  else if (nfNEW == 6. && nfORG == 5.)
1242  C2 = -16. * (log_mu2_mf2 * log_mu2_mf2 / 36. + 19. / 24. * log_mu2_mf2 + 7. / 24.);
1243  else {
1244  if (nfNEW < nfORG)
1245  C2 = -16. * (log_mu2_mf2 * log_mu2_mf2 / 36. - 19. / 24. * log_mu2_mf2 + 11. / 72.);
1246  else
1247  C2 = -16. * (log_mu2_mf2 * log_mu2_mf2 / 36. + 19. / 24. * log_mu2_mf2 - 11. / 72.);
1248  }
1249  logLambdaNEW += 1. / 2. / Beta0(nfNEW) / Beta0(nfORG) / L
1250  * (rORG * (rNEW - rORG) * log_L + rNEW * rNEW - rORG * rORG
1251  - Beta2(nfNEW) / Beta0(nfNEW) + Beta2(nfORG) / Beta0(nfORG)
1252  + rNEW * C1 - C1 * C1 - C2);
1253  }
1254 
1255  return logLambdaNEW;
1256 }
double Beta1(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:892
complex log(const complex &z)
double Beta0(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:887
double Beta2(const double nf) const
The coefficient for a certain number of flavours .
Definition: QCD.cpp:897
double QCD::Nf ( const double  mu) const

The number of active flavour at scale \(\mu\).

Parameters
[in]mua scale \(\mu\) in GeV
Returns
active N_f

Definition at line 867 of file QCD.cpp.

868 {
869  int i;
870  for (i = 1; i < 5; i++)
871  if (mu >= Thresholds(i))
872  return (7. - (double) i);
873 
874  throw std::runtime_error("Error in QCD::Nf()");
875 }
double Thresholds(const int i) const
For accessing the active flavour threshold scales.
Definition: QCD.cpp:835
std::string QCD::orderToString ( const orders  order) const

Converts an object of the enum type "orders" to the corresponding string.

Parameters
[in]orderan object of the enum type "orders"
Returns
the string of the given "order"

Definition at line 267 of file QCD.cpp.

268 {
269  switch (order) {
270  case LO:
271  return "LO";
272  case NLO:
273  return "NLO";
274  case FULLNLO:
275  return "FULLNLO";
276  case NNLO:
277  return "NNLO";
278  case FULLNNLO:
279  return "FULLNNLO";
280  default:
281  throw std::runtime_error(orderToString(order) + " is not implemented in QCD::orderToString().");
282  }
283 }
Definition: OrderScheme.h:33
std::string orderToString(const orders order) const
Converts an object of the enum type "orders" to the corresponding string.
Definition: QCD.cpp:267
bool QCD::PostUpdate ( )
virtual

The post-update method for QCD.

This method runs all the procedures that are need to be executed after the model is successfully updated. This includes

  • computing the decay constant \(F_{B_D}\) from \(F_{B_s}\)
  • computing the bag parameters \(B_{B_d}\) from \(B_{B_s}\)
  • computing the \(\overline{\rm MS}\) mass of the top quark at the \(\overline{\rm MS}\) mass, \(m_t^{\overline{\rm MS}}(m_t^{\overline{\rm MS}})\) and setting the scale at the same value.
    Returns
    a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, NPEffectiveGIMR, NPZbbbarLinearized, NPZbbbar, NPEpsilons, THDM, NPbase, SUSY, SUSYMassInsertion, FlavourWilsonCoefficient, myModel, and GeneralSUSY.

Definition at line 323 of file QCD.cpp.

324 {
325  if (computeFBd)
326  mesons[B_D].setDecayconst(mesons[B_S].getDecayconst() / FBsoFBd);
327  if (computeBd)
328  BBd.setBpars(0, BBs.getBpars()(0) / BBsoBBd);
329  if (computemt) {
330  quarks[TOP].setMass(Mp2Mbar(mtpole, FULLNNLO));
331  quarks[TOP].setMass_scale(quarks[TOP].getMass());
332  }
333 
334  myh_0 = gslpp::complex(absh_0,argh_0,true);
335  myh_p = gslpp::complex(absh_p,argh_p,true);
336  myh_m = gslpp::complex(absh_m,argh_m,true);
337  myh_0_1 = gslpp::complex(absh_0_1,argh_0_1,true);
338  myh_p_1 = gslpp::complex(absh_p_1,argh_p_1,true);
339  myh_m_1 = gslpp::complex(absh_m_1,argh_m_1,true);
340  myh_0_2 = gslpp::complex(absh_0_2,argh_0_2,true);
341  myh_p_2 = gslpp::complex(absh_p_2,argh_p_2,true);
342  myh_m_2 = gslpp::complex(absh_m_2,argh_m_2,true);
343 
344  return (true);
345 }
Definition: QCD.h:735
Definition: QCD.h:721
Definition: QCD.h:719
A class for defining operations on and functions of complex numbers.
Definition: gslpp_complex.h:35
bool QCD::PreUpdate ( )
virtual

The pre-update method for QCD.

This method resets the internal flags #requireYu, #requireYd, #computeBd, #computeFBd and #computemt before updating the model parameters with the method Update().

Returns
a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, THDM, SUSY, SUSYMassInsertion, FlavourWilsonCoefficient, myModel, and GeneralSUSY.

Definition at line 296 of file QCD.cpp.

297 {
298  requireYu = false;
299  requireYd = false;
300  computeBd = false;
301  computeFBd = false;
302  computemt = false;
303 
304  return (true);
305 }
bool QCD::setFlag ( const std::string  name,
const bool  value 
)
virtual

A method to set a flag of QCD.

Parameters
[in]namename of a model flag
[in]valuethe boolean to be assigned to the flag specified by name
Returns
a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, NPEffectiveGIMR, NPEpsilons, SUSY, FlavourWilsonCoefficient, and myModel.

Definition at line 816 of file QCD.cpp.

817 {
818  std::cout << "WARNING: wrong name or value for ModelFlag " << name << std::endl;
819  return (false);
820 }
std::string name
The name of the model.
Definition: Model.h:203
bool QCD::setFlagStr ( const std::string  name,
const std::string  value 
)
virtual

A method to set a flag of QCD.

Parameters
[in]namename of a model flag
[in]valuethe string to be assigned to the flag specified by name
Returns
a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, and THDM.

Definition at line 822 of file QCD.cpp.

823 {
824  std::cout << "WARNING: wrong name or value for ModelFlag " << name << std::endl;
825  return (false);
826 }
std::string name
The name of the model.
Definition: Model.h:203
void QCD::setNc ( double  Nc)
inline

A set method to change the number of colours \(N_c\).

Parameters
[in]Ncthe number of colours

Definition at line 849 of file QCD.h.

850  {
851  this->Nc = Nc;
852  }
double QCD::Thresholds ( const int  i) const

For accessing the active flavour threshold scales.

Parameters
[in]ithe index referring to active flavour thresholds.
Returns
the threshold scale: 1.0E10 (i = 0), \(\mu_t\) (i = 1), \(\mu_b\) (i = 2), \(\mu_c\) (i = 3) and 0. (default)

Definition at line 835 of file QCD.cpp.

836 {
837  if (!(mut > mub && mub > muc))
838  throw std::runtime_error("inverted thresholds in QCD::Thresholds()!");
839 
840  switch (i) {
841  case 0: return (1.0E10);
842  case 1: return (mut);
843  case 2: return (mub);
844  case 3: return (muc);
845  default: return (0.);
846  }
847 }
bool QCD::Update ( const std::map< std::string, double > &  DPars)
virtual

The update method for QCD.

This method updates all the model parameters with given DPars.

Parameters
[in]DParsa map of the parameters that are being updated in the Monte Carlo run (including parameters that are varied and those that are held constant)
Returns
a boolean that is true if the execution is successful

Implements Model.

Reimplemented in StandardModel, THDM, SUSY, SUSYMassInsertion, FlavourWilsonCoefficient, myModel, and GeneralSUSY.

Definition at line 307 of file QCD.cpp.

308 {
309  if (!PreUpdate()) return (false);
310 
311  UpdateError = false;
312 
313  for (std::map<std::string, double>::const_iterator it = DPars.begin(); it != DPars.end(); it++)
314  setParameter(it->first, it->second);
315 
316  if (UpdateError) return (false);
317 
318  if (!PostUpdate()) return (false);
319 
320  return (true);
321 }
virtual bool PostUpdate()
The post-update method for QCD.
Definition: QCD.cpp:323
bool UpdateError
A boolean set to false if update is successful.
Definition: Model.h:192
virtual bool PreUpdate()
The pre-update method for QCD.
Definition: QCD.cpp:296
virtual void setParameter(const std::string name, const double &value)=0
A method to set the value of a parameter of the model.

Member Data Documentation

const int QCD::NQCDvars = 186
static

The number of model parameters in QCD.

Definition at line 739 of file QCD.h.

const std::string QCD::QCDvars
static

An array containing the labels under which all QCD parameters are stored in a vector of ModelParameter via InputParser::ReadParameters().

Definition at line 745 of file QCD.h.


The documentation for this class was generated from the following files: