1 Introduction

Amongst all the hadrons that appear in Nature’s catalogue of bound states, pions would seem to be the simplest. There are three such states (\(\pi ^\pm \), \(\pi ^0\)) and within the Standard Model (SM) they are described as systems built from a valence-quark and a valence antiquark, e.g., as illustrated in Fig. 1, the positively charged pion, \(\pi ^+\), is a \(u {\bar{d}}\) state: one valence u-quark and one valence \({\bar{d}}\)-quark. There are complications, however. In the absence of Higgs boson couplings into quantum chromodynamics (QCD) – the SM’s strong interaction part, pions are massless Nambu-Goldstone bosons. Even with Higgs couplings restored, u, d quarks remain light; and so do pions [1], possessing masses, \(m_{\pi ^\pm }\approx m_{\pi ^0}\), far less than that of the proton, \(m_p\). Amongst other things, these features are critical to the stability of all nuclei [2] and, therefore, to the emergence of the known Universe.

Fig. 1
figure 1

Referring to QCD’s Lagrangian degrees-of-freedom, the \(\pi ^+\) is built from one valence u-quark, one valence \({\bar{d}}\)-quark, and, owing to the properties of QCD, infinitely many gluons and sea quarks, drawn here as “springs” and closed loops, respectively. The \(\pi ^-\) is \(d{\bar{u}}\) and the \(\pi ^0\) is \(u {\bar{u}} - d{\bar{d}}\)

Full size image

Indeed, Nature has no more fundamental Nambu-Goldstone boson than the pion. If it is a QCD bound state, then it must have a complex structure – Fig. 1, because although QCD has not been solved, many remarkable properties have been exposed. For instance, gluons, which appear as massless gauge-boson degrees-of-freedom in the QCD Lagrangian, acquire a momentum dependent mass function, \(m_g(k^2)\), whose value on \(k^2\simeq 0\) is characterised by a renormalisation group invariant mass \(m_0 \approx m_p/2\) [3,4,5,6,7,8,9,10]. This is a primary sign of the dynamical violation of scale invariance in QCD [11], whose origin is strong gluon self-interactions, and the consequent emergence of hadron mass (EHM); namely, the appearance of a nuclear-size mass-scale in a theory that, a priori, contains no masses at all.

In connection with the Nambu–Goldstone pions, two key corollaries of EHM are the appearance of a mass-function for quarks, whose value on \(k^2\simeq 0\) is roughly \(m_p/3\) [12], and the intimate relationship [13,14,15,16,17,18] between that mass function and the pion’s Poincaré-covariant bound-state wave function, \(\chi _\pi (k,k-P)\), where P is the pion’s total momentum and k is the momentum of the valence-quark. Combined, these two properties entail that although it can be exposed using a wide variety of means [19,20,21,22], the cleanest observable manifestations of EHM are to be found in pion properties; so imparting utmost importance to the goals of revealing and explaining pion structure [23,24,25,26,27,28]. In this connection, much can be learnt from data that may be used to extract pion parton distribution functions (DFs).

Section 2 explains QCD constraints on the behaviour of the pion valence-quark DF at large values of the light-front longitudinal momentum fraction, x; and Sects. 3, 4 review aspects of existing experiments [29,30,31,32] and associated phenomenological fits [33,34,35,36] which relate to the determination of this behaviour. Section 5 presents a novel perspective on the data in Ref. [32, E615], referring it to the character of the quark+antiquark interaction; and also describes the importance and expression of key QCD global symmetries in calculations of , along with a natural definition for the hadron scale, \(\zeta _\mathcal{H}\). The principle and practical implications of all-orders evolution are detailed in Sect. 6, enabling an assessment of contemporary fits to data and the methods they employ in Sects. 7, 8. Section 9 provides a summary and an outlook.

2 Pion valence-quark distribution at large \({\mathbf {x}}\)

In quantum mechanics, a bound-state’s Schrödinger wave function contains a great deal of information about the system. Owing to the complexities of quantum field theory, the nearest analogue for the pion is its light-front wave function (LFWF) [37, 38], \(\psi _\pi (x,\mathbf {k}_\perp ;P)\), which can be obtained via a light-front projection of \(\chi _\pi (k,k-P)\) [39, 40] and permits interpretation as a probability amplitude, unlike the Bethe-Salpeter wave function. Here, using linearly independent light-like four-vectors n, \({\bar{n}}\), with \(n^2=0={\bar{n}}^2\), \(n\cdot {\bar{n}}=-1\): \(x=n\cdot k/n\cdot P\), i.e., the light-front fraction of the pion’s total momentum carried by the valence-quark; and \(k_\perp \) is the two-vector built from the in-general nonzero components of \(k_\mu = O^\perp _{\mu \nu }k_\nu \), \( O^\perp _{\mu \nu } = \delta _{\mu \nu } + n_\mu {\bar{n}}_\nu + {\bar{n}}_\mu n_\nu \), viz. that part of the valence-quark’s momentum which lies in the light-front transverse plane. When discussing \(\psi _\pi \), it is common, although not mandatory, to consider its expansion over an infinite dimensional Fock space of basis states built using non-interacting gluon and quark partons [41]. Then the expansion coefficient associated with any particular configuration of gluon and quark partons is the probability amplitude for finding the hadron with the identified partonic configuration carrying and sharing all its properties.

A hadron’s LFWF provides direct access to one of the first quantities used to reveal the presence of quarks within hadrons [42,43,44,45]; namely, the hadron’s parton DFs [46]. As an example, describes the probability density for finding a valence u-quark with light-front momentum fraction x when the pion is resolved at scale \(\zeta \). Notably, at any finite scale, this valence-quark will not be equivalent to a valence-quark-parton; instead, the valence-quark will be linked to the quark-parton as an object dressed by QCD interactions in the manner described by the quark gap equation [47].

