pQCD ingredients

Strong Coupling

Implementation: StrongCoupling.

We use perturbative QCD with the running coupling \(a_s(\mu_R^2) = \alpha_s(\mu_R^2)/(4\pi)\) given at 5-loop by [HRU+17] [LMMS16] [BCK17]

\[\frac{da_s(\mu_R^2)}{d\ln\mu_R^2} = \beta(a_s(\mu_R^2)) \ = - \sum\limits_{n=0} \beta_n a_s^{n+2}(\mu_R^2)\]

It is useful to define in addition \(b_k = \beta_k/\beta_0, k>0\).

We implement two different strategies to solve the RGE:

method="exact": Solve using scipy.integrate.solve_ivp(). In LO we fall back to the expanded solution as this is already the true solution.
method="expanded": using approximate solutions:

\[\begin{split}a^{\text{LO}}_s(\mu_R^2) &= \frac{a_s(\mu_0^2)}{1 + a_s(\mu_0^2) \beta_0 \ln(\mu_R^2/\mu_0^2)} \\ a^{\text{NLO}}_{s,\text{exp}}(\mu_R^2) &= a^{\text{LO}}_s(\mu_R^2)-b_1 \left[a^{\text{LO}}_s(\mu_R^2)\right]^2 \ln\left(1+a_s(\mu_0^2) \beta_0 \ln(\mu_R^2/\mu_0^2)\right) \\ a^{\text{NNLO}}_{s,\text{exp}}(\mu_R^2) &= a^{\text{LO}}_s(\mu_R^2)\left[1 + a^{\text{LO}}_s(\mu_R^2)\left(a^{\text{LO}}_s(\mu_R^2) - a_s(\mu_0^2)\right)(b_2 - b_1^2) \right.\\ & \hspace{60pt} \left. + a^{\text{NLO}}_{s,\text{exp}}(\mu_R^2) b_1 \ln\left(a^{\text{NLO}}_{s,\text{exp}}(\mu_R^2)/a_s(\mu_0^2)\right)\right]\end{split}\]

When the renormalization scale crosses a flavor threshold matching conditions have to be applied [CKS06, SS06].

Splitting Functions

The Altarelli-Parisi splitting kernels can be expanded in powers of the strong coupling \(a_s(\mu^2)\) and are given by [MVV04, VMV04]

\[\begin{split}\mathbf{P}(x,a_s(\mu^2)) &= \sum\limits_{j=0} a_s^{j+1}(\mu^2) \mathbf{P}^{(j)}(x) \\ {\gamma}^{(j)}(N) &= -\mathcal{M}[\mathbf{P}^{(j)}(x)](N)\end{split}\]

Note the additional minus in the definition of \(\gamma\).

Scale Variations

The usual procedure in solving DGLAP that is also applied here is to rewrite the equations in term of the running coupling \(a_s\) assuming the factorization scale \(\mu_F^2\) (the inherit scale of the PDF) and the renormalization scale \(\mu_R^2\) (the inherit scale for the strong coupling) to be equal. This constraint can however be lifted by a suitable redefinition of the splitting kernels [Vog05]:

\[\begin{split}\gamma^{(1)}(N) &\to \gamma^{(1)}(N) - \beta_0 \ln(\mu_F^2/\mu_R^2) \gamma^{(0)} \\ \gamma^{(2)}(N) &\to \gamma^{(2)}(N) - 2 \beta_0 \ln(\mu_F^2/\mu_R^2) \gamma^{(1)} - ( \beta_1 \ln(\mu_F^2/\mu_R^2) - \beta_0^2 \ln^2(\mu_F^2/\mu_R^2) ) \gamma^{(0)}\end{split}\]

while keeping the evaluation of the strong coupling always at \(\mu_R^2\). Estimating the theoretical uncertainties imposed on PDF determination due to missing higher order corrections using scale variation in the evolution corresponds to schemes A and B in [AK+19].

Heavy Quark Masses

In QCD also the heavy quark masses (\(m_{c}, m_{b}, m_{t}\)) follow a RGE and their values depend on the energy scale at which the quark is probed. Masses do not play any role in a single flavour patch, but are important in VFNS when more flavour schemes need to be joined (see matching conditions).

EKO implements two strategies for dealing with the heavy quark masses, managed by the theory card parameter HQ. The easiest and more common option for PDFs evolution is POLE mass, where the physical quark masses are specified as input.

On contrary selecting the option MSBAR the user can activate the mass running in the \(\overline{MS}\) scheme, as described in the following paragraph.

