{'fraction_non_alphanumeric': 0.03974970834659031, 'fraction_numerical': 0.012217626471524022, 'mean_word_length': 4.717438758185787, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Universality of QCD traveling-waves with running coupling

The Balitsky-Kovchegov QCD equation for rapidity evolution describing saturation effects at high energy admits universal asymptotic traveling-wave solutions when the nonlinear damping becomes effective. The asymptotic solutions fall in universality classes depending only on some specific properties of the solution of the associated linear equation. We derive these solutions for the recent QCD formulations of the Balitsky-Kovchegov equation with running coupling constant obtained from quark-loop calculation. While the associated linear solutions depend in different ways with observables and higher-order effects, we show that the asymptotic traveling-wave solutions all belong to the same universality class whose solutions are given. Hence the influence of saturation stabilizes the QCD evolution with respect to higher order effects and leads to universal features at high enough rapidity, such as the form of the traveling waves, the intercept of the saturation scale and geometric scaling in square-root of the rapidity.


I. INTRODUCTION
In the leading-logarithm (LL) approximation of high-energy (high-density) QCD, the evolution with rapidity Y of deep-inelastic scattering cross-sections is driven by the nonlinear Balitsky-Kovchegov (BK) equation [1][2][3]. This equation is supposed to capture essential features of saturation effects in the "mean-field" approximation where fluctuations (or "Pomeron-loop" effects [4]) can be neglected. More specifically, for the dipole-target amplitude in 2-dimensional transverse position space N (x, y, Y ), it reads where |xy| 2 ≡ (x−y) 2 , x and y are the transverse space positions of the quark and antiquark constituting the QCD dipole. In Eq.(1), the couplingᾱ = α s N c /π is kept fixed, since the derivation of the equation is made at leading logaritmic approximation, and thus zeroeth perturbative order in the coupling. Indeed, this nonlinear equation corresponds to resumming QCD "fan diagrams" in the LL approximation [2,3].
In the 1-dimensional case, without impact-parameter dependence (e.g. for large targets), N (x, y) ≡ N (|xy|), the equation acquires a simple form in the Fourier-transformed momentum space, namely [3] where L(k) = log(k 2 /Λ 2 ), and Λ is here an arbitrary scale. The Fourier-transformed amplitude can be related by convolution (see e.g. [5]) to the unintegrated gluon distribution in the 1-dimensional space in transverse momentum k ≡ |k|. In this LL approximation the characteristic function of the kernel has the standard Balitsky-Fadin-Kuraev-Lipatov (BFKL) form [6][7][8], namely Quite recently, the extension of the BK equation to running coupling has been the subject of interesting theoretical studies. More generally, the extension of the saturation formalism beyond LL order is of present concern. The goal of the present paper is to apply the traveling-wave method to the problem of solving the nonlinear BK equation with running coupling. This method allows to find asymptotic solutions to nonlinear equations of the BK type and to discuss their universality properties.
Indeed, we will consider the recent theoretical advances towards a full QCD formulation of the BK equation at next-to-leading logarithmic (NLL) order in the coupling constant and beyond. This implies in particular taking into account in a proper way the running QCD coupling and NLL and higher order contributions to the kernel. The first part of our study is devoted to the equations incorporating the terms determined through the quark-loop contributions to the BK equation [9][10][11] (see also related older [12] and recent [13][14][15] references). Then, we consider the full NLL contributions to the kernel [16,17] and the higher-order contributions implied by the renormalization-group constraints [18][19][20][21].
Since a few years, the solutions of saturation with running coupling have been considered in the literature. Already in the papers [22,23], the dominant asymptotic term of the solution has been obtained from the knowledge of the BFKL equation. Using this equation with absorbing boundary conditions, the subasymptotic term could also be derived [24] and extended to NLL BFKL kernels [25].
Our tools for the analysis rely on the powerful traveling-wave method. Indeed, it has been remarked that, contrary perhaps to naive expectation, the nonlinear character of the BK equation leads to some simplification w.r.t. the linear BFKL equations. Interesting universal properties, i.e. asymptotic solutions which do not depend either on initial conditions, on the precise form of the nonlinear terms or on details of the linear kernel. This is tightly linked to a mathematical property of a well-known class of nonlinear equations [26,27] which admits asymptotic solutions in terms of traveling-waves [28][29][30]. The relevance of this class of equations and traveling-waves for the BK equation has been raised in [31,32].
There exist previous analytical results relevant for our study on the BK equation with running coupling. Heuristically, the running coupling has been considered together with the 1-dimensional BK equation (2) by substituting the fixed coupling in front of the equation with a running couplingᾱ(L) depending on the gluon momentum scale [24,32]. One writesᾱ We define L(k) = log(k 2 /Λ 2 ), where Λ is now the QCD scale. One considers the LL kernel function as in Eq.(4) and the running coupling readsᾱ In this case, the asymptotic universal solution has been found, either by approximation of the linear regime with absorbing conditions [24,25] or directly [31,32] from the traveling wave method of solving equation (5). In Ref. [9], the quark-loop contributions to the small-x evolution has been computed for the dipole-target amplitude in transverse-position space. Through a separation scheme between the running coupling at next-to-leading order level and the kernel at NLL level, an evolution equation has been derived. For the 1-dimensional amplitude N (|xy|, Y ) the correspondingly modified BK equation (see [9], formula (52) where the equation is in fact written using traces over Wilson-line operators, easily translated in terms of N ) reads whereᾱ Note that the expression (7) needs some regularization procedure, since Landau poles appear in the z-integration over α(|xz|) andᾱ(|yz|). We will consider this important point in detail later on for our study. In the same time, another consistent derivation of the BK equation with running coupling constant has been performed both for the dipole amplitude in position space and for the unintegrated gluon distribution in momentum space [10,11]. The same quark-loop contribution as in [9] has been found, but using a different separation scheme. The quark-loop contributions have been splitted differently into a contribution to the running coupling and a contribution to the kernel. The equation found in position space is thus different from (7). The involved running couplings appear to form a "triumvirate" structure in terms of three different scales. More precisely (see [11], where formula (48) the equation is written with additional factors e −5/3 absorbed here in Λ) one writes φ(k, Y ) being the unintegrated gluon distribution in tranverse momentum space. The nonlinear terms in Eq.(9) are not explicitely written, but as already known from other studies of the BK equation in momentum space [33], their precise form is not needed for the derivation of the traveling-wave solutions in the transition region towards saturation, provided they ensure the necessary damping. This "triumvirate" of coupling constants has been advocated long ago from bootstrap properties [12], while the definition (10) has been introduced more recently [11]. Quite interestingly, the consistency of the equation for the dipole amplitude (7) with the one for the unintegrated gluon distribution has been established in [11]. However, the corresponding Pomeron intercept driving the leading rapidity behaviour of the linear equation is different in both cases at NLL order, by contrast with the identical rapidity behaviour at LL level. The question of the universality properties of the asymptotic solutions of these equations (5,7,9) when taking into account nonlinear damping terms is thus interesting to investigate. Indeed, the dependence on the separation schemes and the dependence of the Pomeron intercept on the observable valid for the linear regime at NLL order and beyond have to be reexamined when including saturation effects.
Concerning the scheme dependence related to the renormalisation-group constraints, one starts with the expectation that the linear regime of the BK equation is driven by the full NLL BFKL kernel already available from the calculation of [16,17]. It is known that one introduces a "renormalization-group improved" (RG-improved) kernel at all orders of perturbation. One has to get rid of unwanted spurious singularities brought together with the NLL calculation. This boils down to the existence of different RG-improved schemes which are equivalent at NLL accuracy but with different resummations at higher orders (see, e.g., the S3 and S4 schemes [18] and the CCS scheme [20]).
The traveling-wave method for the BK equation with the NLL and higher-order kernels have already been used in Refs. [34,35]. In Ref. [34], the NLL kernels have been considered together with a fixed value of the QCD coupling. In this case also the method works for finding suitable asymptotic solutions, and they depend on both the value of the coupling and on the NLL framework. In some sense, our study is a continuation of Ref. [35] where the RG-improved scheme dependence is studied for the heuristic equation (5). Note that using the fixed coupling in the definition of the RG-improved scheme leads to a scheme dependence of the traveling waves [35]. We want to reanalyze the non linear problem when using a running coupling for the scheme definition, as done for the corresponding linear problem in [20,21].
It is important to note that, for the Balitsky scheme, Eq.(7), as well as for the "triumvirate case", Eq.(9), the question of spurious singularities and of renormalization-group improved kernels have not yet been studied. Moreover, the regularization of Landau poles may introduce some extra subtelties in this problem. We will show, by a careful study of singularities in the NLL kernels, that such collinear singularities do not change the determination of the universal terms in the asymptotic expansion but they are responsible of a slowing down of the convergence of the asymptotic expansion and probably introduce stronger nonuniversal subasymptotic contributions. Such subasymptotic contributions have been discussed e.g. in Ref. [25]. They certainly deserve a separate further study in the framework of traveling waves.
The plan of the paper is as follows. In section II we consider the BK equation (5) with momentum-space running coupling and recall the known results on its traveling-wave solutions [32,35], including a new discussion of the effect of NLL kernel singularities on the universality properties.In section III we turn to position space and derive the corresponding linear-regime and traveling-wave solutions for the Balitsky scheme of [9]. In section IV we consider the "triumvirate " case [10,11], and exhibits its solutions in the same way. In section V we discuss the RG-improved scheme dependence. In section VI we summarize our results by concluding that all corresponding traveling-wave solutions fall into the same universality class and indicate the corresponding asymptotic predictions.

A. Solution of the linear equation
Let us first remind briefly the method used in Ref. [32] to obtain the traveling-wave solutions of the Balitsky equation (5). Following the general method [32], we first write the solution to the linearized version of the equation. In the saddle-point approximation, it has the form of a double Mellin transform [20] with the kernel dependence appearing through the function γ being an unspecified constant. Indeed, using the saddle-point method for the integration over γ at large enough L, one gets the equation or equivalently in the form of an operator acting on the amplitude Hence the solution (11) verifies the restriction to linear terms of (5), within the saddle-point approximation.
As a next step one performs the saddle point integration over ω in the limit of large Y . The saddle point ω s is given by behaving like ω s ∼ Y − 1 2 . With this form of ω s the gluon density is given by where the time variable is interpreted as t = √ Y and B. Traveling-wave solutions A critical group velocity (defined as the minimum of the phase velocity in the wave language) is obtained as However, γ c determined in such a way still depends on the arbitrary constantγ. Thus, requiring υ g to be independent on the choice ofγ means dυ g (γ)/dγ = 0, which in turn gives dυ g (γ c )/dγ c = 0 since the dependence of the velocity on γ comes through γ c only. Applying this condition to Eq. (17) one gets eliminating all dependence on the arbitrary constantγ in the definition (12). As a consequence, one finds the equation for the critical exponent which is the well-known value (γ c = .6275) associated to the saturation solutions using the LL BFKL kernel (4). Finally, using an ansatz technique borrowed from statistical physics [29], one finds [32] the result for the gluon density in agreement with [24]. ξ 1 = −2.338 is the first zero of the Airy function Ai(ξ) and the saturation scale has the form up to a non universal multiplicative constant. Hence the two first terms of the expansion a large t of d log(Q 2 s )/dt (and thus of the saturation intercept d log(Q 2 s )/dY ) are completely specified. Note that the result (21) is known to be valid in the transition region towards saturation but not in the deep infrared saturation domain. However, the infra-red region is also universally constrained by unitarity, namely N ∼ log(Qs/k).
The important output of this traveling-wave solution with running coupling constant α(L) is that the amplitude (21) depends only on the ratio k 2 /Q 2 s (t) at asymptotic values of t. It thus corresponds to a geometric-scaling property [36] of the gluon density as a function of t ≡ √ Y [37], instead of Y in the non-running case. It is useful to note that the ansatz leads to the solution (21,22) which do not depend on the saddle-point approximation used initially in (11). Indeed, all parameters of the ansatz are determined self-consistently from the asymptotic analysis of the original nonlinear equation. The saddle-point approximation has been useful to determine the appropriate form of the ansatz [32].

C. Extension of traveling waves to a ω-dependent kernel
In the following, and throughout the rest of the paper, we introduce the possibility of dealing with a more general kernel χ(−∂ L , ∂ Y ). This generalization comes naturally from the NLL studies. For instance it has been introduced in the ω-expansion method of Ref. [19], in order to define RG-improved kernels. More generally, it is useful in order to perform a change of variableᾱ(L) → ω in the linearized equation, which is the starting point of the traveling wave method. It allows to express the solution of the linear equation in terms of a linear superposition of waves. The nonlinear term will then select the critical wave as the traveling wave solution.
Let us recall for further use the extension [35] of the traveling wave method in the case of a kernel In double Mellin space, the linear part of the equation reads bLω = χ(γ, ω). Introducing the saddle-point integration over ω leads to the asymptotic solution in Y where ω s is given by the saddle-point equation and the "dot" means the derivative with respect to ω. The solution ω s (γ, Y ) of Eq.(26) defines the asymptotic properties of the amplitude (25). In order to solve it, one has to take into account the analytic properties of the kernel integral (24), which posesses singularities near γ = 0, 1. The solution will depend on the position of the critical value γ = γ c for the saturation solution with respect to the singularities. Indeed, as we shall now prove, one has to avoid some well-defined regions around the singularities in order to find the universality properties.
Let us first consider the situation where γ c is not in the neighbourhood of the singularities of the kernel. Then, the expansion of the integral of the kernel (24) and its derivative at ω = 0 is analytic near γ ∼ γ c . One can use the Taylor expansion around ω = 0 and by collecting terms with the same powers of ω s , Eq. (26) writes Following [35], the term in braces will contribute only to subleading nonuniversal terms in the traveling wave expansion while the left-hand term with the second derivative w.r.t. ω can be absorbed by a translation in Y. This term (and the higher derivative terms) depends thus on initial conditions, and we have the universal asymptotic behaviour at large rapidity Hence the previous relations (15)(16)(17)(18)(19)(20) remains valid once substituting χ(γ) → χ(γ, ω = 0). Now, we have to take into account the singularity regions, in order to discuss the validity range of the previous derivation. The singularities of the function X(γ, ω) depend crucially of the use of a renormalization-group improved kernel or not, due to the presence of spurious singularities in NLL contributions of the kernel, which have to be cancelled by the renormalization-group improvement scheme. We shall however consider also the case when spurious singularities are still present, by comparison.
• Renormalization-group improved kernels. The behaviour of the kernels χ near the singularities are simple poles (apart possibly from the quark contributions, which we shall not consider in this derivation). By integration, see (24), this leads to logarithmic singularities in the function X(γ, ω) and thus in Eq. (26). Hence, the singularity Inside this region, the expansion (28) and thus the universal hierarchy in ω s is spoiled.
• NLL kernels with spurious singularities. The NLL kernels χ contain singularities up to triple poles. By the integration in (24), new single and double-pole singularities appear at next-to-leading order. Hence, the singularity region to be avoided becomes Inside this region, the expansion (28) and thus the universal hierarchy in ω s is spoiled. Thus the useful region of analyticity requires a much higher rapidity to be valid.
This discussion shows that the universality properties will show up at different scales in rapidity, depending on the singularity structure of the NLL kernels. In fact, if one considers NLL contributions keeping spurious singularities, subasymptotic nonuniversal terms may be dominant up to high rapidities. Hence, the singularity structure of the NLL kernels will affect the subasymptotic behaviour.
We note that the running coupling case, introducing the integration (24), softens the singularities. It may also happen that the ω-expansion itself softens the singularities (even the spurious ones, if present), hence possibly improving the rapidity range of universality properties. This, and the question of quark contributions deserve a more detailed study. In the present work we will focus on the case where the critical γ c is well inside the analyticity region. Since γ c takes the same universal value γ c = .6275 as for the LL kernel, the universality properties will appear when the region of singularities is avoided and thus for large enough rapidities. The price to pay for the eventual presence of spurious singularities is that the convergence to the universal features may be sizeably pushed away in rapidity.
In order to finally get the universal asymptotic terms of traveling-wave solutions, one has to express the linearized version of the BK equation (23) retaining only terms in the kernel relevant for the asymptotc analysis. By expanding around ω = 0, γ = γ c up to the second order, one gets where the "prime" is the derivative with respect to γ. The terms in the first line of (30) correspond to the same expansion as for the ω-independent case and they contribute up to order Y −1/6 , which is enough to determine the two first asymptotic terms of the amplitude and the saturation intercept. The second line contains new terms, with derivatives in ω, which contribute only to higher order.
It is interesting to note however that the first dominant term in the second line of (30) contributes with order Y −1/3 . The nonuniversal contributions are expected to start at order Y −1/2 , since they correspond to a shift t → t + t 0 which results in nonuniversal terms of order Y −1/2 . Hence the first term in the second line of (30) gives a contribution related to the ω-dependence of the kernel, since it depends onχ which could be the remaining track of NLL effects in the universal traveling-wave solutions.

A. Balitsky's formalism and regularization
Let us now derive the traveling-wave solutions for Eq. (7). Following the abovementionned method, one first look for the solution of its linear part. However, in order to give mathematical and physical consistency, one has to introduce some regularization procedure of the integration over z in position space in order to avoid the Landau poles atL(|xz|) = 0 andL(|yz|) = 0.
Starting with solving the linear part of (7) bL where the subscript R denotes the regularization procedure, we insert the Mellin transformsN (γ, Y ) of the dipole amplitude in position space.
Applying the Mellin-transform (32) to the right-hand side of (31), one obtains the diagonalized form of the linear kernel The first term in braces corresponds to the LL kernel (4). The two other terms give an additional contribution to higher orders of perturbation which by symmetry reads In order to regularize (34), we choose to introduce a truncation of the perturbative expansion as a function of α(|xy|) ≡ [−b log(|xy| 2 Λ 2 )] −1 . We will discuss later the regularization dependence. Using a rescaling variable λ = (x−z)/|xy| and going to the complex λ-plane one writes (bᾱ(|xy|)) n log(λλ)) n−1 , where N defines the level of truncation. For instance, N = 1 corresponds to the contribution of the running coupling to the NLL term of the kernel. Note that the only remaining scale dependence comes fromᾱ(|xy|).
It is then possible to express analytically the expansion (35) term-by-term. One writes where the integral I 1 (γ, ǫ, δ) is given in Appendix A.
In particular, the contribution (N = 1) to the NLL kernel reads which agrees with the result of [9], confirmed by [11] (up to the inclusion of the factor 5 3 in a redefinition of Λ and a change γ → 1−γ w.r.t. [11]).
Finally, the saddle-point (7) gives, in double-Mellin in position space at NLL order where As well-known, the NLL contribution (39) to the kernel has a spurious singularity in 1/γ 2 . As we shall now see, this will not influence the determination of the universal terms in the traveling-wave solution. However, as explained in the introduction, it may delay the convergence of the asymptotic expansion.