Using contemporary tools:

(1)

where is the valence u-quark generalised parton distribution (GPD), which, on the kinematic domain of interest herein (zero skewness), may be expressed in terms of the pion’s valence-quark LFWF [48]:

(2)

where \(k_{\perp \pm }=k_\perp \pm (1-x)t/2\), with t the squared momentum-transfer to the pion in processes designed to measure the GPD. In the \(\mathcal{G}\)-parity symmetry limit, which is an accurate reflection of Nature, .

Since calculation of any bound-state wave function is a nonperturbative problem, perturbative QCD (pQCD) cannot be used to calculate the x-dependence of . However, it can be employed to predict the behaviour on \(x\simeq 1\). Owing to implicit kinematic constraints on this domain, only the leading two-body Fock-space component of the pion’s complete LFWF is sampled, enabling simplifications that lead to the result

(3)

where \(\gamma (\zeta ) \ge 0\) grows logarithmically with \(\zeta \), expressing the impact of gluon radiation from the struck quark that is realised in QCD evolution equations (DGLAP) [49,50,51,52]. The hadron scale is that value of \(\zeta =\zeta _\mathcal{H}\) at which begins gluon emission from the valence quarks in the complete wave function [53]: \(\gamma (\zeta _\mathcal{H})=0\).

Equation (3) has been derived in various ways, e.g., Refs.  [53,54,55,56,57]; and as we now indicate, it can also be obtained straightforwardly from Eq. (2). The value of the integral in Eq. (2) is determined by the large-\(k_\perp ^2\) behaviour of \(\psi _{\pi }^{u}\left( x,{k}_\perp ^2;\zeta \right) \), which may be obtained using the operator product expansion [58]:

$$\begin{aligned} \psi _{\pi }^{u}\left( x,{k}_\perp ^2;\zeta \right) {\mathop {=}\limits ^{k_\perp ^2 \gg m_0^2}} d_0(\zeta ) \, x(1-x)\frac{\alpha (k_\perp ^2)}{k_\perp ^2} + \ldots \,, \end{aligned}$$
(4)

where \( d_0(\zeta )\) measures the strength of the leading term at resolving scale \(\zeta \), \(\alpha (k_\perp ^2)\) is the QCD running coupling, and the ellipsis indicates omitted terms that are suppressed by additional powers of \(1/\ln (k_\perp ^2/\varLambda ^2_{\mathrm{QCD}})\), where \(\varLambda _{\mathrm{QCD}}\) is a renormalisation group invariant mass scale whose size is an expression of EHM. Inserting Eq. (4) into Eq. (2) and recognising that no infrared divergence can occur because all (regularised or renormalised) Fock-space elements are independently normalisable, then Eq. (3) is recovered via Eq. (1).

It is now evident that Eq. (3) expresses intrinsic properties of the pion LFWF. This fact is expressed in all existing calculations, as highlighted elsewhere [59,60,61,62,63,64,65,66,67,68,69,70,71,72]. Indeed, working in the chiral limit and assuming that chiral symmetry is dynamically broken, so \(m_\pi ^2=0\) in the presence of a dynamically generated dressed-quark mass [13,14,15,16,17,18], one can show algebraically [10, Sec. 5A] that when quark+antiquark scattering in a four-dimensional theory is mediated by the exchange of a vector boson whose propagator is monotonically decreasing and behaves as \([1/k^2]^{n\ge 0}\) on \(k^2 \gg m_p^2\), thenFootnote 1

(5)

These remarks highlight that experimental validation of Eq. (3) is a stern test of QCD and its role as the theory of strong interactions within the SM. For that reason, we restate the equation in words:

T1: If QCD describes the pion, then at any scale for which an analysis of data using known techniques is valid, the form extracted for the pion’s valence-quark DF must behave as \((1-x)^{\beta }\), \(\beta >2\), on \(x\gtrsim 0.9\) [10, 59, 73, 74]. A result with \(\beta <2\) entails one of the following: [a] the analysis is incomplete, omitting or misrepresenting some aspect or aspects of the processes involved; [b] (some of) the data being considered are not a true expression of a quality intrinsic to the pion; or [c] QCD, as it is currently understood, is not the theory of strong interactions.

3 Experiments and the pion valence-quark distribution

The Drell-Yan process [75, 76] \(\pi + A \rightarrow \ell ^+ \ell ^- + X\), where \(\ell \) is a lepton, A is a nuclear target, and X denotes the debris produced by the deeply inelastic reaction, has been used as the basis for all existing attempts to infer the large-x behaviour of from experiment [29,30,31,32]. Analysing that reaction at leading order in pQCD, two terms appear in the cross-section [55, 77]. One is directly sensitive to the pion LFWF and contributes the Eq. (3) piece to the transverse part of the overall structure function. The other is generated by what may be described as an initial state interaction between the valence antiquark in the pion and a valence quark in the target. It is expressed as a \((1-x)^0\langle k_\perp ^2\rangle / \zeta ^2\) (higher-twist) longitudinal term in the overall cross-section, where \(\langle k_\perp ^2\rangle \ll \zeta ^2\) measures the mean transverse momentum of the annihilating valence-antiquark and -quark. Evidently, in any realisable experiment, there will always be a neighbourhood \(x\simeq 1\) whereupon the longitudinal cross-section exceeds the transverse and access to is obscured.

