Infrared singularities of QCD scattering amplitudes in the Regge limit to all orders

Scattering amplitudes of partons in QCD contain infrared divergences which can be resummed to all orders in terms of an anomalous dimension. Independently, in the limit of high-energy forward scattering, large logarithms of the energy can be resummed using Balitsky-Fadin-Kuraev-Lipatov theory. We use the latter to analyze the infrared-singular part of amplitudes to all orders in perturbation theory and to next-to-leading-logarithm accuracy in the high-energy limit, resumming the two-Reggeon contribution. Remarkably, we find a closed form for the infrared-singular part, predicting the Regge limit of the soft anomalous dimension to any loop order.

On the theoretical front, a separate line of investigation concerns the structure of partonic scattering amplitudes in the high-energy limit [14][15][16][17][18][19][20][21][22][23][24]. Scattering amplitudes of quarks and gluons are dominated at high energies by the t-channel exchange of effective excitations dubbed Reggeized gluons. In this context the BFKL equation and its generalisations provide again a highly-valuable tool: by solving these equations iteratively one can compute high-energy logarithms order-by-order in perturbation theory [23,24].
The real part of a 2 → 2 partonic amplitude (i.e. its signature-odd part, see eq. (2.1)) is governed by an odd number of Reggeized gluons. The leading high-energy logarithms simply exponentiate, dressing the t-channel gluon propagator by a power of s/t. In Regge theory (see e.g. [25]) this behaviour corresponds to a Regge pole in the complex angular momentum plane. QCD amplitudes can thus be factorized in the high-energy limit into a t-channel Reggeized gluon exchange which captures the dependence on the energy, and energy-independent impact factors that depend on the colliding partons. However, this simple picture does not extend beyond next-to-leading logarithms (NLL) due to multiple Reggeized gluon exchange, which form Regge cuts. This was recently demonstrated explicitly in ref. [24], where these effects were computed through three-loops, by constructing an iterative solution of the relevant BFKL or Balitsky-JIMWLK equation, describing the evolution of three Reggeized gluons and their mixing with a single Reggeized gluon.
In the this paper we extend this study, focusing on the imaginary part of 2 → 2 partonic amplitudes, which are governed by the exchange of an even number of Reggeized gluons, which also form Regge cuts. The leading logarithmic corrections to the even amplitude are determined to all orders by a wavefunction of a pair of Reggeized gluons, which solves the celebrated BFKL evolution equation. This iterative solution, which will be central to the present work, can be famously described by ladder graphs, where an additional rung is generated at each order in the loop expansion.
The study of scattering amplitudes in the high-energy limit [14][15][16][17][18][19][20][21][22][23][24] is intimately linked to the study of their infrared singularity structure. Indeed, the gluon Regge trajectory α g (t) is infrared-singular, and its exponentiation along with the energy logarithms, which is a manifestation of Reggeization, is readily consistent with the exponentiation of soft singularities through the relevant renormalization group equation. The latter of course holds also away from the high-energy limit, as guaranteed by infrared factorization theorems. The correspondence between the structure of amplitudes in the high-energy limit, which is governed by rapidity evolution equations, on the one hand, and the structure of infrared singularities on the other, becomes more complicated at subleading orders. While both separately provide means to explore the structure of amplitudes to all orders in perturbation theory, the interplay between the two provides additional insight in either direction, as demonstrated multiple times over the past few years [19][20][21][22][23][24].
Infrared singularities of massless scattering amplitudes are now fully known, for general colour, kinematics and any number of partons, through three loops, owing to an explicit computation of the soft anomalous dimension at this order [26,27]. While through two loops infrared singularities are governed exclusively by a sum over colour dipoles formed by pairs of the hard-scattered partons [28][29][30][31], at three loops one encounters for the first time infrared singularities that are simultaneously sensitive to the colour and kinematics of three and four hard partons. Subsequently, ref. [24] specialised these results to the high-energy limit, and provided a detailed comparison between the singularity structure deduced from the soft anomalous dimension and what has been established there through three loops via computations in the high-energy limit. While full consistency was found, remarkably, it was shown that at three loops (see eq. (4.11) there) the real part of the amplitude is only sensitive to non-dipole corrections starting at N 3 LL accuracy, while for the imaginary part of the amplitude they appear already at NNLL accuracy.
As an application of the interplay between these limits, it was recently demonstrated [32] that the functional form of the three-loop soft anomalous dimension in general kinematics can in fact be fully recovered via a bootstrap procedure using the high-energy limit of 2 → 2 scattering, alongside other information, as input. The bootstrap programme of the soft anomalous dimension can be extended beyond three loops, provided that information from special kinematic limits is available. The imaginary part of 2 → 2 amplitudes is a natural place to start; indeed, already in ref. [23], a non-dipole contribution at four-loops and NLL accuracy could be predicted using BFKL theory.
In the present paper we continue to develop this line of investigation of the highenergy limit of 2 → 2 scattering, focusing on the imaginary (signature-even) part of the amplitude, which is governed, as mentioned above, by the exchange of a pair of Reggeized gluons satisfying the BFKL evolution equation. The leading-order equation is sufficient to determine an infinite tower of high-energy logarithms in the soft anomalous dimension 1 .
Although the BFKL Hamiltonian has been diagonalised in many instances [3], to study partonic amplitudes requires us to use the dimensionally-regulated Hamiltonian, which is comparatively less understood. We will nonetheless find an exact iterative solution! This hinges on the following reasons: the two-Reggeon wavefunction itself turns out to be finite at all orders, so that infrared divergences are controlled by the limit of the wavefunction where a Reggeized gluon becomes soft. The evolution equation then closes within that limit, dramatically simplifying its solution. This will enable us to obtain the soft limit of the two-Reggeon wavefunction to all loop orders and NLL accuracy, and corresponding closed-form expressions for the singular part of the amplitude (see eq. (3.18)) and soft anomalous dimension (see eq. (4.20) with (4.21)), which turns out to be an entire function of the coupling.
The structure of the paper is as follows. In section 2 we recall the basic notions regarding the high-energy limit of 2 → 2 amplitudes and explain how the BFKL evolution equation can be solved iteratively to determine the two Reggeized gluon wavefunction and the imaginary part of the amplitude. In section 2 we also reformulate the equation so as to explicitly display the fact that the evolution retains infared-finiteness, comment on the symmetries displayed by the evolution and recover the four-loop results of ref. [23]. Appendix A completes this review by explaining how the particular form of the evolution equation used here follows from the more general non-linear set up used in refs. [8,23,24]. In section 3 we consider the soft approximation, show that the evolution closes in this limit, and exploit this simplification to derive all-order solutions for the wavefunction and amplitude. Finally in section 4 we study the implications of our results in the high-energy limit regarding the soft anomalous dimension, obtaining a closed-form solution for the latter at next-to-leading order in high-energy logarithms to all orders, and verify the consistency of our BFKL-based result with infrared exponentiation.