If the initial condition for the mass is not given at a scale coinciding with the mass itself (i.e. in the input theory card Qmh≠mh), EKO needs to compute the scale at which the mass running function intersects the identity function, in order to properly initiate the ThresholdAtlas and set the evolution path.

For each heavy quark \(h\) we solve for \(m_h\):

\[m_{\overline{MS},h}(m_h^2) = m_h\]

where the evolved \(\overline{MS}\) mass is given by:

\[m_{\overline{MS},h}(\mu^2) = m_{h,0} \int_{a_s(\mu_{h,0}^2)}^{a_s(\mu^2)} \frac{\gamma(a_s)}{\beta(a_s)} d a_s\]

and \(m_{h,0}\) is the given initial condition at the scale \(\mu_{h,0}\). Here there is a subtle complication since the solution depends on the value \(a_s(\mu_{h,0}^2)\) which is unknown and depends again on the threshold path. To overcome this issue, EKO initialize a temporary instance of the class StrongCoupling with a fixed flavor number scheme, with \(n_{f_{ref}}\) active flavors at the scale \(\mu_{ref}\).

Then we check that, heavy quarks involving a number of active flavors greater than \(n_{f_{ref}}\) are given with initial conditions:

\[m_h (\mu_h) \ge \mu_h\]

while the ones related to fewer active flavors follow:

\[m_h (\mu_h) \le \mu_h\]

So for the former initial condition we will find the intercept between RGE and the identity in the forward direction (\(m_{\overline{MS},h} \ge \mu_h\)) and vice versa for the latter.

In doing so EKO takes advantage of the monotony of the RGE solution \(m_{\overline{MS},h}(\mu^2)\) with a vanishing limit for \(\mu^2 \rightarrow \infty\).

Now, being able to evaluate \(a_s(\mu_{h,0}^2)\), there are two ways of solving the previous integral and finally compute the evolved \(m_{\overline{MS},h}\). In fact, the function \(\gamma(a_s)\) is the anomalous QCD mass dimension and, as the \(\beta\) function, it can be evaluated perturbatively in \(a_s\) up to \(\mathcal{O}(a_s^3)\):

\[\begin{split}\gamma(a_s) &= - \sum\limits_{n=0} \gamma_n a_s^{n+1} \\\end{split}\]

Even here it is useful to define \(c_k = \gamma_k/\beta_0, k>0\).

Therefore the two solution strategies are:

method = "exact": the integral is solved exactly using the expression of \(\beta,\gamma\) up to the specified perturbative order
method = "expanded": the integral is approximate by the following expansion:

\[\begin{split}m_{\overline{MS},h}(\mu^2) & = m_{h,0} \left ( \frac{a_s(\mu^2)}{a_s(\mu_{h,0}^2)} \right )^{c_0} \frac{j_{exp}(a_s(\mu^2))}{j_{exp}(a_s(\mu_{h,0}^2))} \\ j_{exp}(a_s) &= 1 + a_s \left [ c_1 - b_1 c_0 \right ] + \frac{a_s^2}{2} \left [c_2 - c_1 b_1 - b_2 c_0 + b_1^2 c_0 + (c_1 - b_1 c_0)^2 \right]\end{split}\]

The procedure is iterated on all the heavy quarks, updating the temporary instance of StrongCoupling with the computed masses.

To find coherent solutions and perform the mass running in the correct patches it is necessary to always start computing the mass scales closer to \(\mu_{ref}\).

Eventually, to ensure that the threshold values are properly set, we add a consistency check, asserting:

\[m_{\overline{MS},h} (m_h) \leq m_{\overline{MS},h+1} (m_h)\]

We provide the following as an illustrative example of how this procedure works: when the strong coupling is given with boundary condition \(\alpha_s(\mu_{ref}=91, n_{f_{ref}}=5)\) then the heavy quarks initial conditions must satisfy:

\[\begin{split}& \mu_{b} \le \mu_{ref} \le \mu_t \\ & m_c (\mu_c) \le \mu_c \\ & m_b (\mu_b) \le \mu_b \\ & m_t (\mu_t) \ge \mu_t\end{split}\]

and EKO will start solving the equation \(m_{\overline{MS},h}(m_h^2) = m_h\) in the order \(h={t,b,c}\).

Since the charm mass will be computed only when both the top and bottom threshold scales are known, the boundary condition \(m_c(\mu_{c})\) can be evolved safely below the scale \(m_{\overline{MS},b}\) where the solution of \(m_{\overline{MS},c}(m_c^2) = m_c\) is sitting.