B. Derivation of the traveling-wave solution in position space
Keeping in a first stage the truncation at NLL level, namely N = 1 in (35,36), the nonlinear equation in position space can be formally written One should care about the exact meaning of Eq. (40), which remains at this stage formal in two aspects. The nonlinear term N ⊗2 coming from the quadratic part of the integrant in (31) is not explicitely given and the operator χ Bal 1 (−∂L,ᾱ(L)) is to be well-defined, since it contains two non-comutative variables. Concerning the nonlinear terms, we know from general traveling-wave properties that its precise analytic form is irrelevant for the asymptotic solutions. The same remains true with the yet unknown NLL contributions to the LL nonlinear term.
For the action of the kernel in Eq.(40) we consider the kernel obtained in diagonalized form (38,39), namely By solving (38) as a function of ω, one obtains a new ω-dependent kernel κ 1 (γ, ω) which leads to the equation where we have only kept the physical root. Then Eq.(40) takes the operator form The derivation of the traveling-wave solutions comes from a simple extension of the arguments of subsection II C on the extension of the traveling-wave method to a ω-dependent kernel. Indeed, Eq.(43) is similar to (23) in section II, substituting L byL, see (8).
Taking into account the same same substitution L →L and considering formula (42), the traveling-wave solution is driven by the value of the ω-dependent kernel κ 1 (γ, ω) at ω = 0, namely We see from Eq.(44), that for ω small enough, that is for the universal terms dominant at large enough rapidity, the NLL contribution will not modify the saturation saddle-point γ c . From this general argument, we see that the traveling-wave solution is given from the LL kernel with running coupling constant in position space. By the momentum-position substitution property and the independence of the result on the precise form the nonlinear damping (despite the fact that they possess a different form in both representations) one is led to the universal form of the solution where ξ 1 , γ c are the same as in formula (21) (note that (45), as well as (21), is not valid in the infra-red. However, the infra-red region is also universally constrained by unitarity in position space, namely N ∼ 1). The saturation scale is then identical to (22) up to a non universal multiplicative constant, meaning that the saturation intercept is predicted to be the same as (22) at large rapidity.

