On the different forms of the kinematical constraint in BFKL

We perform a detailed analysis of the different forms of the kinematical constraint imposed on the low $x$ evolution that appear in the literature. We find that all of them generate the same leading anti-collinear poles in Mellin space which agree with BFKL up to NLL order and up to NNLL in $N=4$ sYM. The coefficients of subleading poles vanish up to NNLL order for all constraints and we prove that this property should be satisfied to all orders. We then demonstrate that the kinematical constraints differ at further subleading orders of poles. We quantify the differences between the different forms of the constraints by performing numerical analysis both in Mellin space and in momentum space. It can be shown that in all three cases BFKL equation can be recast into the differential form, with the argument of the longitudinal variable shifted by the combination of the transverse coordinates.


Introduction
The sharp increase of the structure function F 2 in Deep Inelastic Scattering process with decreasing value of Bjorken variable x was one of the biggest discoveries at HERA collider [1,2]. This rise of F 2 is being driven in QCD by the increase in the gluon density due to dominance of gluon splitting when the parton's longitudinal momentum fraction is small. Understanding the behavior of the parton densities at small x is not only crucial for the DIS processes but also for the high energy hadronic collisions in accelerators and in the cosmic ray interactions. The LHC opened up the kinematic regime which is particularly sensitive to small x values and thus the region where the gluon density dominates. The increase of the gluon density towards small x translates then into the rise of the hadronic cross sections for variety of observables with the increasing center-of-mass energy. Within the framework of collinear factorization [3] the evolution of the gluon density is governed by the DGLAP [4,5,6] evolution equations which lead to the scaling violations. Although the DIS data are well described by the fits based on the DGLAP evolution, there are hints of novel physics phenomena that can be present at very small values of x [7]. The alternative approach to the description of the parton evolution stems from the analysis of the Regge limit in QCD, where the energy of the interaction is assumed to be large s → ∞, or in the case of the DIS process, when x → 0, while the scale Q 2 is perturbative but fixed. In this limit the cross section is dominated by the exchange of the so-called hard Pomeron which is the solution to the famous BFKL evolution equation [8,9,10], see [11] for a review. The two notable differences between DGLAP and BFKL evolution equations are that in the latter the solution exhibits the power like behaviour, x −λ at small values of x and the transverse momenta in the gluon cascade are not ordered. The solution thus depends on the transverse momenta of the exchanged reggeized gluons. BFKL equation is known up to NLL order in QCD [12,13] and up to NNLL order in N = 4 sYM theory [14,15,16]. A known issue with the BFKL at the NLL order is that the corrections are large compared to the LL result, which renders the result unstable. Thus it was early realized that higher order corrections need to be resummed. In the pioneering works [17,18,19] it was pointed out that there are kinematical constraints which need to be included in the evolution. Formally of higher order, when viewed from the perspective of the Regge limit, nevertheless such corrections are numerically very large as was established in [19]. Kinematical constraint was later included in the formalism that combined DGLAP and BFKL equations, and it was demonstrated that very good description of the HERA data was obtained [20]. Various elaborated resummation schemes for the small x evolution were later developed [21,22,23,24,25,26,27,28,29,30,31,32,33,34,35]. More recently the collinear resummation was applied to the non-linear evolution equations in transverse coordinate space [36,37,38,39].
Given the increased precision of the experimental data and progress in theoretical computations, especially the information about the higher orders in N = 4 sYM theory, it is worth to revisit the assumptions in the resummed schemes. Even though the kinematical constraint is well motivated, different approximations, which are used are usually done based on general grounds (i.e. for example low z approximation) but without detailed quantitative analysis. The differences between different forms of constraints may be of subleading order but it could happen that they are quite large numerically. In context of the BFKL equation, which utilizes the kinematical constraint as an ingredient, usually one uses a particular form of this constraint, namely the one which limits the transverse momenta of the exchanged gluons. Such choice is also dictated by practical applications, since in this form the angular variable in the evolution equation is integrated over which reduces the number of variables and makes it easier to solve the equation. On the other hand, there are different forms of the kinematical constraint that appear in the literature, which put limit on the transverse momenta of the emitted gluons. Actually in the CCFM equation which receives contribution from both large and low x it is a natural choice [40]. 1 The main goal of this paper is to perform a systematic semi-analytic and numerical comparison of the small x evolution equation with different forms of the kinematical constraint, which appear in the literature.
In the analysis performed in this work we found that the constrains all generate the same structure in the Mellin space for the leading and first subleading anti-collinear poles up to NLL order in ln 1/x in QCD and in N = 4 sYM and up to NNLL in N = 4 sYM. As such they are all consistent up to the highest known order of perturbation theory. The calculations were performed both in the Mellin space and also by performing the direct solution of the BFKL with kinematical constraint in the momentum space. To this aim we have developed new numerical method for the solution of the BFKL equation with kinematical constraint which is based on reformulating this equation in the differential form. It turns out that all the constraints can be recast as effective shift of the longitudinal momentum fraction in the argument of the unintegrated gluon distribution. We see that, as expected the solutions are substantially reduced with respect to the LL case, but the differences between the different forms of the kinematical constraints are non-negligible and can reach up to 10-15 % for the value of the intercept. This can have important effect for the phenomenology. Figure 1: Schematic diagram representing gluon emission in the BFKL chain.
x and x/z are the longitudinal momentum fractions of the target's momentum carried by the exchanged gluon. k T , k ′ T and q T denote the two-dimensional transverse momenta of the exchanged gluons and the emitted gluon respectively.

