Large impact parameter behaviour in the CGC/saturation approach: a new non-linear equation

In this paper we propose a solution to the long standing problem in the CGC/saturation approach: the power-like fall off of the scattering amplitudes at large $b$. We propose a new non-linear equation, which takes into account random walks both in transverse momenta of the produced gluons and in their impact parameters. We demonstrate, that this equation is in accord with previous attempts to include the diffusion in impact parameters in the BFKL evolution equation. We show in the paper, that the solution to a new equation results in the exponential decrease of the scattering amplitude at large impact parameter, and in the restoration of the Froissart theorem.


INTRODUCTION
It is well known that perturbative QCD has a fundamental problem: the scattering amplitude decreases at large impact parameters (b) as a power of b. In particular, the CGC/saturation approach [1], which is based on perturbative QCD, is confronted by this problem. At large b the scattering amplitude is small and, therefore, only the linear BFKL equation [2] is determined by the scattering amplitude in perturbative QCD. It is known that the eigenfunction of this equation (the scattering amplitude of two dipoles with sizes r and R) has the following form [3] φ γ (r, R, b) = One can see that at large impact parameter b the amplitude has a power-like decrease, which leads to the violation of the Froissart theorem [4]. The violation of the Froissart theorem stems from the growth of the radius of interaction as a power of the energy. Since in Ref. [3] it was proven that the eigenfunction of any kernel with conformal symmetry has the form of Eq. (1), we can only change the large b behaviour by introducing a new dimensional scale in the kernel of the equation. The problem has been known from the beginning of the saturation approach [5,6] and several ideas have been proposed, of how to introduce a new dimensional scale in the kernel of the BFKL equation (See Refs. [5][6][7][8][9][10]). However, for the high energy community at large, the problem was appreciated only after the papers of Refs. [11,12] were published, where it was demonstrated, that the violation of the Froissart theorem cannot be avoided in the framework of the CGC approach. First, we wish to illustrate why the Froissart theorem is violated for the BFKL equation. The general solution to the BFKL for the dipole scattering amplitude equation, has the following form: where 1 and the amplitude is equal to Using Eq. (4) we can attempt to determine the upper bound from the unitarity constraints where G in describes the contribution of all inelastic processes. Recalling that N is the imaginary part of the scattering amplitude and assuming that the real part of the amplitude is small, as is the case for the BFKL equation, one can see that the unitarity constraint has the solution: where Ω > 0 denotes an arbitrary function. Now we can find the bound for the total cross section following Ref. [4].
We need to solve the following equation to find the value of b 0 , N (r, Y ; b 0 ) = φ 1 2 At large Y this solution gives b 2 0 ∝ e D+ √ D(D+4ω0) Y and therefore, Note, that if for the soft amplitude of the typical form N soft ∝ e ω0Y −µb and b 0 ∝ Y . One can see, that the first integral in Eq. (7) leads to σ ≤ Const y 2 , which is the Froissart theorem. This amplitude violates the Froissart theorem and can be considered only at large values of b where N soft < 1. However, the eikonal solution of Eq. (6) with Ω = N sof t , satisfies the unitarity constraints and leads to the amplitude, that describes the Froissart disc: the amplitude is equal to 1 for b < (1/µ)Y and it has an edge which behaves as N soft .
We hope that this simple estimate shows, that the power-like decrease of the scattering amplitude is the source of the problem. To solve the problem in the framework of the CGC/ saturation approach we need to introduce a new dimensional scale µ, in addition to the saturation momentum. In Ref. [4] it is shown that this new scale is related to the mass of the lightest hadrons. Perturbative QCD cannot reproduce the observed spectrum of hadrons, and attempting to solve this problem, we are doomed to introduce something from non-perturbative QCD estimates. Since the non-perturbative approach is still in an embryonic stage, one can only guess how to introduce this scale, which depends crucially on non-perturbative estimates in lattice QCD, on phenomenology of high energy interactions and on intuition, which comes from considering different theoretical models. We hope that more or less full list of the attempts to solve this problem, and to find and introduce a new dimensional scale appear in Refs. [5][6][7][8][9][10][12][13][14][15][16][17][18][19][20][21][22][23][24].
At first sight the recent papers [25,26] have questioned the need of a new dimensional scale, since they demonstrate that the next-to-leading order BFKL equation generates the exponential type of the impact parameter behavior without the need of a new dimensional scale. However, it turns out [27] that the NLO corrections do not change the power-like decrease of the scattering amplitude at large impact parameter, generating the exponential-type decreases in the large but limited region of the values of the impact parameter. * In this paper we re-visit one of the possible ways of introducing a new dimensional scale: to incorporate in the BFKL equation the diffusion in impact parameters(b). The first such attempt was undertaken in the distant 1990's [5,6] and during three decades we have worked in this area [7-10, 14, 23, 24]. This is the reason that in the first part of the paper we give a brief review of these efforts, while in the second part we propose our generalization of the BFKL equation, which takes into account the diffusion on b in accord with QCD estimates. We will discuss the structure of the scattering amplitude at high energies, which at the moment appears to be a black disc with the radius increasing as a power of energy. This paper is partly motivated by Ref. [24], in which many questions concerning QCD has been formulated on the black disc behaviour at high energy † from the point of boost invariance and the parton model. We hope that we answer some of these questions in this paper. In particular, we will demonstrate that QCD leads to the Froissart disc at high energies, with specific behaviour of the amplitude at the edge of this disc.