Scattering amplitudes by iterated solution of the BFKL equation
The well-known BFKL evolution equation predicts the rapidity dependence of two-parton amplitudes in the high-energy limit [1,2]. In the following we briefly summarise the conclusions from this approach regarding the leading contributions to the signature-even amplitude, or the two-Reggeon cut.

2.1
The even amplitude from the BFKL wavefunction Figure 1. The t-channel exchange dominating the high-energy limit, s −t > 0. The figure also defines our conventions for momenta assignment and Mandelstam invariants. We shall assume that particles 2 and 3 (1 and 4) are of the same type and have the same helicity.
Let us consider a 2 → 2 scattering amplitude M ij→ij , where i, j can be a quark or a gluon. The momenta are assigned as indicated in figure 1. In the following we will suppress the species indices i, j, unless explicitly needed. The high-energy limit corresponds to a configuration of forward scattering, such that the Mandelstam variables satisfy s −t > 0. In analysing this limit it is convenient to decompose the amplitude into its odd and even components with respect to s ↔ u exchange, the so-called signature: where M (+) , M (−) are referred to, respectively, as the even and odd amplitudes. As shown in ref. [24], these have respectively real and imaginary coefficients, when expressed in terms of the natural signature-even combination of logarithms, and have independent factorisation properties in the high-energy limit. The effect we discuss in the following originates from the exchange of two Reggeons, therefore it proves useful 2 to 2 The full advantage of considering the reduced amplitude will become clear in what follows. First, BFKL evolution of the reduced amplitude involves an extra term proportional to T 2 t in (2.17). This term renders the wavefunction finite. Second, upon performing infrared factorization of the reduced amplitude one is able to identify the NLL terms that originate in the soft anomalous dimension -see eq. (4.10).
define a reduced amplitude, as introduced in ref. [24], dividing by the effect of one-Reggeon exchange:M where T 2 t represents the total colour charge exchanged in the t channel (see eq. (2.10) below). The function α g (t) in eq. (2.3) represents the gluon Regge trajectory having the perturbative expansion Given that we work to NLL accuracy, we will only need the gluon Regge trajectory to first order in α s , where in d = 4 − 2 dimensions Here, B 0 ( ) is a ubiquitous loop factor and the first of a class of bubble integrals, cf. eq. (3.6), to become important in section 3. For now, it suffices to know that In the following we will consider the leading contributions to the signature-even amplitude to all orders, corresponding to the two-Reggeon exchange. These corrections -which we denote byM (+) NLL -were studied long ago [1,2] and can be expressed in terms the two-Reggeized-gluon wavefunction Ω(p, k) as follows: where the integration measure is (2.8) with B 0 ≡ B 0 ( ) and the tree amplitude is where λ i for i = 1 through 4 are helicity indices. The colour operator T 2 s−u in eq. (2.7) acts on M (tree) ij→ij and it is defined in terms of the usual basis of Casimirs corresponding to colour flow through the three channels [22,33]: where T i represent the colour charge operator [34] in the representation corresponding to parton i. The wavefunction Ω(p, k) has a perturbative expansion in the strong coupling, taking the form where we set the renormalization scale equal to the momentum transfer, µ 2 = −t = p 2 . The amplitude itself then has the corresponding expansion We emphasise that while these corrections are the leading-logarithmic contributions to the even amplitude, we denote them by NLL to recall that the power of the logarithm L is one less than the loop order. This can be contrasted with the single-Reggeized-gluon contribution to the odd amplitude M In the normalisation used in eq. (2.13), the leading-order wavefunction is simply (2.14) At loop level the wavefunction is then obtained iteratively by applying the BFKL Hamiltonian: 16) and the functionJ(p, k) is Eq. (2.15) is the standard BFKL Hamiltonian (see eq. (17) of the initial reference [1]) written using dimensional regularisation as an infrared regulator.J(p, k) accounts for the Regge trajectories of the individual Reggeized gluons, minus the overall Regge trajectory with colour charge T 2 t which was subtracted in the exponent of the reduced amplitude (2.3). As discussed in refs. [23,24] and briefly reviewed in appendix A, this equation and its higher-order generalisations can be understood by considering the expectation value of Wilson lines associated to the colour flow of the external partons [8], which are described as "target" and "projectile" in the (high-energy) forward scattering configuration of figure 1. The wavefunction then represents the transverse momenta in each of two Wilson lines and the BFKL equation is obtained as an appropriate limit of the more general Balitsky-JIMWLK evolution equation.
A graphical representation of eq. (2.13) is provided in figure 2. As a result of BFKL evolution, the amplitude at NLL accuracy can be represented as a ladder. At order it is obtained by closing the ladder and integrating the wavefunction of order ( − 1) over the resulting loop momentum, according to eq. (2.13). The wavefunction Ω ( −1) (p, k), in turn, is obtained by applying once the leading-order BFKL evolution kernel to the wavefunction of order ( − 2). Graphically, this operation corresponds to adding one rung to the ladder. Figure 2. Graphical representation of the amplitude at NLL accuracy, as obtained through BFKL evolution. The addition of one rung corresponds to applying once the leading-order BFKL evolution onto the projectile wavefunction or impact factor at order ( − 2). This gives the wavefunction at order ( − 1), according to eq. (2.18). Closing the ladder and integrating over the resulting loop momentum gives the reduced amplitude, according to eq. (2.13).

