Helicity at Small x : Oscillations Generated by Bringing Back the Quarks
HHelicity at Small x : Oscillations Generated by Bringing Back the Quarks Yuri V. Kovchegov ∗ and Yossathorn Tawabutr † Department of Physics, The Ohio State University, Columbus, OH 43210, USA
We construct a numerical solution of the recently-derived large- N c & N f small- x helicity evolutionequations [1] with the aim to establish the small- x asymptotics of the quark helicity distributionbeyond the large- N c limit explored previously in the same framework. (Here N c and N f are thenumbers of quark colors and flavors.) While the large- N c helicity evolution involves gluons only, thelarge- N c & N f evolution includes contributions from quarks as well. We find that adding quarks tothe evolution makes quark helicity distribution oscillate as a function of x . Our numerical resultsin the large- N c & N f limit lead to the x -dependence of the flavor-singlet quark helicity distributionwhich is well-approximated by∆Σ( x, Q ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c & N f ∼ (cid:18) x (cid:19) α qh cos (cid:20) ω q ln (cid:18) x (cid:19) + ϕ q (cid:21) . (1)The power α qh exhibits a weak N f -dependence, and, for all N f values considered, remains veryclose to α qh ( N f = 0) = (4 / √ (cid:112) α s N c / (2 π ) obtained earlier in the large- N c limit [2, 3]. The noveloscillation frequency ω q and phase shift ϕ q depend more strongly on the number of flavors N f (with ω q = 0 in the pure-glue large- N c limit). The typical period of oscillations for ∆Σ is rather long,spanning many units of rapidity. We speculate whether the oscillations we find are related to thesign variation with x seen in the strange quark helicity distribution extracted from the data [4–7]. PACS numbers: 12.38.-t, 12.38.Bx, 12.38.Cy
CONTENTS
I. Introduction 1II. Helicity Evolution Equations at Large- N c & N f I. INTRODUCTION
Understanding the partonic (quark and gluon) structure of the proton is an essential part of our understanding ofQuantum Chromodynamics (QCD). One of the big open questions in the studies of proton structure is the protonspin puzzle. Proton is a composite spin-1 / ∗ Email: [email protected] † Email: [email protected] a r X i v : . [ h e p - ph ] J un should be a sum of all the spins and orbital angular momenta (OAM) carried by the quarks and gluons comprisingthe proton. This statement is formalized in terms of helicity sum rules, one due to Jaffe and Manohar [8] and anotherone due to Ji [9]. The former reads S q + L q + S G + L G = 12 (2)where S q and S G are the total spin carried by the quarks and gluons, respectively, whereas L q and L G are the OAMcarried by the quarks and gluons in the proton.One can write S q and S G , which are functions of the momentum scale Q , as integrals of helicity distributionfunctions over the Bjorken x variable, S q ( Q ) = 12 (cid:90) dx ∆Σ( x, Q ) , S G ( Q ) = (cid:90) dx ∆ G ( x, Q ) , (3)where ∆Σ( x, Q ) = (cid:88) f = u,d,s,... (cid:2) ∆ f ( x, Q ) + ∆ ¯ f ( x, Q ) (cid:3) (4)is the flavor-singlet quark helicity distribution function. Here ∆ f and ∆ ¯ f denote the helicity distribution functionsof quarks and antiquarks, respectively. In Eq. (3), ∆ G is the gluon helicity distribution. The quark and gluonOAM, L q and L G , can also be written as x -integrals of the OAM distributions [10–14]. Note that the latter cannotbe expressed as expectation values of twist-2 operators in the proton state, unlike the helicity parton distributionfunctions (hPDFs).The proton spin puzzle started with the discrepancy between the prediction of the (perhaps, na¨ıve) constituent quarkmodel stating that all the proton spin is carried by the constituent quarks (such that S q = 1 /
2) and the experimentalmeasurements by the European Muon Collaboration (EMC) [15, 16] that reported a much lower number for S q . Morerecently, the experiments at the Relativistic Heavy Ion Collider (RHIC) measured a non-zero gluon spin S G [17, 18].The current experimental values for quark and gluon spin are S q ( Q = 10 GeV ) ≈ . ÷ .
20, integrated over0 . < x <
1, and S G ( Q = 10 GeV ) ≈ . ÷ .
26, integrated over 0 . < x < S q and S G do not add up to 1/2 and the proton spin puzzle is still open. (Note alsothat none of the terms in Eq. (2) is positive-definite.) We conclude that the missing part of the proton’s spin mustcome either from the smaller- x regions in the integrals of (3) for the spin terms S q and S G than the x -ranges reportedabove or from the OAM terms L q and L G .No experiment, present or future, can measure helicity PDFs all the way down to x = 0, since this would requireinfinite energy. In addition, at very small (but non-zero) values of x higher-twist corrections become comparable tothe leading-order contribution for many observables, making it hard (if not impossible) to extract twist-2 hPDFs fromthe data. (The way such corrections may come in is described using the parton saturation framework, see [24–29]for reviews.) Therefore, to assess the amount of parton spin and OAM coming from the small- x region, one has todevelop a robust theoretical formalism, which, if able to describe the existing experimental data and predict the valuesof the future longitudinal spin measurements, can be extrapolated down to x = 0 with a (hopefully) good level ofconfidence.First theoretical calculation of helicity distributions at small x was carried out over two decades ago by Bartels,Ermolaev, and Ryskin [30, 31] using infrared evolution equations technique constructed by Kirschner and Lipatov [32](see also [33–35]) to resum powers of the leading parameter α s ln (1 /x ) with α s the strong coupling constant. Resum-mation of α s ln (1 /x ) is referred to as the double logarithmic approximation (DLA). The more recent years have seenrenewed efforts to construct a formalism capable of describing helicity PDFs and OAM distributions at small- x [1, 36–39]. In a series of papers [1–3, 40–43] the helicity evolution at small- x was constructed using the shock-wave formalism[44–48] while employing the so-called polarized Wilson lines, similar to (but not the same as) the formalism employed inderiving the unpolarized Balitsky–Kovchegov (BK) [44–47] and Jalilian-Marian–Iancu–McLerran–Weigert–Leonidov–Kovner (JIMWLK) [48–53] evolution equations, which employed regular light-cone Wilson lines.Helicity evolution equations were derived in [1, 40, 42] using the DLA. Similar to the BK evolution, the equationsclose in the large- N c limit, with N c the number of quark colors. In addition, and different from the unpolarized case,the helicity evolution equations [1, 40, 42] also close in the large- N c & N f limit, where N f is the number of quarkflavors. Owing to the complexity of the large- N c & N f equations [1, 42], only the large- N c equations have been solvedto date, resulting in the following power-law small- x asymptotics for the quark and gluon helicity distributions inDLA [2, 3, 41] ∆Σ( x, Q ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c ∼ (cid:18) x (cid:19) α qh , ∆ G ( x, Q ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c ∼ (cid:18) x (cid:19) α Gh (5)with α qh (cid:12)(cid:12)(cid:12)(cid:12) large- N c = 4 √ (cid:114) α s N c π , α Gh (cid:12)(cid:12)(cid:12)(cid:12) large- N c = 134 √ (cid:114) α s N c π . (6)(When writing equations like (5) we suppress the potentially non-trivial Q dependence of the involved observables,along with the possible sub-dominant x -dependence in the prefactor, and concentrate on the leading x dependence ofthe quantities.) The small- x asymptotics of quark and gluon OAM distributions at large N c were found in [38] usingthe same formalism.The goal of this work is to numerically solve the large- N c & N f helicity evolution equations derived in [1, 42] and usethe solution to assess possible corrections to the large- N c result (5) for the asymptotic behavior of ∆Σ( x, Q ) at small x . The results can be (and are) combined with the values of ∆Σ extracted from the experiment [4–7] to extrapolate∆Σ to smaller values of x and estimate the total quark contribution to the proton spin S q . In addition, our resultscan be employed to find the small- x asymptotics of ∆ G and the OAM distributions for quarks and gluons in thelarge- N c & N f limit, that is, beyond the large- N c expressions for ∆ G in Eq. (5) and for the OAM small- x asymptoticsderived in [38].The paper is structured as follows. In Section II, we rewrite and simplify the large- N c & N f helicity evolutionequations from [1, 42]. The equations are written in terms of the quark and gluon polarized dipole amplitudes, whichare defined below as correlators of the polarized and regular light-cone Wilson lines. Our numerical algorithm ispresented in Sec. III, where we rewrite the equations in terms of the variables convenient for our numerical approach,discretize the equations and cast them in the form suitable for implementing our algorithm. The resulting numericalsolutions of the large- N c & N f helicity evolution equations for N f = 2 , , oscillate as functions of energy (or of ln(1 /x )), changing the sign back and forth. These oscillations of dipole amplitudes, inturn, result in oscillations in ∆Σ( x, Q ) as a function of ln(1 /x ), as demonstrated in Eq. (32) or, equivalently, Eq. (1)above, and depicted in Fig. 4. Such oscillating behavior was absent in the large- N c (and N f = 0) analyses carriedout in the previous works [2, 3, 38, 41], which only saw a power-of-1 /x growth (5) of helicity distributions at small x .The oscillations of ∆Σ( x, Q ) with ln(1 /x ) are the main qualitative result of this paper.The parameters α qh , ω q , and ϕ q from Eq. (1), describing ∆Σ( x, Q ) at small x and at large- N c & N f , are extractedfrom our numerical solution in Appendix B and are summarized in Eq. (33). In Sec. IV we also present a fit (29) forthe N f -dependence of the oscillation frequency ω q which increases with N f , meaning that that the oscillation periodgets shorter and, therefore, the oscillations become more pronounced, in the regime where more quark flavors becomedynamically relevant. In Sec. V we follow the method from [2] to construct a preliminary estimate of the impact ofthe asymptotic form for ∆Σ( x, Q ) given by Eq. (1) on the amount of spin carried by the quarks in the proton. Theresults are given in Figs. 5 and 6, similar to [2] indicating a potential for a large correction. We conclude in Sec. VIand speculate whether the oscillations we find in our solution are related to the oscillation of the strange quark hPDFsextracted from the experimental data in [4–7]. II. HELICITY EVOLUTION EQUATIONS AT LARGE- N c & N f At small x , in the DLA, the flavor-singlet quark helicity PDF can be written as [1, 40, 42]∆Σ( x, Q ) ≡ (cid:88) f [∆ q f ( x, Q ) + ∆¯ q f ( x, Q )] = N c N f π (cid:90) Λ /s dzz z Q (cid:90) z s dx x Q ( x , z ) , (7)where s ≈ Q /x is the center-of-mass energy squared, Λ is a transverse momentum scale characterizing the targetbefore small- x evolution, and the impact parameter-integrated polarized quark dipole amplitude is defined by Q ( x , z ) = (cid:90) d (cid:18) x + x (cid:19) Q ( z ) , (8)where the polarized quark dipole amplitude is Q ( z ) ≡ N c Re (cid:28) T tr (cid:20) V (cid:16) V pol (cid:17) † (cid:21) + T tr (cid:104) V pol V † (cid:105)(cid:29) . (9)The dipole consists of a quark and an anti-quark located at transverse positions x and x depending on the term inEq. (9). In our notation the transverse vectors are denoted by x = ( x , x ) with x ij = x i − x j , such that x = x − x .In addition, x ij denotes the magnitude of the vector x ij . The light-cone variables are defined by x ± = ( x ± x ) / √ x + direction. The quantity z is the fraction of some original ( x − -direction-moving) probe’s momentum carried by the softest (anti-)quark in the dipole. The angle brackets in Eq. (9)denote the averaging in the proton wave function in the small- x /saturation sense [24–29], albeit now taking intoaccount that the proton is longitudinally polarized.The polarized dipole amplitude (9) is defined in terms of the light-cone fundamental Wilson lines V x [ b − , a − ] = P exp ig b − (cid:90) a − dx − A + ( x + = 0 , x − , x ) (10)and the so-called polarized Wilson lines V pol [1, 40–42], consisting of regular light-cone Wilson lines along with thesub-eikonal helicity-dependent local operator insertion(s), V polx = ig p +1 ∞ (cid:90) −∞ dx − V x [+ ∞ , x − ] F ( x − , x ) V x [ x − , −∞ ] (11) − g p +1 ∞ (cid:90) −∞ dx − ∞ (cid:90) x − dx − V x [+ ∞ , x − ] t b ψ β ( x − , x ) U bax [ x − , x − ] (cid:20) γ + γ (cid:21) αβ ¯ ψ α ( x − , x ) t a V x [ x − , −∞ ] . (Note that the definition (11) is different by an extra factor of s in overall normalization as compared to that usedin earlier works [41, 42].) We also use an abbreviated notation where Q = Q x ,x , V pol = V polx , etc. As usual, T in Eq. (9) denotes time ordering of the operators. In equations (10) and (11), A + is the eikonal gluon field of theshock wave (often referred to as the target, or the target proton), F is the helicity-dependent sub-eikonal field-strength tensor component of the target gluon field (both in the fundamental representation), while ψ and ¯ ψ are the(sub-eikonal) quark and anti-quark fields of the shock wave. In addition, p +1 is the large light-cone momentum of theparton in the shock wave generating those sub-eikonal fields. In Eq. (11) we have also employed the adjoint light-coneWilson line U x [ b − , a − ] = P exp ig b − (cid:90) a − dx − A + ( x + = 0 , x − , x ) . (12)The matrices t a and t b are the fundamental generators of SU( N c ) with the indices a, b = 1 , . . . , N c − N c & N f helicity evolution equations we are about to study mix the polarized quark dipole amplitude (9)with the polarized gluon dipole amplitude [41, 42] G ( z ) = 12( N c −
1) Re (cid:68) T Tr (cid:104) U U pol † (cid:105) + T Tr (cid:104) U pol U † (cid:105)(cid:69) (13)defined in terms of the adjoint polarized Wilson lines [42]( U polx ) ab = 2 i g p +1 + ∞ (cid:90) −∞ dx − (cid:0) U x [+ ∞ , x − ] F ( x + = 0 , x − , x ) U x [ x − , −∞ ] (cid:1) ab (14) − g p +1 ∞ (cid:90) −∞ dx − ∞ (cid:90) x − dx − U aa (cid:48) x [+ ∞ , x − ] ¯ ψ ( x − , x ) t a (cid:48) V x [ x − , x − ] 12 γ + γ t b (cid:48) ψ ( x − , x ) U b (cid:48) bx [ x − , −∞ ] − c.c.. We would like to stress that while the polarized quark dipole amplitude Q is related to the flavor-singlet quark hPDF ∆Σ, the polarizedgluon helicity amplitude G is not related to the gluon hPDF ∆ G . The distribution ∆ G at small- x is related to a different type of thepolarized dipole amplitude, as detailed in [41]. Q ( z ) = x x x x G ( z ) = x x FIG. 1: Diagrammatic representation of the quark ( Q ) and gluon ( G ) polarized dipole amplitudes. The shadedrectangle is the shock wave, while the square represents insertions of sub-eikonal operator(s).(Note again an additional s in the normalization of Eq. (14) as compared to [42]; in addition, G ( z ) in Eq. (13) waslabeled G adj in [42].) The impact parameter-integrated polarized gluon dipole amplitude is defined similar to Eq. (8) G ( x , z ) = (cid:90) d (cid:18) x + x (cid:19) G ( z ) . (15)The polarized quark and gluon dipole amplitudes are depicted diagrammatically in Fig. 1. The polarized proton,along with all the partons to be produced by the small- x evolution, are depicted by the rectangular shock wave. Thesquare on one of the lines depicts insertions of sub-eikonal operators present in Eqs. (11) and (14).In the limit of large numbers of flavors, N f , and colors, N c , the small- x DLA evolution equations for the polarizedamplitudes Q and G are [1, 42] G ( x , z ) = G (0) ( x , z )+ α s N c π z (cid:90) max { Λ , /x } /s dz (cid:48) z (cid:48) x (cid:90) / ( z (cid:48) s ) dx x (cid:2) Γ( x , x , z (cid:48) ) + 3 G ( x , z (cid:48) ) (cid:3) (16a) − α s N f π z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) / ( z (cid:48) s ) dx x Γ gen ( x , x , z (cid:48) ) ,Q ( x , z ) = Q (0) ( x , z )+ α s N c π z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) / ( z (cid:48) s ) dx x Q ( x , z (cid:48) ) (16b)+ α s N c π z (cid:90) max { Λ , /x } /s dz (cid:48) z (cid:48) x (cid:90) / ( z (cid:48) s ) dx x (cid:20) G ( x , z (cid:48) ) + 12 Γ( x , x , z (cid:48) ) + Q ( x , z (cid:48) ) − Γ( x , x , z (cid:48) ) (cid:21) , Γ( x , x , z (cid:48) ) = G (0) ( x , z )+ α s N c π z (cid:48) (cid:90) max { Λ , /x } /s dz (cid:48)(cid:48) z (cid:48)(cid:48) min { x ,x z (cid:48) /z (cid:48)(cid:48) } (cid:90) / ( z (cid:48)(cid:48) s ) dx x (cid:2) Γ( x , x , z (cid:48)(cid:48) ) + 3 G ( x , z (cid:48)(cid:48) ) (cid:3) (16c) − α s N f π z (cid:48) (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) / ( z (cid:48)(cid:48) s ) dx x Γ gen ( x , x , z (cid:48)(cid:48) ) , Γ( x , x , z (cid:48) ) = Q (0) ( x , z )+ α s N c π z (cid:48) (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) / ( z (cid:48)(cid:48) s ) dx x Q ( x , z (cid:48)(cid:48) ) (16d)+ α s N c π z (cid:48) (cid:90) max { Λ , /x } /s dz (cid:48)(cid:48) z (cid:48)(cid:48) min { x ,x z (cid:48) /z (cid:48)(cid:48) } (cid:90) / ( z (cid:48)(cid:48) s ) dx x (cid:20) G ( x , z (cid:48)(cid:48) ) + 12 Γ( x , x , z (cid:48)(cid:48) ) + Q ( x , z (cid:48)(cid:48) ) − Γ( x , x , z (cid:48)(cid:48) ) (cid:21) , where the generalized quark dipole amplitude is defined as [41, 42]Γ gen ( x , x , z (cid:48) ) = θ ( x − x )Γ( x , x , z (cid:48) ) + θ ( x − x ) Q ( x , z (cid:48) ) . (17)To properly impose the light-cone time ordering, the evolution equations (16) employ two auxiliary functions, thegluon and quark polarized “neighbor” dipole amplitudes, Γ( x , x , z (cid:48) ) for G ( x , z ) and Γ( x , x , z (cid:48) ) for Q ( x , z )[1, 40–42]. The operator definitions for Γ( x , x , z (cid:48) ) and Γ( x , x , z (cid:48) ) are given by the same Eqs. (9) and (13),respectively (integrated over all impact parameters). However, implicit in those definitions is the x − -light-cone life-time of the dipole, proportional to z x [1, 43, 54, 55]. For the “neighbor” dipole amplitudes Γ( x , x , z (cid:48) ) andΓ( x , x , z (cid:48) ) the life-time variable is proportional to z (cid:48) x , different from z (cid:48) x one would expect based on theirdipole size x .As follows from Eq. (7), we only need to find the amplitude Q ( x , z ) in order to construct the quark helicity PDF.However, the evolution equations (16) mix all four involved amplitudes, and have to be solved for all the amplitudesin order to obtain Q ( x , z ).In this paper we are chiefly interested in the asymptotic behaviors of the dipole amplitudes at high energies. In [2]it was shown that the high-energy asymptotics of helicity amplitudes at large- N c is largely independent of the initialconditions/inhomogeneous terms in the evolution, which only affect the overall prefactor in the asymptotic expression.(In addition, it is well-known that the small- x asymptotics of the Balitsky–Fadin–Kuraev–Lipatov (BFKL) evolution[56, 57] depends on the initial conditions only in its prefactor as well.) Inspired by these examples of independence ofthe initial conditions, we assume it to be the case in the large- N c & N f approximation as well, and, for convenience,we take the inhomogeneous terms in Eqs. (16) to be G (0) ( x , z ) = Q (0) ( x , z ) = 1 , (18)instead of using the proper Born-level initial condition [2].In the regime of our interest, z (cid:48) s (cid:28) x (cid:28) x , the definition (17) implies that z (cid:48) (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) / ( z (cid:48)(cid:48) s ) dx x Γ gen ( x , x , z (cid:48)(cid:48) ) (19a)= z (cid:48) (cid:90) max { Λ , /x } /s dz (cid:48)(cid:48) z (cid:48)(cid:48) min { x ,x z (cid:48) /z (cid:48)(cid:48) } (cid:90) / ( z (cid:48)(cid:48) s ) dx x Γ( x , x , z (cid:48)(cid:48) ) + x z (cid:48) /x (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) max { x , / ( z (cid:48)(cid:48) s ) } dx x Q ( x , z (cid:48)(cid:48) ) , z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) / ( z (cid:48) s ) dx x Γ gen ( x , x , z (cid:48) ) (19b)= z (cid:90) max { Λ , /x } /s dz (cid:48) z (cid:48) x (cid:90) / ( z (cid:48) s ) dx x Γ( x , x , z (cid:48) ) + z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) max { x , / ( z (cid:48) s ) } dx x Q ( x , z (cid:48) ) , allowing us to rewrite the evolution equations (16) as G ( x , z ) = 1+ α s N c π z (cid:90) max { Λ , /x } /s dz (cid:48) z (cid:48) x (cid:90) / ( z (cid:48) s ) dx x (cid:20) Γ( x , x , z (cid:48) ) + 3 G ( x , z (cid:48) ) − N f N c Γ( x , x , z (cid:48) ) (cid:21) (20a) − α s N f π z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) max { x , / ( z (cid:48) s ) } dx x Q ( x , z (cid:48) ) ,Q ( x , z ) = 1+ α s N c π z (cid:90) Λ /s dz (cid:48) z (cid:48) x z/z (cid:48) (cid:90) / ( z (cid:48) s ) dx x Q ( x , z (cid:48) ) (20b)+ α s N c π z (cid:90) max { Λ , /x } /s dz (cid:48) z (cid:48) x (cid:90) / ( z (cid:48) s ) dx x (cid:20) G ( x , z (cid:48) ) + 12 Γ( x , x , z (cid:48) ) + Q ( x , z (cid:48) ) − Γ( x , x , z (cid:48) ) (cid:21) , Γ( x , x , z (cid:48) ) = 1 − α s N f π x z (cid:48) /x (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) max { x , / ( z (cid:48)(cid:48) s ) } dx x Q ( x , z (cid:48)(cid:48) ) (20c)+ α s N c π z (cid:48) (cid:90) max { Λ , /x } /s dz (cid:48)(cid:48) z (cid:48)(cid:48) min { x ,x z (cid:48) /z (cid:48)(cid:48) } (cid:90) / ( z (cid:48)(cid:48) s ) dx x (cid:20) Γ( x , x , z (cid:48)(cid:48) ) + 3 G ( x , z (cid:48)(cid:48) ) − N f N c Γ( x , x , z (cid:48)(cid:48) ) (cid:21) , Γ( x , x , z (cid:48) ) = 1+ α s N c π z (cid:48) (cid:90) Λ /s dz (cid:48)(cid:48) z (cid:48)(cid:48) x z (cid:48) /z (cid:48)(cid:48) (cid:90) / ( z (cid:48)(cid:48) s ) dx x Q ( x , z (cid:48)(cid:48) ) (20d)+ α s N c π z (cid:48) (cid:90) max { Λ , /x } /s dz (cid:48)(cid:48) z (cid:48)(cid:48) min { x ,x z (cid:48) /z (cid:48)(cid:48) } (cid:90) / ( z (cid:48)(cid:48) s ) dx x (cid:20) G ( x , z (cid:48)(cid:48) ) + 12 Γ( x , x , z (cid:48)(cid:48) ) + Q ( x , z (cid:48)(cid:48) ) − Γ( x , x , z (cid:48)(cid:48) ) (cid:21) . Our aim is to solve this system of integral equations numerically. Its solution for Q ( x , z ), via Eq. (7), would give usthe quark helicity distribution at small Bjorken x and in the regime where both the gluon and quark contributions toevolution are relevant (due to the large- N c & N f limit employed in deriving Eqs. (20)). III. NUMERICAL SOLUTION: DISCRETIZATION AND ALGORITHM
To numerically evaluate the integrals in Eqs. (20), we first make the following changes of variables η ( n ) = (cid:114) α s N c π ln z ( n ) s Λ , s kl = (cid:114) α s N c π ln 1 x kl Λ , (21)where z ( n ) can be z , z (cid:48) or z (cid:48)(cid:48) which relates to η , η (cid:48) or η (cid:48)(cid:48) respectively. In addition, kl = 10 ,
21 or 32. In terms of thesenew variables, Eqs. (20) can be re-written as G ( s , η ) = 1 + η (cid:90) max { ,s } dη (cid:48) η (cid:48) (cid:90) s ds (cid:20) Γ( s , s , η (cid:48) ) + 3 G ( s , η (cid:48) ) − N f N c Γ( s , s , η (cid:48) ) (cid:21) (22a) − N f N c η (cid:90) dη (cid:48) min { s ,η (cid:48) } (cid:90) s + η (cid:48) − η ds Q ( s , η (cid:48) ) ,Q ( s , η ) = 1 + η (cid:90) max { ,s } dη (cid:48) η (cid:48) (cid:90) s ds (cid:20) G ( s , η (cid:48) ) + 12 Γ( s , s , η (cid:48) ) + Q ( s , η (cid:48) ) − Γ( s , s , η (cid:48) ) (cid:21) (22b)+ 12 η (cid:90) dη (cid:48) η (cid:48) (cid:90) s + η (cid:48) − η ds Q ( s , η (cid:48) ) , Γ( s , s , η (cid:48) ) = 1 + η (cid:48) (cid:90) max { ,s } dη (cid:48)(cid:48) η (cid:48)(cid:48) (cid:90) max { s ,s + η (cid:48)(cid:48) − η (cid:48) } ds (cid:20) Γ( s , s , η (cid:48)(cid:48) ) + 3 G ( s , η (cid:48)(cid:48) ) − N f N c Γ( s , s , η (cid:48)(cid:48) ) (cid:21) (22c) − N f N c η (cid:48) + s − s (cid:90) dη (cid:48)(cid:48) min { s ,η (cid:48)(cid:48) } (cid:90) s + η (cid:48)(cid:48) − η (cid:48) ds Q ( s , η (cid:48)(cid:48) ) , Γ( s , s , η (cid:48) ) = 1 + η (cid:48) (cid:90) max { ,s } dη (cid:48)(cid:48) η (cid:48)(cid:48) (cid:90) max { s ,s + η (cid:48)(cid:48) − η (cid:48) } ds (cid:20) G ( s , η (cid:48)(cid:48) ) + 12 Γ( s , s , η (cid:48)(cid:48) ) + Q ( s , η (cid:48)(cid:48) ) − Γ( s , s , η (cid:48)(cid:48) ) (cid:21) + 12 η (cid:48) (cid:90) dη (cid:48)(cid:48) η (cid:48)(cid:48) (cid:90) s + η (cid:48)(cid:48) − η (cid:48) ds Q ( s , η (cid:48)(cid:48) ) . (22d)Note that the integrals in equations (22) include the regions with s < , s <
0. This means that Λ is not aninfrared (IR) cutoff, but a perturbative transverse momentum scale characterizing the proton at the start of ourevolution.Next, we discretize the resulting equations (22). In particular, let i ( n ) and j ( n ) be the discretized version of s kl and η ( n ) , with step sizes ∆ s and ∆ η , respectively, such that η ( n ) = j ( n ) ∆ η ≡ η j ( n ) along with s = i ∆ s ≡ s i , s = i (cid:48) ∆ s ≡ s i (cid:48) , and s = i (cid:48)(cid:48) ∆ s ≡ s i (cid:48)(cid:48) . The discretized amplitudes are defined as G ij ≡ G ( s i , η j ), Q ij ≡ Q ( s i , η j ),Γ ikj ≡ Γ( s i , s k , η j ), and Γ ikj ≡ Γ( s i , s k , η j ). The equations (22) are self-contained over the following region in ( s , η )-plane: η ∈ [0 , η max ], η − η max ≤ s ≤ η , where η max is some arbitrary positive upper value of the η -range. This isthe region where we will solve them numerically. Discretized Eqs. (22) are G ij = 1 + ∆ η ∆ s j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) = i (cid:20) Γ ii (cid:48) j (cid:48) + 3 G i (cid:48) j (cid:48) − N f N c Γ ii (cid:48) j (cid:48) (cid:21) − N f N c j − (cid:88) j (cid:48) =0 min { i,j (cid:48) } (cid:88) i (cid:48) = i + j (cid:48) − j Q i (cid:48) j (cid:48) , (23a) Q ij = 1 + ∆ η ∆ s j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) = i (cid:20) G i (cid:48) j (cid:48) + 12 Γ ii (cid:48) j (cid:48) + Q i (cid:48) j (cid:48) − Γ ii (cid:48) j (cid:48) (cid:21) + 12 j − (cid:88) j (cid:48) =0 j (cid:48) (cid:88) i (cid:48) = i + j (cid:48) − j Q i (cid:48) j (cid:48) , (23b)Γ ikj = 1 + ∆ η ∆ s j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) =max { i,k + j (cid:48) − j } (cid:20) Γ ii (cid:48) j (cid:48) + 3 G i (cid:48) j (cid:48) − N f N c Γ ii (cid:48) j (cid:48) (cid:21) − N f N c i + j − k − (cid:88) j (cid:48) =0 min { i,j (cid:48) } (cid:88) i (cid:48) = k + j (cid:48) − j Q i (cid:48) j (cid:48) , (23c)Γ ikj = 1 + ∆ η ∆ s j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) =max { i,k + j (cid:48) − j } (cid:20) G i (cid:48) j (cid:48) + 12 Γ ii (cid:48) j (cid:48) + Q i (cid:48) j (cid:48) − Γ ii (cid:48) j (cid:48) (cid:21) + 12 j − (cid:88) j (cid:48) =0 j (cid:48) (cid:88) i (cid:48) = k + j (cid:48) − j Q i (cid:48) j (cid:48) . (23d)In principle, Eqs. (23) can already be solved numerically, similar to [2, 41]. Instead we will simplify these equationswhich will allow us to implement a faster algorithm in the numerical solution. We replace j by j − G ij , this replacement gives G i ( j − = 1 + ∆ η ∆ s j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) = i (cid:20) Γ ii (cid:48) j (cid:48) + 3 G i (cid:48) j (cid:48) − N f N c Γ ii (cid:48) j (cid:48) (cid:21) − N f N c j − (cid:88) j (cid:48) =0 min { i,j (cid:48) } (cid:88) i (cid:48) = i + j (cid:48) − j +1 Q i (cid:48) j (cid:48) . (24)Comparing Eq. (23a) to Eq. (24), we obtain G ij = G i ( j − + ∆ η ∆ s j − (cid:88) i (cid:48) = i (cid:20) Γ ii (cid:48) ( j − + 3 G i (cid:48) ( j − − N f N c Γ ii (cid:48) ( j − (cid:21) − N f N c j − (cid:88) j (cid:48) =0 Q ( i + j (cid:48) − j ) j (cid:48) − N f N c Q i ( j − . (25)Since zs ≥ /x , we have η ≥ s , which leads to j ≥ i . For i = j we have Q j ( j − = 1 in the last term of Eq. (25):this is determined by the initial conditions. Including this additional contribution from a single point in the i, j griddoes not significantly affect the numerical solution.Writing each equation in (23) in this recursive form allows for a numerical evaluation with one fewer layer of loops,resulting in much shorter computation time for smaller step sizes, ∆ η and ∆ s , and for a larger η -range, defined by η ∈ [0 , η max ]. In order to write down the recursive equations similar to Eq (25) for Q ij , Γ ikj and Γ ikj , notice that j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) =max { i,k + j (cid:48) − j } Γ ii (cid:48) j (cid:48) − j − (cid:88) j (cid:48) =max { ,i } j (cid:48) (cid:88) i (cid:48) =max { i,k + j (cid:48) − j +1 } Γ ii (cid:48) j (cid:48) = j − (cid:88) j (cid:48) =max { i + j − k, } Γ i ( k + j (cid:48) − j ) j (cid:48) + j − (cid:88) i (cid:48) = k Γ ii (cid:48) ( j − , (26)where we have used the fact that i < k < j , a consequence of the regime in which Γ (and Γ) are defined and employed,1 /z (cid:48) s (cid:28) x (cid:28) x , or, equivalently, η (cid:48) (cid:29) s (cid:29) s . Equations (23) and (26) imply that (again for i < k < j ) Q ij = Q i ( j − + ∆ η ∆ s j − (cid:88) i (cid:48) = i (cid:20) G i (cid:48) ( j − + 12 Γ ii (cid:48) ( j − + 32 Q i (cid:48) ( j − − Γ ii (cid:48) ( j − (cid:21) + 12 j − (cid:88) j (cid:48) =0 Q ( i + j (cid:48) − j ) j (cid:48) , (27a)Γ ikj = Γ ik ( j − + ∆ η ∆ s (cid:34) j − (cid:88) i (cid:48) = k (cid:20) Γ ii (cid:48) ( j − + 3 G i (cid:48) ( j − − N f N c Γ ii (cid:48) ( j − (cid:21) − N f N c i + j − k − (cid:88) j (cid:48) =0 Q ( k + j (cid:48) − j ) j (cid:48) − N f N c Q i ( i + j − k − + j − (cid:88) j (cid:48) =max { i + j − k, } (cid:20) Γ i ( k + j (cid:48) − j ) j (cid:48) + 3 G ( k + j (cid:48) − j ) j (cid:48) − N f N c Γ i ( k + j (cid:48) − j ) j (cid:48) (cid:21) (cid:35) , (27b)Γ ikj = Γ ik ( j − + ∆ η ∆ s (cid:34) j − (cid:88) i (cid:48) = k (cid:20) G i (cid:48) ( j − + 12 Γ ii (cid:48) ( j − + 32 Q i (cid:48) ( j − − Γ ii (cid:48) ( j − (cid:21) + 12 j − (cid:88) j (cid:48) =0 Q ( k + j (cid:48) − j ) j (cid:48) + j − (cid:88) j (cid:48) =max { ,i + j − k } (cid:20) G ( k + j (cid:48) − j ) j (cid:48) + 12 Γ i ( k + j (cid:48) − j ) j (cid:48) + Q ( k + j (cid:48) − j ) j (cid:48) − Γ i ( k + j (cid:48) − j ) j (cid:48) (cid:21) (cid:35) . (27c)We solve Eqs. (25) and (27) numerically in steps along the η -axis. To obtain the values for 0 ≤ η ≤ η max , westart from η = ∆ η , which is equivalent to j = 1, at which we use Eqs. (25) and (27) to determine G i and Q i for1 − η max ∆ η ≤ i ≤
1, together with Γ ik and Γ ik for 1 − η max ∆ η ≤ i ≤ k ≤
1, assuming that G i = Q i = Γ ik = Γ ik = 1are determined by the inhomogeneous term in Eqs. (23). Afterward, we repeat the same steps by applying Eqs. (25)and (27) for j = 2, then j = 3, and so on, until j max = η max ∆ η . At each j , we compute the dipole amplitudes G ij and Q ij for i in the range j − η max ∆ η ≤ i ≤ j and also find the amplitudes Γ ikj and Γ ikj for all pairs of i and k satisfying j − η max ∆ η ≤ i ≤ k ≤ j . This process determines the values of Q ( s , η ) in the η ∈ [0 , η max ], η − η max ≤ s ≤ η region, which we will use below to determine the high energy asymptotics of the polarized quark dipole amplitude.We always take ∆ s = ∆ η for simplicity. IV. NUMERICAL SOLUTION: RESULTS
In this Section we present the results of our numerical solution of the large- N c & N f evolution equations (22) forthe quark and gluon polarized dipole amplitudes Q ( s , η ) and G ( s , η ), along with the corresponding quark helicityPDF ∆Σ( x, Q ). The plots of ln | G ( s , η ) | and ln | Q ( s , η ) | in the ( s , η )-plane are shown in Fig. 2 for N f = 2 , , N c = 3.The plots in Fig. 2 demonstrate an approximately linear rise of ln | G ( s , η ) | and ln | Q ( s , η ) | with η , similar tothe large- N c case studied previously in [2, 3]. The only difference is that now, in Fig. 2, the rise of those functions isnot monotonic, and appears to be periodically interrupted by lines of sharp local minima.To illustrate the origin of this non-monotonicity, we plot sgn[ G (0 , η )] ln | G (0 , η ) | and sgn[ Q (0 , η )] ln | Q (0 , η ) | asfunctions of η in Fig. 3 for N f = 3. From these plots we see that G (0 , η ) and Q (0 , η ) oscillate with η . The oscillationsexplain the non-monotonic behavior we saw in Fig. 2. These oscillations are the main qualitative difference of thesmall- x asymptotics for the quark helicity distribution in the large- N c & N f limit, as compared to the large- N c case.It appears that introducing the quarks back into the helicity evolution equations generates oscillations. While theabsolute values of Q ( s , η ) and G ( s , η ) still grow exponentially with η , both dipole amplitudes also oscillate.Moreover, from Fig. 2 one can see that the oscillation period appears to get smaller (and the oscillation frequencyappears to get larger) with increasing N f , which is also consistent with the fact that these oscillations are absent inthe gluon-only large- N c equations solved in [2, 3].While an analytic solution of Eqs. (22) is beyond the scope of this work, we can try to find analytic formulasapproximating our numerical results in Fig. 3, at least in the large- η asymptotics. Combining the oscillations withthe exponential growth of the maxima of | Q (0 , η ) | and | G (0 , η ) | with η we propose the following asymptotic forms forthe polarized dipole amplitudes: G (0 , η ) ∼ e α G η cos ( ω G η + ϕ G ) , (28a) Q (0 , η ) ∼ e α Q η cos ( ω Q η + ϕ Q ) . (28b)The oscillation frequencies are denoted by ω G and ω Q , while the initial phases are denoted by ϕ G and ϕ Q . As onecan already see from Fig. 2, both the frequencies and the initial phases depend on N f . This is confirmed by thedetailed analysis of our numerical solution in Appendix A. Furthermore, the amplitudes of oscillations in Q (0 , η ) and G (0 , η ) grow exponentially with η , while the exponents α G and α Q ( the intercepts in Regge terminology) also appearto depend on N f . The results of the analysis carried out in Appendix A, where we fit our numerical solution withthe ansatz (28) and extract the corresponding frequencies, phases, and intercepts, are summarized in Table I.From Table I we see that, for each N f , the two frequencies are equal within the numerical accuracy, ω G = ω Q . Thefrequencies increase with N f : while the exact analytic dependence of ω G = ω Q on N f is yet to be determined, weconstruct a Pad´e approximant to write ω Q = ω G ≈ . N f . N f ≈ . πN f . N f . (29)0 (a) ln | G ( s , η ) | at N f = 2 (b) ln | Q ( s , η ) | at N f = 2(c) ln | G ( s , η ) | at N f = 3 (d) ln | Q ( s , η ) | at N f = 3(e) ln | G ( s , η ) | at N f = 6 (f) ln | Q ( s , η ) | at N f = 6 FIG. 2: Plots of ln | G ( s , η ) | and ln | Q ( s , η ) | for N f = 2 , , N c = 3. All the graphs result from numericalcomputations with the step size ∆ η = 0 .
075 and η max = 30. (a) sgn[ G (0 , η )] ln | G (0 , η ) | (b) sgn[ Q (0 , η )] ln | Q (0 , η ) | FIG. 3: Plots of sgn[ G (0 , η )] ln | G (0 , η ) | and sgn[ Q (0 , η )] ln | Q (0 , η ) | versus η for N f = 3 and N c = 3. Both graphsresult from numerical computations with the step size ∆ η = 0 .
075 and η max = 30.The intercepts, α G and α Q , of the polarized dipole amplitudes’ exponential growth given in Table I, are also equalto each other for each N f with the precision of our numerical solution, α Q = α G . For all N f studied they remainclose to α qh ( N f = 0) = √ ≈ . G (0 , η ) in the large- N c pure-glue limit with N f = 0 (inunits of (cid:112) α s N c / (2 π )), derived analytically in [3]. It appears, though, that α Q = α G is decreasing slowly with N f .The initial phase ϕ in Table I is always between 0 and π for G (0 , η ) and between − π and 0 for Q (0 , η ). The values1 N f G (0 , η ) Q (0 , η ) α G ω G ϕ G α Q ω Q ϕ Q . ± .
003 0 . ± .
001 0 . ± .
011 2 . ± .
001 0 . ± . − . ± . . ± .
004 0 . ± .
001 0 . ± .
008 2 . ± .
006 0 . ± . − . ± . . ± .
002 0 . ± .
004 0 . ± .
009 2 . ± .
002 0 . ± . − . ± . TABLE I: The resulting intercept α , frequency ω , and initial phase ϕ , for the gluon and quark polarized dipoleamplitudes G (0 , η ) and Q (0 , η ) at N f = 2 , , N c = 3.of ϕ Q and ϕ G do not display a clear relation with each other. Their functional dependence on N f is non-monotonic,and its form is also not obvious from Table I. In fact, the initial phase depends greatly on the choice of the initialcondition (the inhomogeneous term) in Eqs. (16). For instance, if one performs a similar computation using theBorn-level-inspired dipole amplitudes as initial conditions (cf. [2]), while still taking the inhomogeneous terms for G and Q to be equal for simplicity, G (0) ( s , η ) = Q (0) ( s , η ) = α s C F π N c ( C F η − η − s )) , (30)for N f = 3, N c = 3, α s = 0 .
35, ∆ η = 0 . η max = 20, and the fundamental Casimir operator C F = ( N c − / N c ofSU( N c ), one finds ϕ G = − .
146 and ϕ Q = 1 . ϕ G = − .
146 and ϕ Q = 1 .
530 and thephases in Table I is too large to be attributed to discretization errors. We thus conclude that the phases ϕ Q and ϕ G are indeed very much dependent on the initial conditions to the evolution, and their values listed in Table I are notuniversal.Next let us determine what our solution implies for the quark hPDF ∆Σ. Rewriting Eq. (7) in terms of η and s from Eq. (21) we arrive at∆Σ( x, Q ) = N f α s π √ αs Nc π ln Q x Λ2 (cid:90) dη η (cid:90) η − √ αs Nc π ln x ds Q ( s , η ) . (31)To determine the small- x asymptotics of ∆Σ we need to evaluate the integral in Eq. (31). Note that for Q ≥ Λ theintegration region in Eq. (31) lies within the area η ∈ [0 , η max ], η − η max ≤ s ≤ η where we found the numericalsolution for Q ( s , η ) if we choose η max = (cid:113) α s N c π ln Q x Λ . Since Q ( s , η ) is known numerically, we perform a numericalintegration to obtain the values of ∆Σ( x, Q ) as a function of Bjorken x at fixed Q = 10 GeV for the case of 3quark flavors, N f = 3, while choosing, for simplicity, Λ = Q , such that η max = (cid:113) α s N c π ln x .FIG. 4: Plot of sgn[∆Σ] ln | ∆Σ | versus x resulting from the numerical integration of Eq. (31) for N f = 3, α s = 0 . Q = Λ , and with the step size ∆ η = 0 . x, Q ), similar to the polarized dipole amplitudes Q and G , is an oscillating function of ln(1 /x ),with the oscillation amplitude growing exponentially with ln(1 /x ). In Fig. 4 we plot sgn[∆Σ( x, Q )] ln | ∆Σ( x, Q ) | as a function of x , demonstrating the oscillations explicitly. Inspired by the success of the ansatz (28) for the dipoleamplitude, we propose the following ansatz for the small- x asymptotics of ∆Σ( x, Q ),∆Σ( x, Q = 10 GeV ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c & N f ∼ (cid:18) x (cid:19) α qh cos (cid:20) ω q ln 1 x + ϕ q (cid:21) . (32)The parameters α qh , ω q and ϕ q are extracted from our numerical results in Appendix B by applying the parameterfitting process outlined in Appendix A. We use α s (10 GeV ) ≈ .
25. This gives (again, for N f = 3) α qh = (2 . ± . (cid:114) α s N c π , ω q = (0 . ± . (cid:114) α s N c π , and ϕ q = − . ± . . (33)The intercept and frequency are within the margins of error from those for the dipole amplitudes Q and G at N f = 3,multiplied by the factor of (cid:113) α s N c π : α qh = α Q (cid:113) α s N c π = α G (cid:113) α s N c π , ω q = ω Q (cid:113) α s N c π = ω G (cid:113) α s N c π . This shows that theintercept and frequency of the dipole amplitude Q determine the intercept and frequency of ∆Σ: both are uniquelydetermined by the evolution. (Note that one cannot simply substitute Eq. (28b) into Eq. (31) to obtain this resultanalytically, since the former is valid only at s = 0, while the latter has an integral over a range of non-zero valuesof s .) However, as remarked previously for the dipole amplitudes, the initial phase, ϕ q , is, in fact, a by-product ofour choice of initial conditions, G (0) and Q (0) . In practical phenomenological applications, G (0) and Q (0) , and, hence,the initial phase ϕ q , have to be determined from the data. V. QUARK HELICITY: AN ESTIMATE
In this Section we follow the strategy employed in [2] to estimate the possible impact of our new functional formfor the small- x asymptotics of ∆Σ( x, Q ) (32) on the amount of the proton spin carried by the small- x quarks. Ouranalysis below should be understood as a rough estimate, with the anticipation of a more detailed phenomenology tobe done in the future.As the asymptotic form (32) only holds for small x , we use that expression to extrapolate the quark helicitydistribution from the DSSV14 result in [6], starting at a particular connecting point, x (cid:28)
1, into the small- x region, x min ≤ x ≤ x . In particular, we write the small- x quark helicity distribution for x < x as∆Σ( x, Q = 10 GeV ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c & N f = K (cid:16) x x (cid:17) α qh cos (cid:104) ω q ln x x + ϕ q (cid:105) , (34)with α qh and ω q given in Eq. (33). For x > x we take the quark hPDF to be given by the DSSV14 parameterization[6], ∆Σ DSSV ( x, Q ). We, therefore, need to match our ∆Σ( x, Q = 10 GeV ) from Eq. (34) onto ∆Σ DSSV ( x, Q ) at x = x .Unlike [2], where the ansatz for ∆Σ at small x contained only the power of 1 /x , we now have the function inEq. (34) with two unknown parameters, the overall normalization factor K and the phase ϕ q . We assume thata more complete phenomenological approach would be able to uniquely determine ϕ q from the large- x ( x > x )data. We further assume that the value of ϕ q obtained from a more complete approach would generate a smoothmatching of ∆Σ at x = x , ensuring continuity of both ∆Σ( x, Q ) and its derivative ∂ ∆Σ( x, Q ) /∂x at x = x .This assumption is not very reliable for our purposes, since the functional form (34) is asymptotic, and is, therefore,valid only for x (cid:28) x : using it to ensure the continuity of ∆Σ( x, Q ) and ∂ ∆Σ( x, Q ) /∂x at x = x is somewhatquestionable. Strictly-speaking we have to admit that ϕ q is an almost arbitrary parameter, whose value appears tobe very important for assessing the amount of quark helicity at small x . While the more detailed phenomenologyshould better constraint the allowed ranges of K and ϕ q , we will fix them here by simply requiring the continuity of∆Σ( x, Q ) and ∂ ∆Σ( x, Q ) /∂x at x = x between our asymptotics (34) and ∆Σ DSSV ( x, Q ) from [6].We perform this computation for x = 0 .
01 and x = 0 . x . The resulting x ∆Σ( x, Q = 10 GeV ) is plotted versus x in Fig. 5, which depicts the original DSSV14 [6] curve, extrapolated by apower-law in x to very small x , along with the two curves resulting from matching our asymptotics (34) to DSSV14at x = 0 .
01 and x = 0 . x ∆Σ( x, Q = 10 GeV ) versus x . The thin (red) line depicts the results from DSSV14 [6], whilethe thick (green) and medium-thick (blue) lines show the results of matching our asymptotics (34) to DSSV14 at x = 0 .
01 and x = 0 . S q ( Q = 10 GeV ), we consider the integralsimilar to Eq. (3) but with the lower limit at some small finite x min , such that 0 < x min < x (cf. [2])∆Σ [ x min ] ( Q ) = (cid:90) x min dx ∆Σ( x, Q ) . (35)From Eqs. (3) and (35) we see that S q ( Q ) = (1 /
2) ∆Σ [ x min ] ( Q ) in the x min → [ x min ] ( Q = 10 GeV ) versus x min for the DSSV14 parameterization [6] along with the two curves resulting frommatching our (34) to DSSV14 at x = 0 .
01 and x = 0 . x ∆Σ obtained from thelarge- N c pure-glue evolution [2, 3], and adjusting its overall normalization to ensure the continuity of ∆Σ( x, Q ) withDSSV14 at x = 0 .
01 and x = 0 . [ x min ] ( Q = 10 GeV ) defined in Eq. (35) versus x min . The thin (red) line is DSSV14 [6], whilethe thick (green) and medium-thick (blue) solid lines show the results of matching our asymptotics (34) to DSSV14at x = 0 .
01 and x = 0 . x extrapolation done in [2] usinga power-law expression for ∆Σ( x, Q ) at small x resulting from the large- N c evolution, also matched onto DSSV14 at x = 0 .
01 (thick dashed, green) and x = 0 .
001 (medium-thick dashed, blue).Figure 6 illustrates the potential of small- x evolution to significantly affect the amount of the quark spin content inthe proton, which was already observed in [2]. We conclude from Fig. 6 that the small- x extrapolation of ∆Σ( x, Q )at Q = 10 GeV varies significantly with the matching point, x . For the larger x -value, x = 0 .
01, we also see astrong variation between using large- N c (KPS16) and large- N c & N f small- x evolution for helicity. It is also interesting4to see in Fig. 6 that for x = 0 .
001 both the KPS16 curve and the curve based on Eq. (34) are rather close to theDSSV14 curve and to each other. Again, let us stress that our use of the asymptotic expression (34), valid for x (cid:28) x ,to perform the matching of both ∆Σ( x, Q ) and ∂ ∆Σ( x, Q ) /∂x at x = x , which is probably outside of its region ofapplicability, gives us only a rough estimate of quark hPDF at small x . We expect that future more detailed studieswould place this matching under a more solid theoretical control. VI. CONCLUSIONS AND DISCUSSION
In this work, we have numerically computed the asymptotic high-energy behavior of the polarized quark dipoleamplitude resulting from the double-logarithmic small- x helicity evolution [1, 42] in the limit of large number of quarkcolors and flavors. The obtained amplitude Q is plotted in Figs. 2 and 3 as a function of η and s defined in Eq. (21).The amplitude Q in the large- N c & N f limit displays an oscillatory pattern as a function of the center-of-mass energyand of the dipole transverse size, on top of the exponential growth seen before in the large- N c limit with N f = 0 [2, 3],when all the evolution was gluon-driven. We conclude that these oscillations appear after including quarks back intothe small- x helicity evolution of [1, 42].Our numerical results for Q are well-approximated by Eq. (28b). The frequency ω Q and initial phase ϕ Q of theoscillations depend on the number of flavors, with this dependence summarized in Table I and Eq. (29). The intercept α Q , also given in Table I, exhibits weak dependence on N f , largely staying close to the large- N c , N f = 0 value α qh ( N f = 0) = √ ≈ . η ). Let us also note that, while the periodic oscillations withenergy we describe here are a new result for helicity distributions, a single sign reversal with decreasing x was predictedfor these quantities in [31] using a different approach.The polarized quark dipole amplitude oscillations result in the similar oscillating behavior of ∆Σ as a function ofln(1 /x ) with the oscillation amplitude growing as a power of 1 /x , as shown in Eq. (1), which we reproduce here∆Σ( x, Q ) (cid:12)(cid:12)(cid:12)(cid:12) large- N c & N f ∼ (cid:18) x (cid:19) α qh cos (cid:20) ω q ln (cid:18) x (cid:19) + ϕ q (cid:21) . (36)This is the main result of this work. It is illustrated in Fig. 4. The intercept, α qh , and the frequency, ω q , given byEq. (33), are within the margin of error from the corresponding parameters, α Q and ω Q , for the dipole amplitude, Q .The initial phase, ϕ q , on the other hand, differs greatly from ϕ Q , due to the relation (31) between the two quantities.In general, the phase ϕ Q , and, therefore, the phase ϕ q , exhibit a very strong dependence on the initial conditions(inhomogeneous terms) for our large- N c & N f evolution. This is in contrast to α qh ≈ α Q (cid:113) α s N c π and ω q ≈ ω Q (cid:113) α s N c π which are independent of the initial conditions and are universal properties of the evolution. This ambiguity in theinitial phase unfortunately trickles down to the ambiguity in our prediction for the quark helicity inside a proton.Fixing the phase should be done by determining the initial conditions for the evolution. Resolving this issue wouldrequire a more detailed phenomenological work, which is beyond the scope of this paper.Note again that the oscillation frequency ω q vanishes in the N f = 0 limit, such that the oscillation is a property ofhaving quarks in the evolution. The N f = 0, N c → ∞ limit itself warrants a little further discussion. In [1–3, 40–42]the large- N c limit was understood as N f = 0 gluons-only evolution. In this regime the polarized Wilson lines aregiven only by the first ( ∼ F ) terms in Eqs. (11) and (14), such that Q ( x , z ) ≈ G ( x , z ) / x asymptotics of the quark and gluon amplitudes Q and G were, therefore, given by the same power of 1 /x . However,strictly-speaking one needs to clarify what is implied by the quark dipole amplitude at N f = 0. The quark dipoleamplitude obtained as a part of the gluon dipole amplitude in the large- N c limit, as employed in [1–3, 40–42] andgiven by Q ( x , z ) ≈ G ( x , z ) /
4, is not exactly the right object to determine the quark helicity PDF ∆Σ in Eq. (7).To properly define ∆Σ as the quark helicity distribution in the large- N c limit, one needs to have a non-zero N f . Onemay assume that the corresponding limit of small finite N f > N c → ∞ , with α s N c = const, can be imposedin equations (16) by dropping all the N f -terms, that is, formally by putting N f = 0 in the equations. We reporthere that the solution of the resulting equations performed as part of this work leads to the intercept α Q = 2 .
39. Ifwe compare this intercept to those listed in Table I we see that this result appears to support our earlier conclusionabout mild N f -dependence of the intercept.However, the large- N c finite- N f limit may need to be taken more carefully. Note that the small N f > N c → ∞ , α s N c = const limit is further complicated by the fact that in Q ( x , z ) the Born-level interaction with the quark5target, which should be used for calculating Q (0) ( x , z ), consists of a sum of two terms with different N c -scaling[40] (see Eq. (30) above): the t -channel quarks exchange is N c -enhanced compared to the t -channel gluon exchangecontribution. In general, the N c scaling of the terms in Q (0) and G (0) depends on whether the interaction happenswith a quark or with a gluon in the target. We see that the small N f > N c → ∞ , α s N c = const limit shouldbe taken by systematically performing the 1 /N c expansion in both the evolution equations (16) (or, more precisely,in the equations for Wilson lines from which Eqs. (16) were obtained in [1, 42]) and in the inhomogeneous terms Q (0) and G (0) . A careful implementation of this expansion shows that imposing the aforementioned limit by putting N f = 0 in equations (16) is correct only for a pure-glue shock wave, that is, for probing a polarized glueball stateinstead of the proton. For the scattering on the actual proton, or on a single polarized quark, the limit is more subtle,and requires a dedicated study. (It appears that, in this case, at leading- N c (and finite small N f ) one has to discardthe gluon dipole amplitudes G and Γ from Eqs. (16), which is likely to result in a solution for Q which is ratherslowly-growing with energy. At the first subleading- N c order, one would get a sum of the solution of the Eqs. (16)with N f = 0 but with N c -suppressed part of the inhomogeneous terms, and the solution of Eqs. (16) with leading- N c initial conditions and exactly one iteration of the N f -terms.) This regime has not been explored in this work due tothe higher phenomenological relevance of the large- N c & N f approximation considered here.Finally, let us comment on the potential phenomenological implications of our main qualitative result: at small- x the flavor-singlet quark helicity distribution ∆Σ( x, Q ) oscillates in ln(1 /x ). Similar oscillation has been found in thestrange quark helicity distribution ∆ s extracted from the experimental data by the PDF collaborations [4–7]. Thisoscillation is the driving force behind the sign change of ∆Σ DSSV ( x, Q ) in Fig. 5 above. If the ∆ s oscillation in x is confirmed by the future data extractions, it appears reasonable to ask a question whether it is related to ouroscillating result (36) for ∆Σ( x, Q ). While we do not separately consider individual quark flavors, the frequencyof ∆Σ oscillations increases with N f , and should be more pronounced if more flavors are included in the helicityevolution. This may be related to the oscillation observed in ∆ s , and not in ∆ u or ∆ d . On the other hand the periodof these oscillations in ln(1 /x ) is T = 2 πω q . (37)Using ω q = (0 . (cid:113) α s N c π from Eq. (33) with α s = 0 . N c = 3, we obtain T ( N f = 3) ≈
35, which is a verylarge number for rapidity or for ln(1 /x ). However, the first sign flip one encounters depends on the initial phase ofthe oscillation, which is hard to determine. If we start from the maximum of the cosine function in Eq. (36) at some x , then the first sign flip would happen at x/x ≈ e − T/ ≈ − , which is a more phenomenologically-reasonablenumber, but is still rather low to be relevant for the upcoming Electron-Ion Collider (EIC) experimental program[19, 23, 58]. On yet another hand, the period of oscillations found above may be significantly affected by the higher-order corrections in α s , even by the running of the coupling. In addition, we are encouraged by the similarity of theshapes between our curves and the DSSV14 line in Fig. 5, indicating that some of the physics behind the DSSV14line might be accurately described by our small- x evolution. We, therefore, leave the final verdict on the issue of thephenomenological relevance of the ∆Σ oscillations found in this work for the future investigations. VII. ACKNOWLEDGMENT
The authors would like to thank Mr. Daniel Adamiak for providing his code to help us construct Figs. 5 and6, which, in turn, was based on the code provided by Prof. Daniel Pitonyak, to whom we are also grateful. Wealso thank Prof. Pitonyak for a discussion of numerical simulations for helicity at small x . YK would like to thankMr. Mohammed Karaki for his work on this project in its very early stages.This material is based upon work supported by the U.S. Department of Energy, Office of Science, Office of NuclearPhysics under Award Number de-sc0004286. Appendix A: Analysis of the solution for the polarized dipole amplitudes
This Appendix describes how we fit our numerical solution for the polarized dipole amplitudes Q (0 , η ) and G (0 , η )using the ansatz (28) and extract the corresponding intercepts, frequencies, and phases.Given the numerical values of G (0 , η ) and Q (0 , η ) for a set of η values with discrete spacing, we utilize the followingmethod to deduce the intercepts ( α G and α Q ), frequencies ( ω G and ω Q ), and initial phases ( ϕ G and ϕ Q ) of theoscillations. Consider a numerical simulation for a function of the form f ( η ) = Ke αη cos ( ω η + ϕ ) (A1)6with some constants, α , ω , ϕ , and K . This is the asymptotic form assumed in Eqs. (28) for G (0 , η ) and Q (0 , η ). Thesecond derivative of the logarithm of | f ( η ) | is d dη ln | f ( η ) | = ddη [ α − ω tan ( ω η + ϕ )] = − ω cos ( ω η + ϕ ) . (A2)A local maximum of this second derivative occurs when cos ( ω η + ϕ ) = ±
1, and hence the frequency ω can be foundfrom the value of the numerically-obtained second derivative at the maximum,max (cid:20) d dη ln | f ( η ) | (cid:21) = − ω , (A3)where we also adopt a convention in which ω >
0. We use the largest- η maximum available in our simulation to extract ω using Eq. (A3). (Indeed the extracted value of ω can be cross-checked by comparing π/ω to the spacing between thepositions of the local maxima along the η -axis in the numerical solution.) The phase ϕ can then be determined fromthe second derivative maximum condition ω η ∗ + ϕ = π n , where η ∗ is the numerically-extracted position of the samelargest- η maximum of the second derivative (A2) and n is an integer. The value of n is adjusted so that ϕ ∈ ( − π, π ],where the choice between ϕ ∈ (0 , π ] and ϕ ∈ ( − π,
0] is done by making sure that the corresponding f ( η ∗ ) given byEq. (A1) is positive or negative (assuming that K >
0) in agreement with the numerical value of the function f ( η ) at η = η ∗ . (a) d dη ln | G (0 , η ) | (b) d dη ln | Q (0 , η ) | FIG. 7: Plots of d dη ln | G (0 , η ) | and d dη ln | Q (0 , η ) | for N f = 3 and N c = 3. Both graphs result from our numericalcomputation with step size ∆ η = 0 .
075 and η max = 30.Finally, a linear regression on ln (cid:12)(cid:12)(cid:12)(cid:12) f ( η )cos ( ω η + ϕ ) (cid:12)(cid:12)(cid:12)(cid:12) = α η + ln K (A4)allows us to extract α from the slope of this function. For the numerical values of G (0 , η ) and Q (0 , η ) found in therange η ∈ [0 , η max ], we only use η ∈ [0 . η max , η max ] to extract the intercepts α G and α Q using Eq. (A4), in additionavoiding the values of η close to the nearest cosine zero, η n = ω (cid:2) π − ϕ + π n (cid:3) , by at least 5% of the cosine’s period,(0 . T = π ω . This is done in order to obtain the intercept as close as possible to the asymptotic value, minimizingthe errors due to numerical artifacts and oscillation.Following the above prescription, in Fig. 7 we plot d dη ln | G (0 , η ) | and d dη ln | Q (0 , η ) | as functions of η for N f = 3.The corresponding plots for other values of N f also display the same qualitative behavior. For large η , the shapesof the graphs approach that of the function in Eq. (A2), displaying periodic local maxima below the η -axis. Thisprovides another justification for the proposed asymptotic form, (28), for G (0 , η ) and Q (0 , η ).The method outlined above is employed to determine the asymptotic forms for G (0 , η ) and Q (0 , η ) at N f = 2 , , N f , ∆ η and η max , we find the values of ω G , ω Q , ϕ G and ϕ Q from the largest maximum ( η ∗ )of the graphs in Fig. 7, which, in turn, correspond to the function in (A2), in order to get as close as possible to theasymptotic behavior at large η . The frequencies are found by using Eq. (A3), while the phases are extracted using ω η ∗ + ϕ = π n , with the integer n adjusted as described above. The parameters α G and α Q can be deduced from theslope of the function in Eq. (A4).7As expected for any numerical and asymptotic solution, the resulting intercepts, α G and α Q , frequencies, ω G and ω Q , and initial phases, ϕ G and ϕ Q , all vary slightly with the step size ∆ η and the maximum η max of the computationrange η ∈ [0 , η max ]. Since the exact continuum and asymptotic solution corresponds to the limit where ∆ η → η max → + ∞ , we perform the computation with various ∆ η ’s and η max ’s, obtaining these parameters for eachcomputation. Then, these values of the parameters are fitted with second-order polynomials in ∆ η and 1 /η max (cf.[2]). For each parameter, the value of its best-fit quadratic surface at ∆ η = 0 , η max = 0 is taken to be our numericalestimate for the parameter. This technique is employed to obtain the estimates of α G , α A , ω G , ω A , ϕ G and ϕ A for N f = 2 , , η → , η max →
0. Fig. 8 displays by the dots the values of these parameters at N f = 3 foreach pair of ∆ η and η max that we performed the computation for, together with the best-fit quadratic surfaces. Forcompleteness, let us list the equations describing the best-fit quadratic surfaces (for N f = 3): α G (∆ η, η max ) = 2 . − .
086 (∆ η ) − .
608 (1 /η max ) − .
472 (∆ η ) + 1 .
350 (∆ η/η max ) − .
108 (1 /η max ) , (A5a) ω G (∆ η, η max ) = 0 .
470 + 0 .
098 (∆ η ) + 0 .
018 (1 /η max ) − .
388 (∆ η ) + 0 .
018 (∆ η/η max ) − .
118 (1 /η max ) , (A5b) ϕ G (∆ η, η max ) = 0 .
327 + 0 .
516 (∆ η ) − .
403 (1 /η max ) − .
026 (∆ η ) − .
618 (∆ η/η max ) + 11 .
143 (1 /η max ) , (A5c) α Q (∆ η, η max ) = 2 . − .
014 (∆ η ) − .
008 (1 /η max ) − .
331 (∆ η ) − .
206 (∆ η/η max ) + 9 .
654 (1 /η max ) , (A5d) ω Q (∆ η, η max ) = 0 .
469 + 0 .
094 (∆ η ) + 0 .
012 (1 /η max ) − .
385 (∆ η ) + 0 .
077 (∆ η/η max ) − .
052 (1 /η max ) , (A5e) ϕ Q (∆ η, η max ) = − .
409 + 0 .
386 (∆ η ) − .
230 (1 /η max ) + 0 .
032 (∆ η ) − .
795 (∆ η/η max ) + 22 .
008 (1 /η max ) . (A5f)The qualitative features of these plots and quadratic fit functions are similar for N f = 2 ,
6, but we omit them forbrevity. The resulting values of the parameters extracted with the quadratic fit are given in Table I of the main text.For each of the parameters, the quadratic fit gives the value for the coefficient of determination, R , of at least 0 . η → η max → + ∞ extrapolations of the quadratic and linear fits in ∆ η and 1 /η max to the data points in Fig. 8. (a) α G for N f = 3 (b) ω G for N f = 3 (c) ϕ G for N f = 3(d) α Q for N f = 3 (e) ω Q for N f = 3 (f) ϕ Q for N f = 3 FIG. 8: Plots of parameters α G , α Q , ω G , ω Q , ϕ G and ϕ Q as functions of ∆ η and 1 /η max for N f = 3 and N c = 3. Thedots represent our numerical evaluation, while the solid surfaces depict the best fits (quadratic in ∆ η and 1 /η max )used for extrapolating to the continuum asymptotic values at ∆ η = 0 and η max = 0.One may wonder why α G and α Q in Fig. 8 approach ∆ η = 0 and η max = 0 from below, while in [2] the inter-cept approached the same limit from above. To better understand the difference, we re-ran the large- N c evolutionsimulations done in [2] with the initial conditions different from those used in [2]. We used G (0) = 1, while in [2]Born-level initial conditions (30) were employed. Using the trivial initial condition G (0) = 1 resulted in the intercept8approaching the ∆ η = 0 and η max = 0 limit from below for the large- N c evolution. We thus conclude that, while theasymptotic and continuum value of the intercept appears to be independent of the initial conditions, the approach tothis intercept in finite-step-size and finite- η -range numerical simulations appears to depend on the initial conditions.As a final cross check for the asymptotic form (28), we plot e − α G η G (0 , η ) and e − α Q η Q (0 , η ) versus η in Fig. 9. Wesee that the functions display clear sinusoidal pattern for η (cid:38)
10, demonstrating sinusoidal oscillation in the large- η asymptotics, as expected from the ansatz (28). (a) e − α G η G (0 , η ) at N f = 2 (b) e − α G η G (0 , η ) at N f = 3 (c) e − α G η G (0 , η ) at N f = 6(d) e − α Q η Q (0 , η ) at N f = 2 (e) e − α Q η Q (0 , η ) at N f = 3 (f) e − α Q η Q (0 , η ) at N f = 6 FIG. 9: Plots of e − α G η G (0 , η ) and e − α Q η Q (0 , η ) at N f = 2 , , N c = 3. All the graphs are numerically computedwith step size ∆ η = 0 .
075 and η max = 30. Appendix B: Analysis of the numerical results for ∆Σ(a) α q for N f = 3 (b) ω q for N f = 3 (c) ϕ q for N f = 3 FIG. 10: Plots of parameters α q , ω q and ϕ q as functions of ∆ η and 1 /η max for N f = 3 and N c = 3. The dotsrepresent our numerical evaluation, while the solid surfaces depict the best fits (quadratic in ∆ η and 1 /η max ) usedfor extrapolating to the continuum asymptotic values at ∆ η = 0 and η max = 0.Here we use the fitting method outlined in Appendix A to compute the parameters given in Eq. (33) describingthe small- x asymptotics of ∆Σ( x, Q ) in Eq. (32). The only difference is that the variable η in Eqs. (A1)-(A4) nowbecomes (cid:113) α s N c π ln x . At the end, the parameters α q , ω q , and ϕ q are extracted by using the following quadratic9best-fit surfaces for α s = 0 .
25, resulting in the values listed in Eq. (33): α q (∆ η, η max ) = 0 . − .
016 (∆ η ) − .
362 (1 /η max ) − .
192 (∆ η ) + 0 .
216 (∆ η/η max ) + 2 .
852 (1 /η max ) , (B1a) ω q (∆ η, η max ) = 0 .
162 + 0 .
032 (∆ η ) + 0 .
052 (1 /η max ) − .
132 (∆ η ) + 0 .
024 (∆ η/η max ) − .
444 (1 /η max ) , (B1b) ϕ q (∆ η, η max ) = − .
25 + 0 .
96 (∆ η ) − .
36 (1 /η max ) + 0 .
32 (∆ η ) − .
66 (∆ η/η max ) + 45 .
43 (1 /η max ) . (B1c)Similarly, for all the parameters, we obtain the coefficient of determination R of at least 0 . η and η max that we performed the computation for, togetherwith the best-fit quadratic surfaces. The error bars in Eq. (33) are also calculated by the difference between the linearand quadratic extrapolations to ∆ η = 0 and η max = 0. [1] Y. V. Kovchegov, D. Pitonyak and M. D. Sievert, Helicity Evolution at Small-x , JHEP (2016) 072, [ ].[2] Y. V. Kovchegov, D. Pitonyak and M. D. Sievert, Small- x asymptotics of the quark helicity distribution , Phys. Rev. Lett. (2017) 052001, [ ].[3] Y. V. Kovchegov, D. Pitonyak and M. D. Sievert,
Small- x Asymptotics of the Quark Helicity Distribution: AnalyticResults , Phys. Lett.
B772 (2017) 136–140, [ ].[4] D. de Florian, R. Sassot, M. Stratmann and W. Vogelsang,
Extraction of Spin-Dependent Parton Densities and TheirUncertainties , Phys. Rev.
D80 (2009) 034030, [ ].[5]
NNPDF collaboration, R. D. Ball, S. Forte, A. Guffanti, E. R. Nocera, G. Ridolfi and J. Rojo,
Polarized PartonDistributions at an Electron-Ion Collider , Phys. Lett. B (2014) 524–531, [ ].[6] D. de Florian, R. Sassot, M. Stratmann and W. Vogelsang,
Evidence for polarization of gluons in the proton , Phys. Rev.Lett. (2014) 012001, [ ].[7] J. Ethier, N. Sato and W. Melnitchouk,
First simultaneous extraction of spin-dependent parton distributions andfragmentation functions from a global QCD analysis , Phys. Rev. Lett. (2017) 132001, [ ].[8] R. L. Jaffe and A. Manohar,
The G(1) Problem: Fact and Fantasy on the Spin of the Proton , Nucl. Phys.
B337 (1990)509–546.[9] X.-D. Ji,
Gauge-Invariant Decomposition of Nucleon Spin , Phys. Rev. Lett. (1997) 610–613, [ hep-ph/9603249 ].[10] S. Bashinsky and R. L. Jaffe, Quark and gluon orbital angular momentum and spin in hard processes , Nucl. Phys.
B536 (1998) 303–317, [ hep-ph/9804397 ].[11] P. Hagler and A. Schafer,
Evolution equations for higher moments of angular momentum distributions , Phys. Lett.
B430 (1998) 179–185, [ hep-ph/9802362 ].[12] A. Harindranath and R. Kundu,
On Orbital angular momentum in deep inelastic scattering , Phys. Rev.
D59 (1999)116013, [ hep-ph/9802406 ].[13] Y. Hatta and S. Yoshida,
Twist analysis of the nucleon spin in QCD , JHEP (2012) 080, [ ].[14] X. Ji, X. Xiong and F. Yuan, Probing Parton Orbital Angular Momentum in Longitudinally Polarized Nucleon , Phys.Rev.
D88 (2013) 014041, [ ].[15]
European Muon collaboration, J. Ashman et al.,
A Measurement of the Spin Asymmetry and Determination of theStructure Function g(1) in Deep Inelastic Muon-Proton Scattering , Phys. Lett.
B206 (1988) 364.[16]
European Muon collaboration, J. Ashman et al.,
An Investigation of the Spin Structure of the Proton in Deep InelasticScattering of Polarized Muons on Polarized Protons , Nucl. Phys.
B328 (1989) 1.[17]
STAR collaboration, L. Adamczyk et al.,
Precision Measurement of the Longitudinal Double-spin Asymmetry forInclusive Jet Production in Polarized Proton Collisions at √ s = 200 GeV , Phys. Rev. Lett. (2015) 092002,[ ].[18]
PHENIX collaboration, A. Adare et al.,
Inclusive cross section and double-helicity asymmetry for π production atmidrapidity in p + p collisions at √ s = 510 GeV , Phys. Rev.
D93 (2016) 011501, [ ].[19] A. Accardi et al.,
Electron Ion Collider: The Next QCD Frontier , Eur. Phys. J.
A52 (2016) 268, [ ].[20] E. Leader and C. Lorc,
The angular momentum controversy: What’s it all about and does it matter? , Phys. Rept. (2014) 163–248, [ ].[21] E. C. Aschenauer et al.,
The RHIC Spin Program: Achievements and Future Opportunities , .[22] E.-C. Aschenauer et al., The RHIC SPIN Program: Achievements and Future Opportunities , .[23] C. A. Aidala et al., Probing Nucleons and Nuclei in High Energy Collisions , .[24] E. Iancu and R. Venugopalan, The color glass condensate and high energy scattering in QCD , hep-ph/0303204 .[25] H. Weigert, Evolution at small x bj : The Color Glass Condensate , Prog. Part. Nucl. Phys. (2005) 461–565,[ hep-ph/0501087 ].[26] J. Jalilian-Marian and Y. V. Kovchegov, Saturation physics and deuteron-Gold collisions at RHIC , Prog. Part. Nucl.Phys. (2006) 104–231, [ hep-ph/0505052 ].[27] F. Gelis, E. Iancu, J. Jalilian-Marian and R. Venugopalan, The Color Glass Condensate , Ann.Rev.Nucl.Part.Sci. (2010) 463–489, [ ]. [28] J. L. Albacete and C. Marquet, Gluon saturation and initial conditions for relativistic heavy ion collisions , Prog.Part.Nucl.Phys. (2014) 1–42, [ ].[29] Y. V. Kovchegov and E. Levin, Quantum chromodynamics at high energy , vol. 33. Cambridge University Press, 2012.[30] J. Bartels, B. Ermolaev and M. Ryskin,
Nonsinglet contributions to the structure function g1 at small x , Z.Phys.
C70 (1996) 273–280, [ hep-ph/9507271 ].[31] J. Bartels, B. Ermolaev and M. Ryskin,
Flavor singlet contribution to the structure function G(1) at small x , Z.Phys.
C72 (1996) 627–635, [ hep-ph/9603204 ].[32] R. Kirschner and L. Lipatov,
Double Logarithmic Asymptotics and Regge Singularities of Quark Amplitudes with FlavorExchange , Nucl.Phys.
B213 (1983) 122–148.[33] R. Kirschner,
Reggeon interactions in perturbative QCD , Z.Phys.
C65 (1995) 505–510, [ hep-th/9407085 ].[34] R. Kirschner,
Regge asymptotics of scattering with flavor exchange in QCD , Z.Phys.
C67 (1995) 459–466,[ hep-th/9404158 ].[35] S. Griffiths and D. Ross,
Studying the perturbative Reggeon , Eur.Phys.J.
C12 (2000) 277–286, [ hep-ph/9906550 ].[36] Y. Hatta, Y. Nakagawa, F. Yuan, Y. Zhao and B. Xiao,
Gluon orbital angular momentum at small- x , Phys. Rev.
D95 (2017) 114032, [ ].[37] G. A. Chirilli,
Sub-eikonal corrections to scattering amplitudes at high energy , JHEP (2019) 118, [ ].[38] Y. V. Kovchegov, Orbital Angular Momentum at Small x , JHEP (2019) 174, [ ].[39] R. Boussarie, Y. Hatta and F. Yuan, Proton Spin Structure at Small- x , Phys. Lett.
B797 (2019) 134817, [ ].[40] Y. V. Kovchegov, D. Pitonyak and M. D. Sievert,
Helicity Evolution at Small x : Flavor Singlet and Non-SingletObservables , Phys. Rev.
D95 (2017) 014033, [ ].[41] Y. V. Kovchegov, D. Pitonyak and M. D. Sievert,
Small- x Asymptotics of the Gluon Helicity Distribution , JHEP (2017) 198, [ ].[42] Y. V. Kovchegov and M. D. Sievert, Small- x Helicity Evolution: an Operator Treatment , Phys. Rev.
D99 (2019) 054032,[ ].[43] F. Cougoulic and Y. V. Kovchegov,
Helicity-dependent generalization of the JIMWLK evolution , Phys. Rev.
D100 (2019)114020, [ ].[44] I. Balitsky,
Operator expansion for high-energy scattering , Nucl. Phys.
B463 (1996) 99–160, [ hep-ph/9509348 ].[45] I. Balitsky,
Factorization and high-energy effective action , Phys. Rev.
D60 (1999) 014020, [ hep-ph/9812311 ].[46] Y. V. Kovchegov,
Small-x F structure function of a nucleus including multiple pomeron exchanges , Phys. Rev.
D60 (1999) 034008, [ hep-ph/9901281 ].[47] Y. V. Kovchegov,
Unitarization of the BFKL pomeron on a nucleus , Phys. Rev.
D61 (2000) 074018, [ hep-ph/9905214 ].[48] H. Weigert,
Unitarity at small Bjorken x , Nucl. Phys.
A703 (2002) 823–860, [ hep-ph/0004044 ].[49] J. Jalilian-Marian, A. Kovner and H. Weigert,
The Wilson renormalization group for low x physics: Gluon evolution atfinite parton density , Phys. Rev.
D59 (1998) 014015, [ hep-ph/9709432 ].[50] J. Jalilian-Marian, A. Kovner, A. Leonidov and H. Weigert,
The Wilson renormalization group for low x physics:Towards the high density regime , Phys. Rev.
D59 (1998) 014014, [ hep-ph/9706377 ].[51] E. Iancu, A. Leonidov and L. D. McLerran,
The renormalization group equation for the color glass condensate , Phys.Lett.
B510 (2001) 133–144.[52] E. Iancu, A. Leonidov and L. D. McLerran,
Nonlinear gluon evolution in the color glass condensate. I , Nucl. Phys.
A692 (2001) 583–645, [ hep-ph/0011241 ].[53] E. Ferreiro, E. Iancu, A. Leonidov and L. McLerran,
Nonlinear gluon evolution in the color glass condensate. II , Nucl.Phys.
A703 (2002) 489–538, [ hep-ph/0109115 ].[54] E. Iancu, J. D. Madrigal, A. H. Mueller, G. Soyez and D. N. Triantafyllopoulos,
Resumming double logarithms in theQCD evolution of color dipoles , Phys. Lett.
B744 (2015) 293–302, [ ].[55] B. Duclou´e, E. Iancu, A. H. Mueller, G. Soyez and D. N. Triantafyllopoulos,
Non-linear evolution in QCD at high-energybeyond leading order , JHEP (2019) 081, [ ].[56] E. A. Kuraev, L. N. Lipatov and V. S. Fadin, The Pomeranchuk singlularity in non-Abelian gauge theories , Sov. Phys.JETP (1977) 199–204.[57] I. Balitsky and L. Lipatov, The Pomeranchuk Singularity in Quantum Chromodynamics , Sov.J.Nucl.Phys. (1978)822–829.[58] D. Boer et al., Gluons and the quark sea at high energies: Distributions, polarization, tomography ,1108.1713