Kinematical constraints
In this section we shall review the origin and different forms of the kinematical constraints that appear in the literature. The kinematical constraint in the initial state cascade at low x was first considered in works [17,41,18].
We shall follow here closely the derivation presented in [19], where the kinematical constraint was implemented both in BFKL and the CCFM equations [17,41,42,43,44]. The latter evolution equation is based on the idea of coherence [45] that leads to the angular ordering of the emissions in the cascade. The BFKL equation for the unintegrated gluon density in the leading logarithmic approximation [10,9,8] can be written as where function F(x, k 2 T ) is the small x unintegrated gluon density, k T and k ′ T = k T + q T are the transverse momenta of the exchanged gluons and q T is the transverse momentum of the gluon emitted. The longitudinal momentum fractions of the exchanged gluons are x and x z respectively. The flow of the momenta in the BFKL cascade is illustrated in diagram Fig. 1. We will also use the notation k 2 T ≡ k 2 T for the squared transverse momenta for the rest of the paper. The rescaled coupling constant is defined asᾱ s ≡ αsNc π . In the leading logarithmic approximation the integration over q T in the Eq. (1), is not constrained by an upper limit thus violating the energy-momentum conservation. Thus for example in the context of the DIS process in principle there should be a limit on the integration over the transverse momentum where Q 2 is the hard scale of the DIS process, x is the Bjorken variable and W 2 is the c.m.s energy squared of the photon-proton system in DIS.
There is a stronger constraint arising however from the requirement that in the low x formalism the exchanged gluons have off-shellness dominated by the transverse components, i.e. one keeps terms that obey |k 2 | ≃ k 2 T . The derivation of the kinematical constraint here follows [19].
The gluon four momentum is as usual decomposed into light cone and transverse components where k ± = k 0 ± k 3 . Now the exchanged gluon virtuality in these variables is equal to The condition |k 2 | ≃ k 2 T translates approximately to The emitted gluon is on-shell q 2 = q + q − − q 2 T = 0, and therefore we can express this fact as On the other hand in the multi-Regge kinematics there is a strong ordering of longitudinal momenta so Using the condition (6) and inserting it into (5) one finally obtains or as limit on q 2 T integration Now there are several approximations that can be made to this constraint. In the small z limit (9) can be approximated to This form of the approximation was used in [17] and also studied in the context of small x approximation to the CCFM evolution in [19]. The lower bound on z, i.e. z > x results in the upper bound on q 2 T < k 2 T /x providing local condition for energy-momentum conservation. Finally, (9) can be further rewritten as a condition on the transverse momentum of the exchanged gluon k ′ T . For a given value of k T a high value of k ′ T means also high value of q T . Rewriting it as and averaging over angle between k T and k ′ T and taking large k ′ T limit we get: This form of the constraint was used for example in [18] and in [20,46,29,47,48]. The nice feature of (12) is the fact that the kernel with kinematical constraint has a Mellin representation which results in a simple shift of poles in the Mellin space. In the rest of the paper we shall analyze in detail all the forms of the constraints and quantify the differences between them both in Mellin space and through direct numerical solution ot the BFKL equation in momentum space.

