The Polarized Transition Matrix Element $A_{gq}(N)$ of the Variable Flavor Number Scheme at $O(\alpha_s^3)$

We calculate the polarized massive operator matrix element $A_{gq}^{(3)}(N)$ to 3-loop order in Quantum Chromodynamics analytically at general values of the Mellin variable $N$ both in the single- and double-mass case in the Larin scheme. It is a transition function required in the variable flavor number scheme at $O(\alpha_s^3)$. We also present the results in momentum fraction space.


Introduction
The variable flavor number scheme (VFNS) can be used to translate twist-2 parton distributions from a scheme with N F light flavors to a scheme with N F + 1 light flavors at a scale µ 2 . Thus, it allows for a process-independent description of the transition from a massive quark to a massless quark. In the single heavy mass case this has been worked out to 2-loop order in [1] and to 3-loop order in [2]. In terms of the N F -flavor distributions the new (N F + 1)-flavor massless parton densities are given by , The quark and antiquark parton densities are denoted by f k and fk, respectively, and G(N F , µ 2 ) is the gluon density. We write Σ(N F , µ 2 ) = N F k=1 (f k + fk) for the singlet-quark density. The massive operator matrix elements (OMEs) A ij (N F , m 2 /µ 2 ) are process-independent quantities and have an expansion in the strong coupling constant a s = α s /(4π), Here, µ denotes the decoupling scale and m is the mass of the decoupling heavy-quark flavor Q. The OMEs explicitly depend on the mass through logarithms. In total, there are seven different OMEs contributing to the matching relations. Moreover, we introduce the shorthand notations The VFNS is important for obtaining parton distribution functions at very large virtualities as they are required for, e.g., scattering processes at the Large Hadron Collider (LHC), the Tevatron, RHIC and the EIC, in particular for precision measurements of observables in Quantum Chromodynamics (QCD) and the determination of the strong coupling constant α s (M Z ) [3]. The relations (1.1) apply structurally to both the unpolarized and polarized case. The corresponding relations have to be generalized in the case of two heavy-quark contributions, as for charm and bottom quarks, the masses of which are very similar, cf. [4][5][6]. We also provide the 2-mass 3-loop OME A (3),two−mass gq , for which a closed form expression can be obtained [7].
The individual OMEs start contributing at different orders: The OMEÃ S Qg starts already at 1-loop order, whileÃ PS Qq , A NS qq,Q , A S gg,Q and A S gq,Q only contribute from 2-loop order onward. At 3-loop order, alsoÃ PS qq,Q andÃ S qg,Q appear. The complete 2-loop corrections were calculated in Refs. [1,6,[8][9][10][11][12][13] for the unpolarized and polarized cases. Mellin moments for the OMEs at 3-loop order were calculated in [2,14] and for transversity A NS,TR qq,Q in [15]. Depending on the process, the moments reached up to N = 13 at most.
In [16], the unpolarized OMEsÃ PS qq,Q (N) andÃ S qg,Q (N), as well as the have been calculated at 3-loop order. Here, denote the color factors for the gauge group SU(N c ). For A S gg,Q and A S gq,Q the corresponding contributions were calculated in [17]. Furthermore, the unpolarized 3-loop OMEs A NS qq,Q ,Ã PS Qq , A S gq,Q and A S gg,Q were obtained in [18][19][20][21][22] and the logarithmic contributions to 3-loop order in Ref. [23]. Partial results forÃ (3) Qg were calculated in [24]. The massive 3-loop polarized OMEs are known in the non-singlet and pure-singlet cases [18,25]. The two-mass contributions up to 3-loop order are known in the unpolarized case [4,5,26,27] and in the polarized case in [28,29], in both cases up toÃ (3) Qg . Heavy-flavor contributions to charged current processes up to 3-loop order and the polarized structure function g NS 1 (x, Q 2 ) were dealt with in [30][31][32][33].
In the present paper we compute the complete polarized OME A S, (3) gq,Q (N) for general values of N in both the single-and two-mass cases. Note that we will drop the superscript S in the following, since for this OME only the singlet part contributes. Due to the crossing relations, cf. [34], only the odd moments contribute and they are used to construct the analytic continuation to complex values of N or the Bjorken x-space, respectively. From the O(1/ε) pole term one obtains the contribution to the 3-loop polarized anomalous dimension γ (2) gq ∝ T F , cf. [35,36]. We perform the calculation using the Larin scheme [37], which is a consistent scheme w.r.t. the γ 5 problem, see also [35]. It is convenient to work in this scheme also for the description of the observables in polarized deep-inelastic scattering. This requires that also the polarized parton distribution functions are evolved using this scheme and that one calculates the massless Wilson coefficients in this scheme. Observables, such as the deep-inelastic structure functions, are then scheme-independent. The OME A (3) gq,Q (N) requires to use a special projector for the external quark lines, which has first been derived in Ref. [35], Eq. (11).
The paper is organized as follows. We discuss technical details of the calculation in Section 2. In Section 3, we present the constant part of the unrenormalized single-mass 3-loop OME A (3) gq,Q in Mellin-N space and discuss its small-and large-x behavior. Section 4 is devoted to the analytic calculation of the OME A . Section 5 contains the conclusions. In the Appendix, we present the polarized single-and two-mass OME A (3) gq,Q in Mellin-and momentum-fraction space, treating the heavy-quark mass in both in the on-shell and MS scheme.

The Formalism
Concerning the formalism, we follow closely Ref. [21], in which the corresponding result in the unpolarized case has been calculated. Renormalizing the heavy-quark mass in the on-shell scheme and the coupling constant in the MS scheme, the massive operator matrix element A (3) gq,Q has the gq,Q . ij up to three-loop order (i.e. k = 0, 1, 2) [35,36], the expansion coefficients of the QCD β-function [38][39][40][41][42][43][44], terms from mass renormalization [45] and the constant parts of the unrenormalized massive OMEs a (k) ij in the polarized case [6,[11][12][13], again up to three-loop order. The expansion coefficients related to the QCD β-function are given by 3) and the terms stemming from mass renormalization read From the logarithmic terms one can extract the complete 2-loop anomalous dimension γ (1) gq [35,36,46,47] and the contributions ∝ T F of γ (2) gq (N), one of the anomalous dimensions at 3-loop order [35,36].
We use a well established approach to calculate the 86 contributing Feynman diagrams. First the diagrams are generated with an extension of QGRAF [48] which can deal with local operator insertions [2]. The Feynman rules are then inserted in TFORM [49] where also the Dirac-and colortraces are calculated. The local operator insertions are resummed into generating functions using the auxiliary variable t, cf. [50]. This introduces on top of the usual denominators from particle propagators denominators which depend linearly on the loop momenta and the variable t. The scalar integrals are subsequently reduced to a minimal set of master integrals using the implementation of integration-by-parts reduction [51] in Reduze2 [52], which can also deal with linear propagators. The solutions of the master integrals are obtained using standard techniques. 1 This includes methods based on hypergeometric functions [9,50,[54][55][56][57], Mellin-Barnes representations [58] and differential equations [59,60]. For the analytic continuation of Mellin-Barnes integrals, the packages MB [61] and MBresolve [62] were used. When applying methods based on direct integration and Mellin-Barnes representations the results are typically given by multiple sums over hypergeometric expressions which can still depend on the dimensional parameter ε = D −4. These expressions can be expanded in ε and the resulting sums can afterwards be performed utilizing modern summation technology [63][64][65][66][67][68][69][70][71] as encoded in the packages Sigma [72,73], HarmonicSums [74][75][76], EvaluateMultiSums, SumProduction [77], and ρ-Sum [78]. For one of the master integrals it was essential to apply the multivariate Almkvist-Zeilberger algorithm [79] as implemented in the package MultiIntegrate [59,75] on the Mellin-space representation of the master integral. This way we were able to directly compute a difference equation for the Mellin-space result which we solved using the same summation technology cited before. The final results for individual master integrals, diagrams and the full final result have been checked by computing a number of integer moments with MATAD [80].
As in the unpolarized case, the polarized OME A gq,Q can be completely expressed by harmonic sums S a (N) and ζ-values [81] in Mellin-space and harmonic polylogarithms H a (x) and ζ-values in Bjorken x-space. The definitions of harmonic sums and harmonic polylogarithms are given by the iterative formulas [82] and [83] with a ∈ {−1, 0, 1} and The full renormalization and mass factorization of all massive operator matrix elements including A (3) gq,Q up to 3-loop order has been presented in Ref. [2] for the single mass case and in Ref. [4] for the two-mass case. The necessary steps are the renormalization of the masses, the coupling constant and the twist-2 light cone operators. Furthermore, collinear singularities have to be removed by mass factorization. Contrary to the massless case, the Z-factors related to the ultraviolet renormalization in the massive case are not inverse to those describing the collinear singularities. Moreover, the coupling constant is first renormalized in a MOM scheme using the background-field method [84] and afterwards translated to the usual MS scheme in order to fulfill the on-shell condition of the external partonic states.

The Single-Mass Correction
In the single-mass case, the renormalized OME (2.1) can be expressed in terms of lower-order terms as well as the newly evaluated constant part a (3) gq (N) of the unrenormalized OME. We definep and use the shorthand notation S k (N) ≡ S k . One obtains with the polynomials 10) (3.14)

15)
The color factors in QCD take on the values C A = 3, C F = 4/3, T F = 1/2. The nested sums and constants in a appears. Here, we write Li k (x) = ∞ l=1 x l /l k , |x| ≤ 1 for the classical polylogarithm. As a cross-check, we compared this result to an independent calculation of the moments N = 3, 5, 7 using MATAD [80] and we find agreement.
The nested sums can be mapped to a basis using algebraic reduction [85] after which only the sums appear. In addition, structural relations, such as multiple argument relations and differentiation, [86], can be applied, which leaves us only with as basic sums. At N = 1, the OME a gq (N) has a removable singularity. In this limit the expression becomes This agrees with the expectation that the rightmost singularity for gluonic OMEs occurs at N = 0. As was observed in [87], removable singularities can also appear at rational values of N > 0 for massive OMEs. The OME A gq (N) is a meromorphic function [86] since it can be expressed in terms of harmonic sums [82] over Q(N) with rational weights whose denominators factorize, with factors of the form (N − k) l , k ∈ Z, l ∈ N. The poles of this OME are located at negative integers, N ≤ 0.
We now turn our discussion to the behavior of a Around N → ∞ the constant part of the OME the asymptotic behavior is given by where we abbreviate L(N) = ln(N) + γ E and we write γ E for the Euler-Mascheroni constant.
The N → ∞ limit in N-space corresponds to the x → 1 limit in x-space. This allows us to deduce that in this limit the leading singular term is ∝ α 3 s ln 4 (1 − x). The position of the so-called 'leading-poles' can be inferred from an analysis of the anomalous dimensions of different scattering processes in fixed-order perturbation theory. For massless vector operators they are located at N = 1 [88], for massless quark operators at N = 0 [89,90] and for massless scalar operators at N = −1 [91]. The leading term of a The leading behavior in x-space is ∝ α 3 s ln 4 (1/x), but the coefficients of the sub-leading terms have an oscillating sign while their magnitude increases with increasing logarithmic order, which strongly compensates the leading term in the physical region which is relevant, for example, at the EIC [92]. This is in line with earlier observations in other cases, cf. Refs. [90,91,93,94]. For the complete OME A