A. Regge approach and Gribov's diffusion in impact parameters
In the framework of the Regge approach the high energy amplitude is given by the exchange of the Pomeron, and has the following form [28][29][30][31]: where g 1 , g 2 and trajectory α I P Q 2 T ≡ 1 + ∆ I P (Q T ) = 1 + ∆ I P − α I P Q 2 T + O Q 4 T are the functions that have to be taken from the phenomenology. Y = ln (s) and Q T is the momentum transferred by the Pomeron ‡ . Eq. (10) can be viewed as the solution to the following equation: with the initial condition In the impact parameter representation solution of Eq. (10) takes the form: In Eq. (13) we have neglected the Q T dependence of g 1 and g 2 which does not contribute at high energies. In Ref. [30] the simple fact is noted i.e. n (Y, b) is the solution of the diffusion equation: Eq. (14) together with Eq. (10) for the total cross section: have very simple interpretations in the parton model. In the parton model [30][31][32] it is assumed, that we can describe the interaction by a field theory, in which all integrals over transverse momenta are convergent, and they lead to the mean transverse momentum, which does not depend on energy. In such a theory, the contribution to the total cross section of the scattering amplitude for production of n partons in each order of perturbation approach, can be viewed as In Eq. (16) we assume that in the proposed theory the amplitude is not equal to zero, when rapidities of emitted partons are equal to zero, and choose the largest contribution which comes from the ordering 0 < y 1 < y 2 < · · · < y i < · · · < y n < Y One can see that in Eq. (15) M 2→2+n ({p i,T }) n i=0 d 2 p i,T = ∆ n I P § and from this equation we can conclude that the number of emitted partons n = ∆ I P Y . The Gribov idea that the emission of partons has no other correlations except the fixed transverse momentum, can be viewed as a random walk in two dimensional space. For each emission due to uncertainty principle Therefore, after each emission the position of the parton will be shifted by an amount ∆b 2 from Eq. (18), which is on average the same. After n emissions, we have the picture shown in Fig. 1, and the total shift in b is equal to In the case of the deep inelastic processes Y = ln (1/x), where x is the Bjorken variable. § In this estimate we assume that σ 0 in Eq. (13) does not depend on energy. In the field theories which can be a realization of the parton model, usually σ 0 ∝ 1/s 2 but the first result from QCD was the understanding that in this approach σ 0 does not depend on energy [33].
Therefore, this diffusion reproduces the shrinkage of the diffraction peak. Indeed, Comparing Eq. (19) and Eq. (20) one can see that Eq. (11) can be re-written in the impact parameter representation for Let us write the equation for the radius of interaction (see Eq. (20)) . First, we see that for we have the following equation: The second term vanishes due to integration over the angle and finally we have the following equation with the solution We can obtain Eq. (24) using the Mueller diagrams [34] of Fig. 2. Indeed, one can see that Therefore, for R 2 we obtain R 2 = 4α I P Y e ∆ I P Y /e ∆ I P Y = 4α I P Y in accord with Eq. (24).
The amplitude of Eq. (13) increases with energy and violates the unitary constraints (see Eq. (5)). The eikonal unitarization of Eq. (6) leads to the following amplitude One can see that this amplitude tends to unity at b > 2 ∆ I P α I P Y leading to the total cross section σ tot ∝ Y 2 in accord with the Froissart theorem [4].
In spite of the primitive level of calculations, especially if you compare them with typical QCD calculations in DIS, the parton model was a good guide for the Pomeron structure for years and, it is still the model where we can see all typical features of the soft Pomeron. It turns out that the time structure of the parton cascade is preserved for QCD and simple parton estimates can help develop our intuition regarding the solution of the QCD problems.
B. BFKL approach and diffusion in ln (pT ).