Iterative solution for the wavefunction and amplitude
Eq. (2.13) shows that the -th order amplitude is obtained in terms of iterated integrals, which arise upon evaluating the wavefunction Ω ( −1) (p, k) to order ( − 1). It is straightforward to compute the first few orders, which gives us an opportunity to revisit the findings of ref. [23]. We will be able to explain why a new colour structure emerges for the first time at four loops, and explore the general structure of the relevant iterated integrals. A useful fact is that the evolution admits one well-known solution in the case where the exchanged state is colour-adjoint and Ω(p, k) is constant (independent of k) [1,2], which gives a positive-signature state with the same leading-order trajectory as the Reggeized gluon. This enables one to rewrite the Hamiltonian (2.15) as a part which vanishes when Ω(p, k) is constant, plus a part proportional to (C A − T 2 t ): where the function J(p, k) is defined by (2.20) The first interesting feature to note is that theĤ i operator in eq. (2.18) vanishes when acting on Ω (0) (p, k) = 1. Therefore the wavefunction to one-loop involves a single colour structure: The second colour structure appears for the first time at the second order: (2.22) Inserting the explicit form of J(p, k) from eq. (2.20) into eq. (2.22), one finds that it involves bubble integrals, as well as three-mass triangle integrals with massless propagators, such as which is represented in figure 3. The wavefunction at higher orders can be expressed p k p − k formally by introducing a class of functions , where Ω ∅ (p, k) ≡ 1, and each of the indices a j can take the value "i" or "m", which stand for integration and multiplication, respectively, according to the action of the two Hamiltonian operators in eq. (2.19). In this notation, the one-and two-loop wavefunctions read, respectively, and it is also easy to write the wavefunctions at higher loops, for example: The wavefunctions written thus far are sufficient to evaluate the reduced amplitude up to four loops. At one and two loops, inserting respectively Ω (0) (p, k) = 1 and eq. (2.21) into eq. (2.13) and performing bubble integrals one gets immediatelŷ We notice that the amplitude depends solely on the colour structure (C A − T 2 t ), and this in turn is a consequence of the fact that the wavefunctions Ω (0) and Ω (1) have only one colour component. Based on this consideration alone, one would expect the second colour structure, (2C A − T 2 t ), to contribute to the amplitude starting at three loops, given that it appears in Ω (2) (p, k) of eq. (2.25). However, this contribution of Ω (2) (p, k) to the amplitudê M (+,3) NLL cancels by symmetry: where in the last line we used the property which makes evident that eq. (2.28) vanishes by antisymmetry with respect to k ↔ k . Because of this, the amplitude at three loops has again a single colour component, proportional Figure 4. Graphical representation of the BFKL ladder at four loops. The fact that in conjunction with the target-projectile symmetry imply that the first rungs on either side can only give rise to contributions proportional to (C A −T 2 t ). As a consequence, distinct colour structures can appear for the first time at four loops.
This symmetry relation generalises to higher orders, i.e. one has for any a 1 . . . a n . While this symmetry ensures that there is only one colour structure at three loops, this is no longer the case starting at four loops. There, one obtains [23] M (+,4) One sees that the integrated result involves two colour structures, and in the final expression in eq. (2.32) we rearranged them so as to single out a factor of C A . In section 4 below we will see that in this form it is easy to compare the amplitude with the structure of infrared divergences. Specifically, we will see that corrections involving the colour structure (C A − T 2 t ) −1 at loop order emerge directly from the simplest "dipole" formula of the soft anomalous dimension, while other colour structures, namely C j A (C A − T 2 t ) −j−1 with j ≥ 1, identify deviations from the dipole formula, as was first observed in ref. [23] for = 4.
Inspecting the diagrammatic representation of BFKL evolution in figure 2, one can interpret the delayed appearance of a new colour structure to four loops, as a consequence of the target-projectile symmetry. Recall that for the first rung of the ladder, only the second termĤ m in eq. (2.18) contributes, so the wavefunction has a single colour structure (C A − T 2 t ). Considering more rungs, using target-projectile symmetry one can deduce that the same is true for the first rung on the opposite side of the ladder. As a consequence, despite the fact that Ω (2) (p, k) contains two structures (see eq. (2.25)), the effect of the second one, (2C A − T 2 t ), cannot appear in the three-loop amplitude, where each of the two rungs contribute a factor of (C A − T 2 t ). As shown in figure 4, distinct colour structures can only appear in the amplitude starting at four loops, where the middle rung -and only that rung -gives rise to both colour factors.

