Parton distribution function for the gluon condensate

Motivated by the desire to understand the nucleon mass structure in terms of light-cone distributions, we introduce the twist-four parton distribution function $F(x)$ whose first moment is the gluon condensate in the nucleon. We present the equation of motion relations for $F(x)$ and discuss the possible existence of the delta function (`zero mode') contribution at $x=0$. We also perform one-loop calculations for quark and gluon targets.


I. INTRODUCTION
The hadronic matrix element of the dimension-four scalar gluonic operator, or the 'gluon condensate' is fundamentally important in hadron physics and beyond. This is primarily because the trace anomaly of QCD imparts mass to the nucleons and nuclei, hence to the visible universe, through the matrix element in Eq. (1) [1]. It thus plays a pivotal role in understanding the origin of the nucleon mass, a problem recently proclaimed by the National Academy of Science [2] as one of the main scientific goals of the future Electron-Ion Collider (EIC) [3]. However, the precise determination of Eq. (1) turns out to be an extremely challenging task. A direct calculation from lattice QCD is notoriously difficult due to the vacuum quantum numbers of the operator involved (see a recent attempt [4]). Another possibility is that the matrix element can be probed experimentally in near-threshold quarkonium production [5][6][7].
In this paper, we propose to study the partonic structure of the gluon condensate in Eq. (1) as a novel direction in the research of nucleon mass structure. Since this is a rather unusual proposal, to motivate the reader let us first draw an analogy to the study of nucleon spin structure. The Jaffe-Manohar sum rule [1] tells how the total nucleon spin of 1/2 is distributed among the helicity ∆Σ, ∆G and orbital angular momentum L q,g of quarks and gluons. Each of these components can be expressed by the first moment of the corresponding parton distribution ∆Σ = dx∆q(x), L q = dxL q (x), etc [8], where x is the longitudinal momentum fraction. Such distributions are not only useful for extracting the moments from experiments, but also interesting in their own right, as they provide a more detailed, higher-dimensional description on the spin structure.
Returning to the problem of mass, similarly to Eq. (2), one can decompose the nucleon mass M as [9] where M q,g kin are the kinetic energies carried by quarks and gluons, M m ∼ P |ψψ|P is the contribution from the nucleon sigma term, and M a is from the gluon condensate Eq. (1). As in the case of spin decomposition, one naturally asks how partons with a given momentum fraction x contribute to the four components in Eq. (3). For kinetic energy, this can be quantified by noticing that M q,g kin are related to the second moment of the ordinary parton distribution functions (PDFs). One then sees that M q kin is dominated by valence quarks at large-x. Gluons tend to have smaller x values, but because there are so many of them, M g kin can become sizable. On the other hand, regarding the remaining two entries M m,a , most of the work done so far has been limited to 'zero-dimensional' physics. While the parton distribution function for M m does exist in the literature, called e(x), its connection to hadron masses is not often emphasized. For M a , the corresponding x-distribution was almost nonexistent until very recently when related distributions were briefly mentioned in [15]. In principle, it is a simple matter to write down the twist-four distribution Together with the twist-two PDFs and e(x), this provides a complete set of parton distributions for the nucleon mass structure.
In this paper, we present the first analysis of F (x). We use the QCD equation of motion to reveal its multi-partonic nature of the distribution. We then present one-loop calculations of F (x) for quark and gluon targets. Particular attention is given to the question of whether F (x) contains the delta function δ(x). The (non)existence of δ(x) in e(x) has been a subject of debate in the literature. We shall see that the discussion is entirely analogous for F (x). We shall present both model-independent and model-dependent arguments in favor of the existence of the delta function.
Before introducing the twist-four gluon distribution F (x), we first give a review of the twist-three, chiral-odd quark distribution e(x). Our purpose is mostly to emphasize the similarity to F (x) studied in the next section, but the present section also contains some original discussions.
e(x) is defined by where M is the proton mass and W is the straight Wilson line along the light-cone which makes the nonlocal operator gauge invariant. The distribution is defined for each quark flavor q with mass m. The first and second moments are proportional to the nucleon sigma term and the number of valence quarks N q , respectively In what follows, we shall omit the subscript q for simplicity. By using the equation of motion and Lorentz invariant relation one can write [10][11][12][13] e(x) = e sing (x) + e tw3 (x) + e mass (x), where e sing is proportional to the delta function at x = 0, e mass (x) is related to the twist-two quark distribution q(x) as Clearly, dx e mass (x) = 0. The 'genuine twist-three' distribution e tw3 (x) also has a delta function at x = 0, where is the quark-gluon-quark mixed distribution. [Our sign convention for the QCD coupling is such that the covariant derivative reads D µ = ∂ µ + igA µ .] It is easy to see that dx e tw3 (x) = 0, dx xe tw3 (x) = 0.
This matrix element is related to the electric dipole moment of the nucleon [14]. One thus arrives at the relation There have been discussions about the nature of the delta function terms, or 'zero modes', in Eq. (14). Ref. [12] argues that the sum rule dxe(x) ∝ P |ψψ|P is of 'no practical use' because the only contribution comes from the delta function at x = 0 which experiments cannot access. [Remember that dx e mass (x) = dx e tw3 (x) = 0.] The presence of 'zero modes' signifies the nonperturbative dynamics of QCD which leads to confinement and the generation of hadron masses. On the other hand, one can make an argument that the delta function may actually be absent. This is indeed the case in the naive parton model owing to the Weisberger relation [16] which in the modern notation reads [17] Since the genuine twist-three physics is absent in the parton model, the expression inside the square brackets in Eq. (14) vanishes. However, Eq. (16) is obviously problematic because the x-integral does not converge in real QCD.
Going beyond the parton model, very recently the authors of [18] claim to have shown that the coefficient of the delta function vanishes exactly in full QCD. Their proof starts by writingψψ =ψ + ψ − +ψ − ψ + where ψ ± = 1 2 γ ∓ γ ± ψ are the so-called 'good' and 'bad' components of the quark field, respectively. It is often stated in the literature that ψ − is not an independent field. Using the equation of motion one can write where i = 1, 2. The general solution to Eq. (17) is where K(z − ) is the Green function subject to the boundary condition. Common choices are x−i and P 1 x , respectively. (P denotes the principal value.) ψ 0 − is not constrained by the equation of motion and should be treated as an independent field. It is essentially the zero mode as it involves an unconstrained integration over z − (up to a gauge rotation). In the literature, this term is routinely neglected when one works in the light-cone gauge A + = 0 and specifies the boundary condition at z − = ±∞ in order to quantize the theory. Often the antisymmetric boundary condition, corresponding to P 1 x , is employed (see, e.g., [19]), but this implicitly assumes the subtraction of the zero mode. While such a procedure may be justified for most purposes, like doing perturbation theory and computing the S-matrix, it may not capture the long-distance physics responsible for the generation of hadron mass.
Ref. [18] only kept the first term of Eq. (18) with the advanced boundary condition K(z − ) = −θ(−z − ) and showed that the coefficient of the delta function in Eq. (14) vanishes exactly. Actually, it does not matter which boundary condition is adopted, because in the end only the combination K(x) + K * (x) = 2P 1 x appears in the sum 1 Here is a quick derivation of the Weisberger relation in the parton model.
The factor 1/x comes from the relativistic normalization of states.
However, the ψ 0 − term does not cancel and leads to a nonvanishing coefficient There is vast literature on the zero mode problem in light-front quantization (see, e.g., [20] and reviews [21,22]). One might argue that in continuum theory the zero mode has no effect on physical observables because it has measure zero in the path integral sense. On the other hand, entirely neglecting the zero mode causes serious inconsistencies such as the lack of Lorentz invariance [20]. This is still an open problem, and discussions of the quark and gluon condensates cannot be complete without a full consideration of the zero mode. For the moment, it seems to us that the coefficient of the delta function is likely nonvanishing, and can be determined only nonperturbatively possibly along the line recently suggested in [15]. A.
e(x) to one-loop In Ref. [23], the authors have shown in the massive quark model to one-loop that e(x) indeed contains the delta function δ(x). This is consistent with the above observation that the delta function is nonvanishing in general. In the massive quark model where |p is a single quark state, it is appropriate to employ the scale invariant mass for the 'hadron' mass M in Eq. (5), where where µ is the renormalization scale. As observed in [23], without the delta function the sum rule cannot be satisfied. Eq. (21) is derived from the following one-loop integral in the light-cone gauge n · A = A + = 0 in where x = k + /p + . [We use the same letter for the small dimension in dimensional regularization and in the i prescription of the propagator, but the distinction should be obvious.] The first term in the numerator is proportional to which is the origin of the delta function δ(x) in Eq. (21).
Let us consider the same matrix element but now |p is an on-shell gluon p 2 = 0 with transverse polarization ε · p = ε · n = 0 and ε · ε * = −ε ⊥ · ε * ⊥ = −1. The one-loop diagrams give where T F = 1/2 and Θ(1 > x > 0) denotes a step function which has support on 1 > x > 0. The delta function δ(x) arises from the same integral in Eq. (24). Integrating over x, we get zero. This is consistent with the fact that the local operatorψψ does not mix with F µν F µν , and the delta function δ(x) is crucial to ensure this property. We also see that the mixing does occur at the level of the x-distributions.

III. GLUON CONDENSATE DISTRIBUTION
Let us now come to the main object of interest. With the motivation stated in the introduction, we consider the twist-four distribution The second moment vanishes dx xF (x) = 0 because F (x) is an even function in x. Similarly to e(x), and as conjectured in [15], we expect that F (x) also has a delta function piece To obtain insights into the structure of F (x), consider the following operator relation where we used the Bianchi identity and D µ represents the translation operator: O(0, z)). Further using the equation of motion, we immediately obtain We shall interpret 1 x as the principal value P 1 x to be consistent with the property F (x) = F (−x). Notice that because gψγ ν F ν+ ψ = −∂ ν T ν+ q is a total derivative. (T µν q is the quark part of the energy momentum tensor.) Thus the coefficient of the delta function is However, the recent work [18] suggests that C may actually be zero, or at least there is a significant cancellation among the three terms in C. From the equation of motion one can formally write Therefore, The first term on the right hand side can be written as, after taking the forward matrix element P |...|P and using translational symmetry, and similarly for F +i . The zero modes F 0 +− , F 0 +i are not constrained by the equation of motion and should be regarded as independent degrees of freedom. We thus expect that, in general, the cancellation is incomplete and there exists a delta function δ(x) in F (x).

IV. ONE-LOOP COMPUTATION OF F (x)
In order to gain insight into the x-dependence of F (x), in this section we perform one-loop calculations for quark and gluon targets. We shall be particularly interested in whether F (x) contains the delta function δ(x) or not.

A. Quark target
We use the light-cone gauge n · A = A + = 0 to eliminate the Wilson line. The gluon propagator is proportional to the tensor We specify the prescription for the pole 1/k · n when need arises. For an on-shell quark external state p 2 = m 2 , we find where k + = xp + . Note that the pole 1/n · k has canceled between the two diagrams. The first term on the last line of Eq. (39) is nonvanishing when 1 > x > 0 and can be evaluated in a standard manner. The second term is proportional to the delta function at x = 1, The third term vanishes. We thus find, for 1 ≥ x ≥ −1, and consequently, Eq. (42) is the expected result consistent with the known operator relation where the left hand side is the bare operator. β 0 = 11Nc and γ m0 = 3C F 2π is the first coefficient of the mass anomalous dimension γ m = γ m0 α s + · · · . Our result gives an interesting new perspective on this well-known result in Eq. (43). The one-loop anomalous dimension γ m0 originates from the delta function spikes at x = ±1 (meaning that the gluon carries away all the quark's energy) and an almost flat distribution for 1 > x > −1. Curiously, the delta function δ(x) is absent, in contrast to e(x) in the same model. In the next subsection we perform the same analysis for the coefficient of F 2 in Eq. (43).

B. Gluon target
Next we consider the case where the target |p is a single gluon with transverse polarization. To regularize the infrared divergence, the gluon is assumed to be off-shell with spacelike momentum p 2 = 2p + p − < 0. Accordingly, we take M 2 = −p 2 . To zeroth order To one-loop, the diagrams which give nonvanishing contributions are listed in Fig. 1. There are also the self-energy diagrams to be considered later. After straightforward calculations we find, for 1 > x = k + /p + > 0, Fig. 1(a)+(b): Fig. 1(c): where x = k + /p + . The result for x < 0 is simply obtained by x → −x, k µ → −k µ . At this point we must specify the prescription for the spurious poles 1/k + ∼ 1/x and 1/(p + − k + ) ∼ 1/(1 − x). If one uses the principal value (pv) prescription the k − integral does not interfere with the poles. Then the terms proportional to (p − k) 2 and k 2 in the numerator can be dropped. However, the remaining integrals contain frame-dependent divergences ∼ ln p + /δ whose cancellation is nontrivial. This is a well-known symptom of the principal value prescription. Here we instead adopt the Mandelstam-Leibbrandt (ML) prescription [25], With this choice, one can write and use the master integrals collected in Appendix B. The result for the total contribution from the three diagrams is where the plus-prescription 1/[x] + is defined as We thus see that, similarly to e(x), F (x) also contains the delta function at x = 0. However, the way it appears is somewhat unexpected. The coefficient of δ(x) is divergent, and its only role is to cancel the familiar soft gluon singularity 1/x in the first moment. This is a potentially important observation that may find other applications. Notice that the x-integral of Eq. (50) vanishes exactly including the finite terms 1 0 dx ((a) + (b) + (c)) = 0. This is a special feature of the ML prescription which is not shared by the principal value prescription. It actually agrees with the result obtained in the background field gauge [24] (for the divergent part), namely, the renormalization of the local operator F µν F µν solely comes from the self-energy insertion into the external legs. However, in the x-space we find an interesting redistribution of partons. The finite part (obtained after removing the 1/ pole and setting µ 2 = −p 2 ) is plotted in Fig. 2. The density of F 2 is negative in the large-x region 1 > x > 0.3, and this depletion is exactly compensated by the positive region at small-x and the delta functions at x = 1. x , see Eq. (50), for 1 > x > 0.
The self-energy diagrams modify the leading term as, again in the ML prescription [26], Adding all contributions, we arrive at, for 1 > x > 0, The result for 0 > x > −1 is simply given by F (x) = F (−x). The first moment reads in agreement with Eq. (43). Incidentally, the n-th moment is given by, for even n ≥ 2, It is tempting to relate this result to the anomalous dimension of the operator F µν (D + ) n F µν . However, this is nontrivial because for high-dimension operators one has to compute multi-point Greens' function, not just the twopoint function, in order to disentangle the mixing with other operators. A proper treatment in the case of e(x) has been given in [23]. Yet, very little is known about the anomalous dimension of high-dimensional, higher-twist gluonic operators [27,28]. We leave this to future work.

V. CONCLUSIONS
In this paper we have introduced the twist-four parton distribution function F (x) which integrates to the gluon condensate P |F 2 |P and studied its properties based on the equation of motion relations and one-loop calculations. Our work literally adds a new dimension-momentum fraction x-to the study of nucleon mass structure. In the future, it would be interesting to further include the dependence on the transverse momentum F (x, k ⊥ ) as was done for the quark distribution e(x, k ⊥ ) (see e.g., Ref. [13]). However, at the moment, all this is highly formal and mostly of conceptual interest. The first moment P |F 2 |P can be probed in near-threshold quarkonium production [5][6][7], but identifying experimental processes that are sensitive to the x-dependence will be more challenging.
Both the operator analysis and one-loop calculations suggest that F (x) contains the delta function δ(x). After all, this is physically reasonable and could have been anticipated since the zero mode x = 0 is the genuine nonperturbative sector of light-front quantization [20][21][22], and therefore it has to do with the generation of hadron masses. In perturbation theory, there is of course no issue of mass generation. Still, the delta function is necessary for the consistency of the calculation, like reproducing the correct anomalous dimension as we have seen and restoring Lorentz invariance as emphasized elsewhere (see, e.g., [29] for a recent discussion). Finally, we emphasize that the structure at finite 1 > x > 0 is equally interesting and has a better chance to be explored either experimentally or in lattice QCD. In particular, we predict the enhancement at small-x due to the familiar soft gluon divergence in QCD. It would be interesting to study higher order corrections to this behavior (for example along the line of [30,31]) and also the possible impact of the gluon saturation.