These remarks are important because existing attempts to extract the large-x behaviour of the pion valence-quark DF [33,34,35,36] are dominated by the data reported in Ref.  [32, E615]. The angular distributions associated with that set and with the data in Ref.  [78, NA10] show signs of the higher-twist contribution [79]. However, extant analyses aimed at determining do not remove this contamination. In addition, the binning of E615 data is not described in Ref.  [32]; and this may introduce substantial uncertainty at large x. For these reasons and also simply because experimental capabilities and techniques have improved in the thirty years since the E615 experiment, new data with much improved precision are essential if Eq. (3) is to be properly tested. Such measurements are in train [23,24,25,26,27,28].

4 Analyses of Drell-Yan data

Data described as [32, E615] “Measured values for the pion valence structure function”, obtained via a leading-order (LO) pQCD analysis of Drell-Yan measurements and associated with a scale \(\zeta =\zeta _5=5.2\,\)GeV, are presented in Fig. 2. The set was used in Ref.  [32] to argue that . This is a marked contradiction of Eq. (3), a fact highlighted in Ref.  [59].

Fig. 2
figure 2

SCI Pion valence-quark DF, , obtained via evolution of the hadron scale DF in Eq. (6): is drawn as the dotted grey curve. Evolution from: \(\zeta =\zeta _H=m_G\), Eq. (8) – dot-dashed black curve, using which \(\chi ^2/\mathrm{datum}=2.7\); \(\zeta =1.05 \zeta _H\) – solid blue curve, \(\chi ^2/\mathrm{datum}=1.0\); \(\zeta =1.10 \zeta _H\) – long-dashed blue curve, \(\chi ^2/\mathrm{datum}=1.0\); and \(\zeta =1.15 \zeta _H\) – upper edge of shaded blue band, \(\chi ^2/\mathrm{datum}=2.2\); Data – black up-triangles: described in Ref.  [32, E615] as measured values for the pion valence structure function

Full size image

Working with the same data, a subsequent next-to-leading-order (NLO) pQCD analysis [33], which also included next-to-leading logarithm (NLL) threshold resummation effects (soft gluon resummation) [80, 81], arrived at the result , a markedly different outcome that is consistent with Eq. (3).

More recently, the effects of threshold resummation in analyses of the E615 data were reconsidered [36], with comparison of three different methods for implementing such effects. Two of the approaches, which lie within the class of Mellin-Fourier schemes [80, 81] and may be called Mellin-Fourier-cosine (MFc) and Mellin-Fourier-expansion (MFe), confirmed the outcome in Ref.  [33], returning a large-x exponent \(\beta (\zeta _c)>2\) at the scale \(\zeta _c=1.27\,\)GeV. The third procedure, described as double-Mellin (dM) [82], produced an effective large-x exponent [36, Eq. (5)] \(\beta _{\mathrm{eff}}(\zeta _c) \sim 1.2\), albeit with large uncertainty. This disagreement admits the interpretation that E615 data and available methods for its analysis are unable to deliver an unambiguous result for \(\beta \); hence, cannot be used to test Eq. (3). This, again, is a call for new, precise data.

5 Content of Drell-Yan data

Given the situation described above, we judge it worth considering the E615 data from a different perspective; hence, ask the following question:

Q1: Suppose one has in hand a body of data which is truly an expression of pion structure and that a complete analysis of that data returns a valence-quark DF whose character is fairly represented by the points drawn in Fig. 2, then what is the associated pion LFWF and what quark+antiquark interaction produced it?

Regarding the data in Fig. 2, whose large-x behaviour is known to be approximately linear [32, E615], the beginning of an answer is suggested by the discussion associated with Eq. (5): this DF can be connected with a vector\(\,\times \,\)vector quark+antiquark interaction that is momentum-independent, viz. \([1/k^2]^{n\simeq 0}\). In fact, this has been known since the first model calculations of Nambu–Goldstone boson valence DFs [83].

To elucidate, we consider the valence-quark DF for a \(m_\pi =0.14\,\)GeV pion obtained using a symmetry-preserving regularisation of a vector\(\,\times \,\)vector contact interaction (SCI) [68]:

(6)

Since it provides a largely algebraic framework, the SCI has been widely used [68, 84,85,86,87] to set benchmarks for more complex studies using QCD-kindred interactions. Comparisons with such calculations [86, 88,89,90] and also data [91,92,93] show that SCI results are typically a reliable guide for long-wavelength properties of hadrons, such as masses and decay constants, but, unsurprisingly, fail in applications that are sensitive to and/or reveal structural properties.

The hadron-scale DF in Eq. (6) is symmetric about \(x=1/2\) because, by construction and implementation, dressed-quark degrees-of-freedom carry all properties of the bound-state at this scale: before gluon radiation begins, all glue and sea partons are sublimated into the dressed quarks. As a consequence, one obtains

(7)

viz. each of the valence-quarks carries half the pion’s light-front momentum at this scale. (Isospin symmetry is assumed.)

The same qualities are expressed in every hadron-scale pseudoscalar-meson valence-quark distribution computed using a framework that respects Poincaré covariance and QCD’s vector and axial-vector Ward–Green–Takahashi identities, especially the pattern and consequences of dynamical chiral symmetry breaking (DCSB) [63,64,65]. Any calculation that delivers an asymmetric DF at \(\zeta _\mathcal{H}\) corresponds, implicitly or explicitly, to a treatment of the matrix element that defines the valence-quark distribution which violates at least one of these key QCD symmetry constraints. (If the meson is built using nondegenerate valence degrees-of-freedom, \(q_1\), \({\bar{q}}_2\), the only changes are and .)

In model calculations, following Ref.  [94], the value of \(\zeta _\mathcal{H}\) is held to be a free parameter. Its value is determined a posteriori by requiring that, after using pQCD’s evolution equations [49,50,51,52] to map the DF to another scale \(\zeta =\zeta _E > \zeta _\mathcal{H}\), where “empirical” information about the DF is available, agreement is obtained with some chosen piece of that \(\zeta =\zeta _E\) information. There is an obvious yet typically overlooked objection to this procedure; namely,