The soft approximation
While it would be possible to calculate the wavefunction and amplitude to higher loop orders, in this paper we focus on the infrared divergent part of the latter. We strive to compare its singularities with the predictions made by the infrared factorisation theorem and, consequently, deduce higher-order corrections to the high-energy soft anomalous dimension. With this goal in mind, we highlight at this point another important property of Ω(p, k), which can be verified when inspecting eq. (2.18) more carefully (see below): the wavefunction Ω ( −1) (p, k) is finite for → 0 to all orders in perturbation theory! This is a non-trivial statement, which becomes evident only after the evolution evolution equation is brought from the form in eq. (2.15) to eq. (2.18). A practical implication is that all divergences in the amplitude must originate in the final integration, namely going from the wavefunction to the amplitude as in eq. (2.7). Inspecting the latter equation, we see that divergences arise only in the k → 0 and k → p limits (and ultraviolet power counting in eq. (2.19) using (2.16) excludes divergences from k p, k). Due to the symmetry of the integrand, all divergences of the amplitude can therefore be obtained by evaluating it in one of these two limits, and multiplying the result by two.
Let us now examine more carefully the evolution of the wavefunction according to eqs. (2.18) and (2.19), verify that the wavefunction is indeed finite, and derive a simplified version of the evolution, valid in the small-k or soft approximation: k p. The loop integral in eq. (2.19) can in principle receive contributions from two regions; k k ∼ p and k ∼ k p. Inspecting the form of f (p, k, k ) in the two regions, it is easy to check that only the second region contributes: This means that the soft approximation closes under evolution! In the following, we will identify the region k ∼ k p as soft and add a subscript s to quantities calculated in this limit. From J(p, k) in eq. (2.20) one gets 2) and the evolution in eq. (2.18) becomes where [Dk ] is the previously defined integration measure (2.8). Eq. (3.3) confirms that it is consistent to truncate the Regge evolution to the soft approximation: using the power counting Ψ(p, k) ∼ 1, we see that the k integral is saturated by the soft region k ∼ k, with no sensitivity to larger scales. Inserting the wavefunction Ω ( −1) s (p, k) into eq. (2.13), we get the amplitude in the soft limit at the -th order. In this approximation the last integral becomes divergent and needs an ultraviolet cutoff, which we fix by requiring k 2 < p 2 , based on dimensional analysis and consistency with the soft limit (any cutoff would be consistent, and would not affect the infrared singularities). The integration measure for the last integral therefore reads where we multiplied by a factor of two, in order to take into account the fact that there is an identical contribution from the region where the Reggeized gluon carrying momentum (p − k) is soft. Inserting this result into eq. (2.13), we getM NLL in the soft approximation: We stress that this approximation gives correct results only as far as infrared singularities are concerned. All poles in are exact, since the integrand is finite and divergences arise only from the k → 0 limit of integration. The reduced amplitude in eq. (3.5) ceases to be correct at finite O( 0 ) order, as indicated. The most significant advantage of the soft approximation is that the evolution equation greatly simplifies, and this allows us to obtain closed-form expressions for the wavefunction Ω ( −1) s (p, k) and the amplitudeM (+, ) NLL | s , as we are going to detail in the following.

The wavefunction at NLL to all orders
In analogy to the exercise done in section 2.2, we start by calculating explicitly the wavefunction at the first few orders in perturbation theory, this time in the soft approximation. The initial condition is still given by eq. (2.14), and the evolution obeys eq. (3.3). This equation has a much simpler structure compared to the original one, eq. (2.19), because the soft approximation turns a two-scale problem into a one-scale problem. It is easy to check that the wavefunction reduces to a polynomial in ξ = p 2 /k 2 , which implies that the integrals involved in eq. (3.3) are simple bubble integrals of the type where the integration measure is given in eq. (2.8). This defines a class of one-loop functions mentioned above eq. (2.6), namely Using this we can write the action of the soft Hamiltonian (3.3) on any monomial (m ≥ 0): where we have introduced the loop functionŝ The second line will be useful in what follows for determining the all-order structure of the wavefunction. Applying eq. (3.6) repeatedly up to three loops (which is sufficient to determine the amplitude at four loops) we find The evaluation of a few additional orders allows us to obtain an ansatz for the ( − 1)-th order wavefunction: (3.11) The validity of this all-order formula can be proved directly using the action of the Hamiltonian in the second line of eq. (3.8) by noticing first that, independently of the loop order, the term ξ n can only be generated by acting n times with the second term of eq. (3.8), each of which raises the power of ξ by one. Hence ξ n will always be accompanied by the product (−1) n n−1 . Furthermore, the combinatorial factor −1 n associated with ξ n simply counts the number of different ways of acting ( − 1) times with the Hamiltonian, out of which n times with the second term and − 1 − n times with the first.

The all-order structure of two-parton scattering amplitudes at NLL
The main result of the previous section is that, in the soft approximation, the wavefunction reduces to a polynomial in p 2 /k 2 , given by eq. (3.11). As a consequence, the calculation of the amplitude (3.5) becomes straightforward, because it involves only integrals of the type which allows us to obtain where the factor (1 −B −1 ) follows from rewriting the factor e γ E /Γ(1 − ) = B −1 ( ): Eq. (3.13) looks rather involved but one must keep in mind that, upon expansion in , it contains many finite terms which do not represent the actual amplitude since we are working in the soft approximation. Given the overall factor of 1/(2 ) in eq. (3.13), all the singularities are obtained by retaining only contributions up to −1 in the subsequent factors. When this is taken into account a great simplification arises: indeed, as shown in appendix B, it is possible to prove that eq. (3.13) is equivalent tô It is remarkable that the complicated sum of products of bubble integrals weighed by a binomial factor collapses to a single factor which depends only on one bubble integral, namelyB −1 ( ). The main ingredient of the proof is the fact that the wavefunction itself is finite. Eq. (3.15) constitutes the main result of this section: by iterating the BFKL equation (which was not diagonalised before in d = 4 − 2 dimensions) we obtained the singular part of the even amplitude at NLL accuracy, to all orders in the strong coupling constant. Anticipating comparison with the structure of infrared divergences dictated by the soft anomalous dimension, it proves useful to rearrange eq. (3.15) in such a way to single out the colour structures C A and (C A −T 2 t ). Indeed, as discussed at the end of section 2.2, we know that the dipole formula of infrared divergencies fixes the singularities of the even amplitude in the high-energy limit to be proportional to the colour structure (C A − T 2 t ) −1 T 2 s−u at loops. From eq. (3.15) we obtain Furthermore, by resumming eq. (3.16) according to eq. (2.12) we get the all-order amplitude: This result will be used in the next section to extract the soft anomalous dimension. Before addressing this topic, however, it proves useful to explore in more detail the implications of eq. (3.16) by writing explicitly a few orders in perturbation theory. Up to three loops eq. (3.16) reduces tô i.e. only one colour structure contributes to the amplitude up to three loops, and the singularities are correctly reproduced by the dipole formula of infrared divergences. Starting at four loops, and for the subsequent three orders, one gets an additional contribution proportional to a new colour structure: 20) which matches with the infrared-divergent part of the result reported earlier in eq. (2.32). It can be easily verified (see the next section) that the infrared divergences associated with the first colour structure are predicted by the dipole formula, while the ones associated with the second are not. Next, starting at seven loops, and for the subsequent three orders, yet another colour structure arises: Expanding eq. (3.16) for the next three orders in α s we get It is now easy to understand the pattern singularities implied by eq. (3.16): at each order the first colour structure, proportional to (C A −T 2 t ) −1 , describes the singularities predicted by the dipole formula. Additional colour structures are generated by the expansion of the geometric series 16), such that every three loops a new colour structure arises with an increasing power of C A , replacing one of the factors of (C A − T 2 t ). All these new structures introduce infrared divergences, which are not accounted for by the dipole formula. Now that we understand the result implied by the BFKL evolution equation, we are in the position to investigate how the infrared divergences not accounted for by the dipole formula can be included in the soft anomalous dimension. This will be the subject of the following section.

