Exploring origin of small x saturation in collinear approach

A modification of the collinear evolution equations as an appropriate approach to improve the behavior of parton distribution functions in the region of small longitudinal momentum fractions, and to find more theoretical arguments to clarify the possible appearence of saturation regime is suggested. It is argued that parton diffusion in the rapidity space at large parton densities along the space-time evolution could result in the emergence of a natural saturation scale on which freezing actually occurs.


I. INTRODUCTION
Nowadays it is widely recognized the hadron interactions at very high energies are driven by the states with very high densities of partons (quarks and gluons), in particular, with small longitudinal momentum fractions x. The routine theoretical framework for analyzing such systems is essentially grounded on the QCD collinear factorization where the calculated cross sections are decomposed in the perturbative coefficient functions and nonperturbative parton densities of which evolvement is treated according the DGLAP equations [1][2][3][4]. Already these linear equations capture qualitatively the traits associated with an increase in the gluon densities at small x with extremely large Q 2 values. The latter turn out also quite instrumental, for example, to justify by neglecting any type of higher twist corrections and some perturbative resummation contributions. An idea to follow the evolution within the perturbative paradigm and to evaluate the leading contributions at small x for not very large Q 2 led to the development of the BFKL approach associated to so-called high energy factorization. However, resolving the corresponding BFKL equations [5][6][7] exhibits a very strong raise (power-like) of the gluon density at small x that is stronger than the experimental data analysis demonstrates and leads to apparent violation of unitarity at very small x. It signals some theoretical problems generated by appearance of an infrared instability related to a diffusion with the rapidity evolution and the consistent description of QCD coupling constant α s behavior that should reflect a very sophisticated interplay of perturbative and nonperturbative QCD physics. Apparently, both look like an ensuing result of taking into account the linear evolution only with resummation in these approaches. From phenomenological point of view an observation of a scaling law in wide range of small x and Q 2 was done [8] thereby demonstrating an onset of saturation scale. This fact is quite interesting because it may provide a perturbative scale in high density region of small x where linear evolution approximtion works and provides, in a sense, a boundary condition to the linear evolution equations. In fairness, remember it has been long time ago [9,10] argued that eventually the system under consideration should enter a new regime, where the rate of growing gluon density slows down and saturates possibly curing, thus, a potential conflict with unitarity of the underlying scattering. Actually, the restoration of the unitarity in high energy limit of QCD remains a challenging problem, although several approaches, drawing a scenario with nonlinear behavior, are being explored in past years (see, for example, [11][12][13][14][15][16][17][18] and references therein) but those allow us to conclude only that we are still no essentially closer to knowing where the problem solution lies besides of very general claim about the nonperturbative finite density effects which are left out entirely from the BFKL evolution. The interest in physics of high density regime of small x QCD is greatly increasinng and dictated by an avalanche of experimental data on collisions of relativistic heavy ions overwhelming this area of research in the last decades.
Meanwhile, there is another opportunity to address the problem in the framework of well-known DGLAP approach that we would like to draw attention to in this letter. It concerns one possible modificaion of collinear time-like evolution equations that has been also discussed long time ago [19,20] as well in the context of increasing parton multiplicity in electron-positron annihilation into hadrons. We adapt this modification for the space-like evolution of parton distribution functions and demonstrate it develops a saturated regime of color glass condensate [21]. This is a regime of strong color fields in which nonlinear dynamics come to the perceptible play and signals thereby an appearance of the natural saturation scale on which the evolution is, in fact, frozen, thus indicating also the universality of both phenomena. Actually, such an approach is treated as an effective theory of high energy scattering successefuly describing the data measured in experiments.
The paper is organized as follows. In Sec. II we review briefly the principal features of the DGLAP evolution necessary in further. The particular modified QCD evolution is discussed in Sec. III. In Sec. IV the extension to the double parton distribution functions is considered. The possible phenomenological issues are discussed in Sec. V, together with some conclusions.

