Matching Quasi Generalized Parton Distributions in the RI/MOM scheme

Within the framework of large momentum effective theory (LaMET), genenaralized parton distributions (GPDs) can be extracted from lattice calculations of quasi-GPDs through a perturbative matching relation, up to power corrections that are suppressed by the hadron momentum. In this paper, we focus on isovector quark GPDs, including the unpolarized, longitudinally and transversely polarized cases, and present the one-loop matching that connects the quasi-GPDs renormalized in a regularization-independent momentum subtraction (RI/MOM) scheme to the GPDs in MS scheme. We find that the matching coefficient is independent of the momentum transfer squared. As a consequence, the matching for the quasi-GPD with zero skewness is the same as that for the quasi-PDF. Our results provide a crucial input for the determination of quark GPDs from lattice QCD using LaMET.


I. INTRODUCTION
Understanding the internal structure of nucleons has been an important goal of hadron physics. For many decades, our knowledge on the structure of nucleons has been mainly relying on experimental measurements of their form factors (FFs) and parton distribution functions (PDFs). The FFs describe the spatial distribution of charge and current within the nucleon and can be probed in elastic lepton-nucleon scattering, while the PDFs characterize the longitudinal momentum distribution of quarks and gluons in the nucleon and can be measured in deep-inelastic scattering processes.
The proposal of generalized parton distributions (GPDs) (for a review, see e.g. [1][2][3]) provides a novel opportunity to characterize the partonic structure of nucleons. As a generalization of the PDFs to off-forward kinematics, the GPDs contain a wealth of new information on nucleon structure. They naturally encompass the FFs, PDFs as well as the distribution amplitudes (DAs), and offer a description of the correlations between the transverse position and longitudinal momentum of quarks and gluons inside the nucleon, thereby giving access to quark and gluon orbital angular momentum contributions to the nucleon spin. Experimentally, the GPDs can be accessed through hard exclusive processes like deeply virtual Compton scattering or meson production. Much effort has been devoted to measuring such processes at completed and on-going experiments (HERA [4][5][6][7][8][9], COMPASS [10], JLab [11][12][13][14]), and will be continued at planned future facilities such as EIC [15,16] and EICC [17]. Given the complicated kinematic dependence of GPDs, extracting them from the accumulated experimental data is in general rather difficult, and one usually needs to resort to certain models that allow for an extrapolation to kinematic regions that are not accessible directly [18].
On the other hand, lattice effort of studying GPDs has been mainly focused on the computation of their moments [19][20][21][22][23][24][25]. The full distribution can be reconstructed in principle if all their moments are known. However, the number of moments that are calculable on lattice is very limited, owing to power divergent mixing between different moments operators and increasing stochastic noise for high moments operators.
In the past few years, a new theoretical framework has been developed to circumvent the above difficulties, which is now known as the large momentum effective theory (LaMET) [26,27]. According to LaMET, the GPDs can be extracted from lattice QCD calculations of appropriately constructed static-operator matrix elements, which are named the quasi-GPDs. The quasi-GPDs are usually hadron-momentum dependent but time independent, and thus can be readily computed on the lattice. After being renormalized non-perturbatively in an appropriate scheme, the renormalized quasi-GPDs can then be matched onto the usual GPDs through a factorization formula accurate up to power corrections that are suppressed by the hadron momentum [28,29].
As for the GPDs, there have been studies on the perturbative matching of the isovector quark quasi-GPDs [28,29]. Such studies were performed in a transverse momentum cutoff scheme and therefore not well-suited for the lattice implementation. In this paper, we reconsider the one-loop matching for isovector quark quasi-GPDs in the RI/MOM scheme. The results can be used to match the quasi-GPDs calculated in lattice QCD and renormalized in the RI/MOM scheme onto the GPDs in MS scheme.
The rest of the paper is organized as follows: In Sec. II, we establish our definitions and conventions. In Sec. III, we present a rigorous derivation of the factorization formula for the isovector quark quasi-GPD based on operator product expansion (OPE). Sec. IV and V are devoted to the RI/MOM renormalization and matching procedure, respectively. We also explain how to obtain the matching coefficients of DAs from the one-loop results of GPDs in Sec. V. Our summary is given in Sec. VI.