C. Truncation and regularization independence
Let us sketch the general argument leading to the independence of the traveling-wave solution w.r.t. the truncation and the regularization scheme. For this sake, we shall use a method of "ω-expansion" introduced in the derivation of the linear BFKL equations at NLL order with renormalization-group constraints [19][20][21]. Briefly, this amounts to transform the perturbative expansion of the kernel into an equivalent expansion in ω.
The asymptotic universality property can be extended to a given finite value of N in a more general truncation of the expansions (35,36). Indeed one can write the corresponding saddle-point equation as where χ Bal N (γ, α(L)) is a polynomial of degree N in α(L). This can be considered as as an implicit equation forL as a function of ω namelyL where the function κ N verifies thus Then, coming back to the truncation equation (36), one finds that χ Bal N (γ,ᾱ = 0) = χ(γ). Then the implicit equation (48) leads to κ N (γ, ω = 0) ≡ χ(γ). One again falls in the same universality class. Indeed, this is expected from the remark that ω and α(L) are of the same order and jointly go to zero at the limit.
The same result should apply for any consistent scheme for the regularization of the Landau poles, e.g. by freezing of the coupling constants beyond some scale |xz|, |yz| > r 0 in the integral (31). The perturbative consistency of the regularization implies the existence of a convergence domain inᾱ, ω near zero. In physical terms, the traveling-wave solutions are valid in the domain of high enough rapidity to have a perturbative saturation scale and driven by the ultra-violet behaviour. Hence the nonperturbative effects related to the regularization should not be relevant. This is a natural constraint on the regularization scheme.
The regularization, which corresponds to a given kernel χ Bal R will then give rise to an implicit equation defining a new ω-dependent kernel κ R . By the same arguments, one again expects κ R (γ, ω = 0) ≡ χ(γ).

