Next-to-leading non-global logarithms in QCD

Non-global logarithms arise from the sensitivity of collider observables to soft radiation in limited angular regions of phase space. Their resummation to next-to-leading logarithmic (NLL) order has been a long standing problem and its solution is relevant in the context of precision all-order calculations in a wide variety of collider processes and observables. In this article, we consider observables sensitive only to soft radiation, characterised by the absence of Sudakov double logarithms, and we derive a set of integro-differential equations that describes the resummation of NLL soft corrections in the planar, large-$N_c$ limit. The resulting set of evolution equations is derived in dimensional regularisation and we additionally provide a formulation that is manifestly finite in four space-time dimensions. The latter is suitable for a numerical integration and can be generalised to treat other infrared-safe observables sensitive solely to soft wide-angle radiation. We use the developed formalism to carry out a fixed-order calculation to ${\cal O}(\alpha_s^2)$ in full colour for both the transverse energy and energy distribution in the interjet region between two cone jets in $e^+e^-$ collisions. We find that the expansion of the resummed cross section correctly reproduces the logarithmic structure of the full QCD result.


Introduction
The accurate theoretical description of non-global QCD observables [1][2][3] is among the main obstacles on the path towards precision collider phenomenology. Non-global observables are commonly characterised by kinematic constraints on limited angular regions of the radiation phase space, and occur ubiquitously at colliders, for instance via the use of jets or often when specific fiducial cuts are applied on experimental measurements. Such observables are sensitive to the coherent pattern of soft radiation outside the measured region of phase space, and their description therefore requires the calculation of such effects at all perturbative orders in the strong coupling. To achieve this resummation one meets two types of theoretical challenges. Firstly, the presence of angular cuts in the definition of the observable makes it impossible to handle the geometry of the problem analytically. This is because the distribution of the soft radiation in the angular region where the observable is defined (e.g. a jet cone or a rapidity interval) depends on the full angular pattern of the radiation in the event after the evolution from the hard scattering scale down to the low scales at which the measurement is performed. Secondly, the evolution itself is complicated because the colour structure of the squared amplitude grows drastically with each new soft emission. The most striking consequence of such a peculiar structure is that already at leading logarithmic (LL) accuracy, the standard exponentiation of soft LL corrections does not hold, and one has to solve a non-linear evolution equation, the Banfi-Marchesini-Smye (BMS) equation [1][2][3], to carry out the resummation of logarithmically enhanced corrections.
Besides the relevance of non-global resummations for collider phenomenology, a theoretical understanding of their dynamics is instrumental in the context of developing more accurate parton-shower algorithms (see e.g. Refs. [4][5][6][7][8][9][10][11][12][13] for recent work). Specifically, the resummation of next-to-leading logarithmic (NLL) non-global logarithms is a crucial ingredient for the development of NNLL algorithms, that are necessary to achieve sufficiently accurate event simulation both at present and future colliders. Moreover, their study is also motivated by purely theoretical interests, related to the discovery of a connection between the evolution of non-global dynamics and the Balitsky-Kovchegov (BK) equation [14][15][16].
The resummation of LL non-global corrections in the large-N c limit was formulated about 20 years ago in the seminal work of Refs. [1][2][3], and has seen substantial interest in recent years [17][18][19][20][21][22][23][24][25][26][27][28][29][30][31][32]. Furthermore, the authors of Refs. [33][34][35] have extended the LL resummation to include finite-N c corrections, finding that subleading-colour corrections are numerically small in common applications. Their study is however of paramount importance for the understanding of the structure of super-leading logarithmic corrections in non-global observables at hadron colliders [36,37]. The calculation of NLL corrections has inspired a considerable amount of remarkable theoretical work along the years, and formulations of the NLL resummation have been achieved in different theoretical formalisms [16,20,22]. However, a full resummation of NLL corrections for a physical observable has not yet been achieved.
In this article, we develop a formalism to resum non-global logarithms at NLL accuracy in the large-N c limit. The resummation relies on a set of non-linear evolution equations that can be solved numerically, for instance by means of Monte Carlo (MC) simulations. We apply the developed formalism to the fixed order calculation of the energy and transverse energy distribution in the region between two cone jets in the process e + e − → 2 jets, and compare our findings to an exact fixed-order prediction. The evolution equations are derived in dimensional regularisation, and later recast in a form that is manifestly finite in four dimensions, hence making it very suitable for a numerical integration. The numerical solution of the proposed equations, as well as the corresponding resummation of NLL corrections, will be presented in a forthcoming publication. The paper is structured as follows: Section 2 introduces the formalism used throughout the paper, and Section 3 presents the strategy used to derive the evolution equations. The detailed derivation of the NLL evolution equation is discussed in Section 4, while Section 5 presents the calculation of the one-loop hard matching coefficients necessary for a NLL calculation of the physical cross section for the observables considered here. Finally, in Section 6 we perform a fixedorder expansion up to O(α 2 s ) and compare our findings to the exact calculation obtained with the program Event2 [38]. Section 7 contains our concluding remarks and outlook.