The BFKL equation.
The BFKL equation was derived in momentum representation [2] and has the following form: whereᾱ S = (N c /π) α S . Kernel K em describes the emission of a gluon, while kernel K reg is responsible for the reggeization of gluons in t-channel. They have the forms: This equation is rewritten in the coordinate representation for the scattering amplitude of a dipole with the size r at impact parameter b [3,38]: with K (r , r − r ; r) = r 2 In Eq. (29) 2. Solutions and random walk in ln (pT ).
We have discussed solutions in the coordinate representation (see Eq. (1), Eq. (2) and Eq. (4)). These solutions are very useful for discussions of the non-linear corrections, since the unitarity constraints are diagonal for the dipole scattering amplitude (see Eq. (5)). However, discussing the random walk in ln (p T ) we need a solution in the momentum representation in which p i,T are the momenta of produced gluons. Note, that in Eq. (28) q − q = p T . For Q T the solution can be easily obtained, since the eigenfunction have the following form [2] Comparing Eq. (32) with Eq. (1) one can see that q 0 ∼ 1/R, where R denotes the size of the target. Repeating all the steps, that are given by Eq.
withξ = ln q 2 /q 2 0 . One can recognize that for the function n we have the diffusion equation in the form: Therefore, the BFKL equation describes that at each emission, ln q 2 /q 2 changes its value by a constant (ln q 2 /q 2 ) 2 = δ, and as a result after n emissions we obtain ln 2 p 2 T /q 2 0 = δ n. Since n = ω 0 Y (see the previous section) we obtain (ln q 2 /q 2 0 ) 2 δω 0 Y . This estimate shows that 4D = δω 0 . From this estimate we see that after n emissions the typical transverse momenta increase as < |p 2 . Therefore, only a small number of steps at the beginning could participate in the increase of b (see red lines in Fig. 1-b ).

The Green function of the BFKL Pomeron.
The solution of Eq. (33) can be re-written in the following form: We introduce the Green function of the BFKL Pomeron as follows: The Green function in Y representation can be calculated as Integrating over κ 0 we obtain the solution of Eq. (33).
One can see that G BFKL Y = 0,ξ; Q T = 0 = δ ξ and therefore, the scattering amplitude can be found where Φ is the initial condition for the scattering amplitude. It should also be mentioned that factors 1/(q q 0 ) are absorbed in integration ofξ in the diagrams for Pomeron interactions. 4. The BFKL approach: random walk in b.
As one can see from Eq. (33) we have introduced a new dimensional scale: q 0 . It was introduced as the nonperturbative size of the target q 0 ∼ 1/R, but its actual meaning is the separation scale: in perturbative QCD we can calculate only for p T > q 0 , while smaller transverse momenta have to be treated in non-perturbative approaches. From the solution of Eq. (33) one can see that the probabilty to have p T = q 0 is small ∝ The gluons with transverse momenta of the order of q 0 could participate in the random walk in b, leading to [5,6] < |b 2 | > n = 1 where P n q0 denotes the probability to have the minimum momentum after n emissions. In Refs. [5,6] the numerical coefficient in Eq. (38) is evaluated. Since the average < |b 2 | > n diverges at small q, in some sense the value of the coefficient is not important. However, we feel it is instructive to understand two qualitative features: the infrared divergency and the energy dependence.
First, we calculate the main ingredients of the calculation that we have discussed in section II-A : As we have discussed in section II-A (see Eq. (23)), we can write the equation In doing so, we obtain: To obtain the complete system of equations we need to add the equation for ∇ Q T N (Y ; q , Q T ) | Q T =0 . It takes the form: Eq. (41) can be re-written as In the double Mellin transform of Eq. (36) the solution to Eq. (42) has the form: Plugging the solution of Eq. (43) into Eq. (40), we reduce this equation to the following one: The main contributions to the integrals over q stems from the region q → 0, Taking this into account the solution in the ω-representation takes the form: Therefore, we see from Eq.
As was shown in Refs [5,6], the second term does not have the suppression of the order ofᾱ S , as we can see from Eq. (45). Such an enhancement comes from the replacement in Eq. (40) stands for all needed integrations. The eigenfunction φ 1 [3] depends on the angles between q" and q and has an intercept which is negative and ∝ᾱ S . Integration over Y leads to the factor 1/ᾱ S which compensates the small value ofᾱ S . This contribution actually leads to the change of the numerical coefficient, which does not have much meaning, since b 2 is infrared unstable and q → q 0 gives the main contribution. The value of the scale q 0 is not determined in our approach. Hence, we introduce a new dimensional scale α eff which, we hope, will heal the large b behaviour of the scattering amplitude.

III. SEARCHING FOR A NEW DIMENSIONAL SCALE IN THE CGC/SATURATION APPROACH
In this section we would like to give a brief review of the attempts to introduce both the random walk in ln (p T ), and in b, in the effective theory describing QCD at high energies: the CGC/saturation approach. These efforts cover a long span in time, from the first GLR equation [35](see also Refs. [36][37][38]) till recent papers. Our goal is to introduce the random walk in b in the framework of the non-linear evolution equations, which guarantee the unitarity of the scattering amplitude.