D. Traveling-wave solution in momentum space
Now we shall demonstrate that the traveling-wave solution obtained in position space in the preceeding subsection is in fact in the same universality class as the solution found for the initial BK equation (5) with running coupling in momentum space. Since a Fourier transform of the asymptotic position-space solution (21) is not mathematically justified, it is convenient to start directly by transforming the equations from position to momentum space and then find the asymptotic solutions.
Due to the previous universality properties in this space, the problem reduces to the solution of the initial equation (7) when keeping only the LL kernel: where the coupling depends only on the parent-dipole size |xy|. Then, performing the Fourier transform (3) on the two sides of the equation gives Using the identity one arrives at the equation: where an extra term appears with respect to (5), depending still on the derivative of the position-space dipole amplitude N . DenotingN the Mellin transform of N (|xy|, Y ) andN , the Mellin transform of N (k, Y ), they can be simply related (see appendix B) byN For the derivative term, by insertion of Mellin-transform, one writes the relations Using now the relation (54) between Mellin transforms, one finally finds the simple expression Following the same derivation as previously, the linear solution of the linear regime writes where X(γ) is the function considered in the LL case, see (12). This solution can thus be identified as a change of the kernel eigenvalue χ(γ) → χ(γ, ω) ≡ χ(γ) + bωϕ(γ). Hence, the universality class of the non-linear equation (57) is determined by χ(γ, ω = 0) = χ(γ). This result can also be obtained by noting that the term γ γ dγ ′ ϕ(γ ′ ) = log 2 2γ−2γ Γ(γ)Γ(1−γ) Γ(γ)Γ(1−γ) can be absorbed in a redefinition of the impact factor N 0 (γ, ω). As easy to realize from the initial saddle point equation (13), a contribution of the impact factor may only consist in a change of the non universal reference scale Q 0 .
The final conclusion is thus that the traveling-wave solutions of the Balitsky equation (7) are in the same universality class as the original BK equation (5) with LL kernel both in position and in momentum space. Note that this result extends in a non trivial way to running coupling a property which was easy to derive when the coupling is fixed and thus factorized from the Fourier transform.