O1. DGLAP evolution expresses properties intrinsic to four-dimensional QCD. It is logically inconsistent to impose QCD-specific gluon radiation and splitting on a DF obtained from a quark+antiquark interaction that has no discernible link with QCD dynamics.

If this caveat is ignored, then one merely arrives at a practitioner-preferred DF fit to the empirical \(\zeta =\zeta _E\) information which supersedes its \(\zeta =\zeta _H\) origin.

It is now worth observing that when considering model calculations which produce a pion-like pseudoscalar-meson valence-quark DF, , that is not symmetric around \(x=1/2\), so that , it is logically consistent to interpret the result as corresponding to a model scale , so long as the model ensures that physical observables are independent of this choice of scale. (If the latter is not true, then the model can be improved so this weakness is remedied.) Then, since the model’s proponent is usually prepared to use pQCD DGLAP evolution to map the DF to a higher scale, it is equally acceptable to use the same procedure and work backwards to determine the scale, \(\zeta _\mathcal{H}\), at which the model’s valence-quarks saturate the momentum sum rule: . (So evolved, the DF might also be symmetric around \(x=1/2\), in which case the procedure delivers a \(\zeta =\zeta _\mathcal{H}\) DF that is implicitly corrected for the symmetry violation(s) in the model’s formulation.) It is now plain that all model studies support the same definition of \(\zeta _\mathcal{H}\).

In the QCD context, distinct from both modelling and data fitting, it is possible to predict a unique value of \(\zeta _H\) based on the behaviour of QCD’s process-independent (PI) effective charge [5, 9], \({\hat{\alpha }}(k^2)\). This charge agrees to better than 0.1% with pQCD’s one-loop coupling on \(k^2 \gtrsim 4 m_p^2\). However, as \(k^2\) continues to run toward zero, \({\hat{\alpha }}(k^2)\) exhibits a qualitative change; and its behaviour on

$$\begin{aligned} k^2 < m_G^2 := (0.331(2)\,\mathrm{GeV})^2\,, \end{aligned}$$
(8)

indicates that gauge sector modes with \(k^2 \lesssim m_G^2\) are screened from interactions and QCD has entered a practically conformal domain. The value of \(m_G\) is determined by \(m_0\), QCD’s renormalisation group invariant gluon mass scale [5, 9, 10]. It follows that the line \(k^2=m_G^2\) separates long- and short-wavelength physics; hence, serves as the natural definition for the hadron scale, viz. \(\zeta _H=m_G\), because, below such momenta, gauge sector modes have effectively decoupled and gluon emission is frozen out. This is the position introduced in Refs.  [9, 64, 65, 69,70,71], which further argue that evolution should be implemented by using \({\hat{\alpha }}(k^2)\) to integrate the one-loop DGLAP equations. Additional features of this all-orders evolution scheme are described elsewhere [88, Sec. VII]. Naturally, given the properties of \({\hat{\alpha }}(k^2)\), the approach is equivalent to standard DGLAP evolution on any domain upon which pQCD is applicable.

Neglecting O1 for the moment, then the all-orders evolution scheme yields the dot-dashed black curve in Fig. 2. This prediction is made without reference to E615 data so the \(\chi ^2/\mathrm{datum}=2.7\) outcome is striking. As shown in Fig. 2, supposing only that the prediction for \(\zeta _\mathcal{H}\) in Eq. (8) has a 5% uncertainty, leads to a description of the data with \(\chi ^2/\mathrm{datum}=1.0\).

Returning to the question posed at the beginning of this section, Q1, the outcomes just described invite the following answer: A valence-quark DF whose character is fairly represented by the points drawn in Fig. 2 is derived from a pion LFWF that is produced by a momentum-independent quark+antiquark interaction. However, this position is untenable because DGLAP evolution cannot be derived from such an interaction. Logically, therefore, the data in Fig. 2 cannot be a true representation of the valence-quark DF of a pion generated by QCD interactions. Thus, either the analysis leading to such data was incomplete or the experimental data on which the analysis was based is not a veracious expression of an intrinsic property of the pion, or both issues have played a role.

6 Character and consequences of all orders evolution

A dM approach to including NLL threshold resummation in the analysis of Drell-Yan data is described in Ref.  [36]. Based on the information therein, we have developed the following approximate form for the extracted valence-quark DF:

(9a)
(9b)

where ensures unit normalisation. Figure 3 depicts the E615 data from Fig. 2 projected onto the central curve in Eq. (9) after its evolution \(\zeta _c \rightarrow \zeta _5\). The original and new analyses are qualitatively equivalent.

On the other hand, however, there is a key quantitative difference between the original and new analysis, viz. the Ref.  [36] fits lodge less of the pion’s light-front momentum with the valence quarks: at \(\zeta _c\), the fraction is 0.46(3). This value is significantly smaller than the prediction 0.52(2) that is typical of contemporary continuum and lattice QCD calculations at this scale [70, 71].

Fig. 3
figure 3

Pion valence-quark DF (dot-dashed red curve) inferred from the fit to data in Ref.  [36] that was obtained with NLL threshold resummation performed using the dM method, evolved \(\zeta _c \rightarrow \zeta _5\). [The encompassing band expresses the uncertainty in Eq. (9).] At this scale, the large-x exponent is \(\beta \approx 1.5\). Data recorded in Ref.  [32, E615] – black open up-triangles; and red down-triangles – projection of that data onto the \(\zeta _c \rightarrow \zeta _5\) evolution of Eq. (9), which we denote as E615\(_{\mathrm{dM}}\)

Full size image