Formalism and notation
Let us consider the production of two jets in e + e − annihilation at a centre-of-mass energy √ s. Considering the thrust axis as a reference axis, we define two jets in the two opposite hemispheres of the event by considering two cones of opening angle θ jet around the thrust axis. We focus on the rapidity region between the two cone jets, of total width ∆η := ln 1 + c 1 − c , c = cos θ jet . (2.1) We will informally refer to this region as the rapidity slice centred at η = 0 with respect to the thrust axis (see Fig. 1). In this paper, we study the distributions of both the energy θ jet ∆η Figure 1: The rapidity slice where the measurement is performed.
(E) and the transverse energy (E t ) in such a slice, and denote by Σ(v) the cumulative distribution for either observable to be less than v, defined as follows where σ 0 is the Born cross section for e + e − → hadrons. In the limit in which v = {E, E t } is small, large logarithms L = ln( √ s/v) spoil the convergence of fixed-order perturbative expansions and must be resummed at all perturbative orders. Since this observable is affected only by soft emissions at wide angles, the largest logarithms in Σ(v), the leading logarithms, are of the form α n s L n . For α s L ∼ 1, all terms suppressed by an extra power of α s give next-to-leading logarithmic (NLL, α n s L n−1 ) contributions, and so on. The phasespace constraint for the observables we consider admits a factorised expression in Laplace space of the type Here k ti is the transverse momentum of particle k i with respect to the thrust axis, ω i is its energy in the lab frame, and we defined u(k i ) as the source corresponding to the measurement. It takes the form where the trigger function Θ in (k) (Θ out (k)) is 1 if the particle k is inside (outside) a rapidity slice of total width ∆η, and zero otherwise. The contour γ lies parallel to the imaginary axis to the right of all singularities of the integrand. Without any emissions, at the lowest order in perturbation theory, Σ(v) = 1, and the event is made up of a quark of momentum p 1 and an antiquark of momentum p 2 , back-toback and aligned along the thrust axis. When extra radiation is considered, Σ(v) can be expressed as where the hard factors H n := H 1...n (2.6) describe configurations with n hard QCD partons along the light-like directions n 1 , . . . , n n (with n 2 i = 0 and | n| 2 = 1), while the soft factors describe the emission of soft radiation off a hard system with n hard emitters along the same directions. The convolutions in Eq. (2.5) are meant to indicate that the directions of the hard emitters in the hard and soft factors are the same, namely where Ω i indicates the solid angle of the i-th hard emitter, namely the direction of the n i vector, specified by a longitudinal (θ) and an azimuthal (φ) angle. Each of the above ingredients admits a perturbative expansion in the strong coupling constant where H (0) 2 = δ(cos θ 1 − 1)δ(cos θ 2 + 1)δ(φ 1 )δ(φ 2 ) and S (0) n = 1. At LL accuracy, the resummation for each observable requires only the first term in the r.h.s. of Eq. (2.5), with the leading order H 2 and the LL soft factor S 2 whose evolution is governed by the BMS equation [3]. At NLL one needs to include both the first (H 2 ⊗ S 2 (v)) and second (H 3 ⊗ S 3 (v)) term on the r.h.s. of Eq. (2.5), where H 2 is to be computed at next-to-leading order and H 3 at leading order. Each of the latter hard factors is individually infrared finite. The large logarithms of v are entirely contained in the soft factors S n , and all convolutions defined in Eq. (2.8) can now be carried out in 2 dimensions. At this point we can achieve NLL accuracy by including the O(α s ) corrections to H 2 and H 3 , as well as the soft factor S 3 at LL and the soft factor S 2 at NLL. In the following sections we will define each of the above contributions in detail.

Evolution equations for the soft factors in dimensional regularisation
We start by deriving the evolution equation for the soft factors S n that appear in Eq. (2.5). Given the complexity of the problem due to the fast growth of the colour space, in the following we will work in the widely used large-N c limit, where all real and virtual graphs are considered in the planar limit. Beyond this limit, the LL calculation has been performed only in a few selected cases [33,34]. The crucial advantage of using the large-N c limit is that it makes it possible to write closed evolution equations in terms of colour dipoles. This in turn makes them suitable for Monte-Carlo integration. We will first re-derive the LL evolution equation formulated in Ref. [3], and we then move on to derive the new NLL evolution equations that constitutes one of the main results of this paper.