IV. TRIUMVIRATE OF RUNNING COUPLINGS
Let us now apply the traveling-wave method to the evolution equation (9). As previously, one has to consider first its linear part. As in section III, one has to introduce a regularization procedure to avoid the Landau poles, except that the integration is now done in momentum space. Let us choose the same type of truncation R as for the previous equation: we expand the Landau poles denominators and truncate the resulting series.
Starting with the linear part of Eq.(9) we rewrite the triumvirate of running couplinḡ Inserting (60) and the Mellin representationφ(k, Y ) = dγ 2πi e −γLφ (γ, Y ) in the right-hand side of (59), one obtains the diagonalized form of the regularized kernel χ T ri Using a rescaling λ = q/k and going to the complex λ plane, one writes where N defines the level of truncation. For instance, N = 0 corresponds to the LL BFKL kernel, and N = 1 to the NLL term of the kernel. It is then possible to express the expansion (62) analytically term-by-term. One formally writes where the generating functional I 2 (γ, ǫ, δ) is given in Appendix C. At the NLL level one finds Once going from the functionφ (see (10)) to the unintegrated gluon distribution φ, which results in changing the factor 3/2 into 1/2, this result gives back the one found in [11]. Eq.(9) then reads in double Mellin space at NLL order where We note that this result is different than for the dipole amplitude, (39), as discussed in [11].
Keeping in a first stage the truncation at NLL level, namely N = 1 in (62,63), the nonlinear equation in momentum space can be formally written where with the nonlinear contribution in momentum space denoted by N ⊗2 (L, Y ). Again and also in previous QCD traveling wave studies at nonzero transverse momentum [33] the precise form of N ⊗2 will not matter for the universal behaviour of asymptotic solutions and thus need not be explicitely derived. As in the previous section, one has to solve (65) as a second-order equation for L. One finds bLω = χ(γ) + χ(γ) 2 + 4bωχ T ri (1) (γ) 2 .
As explained in the introduction, and valid already for the Balitsky scheme, the spurious collinear pole contained in the NLL term of (68) will not change the universal terms, while it may delay the convergence of the asymptotic expansion.
The whole derivation of traveling wave solutions parallels the one described in the previous section, except for the transformation from position to momentum space. Hence finally the traveling wave solutions are directly expressed as in (21,22).
Once again, the generalisation to a fixed perturbative truncation is expected to hold. The universality class of a regularized BK equation with the triumvirate of running couplings is the same as the equation including the running coupling with the "external" gluon momentum (5). Thus, including the effects of the transition to saturation leads to an unification of the asymptotic solutions, irrespective of the differences at the level of the solutions of the linear regime.

