Parton showering with higher-logarithmic accuracy for soft emissions

The accuracy of parton-shower simulations is often a limiting factor in the interpretation of data from high-energy colliders. We present the first formulation of parton showers with accuracy one order beyond state-of-the-art next-to-leading logarithms, for classes of observable that are dominantly sensitive to low-energy (soft) emissions, specifically non-global observables and subjet multiplicities. This represents a major step towards general next-to-next-to-leading logarithmic accuracy for parton showers.

Parton showers simulate the repeated branching of quarks and gluons (partons) from a high momentum scale down to the non-perturbative scale of Quantum Chromodynamics (QCD).They are one of the core components of the general-purpose Monte Carlo event-simulation programs that are used in almost every experimental and phenomenological study involving high-energy particle colliders, such as CERN's Large Hadron Collider (LHC).Parton-shower accuracy is critical at colliders, both because it limits the interpretation of data and because of the increasing importance of showers in training powerful machine-learning based data-analysis methods.
In the past few years it has become clear that it is instructive to relate the question of parton-shower accuracy to a shower's ability to reproduce results from the field of resummation, which sums dominant (logarithmically enhanced) terms in perturbation theory to all orders in the strong coupling, α s .Given a logarithm L of some large ratio of momentum scales, resummation accounts for terms α n s L n+1−p , N p LL in a leadinglogarithmic counting for L ∼ 1/α s , or α n s L 2n−p , N p DL in a double-logarithmic counting, for L ∼ 1/ √ α s .
In this letter we shall demonstrate a first major step towards the next order in resummation in a full parton shower, concentrating on the sector of phase space involving soft partons.This sector is connected with two important aspects of LHC simulations, namely the total number of particles produced, and the presence of soft QCD radiation around leptons and photons ("isolation"), which is critical in their experimental identification in a wide range of LHC analyses.The corresponding areas of resummation theory, for subjet multiplicity [13][14][15] and so-called non-global logarithms , have seen extensive recent developments towards higher accuracy in their own right, with several groups working either on next-to-next-to-double logarithmic (NNDL) accuracy, α n s L 2n−2 , for multiplicity [43,44] or next-to-single logarithmic (NSL) accuracy, α n s L n−1 , for non-global logarithms [45][46][47][48].
To achieve NSL/NNDL accuracy for soft-dominated observables, a crucial new ingredient is that the shower should obtain the correct matrix element even when there are pairs of soft particles that are commensurate in energy and in angle with respect to their emitter.Several groups have worked on incorporating higher-order soft/collinear matrix elements into parton showers [49][50][51][52][53][54][55][56][57][58].Our approach will be distinct in two respects: firstly, that it is in the context of a full shower that is already NLL accurate, which is crucial to ensure that the correctness of any higher-order matrix element is not broken by recoil effects from subsequent shower emissions; and secondly in that we will be able to demonstrate the logarithmic accuracy for concrete observables through comparisons to known resummations.
We will work in the context of the "PanGlobal" family of parton showers, concentrating on the final-state case [2].As is common for parton showers, it organises particles into colour dipoles [59], a picture based on the limit of a large number of colours N c .Such showers iterate 2 → 3 splitting of colour dipoles, each splitting thus adding one particle to the ensemble, and typically breaking the original dipole into two dipoles.The splittings are performed sequentially in some ordering variable, v, for example in decreasing transverse momentum k t .Given a dipole composed of particles with momenta pi and pj , the basic kinematic map for producing a new particle k is followed by a readjustment involving all particles so as to conserve momentum [60], § 1.For the original PanGlobal NLL shower, the splitting probability was given by Here P ĩ→ik (a k ) is a leading-order QCD splitting function, η = 1 2 ln a k /b k + const., with the constant arranged so that η = 0 when the emission bisects the dipole in the event centre-of-mass frame, and f (η) = 1/(1 + e −2η ) is a partitioning function.Additionally, the MS coupling, α s (k t ), uses at least 2-loop running, and K cmw = 67/18 − π 2 /6 C A − 5/9 n f [61].
In moving towards higher accuracy, the two relevant elements are the analogues of the real and virtual corrections in a fixed-order calculation.We focus first on the real term, where we require the shower to generate the correct double-soft matrix element when two particles are produced at commensurate angles and (small) energies, while well-separated from all other particles.
Our approach is illustrated in Fig. 1.Consider the case where a dipole ab first emits a soft gluon 1, followed by a splitting of the dipole 1b whereby a new particle 2 is emitted, and 1 becomes 1 after recoil.When the branching from Eq. (1) produces a particle 2 from the 1b dipole, if p 1 .p 2 < p 2 .pb , we select the {1, 2} pair as the one whose double-soft effective matrix element needs correcting.To evaluate the double-soft correction to this configuration, we first identify all shower histories that could have produced the same nearby {1, 2} pair.This includes the history actually followed by the shower, as well as the case where 2 was emitted from the a 1 dipole, and two extra configurations where the shower produced a particle 2 before 1, i.e. where, in the second splitting, gluon 1 was radiated with 2 taking the recoil.
Each history h is associated with an effective squared shower matrix element |M shower,h | 2 , reflecting the probability that the shower, starting from the ab system, would produce the {1, 2} pair in that order and colour configuration (we address the question of the flavour configura-tion below).|M shower,h | 2 is evaluated in the double-soft limit ([60], § 2 a).In principle, emission 2 should be accepted with probability ( where |M ds | 2 is the known double-soft matrix element for emitting the {1, 2} soft pair from the ab dipole [62][63][64].In practice, however, there are regions where the shower underestimates the true matrix element, leading to P accept > 1.Nevertheless, we find that P accept always remains smaller than some finite value Ω.We therefore enhance the splitting probability Eq. ( 2) by an overhead factor Ω, and accept the emission with probability P accept /Ω.The numerator and denominator in Eq. ( 3) are evaluated in the same double-soft limit, defined by rescaling p 1 → λp 1 , p 2 → λp 2 and taking the limit λ → 0. This ensures that P accept = 1 when 1 and 2 are well separated, thus not affecting regions where the shower was already correct.
The acceptance procedure is sufficient to ensure the proper generation of the {1, 2} kinematics, but not the relative weights of the a12b and a21b colour connections, which is crucial to reproduce the pattern of subsequent much softer radiation from the {a, 1, 2, b} system, as required for NSL accuracy.To address this problem, we evaluate F (12) shower , the fraction of the shower effective double-soft matrix element associated with the a12b colour connection, and similarly F (12) ds for the full double-soft matrix element, in its large-N c limit [63,64].If the shower has generated the a12b colour connection and F (12) shower > F (12) ds , then we swap the colour connection with probability shower − F We apply a similar procedure when the shower generates the a21b connection.In practice, we precede the colour swap with an analogous procedure for adjusting the relative weights of gg and q q flavours for the {1, 2} pair.An alternative would have been to apply P accept separately for each colour ordering and flavour combination, however when we investigated that option for the PanGlobal class of showers, we encountered regions of phase space where the acceptance probability was unbounded.Illustrative plots of the shower matrix element and corrections are given in the supplemental material [60], § 2 b.Next, we address the question of virtual corrections.When 1 is produced in the deep soft-collinear region of the ab dipole, i.e. θ a 1 ≪ θ ab or θ 1b ≪ θ ab , the inclusion of K cmw in Eq. ( 2) already accounts for second order contributions to the branching probability in the softcollinear region, as required for NLL accuracy for global event shapes.However, in general, K cmw alone is not sufficient when θ a1 ∼ θ 1b ∼ θ ab , notably because of the non-trivial η dependence in Eq. ( 2) and the way in which it connects with the overall event momentum Q. Therefore we need to generalise K cmw → K(Φ 1,ab ), where the full K is a function of the kinematics of 1 and of the opening angle of the ab dipole.In the same vein as the MC@NLO [65] and POWHEG [66,67] methods and their MINLO [68,69] extension, the correct next-to-leadingorder (NLO) normalisation for the emission is given by Here, V is the exact QCD 1-loop contribution for a single soft emission, renormalised at scale µ = k t, 1; 2 is the product of shower phase space and matrix element associated with real 1 → 12 branching, including double-soft corrections; and ∆ (ps,1) 1 is the coefficient of α s /(2π) in the fixed-order expansion of the shower Sudakov factor.To aid in the evaluation of K(Φ 1,ab ) we make use of two main elements: firstly, in the soft-collinear limit, K(Φ 1,ab ) → K cmw ; secondly, both V (Φ 1,ab ) and ∆ (ps,1) 1 are independent of the rapidity of 1, as long as 1 is soft and (for ∆ (ps,1) 1 ) kept at some fixed value of the evolution scale.We can therefore reformulate Eq. ( 5) as In the second term, 1sc is at the same shower scale v as 1, but shifted by a constant in rapidity with respect to ab so as to be in the soft-collinear region, wherein K(Φ 1sc,ab ) → K cmw .The labels r and r sc indicate a regularisation of the phase space, which should be equivalent between the two terms.Specifically, we separate M DS in Eq. ( 3) into correlated and uncorrelated parts, respectively those involving C F C A versus C 2 F colour factors for the qggq matrix element.For the correlated part, we cut on the relative transverse momentum of 1 and 2, while for the uncorrelated part, we cut on the transverse momentum with respect to the ab dipole and impose |∆y 12 | < ∆y max .In practice we tabulate ∆K as a function of θ ab , η1 , and ϕ 1, though one could also envisage on-thefly evaluation.We incorporate ∆K in Eq. ( 2), through a multiplicative factor 1 + tanh . This form keeps the correction positive and bounded.It also leaves the shower unmodified in the hard-collinear region.
We study the above approach with several variants of the PanGlobal shower.All have been adapted relative to Ref. [2] with regards to the precise way in which they restore momentum conservation after the map of Eq. ( 1).This was motivated by the discovery that in higher-order shower configurations involving three similarly collinear hard particles, the original recoil prescription could lead to unwanted long-distance kinematic side effects.Details are given in the supplemental material [60], § 1, and tests were carried out using the method of Ref. [70].We will consider three variants of the PanGlobal shower: two choices of the ordering variable, ∼ k t θ β with β = 0 (PG β=0 ) and 1/2 (PG β=1/2 ), and also a "splitdipole-frame" β = 0 variant (PG sdf β=0 ), which replaces f (±η) → f (±η) in Eq. ( 2), with η = 1 2 log a k /b k .The η = 0 transition region bisects the dipole in its rest frame rather than the event frame.This makes the 1 → 12 branching independent of the 1 rapidity in the dipole frame, resulting in ∆K = 0. Illustrative plots of ∆K and its impact are given in Ref. [60], § 2 c.For the three shower variants, the overhead factors Ω associated with Eq. ( 3) are respectively taken equal to 3.1, 20 and 4, independently of the dipole kinematics.
All results, both with and without double-soft corrections, include NLO 2-jet matching [71], which is required for the NNDL/NSL accuracy that we aim for.Spin correlations [72,73] are turned off, because we have yet to integrate them with the double-soft corrections.The double-soft corrections are implemented at large-N c , in such a way as to preserve the full-N c NLL/NDL accuracies obtained in Ref. [74] for global observables and multiplicities.All events have (positive) unit weight.
To test the enhanced logarithmic accuracy of the shower, the first observable that we consider is the Lund subjet multiplicity [43] in e + e − → q q events.This is a perturbatively calculable observable that is conceptually close to the experimentally important total chargedparticle multiplicity.For a centre-of-mass energy Q and a transverse momentum cutoff k t , the subjet multiplicity has a double-logarithmic resummation structure α n s L 2n with L = ln k t /Q.The PanGlobal showers already reproduce terms up to NDL α n s L 2n−1 .The addition of the double-soft corrections and matching [71] is expected to bring NNDL accuracy, α n s L 2n−2 .To test this, in Fig. 2, we examine NSL accuracy tests: energy in a slice where N ps is the parton-shower result and N nndl (N dl ) is the known analytic NNDL (DL) result [43].The α s → 0 limit follows the procedure from earlier work [2].Eq. ( 7) is expected to be zero if the parton shower is NNDL accurate.The original showers, without double-soft corrections (left), clearly differ from each other and from zero, by up to ∼ 100%.With double-soft corrections turned on (right), all three PanGlobal variants are consistent with zero, i.e. with NNDL accuracy, to within ∼ 1%.
Next we turn to the study of non-global logarithms at leading colour.These were recently calculated at NSL accuracy [45,46,48], α n s L n−1 , and are available in the corresponding "Gnole" code [46].We again consider e + e − events, and sum the transverse energies (E t ) of particles with |y| < 1, where y is the rapidity with respect to an axis determined by clustering the event to two jets with the Cambridge algorithm [75].The fraction of events where the sum is below some E t,max is denoted by Σ and for a given shower we define Fig. 3 (left) shows Σ (ps) nsl /Σ sl for our three PanGlobal variants without double-soft corrections.As expected, they differ.
Fig. 3 (middle) compares our PG sdf β=0 shower with double-soft corrections to the NSL Gnole code, showing good agreement, within < 1%.Gnole has n f = 0 in the real contribution and counterterm, but keeps the full n f = 5 in the running of the coupling and inclusive K cmw ("n real f = 0").Among our showers it is relatively straightforward to make the same choice with PG sdf β=0 , in particular because ∆K = 0. Also, Gnole uses the thrust axis, while we use the jet axis; this is beyond NSL as the two axes coincide for hard three-parton events.Fig. 3 (right) shows the results from our three Pan-Global showers with complete (full-n f ) double-soft corrections included.They agree with each other to within 1% of the NSL contribution, providing a powerful test of the consistency of the full combination of the double-soft matrix element and ∆K across the variants.That plot also provides the first NSL calculation of non-global logarithms to include the full n f dependence.An extended selection of results and comparisons is provided in § 3 of Ref. [60].
We close with a brief examination of the phenomenological implications of the advances presented here.We consider e + e − → Z * → jets at Q = 2 TeV.This choice is intended to help gauge the size of non-global effects at the energies being probed today at the LHC.Fig. 4 shows results for the distribution of energy flow in a rapidity slice, defined with respect to the 2-jet axis, without double-soft corrections (left) and with them, i.e. at NSL accuracy (right).It uses the nested ordered doublesoft (NODS) colour scheme, which while not full-N c accurate for non-global logarithms, numerically coincides with the full-N c single logarithmic results of Refs.[38][39][40], to within their percent-level numerical accuracy [74].With a central scale choice (solid lines), the impact of the NSL corrections is modest.This is consistent with the observation from Fig. 3 that the NLL PanGlobal showers are numerically not so far from NSL accurate.However, the NSL double-soft corrections do bring a substantial reduction in the renormalisation scale uncertainty, from about 10% to just a few percent.Conclusions are similar for H * → gg.
The results here provide the first demonstration that it is possible to augment parton-shower accuracy beyond NDL/NLL.Specifically, our inclusion of real and virtual double-soft effects has simultaneously brought NNDL/NSL accuracy for two phenomenologically important classes of observable: multiplicities, and energy flows as relevant for isolation.It has also enabled the first leading-colour, full-n f predictions for NSL non-global logarithms.Overall, our methods and results represent a significant step towards a broader future goal of general next-to-next-to-leading logarithmic accuracy in parton showers.

SUPPLEMENTAL MATERIAL 1. Event-wide momentum conservation for the PanGlobal shower
Compared to the original formulation in Ref. [2], this study introduces a revised rescaling procedure for the kinematic map of the PanGlobal parton shower.Here we outline an important flaw identified in the original PanGlobal map, before detailing two amended versions which rectify it, one of which has been used for all results in this article.

a. Unsafety of the original formulation
Given a set of splitting variables {v, η, ϕ} generated according to Eq. ( 2), the original PanGlobal kinematic map of Ref. [2] proceeds by first constructing a set of intermediate post-branching momenta according to Eq. ( 1), conserving only longitudinal momentum within the emitting dipole system.All momenta in the event are subsequently rescaled by a factor r = Q 2 /(Q + k ⊥ ) 2 , which serves to bring the total invariant mass of the event back to Q 2 .Finally, a Lorentz boost is applied to all particles restoring the total four-momentum in the event, from r(Q + k ⊥ ) to Q.
In all kinematic limits where shower emissions are well-separated in logarithmic phase-space, i.e. all limits pertaining to NLL accuracy, the rescaling factor, r, tends to one.While investigating kinematic configurations relevant beyond NLL accuracy, we found a specific sequence of emissions, involving three similarly-collinear hard particles, that leads to r differing from one by an O(1) amount.
The fact that a double-unresolved configuration produces long-range (event-wide) O(1) momentum shifts is in contradiction with the fundamental construction principle of the PanScales showers.This spurious effect evades NLLaccuracy tests, which are by definition insensitive to triple-collinear corrections.Nevertheless, since configurations with three similarly-collinear hard particles can be generated at arbitrarily small angles, the resulting large rescalings technically mean that the original PanGlobal rescaling prescription violates infrared-and-collinear safety.
The origin of large rescaling coefficients can be simply illustrated as follows.Consider the situation where, in the event frame, one has a highly-energetic dipole ĩȷ, with E ĩ ∼ Q, E ȷ ∼ Q, and an opening angle θ ĩȷ ≪ 1.We then radiate a hard-collinear emission k from the ĩȷ dipole using the kinematic map of Eq. ( 1).The transverse momentum of the emission (with respect to the dipole) is given by k Consequently, in the ĩȷ dipole rest frame, the original dipole momenta and the transverse contribution are all commensurate, and of order θ ĩȷ Q. Boosting to the event frame, all momenta receive a large Lorentz boost factor γ ∼ (E ĩ + E ȷ)/m ĩȷ ∼ 1/θ ĩȷ .In particular, this results in the transverse component having an energy of order Q in the event frame, or equivalently k ⊥ .Q ∼ (1/θ ĩȷ )(θ ĩȷ Q).Q ∼ Q 2 .Since transverse components are not conserved by the map of Eq. ( 1), one obtains a rescaling factor r that significantly departs from one.
Below, we introduce two small adaptations of the original kinematic rescaling that each cure this issue.

b. Updated prescriptions with local rescaling
The essential idea of the updated prescriptions that we introduce below is to avoid spuriously large rescalings of the event as a whole by absorbing the potentially large energy component of the transverse momentum k ⊥ via a dipole-local rescaling.We start with the PanGlobal momentum map from Eq. ( 1).
In the first of our two updated prescriptions, instead of applying a global rescaling to all the momenta in the event, we rescale pi , pj and pk by a local rescaling factor r L .As with the original variant, we fix the rescaling by imposing that the invariant mass of the event is returned to Q 2 .Defining time-like four-momenta pijk = pi + pj + pk and pm = Q − pi − pj , we find After the (local) rescaling the whole event undergoes a Lorentz boost so as to restore the total four-momentum of the event to its original value Q (as in the original PanGlobal prescription).This is the prescription used throughout this paper.
A second possibility is as follows.Noting that the origin of the issue with global rescaling lies in the large energy component of k ⊥ , one can proceed with a hybrid approach whereby one first applies a local rescaling r ′ L to restore the dipole energy in the original event frame, i.e. imposing r ′ L pijk .Q = (p i + pj ).Q.One subsequently applies a global rescaling r ′ G to restore the invariant mass Q 2 of the whole event, followed by a Lorentz boost, as in the original formulation.We find Note that similar rescaling strategies can be adopted in the formulation of the PanGlobal shower for pp collisions, deep-inelastic scattering and vector-boson fusion (see Ref. [10] for details).The specific choice between the two updated schemes does not affect the tests carried out in this paper.This owes to the following two facts: firstly, in the double-soft region all rescalings go to one (both in the original formulation and in the two amended versions introduced here); secondly, for the special case of the first emission, where there are only three post-branching momenta, all three maps are identical, meaning that the matching procedure is also unaffected [71].
To validate the new rescaling prescription, we have performed numerical tests similar to the infrared-and-collinearsafety tests carried out for flavoured-jet algorithms in Ref. [70].While the original prescription fails these tests, the new approaches described here pass them.Additionally, we have extended those tests to automate the verification of the more stringent PanScales criteria, proceeding as follows.First, we generate an initial set of emissions.Second, we add new emissions widely separated from the initial ones by at least a distance ∆ in a logarithmic phase-space (equivalently, the Lund plane).We then check that the recoil induced on the initial emissions decays exponentially with increasing ∆.While the original global rescaling approach fails these tests (as do standard dipole showers), the new variants pass.Note that if we further impose that the additional emissions are themselves widely separated in logarithmic phase-space, both the original and new rescaling prescriptions pass the test, in accordance with the expectations from the NLL matrix-element tests performed in earlier work.

Matrix element tests and ∆K
This section records expressions for the effective shower matrix elements, before giving illustrative results regarding the validation of the two key elements of our double-soft corrections, namely the real matrix-element corrections and the virtual ∆K correction.

a. Effective shower matrix elements
The total effective squared shower matrix element in Eq. ( 3), at large-N c , describing the radiation of a double-soft gluon or quark pair, {1, 2}, from a dipole ab, is most conveniently expressed as a sum of contributions associated to the two contributing colour connections, a12b and a21b : The contribution associated to the a21b colour ordering has the form shower,h wherein J , η2 , and z 2 are given in terms of invariants s ij = 2p i .pj and s i = 2p i .Q as follows The angular dipole partitioning functions, f (±η 2 ), and the splitting functions, P 1→12 (z 2 ), are as given in Ref. [2]; for PG β=0 and PG β=1/2 we have η2 = η2 , as in Eq. ( 2), while for PG sdf β=0 this is modified to η2 = η 2 .The δ a→a2 factor is equal to one when the radiated pair comprises of gluons, and zero when it consists of quarks.Colour factors C 1 and C 2 are both equal to C F , with C A = 2C F in the large-N c limit.
The final theta function in Eq. ( 12) reflects the ordering of the emissions with respect to the shower evolution variable v.The evolution variable was defined in terms of k t and η in Ref. [2] through k t = ρve β|η| .Θ(v 1 > v 2 ) in Eq. ( 12) is then fully specified given Eq. ( 13) for η2 together with the following additional components colour flow and flavour separation FIG. 5. Left: rate of emission before (dashed-dot red) and after (solid blue) double-soft matrix-element corrections at O(α 2 s ).Right: decomposition in terms of colour flows and flavour channels.
Finally, the contribution to the total effective squared shower matrix element arising from the a12b colour ordering, h∈a12b |M (12) shower,h | 2 , is given by swapping a ↔ b everywhere on the right-hand side of Eq. ( 12).

b. Real double-soft matrix-element corrections
In Fig. 5 (left) we provide an illustration of one of the matrix element tests that we have carried out.It is shown for the PG sdf β=0 shower.We start with a dipole ab with an opening angle of θ ab = π − 2 in the event centre-of-mass frame.From that system, we generate two soft emissions, and select those configurations where the higher-k t emission (1) is in a window −12 < ln k t1 /Q < −11 and 1 < y 1 < 3, while the lower-k t emission (2) satisfies ln k t2 /k t1 > −1.We determine transverse momenta and rapidities in the dipole centre-of-mass frame and, for the purpose of Fig. 5, restrict our attention to configurations for which, in that frame, the two emissions are both in the ab dipole's primary Lund plane [12].The upper panel shows the differential distribution of the rapidity difference between the two emissions, ∆y 21 = y 2 − y 1 .
The results in Fig. 5 have been normalised to 2(α s C A /π) 2 which is the expected result for large ∆y 21 , at large-N c , as long as particle 2 is still soft.The red (dot-dashed) curve shows the default shower, without any double-soft correction, while the black curve shows the actual double-soft matrix element.Both are shown averaged over ϕ 2 , the azimuth of particle 2. The shower and exact double-soft matrix element differ for ∆y 21 in the vicinity of zero.The red dot-dashed curve in the lower panel shows the ratio of the two curves, illustrating the need for O (1) corrections at small and moderate ∆y 12 values.As |∆y 21 | becomes larger, all curves tend to the same limit corresponding to independent emission.Once particle 2 starts to become hard, and the physical phase space boundary is approached, ∆y 21 ≳ 6, all three predictions begin to depart from that of the independent emission picture, with small technical kinematic cuts also playing a role in that region.Crucially, however, throughout this hard-collinear region the shower predictions with and without double-soft corrections are seen to be in perfect agreement.
When the shower is run with the (fully-differential) double-soft correction factor (upper panel, blue solid line), one sees that it agrees perfectly with the double-soft matrix element for moderate ∆y 21 .One important point is that at large positive ∆y 21 , the correction factor does not modify the shower, even though the shower and the double-soft matrix element differ: in that limit, where the hard-collinear splitting function corrections are relevant, the shower already provides the correct answer, and it is important to maintain that correct answer.
The right-hand plot shows the same differential distribution but broken into flavour and colour channels.Let us start with the upper left panel, which shows the qg 1 g 2 q channel, where the particle labelled 1 is always the one with larger transverse momentum, and the order of the particles corresponds to the order of the colour connections.Of particular interest is the region of negative ∆y 21 , i.e.where the rapidity ordering is opposite to the colour ordering.
In this region the true double-soft matrix element is strongly suppressed, as one would expect.However, the shower's suppression is parametrically stronger.The pattern is similar in the top-right panel for the opposite qg 2 g 1 q colour ordering at positive ∆y 21 .Had we attempted to correct the shower for each colour-channel separately, there would have been regions where the acceptance probability in Eq. ( 3) would have become arbitrarily large.Instead the approach of Eq. ( 4) ensures that we only have to make an occasional swap of the colour ordering.The lower panels show the analogous curves for double-soft quark production.

c. ∆K and evaluation of its impact
Recall that for a soft emission probability (from a qq dipole) as given in Eq. ( 2), NSL accuracy requires an extra (1 + ∆Kα s /2π) correction factor.Fig. 6 (left) shows the size of the ∆K contribution, Eq. ( 6), for our three PanGlobal shower variants.It is plotted as a function of the rapidity, η1 of the soft emission, in the case of a back-to-back parent dipole.The shower with the largest correction is PG β= 1  2 , but for the configuration shown here, that correction remains relatively modest, at most a factor of about (1 + α s ) for η1 = 0.The correction for PG β=0 is much smaller.The PG sdf β=0 variant has the property that ∆K is identically zero, a consequence of the fact that the shower's second emission probability is independent of the rapidity of the first emission, causing the two terms in Eq. ( 6) to exactly cancel.
Fig. 6 (centre and right) illustrates the separate impact of the double-soft real matrix element and ∆K corrections on the slice observable of Fig. 3, for PG β=0 (centre) and PG β= 1 2 (right).It shows the difference in NSL contributions between the PG β result and an NSL-accurate reference, which is taken to be the PG sdf β=0 shower including the full double-soft corrections.The red curve shows the difference with no double soft corrections at all, illustrating e.g. the fortuitous agreement with the full NSL result for PG β= 1  2 .Turning on the real double-soft corrections (blue curve) introduces a highly visible effect, bringing the PG β=0 result in better agreement with the full NSL but causing a significant departure from NSL in the PG β= 1  2 case.Including also the ∆K correction (green curve) results in agreement with the NSL result for both showers.The sign of the ∆K effect is consistent with the left-hand plot: ∆K is always positive, and the resulting higher emission probability reduces the value of Σ.
Finally, let us comment on the numerical accuracy of our results.For λ = −0.35,we find Σ nsl /Σ sl = 4.832 ± 0.004 (PG sdf β=0 ), 4.817 ± 0.010 (PG β=0 ) and 4.787 ± 0.014 (PG β= 1 2 ), where the quoted uncertainties are purely statistical, as obtained from a cubic polynomial extrapolation α s → 0. These numbers are roughly within 2σ of each other.Note however that for PG β= 1  2 , we found the convergence with α s to be slower, making the extraction numerically more challenging.Accordingly, one should also keep in mind that this comes with additional systematic effects.For example, we observed that varying the set of α s values yields variations in Σ nsl /Σ sl of the order of 0.01.We also estimated the effect of varying ∆K within its numerical uncertainty to be of order 0.005.In all cases, we see a convincing agreement to within 1% relative to the size of the NSL correction.

Reference NSL results for non-global logarithms
In this last section, we provide additional results for non-global observables.We consider the transverse energy in a square patch of fixed extent in rapidity and azimuthal angle.In each case, we study the next-to-single logarithmic contribution, normalised to the single-logarithmic result, Σ nsl /Σ sl .We have extracted Σ nsl using the same variants of the PanGlobal shower as in the main text.Our results are presented in Fig. 7, showing an excellent degree of agreement, at the 1-2% level, between the showers across the whole set of observables.These can also serve as reference results for future studies of non-global logarithms at NSL accuracy.

FIG. 2 .
FIG.2.The result of Eq. (7) for three variants of the Pan-Global shower without double-soft corrections (left) and with them (right).The latter are consistent with NNDL accuracy.The bands represent statistical errors in an αs → 0 extrapolation based on four finite αs values.

FIG. 3 .
FIG. 3. Determinations of Σ (ps)nsl /Σsl for the transverse energy in a slice.Left: parton showers without double-soft corrections illustrating NSL differences between them.Middle: with double-soft corrections but n real f = 0 (cf.text for details), for comparison with the Gnole NSL code.Right: with full double-soft corrections, showing NSL agreement between the three PanGlobal showers.

FIG. 4 .
FIG. 4. Distribution of energy in a slice |y| < 0.5 for the PanGlobal shower without double-soft corrections (left) and with them (right).The bands represent renormalisation scale variation, with NLO scale-compensation enabled only for the results with double-soft corrections.

FIG. 6 . 2 (
FIG.6.Left: Plot of the NLO ∆K correction, Eq. (6), for three variants of the PanGlobal shower, as a function of the rapidity η1 of a soft emission from a back-to-back dipole.Kcmw is given for reference.Centre/right: impact of different parts of the double-soft correction on the NSL contributions for the transverse energy in a slice, showing the difference between PG β=0 (centre) or PG β= 1 2 (right) and a reference NSL-accurate shower.

2 NSL
FIG. 7. Next-to-single non-global logarithms contributions for the transverse energy in patches of given rapidity and azimuthal angle range.Each panel corresponds to a fixed rapidity acceptance and shows different ranges in azimuthal angle (different colours) as well as each of the three PanGlobal showers including double-soft corrections (different line styles).