Comparison of improved kernel with results in N = 4 sYM
Since the works of [18,19] it is known that the kinematical constraint provides very important contribution to the NLL order and beyond in BFKL. Even though formally of subleading order, numerically this correction is very large and leads to the strong reduction of the intercept, aligning the BFKL equation more with phenomenology, for selected early works which implement kinematical corrections in BFKL see for example [20,49,50,51]. At the NLL order the kinematical constraint, when viewed as a correction in the Mellin space gives rise to the cubic poles in the Mellin conjugated variable to the transverse momentum in the BFKL eigenvalue. The N = 4 sYM eigenvalue at LL is identical to the QCD case, and at NLL the same leading cubic poles appear in both theories [52]. Since the constraint is originating from the kinematics specific to the Regge limit of the cascade, one could expect it to be universal in both theories. It is useful to explore whether or not the kinematical constraint correctly reproduces the terms in the higher order NNLL known only in the supersymmetric case [14,15,16] and what are the differences between the different forms of the kinematical constraint discussed in the literature.
In this section we shall perform a detailed comparison of leading and subleading poles in the Mellin space originating from the kinematical constraint with the results in N = 4 sYM case up to NNLL. First, we shall perform the comparison on the example of the constraint (12). Later on, in Sec. 4, we shall check that the results are consistent for different forms of the kinematical constraints up to certain order of poles, and we shall identify this order. In that way we can quantify the level of differences between different forms of the constraints.

NLL and NNLL of N=4 sYM and ω expansion
Let us begin with recalling the definition of the BFKL eigenvalue and the Mellin transform.
The Mellin transformation into the ω space is defined as and into the γ space asF Applying it to the BFKL equation we get the algebraic form of the BFKL equatioñ where χ(γ, ω) is the kernel transformed to the double Mellin space. The ω dependence stems from the z dependence of the kernel K due to the kinematical constraint. Of course the fixed order LL and NLL kernels do not have the ω dependence, it is only the resummed kernel, in this case the kernel with the kinematical constraint that obtains this additional dependence. As is well known, see for example [19,20,21,46,29], such dependence in turn generates the series in α s and therefore resums the contributions from all orders in the strong coupling. The result with the constraint (12) is well known [18] where ψ is a polygamma function. When ω = 0 the above result coincides with the LL BFKL eigenvalue in QCD and in N=4 sYM We shall investigate Mellin forms of the kernels with other constraints (9) and (10) later on in Sec. 4. In order to compare the results with NLL and NNLL calculations one needs to specify the correct scale choice. The above result leads to an asymmetric kernel which is valid in the so-called asymmetric scale choice. That means it is valid when considering the DIS process, in which the x Bjorken is defined as x = Q 2 /s with Q 2 the (minus) virtuality of the photon. For the case of the different process, like for example Mueller-Navelet jets with comparable transverse momenta, the appropriate variable would be QQ 0 /s, where Q ≃ Q 0 are the scales of the order of transverse momenta of the jets. The eigenvalue in this case would be different from (16) as it should correspond to the symmetric choice of scales. Therefore one needs to perform the scale changing transformation which generates terms starting at an appropriate order. Up to NLL order this was discussed in [13,12,21]. We shall recall in detail the scale changing transformation, and what terms it generates up to NNLL in the next subsection, here we shall directly start from the symmetric counterpart of the eigenvalue (16) which has the following form The result for the NNLL eigenvalue in the case of the N = 4 sYM was derived originally in [14,15] for the symmetric case. It was later on rederived in [16] by exploiting the correspondence between the soft-gluon wide-angle radiation in jet physics and the BFKL physics. To have direct relation to these results we will focus below on the symmetric case. The poles around γ = 0 of NLL [52] and NNLL kernel in N=4 sYM are given by [14,15,16] Since it is symmetric case the coefficients of the poles around γ = 0 and γ = 1 are identical. For the purpose of simplification therefore we only focus on the expansion around γ = 0. We can retrieve the leading and the vanishing subleading poles of the N = 4 sYM case by doing ω expansion of the shifted eigenvalue (18) where the χ (i) is the i-th derivative of χ ω with respect to ω. The term lowest in order inᾱ s in the expansion is simply χ 0 which is the LL eigenvalue (17). In order to retrieve the term contributing to the NLL order we use the solution to the equation for the intercept in the lowest order of the coupling, i.e. ω 0 =ᾱ s χ 0 , and substitute it into (21) and keep terms up to first power inᾱ s . One obtains This gives the contribution to the NLL order Expanding around γ = 0 one obtains the following pole structure Similarly, we can derive the term contributing at NNLL level from the ω expansion to the second order by taking 2 ω 1 =ᾱ s χ 0 +ᾱ 2 s χ 1 , substituting again into (21) and extracting terms proportional toᾱ 2 s which results in In this case the pole structure is We observe that the leading poles −1/2γ 3 (24) and 1/2γ 5 (26) coincide with the exact result at NLL in N = 4 SYM and in QCD and NNLL order in N = 4 SYM (20). This structure is consistent with the principle of maximal transcendentality (complexity) [52] meaning that all special functions at the NNLO correction contain only sums of the terms 1/γ i (i= 3,5). What is also interesting is the absence of subleading poles, i.e. 1/γ 2 and 1/γ 4 at NLL and NNLL respectively in the ω expansion of kinematical constraint which again coincides with the exact results in N = 4 sYM (20). We shall see later in this Section that one expects this pattern is expected to hold to all orders inᾱ s from the kinematical constraint. Of course in the case of QCD at NLL there are double poles 1/γ 2 , but they originate from the non-singular parts of the QCD DGLAP anomalous dimension and the running coupling [13,12].