The soft anomalous dimension in the high-energy limit to all orders
It is well known that infrared divergences in gauge-theory scattering amplitudes are multiplicatively "renormalizable": finite hard-scattering amplitudes may be obtained by multiplying the original infrared-divergent amplitude by a renormalization factor Z({p i }, µ, α s (µ)), which is matrix-valued in colour-flow space. This factor solves a renormalization group equation, and hence can be written as a path-ordered exponential of a soft anomalous dimension Γ({p i }, µ, α s (µ)), integrated over the scale µ. As such, the soft anomalous dimension constitutes a fundamental ingredient for the calculation of scattering processes at any given order in perturbation theory, and much effort has been devoted to its determination. It has been shown that the soft anomalous dimension has a simple dipole structure up to two loops [28]. Corrections involving three and four partons arise starting at three loops, and a series of analyses has been performed in order to constrain their structure at three loops and beyond [29][30][31][35][36][37]; the complete correction at three loops was calculated recently [26,27].
The general structure of the soft anomalous dimension is fixed by the factorisation properties of soft and collinear radiation, along with symmetry properties, such as rescaling invariance of soft corrections with respect to the momenta of the hard partons. The latter properties link dipole terms to the cusp anomalous dimension and dictate the structure of corrections to the soft anomalous dimension that correlate more than two partons [29][30][31]35]. In particular, they imply that at three loops, non-dipole corrections can only depend on the kinematics via rescaling-invariant cross ratios. The soft anomalous dimension can be further constrained by the behaviour of scattering amplitudes in special kinematic limits, such as the Regge limit [21,22,24] and collinear limits [30,36]. Furthermore, it was recently shown [32] that the space of functions in terms of which the non-dipole correction is expressed (single-valued multiple polylogarithms) can, in fact, be deduced from general considerations. A bootstrap procedure was then set up, which remarkably completely fixes the functional form of the non-dipole correction at three loops (up to an overall rational numerical factor) based on known information from the kinematic limits mentioned above, reproducing the result of the Feynman-diagram computation of ref. [26,27]. The prospects of extending this bootstrap procedure to higher loops provides an additional motivation to determining the soft anomalous dimension in the high-energy limit.
As discussed above, ref. [23] determined the next-to-leading high-energy logarithms (NLL) of 2 → 2 scattering amplitudes at four loops. In this paper we have been able to extend this and computed the infrared singularities at NLL in the high-energy limit to all order in perturbation theory. We are therefore able to determine the soft anomalous dimension in this approximation to all orders.
We start this section by briefly reviewing the structure of the soft anomalous dimension in the high-energy limit, and then determine it to all orders by extracting the O(1/ ) coefficient from the amplitude obtained in section 3.2, which we then analyze numerically in detail. Finally we show that the singularity structure we deduced from the high-energy limit computation, consisting of poles of O(1/ ) through to O(1/ ) at loops, is consistent with infrared factorisation, namely it is exactly reproduced by the expansion of the pathordered exponential of the integral of the soft anomalous dimension.

