Resummation of Next-to-Leading Non-Global Logarithms at the LHC

In cross sections with angular cuts, an intricate pattern of enhanced higher-order corrections known as non-global logarithms arises. The leading logarithmic terms were computed numerically two decades ago, but the resummation of subleading non-global logarithms remained a challenge that we solve in this Letter using renormalization group methods in effective field theory. To achieve next-to-leading logarithmic accuracy, we implement the two-loop anomalous dimension governing the resummation of non-global logarithms into a large-$N_c$ parton shower framework, together with one-loop matching corrections. As a first application, we study the interjet energy flow in $e^+e^-$ annihilation into two jets. We then present, for the first time, resummed predictions at next-to-leading logarithmic accuracy for a gap-between-jets observable at hadron colliders.

Introduction.-There has been impressive progress in the perturbative calculation of processes at the Large Hadron Collider (LHC).However, for observables involving disparate scales, computations beyond fixed perturbative order are necessary.These include cross sections involving a hard scale Q but with sensitivity to a soft scale Q 0 .Such cross sections involve large logarithms in the scale ratio L = ln(Q/Q 0 ) that degrade the perturbative expansion and should be resummed to all orders to obtain reliable predictions.
A generic set of observables involving scale hierarchies are cross sections where hard radiation is vetoed in certain angular regions.Prime examples are exclusive jet cross sections which require a veto on additional hard jets.While ubiquitous, the all-order resummation of such observables is challenging, since they involve a complicated pattern of enhanced higher-order corrections known as non-global logarithms, which arises due to secondary emissions off hard partons [1][2][3].At leading-logarithmic (LL) ∼ (α s L) n accuracy, resummed results both at large [1][2][3] and finite N c [4][5][6][7] are available.Despite continued progress in the understanding of non-global observables over the past 20 years , a full resummation of next-to-leading logarithmic (NLL) ∼ α s (α s L) n corrections remained elusive.In this Letter we solve this problem based on a factorization theorem [13,14] obtained in soft-collinear effective field theory [36][37][38].The factorization theorem splits the cross section into hard and soft functions.To resum the large logarithms, one solves the renormalization group (RG) equations of the hard functions to evolve them from a scale µ ∼ Q down to µ ∼ Q 0 .Since the associated anomalous dimension is a matrix in the (infinite) space of particle multiplicities, we resort to Monte Carlo (MC) methods to solve the RG equations.A key ingredient for NLL resummation is the recently extracted two-loop anomalous dimension [39] that we implement into a parton shower framework, which iteratively generates additional emissions to solve the RG equations.Combined with the one-loop corrections to the hard and soft functions we obtain in this Letter the full set of NLL contri- 1. Representation of the factorization formula (1).The blue lines depict hard radiation associated with the energy scale Q, which is constrained inside the jet cones, while the red lines represent the soft radiation at lower energies Q0.
The soft radiation can cover the entire phase space.
butions for gap-between-jets cross sections at lepton and hadron colliders.For the lepton-collider case NLL results were first presented in [40], based on a very different formalism [41], and we find full agreement within numerical uncertainties.
Methodology.-The basis for our resummation are factorization theorems for jet production in the presence of a veto on radiation in certain angular regions of the phase space.The simplest case is two-jet production in e + e − collisions, which factorizes as [13,14] where Q is the center-of-mass energy and Q 0 is the energy scale above which we veto radiation in the gap outside the jet cones.We impose the veto by demanding that the transverse energy E T of the particles in the gap is below Q 0 .At the order we are working, our constraint is equivalent to imposing that the transverse momentum of the leading jet in the gap region is below Q 0 .Figure 1 shows a pictorial representation of the factorization theorem (1).The hard functions H m describe m hard partons, which we treat as massless, inside the jet cones.To obtain H m , one integrates the squared amplitudes over the energies of the m hard partons while keeping their directions {n} = {n 1 , . . ., n m } fixed.The with ⃗ p tot the total momentum of the final state particles.The constraint Θ in ({n}) = θ in (n 1 ) θ in (n 2 ) . . .θ in (n m ) prevents the hard radiation from entering the veto region, i.e., it forces the hard partons to be inside the jet region.The constant c = e γ E /π was introduced in [39].Each hard particle may emit soft radiation, which we describe through a Wilson line S i (n i ) along the corresponding direction n i .The soft function S m is defined as the matrix element of these Wilson lines where the constraint is imposed on the transverse energy in the veto region outside the jets.Since both hard and soft functions are matrices in the color space [42] of the m partons, we take the color trace, as indicated by Eventually, in order to obtain the cross section, we integrate over all directions {n} which is indicated by the symbol ⊗.
The bare hard and soft functions suffer from divergences that can be removed by renormalization.The cross section (1) must be independent of the associated renormalization scale µ which leads to RG equations To carry out the resummation of large logarithms, we compute the functions H m at the hard scale µ h ∼ Q and use the RG equations (4) to evolve them to the soft scale µ s ∼ Q 0 , where we evaluate the functions S m .After integrating over the directions {n} of the hard partons, we obtain a resummed prediction for the cross section.An advantage of the RG approach is that it provides clear rules for the required ingredients to reach a given accuracy.At LL, we need the lowest order hard function H 2 which in our case involves two hard partons, while the soft function is trivial, i.e. S m = 1.Hard functions with higher multiplicities contribute at higher order.The running at LL is governed by the one-loop anomalous dimension Γ (1) which is the leading-order contribution in the perturbative expansion of the anomalous dimension The one-loop anomalous dimension takes the form This anomalous dimension matrix has a simple physical interpretation.The V ij m terms on the diagonal are associated with soft singularities of virtual corrections between legs i and j, while the R ij m entries are related to real emissions along the direction n q .In the large-N c limit, the sum only includes neighboring dipoles, and the color structure of Γ (1) becomes trivial, which drastically simplifies the solution of (4).In the MS-scheme we obtain where n ab = n a • n b throughout this Letter.The dipole radiator W q ij is the product of the two eikonal factors.To reach NLL accuracy it is necessary to include virtual as well as real corrections to the hard H m and soft S m functions along with the two-loop anomalous dimension Γ (2) which governs the running at subleading accuracy.A result for Γ (2) was presented in [16], but in an unconventional renormalization scheme.The anomalous dimension in the standard MS scheme was extracted in [39] and contains additional terms from expanding angular integrals in d = 4 − 2ϵ around ϵ = 0.In the large-N c limit it reduces to a sum over dipoles Γ (2)  mn = while the general expression also includes three-leg terms [16,39].The three entries describe double virtual v ij m , real-virtual r ij m and double-real d ij m contributions.The result for these reads The directions of the two real or virtual soft gluons are denoted by n q , n r .The angular function K ij;qr can be found in [16,39].The remaining functions are X q ij = W q ij ln(4s 2 q ) and where s q denotes the sine of the azimuthal angle of n q in the rest frame of the emitting dipole, and s qr is the sine of the azimuthal angle difference.These functions emerge after taking into account the extra terms derived in [39].
The resulting angular integrals are Lorentz invariant up to the gap constraints.Manifestly invariant expressions are provided in the Supplemental Material and in [43].
The invariance allows us to generate the emissions in the back-to-back frame of the emitting dipole and was instrumental in finding efficient parametrizations for sampling the integrals, which is crucial for obtaining reliable MC predictions.We note that the angular functions become singular when additional emissions along n q and n r are either collinear to one of the parents n i , n j or to each other.To regularize these singularities we impose a cut tan(ϑ/2) > e −ηcut on all angles ϑ between two directions in the lab frame, with η cut = 5.Once we combine real and virtual contributions the singular regions cancel.We have verified that the remaining cutoff effects are negligible for the values of Q 0 we consider.We observe that M ij;qr is collinear finite -the angular integration over the region where n r is inside the jet vanishes.
To perform the resummation, it is necessary to solve the renormalization group equations iteratively to evolve the hard functions from µ h ∼ Q to the scale µ s ∼ Q 0 associated with soft emissions.We thus calculate where we introduced the evolution time and denoted the LL evolution from k to l partons by U kl (t 0 , t) and the NLL correction through ∆U kl (t 0 , t).We point out that once µ h has been fixed, for instance to µ h = M Z , the mapping µ s → t is unambiguous.The LL evolution factor U kl (t 0 , t) = P exp (t − t 0 )Γ (1)   kl (12) has been implemented in a large-N c parton shower framework [25] following [1].To reach next-to-leading logarithmic accuracy we have augmented the LL implementation with exactly one insertion of the two-loop anomalous dimension Γ (2) .This is sufficient to capture the effects of RG running at two loops.We denote this insertion by In short, this means that we start a LL shower U kl (t 0 , t ′ ) which runs from t 0 until t ′ , where we evaluate all contributions in ( 9), together with a β 1 -correction to the running of the coupling, and then restart a LL shower.A pictorial representation of our NLL shower can be seen in Figure 2. The figure shows an insertion of the twoloop anomalous dimension after the first emission, but in the shower, we insert Γ (2) after after any number of emissions.We provide a detailed description of our MC algorithm to calculate (13) in the Supplemental Material.
To obtain results at full NLL accuracy, it is necessary to include matching terms as well.More precisely, we need one-loop corrections to H 2 , the tree-level result for H 3 and the one-loop soft functions S m .We expand these perturbatively in α s where σ 0 is the leading-order cross section.The hard functions are given by standard QCD amplitudes squared with their infrared singularities subtracted in the MS scheme.The one-loop soft functions S (1)  m , on the other hand, can be calculated by the shower MC code from the final emission into the gap [27].
Gap fraction in e + e − .-Analytical expressions for ( 14)-( 16) have been calculated in [14] and subsequently implemented in a MC framework [27].Here we extend results from [27] from LL ′ to NLL by including corrections to the RG running due to the two-loop anomalous dimension following the aforementioned approach.We work at center-of-mass energy √ s = M Z with α s (M Z ) = 0.119 and use two-loop running for α s .Particles are inside the veto region if their angle relative to the thrust axis is larger than α = π 3 , which corresponds to a rapidity gap of ∆Y = ln(3) ≈ 1.1, see Figure 1.By working with a large angle α ∼ 1 we avoid collinear logarithms, but the underlying formalism is also available for small angles [13,14].We calculate the gap fraction which is the fraction of events with transverse energy E T in the gap below Q 0 .To include power corrections beyond the factorization theorem (1), we match to the fixed-order prediction at order α s .The power-suppressed matching corrections ∆R are included through additive matching.We use profile functions [44] to continuously switch off resummation once the power-suppressed terms become relevant.The shape of the curve at larger Q 0 values is affected by this choice and we use the functions introduced in [27].On the left side of Figure 3 we show our numerical results for the gap fraction of the interjet energy flow in e + e − .The bands arise from varying the soft scale µ s since this effect dominates over the µ h variation.The LL (red curve) in Figure 3 is taken from [4] at finite N c = 3.After adding the NLL corrections computed from our shower (blue curve), we obtain a resummed result that is only missing subleading color contributions at NLL.These effects are expected to contribute at the percent level or below.We would like to emphasize two things.First, the large logarithms significantly reduce the gap fraction R(Q 0 ) in the low-energy regime in comparison to the fixed-order prediction (black line).Secondly, after including the NLL corrections, the width of the scale uncertainty band decreases by a factor of two.We have compared individual ingredients of our computation as well as the full NLL correction to the results of [40,41] which are based on a generating functional formalism implemented in the computer code Gnole.We find agreement within numerical uncertainties, as detailed in the Supplemental Material.Gap fraction in Z production.-As a first hadroncollider application, we consider pp → Z → ℓ + ℓ − [45] and compute the cross section for a gap around the incoming beam, centered at the rapidity of the electroweak boson, as considered in [34].The hadron-collider analog of factorization formula (1) was given in [25,32].In the large-N c limit, most ingredients for the resummation carry over from the e + e − case, in particular, the RGevolution and the one-loop soft function can be directly obtained from our MC code.The only new elements are the hard real emission corrections H (1) 3 .The hard functions are simply the partonic cross sections with their infrared singularities subtracted in the MS scheme which are well known [46][47][48].However, the literature assumes that an infrared safe observable is computed and simplifies the cross section results by integrating soft terms over the full angular region.These simplifications are not applicable to the hard functions because the angular constraint in (2) restricts hard partons to the jet region, while contributions with partons inside the gap are part of the soft function.To avoid a double counting we need to restore the terms which were dropped in [46][47][48].A detailed discussion of this point and explicit expressions for the hard functions in all partonic channels are provided in the Supplemental Material, which includes references [49][50][51][52][53][54][55][56][57].
In the right plot of Figure 3 we show our numerical results for the gap fraction (17).We work at √ s = 13.6 TeV and use the NNPDF23_nlo_as_0119_qed [58] set and associated α s .We again use the LL results of [4] at full color and match to fixed order.We note that the result of [4] does not account for complex phase terms in the anomalous dimension (4).These cancel for e + e − but are present for hadron colliders and lead to double logarithms at higher orders [9].However, for Z production the effect of these so-called super-leading logarithms is strongly suppressed and numerically below the percent level [32].For the same gap size ∆Y = ln(3) ≈ 1.1, we see a similar behavior as in the e + e − case, see Figure 3.To study the effect of large logarithms on the fixed-order expansion, we increase the center-of-mass energy Q to 1 TeV and the gap size ∆Y = 2. Once the hierarchy between Q and Q 0 is large, it is unclear how the renormalization and factorization scales µ r and µ f should be chosen in the fixed-order result.In the left plot of Figure 4 we compare the perturbative predictions for a high value µ = µ f = µ r = 1 TeV (solid black line) and for the choice µ = Q 0 (dashed black line).The difference in the two predictions is an inherent uncertainty in the fixed-order result.The dynamical scale µ = Q 0 leads to better overall agreement with the resummed result, but also to an unphysical negative cross section at low Q 0 .In our effective field theory approach, each part is evaluated at its natural scale.On the right side of Figure 4 we show the LL as well as the NLL results, both obtained in the large-N c limit.In all cases, we observe that at low Q 0 the NLL effects are around 20% and are within the large LL uncertainty band.As expected, the scale uncertainty is greatly reduced at NLL.