Scale changing transformation at NNLL
In this section we shall recall the scale changing transformation and compute it to NNLL order. The scale changing transformations were discussed in [13], see also [46,29]. Let us start with the general formula in the high energy factorization for the process with two scales Q, Q 0 Here, h A and h B are impact factors, that depend on the process in question and G ω (k T , k T 0 ) is the BFKL gluon Green's function. The unintegrated gluon density introduced in the previous section can be interpreted ad the convolution of the Green's function with one of the impact factors, i.e.F In Eq.(27) Q,Q 0 are hard scales for this process. For example in Mueller-Navelet process of production of two jets with comparable transverse momenta Q ≈ Q 0 . On the other hand in DIS Q ≫ Q 0 . The equation for G as a function of rapidity is where Y = ln(ν/k T k T 0 ), ν is the energy available for the BFKL evolution in a given process.
Here we see that the choice of the energy scale is symmetric ν/k T k T 0 , but it could also be asymmetric ν/k 2 T or ν/k 2 T 0 . This will imply that where k > = max(k T , k T 0 ) and k < = min(k T , k T 0 ). For example: This means that depending on the scale choice the kernel will undergo a scale changing transformations as well.
In the Mellin space this will imply the shift of poles in γ by ω for the function χ ω (γ). This can be expressed as [13] So for example in the case of kernel with kinematic constraint k ′2 < k 2 /z, we have in the asymmetric case and in the symmetric case, Let us now recall how the scale changing works at NLL before proceeding to NNLL. We suppose that in the following an additional dependence on ω might be present in the χ function. Going from symmetric case to an asymmetric case we have at this order Expanding in ω and keeping terms up to NLL we get The scale changing part is given by and at this order we need to take ω =ᾱ s χ 0 (γ), therefore scale changing part is see [52,12,13]. In the following, we will focus only on the pole structure and keep leading poles in γ which gives for the scale changing part at NLL This scale change (41) is in agreement with the difference between the leading poles of symmetric and asymmetric NLL kernels (36) and (35) In NNLL order, we have expand (37) to the second order in ω of χ 0 and to the first order in ω in χ 1 , which gives Now in order to find scale changing terms at NNLL accuracy, we need to keep terms in the solution for the ω at least up to NLL i.e.
Plugging this expansion of ω into (44) and keeping terms up to NNLL, we get for the scale changing transformation at this order (see also [53]) The pole structure of this scale changing transformation at γ = 0 and γ = 1 is This result is the same whenever we are using χ 1 either from exact N=4 sYM NLL kernel or the NLL kernel derived from ω expansion. We see that the scale changing does not introduce any subleading poles 1/γ 4 , 1/(1 − γ) 4 . It can be easily verified that the scale changing at NNLL (48) is consistent with the difference of the leading poles of NNLL kernels obtained from expanding (35) and (36)