The infrared factorisation formula in the Regge limit
The infrared divergences of scattering amplitudes can be factorised as where H is a finite hard-scattering amplitude while Z captures all singularities. Z admits a renormalization group equation whose solution (in the minimal-subtraction scheme) can be written as a path-ordered exponential of the soft anomalous dimension: The scale dependence of the soft anomalous dimension Γ ({p i }, λ, α s ) for massless-parton (p 2 i = 0) scattering is both explicit and via the 4 − 2 dimensional coupling. In QCD (with n f light quark flavours) the latter obeys the renormalization group equation For our purposes only the zeroth order solution will be needed: α s (µ) = α s (p) p 2 /µ 2 . The explicit dependence on the scale (Γ is linear in log λ) reflects the presence of double poles due to overlapping soft and collinear divergences. The soft anomalous dimension in multileg scattering of massless partons is an operator in colour space given by [26,[29][30][31]35] where Γ dip. involves only pairwise interactions amongst the hard partons, and is therefore referred to as the "dipole formula". The kinematic variables are −s ij = 2|p i · p j |e −iπλ ij with λ ij = 1 if partons i and j both belong to either the initial or the final state and λ ij = 0 otherwise. The function γ K (α s ) in eq. (4.4) is the (lightlike) cusp anomalous dimension [38][39][40], divided by the quadratic Casimir of the corresponding Wilson lines. The functions γ i (α s ) represent the field anomalous dimension corresponding to the parton i, which governs hard collinear singularities. Both γ K (α s ) and γ i (α s ) are known through three-loop in QCD and their values are summarised in Appendix A of ref. [24]. In eq. (4.4) ∆ (n) for n ≥ 3 accounts for multi-parton correlations. The three-loop correction ∆ (3) , correlating up to four hard partons, was calculated recently [26,27] for any number of partons in general kinematics. Specializing to 2 → 2 parton scattering in the high-energy limit, ref. [24] showed that ∆ (3) contributes starting from NNLL accuracy in the imaginary (even) part of the amplitude, and starting from N 3 LL accuracy in the real (odd) part; we refer the interested reader to eq. (4.11) in ref. [24] for an expression for ∆ (3) in this limit. Given our focus here on NLL accuracy, we shall not discuss it further. While it is known that Γ dip. fully describes the infrared singularities associated with Regge pole factorisation [21,22] -meaning it is exact at leading and NLL accuracy for the real part the amplitude -it does not fully capture the structure of the two-Reggeon cut [23] at NLL accuracy, where ∆ (n) at four loops and beyond, are relevant. To identify the contribution of the soft anomalous dimension in two-parton scattering, ij → ij, at increasing logarithmic accuracy, let us expand Γ in powers of α s , keeping the product α s L fixed, as follows: Γ (α s (λ)) = Γ LL (α s (λ), L) + Γ NLL (α s (λ), L) + Γ NNLL (α s (λ), L) + . . . . (4.5) The N k LL term in eq. (4.5) can be written as an expansion in α m s L m−k for m ≥ 1. Using Regge-pole factorisation it can be shown [21,22] that the leading logarithmic contribution Γ LL takes the one-loop exact form, This exactly corresponds to the infrared-divergent part of the one-loop gluon Regge trajectory in eq. (2.5). Note that the LL anomalous dimension has even signature Γ LL = Γ (+) LL . At NLL the anomalous dimension can be divided into signature-even and odd parts: (4.7) The even part 3 , which is governed by the Regge pole, is two-loop exact. Referring to eq. (4.4), it contains the terms in the one-loop anomalous dimension that are not enhanced by L, as well as the infrared-divergent part of the two-loop gluon Regge trajectory: The odd part is however sensitive the the two-Reggeon cut. At one-loop it can be obtained from the dipole formula [21,22], while higher-order terms have so far been unknown. The reduced amplitude obtained in section 3 contains information on the infrared divergences of next-to-leading high-energy logarithms to all orders in α s , and hence allows us to determine Γ (−) NLL to all orders. In order to make contact with section 3 we need to express the reduced amplitude defined in eq. (2.3) in its infrared-factorised form. Focusing on the even component, we substitute eq. (4.1) there and expand it to NLL accuracy: where we have written the Regge trajectory explicitly according to eq. (2.5). Substituting Γ LL of eq. (4.6) into eq. (4.2) and integrating over the scale (using the zeroth-order scale dependence of α s ) we obtain: Considering the second term in the square brackets of eq. (4.10) we note that Z NLL . This implies that the infrared-singular part of the reduced amplitude is insensitive to H (+) NLL [23] and is given by: Equation (4.12) can be further simplified by noticing that the hard function at LL accuracy is fixed by Regge factorisation: it is simply the exponential of the finite part of the gluon Regge trajectory, i.e. we have where we used the fact that T 2 t = C A when acting on the Regge limit of the tree level amplitude. Moving this (finite) exponential to the left, this result allows us to write eq. (4.12) more explicitly as where it is understood that both sides of this equality are to be projected onto even signature. Below we will abbreviate the l.h.s. asM NLL . The NLL contribution to the pathordered exponential on the second line can be written out fully as (4.15) Finally, integrating the exponents in each of the two brackets as in eq. (4.11) and using again that T 2 t = C A in the right factor upon acting on M (tree) , we obtain, projecting onto the even amplitude: This expression for the even amplitude may be compared directly with the one obtained in eq. (3.18) using the BFKL analysis; exploiting the fact that the exponential on the l.h.s. of eq. (4.14) is finite (and that R( ) is finite), the BFKL prediction can be written as with R( ) defined in eq. (3.17). We now have two expressions for the infrared singularities of the reduced amplitude -an expression in terms of the soft anomalous dimension, eq. (4. 16), and the all-order result of BFKL evolution in the soft approximation, eq. (4.17). In the next section we equate them and extract Γ (−) NLL .

Extracting the soft anomalous dimension at NLL
In minimal subtraction schemes, anomalous dimensions can be extracted by taking the coefficient of pure 1/ single poles. Indeed, to get the coefficient of the single poles in eq. (4.16) we can drop the exponentials to get This result must be set equal to the single poles obtained from eq. (4.17), whose -loop coefficient is Comparing with eq. (4.18) then gives where the subscript indicates that one should extract the coefficient of −1 . Although the notation does not manifest this, the end result is always a polynomial in colour operators C A and T 2 t , since R( ) has a regular series as → 0. Rescaling , this can also be written as where the function R( ) = −2ζ 3 3 + . . . is defined in eq. (3.17).
Equation (4.22) is the main result of this paper: it gives the soft anomalous dimension in the Regge limit to any loop order at next-to-leading logarithmic accuracy (i.e. all terms of the form α s L −1 ); the even contribution Γ (+, ) NLL was given in eqs. (4.6) and (4.8). In other words, we now know eq. (4.9) to all orders: Expanding the above formula explicitly to eight loops: (4.24) These results are valid in any gauge theory, and hold modulo colour operators which vanish when acting on the Regge limit of the tree amplitude (which is given by the t-channel gluon exchange diagram).