Regarding the E615\(_{\mathrm{dM}}\) data in Fig. 3, it is plain that a repetition of the analysis in Sect. 5 would yield a qualitatively equivalent outcome to that exposed by our discussion of the original E615 data. We will therefore follow a more general approach based on a single proposition:

P1: There exists an effective charge, \(\alpha _{1\ell }(k^2)\), that, when used to integrate the one-loop pQCD DGLAP equations, defines an evolution scheme for parton DFs that is all-orders exact. \(\alpha _{1\ell }(k^2)\) need not be unique.

It is here worth recording some background and corollaries of this proposition.

A context for P1 is provided by the discussion of process-dependent charges in Refs.  [95,96,97]; so, this charge need not be process independent, e.g., different charges may be needed for distinct observables. On the other hand, a process-independent charge with this character is not excluded, as evidenced by the efficacy of that discussed in Refs.  [5, 9, 69,70,71].

As detailed elsewhere [98], in being defined by an observable (here, a structure function) \(\alpha _{1\ell }(k^2)\) is consistent with the renormalisation group. It is also renormalisation scheme independent; everywhere finite and analytic; and supplies an infrared completion of any standard running coupling.

As a consequence of finiteness, P1 entails that there is a scale, which we have denoted \(\zeta _\mathcal{H}\) above, such that

(10)

viz. whereat the \(\zeta =\zeta _\mathcal{H}\) dressed valence quarks carry all the pion’s light-front momentum. Momentum conservation now demands that the glue and sea light-front momentum fractions vanish at \(\zeta _\mathcal{H}\); and since all physical DFs are nonnegative on \(x\in [0,1]\), it follows that the glue and sea DFs vanish identically at this scale.

Eq. (10) is guaranteed for any symmetric function, i.e., when

(11)

This is a sufficient but not necessary condition. However, as noted after Eq. (7), Eq. (11) is satisfied by the result of any calculation which respects Poincaré covariance and QCD’s vector and axial-vector Ward–Green–Takahashi identities, particularly the pattern and consequences of DCSB [63,64,65]. A violation of Eq. (11) by a valence-quark DF obtained by whatever means links this DF to a treatment of the underlying matrix element that breaks basic QCD symmetries.

As noted above, when considering a pseudoscalar meson built using nondegenerate valence degrees-of-freedom, \(q_1\), \({\bar{q}}_2\), the only changes are and .

Now, explicitly describing the isospin symmetry limit, consider the Mellin moments:

(12)

As highlighted in Ref.  [65, Eq. (29)], Eq.  (11) entails that for \(n\in {{\mathbb {Z}}}^{\ge }\), is linearly dependent on (completely determined by) the set of even moments with \(m\le n\). For instance,

(13a)
(13b)

etc. These constraints can most straightforwardly be exploited using the following recursion relation for the odd-order Mellin moments:

(14)

Furthermore, given P1, then [88, Sec. VII]:

(15)

where, for \(n_f=4\) quark flavours,

$$\begin{aligned} \gamma _0^n = -\frac{4}{3} \left( 3 + \frac{2}{(n+1)(n+2)} - 4 \sum _{j=1}^{n+1} \frac{1}{j} \right) \,. \end{aligned}$$
(16)

Namely, given the pion valence-quark DF at one scale, here identified as \(\zeta _\mathcal{H}\), then its full x-dependence at any other scale \(\zeta \) is completely determined by the value of its first Mellin moment at \(\zeta \). No other knowledge is needed, including and especially not information on the form of \(\alpha _{1\ell }(k^2)\). Inserting Eq. (15) into Eq. (14), one finds

(17)

Conversely, any DF whose moments fulfill this recursion relation is related by evolution to a symmetric function at \(\zeta _\mathcal{H}\).

This corollary can be examined in connection with any inferred or calculated pion DF. Hence, it can be tested with Eq. (9). Figure 4 displays the first 15 odd moments computed directly from the central curve in Eq. (9) compared with the values for these moments predicted by Eq. (17): there is precise agreement. (Given the factorial growth of the coefficients and exponential decrease in the moments, one must be careful to eliminate round-off error when evaluating the sum in Eq. (17).)

Fig. 4
figure 4

Odd moments computed from the DF in Eq. (9): directly, using Eq. (12) (blue asterisks); and using the recursion relation in Eq. (17) (gold open diamonds). Also, odd moments computed from the DF in Eq. (23): directly, using Eq. (12) (green stars); and using the recursion relation in Eq. (17) (navy-blue open circles). Recall: the recursion relation is satisfied if, and only if, the source function is related by evolution to a DF that is symmetric at some scale \(\zeta _\mathcal{H}\); evidently, this applies to Eq. (9) but not Eq. (23)

Full size image

The agreement verified by Fig. 4 entails that the DF in Eq. (9), which fairly represents that extracted from a particular analysis of data relevant to the pion valence-quark DF [36], is related by evolution to a symmetric DF at a scale \(\zeta _\mathcal{H}<\zeta _c\). As noted above, momentum conservation now demands that the glue and sea DFs must both vanish identically at this scale. Consequently, when analysing data, one should not independently fit valence-quark, glue, and sea distributions at any scale \(\zeta >\zeta _\mathcal{H}\) because the glue and sea distributions at \(\zeta > \zeta _\mathcal{H}\) are completely determined by and the evolution equations: the glue and sea distributions are not independent functions.

It is also instructive to test Eq. (17) using a calculated pion valence-quark DF instead of a fit to data. In this connection, consider that Ref.  [99] used lattice-QCD output to build a x-dependent pion valence-quark DF at \(\zeta =\zeta _5\). (Its large-x behaviour is consistent with Eq. (3) [99, Fig. 12].) The first six nontrivial Mellin moments of this DF are listed in the first column here, with the indicated statistical and systematic errors:

(18)