II. COLLINEAR EVOLUTION
One may take the value of hard scale as the evolution variable in the DGLAP approach. The most popular choice is the transfer momentum squared Q 2 , or its logarithm ξ = ln(Q 2 /Q 2 0 ). The double logarithm that takes into account explicitly the behavior of the effective coupling constant in the leading logarithm approximation proves very instrumental as well where β = (11N c − 2n f )/3 in QCD, Q 0 is the some characteristic scale above which the perturbative theory is applicable, n f is the number of active flavors, Λ is the QCD dimensional parameter and N c = 3 is the color number. In Eq. (1) the one loop running QCD coupling was used.
The DGLAP evolution equations [1][2][3][4] assume the simplest form if we use the natural dimensionless evolution variable t; that is, These equations describe the evolution of single distributions D j h (x, t) of bare quarks, antiquarks and gluons (j = q,q, g) within a hadron h in response to the change of evolution variable t. The kernels, P , of these equations include a regularization at x → x ′ and are known in their appropriate forms.
Equations (3) are explicitly solved by introducing the Mellin transforms which reduce those to a system of ordinary lineardifferential equations at the first order: where In order to obtain the distributions in x representation the inverse Mellin transformation should be performed where the integration runs along the imaginary axis to the right from all n singularities. It can be done in a general form numerically only. However, the asymptotic behavior can be estimated in some interesting and simple enough limits with the technique under consideration. The solutions of the DGLAP equations with the given initial conditions D j h (x, 0) at the reference scale Q 0 (t = 0) can be expressed by the Green functions D j i (z, t) in the following way: These Green functions (gluon distributions at the parton level) D j i (z, t) are the solutions of Eqs.
(3) at the parton level with the singular initial conditions D j i (z, t = 0) = δ(x − 1)δ ij and in the double logarithm approximation (see, for instance, [3,9]) look like where v = 8N c t ln (1/x), a = 11 and I 1 is the standard modified Bessel function. This result just illustrates the unitarity violation at very small x.
In addition, one should also note that the mean number of partons of type j in a parton of type i can not be correctly determined in the collinear approach because the kernels P g→g (0) and P q→g (0) are divergent and some improvements are necessary to be done at very small x.

III. COLLINEAR EVOLUTION WITH DISSIPATION
The modification of collinear time-like evolution equations was discussed in Refs. [19,20] to take into account the formation (so-called pionization) of soft quarkantiquark pairs at a hard quark (gluon) propagating. In analogy with the electron-photon showers the energy outflow was phenomenologically simulated by the dissipative terms in the evolution equations with rather interesting income. Such a modification for the space-like evolution has, of course, another physical motivation in our case due to the parton diffusion in the rapidity space at large parton densities, and the following evolution equations are suggested: with γ j as some parameters characterizing the process of the energy outflow.
In the situation of small dissipation, γ j ≪ 1, the mean number of partons can be calculated [20] by using the Mellin technique. For credibility we bring here the result for gluon multiplicity at the early evolution stage (t ≪ 1) only referring to the transparent, but laborious, calculations performed in [20]: where V = 8N c t ln (1/γ g ), and I 0 is another modified Bessel function. This result (13) reproduces exactly the mean number of gluons with the longitudinal momentum fractions larger than x 0 = γ g as calculated in the DGLAP unmodified approach. The exercise above makes transparent the physical meaning of the dissipative term. It establishes the scale of energy drift because gluons (partons) with the longitudinal momentum fractions less than γ g are simply withdrawn from consideration. Moreover, the evolution is, in fact, frozen at the scale [19] The origin of this freezing scale is similar to the saturation scale in the color glass condensate (CGC) approach [13,[15][16][17][18]21].