Properties of the soft anomalous dimension in the Regge limit
In the previous section we computed Γ (−) NLL , the imaginary part of the soft anomalous dimension in the Regge limit, to all orders. Let us briefly explore its properties addressing the colour structure, the convergence of the expansion, and finally its asymptotic high-energy behaviour.
Considering eq. (4.24), our first observation is that colour structures of increasing complexity emerge every three loops, as dictated by the expansion of R( ) in eq. (3.17): corrections going beyond the dipole formula start at four loops, where the colour structure is proportional to C A to a single power. This correction reproduces precisely that found previously in ref. [23]. Proceeding to five and six loops Γ NLL only incurs extra powers of (C A − T 2 t ). Starting at seven loops, however terms with two powers of C A appear as well. Similarly, a cubic power of C A would emerge at ten loops, and so on. We also note that the zeta values appearing in Γ NLL are of uniform weight, which is, of course, again a mere consequence of the Taylor series of R( ).
To proceed it would be useful to specify the relevant colour charge exchanged in the t channel, T 2 t . To this end consider for example gluon-gluon scattering, where the t channel colour flow can be any of the SU(N c ) representations appearing in the decomposition 4 where the labels refer to their dimensions for N c = 3. Because of Bose symmetry, the symmetry of the colour structure mirrors the signature of the corresponding amplitudes under s ↔ u exchange. Thus, only even representations are relevant for the two-Reggeon amplitude discussed here; these are the singlet, where T 2 t = 0, the symmetric octet with T 2 t = C A = N c , the 27 representation with T 2 t = 2(N c + 1), and the "0" representation, where T 2 t = 2(N c − 1). In the following we restrict the discussion to the first three cases, which are all relevant for QCD with N c = 3 (the latter has a vanishing dimension, and hence it does not contribute).
The next observation, already mentioned in section 2.2, is that the symmetric octet representation with T 2 t = C A , corresponds to a constant wavefunction, and thus a trivial solution to eq. (2.18), with no corrections to the reduced amplitude beyond one loop (as can be verified for example in the explicit results in eqs. (3.19) through (3.22) upon considering T 2 t = C A ). The reduced amplitude for the symmetric octet state is thus one-loop exact, corresponding to a simple Regge-pole behaviour with a gluon Regge trajectory for the original amplitude according to eq. (2.3). This of course reproduces the known behaviour of the symmetric-octet exchange used in the original derivation of the BFKL equation. In turn, for the singlet -the famous Pomeron -and 27 representation, we find non-trivial radiative corrections associated with a Regge cut. We will thus use these two examples in the discussion that follows.
Next let us consider the convergence properties of the perturbative series representing the soft anomalous dimension in eq. (4.20). One immediately notes that this series is highly convergent due to the 1/( − 1)! prefactor in eq. (4.21). Figure 5 illustrates this factorial suppression of the coefficients G ( ) as a function of the order for C A = N c = 3 and for the two relevant representations, the singlet and the 27. Furthermore, we can establish that the anomalous dimension (4.22) has an infinite radius of convergence as a function of x ≡ Lα s /π. To see this we write the resummed soft anomalous dimension as: Γ where the generating function for the expansion coefficients is defined by It is convenient to further identify G(x) as the Borel transform of some function which upon using eq. (4.21), simply evaluates to We may now recover the original G(x) via the integral where the integration contour runs parallel to the imaginary axis, to the right of all singularities of the integrand.
The function g(y) in eq. (4.28) only has isolated poles away from the origin and has a finite radius of convergence: it is well-defined in a disc around the origin. It then follows that G(x) has an infinite radius of convergence, hence this function -and the soft anomalous dimension Γ (−) NLL in eq. (4.26) -is an entire function, free of any singularities for any finite x = Lα s /π.
We stress that our use of the Borel transform is opposite to the usual application of Borel summation (which is ordinarily used to sum asymptotic series): the function G(x), in which we are interested, is an entire function; we make use of its inverse Borel transform, g(y), which has worse behaviour by having merely a finite radius of convergence. Nonetheless we find that numerically integrating eq. (4.30) is a particularly convenient way to evaluate the anomalous dimension. This numerical integration is compared to the partial sums (4.31) in figure 6, where we find good agreement for the given values of x. While it becomes challenging to efficiently compute the coefficients G ( ) at high orders (here we only evaluated them for ≤ 22), we find the numerical integration of eq. (4.30) to be very stable, even for larger values of x. Thus, the remarkable convergence properties of G(x) along with the Borel technique, presents us with the possibility of computing Γ NLL for x = Lα s /π 1, i.e. at asymptotically high energies. This is a rather unique situation in a perturbative settingin other circumstances resummation techniques are limited to the region x = Lα s /π 1. Evaluating the integral (4.30) and plotting G(x) for larger values of x reveals oscillations with a constant period and an exponentially growing amplitude. Since this behaviour is difficult to capture graphically we instead show the logarithm of |G(x)| weighted by the sign of G(x) in figure 7. This observation suggests to approximate (4.30) by G(x) → c e ax cos (bx + d) , (4.32) for sufficiently large values of x. By means of eq. (4.28), this model is equivalent to which is to be integrated as in (4.30) with a contour to the right of the poles. We thus find that to capture the behaviour G(x) at large x it is sufficient to simply consider g 1 η as a pair of complex-conjugated poles at η = a ± ib. Indeed, numerically extracting the rightmost poles of g 1 η of eq. (4.29) to identify the parameters a and b in eq. (4.33), and dividing the full, numerically-evaluated, G(x) by e ax leaves us with almost pure cosine-like behaviour for any x 1, as can be seen in figure 8. For reference, we quote our numerical results for a, b, c and d in

Exponentiation check for higher-order infrared poles
As a final step we confirm the agreement between the BFKL prediction and the soft factorisation theorem. Thus far we have only used the single poles as predicted by the BFKL evolution to extract the NLL soft anomalous dimension Γ (−) NLL . As explained in section 4.1, higher-order poles of the amplitude are generated upon expansion of the path-ordered exponential in eq. (4.16). They have to match the BFKL computation and therefore provide an independent and non-trivial check of our results.
To see how this works, let us expand the BFKL result (4.17) to the first few orders, namelyM This shows that infrared exponentiation works out if, and only if, all the poles in the NLL amplitude can be written as a function of only (i.e. independent of α s ), times the quantity in the square bracket. With hindsight, infrared exponentiation thus explains the compact form of the BFKL result in eq. (4.17). Finally, it is straightforward to substitute in the definition of G (k) from eq. (4.21) and sum up the series over k, recovering the full result for the singularities of the amplitudes in eq. (4.17). This completes the proof that the BFKL result we obtained is consistent with infrared factorisation.