The second column reports the results computed from the first column using Eq. (17). The underlined entries are the linearly independent moments, used as input, and the remaining entries are the Eq. (17) predictions. Notably, the linearly dependent \(n=7\) moment is not reported in Ref.  [99]; so, our result is a prediction. Working with Ref.  [99, Eq. (46) and Table X], one obtains . The agreement between lattice calculation and the even-function recursion relation is striking. It suggests that this lattice result, too, is linked via evolution to a symmetric DF at a scale \(\zeta _\mathcal{H}<\zeta _5\).

These two examples illustrate the following general conclusion:

C1: Given a pion valence-quark DF, , either fitted to data or calculated at some scale \(\zeta _\mathcal{E}\), whose moments satisfy Eq. (17), then there is always a \(\zeta _\mathcal{H}\in [0,\zeta _\mathcal{E})\) such that evolution \(\zeta _\mathcal{E} \rightarrow \zeta _\mathcal{H}\) maps , with . At \(\zeta _\mathcal{H}\), valence-quarks carry all the pion’s light-front momentum and the glue and sea distributions vanish.

7 Reviewing E615\(_{\mathrm{\mathbf{dM}}}\) data

Armed with P1 and its manifold consequences, we now consider what can additionally be revealed about the content of the E615\(_{\mathrm{dM}}\) data in Fig. 3. (Sect. 6 has already shown that the data fit is linked via evolution to a symmetric DF at some scale \(\zeta _\mathcal{H}<\zeta _c\).) We do not constrain ourselves to solely this set, however. Instead, we consider an array of possibilities with the same qualitative character, viz. a set representative of the outcome of any analysis of a body of data that returns a pion valence-quark DF for which the effective large-x exponent is \(\beta (\zeta _c) \approx 1.2\).

To achieve this generalisation, we first suppose that a fair approximation to any such DF is provided by

(19)

where ensures unit normalisation. This is a weak assumption: any of the forms commonly used in fitting data would serve equally well. Then, we proceed as follows. (i) Determine central values of \(\{\alpha _i^\zeta |i=1,2,3\}\) via a least-squares fit to the \(\zeta =\zeta _5\) E615\(_{\mathrm{dM}}\) data. (ii) Generate a new vector \(\{\alpha _i^\zeta |i=1,2,3\}\), each element of which is distributed randomly around its best-fit value. (iii) Using the DF obtained therewith, evaluate

(20)

where \(\{(x_l,u_l\pm \delta _l) | N=1,\ldots \,40\}\) is the E615\(_{\mathrm{dM}}\) data set. This \(\{\alpha _i\}\) configuration is accepted with probability

(21)

where \(d=N-3\) and \(\chi _0^2 \approx d\) locates the maximum of the \(\chi ^2\)-probability density, \(P(\chi ^2;d)\). (iv) Repeat (ii) and (iii) until one has a \(K \gtrsim 1000\)-member ensemble of DFs. This procedure yields the collection of DFs drawn in Fig. 5A.

With any collection of pion DFs in hand, one also has associated ensembles of Mellin moments: \(\{ \mathcal{M}_j(n;\zeta _5) | j = 1, \ldots , K\}\),

(22)

where \(n \in {{\mathbb {Z}}}^{\ge }\) and \(\{ [\alpha _i]^j | j=1,\ldots ,K\}\) are the coefficient vectors associated with the accepted DFs. Equation (15) can be used to evolve each of these moments to any other scale, \(\zeta \); and given that a continuous function of compact support is uniquely defined by its Mellin moments, one can therefrom reconstruct each of the distributions at \(\zeta \). Plainly, the mean and standard deviation of the \(n^{\mathrm{th}}\) moment can be computed at the new scale.

Using Eq. (15), we evolved each one of the K replicas of the E615\(_{\mathrm{dM}}\) fit in Fig. 5 from \(\zeta _5\rightarrow \zeta _\mathcal{H}\) with the results drawn in Fig. 5B: all the evolved curves are consistent with C1. Furthermore, comparison with the SCI result in Eq. (6) shows that every one of the evolved curves is contained within the class of DFs linked to a momentum-independent quark+antiquark interaction.Footnote 2

These outcomes add weight to the discussion in Sect. 5. Namely, the replicas in Fig. 5 serve as a representative body of extracted forms for the pion’s valence-quark DF whose large-x behaviour is determined by the E615 Drell-Yan data. Accepting P1, then that data cannot be connected with a pion LFWF which is derived from a momentum-dependent quark+antiquark interaction; hence, the data cannot be a fair representation of the valence-quark DF of a pion generated by QCD interactions. Consequently, either the analysis of data that formed the basis for the replicas was incomplete or the source Drell-Yan data itself does not truly express an intrinsic property of the pion, or both problems have interfered together.

Fig. 5
figure 5

Upper panel – A. Randomly distributed ensemble of pion valence-quark DFs (red curves) constructed from the E615\(_{\mathrm{dM}}\) data (red down-triangles) using the procedure described in connection with Eq. (21). Lower panel – B. \(\zeta _5\rightarrow \zeta _\mathcal{H}\) evolution of each curve in Panel A (magenta curves). The SCI result in Eq. (6) is drawn as the dotted black curve

Full size image

8 Threshold resummation using Mellin-Fourier methods

Threshold resummation is also discussed in Refs.  [80, 81]. That scheme involves a Mellin transform and Fourier transform, with the latter operation inviting two differing treatments: MFc, in which the cosine is retained, and MFe, where the cosine is approximated by unity because the argument introduces subleading corrections. Ref.  [33] used MFc.

Ref.  [36] adapted these approaches to produce independent fits whose large-x behaviour is constrained by E615 data. A fair representation of those \(\zeta =\zeta _c\) fits is provided by the following parametrisation:

(23a)
(23b)

which represent an average of the MFc and MFe results and express an effective large-x exponent \(\beta (\zeta _c) = 2.45(11)\).

Matching the outcome of the dM analysis, the Ref.  [36] MF fits represented by Eq. (23) produce , i.e., lodge a markedly smaller fraction of the pion’s light-front momentum with the valence quarks than existing calculations predict. This is highlighted in Fig. 6 by comparison with the parameter-free prediction of the pion valence-quark DF in Refs.  [69,70,71], evolved from the \(\zeta =\zeta _\mathcal{H}\) result:

(24)

whose dilated (hardened) profile owes to EHM. Plainly, the area under the dot-dashed red curve is less than that under the solid blue curve.

Fig. 6
figure 6

Valence-quark DF at \(\zeta =\zeta _c=1.27\,\)GeV. Dot-dashed red curve within like-coloured band – Eq. (23), derived from Ref.  [36], with \(\beta _{\mathrm{eff}}(\zeta _c) = 2.45(11)\); solid blue curve and band – prediction from Refs.  [69,70,71], which express \(\beta _{\mathrm{eff}}(\zeta _c) = 2.52(5)\); and dashed purple curve and band – Eq. (26), \(\beta _{\mathrm{eff}}(\zeta _c) = 2.06(2)\)

Full size image

Figure 6 exposes a curious feature of the MF fits in Ref.  [36]; namely, the DFs are not of uniform concavity on \(x\lesssim 0.5\). Potential causes of such “wiggles” in the Ref.  [36] fits are limitations introduced by the simple DF fitting Ansatz employed, written in Eq. (23), and/or choosing to treat valence, glue, and sea DFs as uncorrelated at \(\zeta _c\). These particular features are not evident in the dM results, Eq. (9), possibly because DFs consistent with T1 have a rich structure, whose reliable extraction requires a more nuanced approach to fitting.

We have checked the DF in Eq. (23) against the P1 corollary in Eq. (17), with the result displayed in Fig. 4. Evidently, the recursion relation is not satisfied. This failure can be traced to the wiggles on \(x\lesssim 0.5\). They can be eliminated as follows. (i) Using Eq. (15), evolve the DF in Eq. (23) to \(\zeta _\mathcal{H}\), whereat . (ii) Working with the first 101 Mellin moments at \(\zeta _\mathcal{H}\), determine the best least-squares fit achievable with a symmetric function. The form

(25)

is efficacious: \(\alpha = 1.28(9)\) reproduces the moments with mean absolute relative error 6(5)%. (iii) Evolve the symmetric functions back to \(\zeta _c\) using Eq. (15). This procedure delivers DFs that may be approximated as

(26a)
(26b)

yield , express an effective large-x exponent \(\beta (\zeta _c) = 2.06(2)\), and are represented in Fig. 6 by the dashed purple curve and linked band. Evidently, requiring elimination of the wiggles forces the peak of to shift toward \(x=0\) whilst maintaining the area under the curve.

The DFs approximated by Eq. (26) reproduce the first 101 moments of the DFs in Eq. (23) with mean absolute relative error 7(6)%; and in being uniformly concave-down on \(x\in (0,1)\), we judge their pointwise behaviour to be a more realistic portrayal of pion structure than Eq. (23) – modern calculations of ground-state pseudoscalar-meson LFWFs, e.g., Refs.  [40, 85, 88, 101,102,103,104,105,106,107] do not exhibit features that could produce wiggles in the associated valence-quark DFs. Consequently, we hereafter represent Eq. (26) as a fair and physical sketch of the MF results in Ref.  [36]. Expressed another way, we interpret Eq. (26) as a smoothed approximation to Eq. (23), equivalent in (most) practically measurable respects. It is conceivable that Eq. (26) could emerge following refinements of the analysis scheme in Ref.  [36]; but given that C1 applies to Eq. (26), that would require introduction of intimate correlations between the valence, glue, and sea DFs.

Evolved to \(\zeta =\zeta _5=5.2\,\)GeV, the DFs associated with Eq. (26) are drawn in Fig. 7A. The pion valence-quark DF inferred in Ref.  [33] from E615 data is also drawn in Fig. 7A: on \(x\gtrsim 0.3\) it agrees with the prediction in Refs.  [69,70,71] and yields a compatible valence-quark light-front momentum fraction. The treatment of glue and sea distributions in Ref.  [33] is very different from that in Refs.  [69,70,71], explaining the \(x\lesssim 0.3\) discrepancy.

Fig. 7
figure 7

Upper panel – A. Pion valence-quark DF (dot-dashed red curve) inferred from fits to data in Ref.  [36] obtained with NLL threshold resummation performed using the MF method, evolved \(\zeta _c \rightarrow \zeta _5\): [The encompassing band expresses the uncertainty in Eq. (26).] At this scale, the large-x exponent is \(\beta = 2.24(7)\). Data recorded in Ref.  [32, E615] – black open up-triangles; and red squares – projection of that data onto the \(\zeta _c \rightarrow \zeta _5\) evolution of Eq. (26), which we denote as E615\(_\mathrm{MF}\). Parameter-free prediction from Refs.  [69,70,71]: blue curve within like-coloured band, which yields . MFc fit to data [32, E615] in Ref.  [33]: dashed green curve within like-coloured band, . Lower panel – B. As in A, except: data from Ref.  [32, E615] and fit from Ref.  [33] are removed; and lattice-QCD prediction from Ref.  [100] is included, drawn as dashed grey curve within like-coloured band

Full size image
Fig. 8
figure 8

