Multiplicative-accumulative matching of NLO calculations with parton showers

We propose a new approach for combining next-to-leading order (NLO) and parton shower (PS) calculations so as to obtain three core features: (a) applicability to general showers, as with the MC@NLO and POWHEG methods; (b) positive-weight events, as with the KrkNLO and POWHEG methods; and (c) all showering attributed to the parton shower code, as with the MC@NLO and KrkNLO methods. This is achieved by using multiplicative matching in phase space regions where the shower overestimates the matrix element and accumulative (additive) matching in regions where the shower underestimates the matrix element, an approach that can be viewed as a combination of the MC@NLO and KrkNLO methods.


Introduction
Next-to-leading (NLO) order event generators interfaced to Parton Showers (NLOPS from now on) have, in the past decades, become the state of the art for the simulation of hard collider processes.The MC@NLO algorithm was the first one to be proposed [1].Soon after, the POWHEG method was proposed [2] in order to overcome the problem of negative weights inherent to MC@NLO.Novel methods were proposed later [3][4][5].Among them, the KrkNLO method, like the POWHEG method, has the characteristic of being free from negative weights, although its applicability is at the moment restricted to relatively simple processes.
Both the MC@NLO and the KrkNLO methods are intimately intertwined with the parton shower generator that is adopted.In fact, the two methods provide the corrections one must apply to the parton shower result in order to achieve next-to-leading order (NLO) accuracy.This is at variance with POWHEG, which is largely independent from the shower generator.In fact, POWHEG takes care of the generation of the hardest event in such a way that NLO accuracy is preserved, while the PS takes care of the remaining (less hard) radiation.This feature is convenient, since one can generate events in the Les Houches format [6] and then shower them with any available PS generator that complies with the Les Houches format specifications.
In view of recent developments aimed at the improvement of the logarithmic accuracy of the shower (see Refs. [7][8][9][10], [11][12][13] and [14,15]) it is legitimate to ask whether POWHEG will maintain the same advantage, i.e. whether minor modifications to the POWHEG algorithm will be sufficient to guarantee that, when interfaced to an already NLL accurate shower, NLL accuracy will be maintained, or whether it will be more feasible to correct an already NLL accurate shower to achieve also NLO accuracy, or rather whether both alternatives will be feasible, and will productively compete among each other.It is therefore important to explore the full range of options for matching NLO calculations and parton showers.
In this note, we will address certain limitations of MC@NLO and KrkNLO, and show how they can be overcome.The MC@NLO approach is widely applicable, and widely used, but has the undesirable feature of negative weights. 1 On the other hand, the KrkNLO method has positive weights, but it is difficult to extend it to generic processes.We will show that combining the two methods one could achieve both positive weights and unrestricted applicability.
The paper is organised as follows, in section 2 we will compare the POWHEG, MC@NLO and KrkNLO methods by formulating them in a common language.We will do this by extending the formulation of the POWHEG and MC@NLO methods of ref. [19] also to the KrkNLO method.In section 3 we will present the combined method, and in section 4 we will present some of its possible variants.Finally, in section 5 we present our conclusions.