Conclusions
We considered the even signature component of two-to-two parton scattering amplitudes in the high-energy limit. This amplitude is dominated by the t-channel exchange of a state consisting of two Reggeized gluons, corresponding to the simplest example of a Regge cut in QCD. The amplitude can be evaluated in QCD perturbation theory by iteratively solving the BFKL equation. Each order in perturbation theory corresponds to one additional rung in the BFKL ladder, building up a tower of so-called next-to-leading logarithms, O(α s L −1 ). Although the BFKL Hamiltonian has been diagonalised in many cases [3], the dimensionally-regulated Hamiltonian relevant for partonic amplitudes has remained more difficult to handle. Our first observation was that the wavefunction describing the two Reggeized gluons remains finite through BFKL evolution for any number of rungs, while the corresponding amplitude develops infrared singularities due to the soft limit of the wavefunction. We further observed that the evolution of a state in which one of the two Reggeized gluons is much softer than the other, k p − k, yields again a similar state. In other words, the soft approximation is consistent with BFKL evolution, and as a consequence, one can systematically solve the equation to any loop order within this approximation. We found that the soft approximation leads to a major simplification, where all integrals reduce to products of bubbles, and the wavefunction at any given order is simply a polynomial of that order in p 2 /k 2 . This eventually allowed us to determine the singularities of the amplitude in a closed form to any order, as given in eq. (3.18).
At the next step we contrasted the singularity structure we obtained though BFKL evolution with the known exponentiation properties of infrared singularities. As expected, we found that the two are consistent, and this provides a highly non-trivial check of the calculation. The leading singularity at each order, O(α s L −1 / ), is simply related to the one-loop soft anomalous dimension, and has a colour structure proportional to (C A −T 2 t ) −1 . New singularities, with fewer powers of 1/ and different colour structures, appear starting from four loop. These correspond to new terms in the imaginary part of the soft anomalous dimension, eq. (4.24). We were thus able to determine the soft anomalous dimension at next-to-leading logarithmic accuracy in the high-energy limit to all orders. These results also provide a valuable input for determining the structure of long-distance singularities for general kinematics using a bootstrap approach, as done at the three-loop order in ref. [32].
We point out that the -loop coefficient of the soft anomalous dimension we computed is a linear combination of zeta values of weight ( − 1), which coincides with the maximal (transcendental) weight. This is not surprising given that these corrections are independent of the matter content nor the amount of supersymmetry of the theory, and are thus common for example to QCD and N = 4 super Yang-Mills. We further showed that these corrections to the soft anomalous dimension can be resummed, as in eq. (4.26), into an entire function of x = Lα s /π. Remarkably, this gives us means to determine the asymptotic high-energy behaviour of this anomalous dimension, corresponding to x 1, a regime which is usually inaccessible to perturbation theory. We find that at large x the imaginary part of the anomalous dimension in the Regge limit, in any colour representation, becomes an oscillating function with an exponentially growing amplitude.
While our analysis in this paper was focused on infrared singularities, for which the soft approximation is sufficient, the formulation of the evolution in eq. (2.19) along with the observation that the wavefunction is finite, pave the way to determining the wavefunction beyond the soft approximation, thus evaluating the finite contributions to Regge-cut of two-to-two amplitudes. It would also be interesting to extend the present analysis to the next order, using the known next-to-leading order Hamiltonian; again we expect that a suitable wavefunction will remain finite to all orders, facilitating a direct determination of the infrared singularities.

A The even amplitude at NLL accuracy within the shockwave formalism
In this appendix we briefly review how eq. (2.13) can be derived within the shockwave formalism refs. [23,24]. Amplitudes in the high-energy limit are calculated as expectation values of null Wilson lines: The latter follows the path of colliding partons from the projectile or target (with x + and x − interchanged), and are labelled by transverse coordinates z ⊥ (below we shall omit the subscript ⊥ for lighter notation). The full transverse structure needs to be retained, because the high-energy limit is taken with fixed momentum transfer. Importantly, the number of Wilson lines cannot be held fixed, because the projectile and target contain an arbitrary number of virtual partons. However, in perturbation theory, the unitary matrices U (z) are close to the identity and can therefore be usefully parametrised by a field W : Physically, the colour-adjoint field W a , which is propagating in the transverse space, is interpreted as source for a BFKL Reggeized gluon [23]. At weak coupling a generic projectile is thus formed by a superposition of W states. Up to NLL accuracy one needs to consider up to two Reggeons. In this approximation, a projectile, created with four-momentum p 1 and absorbed with p 4 , is parameterised in momentum space as where the ellipses stand for wavefunction components with three or more Reggeized gluons, which are not relevant at NLL accuracy. We next note that states with an even (odd) number of Reggeized gluons have an even (odd) signature, so where D (1) i (p) is an impact factor which parameterises the dependence of the coefficient on the (transverse) momentum transfer p = p 4 −p 1 with p 2 = −t. At the leading order, there is only one Wilson line U (z) following the original parton, and the two-Reggeon wavefunction is obtained simply by expanding eq. (A.2), which gives, as in the main text: Ω (0) (p, q) = 1. (A.5) The null Wilson lines acquire energy dependence through rapidity divergences, which must be regulated, leading to the Balitsky-JIMWLK rapidity evolution equation: The scattering amplitude can be obtained by computing the overlap between ψ j | and |ψ i , after evolving them to common rapidity, where the overlap is defined as the vacuum expectation value of left-moving and right-moving W -fields. In terms of the reduced amplitude defined in eq. (2.3) one has i 2sM ij→ij = ψ j |eĤ L |ψ i ,Ĥ ≡ H − T 2 t α g (t). (A.7) Evolution at the desired accuracy is obtained by simply considering the Hamiltonian at leading order in g 2 s in terms of W fields, which, to this order, is diagonal: where "LO" and "NLO" means that all ingredients are needed respectively to leading and next-to-leading nonvanishing order. In this paper we focus on the even amplitude, representing the exchange of a pair of Reggeons, corresponding to the second term in eq. (A.9). It is then convenient to compute the inner product in eq. (A.7) by first evolving the wavefunction: 1 ! α s B 0 ( )L π d 2−2 q (2π) 2−2 Ω ( ) (p, q) W a (q)W b (p−q) .
(A.10) As displayed in eq. (2.15), the wavefunctions Ω ( ) may then be obtained iteratively by applying the HamiltonianĤ 2→2 . This Hamiltonian was discussed at length in terms of Wilson lines in ref. [24], to which we refer for further details (−H k→k is given in eq. (3.13) there; note the overall minus sign between our conventions). Acting withĤ 2→2 on the states in eq. (A.10), reproduces precisely the leading order BFKL Hamiltonian recorded in the main text. Finally, computing the overlap with the target state ψ (+) j,2 | produces the integral which closes the ladder in eq. (2.13).

B Proof of the all-order amplitude
In this appendix we show that the singular terms in eq. (3.15) are equal to those in eq. (3.13). We start by noticing that the statement is equivalent to At this point, we realise that the structure of the sum and product is strikingly similar to that appearing in the target-averaged wavefunction in eq. (3.11). In that case, finitness of the -loop wavefunction implies