V. RENORMALIZATION-GROUP IMPROVED SCHEME DEPENDENCE
Let us recall why the solutions of the BFKL equation with NLL corrections require some special treatment and why this introduces a NLL scheme dependence. Indeed, although the NLL corrections are known [16,17] they turn out to be negative and large. This is due to spurious singularities which are in contradiction with constraints coming from the QCD renormalization-group properties. Indeed the cancellation of these singularities can be obtained by suitable contributions at higher orders. However, since these contributions are not yet calculated, there is a need for convenient schemes ensuring the compatibility with the renormalization group, which are not uniquely defined. If one admits that the linear part of the various nonlinear evolution equations we are studying will be driven by the corresponding BFKL kernels, we have to take into account this scheme dependence in our discussion of universality properties. A first discussion of the scheme dependence can be found in ref. [35], where different schemes were discussed for the kernel of equation (5). The aim of the present section is to extend the discussion to the equations considered in the previous sections and to consider the running coupling also in the scheme definition. Note that the RG-improvement plays already an important role in the convergence aspects of the universality properties, see section II-C.
For our discussion, one could consider at least two different classes of schemes, following Ref. [35]. In a first class, containing e.g. the schemes of Refs. [20,21], the dependence of the scheme on the higher orders of the perturbation expansion has been expressed through the ω-dependence, following the method of [19]. In this case, the higher order resummation appears in the kernel only through the dependence over the two Mellin variables, i.e. χ(γ = −∂ L , ω = ∂ Y ), where ω = O(ᾱ) drives the higher-order corrections. For this class of two-variable NLL kernels, the general argument that the traveling-wave solutions are driven by χ(γ, ω = 0) = χ(γ) applies, see section II-C.
In a second class of schemes such as S3 and S4 [18], the analysis of [35] leads to results depending on the value chosen for the fixed coupling appearing in the definition of the renormalization-group scheme. Indeed, the travelingwave solutions are driven by the kernel χ(γ, ω = 0,ᾱ = 0) = χ(γ). Hence, even if the form of the solutions is similar [35], the critical parameters are different and thus the traveling-wave solutions appear to be scheme-dependent. One practical question, for instance, is to know what are the predictions for the traveling-wave solutions when choosing the coupling defining the scheme varying with the momentum scale, see e.g. [38]. Indeed, one would prefer a unified treatment of the running coupling effects taking into account the running coupling in the RG improved kernel itself. We shall show, that in this case, and on variance with the results at fixed coupling for the regularization, the same universality class defined by the LL kernel (with running coupling) is recovered.
For this sake, we can consider in the same way all previous equations, either that corresponding to the BK equation in transverse momentum, (13), or in transverse position space (38), or for the case of the triumvirate (65), using any of the renormalization-group improved NLL kernels. Then the corresponding saddle-point equations become ω = α(L orL) χ N LL (γ, ω, α(L orL)) . (70) Now, in parallel with the previous discussion, Eq.(70) can be interpreted in all cases as an implicit equation for L or L namely where the function κ N LL verifies the implicit equation This boils down to reformulate schemes, starting from S3 and S4, by using an appropriate ω-expansion [19]. These schemes are now in the same first class discussed above. Then, the condition κ N LL (γ, ω = 0) ≡ χ(γ) is expected to be verified in all cases by perturbative consistency, and thus by the same argument, one again falls in the same universality class. Indeed, if we use the ω-expansion method to redefine the schemes S3 or S4, it is possible to show that κ N LL (γ, ω = 0) = χ(γ, ω = 0, α = 0) = χ(γ). All other schemes verifying this relation will again fall in the same universality class, defined by the traveling-wave solutions (21,22). The question remains for further study whether other types of NLL schemes, such as the one of Ref. [39] can be given a similar treatment.