Conclusions and outlook.
-By performing resummation at subleading logarithmic accuracy, we have reached an important milestone in the effective field theory based resummation of jet observables.The simple non-global observables computed here are of limited phenomenological interest but can serve as benchmark results for the development of general-purpose parton showers at subleading accuracy [29,33,[59][60][61].Our MC framework applies directly to more complicated observables.Phenomenologically interesting examples include H + jets production with a veto on additional jets, and (the background to) Higgs-boson production from vector-boson fusion with a jet veto for central rapidities.The shower evolution is general; the only additional ingredients to the resummation are the relevant NLO hard functions, which should be extracted from existing fixed-order codes.It will also be very interesting to apply our framework to jets with small radii, for example, to compute the jet mass or photon isolation effects [62] at subleading logarithmic accuracy.The factorization theorems for these observables are available [17,62].
Note added.-After completion of this work, the Letter [63] appeared, which includes the full fermionic contribution in the limit where one counts n F ∼ N c .We have adapted our resummation to capture this contribution as well, see the Supplemental Material for details.
In addition, we need the two-particle terms K ij;qr , which are split up as [16,39]  (A.9) These functions are finite when n q or n r become collinear to the legs n i and n j , however, collinear divergences arise when n q becomes collinear to n r .The terms in the double-emission contribution d m proportional to n F describe the splitting of a gluon into a quark-antiquark pair.The large-N c structure of this contribution differs from the remaining terms describing the emission of two gluons.To account for this in our shower, the insertion of d m is split up according to where the lines indicate the dipole structure after the insertion.The evolution after the insertion in ( 13) is then performed separately for the two parts, according to their dipole structure.On a technical level, it is convenient to treat n q and n r as an additional dipole with V qr = 0 in the n F part, i.e., as a dipole which never radiates.

B. Hard matching corrections for Z production
In this section we derive the radiative corrections to the hard functions H 3 , H which can be obtained by a standard next-to-leading order (NLO) QCD calculation [46][47][48] for the Drell-Yan process [45].The production of a (virtual) color-singlet boson with invariant mass Q in hadron-hadron collisions at center-of-mass energy √ s occurs through the annihilation of quarks.More precisely, we study the process N 1 + N 2 → Z + X had , where N 1 and N 2 are the colliding hadrons (protons at the LHC) and X had is an arbitrary hadronic final state.The cross section can be expressed as a convolution of the partonic cross sections dσij dQ 2 with parton distribution functions For convenience, we have introduced a prefactor which absorbs the electromagnetic coupling α and the left-and right-handed couplings g L,R of the Z to fermions.At NLO, each partonic channel only involves a single (anti-)quark flavor q.The explicit form of the couplings for different electroweak bosons can be found in [49]; for simplicity, we only consider Z production in the present work.The PDFs f i/N (x, µ f ) are related to the probability of finding a parton i inside the hadron N with momentum fraction x, and µ f is the factorization scale.In the center-of-mass frame, we parametrize the incoming hadron-momenta as 2 (1, 0, 0, ±1) which translates to, as usual, p 1 = x 1 P 1 and p 2 = x 2 P 2 at the partonic level.We find it useful to introduce the variables (see e.g.[50]) where ŝ = (p 1 + p 2 ) 2 = x 1 x 2 s denotes the center-of-mass energy squared in the partonic system.In addition, we also introduce the Mandelstam variable t = (p 1 − k) 2 , where k is the momentum of the final-state gluon or quark at NLO. Factoring out the electroweak prefactors as in (B.1), the partonic cross section can be written in terms of hard scattering kernels C ij with y = 1+ t ŝ+Q 2 , which is related to the scattering angle θ in the partonic center-of-mass frame through cos θ = 2y−1.We have written the phase space of the final-state particles dΠ f in terms of y, z and Q.The hard scattering kernels C ij are expanded perturbatively in α s .At leading order, i.e. α 0 s , only the q q channel contributes, whereas at NLO the inelastic qg channel opens up (of course, together with the directly related qq, qg, gq and g q channels).After combining real and virtual corrections and performing the collinear factorization, we obtain the LO and NLO hard scattering kernels in the MS scheme The remaining channels can be reconstructed from these results.Using the fact that (B.1) is scale independent, we may derive the dependence on µ f in (B.5) and (B.6), using the fact that the PDFs obey the DGLAP evolution equations [51][52][53].As pointed out in the main text, the expressions for partonic cross sections are usually simplified by integrating soft contributions over all angles.In our case, this concerns the terms in the last two lines of (B.5), which can be written in compact form as These terms were omitted in [47] since they integrate to zero in perturbative predictions for infrared safe observables.However, we use the partonic cross section to construct the hard functions (2) in the factorization theorem and due to the restriction Θ in ({n}) in (2) also the δ(1 − z) terms are only integrated over the jet region.The part of the cross section where the gluon is inside the gap is included in the soft function and integrating the terms (B.7) in the hard function over the full angle would lead to a double counting.The terms in the last line of (B.5) depend on the renormalization scale which we denote by µ h , since the above partonic cross sections will be part of the hard function.
For the same reason, the coupling α s ≡ α s (µ h ) is evaluated at this scale.We now want to use the results (B.5) and (B.6) for the calculation of the angular convolutions of the hard and soft functions.For the hard matching corrections, we need the one-loop corrections to the hard function H 2 , together with the leading-order hard function H 3 in the combination where the LL soft functions include the evolution from the hard to the soft scale and the leading-order soft functions are trivial.The symbol ⊗ indicates the integration over the directions of the additional partons generated during the evolution.
For the hard function H 2 , the angular convolution in (B.8) is trivial since the two incoming partons are along the beam directions.For the function H 3 , a single nontrivial angular integral remains, corresponding to the angle of the emitted final-state particle with respect to the beam axis.We parameterize this angle in the variable y so that we have where we use that the soft function S 3 (y, t) reduces to the two-parton soft function when the emitted parton is along the beam direction, i.e. S 3 (1, t) = S 3 (0, t) = S 2 (t).The soft functions S 3 (y, t) are computed by starting the LL shower with an appropriate three-parton configuration.
The advantage of combining the real and virtual corrections in (B.10) is that the combined hard functions are directly related to the partonic amplitudes given above.For example, for the q q channel, we have δ(y)+δ(1−y) 2 H 2,q q (z) + H 3,q q (y, z) = σ q 0 C q q (z, y) Θ in (y, z) , (B.11 where we express the in-jet constraint Θ in (y, z) through the partonic variables y and z.For a gap around the rapidity of the Z-boson, we have where The quark-gluon channel is not present at the lowest order, so that In our code, we first calculate the soft functions S qg 3 (y, t), S qq 3 (y, t) numerically via LL evolution for a grid of y values.We then interpolate in y and evaluate the convolution with the partonic amplitudes.The integrations in y, z as well as x i are performed using VEGAS [54,55].To simplify the matching to fixed order, we compute the order α s corrections to the hard and soft functions for N c = 3.

C. Shower algorithm
In the following, we describe the NLL shower algorithm.The LL shower was described in detail in [25,28] and is closely related to the one introduced in [1].Reference [27] detailed the computation of the one-loop soft function and the implementation of H 3 into the shower.We will assume that the reader is familiar with the LL shower and will explain the computation of A fundamental ingredient of our shower is a list of m vectors {n} = {n 1 , n i1 , . . ., n im−2 , n 2 } describing an event E at time t with weight w.This list of vectors defines color dipoles with which we associate a virtual correction where we used The time steps in our LL shower are computed according to the distribution after which an additional parton is generated and added to the event.The shower terminates if the additional parton is in the veto region, see Figure 5.
To compute the correction (C.1), we run two showers.At shower time t 0 = 0 we start a LL shower with a list of vectors {n} along the direction of the hard partons in the Born level process and weight w = 1/V E .In our case {n} = {n 1 , n 2 } contains two particles which are back-to-back.The purpose of the first shower is to generate event configurations with different amounts of particles at subsequent times.For instance, on the left of Figure 5 we have This represents U 2k ′ (t 0 , t ′ ) in (C.1).At t ′ = t 1 , t 2 , . . ., we then insert the two-loop anomalous dimension correction, which increases the number of vectors by 0, 1 or 2, depending on whether v m , r m or d m is computed.After this, we start a second LL shower, which corresponds to the factor U l ′ l (t ′ , t) in (C.1).This is depicted on the right-hand side of Figure 5 for the example of t ′ = t 2 .
Let us now explain in detail how the insertion is implemented in our shower code, using, as an example, a term involving a single emission, such as the X q ij term in r m , or the single angular integral terms in v m .1. Pick a dipole in E i with the probability V ij /V Ei .Then generate an additional direction n q and evaluate its weight ∆Γ, which is just the integrand of one of the terms in Γ (2) − β 1 /β 0 Γ (1) under consideration.
2. For real contributions we need to update the list of vectors by inserting the additional emission between its parents, i.e.E ′ i = {n 1 , . . ., n i , n q , n j , • • • , n 2 }, while we leave E ′ i = E i for virtual contributions.

Compute an insertion weight
where the factor V Ei /V ij cancels the one introduced when selecting the dipole.For concreteness, we described the implementation of terms involving a single emission or a single angular integral.The contributions with two directions n q , n r are implemented analogously.By repeating the entire procedure N times and averaging the histograms, we get a numerical estimate for the integral in (C.1).

D. Comparison with finite-Nc calculation at LL
In [4] the interjet energy flow at LL was computed for N c = 3.We have included this result in the plots presented in the main text.To illustrate the size of the finite-N c corrections to the interjet energy flow, we now compare to the LL result obtained using the parton shower in the large-N c limit.In the left panel of Figure 6 we show the LL results as a function of the evolution time t.The gap fraction is normalized to one at t = 0 and the deviations between the two results then become larger with increasing t.On the plot's upper edge, we show the corresponding Q 0 values for µ h = M Z for illustrative purposes.The black dashed lines represent the values Q 0 = 45 GeV (left line) and Q 0 = 2 GeV (right line).We see that the finite-N c corrections within the energy range we consider in the main text reach at most 10%.
A variation of the scale µ s translates into a variation of the value of t and the gray shaded areas in the left plot show the µ s scale variations around 2 GeV and 45 GeV.We compute the scale from the profile function µ s (xQ 0 , µ h ) introduced in [27] and vary 1/2 < x < 2. The width of the band becomes smaller at larger Q 0 as the profile function switches off resummation and enforces µ s (xQ 0 , µ h ) → µ h for Q 0 → µ h .In the right panel of Figure 6, we show the same plot as a function of Q 0 .We see that the scale uncertainty of the LL result by itself is larger than the finite-N c effects so their inclusion does not improve the accuracy of the predictions unless NLL corrections are included as well.

E. Consistency checks and comparison with Gnole
In [40] first results for the full set of NLL corrections to the interjet energy flow were obtained using the computer code Gnole which is based on the framework presented in [41].Our LL results agree within statistical fluctuations to few permille accuracy with Gnole over the range t < 0.06 relevant for our predictions.This is not surprising since a previous version of our code [25] has been cross checked against the results of [1].
The matching corrections H (1) 2 and H (1) 3 , required at subleading accuracy, agree as well.These matching corrections are calculated using the LL shower and do not probe contributions due to the two-loop anomalous dimension.To test the latter, we isolate contributions in (9) proportional to γ cusp .This is the simplest part of Γ (2) , and we use it to perform a series of consistency checks to validate that our algorithm correctly implements (C.1).First, we can exponentiate the γ cusp contributions and subsequently expand numerically.This is delicate since we have to numerically take the limit of α s → 0 in order to suppress NNLL terms.After taking the limit we can compare to the result with one insertion of the anomalous dimension.As expected, we find agreement within our framework.In addition, we compare to the γ cusp -piece computed with Gnole and agree to better than a percent for t < 0.06.
It is more difficult to test the remaining parts of ( 9) which involve the angular functions K ij;qr and M ij;qr .One important test is that our numerical code reproduces the analytical integrals over these expressions which arise in the NNLO correction of the interjet energy flow at small Q 0 given in [39], for example for the combination d2 Ω r d 2 Ω q K 12;qr + K 21;qr θ in (n q ) θ out (n q ) .(E.1) We have calculated this quantity with our shower and find that the integration converges to the analytical result.
Together with the validation of our implementation of the single insertion of Γ (2) using the γ cusp -piece, this provides a strong check on our shower code.It is not obvious how the more complicated parts of Γ (2) map between Gnole and our RG framework.The Gnole approach works with full four-vectors instead of directions and to obtain the mapping analytically one would have to integrate out the energies of the generated partons.It will be interesting to map the two formalisms analytically, but for the moment we perform the comparison numerically.We find agreement with Gnole within numerical uncertainties, at the level of a few per cent. 2

FIG. 3 .
FIG. 3.Gap-between-jets cross sections at lepton (left) and hadron colliders (right) using Q = MZ .The statistical fluctuations visible in the results are due to the Nc = 3 LL result.

FIG. 4 .
FIG. 4. Fixed-order predictions for different scale choices (left) and resummed results at LL and NLL (right) for the value Q = 1 TeV.

FIG. 5 .
FIG.5.Pictorial representations of the LL shower (left) and the NLL shower (right).Blue lines denote hard inside the jets generated by the shower evolution.The red line depicts a soft emission into the veto region, which terminates the shower.The pink blob is an insertion of the two-loop double-real contribution dm at time t ′ .

3 FIG. 6 .
FIG.6.Comparison between the Nc = 3 (red) and the large-Nc (green) results for Z-boson prediction at LL accuracy.

4 .
Start a LL shower with weight w new = w I /V E ′ i and list of vectors E ′ i at t ′ .The time values arising in the second shower correspond to values of t in (C.1) and the weight w new is filled into a histogram at each time t.
Pictorial representation of the NLL shower.Black lines represent the initial hard partons in H2(µ h ).Blue lines denote hard emissions inside the jets generated by the shower evolution.The red line depicts a soft emission into the veto region, which terminates the shower.The pink blob is an insertion of the two-loop double-real contribution dm at time t ′ .
bare hard functions in d = 4 − 2ϵ are defined as