The LL (BMS) evolution equation for S 2 and S 3
To derive the evolution equations, it is convenient to work with the Laplace transform of the soft factors of Eq. (2.5), as the observable takes a factorised form in this space. We thus define the Laplace-space soft factors G 12···n as We then start from an initial state made of the quark-antiquark pair {p 1 , p 2 }, which in a large-N c picture defines the radiating colour dipole {12}. When the emission of extra soft gluons strongly ordered in energy is considered, the initial dipole receives radiative corrections according to the squared amplitude (see e.g. [39,40]) where the sum runs over all n! permutations, and we definedᾱ = N c α s (µ)/π, where α s is the QCD coupling in the MS scheme. This approximation is valid at the leading logarithmic order, and we will consider higher-order corrections in the next section. We will derive the evolution equation using a branching formalism. This is based on identifying a resolution variable Q, such that subsequent evolution steps in Q reproduce the correct squared amplitude at a given logarithmic accuracy. Moreover, it is of course necessary that at fixed Q, the contribution of the soft evolution to the cross section be infrared and collinear (IRC) finite. The choice of the resolution scale Q can vary for different problems (it can for instance coincide with the virtuality [41] of the radiating system or with the relative angle [42] of the radiation w.r.t. the emitter). For the problem at hand, the LL corrections originate from radiation strongly ordered in energy. It is therefore natural to chose Q as the energy of the hardest gluon [3] and introduce the real-emission contribution to G 12 [Q; u] defined by where ω i is the energy of gluon k i . The evolution in Q from low scales up to Q ∼ √ s will achieve the resummation of LL corrections. We also introduce the phase-space measure in where Ω i denotes the solid angle. The inclusion of virtual corrections will be discussed shortly. The definition of Q we just adopted is only collinear safe in the strongly-ordered energy limit relevant for LL. However, our aim is to formulate a NLL evolution equation that is well defined for any value of Q and for different choices of the IRC safe observable's source u. This requires modifying the definition of Q in configurations with two gluons with commensurate energy becoming collinear, crucial for attaining NLL accuracy. We will come back to this point in Section 4.
To predict how the multi-gluon system evolves with the scale Q, we derive an evolution equation that describes the dependence of the soft system on the radiation's energy. We start by considering the large-N c real-emission contribution G (R) 12 [Q; u] defined above. We introduce the tree-level eikonal kernel where k t is the transverse momentum with respect to the emitting dipole. In the following we will work in the MS scheme, defined from the bare coupling as where β 0 is the first coefficient of the QCD β function In the large-N c limit, the factorisation properties of the squared amplitude (3.2) lead to the following evolution equation for the real contribution where G (R) 1a and G (R) a2 are defined according to Eq. (3.2) by just replacing either p 1 or p 2 with k a . In this form, Eq. (3.8) is manifestly collinear unsafe, and one needs to introduce the appropriate virtual corrections. This can be done by imposing unitarity, which enforces G ij [u] = 1 for u(k a ) = 1 and for any value of i and j (in this case i, j = 1, 2, a). We then obtain the final LL evolution equation for the physical G 12 distribution, which reads Eq. (3.10) is well defined in dimensional regularisation as the Landau singularity can be avoided by analytically continuing G 12 to complex values [43] (albeit with ( ) < 0). When taking the four-dimensional limit of Eq. (3.9), some extra considerations are necessary and will be discussed in Section 4.2. When the inverse Laplace transform (3.1) is considered, this provides a resummation of the logarithms L = ln(Q/v) in Σ(v).
We can further manipulate Eq. (3.9) and perform a change of evolution variable from energy to the transverse momentum of the soft radiation. At the leading (single) logarithmic level, we can replace δ(Q − ω a ) with δ(Q − k ta ), where k ta is the transverse momentum of k a with respect to the emitting dipole {12}. This is because for soft radiation emitted with a large angle in the event centre-of-mass frame one has k ta ∼ ω a up to a O(1) function of the pseudo-rapidity of the radiation. The latter function gives only rise to NLL corrections, which are entirely accounted for by the source u(k). This gives (3.11) with boundary condition given in Eq. (3.10).
We now rewrite eq. (3.11) in a physically appealing integral form. Let us first introduce the Sudakov form factor ∆ 12 (Q) where the phase space measure as a function of k t and the rapidity η with respect to the emitting dipole {12} is given by The upper bound of the rapidity integral can be consistently expanded around its soft limit as where we neglect O(k 2 t ) terms as they only give rise to subleading power corrections. Eq. (3.11) can then be written as We can now bring the previous equation into an integral form as The above equation identifies the Sudakov form factor ∆ 12 (Q) with the no-emission probability, and the second term on the right-hand side represents the total probability for the emission of one gluon. This physical interpretation of the Sudakov form factor will be instrumental later when deriving the NLL evolution equation. We note that the k t -ordered formulation we have just introduced is significantly different from the energy-(ω-)ordered case of Ref. [3] in the collinear limit, namely when k a is emitted collinear to either of the dipole ends p 1 or p 2 . Specifically, the difference appears in a situation in which k a is collinear to, say, p 2 , and hence k ta → 0. This will automatically also introduce an exponential suppression for all emissions off the dipole {1a} that still has a significant angular phase space available for further emissions. This suppression is in fact not present in the energy-ordered case, as long as ω a is different from zero. In problems sensitive to soft radiation only, such as the one considered in this article, this limit is entirely irrelevant. This is because when k a is collinear to one of the dipole ends one has u = 1 (and hence G = 1), giving no contribution to the observable. This ensures that for the problem at hand working in terms of dipole transverse momenta rather than energies is equivalent. However, problems with sensitivity to collinear radiation (such as observables with Sudakov double logarithms) would present some extra subtleties and some care must be taken in this respect. As stressed multiple times in this paper, we do not consider this class of observable here. Before moving on with the NLL evolution, we discuss the evolution properties of the soft factor S 3 defined on three light-cone directions n 1 , n 2 , n 3 . According to Eq. (3.1), its Laplace transform is defined by In the planar (large−N c ) limit, the emission of a hard gluon p 3 off the initial dipole defined by the quark and anti-quark momenta p 1 and p 2 creates two adjacent dipoles {13} and {32}. In this limit, these two dipoles radiate incoherently, and therefore one can write The above replacement is exact in the large−N c limit, and therefore valid at all logarithmic orders. The evolution of each of the factors in the r.h.s. of the above equation is described by Eq. (3.16). One final aspect that we need to discuss is the initial scale of the evolution for G 13 and G 32 . In this case, the hard gluon p 3 carries a significant fraction of the centre of mass energy √ s. Therefore, in the transverse momentum ordered picture considered here, any subsequent soft emission will have a dipole transverse momentum smaller than that of p 3 with respect to the {12} dipole. We therefore set the initial scale for the evolution of the {13} and {32} dipoles to p t3 where We now observe that p t3 ∼ Q, since S 3 is convoluted with H 3 which vanishes by definition in the soft limit, that is p t3 Q. Therefore, up to subleading (NNLL) corrections we can expand p t3 about Q and write (3.21)