VI. SUMMARY AND PREDICTIONS
To summarize, we considered three versions of the QCD evolution equation in rapidity in the mean-field approximation, i.e. the Balitsky-Kovchegov [1][2][3] equation (2), extended to take into account the running coupling of QCD. The first one considers the 1-dimensional BK equation for the dipole amplitude in momentum space, with the substitution of the fixed coupling by a running coupling in terms of the gluon transverse momentum, see (5). The two other forms come from recent theoretical calculation of the quark-loop contributions. They differ by a separation scheme for higher orders between the running coupling and the kernel. The second equation we consider is written for the dipole amplitude in transverse position space , see (7) [9]. The last one can also be written for the dipole amplitude in transverse position space, leading to a "triumvirate" of running couplings, and differs from the previous one. We study it in its equivalent formulation as the evolution equation for the unintegrated gluon distribution in momentum space see (9) [10,11]. We have also enlarged the discussion by considering the modified BK equations including the renormalization-group improvements of the NLL BFKL kernels [18,20,21].
Let us summarize our results. The saturation effects which are formulated through the nonlinear damping terms in the BK equation leads to asymptotic traveling-wave solutions for all equations. Their behaviour at high rapidity is highly universal. More precisely the form of the solution at high rapidity and the two first terms of the rapidity expansion of the the saturation intercept d log(Q s )/dY are the same. Hence they are identical to those found and derived in [32] for the initial BK equation with running coupling (5).
For establishing this strong universality property, We have examined the different types of higher order dependence which modify the solutions of the linear regime • Observable dependence: The traveling wave solutions at high rapidity happen to be independent of the considered distribution functions, by contrast with the solutions of the linear part of the equations. For instance the dipole amplitude and the unintegrated gluon distribution which lead to different rapidity dependence in the linear regime for the same scheme [11] have the same saturation intercept and traveling-wave spectrum. This leads to an unified asymptotic NLL predictions for observables, independently of its formulation in terms of dipole or gluon distributions.
• Separation scheme dependence: By computing the traveling-wave solutions for both the Balitsky scheme [9] and the Kovchegov-Weigert one [10,11], we find the same results, belonging to the universality class of the simpler equation (5) .
• Regularization dependence: The definition of equations (7,9) requires a regularization, since their consistent formulation should avoid Landau-pole singularities. Using a finite truncation of the QCD perturbative expansion, we show that the order of truncation, which affects the linear regime, does not change the asymptotic travelingwave solutions. We propose a general argument for general perturbatively consistent regularization schemes.
• RG-improved scheme dependence: The QCD coupling appears also in the definition of appropriate schemes for using the NLL BFKL kernels in order to avoid spurious singularities. The improved NLL kernels with constant coupling were shown to give scheme-dependent traveling wave solutions [34], even if the overall coupling is running [35]. If the scheme is defined using the running coupling through the ω-expansion [19], we recover the strong universality property.
• Position vs. momentum-space dependence: The universality class is the same when the equations are expressed either in position, or momentum-space formalism by Fourier transformation. This shows that the saturation intercept is the same. The form of the front is also invariant except for the infra-red (with respect to the saturation scale) regions which are not bound by universality properties. We have already mentionned that the infra-red region in both representations is also universally constrained, by unitarity, namely N ∼ 1 in position space and N ∼ log(Qs/k) in momentum.
• Rapid vs. slow convergence: By discussing the effects of the singularities contained in the NLL kernels, we have distinguished between the RG-improved kernels and those containing spurious singularities. Indeed, the rapidity region where the universal properties appear is expected to be delayed in the case where such singularities remain. It is thus important to use, whenever known, a RG-improved scheme for exhibiting the universality properties.
From our results, we conclude that the saturation effects on the rapidity evolution of the dipole amplitude or the unintegrated gluon distribution with running coupling give a stabilisation of the solutions with respect to higher order corrections at high enough rapidity. They appear to be independent of higher order contribution either to the linear kernel or to the nonlinear damping terms. however, subasymptotic nonuniversal contributions may be present.
In mathematical terms one expresses this property as the existence of a large universality class of solutions given by (21,22).
At this stage, it is possible to formulate predictions for the behaviour of the exact solutions of the evolution equations which can be obtained through numerical simulations. We can list them as follows: • The form of the traveling-wave fronts (except in the full saturation region) should converge at high rapidity for the solution of the different equations in the same representation space, either in position or momentum.
• The saturation-scale intercept d log(Q s )/dY which can be determined from these fronts should converge to the same function.
• Geometric scaling in √ Y , related to the constant limit of the wave speed should be seen from the solutions, at least in the middle of the front.
These predictions can be tested using numerical simulations of the rapidity evolution corresponding to modified BK equations with running coupling. In Ref. [40], the running coupling was heuristically introduced for models in position space and compared with the fixed coupling case. Naturally, they could not include the recent theoretical advances, but some of the conclusions are still interesting to quote and seem to fit approximately with our predictions, even if the equations are different. The traveling-wave regime is observed at high enough rapidity, the saturation scales for running coupling seem to converge, at variance with the fixed coupling case where they depend on the value of the coupling. The form of the front verifies geometric scaling.
While completing this theoretical analysis we were informed of interesting numerical simulations of the relevant equations [41] using a regularization scheme based on the freezing of the coupling. The simulations [41] seem to indicate that the traveling-wave structure is preserved at moderate rapidities and independent of the initial conditions. However, the wave parameters seem to differ from the critical ones. We predict that at higher enough rapidity this difference will decrease and disappear. Note that arguments based on subasymptotic parametric solutions of the nonlinear equations [42] show that the traveling wave structure can be maintained at moderate rapidity but with modified speed as observed in [41]. Also the form of the wave front is modified, which could give an explanation of the difference already noticed in [40] between the running and non running cases.
On a phenomenological ground, such a theoretical solution is characterized by a geometric scaling [36] in √ Y , either in momentum or position space. By contrast with fixed coupling which leads to geometric scaling in Y [24,31], this property is characteristic of the BK equation with running coupling [25,32]. Geometric scaling in Y or √ Y seem to be consistent with data on the proton structure functions [37]. However, the precise comparison with the parameters of the theoretical traveling-wave solutions show that nonuniversal terms are present, and thus nonasymptotic contributions cannot be neglected. Parametric solutions, as in [42] could help to take into account nonuniversal contributions.
As an outlook for further directions of theoretical study, we mention the extension of our investigation beyond the mean-field approximation leading to the BK equation. For this sake, one has to consider the hierarchy of QCD evolution equations including fluctuations, i.e. Pomeron-loop terms, when the coupling constant is running. Also the role of impact factors, e.g. for the coupling to the virtual photon in deep inelastic scattering has not been yet studied at NLL level. Let us also note that a third (and last) universal term in the asymptotic expansion of the saturation intercept, which is known for the fixed coupling case [31], could be the remaining track of the NLL kernels in the universal traveling-wave solutions as discussed after equation (30). Finally a specific study of subasymptotic terms is deserved. It may lead to traveling-wave properties at more moderate rapidities and thus interesting on a phenomenological point of view.
Let us relate these Mellin-transformed amplitudes by rewriting the initial dipole amplitude in momentum space (3). The integral to be computed in (63) is .