A. Non-linear evolution equations in QCD
We start with the current form of the non-linear equation [39], which has the following form for the scattering amplitude of the dipole with size r and rapidity Y at impact parameter b: The kernel K (r , r − r ; r) is given by Eq. (30). For N (Y ; q, Q T ) and N (q, b, Y ) (see Eq. (31) ) can be re-written in the forms: with the kernels from Eq. (28).
The GLR equation The first non-linear equation [35] addressed the problem of large b behaviour of the scattering amplitude. The non-linear evolution equations of the previous subsections can be viewed as the sum of the 'fan' diagrams (see Fig. 3 for the example of such diagrams). It was shown in Ref. [35] that integrations over all transfer momenta Q T carried by the BFKL Pomerons in such diagrams, (see Fig. 3 for example) turns out to be of the order of Q" T ∼ 1/R h , where R h is the size of the hadron target. Therefore, if the dipole sizes are much smaller than R h , as it is for deep inelastic processes, we can replace N (Y ; q, Q T ) in Eq. (49) is a constant which depends on the unknown non-perturbative structure of the hadrons. S h is the area of the hadron which is populated by gluons. The GLR equation is written for the infrared safe observable: the gluon structure function; and it was suggested that it replace the DGLAP equation in the region of small x( high energies). Considering this equation we do not have a problem of the large b behaviour of the amplitude. The intuition behind this equation is very transparent: in deep inelastic processes the number of gluons increases with the growth of energy, while the size of their distribution remains the same until the hadron becomes black. At ultra high energies the radius of the distribution could increase and we cannot trust the equation. From Eq. (51) we can see that the large b problem is intimately related to the black disc behaviour of the scattering amplitude, and we need to develop a different equation to describe the behaviour of the edge of the disc.

C. The first equation incorporating both diffusions
In Ref. [5] the first equation was proposed that incorporates both diffusions in ln (p T ) and in b. This equation has the form: At large b the shadowing corrections (the non-linear term) in the equation makes a small contribution, and one can see that the new kernel leads to the same intercept of the BFKL Pomeron since it does not contribute to Recently, the differential form of Eq. (52) has been proposed in Ref. [24], which has the form ¶ : Neglecting the non-linear term at large b and using the function N 2 (ω; γ, Q T ) we can re-write Eq. (53) in the form: with α eff = κ/q 2 0 . In Eq. (54) we replace 1/q 2 by 1/q 2 0 where q 0 is the smallest momentum where one can use the perturbative QCD approach.
Hence the solution takes the form: Therefore, we see that indeed we have both diffusions, but the diffusion in b does not have b 2 ∝ √ Y . Eq. (55) leads to the solution with a black disc, whose radius can be found from the equation: As we have demonstrated in section II-A Eq. (57) results in the Froissart bound; The main problem with Eq. (52) and Eq. (53) is, that they do not reproduce the correct b 2 ∝ √ Y which as we have discussed, arises in QCD. ¶ We take the liberty of presenting this equation in a slightly different form than in Ref. [24]. This form coincides with the original one in the diffusion approximation for the BFKL kernel.
At first sight Eq. (52) can lead to a different solution, but for N 2 (ω; γ, Q T ) it takes the following form: The solution to Eq. (58) is In Eq. (59) we summed over n using the method of steepest descend and used , that at Q T = 0, the ω(γ) spectrum is the same as in the BFKL equation.

D.
The field theory realization of the parton model: the BFKL Pomeron in spontaneously broken QCD in 2+1 space-time dimensions In Ref. [8] the BFKL Pomeron is studied in spontaneously broken 2 + 1 QCD, The coupling constant of this theory g has dimension of mass, therefore, this is a super renormalizable theory which is, actually, a theoretical realization of the parton model, that has been discussed in section II-A. In this theory we do not have the ultraviolet divergencies, which results in the absence of the scaling violations, in the behaviour of the deep inelastic structure functions. On the other hand, the infrared singularities in the massless limit are stronger than in 3+1 QCD. Being such, this theory is suitable for discussing the influence of the infrared singularities on the Pomeron structure. In particular, we expect Gribov's diffusion in the impact parameters.
The BFKL equation in this theory has the same form of Eq. (27) with replacement and with the change of the kernels of Eq. (28) due to inclusion of the mass of gluon m. It can be written as where the gluon trajectories are where α S = g 2 N c /(4πm) (N c denotes the number of colours).
(63) Eq. (61) is solved analytically in Ref. [8], both in coordinate and momentum representations. It turns out that the rightmost singularity in ω is the Regge pole with the trajectory α I P (Q T ) = 1 + ∆ (2) − α (2) The dependence on the sizes of the scattering dipoles or on their transverse momenta, only enters the vertices of the interactions of this Regge pole with the particles. Therefore, in this theory we obtain the typical soft Pomeron of the parton model, with the diffusion in b but without random walk in ln (p T ).