The NLL evolution equation for S 2
At the NLL order, the evolution equation (3.11) receives radiative corrections to the kernel for the evolution of the soft factor S 2 . Conversely, as stressed in the previous section, the soft factor S 3 obeys the factorisation (3.21) into two independent colour dipoles, each of which evolves according to the LL evolution equation (3.11). Given the technical nature of the derivation of a NLL evolution equation, we explain its structure in this section, and leave the detailed derivation to Section 4 for the interested reader. We start by considering the LL evolution equation (3.11) where we introduced a short hand notation for the LL evolution kernel (3.23) From the previous equation, we see that the LL evolution kernel is made of a term which describes the real emission of a soft gluon, and the corresponding virtual corrections at one loop. The combination of the two is finite for an infrared safe observable described by the source u(k). We note that derivative of G 12 with respect to ln Q in Eq. (3.22) singles out, by construction, the contribution of the hardest of the soft gluons in real configurations. Conversely, virtual corrections are included by unitarity as discussed above. At NLL, we need to compute the radiative corrections to the r.h.s. of Eq. (3.22). At this order the ln Q derivative will resolve at most two emissions with comparable transverse momentum (or equivalently, energy) as well as the virtual corrections to the emission of the hardest gluon considered in the LL kernel (3.23). Specifically, as shown in detail in Section 4, these corrections can be all computed from known amplitudes. We can then parametrise the NLL evolution equation as follows In the next section we will address the calculation of each of these three ingredients necessary to formulate a complete NLL evolution equation. The derivation will be carried out in d = 4 − 2 dimensions, but we will also present a simple procedure to perform a local subtraction of the IRC divergences and formulate each of the three contributions to the NLL kernel (3.24) in a form that is manifestly finite in d = 4 dimensions. The final results will be given in Eqs.