Leading and subleading poles in expansion of kernel with kinematical constraint
There is an interesting feature of the subleading terms extracted from the expansion in ω of the χ function with ω shift. It is expected that a χ k has the leading pole ∼ 1 γ 2k+1 . We also observed that the subleading pole ∼ 1 γ 2k vanish at NLL and NNLL order both in N = 4 sYM results and for the kernel with kinematical constraint. This pattern can be shown to exist for kernel with kinematical constraint at higher orders and can be proved via mathematical induction. Let us assume that the kernel χ k possesses the above discussed properties, and we will prove that the higher order kernel χ k+1 also has the same properties. Let us write the solution for the BFKL equation in Mellin space up to k + 1'th power inᾱ s (this is indicated by the subscript on ω) Now let us assume that the kinematical constraint introduces some ω dependence into the kernel eigenvalue so that the BFKL eigenvalue equation can be written as and let us expand the eigenvalue in ω and keep terms to find ω k+1 with the χ (i) is the i'th derivative of χ ω on ω set at ω = 0 and of course χ 0 = χ (0) . Let us now extract theᾱ k+2 s term from (53). We note, that one needs to keep terms up to ω k on the r.h.s which will contribute to ω k+1 .
The general expansion of (55) could be really complicated. But ignoring the coefficients, we can express an arbitrary term in (55) in the following form k l=0 with a constraint on powers {j l } k l=0 j l = i.
Since a χ l term would automatically bring a power ofᾱ l s , if we demand (56) to have the term of orderᾱ k+1−i s , we have to put a new constraint on {j l }, which is k l=0 j l l = k + 1 − i .
Thus, for a term (56), given that the leading pole of χ l ∼ 1 γ 2l+1 and under the constraints (57), (58), its leading pole can be given as Also, knowing that every χ l in (56) doesn't have a subleading pole, we can see that the subleading pole of thisᾱ k+1−i s term will also vanish. Now, we can get back to expression (54). Knowing that the derivative χ (i) has the leading pole ∼ γ −(1+i) , and a vanishing subleading pole, we can conclude that the leading pole in χ k+1 is ∼ γ −(2k+3) , and its γ −(2k+2) pole vanishes.

Properties of the kernel in Mellin space
In this section we shall investigate in detail the differences between the different constraints on the level of the Mellin transform.
The BFKL equation in momentum space with the constraint (9) is be given by where the kinematical constraint is implemented onto the real emission term in the form of the Heaviside function. Performing the Mellin transformation of Eq. (60) into the ω space gives Performing another Mellin transformation, this time on k 2 T variable, we get the following algebraic form of the BFKL equatioñ where the kernel is One can further simplify this expression by performing the change of variables u = q 2 T /k 2 T . This gives also |k T + q T | 2 k 2 where φ is defined as the angle between q T and k T vectors. Then the kernel (63) becomes We note that the integration is from 0 to 1 here, and the term u ω+1−γ arises when one integrates over the region 1 < u < ∞ after the transformation of u to 1/u. The angular integral can be performed, resulting in the hypergeometric function, For comparison, the BFKL kernel with the constraint k 2 T > zq 2 T , or Eq.(10) the low z approximation, is given by [19] χ(ω, γ) = The leading pole position is of course given by the solution to This transcendental equation cannot be solved analytically with the complicated ω dependence inside the kernel, nevertheless can be solved formally to give an 'effective' kernel In order to obtain the χ eff we have solved the Eq.  The first feature of the effective kernels in all cases is that they are strongly modified in the large γ region. The results from the constraint (9) and (12) are very close in the γ ≃ 1 region. Comparing Eq. (66) and (16) we found that they are in fact identical at γ = 1 for any value of ω. One can observe that for z → 1 both (12) and (9) lead to a strong suppression of the anti-collinear phase-space, i.e. k 2 T < k ′ 2 T , whereas even for z ≃ 1 the constraint (8) leaves some phase space for integration. On the other hand the constraints (10) and (12) are very close in the small γ region while the (9) is lower in that range. This can be understood from the fact that in the collinear region, i.e. for k 2 T > k ′2 T the latter constraint is stronger, particularly if z is large.
Next, we have checked numerically the coefficients in front of the poles in γ and found agreement between different forms of constraints in the leading ∼ 1/(1 − γ) 2k+1 and subleading poles ∼ 1/(1 − γ) 2k poles when the expressions for the kernels are expanded out to NLL and NNLL .
There is one notable difference between (9) and the other two constraints, in that in the first one the collinear pole in γ is absent. This can be understood as for large values of z > 1/2 constraint (9) implies the restriction of the integration over collinear region as well.
From χ eff one can find the leading behavior as x → 0. The intercept which is found from and is plotted in Fig. 3. All the constraints give similar intercept for values ofᾱ s up to about 0.2. The constraint (10) gives somewhat larger values, whenᾱ s is further increased, while the two constraints (12) and (9) give comparable values for the intercepts even for very large values of the strong coupling. This is consistent with the shape of the kernel discussed before.