POWHEG, MC@NLO and KrkNLO
We assume that the phase space with radiation Φ can be written in terms of an underlying Born phase space Φ B and three radiation variables, indicated collectively as Φ rad .We also assume that the mapping from Φ to Φ B is such that in the singular (collinear or soft) limit, the Born configuration matches the full phase space with the collinear pair merged into a single parton or with the soft particle removed.Such mappings are easy to realise for processes with a single singular region, while for more complex processes one must separate the real cross section into contributions having a single singular region.In the following illustrative discussion we ignore these complications.
We define the following quantity, a function of the underlying Born kinematics, where the full phase space Φ is defined as function of Φ B and Φ rad , B 0 is the Born cross section, V comprises the virtual corrections and, for hadron initiated processes, the collinear counterterms integrated at fixed underlying Born, and R s is a part of the real cross section that includes all soft and collinear singularities.In other words, R − R s is non-singular.Notice that if R s is taken equal to R, Bs is the inclusive cross section at fixed underlying Born.We imagine that renormalisation has been carried out, and that the infrared divergences arising in V and in the dΦ rad integral of R s are regularised in some way.Notice also that Bs is finite, since the singularities present in V cancel those arising integrating the R s term.
We also define the Sudakov form factor associated with R s where Φ is defined in terms of the variables Φ B and Φ rad , and t is some definition of hardness, depending upon the full phase space with radiation.One may think, for example, that t is the relative transverse momentum of the splitting pair.If R s was taken equal to R, the Sudakov form factor would represent the probability for not radiating anything harder than t.
The hardest event cross section can be represented in both POWHEG and MC@NLO as where dΦ = dΦ B dΦ rad , and t cut represents a lower limit for radiation, needed to avoid the Landau-pole singularities.Events are generated with a probability proportional to each term of the cross section formula.The first two terms are generated with a Monte Carlo technique.In fact, they satisfy the shower unitarity equation which follows from the fact that the expression under the integral sign is an exact differential.
In POWHEG, the generation of the hardness t is uniform in the Sudakov form factor, and one can generate events with the standard shower algorithm by equating a random number with the Sudakov form factor, and solving for t.If the t value so obtained is above t cut , the radiation kinematics is generated, and the event with the hardest radiation is fed to a shower generator, which takes care of adding subsequent (less hard) radiation.
In the case of MC@NLO, the implementation of the first term is generally more involved.If the shower is ordered in transverse momentum, the hardest emission is the first, and the radiation mechanism is the same as in POWHEG, except that it is implemented within the shower generator, rather than in the NLO program.In case of an angular ordered shower, large angle soft radiation is generated first, and the hardest radiation occurs somewhere down the shower. 2It was shown in ref. [2] that, in this case, by suitable rearrangement of the Sudakov factors for each emission, one reconstructs the transverse momentum Sudakov form factor in formula (2.3).
The last term in Eq. (2.3) is non-singular, and thus is dominated by hard radiation.It is handled essentially in the same way in MC@NLO and POWHEG.
The NLO accuracy of Eq. ( 2.3) can be demonstrated by computing the expectation value of a generic infrared safe observable O(Φ) as follows or equivalently, In view of the unitarity equation (2.4), the first two terms in the curly bracket of equation ] in the third term in the curly bracket suppresses the singular region, so that up to higher order terms the Sudakov form factor can be omitted, and the lower integration limit can be set to zero up to a power-suppressed correction.We thus get Inserting into the above equation the expression for Bs , Eq. (2.1), we get which is the correct NLO expression for the observable.POWHEG implements Eq. (2.3) directly.R s is defined as where F (Φ) ≤ 1, and approaches 1 as Φ approaches the singular region, so that R s carries all the singularity structure of R.
In MC@NLO, R s is the shower approximation to R, typically given by the Born term times a DGLAP splitting function.Thus, in MC@NLO one computes a Born configuration using the Bs function, and passes it directly to the Shower Monte Carlo.The R − R s terms are generated separately, and fed directly to the shower.Negative weights can arise at this stage, since there is no guarantee that R − R s is positive.
Using the same language adopted here for POWHEG and MC@NLO, the key formula for the KrkNLO method can be written as follows (with an implicit θ(t − t cut ) in dΦ).In the literature dealing with the KrkNLO method [5,[20][21][22][23][24][25], particular emphasis has been put on the use of a specific scheme for the parton densities, which considerably simplifies the expression for the Bs function.Here we are instead interested in a simpler aspect of the method, which is that to generate NLO accurate radiation, it uses a multiplicative correction, rather than the additive correction of the MC@NLO method.The NLO accuracy of the KrkNLO formula can be simply demonstrated by showing that KrkNLO is equivalent to MC@NLO at the NLO level.In fact, we can rewrite formula (2.9) as The middle term is now insensitive to the soft region, because the factor in square brackets vanishes there, so we can drop the Sudakov form factor and, neglecting terms of NNLO size, set Bs /B 0 = 1.By doing this we recover exactly Eq. ( 2.

3).
The KrkNLO method leads to positive weighted events.On the other hand, unlike MC@NLO and POWHEG, its cross section at fixed underlying Born does not exactly match the corresponding fixed order result, but differs from it by NNLO terms. 3The KrkNLO method, however, generates weighted events, so, the unweighting efficiency may constitute a problem if one wants to apply the method to generic processes without having to do process-by-process adjustments.As a related problem, the shower generator may not cover the full radiation phase space.This is the same as saying that R s can become zero in certain regions, in which case the method is not applicable.For these reasons, it seems difficult to apply the method to generic processes in automated framework, something that has been available in POWHEG and MC@NLO since a long time.

The new method
By comparing the MC@NLO and KrkNLO method using a common language, it becomes clear that the two methods have much in common, and can in fact be merged in such a way that the KrkNLO positivity is maintained, unweighted events can be generated on the fly, and no issues arise from the limited coverage of the phase space by the parton shower code.The merged method is defined by the following formula (3.1)It is easy to see by inspection that this formula has the same NLO accuracy as MC@NLO.In fact, the θ function in the first term of eq.(3.1) is regular, and as such the Sudakov form factor multiplying it and the Bs /B 0 ratio can both be set to 1 up to NNLO corrections.
Having done this, the contribution proportional to the θ function in the first term becomes identical to the last term, except for the theta function, which has the opposite sign in the argument.The two terms thus combine, and the theta function disappears, yielding the MC@NLO formula.
The method of Eq. (3.1) is immediately seen to have positive weights, since the second term in eq.(3.1) is forced to be positive by the theta function, while the factor in the curly bracket of the first term is forced by the theta function to be less than 1.Notice that regions of phase space not populated by the Monte Carlo only contribute to the additive term in the square bracket, since they have R s = 0.
Let us consider how the method of Eq. (3.1) might be implemented in a scenario where it is the parton shower that is responsible for generating R s . 4We first consider the case of a parton shower ordered in a genuine hardness variable (for example, transverse momentum).As in the MC@NLO approach, the NLO program is responsible for generating a sample of Born events, with a distribution corresponding to the correct Bs (Φ B )dΦ B weight.The event is communicated to the parton shower, which generates a first emission, i.e. accounting for a remaining factor S(t Φ , Φ B ) Rs(Φ) B 0 (Φ B ) dΦ/dΦ B in Eq. (3.1).Next, the parton shower communicates the event, with its first emission, back to the NLO program, which evaluates the contents of the curly bracket in Eq. (3.1).The result of that evaluation, which is bounded to be between zero and one, is to be used as an acceptance probability for the event.If the event is accepted, the parton shower then continues from that point, generating the remaining emissions.In addition, the NLO program is responsible for generating a second (positive definite) sample of events, corresponding to the θ(R−R s ) term of Eq. (3.1), covering the part of the real phase space where the shower R s underestimates the true real matrix element.These events are passed to the parton shower code for normal showering.
The above scheme requires a more connected workflow between the parton shower and the NLO code than either of the MC@NLO or POWHEG approaches.Recall that those approaches pass a Les Houches Event file to the shower code, and then let the shower take over from there.In contrast, the scheme outlined above requires additional action from the NLO code after the first parton shower emission.However, this is simply a technical consideration, which in our view is a small price to pay for an approach that, like the MC@NLO scheme, leaves responsibility for the first emission with the parton shower, while eliminating negative weights.
A further remark concerns angular-ordered parton showers.For such showers, the application of the acceptance probability (factor in curly brackets in Eq. (3.1)) after the first emission would not be correct.Instead, one might envisage an approach in which one runs the complete shower, uses a jet algorithm to map the full event to the real phase space and then applies the acceptance probability.

Variants of the method
Eq. (3.1) can be viewed as a special case of a family of approaches, parameterised by a constant c ≥ 1, This formula was obtained by replacing R s → c R s in all occurrences where it appears in Eq. (3.1), except for the Sudakov form factor.First of all, we should convince ourselves that this formula has the correct behaviour near the soft limit, and that it is NLO accurate.This is seen immediately by writing (4.1) as The first line is equal to the KrkNLO formula, eq.(2.9) (leaving aside the t cut , for simplicity).In the second line, the factor in square brackets is of order α S and it multiplies a term dominated by the hard region (and thus also of order α S ), since Accordingly, the second line is of NNLO order, and we can conclude that formula (4.1) is equivalent to the KrkNLO formula up to NNLO corrections.For a shower where R s /R is always larger than some value r, taking c ≥ 1/r causes Eq. (4.1) to reduces to the KrkNLO approach over all of phase space.The implementation of formula (4.1) goes as follows.The soft events are oversampled by a factor of c, and accepted with a probability proportional to the expression in the curly bracket, while the hard events are generated in a standard way.
There are several reasons why it may be of use to have such a family of approaches parameterised by c.One is that it is useful to have a parameter to help gauge systematic uncertainties associated with terms that are beyond the accuracy of the method.In fact, the c parameter also gauges the amount of radiation that is multiplied by the inclusive Kfactor, with respect to the amount that is added in as hard radiation.In this sense, it would play a similar role to the hdamp parameter in POWHEG [2,26].Another reason is that the θ-functions in Eqs.(3.1), (4.1) can, conceivably, introduce non-smoothness in phase-space coverage even when the underlying R s and R functions are smooth.Sampling over a range of c values would allow one to address that issue. 5o far, the variants that we have discussed involve the rejection of events.For showers with an ordering variable that involve a genuine hardness scale, such as transverse momentum, it is also possible to envisage a variant where instead of rejecting the event with probability Rs−R Rs θ(R s − R), one rejects the shower's first emission with that same probability.If that first emission is rejected, at a value of the ordering variable that we label t 1 , the shower then continues from that scale t 1 , based on an event without the first emission.As the shower continues, each time the shower again attempts to create a new first emission, that emission continues to be rejected with probability Rs−R Rs θ(R s −R).Once the shower has generated a first emission that is accepted, the shower continues as normal.This ensures that the shower (with first emission rejection) remains unitary. 6The Born normalisation factor that multiplies the shower generation then needs to be modified to read

Conclusion
In this article we have proposed an algorithm for NLOPS generators that, like the MC@NLO approach, relies on the parton shower code for driving the showering steps, but avoids the negative-weight events that complicate the practical use of the MC@NLO approach.The approach can be viewed as a hybrid version of the MC@NLO and KrkNLO methods.The method should be straightforward to implement, requiring no more information than is already used for MC@NLO codes, though it should be noted that it will probably require a closer integration of the showering and fixed-order codes.Following the tradition of associating an acronym to NLOPS methods, a suitable one for this method could be MAcNLOPS, reflecting the combination of Multiplicative and Accumulative (additive) steps.