II. DEFINITIONS AND CONVENTIONS
The parent function for the quark GPDs, which we call parent-GPD for simplicity, is defined from the Fourier transform of the off-forward matrix element of a light-cone correlator, where x ∈ [−1, 1], the light-cone coordinates ζ ± = (ζ 0 ± ζ 3 )/ √ 2, and the hadron state |P , S (|P , S ) is denoted by its momentum and spin. The parent-GPD is defined in the MS scheme and µ is the renormalization scale. The kinematic variables are defined as where without loss of generality we choose a particular Lorentz frame so that the average momentum and only consider the case with 0 < ξ < 1.
To calculate the quark GPDs within LaMET, we consider a quark quasi-parent-GPD defined from an equal-time correlator: where µ is the renormalization scale in a particular scheme, and N is a normalization factor that depends on the choice of Γ. For example, N = 2P z for Γ = γ z . The nonlocal quark bilinear is along the z direction with a spacelike Wilson line The kinematic variables are similar to those in Eq.
(2) except that the "quasi" skewness parameter which is equal to ξ up to power corrections. From now on we will replace ξ with ξ by assuming that the power corrections are small. The quasi-GPDs are defined as form factors of of the quasi-parent-GPD, where H, E, H , and E are the quasi-GPDs with support x ∈ (−∞, ∞). Again, H and E are non-zero only for transversity quasi-GPD. In order to minimize operator mixing on lattice, we choose Γ = γ t , γ z γ 5 , and iσ z⊥ for the unpolarized, helicity, and transversity quasi-GPDs [36,51], respectively, which all correspond to the same normalization factor N = 2P t . According to LaMET [26,27], the quasi-GPDs and GPDs are related through a factorization formula. For example, where O M 2 /P 2 z , Λ 2 QCD /P 2 z are mass and higher-twist corrections. Since the choice of Γ corresponds to a uniqueΓ, we suppress the labelΓ in the matching coefficient C Γ . Similar factorization formulas also exist for H , E, and E . Eq. (12) with its explicit form will be rigorously derived in the next section.

III. OPERATOR PRODUCT EXPANSION AND THE FACTORIZATION FORMULA
In this section, we derive the explicit form of the factorization formula for the quasi-GPDs using the OPE of the nonlocal quark bilinear O(Γ, z). The same method has been used for the "lattice cross section" [91] and quasi-PDF [38], which are both forward matrix elements of a nonlocal gauge-invariant operator. In the case of non-singlet quasi-PDF, O(γ z , z, µ) (e.g., in the MS scheme) can be expanded in terms of local gauge-invariant operators in the |z| → 0 limit [38], where e µ = (0, 0, 0, 1), µ 0 = z, C n = 1 + O(α s ) is the Wilson coefficient, and O µ0µ1···µn (µ) is the only allowed renormalized traceless symmetric twist-2 quark operators at leading power in the OPE, where Here Z ij n+1 (ε, µ) are multiplicative MS renormalization factors and (µ 0 · · · µ n ) stands for the symmetrization of these Lorentz indices.
The multiplicative renormalization shown in Eq. (14) is valid for the forward case only, as it is known that in the off-forward case, O µ0µ1...µn (µ) can mix with other twist-2 operators with overall derivatives i∂ µ according to the renormalization group equation [92] where the anomalous dimension Γ is an upper triangle matrix. In off-forward matrix elements, the overall derivative i∂ µ contributes a factor of the momentum transfer ∆ µ . As a result, the OPE in Eq. (13) cannot maintain its form under evolution in µ, so one has to choose the operator bases to be the eigenvectors of Eq. (15) so that each of them is multiplicatively renormalizable. At leading logarithmic (LL) accuracy, Eq. (15) is diagonalized by the conformal operators [92,93], where n µ is an arbitrary four-vector, and C 3/2 n (η) is the Gegenbauer polynomial. Beyond LL, the conformal operators start mixing with each other, but Eq. (15) can still be diagonalized with the "renormalization group impoved" conformal operators [92,94], where B nn = 1. As a result, the nonlocal operator O(γ z , z, µ) should be generally expanded in terms of these improved conformal operators as with the normalization factor so that the OPE has the correct forward limit as Eq. (13). For µ 0 = µ 1 = · · · = µ n = +, the off-forward matrix element of the conformal operator O µ0µ1...µn (µ) is given by which is also known as the Gegenbauer moments. Using Lorentz covariance, we have for µ 0 = µ 1 = · · · = µ n = z, where M 2 , t P 2 z , and we have used ∆ · P = 0. The power corrections originate from the subtracted traces in the kinematic part of the matrix element, and their exact form will be derived in the future.
Plugging Eq. (21) into Eqs. (17) and (18), we have the leading-twist approximation of the off-forward matrix element of O(γ z , z, µ), where the higher-twist terms contribute to O(z 2 Λ 2 QCD ). The polynomiality of C 3/2 n allows us to define for m ≤ n, where C m is also a polynomial that satisfies If we define the matching coefficients as then we can Fourier transform Eq. (22) from z to xP z to obtain the quasi-GPD and its factorization formula, where the second form in Eq. (27) was postulated in Ref. [28,29]. Since xP z is the Fourier conjugate to z, the highertwist contribution of O(z 2 Λ 2 QCD ) in Eq. (22) should be of O(Λ 2 QCD /(x 2 P 2 z )) in momentum space with an enhancement at small x 1 . Based on Eqs. (26) and (27), we can infer that the matching coefficients for the quasi-GPDs H, H , E, and E must be the same.
The two matching coefficients above are related to each other by The factorization formulas are similar to the evolution equations for the GPD [95,96]. Notably, at zero skewness ξ = 0, we have where the matching kernel C γ z (x/y, 0, µ/(yP z )) is exactly the same matching coefficient for the MS quasi-PDF [38], even when t = 0. Moreover, in the forward limit ξ → 0 and t → 0, Eq. (29) is exactly the factorization formula for the MS quasi-PDF [38].
On the other hand, in the limit ξ → 1 and t → 0, we obtain the factorization formula for the quasi DA, whose explicit form has been postulated in Refs. [42,44,97]. The same procedure described above also applies to the Γ = γ t case. This finishes our derivation of the factorization formula for the isovector quark quasi-GPD, which will enable us to identify the matching coefficients from the one-loop calculation in Sec. V.