Derivation of the NLL evolution equation
We now extend the evolution equation (3.16) to NLL accuracy. We stick to the transverse momentum with respect to the emitting dipole as our evolution variable, although one could alternatively use energy as done originally, which would lead to a slightly different form of the final equation, with solutions identical up to subleading logarithmic terms.
It is instructive to start by performing the first iteration of Eq. (3. 16), and expand in a fixed number of emissions. We obtain where we neglected terms describing the emission of more than two gluons in large-N c limit.
With a slight abuse of notation, we now denoted with k ta the usual transverse momentum of k a with respect to the {12} dipole, and with k tb that of k b with respect to the emitting dipole, either {1a} or {a2}, in the last two terms of Eq. (4.1), respectively. Explicitly, each theta function has to be interpreted according to the following definition: tb is the transverse momentum of k b with respect to the "emitting" dipole {ij}. The first line of Eq. (4.1) encodes the single real emission at tree level, while the last four lines encode the double real correction in the two colour flows that contribute to the large-N c pattern. However, these corrections are only accounted for in the strongly ordered limit in Eq. (4.1), and we must instead use the full unordered limit to account for NLL corrections.
There are two types of contributions that arise when we consider the corrections to the squared amplitude of Eq. (3.2). Eq. (3.2) is obtained in the limit of emissions strongly ordered in energy, or equivalently in dipole transverse momentum (we recall that we work in the limit of soft radiation emitted with wide angle w.r.t. the Born legs). This limit leads to the leading terms (α s L) n resummed by Eq. (3.16). We start by noticing that in the kernel of the LL evolution equation (3.16), the eikonal squared amplitude w (0) 12 (k a ), describes the emission of the soft gluon which carries the largest transverse momentum k ta with respect to the emitting dipole. To gain control over NLL terms of order α s (α s L) n , we need to account for the corrections of relative order O(α s ) to the above kernel. These are discussed in the following.
Real corrections. The squared of the double soft current, which is needed to describe the emission of two soft partons of commensurate energy, can be split into two terms in the large-N c limit asw  The colour-ordered double soft squared amplitude at tree level can be found in Section 5.3 of Ref. [44] (see also Section 8.1 of Ref. [45]) and it is reported below 1 where the Lorentz invariants s i...k indicate the standard Mandelstam variables. For later use, it is also convenient to single out the independent emission contribution and writẽ Note, however, that the separation of the independent contribution is immaterial at the level of the single colour flow, and only makes physical sense at the level of the sum w . We also observe that the correlated termw (gg) 12 is not positive definite. The separation (4.4) is useful since the independent emission contribution is correctly described by the BMS equation as it simply arises from the iteration of the leading order squared amplitude. We can therefore focus on the correlated term of Eq. (4.4) in what follows.
The dipole structure of this double-real correction can be read from the last two terms of Eq. (4.1). In that equation, the first line corresponds to the usual dipole structure of the BMS equation, obtained from the replacement Similarly, the product of three Sudakov factors in the third and fifth line encodes the no-emission probability for the three dipoles created by the emission of k a and k b off the initial dipole {12}. The dipole structure for the double real correction then corresponds to replacing As a last step, we need to upgrade the phase space boundary given by Θ(Q − k ta ) in Eq. (4.1) which emerges from the strongly ordered limit. Conversely, the double real correction given byw (gg) 12 describes the emission of two gluons with commensurate k t (and energy) with respect of the emitting {12} dipole. As a consequence, the definition of the resolution scale Q used in Eq. (3.3) now has to constrain the total energy of k a + k b or equivalently its transverse momentum w.r.t. the {12} dipole, in order for the evolution to be collinear safe for any value of Q. In line with our choice of the dipole transverse momentum as evolution scale, we therefore perform the replacement where we defined k t(ab) = k ta + k tb , (4.8) and k t(ab) = | k t(ab) |. We denoted by k tb the transverse momentum of k b with respect to the {12} dipole and Notice that k tb differs significantly from k tb (in the frame of the emitting dipole of k b ) in the limit in which k b is collinear to k a , and the two have comparable energies. As already stressed, and we will show shortly, in this limit the real correction will cancel against the virtual contributions and therefore this limit is irrelevant for the observable considered here. In this case, one can only receive a NLL correction from the regime where k tb ∼ k tb . Subtracting the double counting with the iteration of the BMS equation (4.1) leads to the following term to be added to the right hand side of Eq. (3.16) where the two emissions are ordered in their dipole k t . In its present form, Eq. (4.10) cannot be readily interpreted as a correction to the kernel of the evolution equation (3.16), in that it contains the product of two ratios of Sudakov factors that already indicate an iteration of the kernel. In order to derive a corresponding term at the level of the integral equation (3.16) we will need to make some considerations. As a next step we discuss virtual corrections at NLL order. The BMS equation already contains virtual corrections in the strongly ordered soft limit. However, these do not cancel the collinear singularity present in Eq. (4.10) when k a is collinear to k b . For this, we need first to introduce the full one-loop corrections to the evolution kernel in the soft limit.
Virtual corrections. The second term to consider is the virtual one-loop correction to the leading order kernel in the r.h.s. of Eq. (3.16). The structure of the virtual corrections can be read off the first line of Eq. (4.1). In the large-N c limit, the emission of k a generates two adjacent colour dipoles that emit further radiation incoherently. Therefore, the singularity structure of virtual corrections to the {1a2} configuration will factorise into the product of virtual corrections to the {1a} and {a2} dipoles in this limit (see e.g. [46]). This factorising structure is already encoded at LL in the first line of Eq. (4.1), specifically in the product ∆ 1a (k ta )∆ a2 (k ta ). However, the LL Sudakov of factors in this product do not contain the correct NLL singular structure (i.e. the single poles), and therefore we need to supplement the first line of Eq. (4.1) with an extra correction factor that accounts for the missing terms. The one loop corrections to the emission of the soft gluon k a off the {12} dipole can be written as [47,48] (given here in the large-N c limit) where V (4.12) The mismatch between Eq. (4.12) and the full virtual correction to the Born process considered here is taken into account in the hard matching coefficient H 2 computed in Sec. 5. Moreover, from Eq. (4.11) we also see that choosing µ = k ta absorbs all the k ta dependence into the running coupling. Renormalising the coupling in the MS scheme (cf. Eq. (3.6)) allows us to write whereβ 0 is defined by the coefficient of N c /π in the large-N c β 0 given in Eq. (3.7). The virtual corrections to the evolution kernel are then obtained from the first line of Eq. (4.1) via the replacement where the function γ(k a , ) determines the matching coefficient necessary to obtain the correct virtual corrections to the emission of a soft gluon. It is obtained by matching Eq (4.13) and the one-loop contribution to the r.h.s. of Eq. (4.14). The expansion of the r.h.s. of Eq. (4.14) gives The coefficient of O(ᾱ) in the above equation has to be matched to the one-loop expansion of Eq. (4.13), from which we obtain The last integral reads Cancellation of collinear singularities. We now combine the real and virtual corrections obtained above into a single integral equation in d = 4 − 2 dimensions. We want to achieve a manifest cancellation of soft and collinear singularities between real and virtual contributions independently of the precise form of the source u(k). For this, we make use of the fact that, as stressed earlier, the double-real corrections of Eq. (4.10) only contribute in the unordered regime where ω a ∼ ω b . For the observables considered in this article, this also implies k ta ∼ k tb ∼ k t(ab) . To see this, we observe that k ta ∼ k tb only if ω a ∼ ω b and the two gluons are not collinear to one another. Away from this configuration, when the two gluons become collinear, one can end up in a situation with ω a ∼ ω b but the dipole transverse momenta are strongly ordered, i.e. k ta k tb . This configuration however only contributes to the observable if both gluons are inside the rapidity slice, hence at wide angles w.r.t. the {12} dipole, and the corresponding collinear singularity exactly cancels against the one in γ( ) (4.18) at all orders inᾱ, leaving behind only finite contribution without a logarithmic enhancement (i.e. NNLL). Conversely, for an observable sensitive also to collinear radiation (e.g. the light-hemisphere mass), this configuration can occur with both soft gluons being simultaneously collinear to one of the {12} dipole ends, and to each other. The singularity structure of this configuration would result in an extra logarithmic enhancement and the argument made above would not hold. This extra logarithmic enhancement can be however resummed by means of standard techniques used for global observables, and the formalism presented here can still be used for the calculation of the non-global contributions. We can therefore perform a first-order Taylor expansion and make the following approximation in the double real corrections given in Eq. (4.10) (and similarly for the second term corresponding to the alternative colour flow) where we systematically neglected corrections of order 20) in the limit of soft-wide-angle radiation. This is because the logarithm of the ratios k ta /k t(ab) and k ta /k tb are not large in the region where Eq. (4.10) is non vanishing, and the quantity (4.20) only gives a contribution of O(1) upon integration over the phase space of k a and k b , that is a NNLL correction. We can write the NLL integral equation in d dimensions as where we have also Taylor expanded the scale of G 12 in the double real correction following the same argument as in Eq. (4.20). We will shortly discuss a way to make the subtraction of IRC divergences manifest and eventually to take the limit → 0.
The NLL Sudakov form factor. As a consistency check of the above equation, we can use Eq. (4.21) to derive the NLL soft anomalous dimension that enters the evolution equation of the soft corrections to the Sudakov form factor given in Ref. [43]. We start by setting the source u to 1 (hence G[1] = 1) in Eq. (4.21) and obtain We now introduce a d-dimensional massless momentum k such that k t := | k t | = k t(ab) and its rapidity η with respect to the emitting dipole {12} is that of the k a + k b system η (ab) , this is defined by the corresponding kinematic map 2 where the three-vectors are taken in the rest frame of the {12} dipole. We recast the above equation as where we defined We now divide Eq. (4.24) by ∆ 12 (Q) and take the derivative with respect to ln Q, obtaining the differential evolution equation 2 An analogous map has been considered in the context of the resummation of global observables in refs. [49][50][51].
In the second line of Eq. (4.26), it is convenient to parametrise the double real phase space as in [52] (see also Ref. [50] for the d-dimensional case) where [dk] is given in Eq. (3.13), m 2 ab = (k a + k b ) 2 ∈ [0, ∞), and z ∈ [0, 1] is the fraction of one of the light cone components of k with respect to the {12} dipole carried by one of the daughters k a or k b . For instance if we adopt the following Sudakov parametrisation for the four momenta with κ i being a space-like transverse vector describing the transverse momentum of k i w.r.t. the {12} dipole, we can define z as Finally, dΩ 2−2 is the angular phase space for the vector q := kta z − k tb 1−z with respect to k t , given by Since k is massless, we can integrate inclusively over m ab in the second line of Eq. (4.26), and find This leads to the well known evolution equation for the soft radiative corrections to the Sudakov form factor (see e.g. [43,53,54]), which can be expressed as [50,53] is the coefficient of N c of the large-N c two-loop cusp anomalous dimension in units of α s /π. One can directly use this result and Eq. (4.21) to derive a differential equation for G 12 in d dimensions. However, since we are ultimately interested in a numerical evaluation of its solution, we first discuss how to take the limit → 0.