The Two-Mass Correction
The 2-loop two-mass OME has been calculated in [6]. The 3-loop two-mass OME A (3),two−mass gq is calculated as follows. Using the projector for the external quark lines given in [35] the following representation is obtained where the function I gq,Q (N) is given by The N-dependence completely factorizes from the dependence of the masses. The integral can be performed with analytic Mellin-Barnes integral techniques. We define the mass ratio by and For the unrenormalized OME one obtainŝ with S ε = exp[ ε 2 (γ E − ln(4π))] and the polynomials Since the color-factor contribution of O(C F T 2 F ) does not receive a finite renormalization, it is directly given in the M scheme [36,95]. We have checked the results using q2e/exp [96] by calculating the moments N = 3, 5 expanding in the first powers of η, cf. also [4]. Note that by expanding the OME in powers of η the root-structures in the above expressions disappear showing a dependence on η only.
We can recover the O(T 2 F ) part of the single mass OMEÂ gq,Q by performing the limit η → 1.
which agrees with the result given before. The renormalized two-mass OME in N and x-space is given in the Appendix.

Conclusions
In this paper we have calculated the single and two-mass contributions to the massive operator matrix element A gq (N), which contributes to the matching relations of the VFNS at 3-loop order. On the technical side of the calculation, we have used the the integration-by-parts program Reduze2 to reduce the scalar integrals with local operator insertions to a minimal set of master integrals. The master integrals have been computed using different techniques based on generating functions. These techniques allowed us to find difference equations for the Mellin space results, which were subsequently solved with the packages Sigma, EvaluateMultiSums, SumProduction, ρ-Sum and HarmonicSums. As in the unpolarized case, the polarized matrix element A (3) gq (N) can be expressed in terms of harmonic sums up to weight w = 4 in Mellin space and harmonic polylogarithms up to weight w = 5 in Bjorken x-space. For the two-mass relation the dependence on the Mellin variable N and the squared mass ratio η factorize. Note that other massive operator matrix elements [19,20,22,29] also depend on more complicated structures like generalized harmonic sums and finite binomial sums. As in the unpolarized case [21], diagrams of the Benz topology contribute. Additionally, we presented the results for the renormalization of the heavy-quark mass in the on-shell and MS scheme.

A Appendix
The massive OME A gq,Q (N) with the strong coupling constant renormalized in the MS scheme and the heavy-quark mass in the on-shell scheme obeys the following expansion up to 3-loop order: gq,Q (N) is given in Ref. [6]. The 3-loop OME reads where we define The polynomials Q i read The methods to obtain the analytic continuation of harmonic sums to complex values of N are presented in Refs. [86,[97][98][99].
For the Bjorken x-space representation it is convenient to define The massive operator matrix element A gq,Q (x) in Bjorken x-space is given in terms of harmonic polylogarithms [83]. To shorten the expressions we use H a (x) ≡ H a as a shorthand notation. The full expression is given by:  For the polarized massive operator matrix element the same set of harmonic polylogarithms contribute as in the unpolarized case. The full set is reads: Methods and programs for the numerical evaluation of harmonic polylogarithms are given in [100]. At 2-loop order the OME is not altered when changing the renormalization scheme of the heavy quark. The additional terms when changing from the on-shell scheme to the MS scheme for m =m are given by: Here we identified the masses in the on-shell and MS scheme m =m to shorten the expressions. It is straight forward to obtain the relation between the two renormalization schemes [45] while keeping also the scale dependence. The corresponding relation has been given e.g. in [101].