Differential form of the evolution equation with kinematical constraints
The evolution equation in the presence of the kinematical constraint becomes an integral equation, over the longitudinal variable z as is evident from Eq. (60). Equations in this form were solved in the momentum space via different techniques in the past. One technique utilized in the past in Ref. [20] was based on the extension of the method introduced in [54]. This method incorporates the interpolation in two variables x and k 2 T with orthogonal polynomials. By doing this the integral equations can be transformed into the set of linear algebraic equations and solved by simply inverting the large matrix. An alternative method was used in [29] and it was checked that it coincides with the one proposed in [20].
In this paper we shall propose and utilize another method which reduces the integral equation to the differential one. The differentiation of the BFKL equation by x is motivated by the reduction of computational complexity of the integral on the right hand side of the equation. After performing the derivative we find, that no matter what is the form of the kinematical constraint the resulting equation can be written in the differential form and has the same form containing two terms -real emission term with shifted argument x and virtual correction termunder the remaining integral over the emitted momentum q. The argument of the real emission term is shifted such, that if the kinematical constraint takes the form θ (k C (q, k) − z), then the argument of the unintegrated gluon density function in the real emission term changes in following way: Such reformulation of the BFKL equation with kinematical constraint was first proposed in [18], and similar ideas were also discussed in [37] in the context of the Balitsky-Kovchegov evolution equation in transverse coordinate space. Below we derive the differential form of the BFKL equation for each case of the kinematical constraint. Let us begin by differentiating in x the equation with kinematical constraint in the low z approximation, Eq. (10). The equation formally differentiated in x reads where we performed the change of variable z → x ′ = x/z. After performing the derivative we get: where the expression x max 1, T emerges from two terms generated by the derivative acting on the integral on the right hand side of the (74). The first one coming from derivative on the boundary of the integral generating term containing F x, |k T + q T | 2 Θ k 2 T − q 2 T and the second one acting on the Θ-function from kinematical constraint generating a term containing F x -which combine giving the result in (74).
Rewriting the derivative as a derivative by ln 1/x instead of x results in The case of constraint k ′2 T < k 2 T /z can be treated in exactly the same way with the result ∂F(x, k 2 Thus we see very intersting result in that the entire effect of the kinematical constraint is to shift the longitudinal variable in the argument of the unintegrated gluon density in the real term by a factor which includes a combination of ratio of transverse momenta. The case with the kinematical constraint q 2 T < 1−z z k 2 T is a bit different, since in this case the term originating from the derivative acting on the integral boundary is equal to 0. This is because of the Θ-function being evaluated at −∞. In other words, this constraint means that which is always stronger than x ′ > x and therefore the boundary term will not contribute here. This results in the following differential equation 6 Numerical method and results

Numerical results
In this section we shall present the numerical results for the solution of the BFKL equation with different forms of the kinematical constraint. The details for the numerical method are described in the Appendix. The equations are solved using the boundary condition with p 0 = p 1 = −0.1 and p 2 = 2 and µ = 1 GeV. Since we are not doing any phenomenology here this boundary condition has not being fitted to any data and is for illustrative purposes only. The particular choice of k T dependence is motivated by analitical solution of the LO BFKL equation and guarantes fast convergence of the solution. The term (1 − x) p 2 guaantees that the gluon is suppressed for large values of x so that the integral over the kernel can be extended to large x values. We introduce the lower cutoff onto the transverse momenta Q 0 = 0.1 GeV. The results were obtained for the fixed value of the coupling constantᾱ s = 0.2.
No k.c.  The results for the solution as a function of the k 2 T for fixed values of x and different forms of the kinematical constraint are shown in plots in Figs. 4 together with the solution to the LL BFKL. The results reflect the trend of the semi-analytical results from the section 4. In plots of k T 2 dependence for x = 10 −2 it is clear that the solution with the full kinematical constraint and the solution of the equation with the kinematical constraint in the k ′2 T < k 2 T /z approximation get close to each other for low k T 2 . On the other hand, for large values of k 2 T it is the two constraints (10) and (9) which come close to each other and (9) has a different behavior. This is related to the different behavior in the collinear limit, i.e. the fact that constraint (9) introduces modification in the collinear limit as well as discussed previously. In any case the solution with the constraint (10) always leads to the results which is largest, consistent with the analysis from the Mellin space. Similar conclusions can be reached from the analysis of the plots with the x dependence for small k T 2 shown in Figs. (5). The small z approximation of the kinematical constraint results in the solution being closer to the solution of the equation with the q T 2 ≃ k ′ T 2 approximation of the kinematical constraint for large k T 2 . From Fig.(5) we see the onset of the power behaviour in x, and in all versions of the kinematical constraint the intercept is strongly reduced with respect to the LL approximation. The differences between the solutions with different versions of the kinematical constraint are much smaller than between resummed and LL result. Nevertheless, the differences between different versions of the kinematical constraint are non-negligible and can reach up to a factor 2-3 depending on the range of x and k T . x No k.c. x ℱ (x,(10 GeV No k.c.

Conclusions
In the paper we have performed a detailed analysis of the three versions of kinematical constraints imposed onto the momentum space BFKL equation. We used the fixed coupling BFKL equation as a basic equation on which kernel's the constrains are imposed. This choice allowed us to investigate semi-analytically the properties of the kinematical constrains in the Mellin space and to compare to known results in sYM version of the BFKL equations where the results are known up to NNLL order. In particular, we observed that the leading poles in Mellin space generated by the kinematical constraint obey maximal transcendentality principle. In addition we proved that the subleading poles vanish at NLL and NNLL order both for BFKL and its sYM version, and that such feature can be expected to hold to all orders. Furthermore we observed that the χ ef f function obtained with the most general version of the constraint i.e. (9) for large values of γ agrees well with the more commonly used version of the constraint i.e (12). However, at small values of γ the χ ef f function with (12) agrees with (10) and both of them differ from the case when constraint (9) is imposed. In particular, it is important to note that the last one changes the collinear structure of the kernel i.e. it limits z to z < 1/2. While on formal grounds this effect might be neglected in the small z limit it has nonnegligible impact on the resulting gluon density. We demonstrate that by numerical solutions of the momentum space BFKL equation for unintegrated gluon density for all versions of the constraint as well as for the unmodified BFKL equation. The understanding of action of the kinematical constraint across the whole region of z will be particularly important for ongoing efforts to unify the BFKL and DGLAP evolution schemes [55,56,57] where effects of kinematical constrains may play a crucial role. Furthermore, the inclusion of the kinematical constraint is known to be relevant for BFKL phenomenology [20,58,59,60,61,62,63]. However, so far only the simplest version of it was used in the context of BFKL. It will be interesting to see how different versions of the constraint perform while a fit to F 2 data is done or what are the differences when the description of processes probing gluon density at very low x like high energy neutrinos or very forward dijet data is attempted. interpolation where k T m are values projected from Chebyshev nodes on the k T -space, whereas ζ (k T i ) is the projection of the momentum k T on the interval of [−1, 1] interval, suitable for Chebyshev interpolation. The integers N x and N k are the numbers of grid points in x and k T ranges respectively. Finally we define the function t (x, n) which is a combination of triangular functions providing the trapezoidal interpolation in ln x: t (x, n) = ln (x/x n+1 )/ ln (x n /x n+1 )Θ (x − x n+1 ) Θ (x n − x) Further rewriting (80) into the matrix form results into this form of the equation We define a matrix and finally write: To solve the equation we can use the discretized form of the differential of the function F(x, k T 2 i ) in ln x: The equation above (84) together with (83) with the boundary condition F x 0 , k T 2 i = F 0 x 0 , k T 2 i can be used to solve the equation (76) numerically using the iteration method. The algorithm can be summarized in the following diagram: Distance between the values of F x j , k T 2 i in two subsequent cycles is evaluated after each cycle and the algorithm is terminated when the distance becomes sufficiently small.
We parametrize the grid points in following ways: with x max and x min being the upper and lower limit of the x-axis of the grid, w and Q 0 being the upper and lower limit of the k T -axis of the grid and ζ i = cos π(i+1/2) N k being the Chebyshev node.