Flavor decomposition for the proton helicity parton distribution functions
C. Alexandrou, M. Constantinou, K. Hadjiyiannakou, K. Jansen, F. Manigrasso
FFlavor decomposition for the proton helicity parton distribution functions
Constantia Alexandrou,
1, 2
Martha Constantinou, KyriakosHadjiyiannakou,
1, 2
Karl Jansen, and Floriano Manigrasso
1, 5, 6 Department of Physics, University of Cyprus, P.O. Box 20537, 1678 Nicosia, Cyprus Computation-based Science and Technology Research Center,The Cyprus Institute, 20 Kavafi Street, Nicosia 2121, Cyprus Department of Physics, Temple University, 1925 N. 12th Street, Philadelphia, PA 19122-1801, USA NIC, DESY, Platanenallee 6, D-15738 Zeuthen, Germany Institut fr Physik, Humboldt-Universitt zu Berlin, Newtonstr. 15, 12489 Berlin, Germany Dipartimento di Fisica, Universit di Roma “Tor Vergata”,Via della Ricerca Scientifica 1, 00133 Rome, Italy
We present, for the first time, an ab initio calculation of the individual up, down and strangequark helicity parton distribution functions for the proton. The calculation is performed within thetwisted mass clover-improved fermion formulation of lattice QCD. The analysis is performed usingone ensemble of dynamical up, down, strange and charm quarks with a pion mass of 260 MeV. Thelattice matrix elements are non-perturbatively renormalized and the final results are presented in theMS scheme at a scale of 2 GeV. We give results on the ∆ u + ( x ) and ∆ d + ( x ), including disconnectedquark loop contributions, as well as on the ∆ s + ( x ). For the latter we achieve unprecedented precisioncompared to the phenomenological estimates. Introduction.
The theory of the strong interaction, Quan-tum Chromodynamics (QCD), successfully explains thestructure of hadrons and their interactions. The funda-mental degrees of freedom in QCD are the quarks andgluons, collectively referred to as partons. Partons areresponsible for the rich internal structure of hadrons.Most of the knowledge of the complex hadronic struc-ture comes from parton distribution functions (PDFs), aset of number densities describing the non-perturbativeQCD dynamics. Distribution functions are universalquantities and, therefore, can be accessed by a vari-ety of high-energy scattering processes. The cross sec-tion of such processes can be factorized into a pertur-bative component calculable in perturbative QCD, anda non-perturbative part expressed in terms of the par-tonic densities. The generalized parton distributions(GPDs) and transverse momentum distribution functions(TMDs) complement PDFs, and are necessary for the 3-dimensional mapping of the hadrons.At leading order within the parton model, the PDFshave a simple interpretation. The unpolarized PDFs areinterpreted as the probability to find an unpolarized par-ton with a longitudinal momentum fraction x within an unpolarized nucleon. The helicity PDFs can be inter-preted as the difference between finding quarks with spinsaligned and opposite to that of a longitudinally polar-ized nucleon. The colinear PDFs are completed with thetransversity PDFs, which have the interpretation of find-ing quarks polarized in the same or opposite directionas a transversely polarized nucleon. PDFs play a cen-tral role in the on-going experimental program of majorfacilities, such as, BNL, CERN, DESY, Fermilab, JLaband SLAC (see, e.g., Refs. [1–3]). These experiments pro-vide a wealth of measurements that are jointly analyzedwithin the framework of phenomenological fits. Based onthe available experimental data, the most well-studiedcolinear distributions are the unpolarized, followed bythe helicity with an order of magnitude less experimen-tal data sets, namely a few hundred data sets [4, 5]. Thetransversity PDFs are even less-known [6]. The acces-sible kinematical region is more limited for the helicityand transversity PDFs as compared to the unpolarizedPDFs, and therefore, the reconstruction of PDFs uses in-put from models. Such input introduces dependence onthe functional forms employed. As a consequence, theextraction of the helicity and transversity PDFs are, tosome extent, driven by the fit functions, due to the lack a r X i v : . [ h e p - l a t ] S e p of experimental data (see, e.g., Ref. [6] for a discussion).The dependence on the analysis procedure is evidence bythe tension among some of the global analyses [5, 7–9].The focus of this work are the helicity PDFs, which aretypically accessed experimentally in deep-inelastic scat-tering (DIS), semi-inclusive DIS, Drell-Yan, and proton-proton scattering processes. Currently, the global anal-yses use next-to-leading order (NLO) corrections in per-turbative QCD (NNPDF POL u ( x ) , ∆ d ( x ) are better constrained in the valence sec-tor, with ∆ u ( x ) being more precise. On the other hand,constraining ∆ s ( x ) is not successful, as the kinematic re-gions of some of the data sets (e.g., the W -boson produc-tion data) are not sensitive to the strangeness [5]. Thesituation somewhat improves with the inclusion of kaonproduction SIDIS data, but it is still unsatisfactory, andinfluenced by theoretical assumptions. In the recent workof the JAM Collaboration [7] the authors used inclusiveand semi-inclusive data, and find, for both sets of data,the strange polarization to be very small and consistentwith zero. More details on the global analyses can befound in the recent reports [6, 10].Based on the current status of phenomenological anal-yses, an extraction of the PDFs from first principles ishighly desirable. Here we present the first extraction ofthe up, down and strange helicity PDFs for the protonusing lattice QCD, the only known ab initio formula-tion of QCD. We study both the valence and sea quarkcontributions that allow one to perform a controlled de-composition of the various distributions. To obtain the x -dependence of PDFs, we implement the quasi-PDFmethod [11]. This approach is based on correlation func-tions that are calculable on a Euclidean lattice. Thematrix elements are between include hadron state witha finite momentum (cid:126)P = (0 , , P ). A non-local opera-tor with fermion fields separated by a distance z con-nected by a Wilson line, is inserted between the protonstates. Note that the Wilson line is in the same spatialdirection as (cid:126)P . Naturally, the matrix elements are de-fined in coordinate space, with z varying from zero up tohalf the spatial extent of the lattice. To extract physi-cal quantities, a Fourier transform is applied on the ma-trix elements to obtain the so-called quasi-PDFs, whichare defined in momentum space, x . For large values of P , the momentum boost in the quasi-PDFs can be in-terpreted as a Lorentz boost, recovering the light-conePDF. The difference between quasi-PDFs and light-conePDFs is O (cid:0) Λ /P , m N /P (cid:1) and is calculable in con-tinuum perturbation theory within the Large MomentumEffective Theory (LaMET) [12]. A successful researchprogram on obtaining the PDFs using the quasi-PDFsmethod was developed since Ji’s proposal, leading to the-oretical and numerical advances [13–50]. Recently, an ex-ploratory study appeared on the strange and charm un-polarized PDFs [51] using ensembles with pion mass 310and 690 MeV. However, the work only presents matrix elements in coordinate space. Other methods on extract-ing the x -dependence of distribution functions have beendiscussed [52–80]. For an extensive review of the lat-tice calculations using the quasi-PDFs method, as wellas other approaches to extract PDFs, see Refs. [81, 82]. Lattice implementation.
Based on the quasi-PDFs ap-proach, the light-cone PDFs are obtained by the con-volution of quasi-PDFs and the corresponding analyticexpression for the matching kernel calculated in contin-uum perturbation theory. The quasi-PDFs are definedin momentum space as (cid:101) ∆ q ( x, µ, P ) = 2 P (cid:90) + ∞−∞ dz π e − ixP z M R ( z, P ) , (1)and are Fourier transform of hadronic matrix elements M R ( z, P , µ ) ≡ Z ( z, µ ) M ( z, P ) , (2) M ( z, P ) ≡ (cid:104) N ( P ) | ¯ ψ ( z ) γ γ W (0 , z ) ψ (0) | N ( P ) (cid:105) . (3)These matrix elements are calculable on a Euclideanlattice. In Eq. (3), the proton initial and final states, | N ( P ) (cid:105) , carry the same momentum P = ( P , , , P ), asthe PDFs are obtained in the forward kinematic limit.Here we focus on the helicity PDFs, ∆ q ≡ g q ( x ), andtherefore we use the axial non-local operator, which con-tains a Wilson line W (0 , z ) of length z that guaranteesgauge invariance. The bare matrix elements M ( z, P )must be renormalized with an appropriate renormaliza-tion function, Z ( z, µ ), to remove divergences. We calcu-late Z ( z, µ ) using the RI (cid:48) -type prescription proposed inRefs. [16, 17] Z ( z, µ )12 Z − q ( µ ) Tr (cid:104) V γ γ ( p, z ) (cid:0) V Born g T ( p, z ) (cid:1) − (cid:105)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) p = µ = 1 , (4)which is applied at each value of z separately. We re-fer the Reader to Ref. [34] for notation. Due to thepresence of the Wilson line in the operator, extractingthe singlet renormalization functions is very challenging,as it involves a disconnected diagram. Here we use thenon-singlet function indicated by Z ( z, µ ). We note thatthe difference between the singlet and non-singlet renor-malization functions is expected to be small, as is thecase of the local axial-vector operator [83]. This smalldifference has its origin to the fact that the differencebetween singlet and non-singlet arises to two loops inperturbation theory [84]. In addition to the logarithmicdivergences and finite renormalization, the definition ofEq. (4) also removes the power-law divergence of the Wil-son line. Z ( z, µ ) is obtained at an RI (cid:48) scale µ . In ouranalysis, we convert to the MS scheme at a scale µ = 2GeV. An additional conversion factor is used to bring Z ( z, µ ) in the modified MS-scheme [34]. Therefore, thescale dependence appears in the renormalized matrix el-ement M R ( z, P , µ ). While the matrix elements of localoperators mix under renormalization [85], the non-localoperators under study do not mix in the renormalizationprocess, as discussed in Refs. [32, 33, 35]. This is becausethere is no additional non-local ultraviolet divergence inthe quasi-PDF, an argument that holds to all orders inperturbation theory. However, the mixing occurs at thematching level and should be eliminated. To disentan-gle the singlet PDFs requires the matrix elements of thegluon PDFs, which is beyond the scope of this work. Thenature of the mixing was also discussed earlier in Ref. [21]using the auxiliary field approach.The most widely-used method to obtain the quasi-PDFs is via the discretized Fourier transform of Eq. (1).More recently, alternative reconstruction techniques arebeing explored [44, 69, 80, 86, 87]. In this work, we com-pare the standard Fourier transform, with the Bayes-Gauss-Fourier transform [87]. We find agreement be-tween the two approaches, indicating that the behaviorof the lattice results at the large- x region are not dueto the discretization of the Fourier transform. We thuspresent results using the discretized Fourier transform.As can be seen in Eq. (1), the quasi-PDFs depend onthe nucleon momentum P , which should be finite butlarge. This dependence is expected to be removed by thematching kernel∆ q ( x, µ ) = (cid:90) ∞−∞ dξ | ξ | C (cid:18) ξ, µxP (cid:19) (cid:101) ∆ q (cid:18) xξ , µ, P (cid:19) , (5)which is calculated to a given order in continuum per-turbation theory. The matching kernel for the quasi-PDFs approach has been extensively studied (see, e.g.,Refs. [22–24, 45–49]. In this work we employ the one-loop matching kernel in the modified MS-scheme, as de-fined in Ref. [34]. Note that we choose the factorizationscale to be the same as the renormalization scale µ . Thefinal step in extracting the light-cone PDFs is the appli-cation of the nucleon mass corrections, that have beencalculated analytically in Ref. [13]. Numerical Methods.
Obtaining M ( z, P ) is the mostcomputationally demanding part of the calculation, asit contains a non-local operator, and must be calculatedin the boosted frame. We perform the calculation includ-ing, for the first time, connected and quark-disconnecteddiagrams, for both the light and strange quark. In thelight sector, we extract the isovector and isoscalar com-binations, which are decomposed into the up and downquark helicity PDFs. The calculation is performed usingan ensemble of two light, a strange and a charm quark( N f = 2 + 1 + 1) within the twisted mass fermion formu-lation with clover term. The lattice spacing is a = 0 . ×
64 ( L ≈ m π L ≈ ρ = 0 .
15, on the Wil-son line entering the operator. Both smearing methods contribute to the reduction of the statistical noise. Werefer to Ref. [34] for the details. We use a total number ofmeasurements N meas = 392 , P = 0 . , .
83 and 1 .
24 GeV, respectively. The source-sink separation is t s = 0 .
94 fm for the lowest momentumand t s = 1 .
13 fm for the other two.The evaluation of the quark-disconnected diagrams in-volves the computation of disconnected quark loops thathave to be combined with the nucleon two-point corre-lators. The disconnected quark loop with Wilson linereads L ( t ins , z ) = (cid:88) (cid:126)x ins Tr (cid:2) D − q ( x ins ; x ins + z ) γ γ W ( x ins , x ins + z ) (cid:3) , (6)where D − q ( x ins ; x ins + z ) is the quark propagator, whoseendpoints are connected by a Wilson line of length z .To reduce the stochastic noise coming from the low-modes [89], we computed the first N ev = 200 eigen-pairsof the squared Dirac twisted-mass operator. From theeigen-pairs, the low-modes contribution to the all-to-allpropagator can be exactly reconstructed and the high-modes contribution can then be evaluated with stochas-tic techniques. In particular, the stochastic evaluationof the disconnected loops is based on well-establishedtechniques developed for local operators, such as hier-archical probing [90]. The latter allows for reductionof the contamination of the off-diagonal terms in theevaluation of the trace of Eq. (6), up to a distance 2 k .This is done using Hadamard vectors as basis vectorsfor the partitioning of the lattice. Here, the hierarchi-cal probing algorithm has been implemented with k = 8in 4-dimensions, leading to 512 Hadamard vectors. Inaddition, for the stochastic evaluation of the discon-nected loops in this work we make use of the one-endtrick [91, 92] and fully dilute spin and color sub-spaces.We have recently employed successfully such methods inother studies of disconnected contributions. For moredetails see Refs. [83, 93–95]. Results for the connected and disconnected contributions.
For each value of the proton momentum, P = 0 . , . .
24 GeV, we compute the two-point correlator for200 source positions to reach a good statistical accuracy.We also take all spatial orientations of P and W , thatis, ± x, ± y, ± z . Moreover, both in the two-point anddisconnected three-point functions we average over theforward and backward contributions. The total num-ber of configurations analyzed is 330 for the two small-est momenta, and 480 for the largest one, bringing thetotal statistics to 66,000 and 96,000, respectively. Mo-mentum smearing is applied for the two largest values P = 0 . , .
24 GeV. The gauge links in the Wilson lineentering the disconnected loop of Eq. (6) undergo 10 it-erations of stout smearing, with parameter ρ = 0 . t s , we compute the disconnected three-point correlatorsat t s = 0 . , . , . , . , .
13 fm, and perform a two-state fit analysis, following the procedure described inRef. [93]. We find that the two-state fit gives results thatare in agreement with those obtained form the plateaumethod analysis for t s = 1 .
13 fm. We will use the re-sults from the plateau method for what follows. In Fig. 1 . . . M u + d c o nn ( z , P ) − . − . . M u + d d i s c ( z , P ) P = 0 .
41 GeV P = 0 .
83 GeV P = 1 .
24 GeV0 3 6 9 12 15 z/a − . − . − . . M s ( z , P ) FIG. 1: Real part of the bare matrix elements of the isoscalar u + d connected (upper panel) and disconnected (middlepanel) contributions and the strange quark (lower panel). Thedata for the disconnected matrix elements are obtained withthe plateau fit performed at t s = 1 .
13 fm. we show the real part of the bare matrix elements us-ing t s = 1 .
13 fm for the disconnected contributions. Thedisconnected part of the isoscalar combination u + d , issmaller than the connected isoscalar contribution as ex-pected. The real strange matrix element is about halfas compared to the disconnected u + d . However, inboth u + d and strange we clearly obtain a non-zero sig-nal with the statistical uncertainties under control. Theimaginary part of the bare disconnected matrix elementis compatible with zero at each z , and is not shown. Thematrix elements smoothly decay to zero and for z/a > z = 8 a in the disconnected part of thematrix element, i due to using hierarchical probing withlength 2 k and k = 3. This is verified by repeating theevaluation of the disconnected diagrams with k = 2, andconfirm that the same behavior occurs at z = 4 a andits multiples, reflecting the limitation of the hierarchi-cal probing technique when dealing with large lengths ofthe Wilson line. In taking the Fourier transform in Eq.(1), we choose the cutoff z max such that the renormal-ized matrix element is compatible with zero. Since forthe isoscalar and isovector matrix elements this occursat different values of the Wilson line length z , we use different cutoffs z max for the two quantities. In partic-ular, for the isoscalar case (the sum of connected anddisconnected contributions) at P = 1 .
24 GeV, we use z max = 7 a , which is below the hierarchical probing lengthof 8 a . While, for the isovector case, the matrix elementis compatible with zero at z max /a = 9.Two additional important issues need to be addressedin order to extract the PDFs, namely the dependenceof the results on the momentum boost and the accuracyof the discrete Fourier transform. We examine these is-sues b considering the x ∆ d + ( x ) ≡ x (cid:0) ∆ d + ∆ ¯ d (cid:1) distribu-tion, since the behavior is similar for the other two. Toextract the ∆ d + ( x ) distribution we apply renormaliza-tion and matching procedures separately on the isovector,isoscalar and strange quasi-PDFs. As mentioned above,we neglect the mixing with the gluon PDFs at the match-ing level.In Fig. 2 we show the momentum dependence of x ∆ d + ( x ). We observe that, while when increasing themomentum from 0 .
41 GeV to 0 .
83 GeV there is a dis-crepancy in particular for large values of x , when wefurther increase the momentum to P = 1 .
24 GeV, theresults become compatible. This suggests that conver-gence has been reached within the limits of our currentprecision. In Fig. 2 we also show the dependence of the x ∆ d + ( x ) distribution on the cutoff z max adopted in thecomputation of the isoscalar and isovector quasi-PDFs.Despite the fact that when increasing z max , the result-ing distribution tends to show more pronounced oscil-lations, the results for differrent z max all agree withinuncertainties. In order to estimate the extent of the sys-tematic effect due to the discretization and truncation ofthe Fourier Transform (FT), we employ the Bayes-Gauss-Fourier Transform (BGFT) [87]. As can be seen in Fig.2, the distribution obtained with the BGFT techniqueis compatible with the standard reconstruction based onthe discrete FT. Flavor decomposition and comparison with phenomenol-ogy.
The aim of this work is to obtain the flavor decom-position of the up, down and strange quark distributions,by combining the total isoscalar and isovector contribu-tions at each P value. In Fig. 3 we show our final re-sults at P = 1 .
24 GeV for | x | ∆ q + ( x ) ≡ | x | (∆ q + ∆¯ q ),for q = u, d, s , and compare with the JAM17 data [7].We find that x ∆ d + ( x ) and x ∆ s + ( x ) nicely decay to zeroat x = 1. While x ∆ u + ( x ) is also zero at x = 1 and inagreement with the JAM17 results for x (cid:46) .
6, it decaysslower than the JAM17-determined distribution. On theother hand, we find a remarkable agreement for x ∆ d + ( x )for the whole x region. The lattice determination of thestrange distribution x ∆ s + ( x ) is much more precise ascompared to the one determined from the global anal-ysis. Although small is non-zero for small values of x .This result provides a valuable input for phenomenolog-ical studies. − . − . . . x ∆ d + ( x ) z max /a = { , } z max /a = { , } z max /a = { , }− . − . . . x ∆ d + ( x ) BGFTDiscrete FT0 . . . . . . x − . . . x ∆ d + ( x ) P = 0 .
41 GeV P = 0 .
83 GeV P = 1 .
24 GeV
FIG. 2: Dependence on the cutoff z max of the ∆ d + distri-bution (upper panel). The first(second) number reported be-tween curly brackets indicates the value of z max adopted withthe isoscalar(isovector) matrix element. Comparison betweenthe ∆ d + distribution obtained with discrete Fourier Trans-form and with the BGFT technique [87] (middle panel). Mo-mentum dependence of the distribution ∆ d + (bottom panel). − . . . . . x ∆ u + ( x ) JAM17 P = 1 .
24 GeV − . − . . . x ∆ d + ( x ) . . . . . . x − . . . . x ∆ s + ( x ) FIG. 3: Comparison of lattice data on the up (upper), down(middle), and strange (bottom) quark helicity PDFs (blue) inthe MS scheme at 2 GeV with the JAM17 phenomenologicaldata [7] (gray).
Conclusions.
Results for the up, down and strange quarkhelicity PDFs of the proton, within lattice QCD are pre-sented for the first time using the quasi-PDFs approach.We compute matrix elements with nucleon states boostedto maximum momentum P = 1 .
24 GeV. We verify thatthe ground state matrix elements are well-determined byusing one- and two-state fits, confirming that t s = 1 . | x | ∆ q + areshown in Fig. 3, and are compared with the global fitsof the JAM Collaboration. We find a remarkable agree-ment for the case of ∆ d + for all values of x and for case of∆ u + for x < .
6. We also obtain ∆ s + much more precisethat the phenomenological determination and show thatis clearly non-zero for small values of x . This work pavesthe way for a determination of these helicity PDFs usingensembles simulated with pion mass, which we plan todo in the near future.In the near future, a number of sources of systematicuncertainties will be explored, using the particular en-semble, along the lines of the analysis of Ref. [34]. Othereffects is the implementation of the mixing matching ma-trix between quark and gluon PDFs, that requires knowl-edge of the gluon matrix elements of non-local operators.Systematic uncertainties requiring more than one ensem-ble include discretization effects, volume effects, and pionmass dependence. We plan to assess a proper determi-nation of all sources of systematic uncertainties for theindividual flavor PDFs in the future. Once systematicuncertainties are addressed and quantified, lattice resultscan provide useful input in the global fits for the strangePDFs, as well as the individual light-quark PDFs. Thiscalculation is a first step towards achieving this goal. Acknowledgments
We would like to thank all members of ETMC for theirconstant and pleasant collaboration. We also thank N.Sato for providing the global fits data and J. Green forhis comments. Finally, our thanks go to A. Scapellato forproviding us the data at 0.43 GeV for the connected ma-trix elements. M.C. acknowledges financial support bythe U.S. Department of Energy Early Career Award un-der Grant No. DE-SC0020405. K.H. is supported by theCyprus Research and Innovation Foundation under grantPOST-DOC/0718/0100. F.M. is supported by the Euro-pean Joint Doctorate program STIMULATE of the Eu-ropean Unions Horizon 2020 research and innovation pro-gramme under grant agreement No 765048. This researchincludes calculations carried out on the HPC resourcesof Temple University, supported in part by the NationalScience Foundation through major research instrumenta-tion grant number 1625061 and by the US Army ResearchLaboratory under contract number W911NF-16-2-0189.Computations for this work were carried out in part onfacilities of the USQCD Collaboration, which are funded by the Office of Science of the U.S. Department of Energy.This research used resources of the Oak Ridge Leader-ship Computing Facility, which is a DOE Office of Sci-ence User Facility supported under Contract DE-AC05-00OR22725. The gauge configurations have been gener-ated by the Extended Twisted Mass Collaboration on theKNL (A2) Partition of Marconi at CINECA, through thePrace project Pra13 3304 ”SIMPHYS”. [1] J. Gao, L. Harland-Lang, and J. Rojo, Phys. Rept. ,1 (2018), 1709.04922.[2] C. A. Aidala, S. D. Bass, D. Hasch, and G. K. Mallot,Rev. Mod. Phys. , 655 (2013), 1209.2803.[3] K. J. Eskola, P. Paakkinen, H. Paukkunen, and C. A.Salgado, Eur. Phys. J. C , 163 (2017), 1612.05741.[4] J. J. Ethier and E. R. Nocera, Ann. Rev. Nucl. Part. Sci.pp. 1–34 (2020), 2001.07722.[5] E. R. Nocera, R. D. Ball, S. Forte, G. Ridolfi, and J. Rojo(NNPDF), Nucl. Phys. B887 , 276 (2014), 1406.5539.[6] M. Constantinou et al. (2020), 2006.08636.[7] J. J. Ethier, N. Sato, and W. Melnitchouk, Phys. Rev.Lett. , 132001 (2017), 1705.05889.[8] D. De Florian, G. A. Lucero, R. Sassot, M. Stratmann,and W. Vogelsang, Phys. Rev.
D100 , 114027 (2019),1902.10548.[9] E. R. Nocera, in (2017), 1702.05077.[10] H.-W. Lin et al., Prog. Part. Nucl. Phys. , 107 (2018),1711.07916.[11] X. Ji, Phys. Rev. Lett. , 262002 (2013), 1305.1539.[12] X. Ji, Sci. China Phys. Mech. Astron. , 1407 (2014),1404.6680.[13] J.-W. Chen, S. D. Cohen, X. Ji, H.-W. Lin, and J.-H.Zhang, Nucl. Phys. B911 , 246 (2016), 1603.06664.[14] C. Alexandrou, K. Cichy, M. Constantinou, K. Had-jiyiannakou, K. Jansen, F. Steffens, and C. Wiese, Phys.Rev.
D96 , 014513 (2017), 1610.03689.[15] R. A. Brice˜no, M. T. Hansen, and C. J. Monahan, Phys.Rev.
D96 , 014502 (2017), 1703.06072.[16] M. Constantinou and H. Panagopoulos, Phys. Rev.
D96 ,054506 (2017), 1705.11193.[17] C. Alexandrou, K. Cichy, M. Constantinou, K. Had-jiyiannakou, K. Jansen, H. Panagopoulos, and F. Stef-fens, Nucl. Phys.
B923 , 394 (2017), 1706.00265.[18] X. Ji, J.-H. Zhang, and Y. Zhao, Nucl. Phys. B , 366(2017), 1706.07416.[19] X. Ji, J.-H. Zhang, and Y. Zhao, Phys. Rev. Lett. ,112001 (2018), 1706.08962.[20] T. Ishikawa, Y.-Q. Ma, J.-W. Qiu, and S. Yoshida, Phys.Rev.
D96 , 094019 (2017), 1707.03107.[21] J. Green, K. Jansen, and F. Steffens, Phys. Rev. Lett. , 022004 (2018), 1707.07152.[22] W. Wang, S. Zhao, and R. Zhu, Eur. Phys. J.
C78 , 147(2018), 1708.02458.[23] I. W. Stewart and Y. Zhao, Phys. Rev.
D97 , 054512(2018), 1709.04933.[24] T. Izubuchi, X. Ji, L. Jin, I. W. Stewart, and Y. Zhao,Phys. Rev.
D98 , 056004 (2018), 1801.03917. [25] C. Alexandrou, K. Cichy, M. Constantinou, K. Jansen,A. Scapellato, and F. Steffens, Phys. Rev. Lett. ,112001 (2018), 1803.02685.[26] J.-H. Zhang, J.-W. Chen, L. Jin, H.-W. Lin, A. Schfer,and Y. Zhao, Phys. Rev.
D100 , 034505 (2019),1804.01483.[27] R. A. Briceo, J. V. Guerrero, M. T. Hansen, and C. J.Monahan, Phys. Rev. D , 014511 (2018), 1805.01034.[28] G. Spanoudes and H. Panagopoulos, Phys. Rev. D ,014509 (2018), 1805.01164.[29] C. Alexandrou, K. Cichy, M. Constantinou, K. Jansen,A. Scapellato, and F. Steffens, Phys. Rev. D98 , 091503(2018), 1807.00232.[30] Y.-S. Liu et al. (2018), 1807.06566.[31] A. Radyushkin, Phys. Lett. B , 380 (2019),1807.07509.[32] J.-H. Zhang, X. Ji, A. Schfer, W. Wang, and S. Zhao,Phys. Rev. Lett. , 142001 (2019), 1808.10824.[33] Z.-Y. Li, Y.-Q. Ma, and J.-W. Qiu, Phys. Rev. Lett. ,062002 (2019), 1809.01836.[34] C. Alexandrou, K. Cichy, M. Constantinou, K. Had-jiyiannakou, K. Jansen, A. Scapellato, and F. Steffens,Phys. Rev.
D99 , 114504 (2019), 1902.00587.[35] W. Wang, J.-H. Zhang, S. Zhao, and R. Zhu, Phys. Rev.D , 074509 (2019), 1904.00978.[36] J.-W. Chen, H.-W. Lin, and J.-H. Zhang, Nucl. Phys. B , 114940 (2020), 1904.12376.[37] T. Izubuchi, L. Jin, C. Kallidonis, N. Karthik, S. Mukher-jee, P. Petreczky, C. Shugert, and S. Syritsyn, Phys. Rev.
D100 , 034516 (2019), 1905.06349.[38] K. Cichy, L. Del Debbio, and T. Giani, JHEP , 137(2019), 1907.06037.[39] W. Wang, Y.-M. Wang, J. Xu, and S. Zhao, Phys. Rev.D , 011502 (2020), 1908.09933.[40] H.-D. Son, A. Tandogan, and M. V. Polyakov, Phys. Lett.B , 135665 (2020), 1911.01955.[41] J. R. Green, K. Jansen, and F. Steffens, Phys. Rev. D , 074509 (2020), 2002.09408.[42] Y. Chai et al., Phys. Rev. D , 014508 (2020),2002.12044.[43] V. Braun, K. Chetyrkin, and B. Kniehl, JHEP , 161(2020), 2004.01043.[44] S. Bhattacharya, K. Cichy, M. Constantinou, A. Metz,A. Scapellato, and F. Steffens (2020), 2004.04130.[45] S. Bhattacharya, K. Cichy, M. Constantinou, A. Metz,A. Scapellato, and F. Steffens, Phys. Rev. D , 034005(2020), 2005.10939.[46] S. Bhattacharya, K. Cichy, M. Constantinou, A. Metz,A. Scapellato, and F. Steffens (2020), 2006.12347.[47] L.-B. Chen, W. Wang, and R. Zhu (2020), 2005.13757. [48] L.-B. Chen, W. Wang, and R. Zhu (2020), 2006.10917.[49] L.-B. Chen, W. Wang, and R. Zhu (2020), 2006.14825.[50] X. Ji, Y. Liu, A. Schfer, W. Wang, Y.-B. Yang, J.-H.Zhang, and Y. Zhao (2020), 2008.03886.[51] R. Zhang, H.-W. Lin, and B. Yoon (2020), 2005.01124.[52] K.-F. Liu and S.-J. Dong, Phys. Rev. Lett. , 1790(1994), hep-ph/9306299.[53] W. Detmold and C. J. D. Lin, Phys. Rev. D73 , 014501(2006), hep-lat/0507007.[54] V. Braun and D. Mueller, Eur. Phys. J.
C55 , 349 (2008),0709.1348.[55] G. S. Bali et al., Eur. Phys. J. C , 217 (2018),1709.04325.[56] G. S. Bali, V. M. Braun, B. Glle, M. Gckeler, M. Gru-ber, F. Hutzler, P. Korcyl, A. Schfer, P. Wein, and J.-H.Zhang, Phys. Rev. D98 , 094507 (2018), 1807.06671.[57] W. Detmold, I. Kanamori, C. J. D. Lin, S. Mondal, andY. Zhao (2018), 1810.12194.[58] J. Liang, T. Draper, K.-F. Liu, A. Rothkopf, and Y.-B.Yang (XQCD) (2019), 1906.05312.[59] Y.-Q. Ma and J.-W. Qiu, Phys. Rev.
D98 , 074021 (2018),1404.6860.[60] Y.-Q. Ma and J.-W. Qiu, Int. J. Mod. Phys. Conf. Ser. , 1560041 (2015), 1412.2688.[61] A. Radyushkin, Phys. Lett. B767 , 314 (2017),1612.05170.[62] A. J. Chambers, R. Horsley, Y. Nakamura, H. Perlt,P. E. L. Rakow, G. Schierholz, A. Schiller, K. Somfleth,R. D. Young, and J. M. Zanotti, Phys. Rev. Lett. ,242001 (2017), 1703.01153.[63] A. V. Radyushkin, Phys. Rev.
D96 , 034025 (2017),1705.01488.[64] K. Orginos, A. Radyushkin, J. Karpie, and S. Zafeiropou-los, Phys. Rev.
D96 , 094503 (2017), 1706.05373.[65] Y.-Q. Ma and J.-W. Qiu, Phys. Rev. Lett. , 022003(2018), 1709.03018.[66] A. V. Radyushkin, Phys. Lett.
B781 , 433 (2018),1710.08813.[67] A. Radyushkin, Phys. Rev.
D98 , 014019 (2018),1801.02427.[68] J.-H. Zhang, J.-W. Chen, and C. Monahan, Phys. Rev.D , 074508 (2018), 1801.03023.[69] J. Karpie, K. Orginos, and S. Zafeiropoulos, JHEP ,178 (2018), 1807.10933.[70] R. S. Sufian, J. Karpie, C. Egerer, K. Orginos, J.-W. Qiu,and D. G. Richards, Phys. Rev. D99 , 074507 (2019),1901.03921.[71] B. Jo, J. Karpie, K. Orginos, A. Radyushkin,D. Richards, and S. Zafeiropoulos, JHEP , 081 (2019),1908.09771.[72] A. V. Radyushkin, Phys. Rev. D , 116011 (2019),1909.08474.[73] B. Jo, J. Karpie, K. Orginos, A. V. Radyushkin, D. G.Richards, R. S. Sufian, and S. Zafeiropoulos, Phys. Rev. D100 , 114512 (2019), 1909.08517.[74] I. Balitsky, W. Morris, and A. Radyushkin, Phys. Lett.B , 135621 (2020), 1910.13963.[75] A. Radyushkin, Int. J. Mod. Phys. A , 2030002 (2020), 1912.04244.[76] R. S. Sufian, C. Egerer, J. Karpie, R. G. Edwards, B. Jo,Y.-Q. Ma, K. Orginos, J.-W. Qiu, and D. G. Richards(2020), 2001.04960.[77] B. Jo, J. Karpie, K. Orginos, A. V. Radyushkin, D. G.Richards, and S. Zafeiropoulos (2020), 2004.01687.[78] M. Bhat, K. Cichy, M. Constantinou, and A. Scapellato(2020), 2005.02102.[79] K. Can et al. (2020), 2007.01523.[80] C. Alexandrou, K. Cichy, M. Constantinou, K. Had-jiyiannakou, K. Jansen, A. Scapellato, and F. Steffens(2020), 2008.10573.[81] K. Cichy and M. Constantinou, Adv. High Energy Phys. , 3036904 (2019), 1811.07248.[82] X. Ji, Y.-S. Liu, Y. Liu, J.-H. Zhang, and Y. Zhao (2020),2004.03543.[83] C. Alexandrou, S. Bacchio, M. Constantinou, J. Finken-rath, K. Hadjiyiannakou, K. Jansen, G. Koutsou, andA. Vaquero Aviles-Casco (2019), 1909.00485.[84] M. Constantinou, M. Hadjiantonis, H. Panagopoulos,and G. Spanoudes, Phys. Rev. D , 114513 (2016),1610.06744.[85] J. C. Collins, Renormalization: An Introduction toRenormalization, the Renormalization Group and theOperator-Product Expansion , Cambridge Monographson Mathematical Physics (Cambridge University Press,1984).[86] J. Karpie, K. Orginos, A. Rothkopf, and S. Zafeiropoulos,JHEP , 057 (2019), 1901.05408.[87] C. Alexandrou, G. Iannelli, K. Jansen, and F. Manigrasso(Extended Twisted Mass) (2020), 2007.13800.[88] G. S. Bali, B. Lang, B. U. Musch, and A. Schfer, Phys.Rev. D93 , 094515 (2016), 1602.05525.[89] A. Abdel-Rehim, C. Alexandrou, M. Constantinou,J. Finkenrath, K. Hadjiyiannakou, K. Jansen, C. Kalli-donis, G. Koutsou, A. V. Avils-Casco, and J. Volmer,PoS