IV. GENERALIZING TO DOUBLE PARTON DISTRIBUTIONS
The extension of basic equations to double parton distribution functions is straightforward: Here, the splitting kernels, which appear in the nonhomogeneous part of the equations, are the nonregularized one-loop well-known DGLAP kernels without the "+" prescription. The unmodified equations were derived first in Refs. [22,23] in framework of the DGLAP approach. The functions D j1j2 h (x 1 , x 2 , t) in question have a specific interpretation in the leading logarithm approximation of perturbative QCD. They are the inclusive probabilities which allow one to find two bare partons of types j 1 and j 2 with the given longitudinal momentum fractions x 1 and x 2 in a hadron h.
The dissipative terms provide the energy outflow and establish the scale of energy drift as well. Gluons (partons) with the longitudinal momentum fractions less than γ j are simply removed again from consideration for each of two parton cascade branches practically independently. In the small x region we can restrict ourselves to homogeneous evolution equations because the solutions of nonhomogeneous unmodified equation are substantial at not parametrically small longitudinal momentum fractions only [24]. Moreover, the homogeneous evolution equations (independent evolution of two branches) admit the factorization of double parton distribution functions: as a good approximate solution, if such a factorization was assumed at the reference scale Q 0 (t = 0). Further we hold the leading exponential terms only if those have the same structure [25] both at the parton level and the hadron level under smooth enough initial conditions at the reference scale. Indeed, Eq. (8) in the double logarithm approximation reads with Y = ln(1/x). The y integration is not as a saddlepoint type, and, therefore, one of the edges, just y → 0 (z → 1), dominates, provided that the initial gluon distribution does not increase too much with z decreasing. Actually, one needs zD g h (z, Let's notice that the parametrization of the initial gluon distributions, usually used, satisfies this condition (e.g., the CTEQ parametrization of Ref. [26]). Thus, as a result we have for the double gluon distributions [25] in this appproximation: (20) with the infinite mean number of such gluons. If the two branches evolves independently then by introducing the dissipative terms slow down the rate of gluon density increase and one gets the finite mean gluon numbers as < n gg > h ∼ exp [ 8N c ( t ln (1/γ g ) + t ln (1/γ g ))] (21) since the gluons with the longitudinal momentum fractions less than γ g are simply excluded.

V. DISCUSSION AND SUMMARY
Clearly, the dissipative parameters above can not be determined within the DGLAP approach. They are treated as the phenomenological parameters in numerical simulations and shoud be estimated in the other models for further applications. The phenomena of saturation and slowing down an increase of gluon density take place also in the CGC scenario [13,[15][16][17][18]. However, the saturation scale is energy dependent in that approach and, nevertheless, it comes about quite predictive. For example, in the Golec-Biernat-Wusthoff (GBW) model [27,28] it is parametrized by three parameters : with Q 0 = 1 GeV, x 0 ≃ 0.0001, and λ ≃ 0.3 which have been used to describe accurately the HERA data [29]. The value of characteristic energy (longitudinal momentum fraction) x 0 in Eq. (22) allows us to estimate the dissipative parameter γ g that has a physical meaning similar to x 0 . In fact, it justifies the assumption of small dissipation used in the previous Sections to obtain the crucial estimates (13) and (21) which are pretty encouraging to investigate the properties of modified collinear equations further as a new alternative insight into the saturation physics extending the initial limits of linear approach.
In summary, the modified collinear evolution equations are suggested to extract information on the properties of hot and dense QCD medium created in the experiments on heavy ion collisions searching the quark-gluon plasma, a thermolized phase, that may exist in very spesific regimes for very short periods of time. Comprehensive phenomenological analysis of proton-proton collisions based on the QCD factorization, as a key instrument, made it possible to extract the universal distribution functions validating such an approximation and open up (quite often) transparent ways for introducing the efficient corrections. Truly, these corrections at leading power of the large momentum transfer are fairly general and easily traceable but the corrections within the factorized forms turn out very complicated and too much sensitive to the process details, as it was shown again many years ago, because of the QCD multiple scaterrings which differ [30] hadronic and heavy ion collisions significantly. The model presented in this letter shows a possibility of forming perturbatively a dynamical regime in particular kinematical configuration which could not be foreseen according to the theoretical dogmas. As argued, it concerns a regime of high parton densities and dynamical interactions described definitely by nonlinear equations. The evolution of hadron scattering amplitudes, at least, in the framework of color dipole picture in such a regime is quite similar [31] to the time evolution of a classical particles system undergoing reaction-diffusion processes. Amazingly, by introducing the dissipative terms results in an origin of natural saturation scale on which the evolution is frozen and the gluons with longitudinal momentum less than γ g are simply excluded. In the phenomenological applications the direct numerical solutions of suggested modified equations may occur simpler than the BFKL treatment of very small x region.