The NLL differential equation
In order to take the → 0 limit we wish to make the cancellation of infrared singularities in Eq. (4.21) manifest. This is crucial to solve the evolution equation numerically for different infrared safe observables. To this end we introduce the counterterm is that the generating functionals G 1(ab) and G (ab)2 now depend on the angle of the momentum k (ab) with respect to the {12} dipole defined by the above projection (4.23). We now deal with the integral of the counterterm (4.35), which is instead combined with the real-virtual corrections in the first two lines of Eq. (4.21). Given the massless nature of k (ab) , we can adopt the parametrisation (4.27) and integrate as in Eq. (4.31). The integrated counterterm then cancels the poles of γ( ) in Eq. (4.21), and will later allow us to take the limit → 0 at the integrand level. One can finally derive the corresponding integro-differential equation, and bring it in the form given in Eq. (3.24). As done for the evolution equation for the NLL Sudakov factor, we start by dividing Eq. (4.21) supplemented with the counterterm introduced above by ∆ 12 (Q) and then we take the derivative with respect to ln Q. Using the evolution equation for the Sudakov factor (4.33) one finds that the dependence on the Sudakov factors entirely drops out and we obtain the NLL evolution equation (4.36) where we have defined as the kernel correction due to the virtual and subtracted real-virtual corrections, for the double real corrections, and finally for the subtraction of the double counting with the iteration of the LL evolution kernel. In taking the four-dimensional limit in Eqs. (4.38) and (4.39), while keeping NLL contributions only, we implicitly neglect NNLL configurations in which k a and k b are both inside the rapidity slice. This is also crucial to guarantee the collinear safety of Eq. (4.39). Finally, we observe that since the observables under consideration are sensitive only to radiation at small rapidities, the exact rapidity bound (3.14) in the phase space integral is irrelevant. Therefore, it can be relaxed and set as − ∞ < η a/b < +∞ (4.40) in the differential equations (3.22) and (4.36), which are both collinear safe and thus well defined since outside the slice there is a complete cancellation between real and virtual corrections. On the other hand, in the corresponding integral equations, virtual corrections are encoded in the Sudakov form factors. Therefore, one cannot use Eq. (4.40) since real and virtual terms require a finite rapidity upper bound in order to be separately well defined. In this case the initial rapidity bound (3.14) must be retained although the insensitivity of the observable to the large rapidity region ensures that the dependence on this cancels out eventually. The boundary conditions to Eq. (4.36) in d dimensions are given in Eq. (3.10).

Limit to d = 4 and boundary conditions
We now briefly comment on taking the → 0 limit of Eqs. (3.22) and (4.36). The righthand side of both equations is now finite and well defined in four dimensions as long as the evolution scale is larger than Λ QCD , so that the limit → 0 can be directly taken at the integrand level. However, the boundary condition given in Eq. (3.10) is not well defined in this limit because of the Landau singularity. One has to supplement the evolution equations for G 12 with some non-perturbative modelling that allows the computation at very small scales. A simple prescription in the context of a numerical calculation is to introduce a freezing of the coupling at some small scale Q 0 > Λ QCD such that An alternative prescription is to require that G 12 is kept constant below the freezing point Q 0 , resulting in the boundary condition while the unitarity condition G 12 [Q; 1] = 1 is unchanged. The physical picture corresponding to Eq. (4.42) is that, below the resolution scale Q 0 , real radiation cancels exactly virtual corrections by virtue of unitarity. This of course requires that the freezing scale Q 0 v, which is the typical scale of soft radiation inside the observed rapidity gap. If this condition is satisfied, the dependence on the infrared scale Q 0 becomes numerically negligible. In their four-dimensional formulation, Eqs. (3.11) and (4.36) can be evaluated numerically with the boundary conditions (4.41), or (4.42), either via discretisation techniques or by Monte Carlo methods. We will address the numerical solution in a forthcoming publication.