IV. RENORMALIZATION
Following the strategy in Ref. [97], the UV divergence of the quasi-GPD only depends on the operator O(Γ, z), not on the external states. We can choose the same renormalization factor as the one for the quasi-PDF [35,67]. For each value of z, the RI/MOM renormalization factor Z is calculated nonperturbatively on lattice by imposing the condition that the quantum corrections of the correlator in an off-shell quark state vanish at scales In a systematic calculation of GPD, we start with the bare matrix element of the nonlocal quark bilinear on lattice, h(Γ, z, ξ, t, P z , a) = 1 N P , S | O(Γ, z, a)|P , S .
After performing RI/MOM renormalization and taking the continuum limit, the renormalized matrix element is which is to be Fourier transformed into the x-space Next, we need to disentangle the terms with different kinematic dependences to extract quasi-GPDs from F . Finally, we match quasi-GPDs in the RI/MOM scheme to GPDs in MS scheme according to Eq. (12). Note that the continuum limit has been taken after the RI/MOM renormalization, we can therefore calculate the matching coefficient in the continuum as the result is regularization independent. For simplicity, we choose dimensional regularization in our calculation. The one-loop result will be presented in the next section.

V. ONE LOOP MATCHING COEFFICIENT
When the hadron momentum P z is much greater than M and Λ QCD , the RI/MOM quasi-GPD can be matched onto the MS GPD through the factorization formula [35,38], where r = µ 2 R /(p z R ) 2 . Here we have chosen the explicit form of factorization in Eq. (27). To obtain the matching coefficient, we calculate their on-shell massless quark matrix element in perturbation theory by replacing the hadron states in Eqs. (1) and (7) with the quark states carrying momentum p + ∆/2 and p − ∆/2 with p µ = (p t , 0, 0, p z ).
At tree level, the GPDs and quasi-GPDs are At one-loop order, H (1) and H (1) are nonzero and not equal, so their next-to-leading order (NLO) matching kernel is non-trivial; since H (1) = H (1) , a two-loop calculation is needed to extract the NLO matching kernel; E (1) , E (1) , E (1) , and E (1) vanish for massless quarks. For massive quarks, E (1) = E (1) = 0 and E (1) = E (1) = 0 according to Ref. [28,29], so the NLO matching kernel for E (1) and E (1) can only be extracted from the two-loop matrix elements in massive quark states. This can be a cross check of the factorization formulas in Eqs. (26) and (27), which, however, is beyond the scope of this work.
In order to combine the "real" and "virtual" contributions (defined in Ref. [35]) in a compact form at one-loop level, we introduce a plus function defined as with two arbitrary functions h(x) and g(x). By taking the limit p t → p z , we obtain the matching kernel for the gauge-invariant bare quasi-GPD and MS GPD in a quark, where the subscript B denotes "bare"; the infrared (IR) divergence is regulated by t and dimensional regularization (D = 4 − 2ε), and canceled in C The results are where Some technical details of the calculation are provided in the Appendix. The above calculation has been carried out in momentum space. In principle, the same result can be obtained from calculations in coordinate space and then taking a Fourier transform. For examples in the case of meson DA and nucleon PDF see Refs. [73,98]. However, as noticed in [38,73], the step of taking Fourier transform is highly non-trivial. We observe that the bare matching coefficients for Γ = γ t , γ z γ 5 , and iσ z⊥ reduce to that for the quasi-PDFs [38,84] when ξ = 0 even if t = 0. One can also obtain the bare matching coefficients of DAs [97] by crossing the external state with the following replacement ξ → 1/(2y − 1), x/ξ → 2x − 1, and the external momentum p z to p z /2 [28].
Next we need the counterterm of the quasi-GPD in RI/MOM scheme. As we argued in Sec. IV, we can use the renormalization factor for the quasi-PDF to renormalize the quasi-GPD, which leads to the one-loop RI/MOM counterterm [35,67] where , r) is the real part of the off-shell quark matrix element of the quasi-PDF calculated at the subtraction point { µ}, and F (1) is the coversion factor between RI/MOM and MS schemes for the local operator O(Γ, 0). We choose Landau gauge, which is convenient for lattice simulation, and project out the coefficient of Γ (also known as the minimal projection according to [67]) to obtain q (1) and F (1) . The results for different spin structures are [67,84], and Finally, combining Eqs. (41) and (50), and making the replacements x → x/y, ξ → ξ/y, and p z → yP z [35,38], we obtain the one-loop matching coefficient C Γ in Eq. (36),