E.
The hierarchy of scales: the interface of soft and hard Pomerons.
In the previous section we saw, that the infrared divergency could generate the soft Pomeron: the Reggeon with the intercept larger than 1 and with the trajectory of the type of Eq. (64). However, such a soft Pomeron should exist simultaneously with the BFKL Pomeron, and we have to find the high energy asymptotic behaviour of the scattering amplitude, taking into account the interface of these two Pomerons. The second problem is to find a new dimensional scale, which in the spontaneously broken 2 + 1 QCD is the coupling constant g, which has the dimension of mass.
In Ref. [9] it is shown that such a scale exists in the QCD vacuum as a consequence of scale anomaly, and it is related to the density of vacuum gluon fields due to semi-classical fluctuations: numerically, M 2 0 4 ÷ 6 GeV 2 . Fig. 4-b and Fig. 4-c illustrate how this new scale enters the structure of the soft Pomeron. Indeed, among the higher order, O(α 2 S ) (α S = g 2 /4π) corrections, to the BFKL kernel (see Fig. 4-b). We isolate a particular class of diagrams, which include the propagation of two gluons in the scalar color singlet channel J P C = 0 ++ . The vertex of their production, generated the four-gluon coupling in the QCD Lagrangian, is ∼ α s F µν a F a µν . We observe, that this vertex is proportional to the trace of the QCD energy-momentum tensor ( θ µ µ ) in the chiral limit of massless quarks: The contribution of Fig. 4-b, therefore, appears proportional to the correlator of the QCD energy -momentum tensor, for which we have   (Fig. 4-a) and to the soft Pomeron ( Fig. 4-b -Fig. 4-d). The helical lines describe the gluons, the dashed lines correspond to pions. Fig. 4-d shows the contribution of the QCD instantons to the soft Pomeron structure.
The disappearance of the coupling constant in the spectral density of the g 2 F 2 operator can be explained, if we assume, that non-perturbative QCD vacuum is dominated by the semi-classical fluctuations of the gluon field. Since the strength of the classical gluon field is inversely proportional to the coupling, F ∼ 1/g, the quark zero modes, and the spectral density of their pionic excitations, appear independent of the coupling constant. Summing the ladder diagrams of Fig. 4-c we obtain the soft Pomeron with the intercept: In Ref. [10] the particular type of these classical fields are considered: we will assume, that the field is given by the instanton solution. It is shown, that these fields lead to the soft Pomeron (see ladder diagram of Fig. 4-d). The new dimensional scale in this particular approach is the typical size (ρ 0 ) of the instanton, which comes from the lattice calculation (ρ 0 ∼ 0.3f m). Summing the diagrams of Fig Bearing in mind the hierarchy of the scales shown in Fig. 5, one can see, that the equation for the resulting Pomeron Green function takes the form [14] (see Fig. 6 In Eq. (70) K BFKL (r , r") is familiar kernel of the BFKL equation (see Eq. (30)), while K S denotes the kernel of the soft Pomeron which we take in the form: with ∆ (Q T ) = ∆ I P − α I P Q 2 T (see Eq. (69) ) and D (r) is the function, which decreases with r for r > r soft = 1/M 0 , where M 0 ≈ 2 GeV , as we have discussed above. In addition we normalize d 2 r D (r) = 1.  Fig. 6-a). G 0 (ω, ξ, QT ) denotes the Green function for the exchange of two gluons while G H (ω, ξ, QT ) is the Green function of the BFKL Pomeron. Fig. 6-b illustrates the solution of Eq. (73).
The solution to this equation, which is clear from Fig. 6-b, has the form (see appendix B of Ref. [14]): The new singularities of the Green function stem from the equation: Integrals over r and r" lead to r = r ≈ r soft , and the Green function of the BFKL Pomeron (see Eq. (36) and Eq. (37)) is equal to . Therefore, Eq. (74) reduces to From Eq. (75) we see, that the resulting intercept increases, which means that the soft Pomeron influences the energy behaviour of the deep inelastic processes. However, for this paper the most significant result is the impact parameter behaviour of the Green function.
Function G I P (ω, Q T ; r 1 , r 2 ) has been found in Ref. [3], and it takes the following form (see Eq.35 in Ref. [3]) From Eq. (76) one can see, that in the kinematic region of small We will restrict ourselves to the region of small Q T , since we are interested in the large impact parameter dependence of the Green function.
Introducing a new variable t = √ δω − ∆ I P we can re-write Eq. (77) in the form: whereξ i = ξ/ √ D andα I P = 2 √ Dα I P . Taking the integral of Eq. (78) using the method of steepest descent, we obtain the following equation for the saddle point value: Substituting t SP into the conditions for large b we obtain: Plugging this t SP into Eq. (78) one can see that . Two lessons that we can learn from Eq. (81): (i) the random walk in b leads to b 2 ∝ 2α I P √ Y and (ii) that the eikonal unitarization leads to the Froissart disc with radius One can see that Eq. (82) satisfy the conditions of Eq. (80). Note, that we obtain the value of b 2 and the radius of the Froissart disc in accord with the general analysis in section II-B-4.

F. Summing pion loops
In Ref. [23] a theoretical approach is proposed, which is based on the assumption, that the BFKL Pomeron, being perturbative in nature, takes into account rather short distances (say of the order of 1/m G , where m G is the mass of the lightest glueball); and the long distance contribution can be described by the exchange of pions (1/2µ 1/m G , where µ is the pion mass). Indeed, the proof of the Froissart theorem indicates, that the exponential fall-off at large b is closely related to the contribution of the exchange of two lightest hadrons: 2 pions, in t-channel (see Fig. 7-a). Fig. 7-b demonstrates both the main idea and the practical method of our approach. The set of diagrams of Fig. 7 can be summed, using the following equation: where G I P denotes the resulting Green's function of the Pomeron with the pion loops for r 1 = r 2 = R π where R π is the size of the pion. G BFKL I P is the Green function of the BFKL Pomeron.
The t-channel unitarity for G I P (ω, Q T ; R π , R π ) is used to calculate Σ, as it was suggested in Ref. [40]. The unitarity constraints for G I P takes the form [28,30] Eq. (85) is written for t > 4µ 2 where µ is the mass of pion. We wish to only sum the two pion contribution, assuming that all other have been taken into account in the BFKL Pomeron. In other words, we believe that Eq. (85) can be trusted for t < m 2 G where m G is the lightest glueball. ρ (ω, t) is equal to (t − 4µ 2 ) where vertex Γ I P ππ describes the interaction of the BFKL Pomeron with π-meson and it has been evaluated in Ref. [23]: Γ I P ππ = = 8 9ᾱ S R π . Experimentally < R 2 π >= 0.438 ± 0.008 f m 2 [42] and we will bear this value in mind for the numerical estimates.
We use dispersion relations to calculate Σ. Since we believe, that non-perturbative corrections will not change the intercept of the BFKL Pomeron (see Refs. [21,22]) we write the dispersion relation with one subtraction. It takes the form Factor 3/2 stems from two states in t-channel: π + π − and π 0 π 0 . The integral over t for ω 1 has been taken in Ref. [40]. For fixed ω where F 1 is the Appell F 1 function ( see Ref. [41] formulae 9.180 -9.184).
The linear term in Q 2 T can be easily found since two arguments in F 1 coincide, and we can use the relation (see Ref. [41] formula 9.182(11)) and reduce Eq. (88) to the form For m G µ Eq. (89) takes the form We evaluate α ef f (ᾱ S ) using R 2 π = 0.436f m 2 and the value for m 2 G = 5 GeV 2 that stems from lattice calculation of the glueball masses [43] for the Pomeron trajectory. In principle we can consider this coefficient as the only nonperturbative input that we need, but it is pleasant to realize that simple estimates give a sizable value of 0.1 GeV −2 forᾱ S = 0.3 (see Ref. [23]).
Comparing Eq. (84) with Eq. (77) we see that the difference is only that in Eq. (84) we put ∆ I P = 0, and the new dimensional scale is determined by the mass of the lightest glueball m G . Bearing this in mind and using Eq. (81) we obtain: with r soft = 1/m G . In Eq. (91) we replace Eq. (83) by a more general Eq. (70) with the solution of Eq. (73). We believe that Eq. (73) gives the solution to the problem of large b behaviour.
As we have mentioned, obtaining this equation we used, based on the results of Refs. [21,22], that the nonperturbative corrections do not change the intercept of the BFKL Pomeron. One way to check this key assumption is to try to describe the experimental data for the soft interaction at high energy, assuming, that we have only the BFKL Pomeron. In Refs. [44] the model for the soft interaction was built, which based on two main ingredients: (i) CGC approach for the scattering amplitude of two dipoles; and (ii) phenomenological description of the hadron structure. In this model, we successfully describe the DIS data from HERA, the total, inelastic, elastic and diffractive cross sections, the t-dependence of these cross sections, as well as the inclusive production and rapidity and angular correlations in a wide range of energies, including that of the LHC.

G. Resume
Concluding this brief review, we believe, that Eq. (91) incorporates correctly the QCD random walk in the impact parameters (b) into the BFKL equation. We would like to emphasize, that two main ingredients formed the basis of this conclusion: the assumption, that the non-perturbative modification of the BFKL Pomeron does not change its intercept; and our result of section II-B-4, that b 2 for the QCD diffusion in b is proportional to √ Y . Insufficiency of Eq. (91) is, that it gives the modification of the Green function of the perturbative BFKL Pomeron, but it is not written as an alternative for the non-linear Balitsky-Kovchegov equation. On the other hand, for discussion of the asymptotic behaviour of the scattering amplitude at high energy, we need a new evolution equation, which will include on the same footing both the shadowing corrections, and the correct behaviour at large impact parameter.

IV. A NEW NON-LINEAR EVOLUTION EQUATION
We proposed the following generalization of the Balitsky-Kovchegov non-linear evolution equation: Where α eff is a new dimensional scale. We believe that Eq. (92) is the correct way to introduce this scale in the non-linear evolution.
We can obtain the solution to this equation using it in the form: is the Green's function of the BFKL equation and Φ (Y ; b) satisfies the following equation In ω-representation: Eq. (93) takes the form The general solution to Eq. (95) can be written (see Ref. [51] formula 9.4.3) as where Φ ± are the solution to the following equations: We restrict ourself by the solution Φ (b, ω) = Φ + (b, ω) for which the Green's function is equal to where φ in + (ω) should be found from the initial condition at Y = 0. Searching for the Green's function of Eq. (98) we impose the initial condition: which results in φ in . The similarities between this equation and Eq. (77) and Eq. (84) are obvious. Taking integral over d 2 Q T we obtain: and For large b we can take the integral over ω plugging the asymptotic behaviour of the K 0 function and using the method of steepest descend. Eq. (101) takes the form: The equation for the saddle point has the following form From Eq. (103) we obtain for the integral over ω: and Φ (Y ; b) is given by Eq. (104). In Eq. (106) we integrated over γ using the diffusion expansion for the kernel, χ(γ) for γ → 1 2 . For the region of not large b we obtain the solution which has a cumbersome form: where 0 F 2 is the hyperbolic function (see formula 9.14 of Ref. [41]). As we have mentioned above the solution of Eq. (104) is the Green's function of Eq. (93), and it preserves the remarkable property of the Green's function of the BFKL Pomeron: In Eq. (108) we use Eq. (98) for Φ + (b, ω) which leads to the following representation for Φ (b, Y ) One can see, that the integration over d 2 b in Eq. (108) generates δ(2) Q T − Q T , which results in Eq. (108). Using Eq. (108) we can prove that for the Green's function of Eq. (106) we have the following equation: where ξ = ln r 2 /R 2 . The integration over ξ and b gives δ (γ − γ ) and δ (2) Q T − Q T , which leads to tequation Eq. (111).
One can see that Eq. (106) has the same main features as Eq. (91). In particular, the unitaristion leads to the Froissart disc with radius The behaviour of the scattering dipole amplitude in the vicinity of the saturation scale can be found from the solution of the linear equation [35,45]. We will use the solution to Eq. (92) in the form of Eq. (104), and will take the integral over γ using the method of steepest descend. The equation for the saddle point in γ has the form: The second condition is that the scattering amplitude is a constant for r 2 Q 2 s = 1, which has the form Plugging Eq. (113) into Eq. (114) we obtain the following equation for γ SP .
For Z = 0 the solution to Eq. (115) is γ SP = γ cr ≈ 0.37. If Z ᾱ S Y we can find the solution to Eq. (115) assuming, that γ SP = γ cr + δγ with δγ γ cr . The value of δγ from Eq. (115) turns out to be equal to Substituting γ SP = γ cr + δγ into Eq. (113) we obtain the expression for the saturation momentum: In the vicinity of the saturation scale the scattering amplitude shows geometric scaling behaviour [46][47][48][49] and has the following form [45] for τ ∼ 1: where N 0 is a constant and we assumed that δγ ln τ 1. Using Eq. (117) for Q s we see that N (r, b, Y ) of Eq. (118) can be re-written in the form: for the kinematic region C. The scattering amplitude in the deep saturation region We use the approach of Ref. [47] for linearizing Eq. (92) deep inside the saturation region (r 2 Q 2 s (Y, b) 1): we assume that N (r, b, Y ) → 1 and neglect the contributions, which are proportional to ∆ 2 (r, b, Y ) 1 where N (r, b, Y ) = 1 − ∆ (r, b, Y ). Eq. (92) reduces to the following linear equation for ∆ (r, b, Y ) In the previous section we found, that the scattering amplitude for τ ∼ 1 shows geometric scaling behaviour, being a function of only one variable τ . The first term in the r.h.s. of Eq. (121) is equal to −ᾱ S ln τ ≡ −ᾱ S ζ. We can re-write ∂ ∂Y as Assuming that 1 4(1−γcr) Y Z 1 we can neglect this contribution and re-write Eq. (121) in the form where The solution to Eq. (123) can be written in the form (see section III-A): This solution violates the geometric scaling behaviour of the scattering amplitude, leading to the additional suppression of the scattering amplitude at large b. However, it should be emphasized, that the main dependence on the impact parameter is included in b dependence of the variable ζ (see Eq. (124)). For the BK equation the asymptotic behaviour of the scattering amplitude is determined by ∆ (ζ, b) = Const exp − ζ 2 2κ , with a constant which we cannot estimate. The solution of Eq. (125) states that instead of Const we can have function Φ (Y, b). However, in the next section we will argue that actually we can replace Φ (Y, b) by 1. The general non-linear evolution that is given by Eq. (92) is difficult to analyze analytically for the full BFKL kernel of Eq. (30) or Eq. (3). This kernel includes the summation over all twist contributions. We would like to start with a simplified version of the kernel in which we restrict ourselves to the leading twist term only, which has the form for τ = rQ s < 1 summing (ln (1/(r Λ QCD ))) n ; 1 1 − γ for τ = rQ s > 1 summing (ln (rQ s )) n ; (126) instead of the full expression of Eq. (30).
As indicated in Eq. (126) we have two types of logs: ᾱ S ln (r Λ QCD ) n in the perturbative QCD kinematic region ) denotes the saturation scale. To sum these logs it is necessary to modify the BFKL kernel in different ways in the two kinematic regions, as shown in Eq. (126).
Inside the saturation region where τ = r 2 Q 2 s (Y, b) > 1 the logs originate from the decay of a large size dipole into one small size dipole and one large size dipole. However, the size of the small dipole is still larger than 1/Q s . This observation can be translated in the following form of the kernel Inside the saturation region Eq. (92) takes the form The advantage of the simplified kernel of Eq. (126) is that, in the Double Log Approximation (DLA) for τ < 1, it provides a matching with the DGLAP evolution equation [50]. Introducing we can re-write Eq. (128) in the form: Assuming that φ (ζ, b) = φ (Y, ξ) +φ (Y, b) we can re-write Eq. (130) as follows: Using Eq. (122) and Eq. (124) with We can reduce Eq. (131), applying ∂ ∂ξ to the both sides of this equation, to the following equation Introducing dφ(ζ) dζ = F (φ) we can re-write Eq. (133) in the form Finally, Therefore, in this simplified version of the non-linear equation we found, that the scattering amplitude in the saturation region has the geometric scaling behaviour being a function of one variable τ = r 2 Q 2 (Y, b) with Hence we find that the b dependence is concentrated in the dependence of the saturation scale, as it has been assumed in the numerious attempts to introduce b-dependence in the saturation models (see Ref. [52,53] and references therein).

V. THE FROISSART DISC IN QCD.
In sections IV-C and IV-D we showed, that (i) the scattering amplitude shows geometric scaling behaviour being a function of one variable ζ(see Eq. (124)) and (ii) 1 − N (ζ) = ∆ (ζ) 1 for ζ ≥ 1. The condition ζ ≥ 1 means that Therefore, for where R F is the radius of the Froissart disc. In Fig. 8 we plot the solution to the leading twist non-linear equation (see previous section) as function of b. The solution demonstrates all features of the Froissart disc with a radius, which is proportional to Y . The edge of the disc decreases as exp −C 1 (b − R F ) 3 where C 1 is a constant that can be calculated. Indeed, we can see this plugging in the expression for the amplitude N (r, Y ; b) = N 0 r 2 Q 2 s (Y, b) γ b = R F + ∆b and expanding the amplitude with respect to ∆b. One can see, that for ∆b R F , the amplitude is equal which describe the numerical solution within the accuracy of 4% (see Fig. 9).

VI. CONCLUSIONS
The main result of this paper is the new non-linear evolution equation which includes the random walk both in the transverse momentum and in impact parameter of the produced gluons (dipoles). We showed, that the solution to new equation results in the exponential decrease of the scattering amplitude at large impact parameter and in the restoration of the Froissart theorem. Therefore, this equation solves the long standing problem of the CGC/saturation approach. We demonstrated, that the new equation generates the amplitude, which approaches 1 for b ≤ Const Y and which decreases as exp (−µb) at b > Const Y .
We gave a brief review of the previous attempts to include the diffusion in b, in the framework of the BFKL equation, and showed that the new equation provides a kind of the generalization of all these attempts.
We found the solution to the equation in the kinematic region, where Z 4γ Y . Hence, a problem for the future is to develop a more general solution, as well as to describe the available experimental data within an approach, based on the new equation.