Computation of the one-loop hard matching factors
We now discuss the hard factors in Eq. (2.5). These account for the contribution of hard radiation that propagates outside of the observed region of phase space (rapidity slice in our case). Below we will give a definition that allows for a numerical implementation of the formulae derived in this article.
To extract the matching coefficients H 2 and H 3 , it is instructive to see how they arise from an explicit NLO calculation in the case of e + e − → jets. We recall that here we use a flavour-based labelling of momenta, where p 1 and p 2 are the momenta of the quark and the antiquark, and p 3 the momentum of an additional hard gluon. The O(α s ) virtual contribution to the cumulative cross section reads (in the MS scheme and with α s = α s (µ)) To obtain the real corrections, and in order to avoid introducing additional abstract notation, we explicitly parametrise the phase space in terms of the energy fraction of the gluon x := 2E g / √ s and the cosine of the angle between the gluon and the quark y := cos θ qg , obtaining The three directions n i (i = 1, 2, 3) in S 3 correspond to the directions of the quark, antiquark and gluon that are entirely specified by the variables x and y. The phase space constraint Θ out (k) ensures that no hard parton is inside the observed rapidity slice, so that the observable is zero in a pure three-parton configuration. It can be entirely parametrised in terms of the kinematics of the gluon k. Configurations with a hard parton inside the slice would just correspond to power corrections and can be accounted for at fixed perturbative order through standard matching procedures. Upon integration over x and y, Eq. (5.2) develops double and single poles in that exactly cancel against those in Eq. (5.1).
In order to derive the expressions for H 2 and H 3 , we now need to subtract the double counting with the real and virtual corrections to the evolution equation at O(α s ). These can be derived from the results of the previous section. The virtual corrections are given in Eq. (4.12) and at O(α s ) read (we set the upper bound of the evolution scale Q = √ s) Similarly, the real corrections can be obtained from the first iteration of the evolution equation (3.16), retaining only the contribution in which the soft gluon k is outside the slice, obtaining indicates that the gluon k defining the third direction n 3 is soft, so no recoil is present in the event kinematics. In the above equation Θ soft out (k) ensures again that no partons are present inside the slice at O(α s ). Its expression differs from that of Θ out (k) in that no kinematic recoil is present in the soft limit and therefore the thrust axis is aligned with the direction of the quark (antiquark).
From there, we can compute the hard contribution to the virtual corrections as which, as expected, contains only a single pole of hard-collinear nature. We now consider the difference between the real corrections given in Eqs. (5.2), (5.4). Before taking the difference, we observe that the phase space for the emission of a soft gluon in Eq. (5.4) is larger than the one imposed by momentum conservation in Eq. (5.2). This difference comes from having expanded consistently the rapidity boundary as in Eq. (3.14) as well as from taking k t ≤ √ s rather than √ s/2 as required by energy conservation. We can therefore recast Eq. (5.4) as the sum of a term with the same phase-space limits as Eq. (5.2) and a remainder as For the first term we can now switch to the same x, y variables used for Eq.
while the integrals in the second line of Eq. (5.6) are now finite and one can take the limit → 0 prior to integration. We can now take the difference between Eq. (5.2) and (5.6) and obtain where we took the → 0 limit of the second line of Eq. (5.6). We expand the integrand in the first three lines in a Laurent series in , obtaining To proceed, we observe that the delta functions δ(1 ± y) act as follows since the gluon is projected onto the collinear limits. Therefore, we can evaluate the integrals in all terms containing δ(1 ± y) and get The hard factors H 2 and H 3 can be directly extracted from the above computation according to the definition given in Eq. (2.8), and they read and where we introduced the projector P soft which maps the gluon momentum k into its soft limit, such that 14) The plus distributions act on the soft factors providing the counterterms necessary to make the integration finite and suitable for a numerical evaluation. They are to be evaluated only within the convolution integral of Eq. (2.8) as they also act on the soft factor S 3 . The expressions for H 2 and H 3 depend on the specific phase-space parametrisation adopted for the three-parton final state as constant terms could be reshuffled between the two contributions, and only the physical combination of the two is invariant. As a check, we can set Θ out (k) = 1 in eq. (5.13) and integrate H 2 and H 3 over all solid angles to obtain which reproduces the well known total cross section for e + e − → hadrons at NLO normalised to the Born result. We finish this section by computing the first-order observableindependent hard contribution to Σ(v), defined by where, to compute H 3 , we need to use the explicit expression for Θ out (k), as follows As pointed out already in Ref. [22], the result depends only on the quantity In particular, it depends only on δ 2 . For c ≥ 1/2, we obtain Eq. (5.16) only accounts for the hard contribution, while there will be an additional constant term at O(α s ) coming from S 2 . In the case of the energy distribution, plugging Eq. (5.16) in eq. (2.5), and expanding at first order with v(k) = 2ω leads to the same result as Ref. [22]. We stress that this comparison must necessarily be carried out at the level of the physical quantity H 1 (c) and not individually for H 2 and H 3 which, as mentioned above, depend on the scheme used to cancel the collinear singularities between the two as well as on the definition of the functions S 2 and S 3 . The region c < 1/2 is not considered in Ref. [22]. We obtain H 1 (c) = −6 + 12δ 2 − 9 2 − 6 ln 2δ 2 1 + δ 2 δ 4 + 2 ln 2 2 − 6 ln 2 − 9 ln(δ 2 ) − ln 2 (δ 2 ) + 4 ln(δ 2 ) ln(1 − δ 2 ) + 6 ln(1 + δ 2 ) − 4 ln 2 ln(1 + δ 2 ) − 4 ln δ 2 ln(1 + δ 2 ) (5.20)