VI. SUMMARY
Within the framework of LaMET, we have derived the one-loop matching coefficients that connect the isovector quark quasi-GPDs renormalized in the RI/MOM scheme to GPDs in the MS scheme. The calculation was performed for the unpolarized, longitudinally and transversely polarized cases defined with Γ = γ t , γ z γ 5 and iσ z⊥ , respectively. We also presented a rigorous derivation of the factorization formula for isovector quark quasi-GPDs based on OPE. The matching coefficient turns out to be independent of the momentum transfer squared t. As a result, for quasi-GPDs with zero skewness the matching coefficient is the same as that for the quasi-PDF. Our results will be used to extract the isovector quark GPDs from lattice calculations of the corresponding quasi-GPDs.

APPENDIX
In this Appendix, we present some technical details in calculating the following dimensionless integral that arises from the vertex diagram in Fig. 1, where N is a function of Feynman parameters y 1 and y 2 and n is the power of the denominator of the integrand; t = −t/(p z ) 2 . In unphysical regions (x < −ξ and x > 1), the integral has no 1/ε pole so that it can be easily calculated by setting ε = 0. However, this is not the case in the ERBL (|x| < ξ) and DGLAP (ξ < x < 1) regions where there are IR divergences.
As an example, we evaluate the integral with N = 1 and n = 3/2. After integration over y 2 , the remaining integrand denoted as F (y 1 , ε) contains hypergeometric functions 2 F 1 . We identify the divergent part of F (y 1 , ε) as A(ε)/y 1+ε 1 in the limit of ε → 0, and then separate it out from the integral, where the first term is convergent so that we can set ε = 0 before the integration. We suppress x, ξ, and t dependences of F and A for simplicity. In Fig. 2, the singularities are shown in integration regions. We use Pfaff transformation to extract the divergent part We obtain Finally, we have More generally, when calculating the vertex diagram in Fig.(1), we encounter integrals similar to Eq.(57) with numerator of the integrand replaced by polynomials of y 1 and y 2 . After integrating out y 2 , we obtain Appell hypergeometric function F 1 . In this case, to separate the divergent part, we need Euler transformation F 1 (α; β, β ; γ; x, y) = (1 − x) −β (1 − y) −β F 1 (γ − α; β, β ; γ; x x − 1 , y y − 1 ) .
In the following, we list integrals used in our calculation: 1. In unphysical region x < −ξ, there is no divergence.