Quark distribution inside a pion in many-flavor (2 + 1)-dimensional QCD using lattice computations: UV listens to IR
JJLAB-THY-21-3258
Quark distribution inside a pion in many-flavor 2+1 dimensional QCD using lattice:UV listens to IR
Nikhil Karthik
1, 2, ∗ Department of Physics, College of William & Mary, 300 Ukrop Way, Williamsburg, VA 23185, USA Thomas Jefferson National Accelerator Facility, Newport News, VA 23606, USA (Dated: January 8, 2021)We study the changes in the short-distance quark structure of the Nambu-Goldstone boson whenthe long-distance symmetry-breaking scales are depleted controllably. We achieve this by studyingthe valence Parton Distribution Function (PDF) of pion in 2+1 dimensional two-color QCD, withthe number N of massless quarks as the tunable parameter that slides the theory from being stronglybroken for N = 0 to the conformal window for N >
4, where the theory is gapped by the fixedfinite volume. We perform our study non-perturbatively using lattice simulations with N = 0 , , , x -dependent PDF itself dramatically changes frombeing broad in the scale-broken sector to being sharply peaked in the near-conformal region, bestreflected in PDF shape observables such as the cumulants. I. INTRODUCTION
QCD offers a unified theoretical description of themass-gapped hadron spectrum in its infrared limit, andof the asymptotically free quarks and gluons in the short-distance limit. While the perturbative facet of QCD hasbeen subject to stringent tests in collider experiments [1],the only first-principle field theoretic description of thelow-energy hadronic physics with a controlled continuumlimit comes from the numerical lattice QCD computa-tions. Even though lattice QCD reproduces the low-energy behavior of QCD precisely (c.f., [2–4] and refer-ences therein), a good theoretical understanding of QCDby abstracting and characterizing the cause of the com-plex nonperturbative low-energy features to few relevantaspects of the theory is sought after. A fascinating as-pect of QCD is the very fact that there is a non-zeromass-gap in this classically conformal field theory. Un-derstanding how length-scale emerges in QCD in termsof the short-distance dynamics of the partons inside theproton and other hadrons is one way to approach thisproblem, and will be a key question that will be stud-ied in the Electron-Ion Collider [5]. A starting point inthis approach is to break-down the energy-momentumtensor [6, 7] into the quark and gluon momentum frac-tions (c.f., [8, 9] for recent results for such a break-downof proton momentum fraction), and the trace-anomalyparts. This approach thus closely ties the understandingof emergence of scale in the infrared to the parton distri-bution functions (PDF), f ( x ) of hadrons and their mo-ments, with x being the momentum fraction of hadronscarried by a parton.One of the low-energy scale generating mechanism isthe spontaneous chiral symmetry breaking (SSB) that ∗ [email protected] leads to a dimensionful quark condensate, with the pionbeing its Nambu-Goldstone (NG) boson. The naturalquestion then is whether we can learn about the SSBphysics by studying the quark and gluon structure ofpion, which has to be constrained in such a way as tomake it exactly massless in the chiral limit. Therefore,not surprisingly, the PDF of pion has been determinedfrom multiple analyses of experimental data with increas-ing levels of sophisticated analysis techniques, processesbeing included and at different perturbative orders [10–18]. Of special theoretical interest has been the valencequark distribution and its puzzling (1 − x ) β large- x behav-ior — whether the value of β ≈ ≥
2. These aspectshave been extensively studied through many model cal-culations [18–30]. With the rapid progress in the leading-twist perturbative matching formalisms (LaMET [31, 32],SDF [33–35], good lattice cross-section using current-current correlators [36, 37], and see reviews [38–42]), thelattice QCD computations of the pion PDF have beenable to weigh in on the large- x behavior [43–49]. Whilethe lattice findings seem to lean closer to β ≈
1, thestudies [43, 49] found that variations in alternate analy-sis methods could make the results consistent with 2 aswell. Thus, our understanding of the pion PDF is stillevolving, and will be guided further by some of the up-coming experiments [50, 51] as well as the future latticecomputations.The aim of this paper is to make use of the leading-twist formalism and extend it to lattice computations ofPDF in a family of QCD-like theories, with the degreeof infrared scale-breaking varying within the family. Bystudying how and which aspects of the quark structureinside the Nambu-Goldstone boson (which we simply callas the pion) change because of variations in the infraredscale, induced by the choice of members in the family oftheories, we aim to understand the correlations betweenthe quark structure of pion and the long-distance vac-uum structure. By such observations on how the PDFs a r X i v : . [ h e p - l a t ] J a n evolve to their functional forms in the strongly-brokentheories as one slowly turns-on the infrared scales fromnear-zero values, it also gives us a new viewpoint of thenonperturbative origin of the parton distribution.The model system we choose to work with is the 2+1dimensional two-color ( N c = 2) QCD coupled to evennumber N of two-component massless Dirac fermions ina parity-invariant manner. In a previous study [52], wefound the global flavor symmetry in this system to bespontaneously broken for N (cid:46)
4, and conformal in the in-frared for N ≥ II. THREE-DIMENSIONALPARITY-INVARIANT SU(2) THEORY AND ITSSYMMETRIES
The zero temperature system is defined on three-dimensional Euclidean torus of physical extents (cid:96) × (cid:96) × (cid:96) with (cid:96) (cid:29) (cid:96) , (cid:96) . Here, we are using the conventionthat µ = 1 , x, y -directions, and µ = 3is the temporal t -direction. We will refer to the the as-pect ratio of the spatial slice as ζ = (cid:96) /(cid:96) . The parity-invariant three dimensional QCD consists of SU(2) val-ued gauge fields coupled to N flavors of two-componentfermions. Writing the action as S = S g + S f , the gaugeaction is S g = 14 g (cid:88) µ,ν =1 (cid:90) d x Tr F µν , (1)with F being the SU(2) algebra valued field strength.The important difference from the 3+1 dimensional QCDis that the gauge coupling g has mass dimension 1, mak-ing the theory super-renormalizable. This makes thescale-setting simpler, as one simply needs to measureall dimensionful quantities in the fundamental units of g . Once UV regulated, the continuum limit is simplyobtained by removing the regulator at the fixed valuesof quantities in units of g . Particularly for this work,it will also greatly simply our computation of the PDFcompared to 3+1 dimensions. The dimensionful natureof g also makes the theory trivially asymptotically free.The SU(2) gauge fields a µ are coupled to an even num-ber N = 2 n flavors of Dirac fermions, which are two-component spinors, in a parity-invariant manner; in or-der to make the action parity-invariant, n of the fermionflavors have mass + m and the other n have mass − m .We will refer to the fermions with positive mass as u andthose with negative mass as d . This is a deliberate choiceto be analogous to the light flavors in 3+1 dimensionalQCD. Throughout this paper, we will simply refer to thetwo-component Dirac fermions as quarks , to make theconnection to the real-world 3 + 1 dimensional QCD eas-ier. The N = 2 n flavor parity-invariant continuum actionis S f = n (cid:88) i =1 u i ( /D + m ) u i + n (cid:88) i =1 d i ( /D − m ) d i ; /D = (cid:88) µ =1 σ µ ( ∂ µ + ia µ )(2)with σ µ being the three 2 × d -quarks couple to as − /D † . This willbe the form of the lattice regulated fermion action. Wewill use the value of n as a tunable knob to control theinfrared fate of the theory.The explicitly massive theory has a globalSp( n ) × Sp( n ) symmetry [67] (with the symplecticgroup being special for SU(2) gauge group due to itbeing pseudo-real. For other generic color, it becomesU( n ) × U( n ) symmetry). At the massless point, thetheory has a larger Sp( N ) symmetry, which gets spon-taneously broken to Sp( n ) × Sp( n ) symmetry when N is smaller than some critical flavor N < N ∗ [68, 69].The conformal window extends for all N above N ∗ .Indications of non-zero value of N ∗ have been seen inprevious large- N Schwinger-Dyson Equation study [70]and in (cid:15) -expansion calculation [71]. In a previous latticestudy [52] to determine N ∗ , we found that it is likelyfor N ∗ to lie somewhere between 4 and 6, with N = 8showing strong evidences in the finite-size scaling oflow-lying Dirac eigen values for being scale-invariant inthe infrared, with mass-anomalous dimension γ m ≈ . N < N ∗ , the theory developsa scalar condensate, Σ ≡ (cid:10) u i u i − d i d i (cid:11) that sets theinfrared scale even after the box size taken to infinity,and sets the scale for the mass-gaps in the theory; thehadronic content in the SU(2) theory are mesons of thetype qq and diquarks (“baryons”) of the type q T τ q with τ being the Pauli matrix in color space. In thescale-broken sector, there will be 4 n Nambu-Goldstone(NG) modes. Of these, the 2 n NG modes will be themesons π + ij = d j u i ; π − ij = u j d i (3)which we simply refer to as pions in this theory, that cou-ple to the conserved flavor currents A µ,ij ( x ) = d j σ µ u i ( x ).Associated with the symmetry-breaking, there is the piondecay constant , (cid:104) |A ∓ µ,ij (0) | π ± ij ; p µ (cid:105) ≡ − ip µ F π , (4)with an on-shell pion at momentum p = ( p , p , E ). Theremaining set of NG modes will be of the diquark type u Ti σ τ d j and d Ti σ τ u j that couple to the correspond-ing conserved currents. Since these conserved currentsand extra NG modes are very special to the SU(2) the-ory, we simply consider the pions π + and π − that ex-ists for any number of color, and roughly belong to theU(2 n ) → U( n ) × U( n ) symmetry breaking part of the en-larged Sp(2 n ) → Sp( n ) × Sp( n ) symmetry-breaking pat-tern for SU(2) gauge theory .Before ending this section dealing with the system inthe continuum, we discuss the subtlety with parity sym-metry in the theory. The spatial parity P acts as x = ( x , x , x ) → x (cid:48) = ( − x , x , x ) , [ a ( x ) , a ( x ) , a ( x )] → [ − a ( x (cid:48) ) , a ( x (cid:48) ) , a ( x (cid:48) )] , [ u i ( x ) , d i ( x )] → [ σ u ( x (cid:48) ) , σ d ( x (cid:48) )] , [ u i ( x ) , d i ( x )] → [ − u i ( x (cid:48) ) σ , − d i ( x (cid:48) ) σ ] . (5) The mass-dimension of F π in d space-time dimensions is ( d − / / d = 3. We will therefore consider F π to be the IR scale inthe subsequent sections At the level of correlators after Wick contraction of fermions, thediquark correlators can be seen to be degenerate with the thatof mesons.
For a single two-component Dirac fermion, N = 1, thissymmetry is broken by the mass term and also becomesanomalous in the massless limit [72, 73]. While it appearsas though this is not a symmetry of Eq. (2) with even N ,in fact, it is a symmetry once the fermions are integratedout, or the symmetry can be made more obvious by per-forming a parity P transformation along with a pairwiseflavor permutation, G : [ u i ( x ) , d i ( x )] → [ d i ( x (cid:48) ) , u i ( x (cid:48) )].This GP operation is usually referred to as the parity inthe literature on parity-invariant theories in 2+1 dimen-sions. Since we are interested in hadron spectroscopy, itis easier to consider the usual notion of spatial parity P above, and whether bilinears are odd or even under it.As in 3+1 dimensions, the pions π and the correspond-ing current A µ in 2+1 dimensions are pseudo-scalars andaxial-vectors under the spatial parity P (However under GP , the bilinears can be linearly combined to becomeeven under it, but this does not play any role in our fur-ther discussions.) III. DESCRIPTION OF METHOD ANDSTATEMENT OF THE PROBLEM
We propose to study the internal quark structure ofNG boson as a function of the varying vacuum struc-ture by varying the number of massless fermion flavors,wherein the theory moves from being scale-broken to be-ing conformal in the infrared. In this section we first dis-cuss how to prepare a well-defined massive valence pionon top of a vacuum containing massless sea quarks, andthen define its valence PDF which we will use to charac-terize the pion quark structure.
A. Setting-up the computation such as to ensurenon-zero mass-gaps for all N In the scale-broken side of small N , the infrared con-tent of the theory are the mass-gapped hadrons, with thetypical gaps, denoted by M H , set by the IR scales suchas the condensate Σ, the decay constant F π , and in thecase of pure-gauge theory, the string tension σ . Thesenon-zero gaps survive the thermodynamic and masslessquark limit. On the other hand, in the conformal side ofthe theory, the eigenstates of the Hamiltonian are gap-pless and continuous in the thermodynamic limit. Weneed to deform the theory to introduce a mass-gappedspectrum in order for us to address the quark structureof a distinct ground state for any number of flavors. Sucha deformation would only be a sub-leading correction inthe scale-broken side, but will be the leading contributionin the conformal side. One can introduce the non-zeromass-gap (1) by studying the theory at finite spatial vol-ume, where the box size (cid:96) sets the infrared scale [74].That is, the mass-gaps become M H ( (cid:96), ζ ) = c ( ζ ) /(cid:96) where ζ is the aspect ratio of the two-dimensional spatial torus,and c is some function of ζ . (2) by introducing finitequark mass m in the theory [75], so that all the massesreceive finite mass corrections with M H ( m ).In this work, we will take a hybrid approach, which iseasier to implement in a lattice calculation, than beingtheoretically pristine. In order to capture the effect of N flavors of fermions without them getting decoupled inthe infrared, we sample the gauge configurations in thetheory coupled to N massless fermions in a finite spatialvolume, (cid:96) with ζ = 1. In the lattice theory terminology,the sea quarks are massless. On the gauge fields sampledthis way, we construct pion states π ij built out of u i and d j quarks which have a finite quark mass m , whichis tuned so that the mass of pion for any flavor N in (cid:96) spatial volume is an arbitrary chosen value. In the latticeterminology, the valence pion is made massive by tuningthe mass of the valence quark masses to non-zero values.This approach has the advantages that (a) the depletionof the infrared scales is preserved due to the presenceof massless fermions, and (b) the pion is massive even inthe scale-broken side which makes computation of matrixelements feasible without large periodicity effects [43],which is a technical boon.We chose the mass of the valence pion, M val π = 0 . g ,and kept this fixed across all N . The reason for thisvalue being that the pion is light enough to have thechiral properties in the scale-broken side, whereas in theconformal side, it will ensure that the mass M val π ( m, (cid:96), ζ ),which is now a function of valence quark mass and spatialvolume, is dominated by the finite m . This preference isbecause we will use hadrons boosted in the x -directionin our computation of PDF, which will cause a Lorentz expansion of extent (cid:96) to γ(cid:96) in the rest frame of thatstate, which effectively will decrease the aspect ratio ζ to ζ/γ in that state’s rest-frame [76]. Our choice of M val π is to minimize the effect of this variation in ζ on theenergy-momentum dispersion for the ground state in theIR conformal sector. One could also use lattices withsmall ζ (i.e., (cid:96) (cid:29) (cid:96) ) to begin with, but we realizedit post facto and intend to improve the calculation with ζ < B. Definition of u − d and valence PDFs Having described the preparation of valence pion stateof mass M val π above, we now specify how to study itsvalence PDF, f v ( x ), which we will use to characterize theUV quark structure of the pion. To define a valence PDF,we should first consider the u − d PDF of the π + ij = u i d j pion which has a well defined operator definition as f u i − d j ( x ) ≡ (cid:90) dξ − π e − ixξ − P + (cid:104) π ij ; P | O σ + | π ij ; P (cid:105) ,O σ + ( ξ ) = n (cid:88) k =1 (cid:18) u k ( ξ − ) σ + W + ( ξ − , u k (0) − d k ( ξ − ) σ + W + ( ξ − , d k (0) (cid:19) . (6) Here, the light-cone coordinates ξ ± = ( x ± x ) / √ σ ± = ( σ ± σ ) / √ W + ( ξ − ,
0) is the straight Wilson-line along the light-cone connecting the quark and anti-quark that are displaced by ξ − . Roughly speaking, thebilocal operator O σ + counts the number of massless u type quark minus the number of d type quark movinglong the light-cone, each carrying a fraction x of the mo-mentum P + . We have written the operator O formallyto be a singlet in the unbroken Sp( n ) × Sp( n ) symmetrybut non-singlet in the full Sp(2 n ) symmetry. For prac-tical purposes, we can simply speak of the operator ofthe type u i u i − d j d j for a pion of type π ij . Since themagnitude of the quark masses are all the same as | m | ,we will drop the indices i, j from π ij and f u i − d j ( x ). The u − d PDF has support from x ∈ [ − , GP symmetry ensures that f u − d ( x ) = f u − d ( − x ). Thus, we can write, f u − d ( x ) = (cid:40) . f v ( x ) for x > , . f v ( − x ) for x < . (7)This defines for us the valence PDF f v ( x ) of pion definedin x ∈ [0 , . Their moments are defined as (cid:104) x n (cid:105) u − d = (cid:90) − x n f u − d ( x ) dx ; (cid:104) x n (cid:105) v = (cid:90) x n f v ( x ) dx, (8)respectively. The even moments (cid:104) x k (cid:105) u − d = (cid:104) x k (cid:105) v , butfor the odd ones, (cid:104) x k − (cid:105) u − d = 0 whereas (cid:104) x k − (cid:105) v (cid:54) = 0.Since we can only determine f u − d via the well-definedoperator definition above, we will be inferring propertiesof f v ( x ) indirectly from f u − d ( x ) in this paper.With the set-up and key quantities defined, the precisequestions we want to address are the following. As weincrease N , the IR scales will vanish, and can be quanti-fied by how F π decreases. Is f v ( x ) of the pion sensitiveto the changes in the symmetry-broken vacuum given itsrole as the NG boson? If so, to what degree the PDFchanges with F π and what aspects of the pion valencePDF and its moments are sensitive to these changes? IV. LEADING-TWIST OPE IN A PLANARWORLD
The defining equation for the PDF involving the quarkand antiquark separated on the light-cone is given in Eq.(6). Instead, one can take the matrix element2 P + M ( ξ − , P + ) ≡ (cid:104) π ; P | O σ + | π ; P (cid:105) , (9)as the defining central object, which is also called as Ioffe-time distribution [77], and one can define the moments By defining antiquark distribution f q ( x ) = − f q ( − x ), and usingthe same symmetry argument, one can see that f v ( x ) = f u ( x ) − f u ( x ) for x > (cid:104) x k (cid:105) u − d of u − d PDF through its expansion as a function ν = P + ξ − , referred to as the Ioffe time in the literature, M ( ξ − , P + ) = M ( ν ) = ∞ (cid:88) k =0 ( − iξ − P + ) k k ! (cid:104) x k (cid:105) u − d with , (cid:104) π ; P | (cid:2) uσ + ( iD + ) k u − ( u ↔ d ) (cid:3) | π ; P (cid:105) ≡ P + ( P + ) k (cid:104) x k (cid:105) u − d . (10)Only even u − d moments are non-vanishing in the aboveequation for the pion. Since, the 2+1 dimensional QCDis super-renormalizable, owing to the dimensionful cou-pling, there are no UV divergences once the theory isregularized. Therefore, there are no additional scales µ entering the matrix elements defining (cid:104) x k (cid:105) in the equa-tion above, unlike in 3+1 dimensions. Hence, one cantalk of the PDF without referencing an MS renormaliza-tion scale of the PDF in 2+1 dimensions, as is the case inthe super-renormalizable 1+1 dimensional QCD as well.A brute-force Monte Carlo computation of M ( ξ − P + )is not possible due to the unequal time separation in theoperator O σ + . Computing the matrix elements of the lo-cal operators [78] defining the PDF moments in Eq. (10)is one possibility. Another recent method [31, 33], whichhas been proven to be very successful in 3+1d, is to com-pute the following equal time bilocal matrix element ofthe pion boosted with a momentum P = ( P , , E ( P )),2 E M B ( z , P ) ≡ (cid:104) π ; P | O σ ( z ) | π ; P (cid:105) ; where ,O σ ( z ) = u (0) σ W ˆ1 (0 , z ) u ( z ) − ( u ↔ d ) , (11)containing a purely spatial displacement z = ( z , , σ ,along the t -direction instead of the σ + present in Eq.(6). The straight Wilson line along the x -direction join-ing the quark and antiquark is denoted as W ˆ1 . This equaltime matrix element in 3+1 dimensions has been calledquasi-PDF matrix element, pseudo-ITD matrix elementor, simply as Ioffe-time Distribution in the literature. Inthis paper, we simply refer to the equal time correla-tion above as the bilocal quark bilinear matrix element(or simply as bilocal matrix elements) due to its centralrole in both quasi- and pseudo- approaches to PDF fromlattice, and the present work can be viewed from anyperspective the reader wants to approach it with. TheOPE of the above equal time bilocal matrix element [79],arranged by twist, gives M B ( z , P ) = ∞ (cid:88) k =0 ( − iz P ) k k ! (cid:104) x k (cid:105) u − d + O (cid:18) ( g z ) , ( F π z ) , ( M val π z ) (cid:19) , (12)with the first term being at leading twist, and the rest arehigher twist corrections due to the non-vanishing P and z present off the light-cone, unlike in the similar expres-sion Eq. (10) on the light-cone. The leading twist term isexactly the same as the one in Eq. (10). In 3+1 dimen-sions, the similar expression [34, 79] will involve a match-ing Wilson-coefficient c n ( z µ ) which is 1 at tree-level and the terms higher order in coupling capture log (cid:0) z µ (cid:1) divergence in the limit of z →
0. In the above expressionfor 2+1 dimensions, the Wilson coefficients take their treelevel value c n = 1, and there are no higher order pertur-bative corrections to this tree-level value at leading twist,making it exact at leading twist. This peculiarity in 2+1dimension arises because the coupling g has mass dimen-sion 1, which means that a perturbative correction to c n increases the twist of the term occurring in the OPE by1. Therefore, we have discarded such higher-order termsas ( g z ) higher-twist corrections. Along with such cor-rections, there could be other genuine higher twist cor-rections coming from higher-dimensional operators thatoccur in the OPE, which we have denoted by a ( F π z ) corrections. Even at leading twist, there will be targetmass corrections [80, 81] coming from the trace termswhich bring factors of P = ( M val π ) . We have denotedthese as the ( M val π z ) corrections above.In the above discussion, we have been a little cav-alier about the Wilson-line. In 3+1 dimensions, theself energy divergence of the Wilson-loop causes a non-perturbative exp( − cz /a ) suppression of Wilson-line as afunction of z [82–84]. The non-perturbative renormaliza-tion of the bilocal operator removes this nonperturbative z dependence. In 2+1 dimensions, there will insteadbe exp (cid:0) − c (cid:48) g z (cid:1) behavior as g is dimensionful; one wayto justify this is to see that the set of 1-loop diagramin real-space that contributes to the α s ( z a ) /a behav-ior of the bare Wilson-line in 3+1 dimensions, now con-tributes g ( z a ) /a ; where, the factor of ( z a ) in boththe dimensions comes from the integral measure whenthe end-points of the gluon loop on the Wilson-line be-come nearly coincident, whereas, the other 1 /a factorin 3+1 dimensions comes from | z | − gauge field prop-agator, and similarly the 1 /a factor in 2+1 dimensionscomes from the corresponding | z | − gauge field propaga-tor. The residual exp (cid:0) − c (cid:48) g z (cid:1) effect of the Wilson-loopinsertion, which is hadron momentum independent, isthen a non-perturbative higher-twist effect, which we re-move by forming ratio as done in 3+1 dimensions [34],namely˜ M ( z , P ) = (cid:18) M B ( z , P ) M B ( z , (cid:19) (cid:18) M B (0 , M B (0 , P ) (cid:19) , (13)which we expect to converge to leading-twist expansionin Eq. (12) better in a moderate range of z and P .The reason for the second factor in the above equationto ensure z = 0 matrix element is 1 by definition, sinceit is the charge of the pion. From the OPE, one canobtain the light-front Ioffe-time Distribution M from theEuclidean construction ˜ M in the limit, M ( ν ) = lim z → ,P →∞ P z = ν ˜ M ( z , P ) . (14)In practice however, we will simply be looking at a setof data that spans a range of z and P . If the leadingtwist expansion works, we expect a scaling ˜ M ( z , P ) =˜ M ( z P ) for all z and the range of P where the scalingviolations from z -type higher twist corrections in Eq.(12) are negligible. Based on fits of Eq. (12) to a subsetof data where the leading twist OPE works the best, wewill be able to infer M ( ν ), and the PDF and its moments. V. LATTICE METHODOLOGY ANDTECHNICAL DETAILS
In this section, we detail the lattice regularization of2+1 dimensional QCD in a parity-invariant manner, thegauge field statistics, and the construction of correlatorsrequired to build the pion bilocal matrix element.We regulate the system defined on (cid:96) × (cid:96) × (cid:96) on aperiodic L × L × L lattice with isotropic lattice spacing a = (cid:96) µ /L µ . In this paper, we will be using 28 × × U µ,x connecting the lattice site x to x + ˆ µ . The gauge action is the lattice regulated Wilsonsingle plaquette action, S g = − β (cid:88) µ>ν =1 (cid:88) x Re Tr P µν ( x ); β = 4 g a , (15)where P µν ( x ) is the SU(2) valued plaquette at lattice site x = ( x , x , x ). Periodic boundary condition is imposedon all three directions (an explicit anti-periodic boundarycondition in the temporal direction is superfluous as − β = 9 . L = L = 28 at this lattice spacing was found to be close tothe thermodynamic limit.The gauge fields are coupled to a system of N =2 n massless fermions, which we regulate by using two-component Wilson-Dirac fermions [52, 86, 87], definedusing the regulated Dirac operator, /D w = /D n + B + m w , (16)where /D n is the naive Dirac operator, /D n = 12 (cid:88) µ =1 σ µ (cid:16) U ( n ) µ,x δ x +ˆ µ,y − U ( n ) † µ,x − ˆ µ δ x − ˆ µ,y (cid:17) , (17) B is the Wilson term, B = − δ x,y + 12 (cid:88) µ =1 (cid:16) U ( n ) µ,x δ x +ˆ µ,y + U ( n ) † µ,x − ˆ µ δ x − ˆ µ,y (cid:17) , (18)and m w is the Wilson fermion mass in lattice units. Thelattice fermion is coupled to the gauge fields via n -stepStout smeared [88] gauge links, U ( n ) µ,x , in order to reducethe lattice artifacts coming from irrelevant UV fluctu-ations [89, 90], with the identification U (0) µ,x = U µ,x . We used 1-step stout smeared links in the Wilson-Dirac oper-ator with optimal value (cid:15) = 0 .
65 for the smearing param-eter. The lattice regulated action that is exactly invariantunder spatial parity is S f = n (cid:88) i =1 u i /D w u i − n (cid:88) i =1 d i /D † w d i , (19)making the partition function, Z = (cid:90) [ dU ] det (cid:16) /D w /D † w (cid:17) n e − S g , (20)with a positive definite measure that can be simulatedwith Monte Carlo algorithms. The theory only has anSp( n ) × Sp( n ) symmetry even when m w is tuned to themassless point, and the full Sp(2 n ) symmetry will be re-covered in the continuum limit. A. Gauge field generation
We studied the theories with N = 0 , , , β = 9 . ×
48 lattices. We generated gauge configurations us-ing the standard Hybrid Monte Carlo algorithm [91] us-ing n copies of Gaussian noise vectors to sample the de-terminant det (cid:16) /D w /D † w (cid:17) . We tuned the value of the Wil-son mass m w to the approximate massless point such thatthe smallest Dirac eigenvalue Λ ( m w ) has a minimum atthe tuned m w in the finite fixed volume. Since the Diraceigenvalues are gapped in finite volume, the eigenvaluesoccurring are not zero at the massless point, and hencemakes the HMC tractable. For N = 2 , ,
8, the valuesof sea quark mass m w = − . , − . , − . N . For thermalization, wediscarded the first 400 trajectories in each stream thatwere started from random configurations. After that,gauge configurations every 5 trajectories were stored andused for correlator measurements. This way, we gen-erated 24.5k, 25.2k, 27.3k and 30.2k configurations for N = 0 , , , ∼
100 trajectories).
B. Choice of valence quark masses to create amassive valence pion
Using the configuration generated with near masslesssea quarks, we tuned the valence Wilson-Dirac quarkmass to produce a valence pion ground state of mass M val π = 0 . g in units of g , or M val π a = 0 . N = 0 , , , m w and in-terpolating the M val π dependence on m w , and zero-in onthe exact tuned mass value. This way, we found thetuned valence quark mass corresponding to 0 . g pionmass to be m val w = − . , − . , − . − . N = 0 , , /D − w . For this choice of valence pion mass, thevalues of e − M val π ( aL − t s ) = 0 . t s = 12 a , implyingonly a small periodicity effect of 0.4% when operators aretemporally separated by 12 lattice units. C. Two point function computations
The first step is to find the ground and the excited statecontributions to the pion two-point function. Since theleading twist formalism demands boosted pion states, weconstruct pion sources that project to definite momentumstates. Namely, we construct the two point functions, C ( t s ; P ) = (cid:68) π S ( x , t s ) π † S ( P , (cid:69) , with ,π S ( P , t s ) = (cid:88) x d ( x , t s ) u ( x , t s ) e − i P · x , (21)using smeared source π S ( P ,
0) and smeared sink π S ( P , t s ) that project to momentum P = ( P , x per configuration. We usethe momenta, aP = 2 πL n , (22)for n = 0 , , , ,
4. At the given fixed β , these mo-menta correspond to P /g = 0 . , . , . . .
09 respec-tively in units of coupling g . In order to suppress thetower of excited states, we use Wuppertal smeared quarksources [94] to construct the pion source and sink. Forthis, we used 10 steps of two-dimensional stout smearedlinks to construct the smearing kernel with smearingparameter (cid:15) = 0 . P = 0,we found the optimal number of steps n wup of Wup-pertal smearing to be 80 with Wuppertal smearing pa-rameter δ = 0 . N = 0 , N = 4 wefound ( n wup = 120 , δ = 0 . N = 8 we found( n wup = 160 , δ = 0 . k = ζ (cid:48) P which then are used to construct theWuppertal sources (with the same smearing radius as at P = 0). We found the optimal boost parameter ζ (cid:48) for n = 1 , , , (cid:104) d ax d by (cid:105) = (cid:16) [ − /D † w ] − (cid:17) abxy and (cid:104) u ax u by (cid:105) = (cid:0) [ /D w ] − (cid:1) abxy .One can then simply use (cid:16) [ − /D † w ] − (cid:17) xy = (cid:0) [ − /D w ] − (cid:1) † cs yx with A † cs meaning a conjugate transpose of A only overcolor-spin space, thereby halving the number of inver-sions, just like in 3+1 dimensions. We used ConjugateGradient algorithm for inversion here (and also in theHMC) with a stopping criterion of 10 − . D. Three point function computations
The next important ingredient in the PDF computa-tion is the three-point function between the pion source,pion sink and the bilocal operator O σ , C ( t s , τ ; P , z ) ≡ (cid:68) π S ( x , t s ) O σ ( z ; τ ) π † S ( P , (cid:69) , (23)with the zero momentum projected bilocal operator thatis inserted at a time-slice τ between the pion source andsink at time-slice t s , O σ ( z ; τ ) = (cid:88) x (cid:18) u ( x ) σ W ˆ1 ( x, x + L ) u ( x + L ) − d ( x ) σ W ˆ1 ( x, x + L ) d ( x + L ) (cid:19) ; x = ( x , τ ) . (24)Here, the quark and antiquark are separated along the x -direction by L = ( z , , W ˆ1 = (cid:81) x (cid:48) ∈ [ x,x + L ] U ( n )1 ,x (cid:48) . We used only 2-level stout smearedlinks U (2)1 ,x (cid:48) for this construction, so as to not risk thesmearing to spoil the UV physics. The pion source andsink are smeared using the same set of parameters usedfor the corresponding two-point function. It should benoted that the u and d quark operators appearing in O σ are simple point operators.The contractions for the above three-point functionwere performed using the sequential-source trick (c.f.,appendix of [45] for details relevant to the bilocal op-erator) to take care of the necessary Fourier summationover two-dimensional time-slice at the sink. It shouldbe noted that, similar to the u − d three-point functionof the pion in 3+1 dimensions, the three-point functionis purely real at all P . Also, there are no fermion-linedisconnected pieces; this comes non-trivially at finite lat-tice spacing, by the parity invariance of the action which . . . . .
20 5 10 15 20 25 n = 0 n = 1 n = 2 n = 3 n = 4 N = 0 E e ff a t s /a . . . . .
20 5 10 15 20 25 n = 0 n = 1 n = 2 n = 3 n = 4 N = 2 E e ff a t s /a . . . . .
20 5 10 15 20 25 n = 0 n = 1 n = 2 n = 3 n = 4 N = 4 E e ff a t s /a . . . . .
20 5 10 15 20 25 n = 0 n = 1 n = 2 n = 3 n = 4 N = 8 E e ff a t s /a . . . . . . . . . . . . . . N = 0 E ( P ) a P aE ( P ) E ( P )from A continuum disp.lattice disp. 0 . . . . . . . . . . N = 2 E ( P ) a P aE ( P ) E ( P )from A continuum disp.lattice disp. 0 . . . . . . . . . . N = 4 E ( P ) a P aE ( P ) E ( P )from A continuum disp.lattice disp. 0 . . . . . . . . . . N = 8 E ( P ) a P aE ( P ) E ( P )from A continuum disp.lattice disp. FIG. 1. The spectral content of pion two-point functions in N = 0 , , , t s /a in lattice units. The bands are the expected effective mass from two-state fits to the two-point function.(Bottom) The dispersion relation between ground state (red circles) and the first excited state (purple triangle) energy onboosted momentum is shown. For comparison, the expected single particle dispersion in the continuum (blue dashed curve)and on the lattice (black curve) are shown. The ground state pion energy extracted from axial-vector A correlator are alsoshown at non-zero momenta. guarantees that (cid:104) Tr( /D − w ) (cid:105) = −(cid:104) Tr( /D † w − ) (cid:105) . If one used2+1 dimensional overlap fermions [96], the cancellationof disconnected piece would have been on each configu-ration. E. Pion decay constant computations
We will quantify the presence of infrared scale in thesystem using the pion decay constant, F π , defined inEq. (4). We extracted this matrix element using theaxialvector-pion two-point function (c.f., [97]), C π −A ( t s ) ≡ (cid:68) A ( x , t s ) π † S ( P = 0 , (cid:69) ; A ( x ) = d ( x ) σ u ( x ) . (25)The pion source was optimally Wuppertal smeared,whereas the current is constructed out of point quarkoperators. We will describe the analysis of the two-pointfunction leading to F π in a subsequent section. VI. ANALYSIS OF CORRELATOR DATA TOOBTAIN THE PION BILOCAL MATRIXELEMENT & F π A. The spectral content of pion correlator
The two point function in Eq. (21) gives informationon the spectrum of states contributing to the pion quan-tum number, that we will use to extract the ground stateboosted bilocal matrix element. In the top panels of Fig.1, we have shown the effective masses for the pion at all flavors as a function of source-sink separation t s , both inlattice units. For each flavor, the effective masses at thefive momenta are shown by the different symbols. First,one can see that the ground state displays a well definedplateau for all N , even for N = 8, thereby demonstratingthe effectiveness of gapping the spectrum by finite valencequark mass and volume even in the otherwise conformalinfrared theories. We can see that the value of the groundstate mass has been tuned well to be ≈ . a in all thetheories, which corresponds to 0 . g physical mass. Thesmallest t s from where one can see a well-defined plateau,at least for the smallest three momenta, increases withmomenta due to the decreasing gap with the excited statewith the boost. However, we have tuned the Wuppertalsmearing and quark boost parameters precisely to reducethe amplitudes of the excited state as best we could, andthe any observed deviation from the plateau at smaller t s was the best we could reduce it to, without compro-mising on the noise at larger t s . Since the range of t s that we will use to analyze the three-point function fallsin the small t s region without the plateau, we need tounderstand the spectral content of C better.We take the spectral decomposition of C , C ( t s ; P ) = N state − (cid:88) i =0 | A i | (cid:16) e − E i ( P ) t s + e − E i ( P )( aL − t s ) (cid:17) , (26)and truncate it at N state = 2, which is referred to as thetwo-state ansatz. We found that this is enough to de-scribe the behavior of C ( t s ) for all the P used, evenstarting from t s = 2 a and be able to reproduce the valueof ground state E , as obtained from one-state fit with . . . . . . . .
640 2 4 6 8 10 12 N = 0 R ( t z , τ ) τ z = 3 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a . . . . . . . . . .
60 2 4 6 8 10 12 N = 2 R ( t z , τ ) τ z = 3 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a . . . . . . . . .
60 2 4 6 8 10 12 N = 4 R ( t z , τ ) τ z = 3 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a . . . . . .
640 2 4 6 8 10 12 N = 8 R ( t z , τ ) τ z = 3 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a − . − . . . . .
20 2 4 6 8 10 12 N = 0 R ( t z , τ ) τ z = 6 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a − . − . − . . .
10 2 4 6 8 10 12 N = 2 R ( t z , τ ) τ z = 6 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a − . − . − . − . .
050 2 4 6 8 10 12 N = 4 R ( t z , τ ) τ z = 6 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a − . − . − . − . − . − . − . − . − . − . − .
020 2 4 6 8 10 12 N = 8 R ( t z , τ ) τ z = 6 a,n z = 3 t s = 4 at s = 6 at s = 8 at s = 10 at s = 12 a FIG. 2. The estimation of matrix element h ( z , P ) by the two-state fits to the ratio, R ( t s , τ ). Some sample fits to theoperator insertion point, τ , and the source-sink separation, t s , dependence of R are shown in the different panels; the top panelsare at z = 3 a and the bottom ones are z = 6 a , at a momentum of n = 3. For each z , the panels from left to right are fromdifferent number of flavors. the minimum t s > a . The uncertainly bands for theeffective mass curves for the different P based on thetwo-state fits over the range t s ∈ [3 a, a ] are also shownin Fig. 1 along with the data. The quality of the fits areseen to be good, which is also reflected in χ / dof ≈ (cid:104)A (0) A ( t s ) (cid:105) correlators also; at P = 0, thisgives the mass of the axial-vector meson, but at non-zeromomentum the lowest mass comes from the pion due tothe non-zero overlap ∼ P F π with the lighter pion state.Thus, at non-zero momenta this gave a cross-check onthe determined ground-state values of the pion.In the bottom panels of Fig. 1, we have shown the bestfit values of the ground state energy E and the first ex-cited state E as a function of boosted momentum P .The different panels are again for the five different N . Wehave compared the data for E ( P ) with the curves forthe single particle dispersion in the continuum, E ( P ) = (cid:113) P + M val π , shown as the blue dashed curves. There isa slight discrepancy which increases with P , as P a ≈ E ( P ) a ) = cosh (cid:0) M val π a (cid:1) + 1 − cos( P a ). Thislattice single particle dispersion curve is shown as blackcontinuous curve in the figures. The nice agreement tellsus that there are possible lattice corrections at the levelof 3 −
4% at the highest momenta, which can be reducedin the future by going to much finer lattices. But thiseffect will persist for all the N , and therefore, we do notexpect this to affect variations as a function of N that weare interested in. As a cross-check, we have also shownthe values of ground-state masses of pion as extractedfrom the axial-vector A correlator at non-zero momenta,which can be seen to agree with the values from the sim- ple pion correlator. While the slight disagreement withthe dispersion curve at higher momenta are understoodas lattice spacing effect, a slight disagreement at the levelof 4% is also seen at the smallest non-zero momentumcorresponding to n = 1 for N = 4 and 8. This tells usthat the valence pion mass for the near-conformal andconformal theories mildly originate from aspect-ratio( ζ )-dependent 1 /(cid:96) effect that we described in Section III, inspite of our effort to use somewhat larger value of valencequark mass. As the pion is boosted, the aspect ratio inthe boosted frame decreases, and causes the observedsmall discrepancy at the smallest non-zero momentum.At the larger momentum, these aspect-ratio variationsare not important as the leading E ∝ P relativistic de-pendence takes over. Therefore, as a precaution, we willavoid using n = 1 momentum in our analysis of three-point function to avoid systematic effects. In a futurecomputation, we aim to rectify this by using lattices with ζ < B. The extraction of the pion bilocal matrixelement from three-point function
The required ground state matrix element of the bilo-cal operator can be obtained from the spectral decompo-sition of the three-point function, C ( t s , τ ; P , z ) = N state − (cid:88) i,j =0 A ∗ i A j h ij ( z , P ) e − E i ( t s − τ ) − E j τ , (27)0 − . . . . . − − − N = 0 ˜ M B ( z , P ) z /a . g . g . g . g − . . . . . − − − N = 2 ˜ M B ( z , P ) z /a . g . g . g . g − . − . . . . . − − − N = 4 ˜ M B ( z , P ) z /a . g . g . g . g − . − . . . . . − − − N = 8 ˜ M B ( z , P ) z /a . g . g . g . g FIG. 3. The matrix elements at N = 0 , , , z /a . The different colored symbols correspond to different fixed momenta P , which are shown in thelegend in units of gauge coupling g . with the amplitudes A i and energies E i , being thesame as obtained from the two-point function anal-ysis. The matrix elements are the terms h ij = (cid:104) E i , P |O σ ( z ) | E j , P (cid:105) . Therefore, the leading term h , is the required ground state matrix matrix element (cid:104) π ; P |O σ ( z ) | π ; P (cid:105) .We extracted this leading term by fitting the t s , τ de-pendencies of the actual C data at various fixed z and P using the above spectral decomposition truncated to N state = 2 (since we found that N state = 2 was able todescribe the corresponding two-point function well evenfrom small t s ≈ − a ). In practice, we constructed theratio, R ( t s , τ ) ≡ C ( t s , τ ) C ( t s ) , (28)with the P and z arguments being the same for bothnumerator and denominator, and hence notationally sup-pressed above. We then fitted R ( t s , τ ) using the ratio ofexpressions in Eq. (27) and Eq. (26), with h ij as the fitparameters. We took the values of the amplitudes | A i | and energies E i from the two-state fit analysis of C with the fit range over t s ∈ [3 a, a ]. We used these( A i , E i ) from the same Jack-knife blocks as used for thethree-point function analysis. We performed these fitsover τ ∈ [2 a, t s − a ] to reduce larger excited state effectsfor insertion closer to source and sink. Further, we usedall the data for t s /a ∈ [6 , , [6 , , [6 , , [6 , , [6 , n = 0 , , , , t s = 10 a, a for n = 0 in order to reduce the 0.4%lattice periodicity effect due to the smaller E , and sim-ilarly we skipped only t s = 12 a for the larger n = 1momentum. We did not used t s = 12 a for n = 4 as thetwo-point function for t s > a was very noisy. While theextrapolated values were insensitive to changes in fittingwindows, we kept the fit systematic fixed for all N sothat even if there is any unnoticed systematic error, it isunlikely to affect the overall variations in the data as afunction of N , which is our interest in this paper. Sucha two-state fit to the ratio R resulted in good fits for all z and P .Some sample data for R along with the results fromthe fits are shown in Fig. 2 for the case of momentum n = 3. The top and the bottom panels are for fixed z = 3 a and 6 a respectively, with the different N shownin the different columns. The fits, shown as bands, agreewith the data well for all N , and the extrapolated valueis shown as the horizontal band. The ground state ma-trix element M B ( z , P ) = h ( z , P ) so extracted, areshown as a function of z /a in Fig. 3; with each panelfor different N , and in each panel, the data for differentmomenta differentiated by color and symbols. The datahas not been symmetrized by hand with respect to z and − z , so as to show that the symmetry emerges automat-ically, which is a simple cross-check on the computation.The local matrix element corresponding to z = 0 shouldbe precisely be 1 if we had used an exact conserved cur-rent on the lattice, which we have not. So, one sees thematrix element at z = 0 to be slightly below 1 at z = 0and this difference with 1 increases with larger momen-tum; we found this lattice spacing effect to be of the kind( aP ) , which was also seen in 3+1 dimensional computa-tion [43]. We will see that such effects are nicely canceledby an overall normalization such that z = 0 matrix ele-ments are exactly 1; this is justified since the informationon the PDF comes from the variations in z and P , andnot by a fixed overall normalization. C. Determination of pion decay constant
We determined the pion decay constant through thespectral decomposition of the correlator C A− π in Eq. (25)as, C A− π ( t s ) = − F π (cid:112) M val π M val π A (cid:16) e − M val π t s − e − M val π ( aL − t s ) (cid:17) + A (cid:48) e − E (cid:48) t s + . . . , (29)which follows from Eq. (4) with µ = 3. The factor of (cid:112) M val π is to convert the lattice normalization of statevectors to the relativistic one used in defining F π . Thefactor A is the amplitude (cid:104) | π S | π (cid:105) , which we take fromthe smeared-smeared pion two-point function at zero mo-mentum; we only determine the magnitude of A , andtherefore we are assuming there is no phase in A . The1 . . . . . . . . . . .
220 5 10 15 20 25 F e ff π a . t s /a N = 0 N = 2 N = 4 N = 8 FIG. 4. The determination of pion decay constant F π . Thedata points are the effective F eff π ( t s ) as determined from the (cid:104)A π (cid:105) ( t s ) correlator, and the curves are fits to F π + Ae − δmt s . correlator is anti-periodic in the t -direction, which can beseen by a rotation in the xt -plane; taking ( x , x , x ) → ( − x , x , − x ) along with u, d → σ u, σ d . The excitedstate contributions captured by e − E (cid:48) t s . We fit the abovefunctional form along with a subleading excited state con-tribution to the C A− π ( t s ) correlator to determine F π .Such fits worked well even starting from t s = 2 a with thefitted value of F π independent of the fit range. In Fig. 4,we show an effective F eff π ( t s ) obtained by inverting right-hand side of Eq. (29) without excited state term to geta t s dependent value of F π . The curves in the plot arethe expectations for F eff π ( t s ) from the excited state fits,which can be seen to perform well. From this analysis,we find F π a . = 0 . , . , . , . N = 0 , , , VII. RESULTS
We will present the results in the following logical or-der; first, we will explain how we measure the presenceof infrared scale, by which we establish that the infraredscales are indeed depleted as a function of number ofmassless fermion flavors N . Then, we will infer the Mellinmoments of the PDF and reconstruct the PDF based ona two-parameter model, and see how the PDF relatedquantities change as a function of N . This induces acorrelation between the infrared scale and the PDF pa-rameters, which we look for. A. The depletion of infrared scales
For the N = 0 pure-gauge theory, the confining in-frared can simply be characterized by the string-tension,which takes a value √ σ = 0 . g for the SU(2) the-ory [85]. For theories with non-zero N , string-tension isnot a good parameter to use; instead we use the conden-sate Σ and the decay constant F π . In a previous studywith R. Narayanan [52], we measured the scalar conden- . . . . . . .
140 0 .
02 0 .
04 0 .
06 0 .
08 0 . √ Σ g − F π g − . . . . . . . .
40 0 .
02 0 .
04 0 .
06 0 .
08 0 . ( M A − M π ) g − F π g − FIG. 5. The changes to the infrared scales when the numberof massless quark flavor is changed. (Top panel) the quarkcondensate as a function of decay constant, both implicitlydepending on the number of quark flavors. (Bottom panel) asimilar plot showing how the mass gap between the pion andthe axial vector changes as a function of decay constant. sate Σ that breaks the Sp( N ) global flavor symmetry,as a function of N in the massless limit of the theory.For this, we compared the finite-size scaling (FSS) ofthe low-lying eigenvalues, λ i ∝ z i Σ − (cid:96) − behavior of theeigenvalues in an (cid:96) box, where the proportionality con-stant z i are the eigenvalues of the random matrix modelfrom the non-chiral Gaussian Orthogonal Ensemble [98],which shares the same symmetries as the Dirac operatorcoupled to SU(2) gauge field in 2+1 dimension. The co-efficient Σ is the condensate in the massless limit. Wefound non-zero Σ /g = 0 . , . , . . − for N = 0 , , , N = 8 and 12 theories were instead likely to be infraredconformal with non-trivial mass anomalous dimensions γ m = 0 . γ m = 0 . λ i ∝ (cid:96) − γ m − rather than an (cid:96) − FSS expected from SSB. Here, weshould remark that we found that it was also possibleto describe the eigenvalue FSS in the N = 4 theory as-suming a rather large γ m = 0 . /(cid:96) FSS corrections; however, in light of theresults on F π in this work, it appears that the N = 42theory indeed is more likely to be scale-broken in the IR.The F π we determined here are at finite volume and finitevalence pion mass, but it cannot change the non-zero F π result for N = 4 because of the following. One shouldnotice that the very small F π = 6 × − g for N = 8is most likely to arise due to finite volume and valencequark mass, and hence it gives the typical correction to F π due to these effects; the value of F π for N = 4 theoryis F π = 0 . g , which is much larger than those typi-cal corrections, and hence justifying our inference aboutthe IR fate of N = 4. Thus, our current understandingabout the IR fate of many-flavor SU(2) gauge theory isthat N = 0 , , N ≥ √ Σ /g , as a function of another scale, F π /g . The two seem to be almost directly proportional.Another infrared scale one could use is the mass-gap, M A − M val π between the pion and the axial-vector meson.In the bottom panel of Fig. 5, we correlate this mass-gapwith F π . Again, it is clear than the mass-splitting alsoshrinks with the other diminishing, perhaps a more fun-damental scale, F π . The mass-splitting does not go tozero even for N = 8 most likely because of the finitefixed volume and the quark mass. The one-to-one posi-tive correlation between the infrared scales also suggeststhat we can now make the number of flavors implicit, andsimply ask for the effect of reducing one infrared scale onanother, as done in Fig. 5. As one would have expected, afactor reduction in F π induces a reduction in other scalesby a similar factor. We now apply this perspective toquark structure of pion, where the effect is not obvious. B. Response of the pion PDF to changes ininfrared
The central object in our analysis of PDF is theequal time bilocal matrix element ˜ M ( z , P ) in Eq.(13), formed by taking ratios of the matrix elements M B ( z , P ) that we obtained directly from the three-point function analysis. We formed these ratios to re-move the presence of exp (cid:0) − c (cid:48) g z (cid:1) behavior due to theWilson-line which is present in the definition of the bilo-cal operator, and hence ensure a better description bythe OPE. In Appendix A, we describe features of M B itself, and here we proceed with using the improved ˜ M .Through its OPE in Eq. (12), ˜ M contains the leadingtwist terms that relate to the pion PDF as well as con-tribution from operators with higher-twist, which couldhave interesting physics in their own right, but for ourpurposes here are contaminations. We can distill theleading-twist PDF terms in a practical lattice compu-tation, when z is small and P is large, so that one has − . − . . . . .
81 0 2 4 6 8 10 12 ˜ M ( z , P ) z /a N = 0; n = 32-stout6-stout FIG. 6. A cross-check that at least an application of smallnumber of Stout smearing to the Wilson-line connectingquark-antiquark is harmless. The pion bilocal matrix elementat n = 3 momentum in N = 0 theory with 2-stout and 6-stout smeared Wilson line insertions are compared and shownto be consistent once the ratio is taken. a finite range of z P which can simply be described thelead-twist part of the OPE in the analysis. Before go-ing further, we need to first make sure that the ratio inEq. (13) indeed cancels any remnant non-perturbative z -dependent factor from the usage of Wilson line in theoperator. For this, we performed the computations of˜ M ( z , P ) with 2-stout and 6-stout smeared Wilson linesfor a sample case with n = 3 momentum in the N = 0theory. The results from the two are compared in Fig.6, where one can see a good agreement between the two,giving confidence that the results are not stout smear-ing dependent, at least for few steps of it. The resultsget less noisier when number of stout smearing steps in-creases, but we use 2-stout in order to be conservative.We perform two types of analysis on ˜ M , namely, amodel-independent determination of the even momentsof the valence quark PDF, and secondly, by model-dependent reconstruction of x -dependent PDF. For boththe ways, the starting point is the working version ofthe leading twist OPE in Eq. (12) with some unknownhigher-twist z corrections, namely,˜ M ( z , P ) = 1+ (cid:34) N max (cid:88) k =1 ( − k ( z P ) k (2 k )! (cid:104) x k (cid:105) v (cid:35) + bz . (30)We have rewritten leading twist part of Eq. (12) in adifferent form above so that is clear that (cid:104) x (cid:105) = 1, ˜ M is purely real and that only even valence PDF moments (cid:104) x k (cid:105) v appear. These are very specific properties of ˜ M for a pion in 2+1 as well as 3+1 dimensions. The upper-cutoff of the sum N max is infinity but for practical imple-mentation, we need to work with smaller N max since thedata is only sensitive to some smaller powers k . Here, themoments of the PDF are the unknowns we are interestedit, but we will also fit the parameter b to effectively takecare any higher twist g z, F π z corrections, and also anytarget mass corrections. We also tried adding lattice cor-rections of the form ( aP ) to the above functional form3 − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 0 P = 1 . g P = 1 . g P = 2 . g − . − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 2 P = 1 . g P = 1 . g P = 2 . g − . − . − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 4 P = 1 . g P = 1 . g P = 2 . g − − . .
51 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 8 P = 1 . g P = 1 . g P = 2 . g FIG. 7. The pion bilocal quark bilinear matrix elements, ˜ M ( z , P ) from N = 0 , , P = 1 . , . , . g are put together and shown as a function of z P .The color and symbols differentiate the data at fixed P . The bands are the expectations for ˜ M ( z , P ) based on the fits tothe leading twist expression in Eq. (30). The bands extend over points included in the fit. of OPE [43], but such terms were found to be unnecessaryand consistent with zero well within errors. Therefore, wedo not present such an analysis here.In any method of analysis, we need to choose the rangeof z and P carefully, since we will not be taking either z → P → ∞ limits, and instead we will simply befitting the data which spans a finite range of z and P .First, we will work with momenta P /g ≥ z , a term like ( P z ) k is larger than a similar order term ( g z ) k . This leavesthe momenta corresponding to n = 2 , ,
4. Through thischoice, we are also guaranteed that M val π /P and F π /P corrections would also be controlled. For the range ofquark-antiquark separation z , we have two choices; wemight want z g < z F π <
1, where the first factoris simply due to the superficial dimensional scale in thesystem and the second is the natural infrared scale. Weassume that the superficial scale will arise simply duethe exp (cid:0) − cg z (cid:1) -type Wilson-line term which we find tobe nicely canceled in the ratio ˜ M . Due to the naturalinfrared scales being at least a factor 10 smaller than g for N = 0, and even smaller for larger N , even a usageof z = 10 a will only lead to F π z = 0 . z ∈ [ a, a ] and changethe maximum z to 6 a and 10 a to check for the robustnessof results. The justifications for the used ranges of z , P , will also bear out in the data.
1. Model-independent inferences
In the model independent analysis [43, 99], we first fitEq. (30) to ˜ M ( z , P ) data over the specified range of z and P with the even moments (cid:104) x k (cid:105) v being the fitparameters. In addition, we also fitted the high-twist pa-rameter b to take care leading higher twist effects; buttheir values were consistent with zero, and when we per-formed the fits without the higher-twist corrections, theresults were consistent with the one including it. Here,we keep this correction nonetheless. Since the valencequark PDF is positive (since the anti- u quark and d quarkarises only radiatively in ud pion , whereas the u quarkis present at tree-level itself), it imposes a set of inequal-ities to be satisfied by the moments as discussed in [43];with the important one being (cid:104) x k (cid:105) v < (cid:104) x m (cid:105) v for k > m .We imposed these constraints in the fit using the meth-ods discussed in [43]. With such constraints, one can addas many moments, N max , in the analysis without over-fitting the data, except that it will result in many of thehigher-moments, which the data is not sensitive, to con-verge to zero. We found that N max = 5 was sufficient todescribe the data in the range we fitted, as we describe4 − . − . − . − . . . . .
81 0 1 2 3 4 5 6 7 8 M ( ν ) ν F π g − . × − . × − . × − . × − FIG. 8. The effect of reduction in F π on the bilocal matrixelement (Ioffe-time distribution) M ( ν ) is shown. The bandsare inferred from the fits to lead-twist expression in Eq. (30)and taking its z → , P → ∞ limit at fixed P z = ν . Thematrix element is shown as a function of ν = P + ξ − (Ioffetime). The different colored bands are at different F π g − . below.In Fig. 7, we show the data for ˜ M ( z , P ) as a func-tion of z P ; the data from different fixed pion momenta P are differentiated by the colored symbols. The fourdifferent panels show the data from N = 0 , , , z P , which demonstrates the dominance ofthe lead-twist part, as is essential for this work. As seenby the early peeling-off of P = 1 . g data from thehigher momenta data beyond P z > N = 4 , N = 0 than for N = 8. This could be because the nat-ural higher-twist scale in the broken phase is F π whichis smaller than the natural scale corrections g z , and thefinite-box scale, z/(cid:96) , which could be important in the con-formal phase. However, for the range of z = 6 a, a, a ,this ensuing higher twist effect is less important even for P = 1 . g , and definitely not important for higher mo-menta. This justifies our choices of fit ranges and the rea-soning we presented before. The data gets increasinglyprecise with increasing N because the fluctuations in thegauge field decreases roughly as 1 / √ N for larger N . Thebands of various colors in Fig. 7 are the expectations for M ( z , P ) from the best fits from the analysis; the colorsmatch the corresponding color for the momenta for thedata. The bands cover the range of P z for each P for afixed range of z up to 8 a , and for P = 1 . g this rangeis within the point where the higher twist effects startbecoming visible at this lower momentum. The qualityof fits are very good with resulting χ / dof ≈ . M ( z , P ) and extracted the light-front bilocal matrix element, M ( ν ), which is the Ioffe-time Distribution (ITD). As we discussed, the variationof infrared scales with N , induces a direct infrared scaledependence of various quantities. Therefore, in Fig. 8, . . . . . . . . .
450 0 .
02 0 .
04 0 .
06 0 .
08 0 . h x n i v F π g − h x ih x ih x ih x ih x i− . − . − . − . . . . . . .
120 0 .
02 0 .
04 0 .
06 0 .
08 0 . κ n F π g − κ κ FIG. 9. (Top) The correlation between decay constant andthe valence PDF moments. The filled symbols are obtainedfrom model-independent fits and the open ones from model-dependent PDF ansatz fits. For model-independent fits onlythe even moments are directly obtained, whereas the odd mo-ments (cid:104) x (cid:105) v and (cid:104) x (cid:105) v were obtained indirectly by definition inEq. (32). The curves are quadratic fits in order to interpolatethe data. (Bottom) The correlation between decay constantand the cumulants κ and κ of u − d PDF. we have shown the ITD dependence on the decay con-stant. This is the main result in this paper, which wewill process further and look at from various angles. Asthe infrared scale-breaking is made stronger, as reflectedin F π , the corresponding valence quark ITD starts peel-ing off from the large- N conformal curve at shorter andshorter ν . In this process, however, the ITD remains al-most universal up until ν ≈
3. This tells us that thelowest non-zero u − d moment, (cid:104) x (cid:105) u − d = (cid:104) x (cid:105) v , mustremain quite insensitive to the scale changes. Thus, theeffect of scale-breaking seems to be encoded in the fall-off rate of the ITD for ν > x behavior, which is typically modeled as aRegge-type x α behavior, and also in the large- x , (1 − x ) β behavior of the underlying valence PDF. This is because,the tail of the ITD typically carries information on thesmall- x asymptotic of the PDF, whereas given the infer-ence that the lowest moment will be almost fixed, will5induce a variation in β as well via the implicit relation (cid:104) x (cid:105) v ( α, β, . . . ), with ( α, β, . . . ) being the parametrizationof the shape of the PDF.In the top panel of Fig. 9, we have plottedthe F π dependence of the first three even moments (cid:104) x (cid:105) v , (cid:104) x (cid:105) v , (cid:104) x (cid:105) v , as obtained directly from the model-independent analysis discussed above, using the closedsymbols. Since we can directly get only the even valencemoments, we infer the odd moments from (cid:104) x (cid:105) v and (cid:104) x (cid:105) v by assuming a two-parameter PDF ansatz, f v ( x ) = N x α (1 − x ) β , (31)with normalization N to ensure (cid:104) x (cid:105) v = 1, and simplysolve for α and β through the two equations,Γ(3 + α )Γ(2 + α + β )Γ(1 + α )Γ(4 + α + β ) = (cid:104) x (cid:105) v , Γ(5 + α )Γ(2 + α + β )Γ(1 + α )Γ(6 + α + β ) = (cid:104) x (cid:105) v . (32)Through this we get (cid:104) x k − (cid:105) v ( α, β ) by this semi-model-dependent analysis. As a cross-check, this procedure alsopredicts the even moment (cid:104) x (cid:105) ( α, β ), which we foundto agree well with the actual value we obtained in themodel-independent analysis. These odd moments (cid:104) x (cid:105) v and (cid:104) x (cid:105) v are also shown in Fig. 9. The inferred valueof (cid:104) x (cid:105) v , the fraction of pion mass carried by a valencequark, seems to be ≈ .
35 in N = 0 theory and increasesto ≈ .
45 as F π decreases to zero. Thus, even in thestrongly confined regime of 2+1 dimensional SU(2) QCD,about 30% of pion mass is carried by gluons and seaquarks, which one might want to contrast with the 3+1dimensional QCD where this fraction is about ≈ QCD . Thus, it might be that thescale-independent value of moments in 2+1 dimensionshas to be compared with PDFs in 3+1 dimensions de-termined at typical non-perturbative hadronic scales toserve as good analogues. We interpolated the data with aquadratic in F π g − , which are shown as the curves in Fig.9. It is quite striking how the individual moments them-selves weakly depend on F π . One should contrast thisbehavior with the commensurate dependence of other IRquantities to this decrease in F π . As we inferred fromthe ITD itself, (cid:104) x (cid:105) v seems to be the least sensitive tochanges in F π .It is the observables that dictate the shape of the full x -dependent PDF that are quite sensitive the infraredrather than the moments themselves. One such observ-able is the log-derivatives [43] of moments β eff ( k ) = − − ∂ log (cid:0) (cid:104) x k (cid:105) v (cid:1) /∂ log( k ) which approaches the large- x exponent β for k → ∞ . As explained in [43], we definethe discretized version of the effective β for any k as β eff ( k ) = (cid:104) x k − (cid:105) v − (cid:104) x k +2 (cid:105) v (cid:104) x k (cid:105) v k − . (33)Using k = 4, we find β eff for N = 0 , , , . , . , . , . u − d PDF, f u − d ( x ) being a positive quantityand having a probabilistic interpretation, also admits thecanonical shape observables, cumulants κ n , κ n ≡ ∂ n ∂s n log (cid:18)(cid:90) − f u − d ( x ) e sx dx (cid:19) (cid:12)(cid:12)(cid:12)(cid:12) s =0 with ,κ = (cid:104) x (cid:105) v − (cid:104) x (cid:105) v ,κ = (cid:104) x (cid:105) v − (cid:104) x (cid:105) v (cid:104) x (cid:105) v + 30 (cid:104) x (cid:105) v . (34)Using the model independent estimates of the evenmoments up to (cid:104) x (cid:105) v , we find the fourth and sixthcumulants, [ κ , κ ], for N = 0 , , , u − d PDF so as to keep the analysis fully model-independent. The aim of this exercise was to point tosome good observables of the pion PDF that seem to besensitive about the IR, and consequently, we were able todeduce simply from the model independent analysis thatthe shape of the PDF will show sharp changes as the the-ory morphs. Next, we will see these inferences concretelyarise in the reconstructed x -dependent valence PDFs.
2. Model dependent analysis: PDF reconstruction
Now we reconstruct the x -dependent valence PDFsthat best describe the real space data for ˜ M ( z , P ). Forthis we use the two-parameter functional form of the PDFin Eq. (31) that was completely sufficient to describe˜ M ( z , P ) at all N and in the range of z and P de-scribed before; in fact, when we tried to make the ansatzmore complex by adding subleading small- x terms x α +0 . and x α +1 , the fits became quite unstable and hence weresort to the simpler two-parameter ansatz. Essentially,the parametrized PDF enters through its correspondingmoments (cid:104) x k (cid:105) v ( α, β ), that is then input into the leadingtwist OPE in Eq. (30) to get the best values of α and β .In the left panels of Fig. 10, we have shown the resultingcurves for ˜ M ( z , P ) from such two-parameter fits super-imposed on the data. The quality of fits are as good asthe one from model-independent fits shown in Fig. 7.The valence PDFs, f v ( x ), corresponding to the bestfits are shown in the middle panels of Fig. 10, and therightmost panels are simply the same data replotted asthe momentum distribution, xf v ( x ). We checked that thereconstructed PDFs were robust against variations in thefit ranges by changing the maximum of the fit range from z = 6 a to z = 10 a . These variations are shown as thebands of different colors in the middle and right panelsof Fig. 10. Since the data points fall on universal curveswell to begin with, the reconstructed PDFs also show al-most no variations; so we simply take the estimate with z = 8 a for further discussions. It should be noted that6 − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 0 P = 1 . g P = 1 . g P = 2 . g . . .
50 0 . . . . f v ( x ) x N = 0 z/a ∈ (1 , z/a ∈ (1 , z/a ∈ (1 ,
10) 00 . . . . . . .
70 0 . . . . . . . . . x f v ( x ) x − . − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 2 P = 1 . g P = 1 . g P = 2 . g . . .
50 0 . . . . f v ( x ) x N = 2 z/a ∈ (1 , z/a ∈ (1 , z/a ∈ (1 ,
10) 00 . . . . . . . . .
90 0 . . . . . . . . . x f v ( x ) x − . − . − . . . . .
81 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 4 P = 1 . g P = 1 . g P = 2 . g . . .
50 0 . . . . f v ( x ) x N = 4 z/a ∈ (1 , z/a ∈ (1 , z/a ∈ (1 ,
10) 00 . . . . .
20 0 . . . . . . . . . x f v ( x ) x − − . .
51 0 1 2 3 4 5 6 7 8 ˜ M ( z , P ) z P N = 8 P = 1 . g P = 1 . g P = 2 . g . . . . f v ( x ) x N = 8 z/a ∈ (1 , z/a ∈ (1 , z/a ∈ (1 ,
10) 00 . . . .
50 0 . . . . . . . . . x f v ( x ) x FIG. 10. Reconstruction of valence PDF, f v ( x ), by fits to two-parameter ansatz. The left panels show the fits (bands) tothe bilocal matrix element ˜ M ( z , P ) (points) via leading-twist expression in Eq. (30). The middle panels show the inferredvalence PDF, f v ( x ). The different colored bands correspond to different fit ranges [0 , z ]. The right panels show xf v ( x ). Topto bottom are N = 0 , , , the scales in the different panels in Fig. 10 are differ-ent, but it is already clear that the PDFs get narroweras N increases. In terms of the exponents [ α, β ] of thePDFs, they change as [0 . , . . , . . . , . . , N = 0 , , , x exponent β for thestrongly broken phase for N = 0 and 2, are the typicalvalue around 1 and 2 as in 3+1 dimensions. The expo-7 . . . . f v ( x ) x F π g − . × − . × − . × − . × − . . .
50 0 . . . . x f v ( x ) x FIG. 11. (left) The reconstructed valence pion PDFs, f v ( x ), and (right) their corresponding momentum distribution, xf v ( x )are shown as a function of pion decay constant that characterizes the vacuum of different N = 0 , , , nent subsequently gets larger as the theory is pushed intothe near-conformal and conformal regimes. Perhaps it isof interest to note that numerically, β ≈ α + 1 for thesePDFs, which makes xf v ( x ) appear almost symmetricalaround their peak positions. As a cross-check, we alsoplot the values of moments from this analysis using PDFansatz in the top panel of Fig. 9 as the open symbols,which nicely agrees with the moments obtained from themodel-independent analysis.We summarize the PDF determination in Fig. 11 byputting together the PDFs from all N , and showing itas a function of the induced dependence on the infraredscale F π . The left and right panels show f v ( x ) and xf v ( x )respectively. The depletion of the IR scales can be seen tohave visible effect on the pion PDF. The effect of strongscale-breaking is to broaden the pion PDF over the entirerange of x ; implying indirectly, the increased importanceof gluons and the sea quarks. This is the case for F π =9 . × − g , . × − g . As the symmetry-breaking ismade about three-times weaker with F π = 3 . × − g ,we start seeing the PDF get sharper around the middlevalues x ≈ . F π ≈ (cid:104) x (cid:105) v ≈ .
44, pointing to a neardominance of the valence quarks. This extreme case canbe seen as a control in this calculation; that is, the quarkstructure of an artificial pion-like state emerging simplybecause of finite mass and volume, being not consistentwith the quark structure of an actual pion state in thescale-broken theories points to the important causal roleof the infrared vacuum structure in shaping the valencequark structure of the Nambu-Goldstone boson.
VIII. CONCLUSIONS AND DISCUSSION
We presented a lattice calculation of the valence quarkstructure of the Nambu-Goldstone boson (which we referto as the pion) of the flavor symmetry breaking in 2+1dimensional SU(2) gauge theory coupled to many mass-less flavors of fermions. The motivation for this workwas to first of all see if the quark structure of the pion issensitive to the long-distance vacuum structure, as onewould expect; and secondly to understand precisely howmuch this dependence is and in what observables thisshows up. For this work, we used N = 0 , , . g for all flavors.We studied the theories at a fixed lattice spacing andfixed finite box size. We used the pion decay constant F π as a measure of the strength of scale-breaking in theinfrared, and correlated its decrease as a function of N with other infrared quantities and to the short-distancequark structure of the pion to F π .We showed that as the strength of the infrared scalebreaking decreases, the pion Ioffe-time distribution (ITD)or bilocal quark bilinear matrix element on the light-conebecomes sensitive to this effect for Ioffe-time (or light-front distance) ν > ν <
3; the effect is seen by a slower fall-offof the ITD at ν > x -dependent valence PDF and equivalently of the u − d PDF; we demonstrated this in terms of the first8few cumulants of the u − d PDF and in terms of thelog-derivative of the moments with respect to the orderof the moment that determine the large- x behavior. Wereconstructed the valence PDF of the pion based on atwo-parameter ansatz. The above behavior of the ITDresulted in a broadening of the valence PDF over smalland large x regions when the value of the F π increased.When the F π was near zero in the near-conformal region,one could see a sharp localization of the PDF around (cid:104) x (cid:105) v ≈ .
43. The key results in this paper are shown inFig. 8 and Fig. 11.As an outcome of this work, we established the 2+1 di-mensional QCD as a good model system to perform com-putational experiments on the nonperturbative aspectsof the internal structure of hadrons using the recent de-velopments in leading-twist matching frameworks. Theshort-coming of the present work is that we do not com-pare and contrast the behavior of the pion PDF with thatof another ordinary non-Goldstone boson, such the axial-vector or the diquark states. We intend to perform thesecomparative studies in future computations, especially byusing lattices which have larger extents in the directionof boost so as to reduce the effect of Lorentz contraction(rather an expansion) of the lattice extent longitudinalto the boost. Another improvement one could do is toextend this calculation to SU(3) theory in 2+1 dimen-sions; this will extend the range of flavor N where thetheory is scale-broken, thereby making the changes tothe PDFs more gradual and easier to study than donehere. Owing to the lower-dimension used, performing anexact massless overlap fermion computation to improveon this work will be feasible. Understanding the obser-vations made in this paper in terms of simplistic modelcalculations will also shed more light on how the UV iscorrelated to the IR. With the availability of many-flavortheory computations in 3+1 dimensions (e.g., [100]), per-formed due to its relevance to composite Higgs models,it would be interesting to use them to understand theevolution of quark structures with scale depletion as N f is changed from 2 to the near-conformal point near 8 or10; especially, ask how does large- x exponent β changefor 3+1 dimensional pion?It would also be amusing to study the properties of thebilocal bilinear matrix element (Ioffe time distribution)in the long distance limit of the quark-antiquark separa-tion when the theory is in the conformal phase for N > higher-twist effects now are actually goingto be due to operators with non-trivial infrared scaling di-mensions, and thereby shed new light into the higher-spinoperators of fixed twist in the infrared CFT and its con-formal blocks, possibly corresponding to scalar-vector-vector-scalar four point function. A recent study [101]of conformal QCD in 4- (cid:15) dimensions might be helpful inthis endeavor, by carefully extrapolating the results to (cid:15) = 1.
ACKNOWLEDGMENTS
N.K. would like to thank Chris Monahan, SwagatoMukherjee, Rajamani Narayanan, Kostas Orginos, Pe-ter Petreczky, Jianwei Qiu and Raza Sufian for helpfulcomments and discussions. N.K. thanks all the membersof the BNL-SBU collaboration and the HadStruc collab-oration for fruitful discussions which greatly helped thepresent work. N.K. is supported by Jefferson ScienceAssociates, LLC under U.S. DOE Contract ). Appendix A: The behavior of P = 0 matrix element M B In Section VI B, we described the extrapolation ofthree-point function to obtain the “bare” matrix ele-ment, M B ( z , P ). The nomenclature “bare” here sim-ply means the matrix element obtained before taking theratio in Eq. (13), as there are no truly divergent behav-iors in 2+1 dimensions due to its super-renormalizability,and even for the Wilson-line, we expect it to contributeonly a benign exp (cid:8) − c (cid:48) g z (cid:9) nonperturbative higher twisteffect. In this appendix, we look at M B ( z , P = 0) it-self.From Fig. 3 in the main text, we can notice that the P = 0 matrix element with pion at rest shows a z de-pendence. To see why it is interesting, for argument-sake, . . . . . . . . . − − − P = 0 M TMC M TMC,WL M B ( z , P = ) g z N = 0 N = 2 N = 4 N = 8 FIG. 12. The matrix element M B ( z , P = 0) (before formingratios) is shown as a function of quark-antiquark separation g z in units of coupling. The different colored symbols arethe data taken from N = 0 , , , P = 0 OPE, then one simply does not expectany z dependence. In Fig. 12, we have put together the P = 0 data (shown as the points) for M B at all flavors N as a function of lattice separation, z /a . For the datashown, the Wilson-line entering the bilocal operator wassmeared using 2-steps of Stout. It is quite surprising thatthe P = 0 matrix element shows absolutely no depen-dence on the flavor or changing F π equivalently. Due tothe finite valence pion mass, there can be z dependencefrom the target mass corrections (TMC) that arises dueto trace terms at leading twist [80, 81]. We expect thisto be described by M TMC ( z ) = 1 − ( M val π z ) (cid:104) x (cid:105) v + O (cid:0) ( M val π z ) (cid:1) . (A1)To see if this arises because of the TMC, we have plottedEq. (A1) as the dashed black curve, using the value of (cid:104) x (cid:105) ≈ . z dependence should, of course, be from the Wilson-linedue to its exp (cid:0) − c (cid:48) g z (cid:1) behavior for larger | z | with some c (cid:48) . Since it should be even with respect to z → − z , we model this behavior as M TMC , WL ( z ) ≡ M TMC ( z )cosh( c (cid:48) g z ) . (A2)In Fig. 12, we plot M TMC , WL as the blue dashed curveusing c (cid:48) = 0 . c (cid:48) will be dependent onthe construction of Wilson-line itself, such as the steps ofsmearing (in our case, the value of c (cid:48) decreases to 0.206when 6-step stout was used). We see that M TMC , WL nicely describes the data at all z and for all N . Thus,through this exercise, we first understand that the non-perturbative screening behavior of Wilson-line is impor-tant in M B , and therefore, it is very important to formratios, like in 3+1 dimensions, to get rid of this trivial z dependence. In Fig. 6, we showed how this cancella-tion works well by using Wilson-lines with two-differentsmearing parameters, thereby justifying the applicationof leading twist framework to the ratio ˜ M . Secondly, dueto universal behavior of M B ( P = 0) at all N , the screen-ing behavior does not care about the infrared physics atall, pointing to the fact that it arises due to the Wilson-line self-interaction at ultraviolet scales. [1] G. Dissertori, I. G. Knowles, and M. Schmelling, Quan-tum chromodynamics: high energy experiments and the-ory , Vol. 115 (Oxford University Press, 2003).[2] C. Liu, PoS
LATTICE2016 , 006 (2017),arXiv:1612.00103 [hep-lat].[3] M. Padmanath, PoS
LATTICE2018 , 013 (2018),arXiv:1905.09651 [hep-lat].[4] R. A. Briceno, J. J. Dudek, and R. D. Young, Rev. Mod.Phys. , 025001 (2018), arXiv:1706.06223 [hep-lat].[5] A. Accardi et al. , Eur. Phys. J. A , 268 (2016),arXiv:1212.1701 [nucl-ex].[6] X.-D. Ji, Phys. Rev. Lett. , 1071 (1995), arXiv:hep-ph/9410274.[7] Y. Hatta and Y. Zhao, Phys. Rev. D , 034004 (2020),arXiv:2006.02798 [hep-ph].[8] C. Alexandrou, M. Constantinou, K. Hadjiyiannakou,K. Jansen, C. Kallidonis, G. Koutsou, A. VaqueroAvil´es-Casco, and C. Wiese, Phys. Rev. Lett. ,142002 (2017), arXiv:1706.02973 [hep-lat].[9] C. Alexandrou, S. Bacchio, M. Constantinou, J. Finken-rath, K. Hadjiyiannakou, K. Jansen, G. Koutsou,H. Panagopoulos, and G. Spanoudes, Phys. Rev. D ,094513 (2020), arXiv:2003.08486 [hep-lat].[10] J. Badier et al. (NA3), Z. Phys. C , 281 (1983).[11] B. Betev et al. (NA10), Z. Phys. C , 9 (1985).[12] J. Conway et al. , Phys. Rev. D , 92 (1989).[13] J. Owens, Phys. Rev. D , 943 (1984).[14] P. Sutton, A. D. Martin, R. Roberts, and W. Stirling,Phys. Rev. D , 2349 (1992).[15] M. Gluck, E. Reya, and A. Vogt, Z. Phys. C , 651(1992).[16] M. Gluck, E. Reya, and I. Schienbein, Eur. Phys. J. C , 313 (1999), arXiv:hep-ph/9903288. [17] K. Wijesooriya, P. Reimer, and R. Holt, Phys. Rev. C , 065203 (2005), arXiv:nucl-ex/0509012.[18] M. Aicher, A. Schafer, and W. Vogelsang, Phys. Rev.Lett. , 252003 (2010), arXiv:1009.2481 [hep-ph].[19] S. J. Brodsky and G. R. Farrar, Phys. Rev. Lett. ,1153 (1973).[20] T. Nguyen, A. Bashir, C. D. Roberts, and P. C. Tandy,Phys. Rev. C , 062201 (2011), arXiv:1102.2448 [nucl-th].[21] C. Chen, L. Chang, C. D. Roberts, S. Wan,and H.-S. Zong, Phys. Rev. D , 074021 (2016),arXiv:1602.01502 [nucl-th].[22] Z.-F. Cui, M. Ding, F. Gao, K. Raya, D. Binosi,L. Chang, C. D. Roberts, J. Rodr´ıguez-Quintero, andS. M. Schmidt, Eur. Phys. J. C , 1064 (2020).[23] C. D. Roberts and S. M. Schmidt (2020)arXiv:2006.08782 [hep-ph].[24] G. F. de Teramond, T. Liu, R. S. Sufian, H. G. Dosch,S. J. Brodsky, and A. Deur (HLFHS), Phys. Rev. Lett. , 182001 (2018), arXiv:1801.09154 [hep-ph].[25] E. Ruiz Arriola, Acta Phys. Polon. B , 4443 (2002),arXiv:hep-ph/0210007.[26] W. Broniowski and E. Ruiz Arriola, Phys. Lett. B ,385 (2017), arXiv:1707.09588 [hep-ph].[27] J. Lan, C. Mondal, S. Jia, X. Zhao, and J. P. Vary,Phys. Rev. D , 034024 (2020), arXiv:1907.01509[nucl-th].[28] P. Barry, N. Sato, W. Melnitchouk, and C.-R. Ji, Phys.Rev. Lett. , 152001 (2018), arXiv:1804.01965 [hep-ph].[29] I. Novikov et al. , Phys. Rev. D , 014040 (2020),arXiv:2002.02902 [hep-ph].[30] K. D. Bednar, I. C. Clo¨et, and P. C. Tandy, Phys. Rev. Lett. , 042002 (2020), arXiv:1811.12310 [nucl-th].[31] X. Ji, Phys. Rev. Lett. , 262002 (2013),arXiv:1305.1539 [hep-ph].[32] X. Ji, Sci. China Phys. Mech. Astron. , 1407 (2014),arXiv:1404.6680 [hep-ph].[33] A. Radyushkin, Phys. Rev. D , 034025 (2017),arXiv:1705.01488 [hep-ph].[34] K. Orginos, A. Radyushkin, J. Karpie, andS. Zafeiropoulos, Phys. Rev. D , 094503 (2017),arXiv:1706.05373 [hep-ph].[35] V. Braun and D. M¨uller, Eur. Phys. J. C , 349 (2008),arXiv:0709.1348 [hep-ph].[36] Y.-Q. Ma and J.-W. Qiu, Phys. Rev. D , 074021(2018), arXiv:1404.6860 [hep-ph].[37] Y.-Q. Ma and J.-W. Qiu, Phys. Rev. Lett. , 022003(2018), arXiv:1709.03018 [hep-ph].[38] M. Constantinou, in (2020) arXiv:2010.02445 [hep-lat].[39] Y. Zhao, Int. J. Mod. Phys. A , 1830033 (2019),arXiv:1812.07192 [hep-ph].[40] K. Cichy and M. Constantinou, Adv. High Energy Phys. , 3036904 (2019), arXiv:1811.07248 [hep-lat].[41] C. Monahan, PoS LATTICE2018 , 018 (2018),arXiv:1811.00678 [hep-lat].[42] X. Ji, Y.-S. Liu, Y. Liu, J.-H. Zhang, and Y. Zhao,(2020), arXiv:2004.03543 [hep-ph].[43] X. Gao, L. Jin, C. Kallidonis, N. Karthik, S. Mukher-jee, P. Petreczky, C. Shugert, S. Syritsyn, and Y. Zhao,Phys. Rev. D , 094513 (2020), arXiv:2007.06590[hep-lat].[44] J.-H. Zhang, J.-W. Chen, L. Jin, H.-W. Lin, A. Sch¨afer,and Y. Zhao, Phys. Rev. D , 034505 (2019),arXiv:1804.01483 [hep-lat].[45] T. Izubuchi, L. Jin, C. Kallidonis, N. Karthik,S. Mukherjee, P. Petreczky, C. Shugert, and S. Syrit-syn, Phys. Rev. D , 034516 (2019), arXiv:1905.06349[hep-lat].[46] B. Jo´o, J. Karpie, K. Orginos, A. V. Radyushkin, D. G.Richards, R. S. Sufian, and S. Zafeiropoulos, Phys. Rev.D , 114512 (2019), arXiv:1909.08517 [hep-lat].[47] H.-W. Lin, J.-W. Chen, Z. Fan, J.-H. Zhang, andR. Zhang, (2020), arXiv:2003.14128 [hep-lat].[48] R. S. Sufian, J. Karpie, C. Egerer, K. Orginos, J.-W.Qiu, and D. G. Richards, Phys. Rev. D , 074507(2019), arXiv:1901.03921 [hep-lat].[49] R. S. Sufian, C. Egerer, J. Karpie, R. G. Edwards,B. Jo´o, Y.-Q. Ma, K. Orginos, J.-W. Qiu, andD. G. Richards, Phys. Rev. D , 054508 (2020),arXiv:2001.04960 [hep-lat].[50] A. C. Aguilar et al. , Eur. Phys. J. A , 190 (2019),arXiv:1907.08218 [nucl-ex].[51] B. Adams et al. , (2018), arXiv:1808.00848 [hep-ex].[52] N. Karthik and R. Narayanan, Phys. Rev. D , 054510(2018), arXiv:1801.02637 [hep-th].[53] N. Seiberg, T. Senthil, C. Wang, and E. Witten, AnnalsPhys. , 395 (2016), arXiv:1606.01989 [hep-th].[54] A. Karch and D. Tong, Phys. Rev. X , 031043 (2016),arXiv:1606.01893 [hep-th].[55] C. Choi, D. Delmastro, J. Gomis, and Z. Komargodski,JHEP , 078 (2020), arXiv:1810.07720 [hep-th].[56] A. Sharon, JHEP , 078 (2018), arXiv:1803.06983[hep-th].[57] S. Gazit, F. F. Assaad, S. Sachdev, A. Vishwanath,and C. Wang, Proc. Nat. Acad. Sci. , E6987 (2018), arXiv:1804.01095 [cond-mat.str-el].[58] X.-Y. Song, Y.-C. He, A. Vishwanath, and C. Wang,Phys. Rev. X , 011033 (2020), arXiv:1811.11182[cond-mat.str-el].[59] R. Ma and Y.-C. He, Phys. Rev. Res. , 033348 (2020),arXiv:2003.05954 [cond-mat.str-el].[60] S. M. Chester and S. S. Pufu, JHEP , 069 (2016),arXiv:1603.05582 [hep-th].[61] S. M. Chester and S. S. Pufu, JHEP , 019 (2016),arXiv:1601.03476 [hep-th].[62] S. S. Pufu, Phys. Rev. D , 065016 (2014),arXiv:1303.6125 [hep-th].[63] X. Ji, Y. Liu, and I. Zahed, Phys. Rev. D , 054008(2019), arXiv:1807.07528 [hep-ph].[64] Y. Jia, S. Liang, X. Xiong, and R. Yu, Phys. Rev. D , 054011 (2018), arXiv:1804.04644 [hep-th].[65] L. Del Debbio, T. Giani, and C. J. Monahan, JHEP , 021 (2020), arXiv:2007.02131 [hep-lat].[66] A. Kock, Y. Liu, and I. Zahed, Phys. Rev. D ,014039 (2020), arXiv:2004.01595 [hep-ph].[67] U. Magnea, Phys. Rev. D , 056005 (2000), arXiv:hep-th/9907096.[68] R. D. Pisarski, Phys. Rev. D , 2423 (1984).[69] C. Vafa and E. Witten, Commun. Math. Phys. , 257(1984).[70] T. Appelquist and D. Nash, Phys. Rev. Lett. , 721(1990).[71] H. Goldman and M. Mulligan, Phys. Rev. D , 065031(2016), arXiv:1606.07067 [hep-th].[72] A. Niemi and G. Semenoff, Phys. Rev. Lett. , 2077(1983).[73] A. Redlich, Phys. Rev. D , 2366 (1984).[74] A. Thomson and S. Sachdev, Phys. Rev. B , 205128(2017), arXiv:1607.05279 [cond-mat.str-el].[75] L. Del Debbio and R. Zwicky, Phys. Rev. D , 014502(2010), arXiv:1005.2371 [hep-ph].[76] K. Rummukainen and S. A. Gottlieb, Nucl. Phys. B , 397 (1995), arXiv:hep-lat/9503028.[77] V. Braun, P. Gornicki, and L. Mankiewicz, Phys. Rev.D , 6036 (1995), arXiv:hep-ph/9410318.[78] G. Martinelli and C. T. Sachrajda, Phys. Lett. B ,184 (1987).[79] T. Izubuchi, X. Ji, L. Jin, I. W. Stewart, and Y. Zhao,Phys. Rev. D , 056004 (2018), arXiv:1801.03917 [hep-ph].[80] J.-W. Chen, S. D. Cohen, X. Ji, H.-W. Lin, and J.-H.Zhang, Nucl. Phys. B , 246 (2016), arXiv:1603.06664[hep-ph].[81] A. Radyushkin, Phys. Lett. B , 514 (2017),arXiv:1702.01726 [hep-ph].[82] V. Dotsenko and S. Vergeles, Nucl. Phys. B , 527(1980).[83] T. Ishikawa, Y.-Q. Ma, J.-W. Qiu, and S. Yoshida,Phys. Rev. D , 094019 (2017), arXiv:1707.03107 [hep-ph].[84] X. Ji, J.-H. Zhang, and Y. Zhao, Phys. Rev. Lett. ,112001 (2018), arXiv:1706.08962 [hep-ph].[85] M. J. Teper, Phys. Rev. D , 014512 (1999), arXiv:hep-lat/9804008.[86] N. Karthik and R. Narayanan, Phys. Rev. D , 045020(2016), arXiv:1512.02993 [hep-lat].[87] N. Karthik and R. Narayanan, Phys. Rev. D , 025003(2015), arXiv:1505.01051 [hep-th].[88] C. Morningstar and M. J. Peardon, Phys. Rev. D , , 028(2006), arXiv:hep-lat/0607006.[90] S. Gupta and N. Karthik, Phys. Rev. D , 094001(2013), arXiv:1302.4917 [hep-lat].[91] S. Duane, A. Kennedy, B. Pendleton, and D. Roweth,Phys. Lett. B , 216 (1987).[92] I. P. Omelyan, I. M. Mryglod, and R. Folk, Phys. Rev.E , 056706 (2002).[93] N. Karthik, (2014), arXiv:1401.1072 [hep-lat].[94] S. Gusken, U. Low, K. Mutter, R. Sommer, A. Patel,and K. Schilling, Phys. Lett. B , 266 (1989).[95] G. S. Bali, B. Lang, B. U. Musch, and A. Sch¨afer, Phys.Rev. D , 094515 (2016), arXiv:1602.05525 [hep-lat]. [96] N. Karthik and R. Narayanan, Phys. Rev. D , 065026(2016), arXiv:1606.04109 [hep-th].[97] E. V. Mastropas and D. G. Richards (Hadron Spec-trum), Phys. Rev. D , 014511 (2014), arXiv:1403.5575[hep-lat].[98] J. Verbaarschot and I. Zahed, Phys. Rev. Lett. , 2288(1994), arXiv:hep-th/9405005.[99] J. Karpie, K. Orginos, and S. Zafeiropoulos, JHEP ,178 (2018), arXiv:1807.10933 [hep-lat].[100] T. Appelquist et al.et al.