NLL corrections up to O(α 2 s ) and fixed order tests
In order to test our result Eq. (2.5), we perform a fixed order expansion of Σ(v) up to O(α 2 s ). Here we set the evolution scale as Q = √ s. Variations of Q around √ s would correspond to standard variations of the resummation scale (used to assess the corresponding theoretical uncertainty), that we do not consider here in the context of a fixed-order expansion. We first start using the explicit expressions for S 2 (v) and S 3 (v) in the large-N c limit. Once the large-N c result is established, we upgrade it to include finite-N c corrections. Working at NLL accuracy, we obtain where Similarly, also H 2 and H 3 can be expanded in powers of α s , with the convention of eq. (2.9). Therefore, using the explicit expression forᾱ = N c α s /π, and keeping all terms up to order α 2 s , we obtain We now compute all inverse Laplace transforms by observing that they affect only the sources, and not the matrix element squared or the phase space. For u(k) as in eq. (2.4), we obtain We use the above information to compute separately each term that depends on the sources in eq. (6.4). Introducing L := ln(Q/v) = {ln Q/E, ln Q/E t }, we obtain: The term containing w gives a contribution of orderᾱ 2 , without any logarithmic enhancement. Therefore, we can neglect that term and obtain, up to NNLL corrections, The term containingw 12 (k b , k a ) corresponds to correlated emissions. This is the term that gives rise to leading non-global logarithms at this perturbative order. The observable constraints can be rearranged as follows Each line in the above equation gives rise to a finite contribution, without soft or collinear divergences. The term in the first line gives rise to leading non-global logarithms. The term in the second line can be further simplified by observing that, when both k a and k b are inside the slice, the parent k (ab) is forced by kinematics to be inside the slice as well. Therefore, the only non-zero contribution to the second line of eq. (6.9) is given by The last line of the above equation corresponds to a collinear and infrared finite contribution, where the phase-space of the parent gluon is integrated over a region where its rapidity is finite and its transverse momentum bounded from above and from below by two quantities of order v, hence giving clearly a finite contribution of order α 2 s . Similarly, to further simplify eq. (6.9), we keep track of whether the parent k (ab) is inside or outside the slice. This gives Again, the contributions in the last two lines of the above equation corresponds to phase space regions in which the rapidity of the parent gluon is bounded (it is in fact inside the gap, and its transverse momentum bounded by two limits of order v). Such regions give contributions of order α 2 s (with no logarithmic enhancement) and can therefore be neglected. Therefore, the only regions of phase space where correlated soft gluon emission can give rise to logarithmically enhanced contributions are those where the parent is outside the gap and only one of its offspring is inside the slice and vetoed, and another where the parent is inside the gap and vetoed and both k a and k b are outside. This gives the final contribution due to correlated emission: Altogether we obtain, up to NNLL corrections dν 2πiν e νv G The last contributions we need to compute at order α 2 s are those involving convolutions of the hard factors H 2 and H 3 with G  To summarise, the expansion of Σ(v) up to order α 2 s can be obtained by inserting Eqs. (6.6), (6.13) and (6.14) into Eq. (6.4).   Comparison with Event2. We now compare the differential distribution in v obtained from the expansion in eq. (6.19) to the full QCD result obtained with the program Event2 [38]. Specifically, we consider both the transverse energy as well as the energy as is clearly confirmed by Fig. 5, which validates our predictions for the next-to-leading non-global corrections, up to O α 2 s . In the case of the energy distribution and c > 1/2, we also compared the result of our calculation to the analytic result of Ref. [22] and found agreement up to O(α 2 s ) for the values of δ adopted there.

Conclusions
In this article we have presented a formalism for the resummation of next-to-leading nonglobal logarithms in QCD. Problems sensitive to non-global logarithms are ubiquitous in collider physics, whenever one considers observables which are sensitive to QCD radiation only in limited angular regions of the phase space. As a case study, we have considered both the transverse energy and energy distribution in a rapidity slice with respect to the thrust axis in e + e − annihilation. We showed that the resummed cumulative cross section (2.5) can be expressed as a sum of convolutions between hard factors (encoding the contribution of hard radiation) and soft factors (that resum the soft radiative corrections). While the hard factors do not contain large logarithms, the soft factors contain logarithmically enhanced corrections that are resummed by a set of non-linear evolution equations in the large-N c limit. We have computed all ingredients necessary for the resummation of the NLL nonglobal corrections, namely the hard factors H 2 , H 3 up to O (α s ) (cf. Eqs. (5.12), (5.13)), as well as the evolution kernel for the resummation of the soft factors S 2 , S 3 up to O α 2 s (cf. Eqs. (3.22), (4.36)). We used these results to carry out a O(α 2 s ) fixed-order calculation at finite N c that is in excellent agreement with the full QCD prediction at O(α 2 s ) for asymptotically small values of the considered observable.
A natural step forward will be the all-order solution of the evolution equations presented here and thus the resummation of NLL non-global corrections with the corresponding study of their phenomenological impact. We envision that the most efficient way to achieve this is by means of Monte Carlo technology. This will be addressed in detail in a forthcoming publication. We finally stress that the formulation derived in this article is not tailored to a specific observable and thus can be applied to other infrared safe observables sensitive to soft radiation emitted away from the large energy flow of the scattering process. In particular, the entire process dependence is encoded in the hard factors while the evolution of the soft factors is universal, which will ultimately allow the resummation of next-to-leading non-global logarithms for hadron collider observables. The precise calculation of these corrections enters the precise description of several collider observables involving jets. An notable example is the fraction of gluon-fusion events in Higgs production in association with two jets with VBF selections cuts [35], or observables defined by means of jet substructure technology [30,55,56].
The application to observables sensitive to collinear radiation (e.g. the light hemisphere mass in e + e − ), on the other hand, requires extra care since our formulation does not immediately apply in these cases. One possible way forward would be the calculation of the global corrections using standard resummation technology supplemented by a non-global correction factor obtained by subtracting the global contributions during the numerical integration of the evolution equation derived in this article. This will require a careful extension of the Monte Carlo method of Ref. [1], where this subtraction is carried out at LL accuracy. Another interesting future direction will concern the comparison of our formalism to the formulations of Refs. [16,22]. This can be done by explicitly computing the evolution kernels in Eqs. (3.22), (4.36) for a specific source corresponding to a given observable. Moreover, from a theoretical point of view it could be interesting to consider the inclusion of subleading-N c corrections to the evolution kernel. Although such corrections have been observed to be usually small in known cases [9,35], their theoretical understanding must be improved in view of high precision phenomenology.