Upper panel – A. Randomly distributed ensemble of pion valence-quark DFs (orange curves) constructed from the E615\(_{\mathrm{MF}}\) data (red squares) using the procedure described in connection with Eq. (21). Lower panel – B. \(\zeta _5\rightarrow \zeta _\mathcal{H}\) evolution of each curve in Panel A (purple curves). Parameter-free \(\zeta _\mathcal{H}\) prediction in Refs.  [69,70,71], which expresses EHM-induced dilation: dashed blue curve. Scale free DF, : dotted black curve

Full size image

Further to these points, working with the E615\(_{\mathrm{MF}}\) data in Fig. 7, we repeated the analysis described in Sect. 7. The \(K\gtrsim 1000\) replicas of this data are drawn in Fig. 8A and their \(\zeta _5\rightarrow \zeta _\mathcal{H}\) evolutions in Fig. 8B: all evolved curves are consistent with C1. Moreover, comparison with the continuum prediction from Refs.  [69,70,71], Eq. (24), whose dilated (hardened) profile owes to EHM, demonstrates that, within mutual uncertainties, the set of evolved curves abuts the class of DFs linked to a quark+antiquark interaction with the \(1/k^2\) ultraviolet behaviour which characterises QCD.Footnote 3 Hence, in contrast to the dM scheme, the MF approaches to NLL resummation in the analysis of Drell-Yan data lead to valence-quark DFs whose large-x behaviour is approximately consistent with Eq. (3) and can thus be connected with pion LFWFs whose momentum dependence matches QCD’s prediction, Eq. (4). In this case, the data may be identified as resulting from measurements which are sensitive to intrinsic properties of the pion.

With the information now at hand, it is possible to follow Ref.  [88, Sec. VIII.A] and record the pion mass-squared fractions carried by the different parton species at \(\zeta =\zeta _c\). The results are listed in Table 1 and illustrated in Fig. 9. Evidently, despite the marked differences between their pointwise behaviour, the dM and MF fits in Ref.  [36] yield practically the same pion mass-squared apportionments; and compared with the predictions associated with Eq. (24), they lodge far more of \(m_\pi ^2\) with the sea-quarks at the cost of the valence-quark fraction. Curiously, the relative stability of the valence-quark fraction indicates that the sea fraction increases as the glue fraction decreases and vice versa. This is counterintuitive, given the character of QCD evolution, with gluon splitting serving to populate the sea. Hence, it is likely an artefact in Ref.  [36] of treating parton DFs as independent at \(\zeta _c\) in the fitting.

Table 1 Pion mass-squared fractions carried by different parton species at \(\zeta =\zeta _c=1.27\,\)GeV
Full size table

9 Summary and outlook

QCD predictions for the behaviour of the pion Bethe-Salpeter- and light-front wave functions at large relative momenta have been known for more than forty years; and we have highlighted [Sect. 2] that all calculations of the pion valence-quark distribution function (DF) which are consistent with these predictions deliver the following result:

(27)

where \(\gamma (\zeta ) \ge 0\) grows logarithmically with \(\zeta \). Equivalently, using any known evolution scheme, no description of pion structure that is consistent with QCD pion wave functions can produce a large-x exponent \(\beta <2\) on \(\zeta \gtrsim m_p\), where \(m_p\) is the proton mass.

Fig. 9
figure 9

Pion mass-squared fractions carried by different parton species at \(\zeta =\zeta _c=1.27\,\)GeV. Upper panel – A. Row 2 in Table 1, inferred in Ref.  [36] Lower panel – B. Row 3 in Table 1, calculated using the DF predictions in Refs.  [69,70,71]

Full size image

Yet some approaches to the analysis of extant data, which should serve to constrain the large-x behaviour of , nevertheless violate Eq. (27). For instance, whereas fits based on next-to-leading-order perturbative QCD treatments of hard scattering kernels that include soft gluon resummation using the Mellin-Fourier (MF) scheme yield DFs in agreement with Eq. (27), those which use the double-Mellin (dM) approach to soft gluon resummation do not. This does not necessarily mean that the MF approaches are the more rigorous, but it does justify continued efforts aimed at understanding their efficacy in delivering DFs that are consistent with QCD predictions. On the other hand, if one chooses to favour the dM approach, for some reason or another, then the associated disagreement with Eq. (27) requires explanation; and these are the only possibilities: [a] the dM scheme is incomplete, omitting or misrepresenting some aspect or aspects of the hard processes involved; [b] (some of) the data being considered in the analysis are not a true expression of a quality intrinsic to the pion; or [c] QCD, as it is currently understood, is not the theory of strong interactions.

In supporting these conclusions, we explained an array of corollaries that follow from a single proposition [P1 in Sect. 6]: there exists an effective charge that defines an evolution scheme for parton DFs that is all-orders exact. A great deal can be concluded from this proposition without specifying the form of the charge. For instance, P1 entails that there is a scale \(\zeta =\zeta _\mathcal{H} \in [0,m_p)\) whereat valence-quarks carry all the pion’s light-front longitudinal momentum; \(\gamma (\zeta _\mathcal{H})\equiv 0\); and, consequently, glue and sea distributions are completely determined by the nonperturbative information embedded in and revealed in evolution from this scale. These outcomes are guaranteed when , a quality expressed in all calculations that respect Poincaré covariance and QCD’s vector and axial-vector Ward–Green–Takahashi identities. Such symmetry entails that any given odd Mellin moment of is linearly dependent upon the even Mellin moments of lower order. This corollary is expressed in a recursion relation that can be used to reveal the character of any DF, whether fitted or calculated.

Forty years after the first experiment to collect data amenable for use in constraining the large-x behaviour of , the answer, for some practitioners, remains uncertain. In large part, the paucity of such data and its imprecision are responsible. Modern and anticipated facilities promise to remedy these issues. New developments in phenomenology and theory are required to eliminate the others.