A first look at two-loop five-gluon scattering in QCD

We compute the leading colour contributions to five-gluon scattering at two loops in massless QCD. The integrands of all independent helicity amplitudes are evaluated using d-dimensional generalised unitarity cuts and finite field reconstruction techniques. Numerical evaluation of the integral basis is performed with sector decomposition methods to obtain the first benchmark results for all helicity configurations of a 2 to 3 scattering process in QCD.


INTRODUCTION
As data continues to pour in from the LHC experiments, the precision of many theoretical predictions for high energy scattering processes are being challenged by experimental measurements. While there has been remarkable progress in Standard Model (SM) predictions for multi-particle final states at next-to-leading-order (NLO) and 2 → 2 scattering processes at next-to-nextto-leading order (NNLO), the computational complexity of 2 → 3 scattering processes at NNLO results in many important measurements being currently (or in the near future) limited by theoretical uncertainties. Pure gluon scattering at two loops in QCD is a key bottleneck in making such predictions which have been known for gg → gg for more than 15 years [1,2]. The one-loop five-gluon amplitudes have been known since 1993 [3] and were among the first results from the on-shell methods that led to the modern unitarity method [4,5].
In this letter we demonstrate how new evaluation techniques based on generalised unitarity [6,7] and integrand reduction [8][9][10][11][12][13][14] can offer a solution to the traditional bottlenecks in these computations and present the first results for a complete set of planar five-gluon helicity amplitudes in QCD. The results extend previous results obtained for 'all-plus' helicity amplitudes [15][16][17][18][19][20][21][22]. These on-shell techniques have also been explored in the context of maximal unitarity [23,24] and numerical unitarity [25][26][27] approaches to QCD amplitudes. Work in this area has received considerable interest due to the phenomenological importance of precision predictions for 2 → 3 scattering. Efforts to complete the unknown two-loop amplitudes for processes such as pp → 3 jets, pp → H + 2 jets or pp → γγ+jet have been further motivated by the recent analytic computations of the planar master integrals (MIs) [18,28] using new differential equation techniques [29,30].
Our approach exploits a parametrisation of the multiparticle kinematics with rational functions combined with numerical evaluation over finite fields [31] to avoid the large intermediate algebraic expressions that traditionally appear. The rational parametrisation of the external kinematics is provided by momentum twistor coordinates [32].

INTEGRAND PARAMETRISATION AND RECONSTRUCTION
We define the unrenormalised leading-colour (planar) five-gluon amplitudes using the simple trace basis: where T a are the fundamental generators of SU (N c ) and S 5 /Z 5 are all noncyclic permutations of the external particles. The overall normalisation is n = m N c α s /(4π) where α s = g 2 s /(4π) is the strong coupling constant and m = i(4π) e − γ E (γ E is the Euler-Mascheroni constant). The L-loop partial amplitude A (L) can be constructed from colour ordered Feynman diagrams. In this article we will compute the pure gluonic contributions to these amplitudes at two loops including the dependence on the spin dimension, d s . Results in the 't Hooft-Veltman (tHV) and four-dimensional-helicity (FDH) schemes can be obtained by setting d s = 4 − 2 and d s = 4 respectively [33].
The integrand of the ordered partial amplitudes can be parametrised in terms of irreducible numerators, ∆, are built from 425 irreducible numerators with 57 distinct topologies. 18 of these 57 can be extracted from the (1-loop) 2 cut configurations as shown in Fig. 1. This means that all topologies with an additional propagator including k 1 + k 2 are computed simultaneously with the (1-loop) 2 cuts. This is more efficient since the parametrisations of the cut loop momentum solutions are much simpler. The remaining 39 can be extracted from a further 31 configurations shown in Fig. 2. The 8 topologies shown in Fig. 3 have divergent maximal cuts and are extracted simultaneously with sub-topologies within the set of 31 2-loop cuts.
The construction of an integrand basis has been discussed before using the language of computational algebraic geometry through polynomial division over a Gröbner basis [10,14]. In this work we took a simpler approach which did not rely on the computation of a Gröbner basis, instead relying on the inversion of a linear system which can be performed efficiently with finite field reconstruction methods. We begin by expanding the loop momenta around a basis of external momenta and transverse directions (similarly to the methods of Van Neerven and Vermaseren [34]), where k lives in the physical space spanned by the external momenta of the topology and k ⊥ lives in the transverse space. We further decompose the transverse space into four dimensional and (−2 ) dimensional spaces, ⊥,i . The size of the 4-d transverse space (which we will call the spurious space) has dimension d ⊥, [4] = 4 − d where d is equal to the number of independent momenta entering the vertices of the topology, up to a maximum value of four. We choose a spanning basis v for the physical space of each topology k µ ,i = d j=1 a ij v µ j and a basis w for the spurious space k µ, [4] ⊥,i = d ⊥, [4] where D are the inverse propagators and k i .q j are the physical space irreducible scalar products (ISPs) for a given topology, where q j are suitable linear combinations of external momenta. The coefficients in the spurious and (−2 )-d spaces are functions of additional ISPs k i .w j and µ ij = −k Having completed this decomposition we find relations between monomials in the ISPs by expanding Eq. (3), From this equation it is easy to obtain a valid basis of monomials for each irreducible numerator of a dimensionally regulated amplitude by using Eq. (4) to remove dependence on the extra dimensional ISPs. This basis is just the most general polynomial in the ISPs k i .q j and k i .w j where the power counting is restricted by the renormalizability constraints [35]. This basis is trivial to obtain without polynomial division but results in high rank tensor integrals with a complicated infrared (IR) pole structure. Instead we prefer to map to a new basis which prefers to keep monomials in µ ij in the numerator and make the → 0 limit easier to perform. The map to the new basis is performed in four steps: 1) write down a complete set of monomials in k i .q j , k i .w j and µ ij obeying the power counting restrictions. 2) Order the monomials with respect to a set of reasonable criteria (for example prefer lower rank monomials or prefer monomials proportional to µ ij ). 3) Map all monomials onto the simple basis and construct a linear system according to the ordering of variables. 4) Solve the linear system for the independent monomials in the new basis. The result of this procedure is a process FIG. 3. The 8 distinct topologies with divergent cuts that must be computed simultaneously with subtopologies independent basis of monomials whose coefficients can be fixed from unitarity cuts in six dimensions. We take a top-down, OPP-like, approach to solving the complete system using information from previously computed cuts to remove known poles from the factorised product of tree amplitudes using the six-dimensional spinor-helicity formalism [36]. The product of tree amplitudes is efficiently evaluated by sewing together Berends-Giele currents [37] as described in Ref. [31].
After completing the integrand level reconstruction, the remaining transverse integration must be performed to obtain a form compatible with traditional integrationby-parts (IBP) relations.
Following a recent approach [38], we have two options in order to achieve this: 1) to integrate the full transverse space to remove k i .w j and µ ij introducing dependence in into the integral coefficients or, 2) integrate only over the spurious space retaining µ ij dependence which can subsequently be removed through dimension shifting identities. In this work we have taken the second approach since it turned out to have better numerical stability to use dimension shifted integrals instead of high rank tensor integrals.
In either case the tensor structure in the transverse space can only involve the metric tensor g µν ⊥ (or g µν ⊥, [4] , g µν ⊥,[−2 ] depending on the particular transverse space being integrated out). This makes the tensor decomposition for non-vanishing integrals in the spurious space rather simple. Further examples of this technique can be found in Ref. [38].
We build integration identities and certain symmetry relations (for example k 1 ↔ k 2 in the 3-propagator sunrise topology) into the integrand basis by using them to create spurious numerators. For example, rather than fitting the coefficient of (k 1 .w 2 ) 2 we replace it with the function [4] k [4] ⊥,1 .k [4] ⊥,1 , which will integrate to zero. In Tab. I we summarise the result of our fit to unitarity cuts listing the number of non-zero coefficients at the integrand level before and after performing the integration over the spurious space. Cuts with scalar loops are required for the reduction from 6 to 4−2 dimensions. We perform the fit taking into account the individual contribution of these scalar loops in order to reconstruct the dependence of the numerator on the spin dimension d s . Setting d s = 2 gives a supersymmetric limit in which the highest rank tensor integrals do not appear in the amplitudes. We use a polynomial expansion of the integrand in (d s − 2) to separate the coefficients into terms of increasing complexity. The fit can be performed efficiently using rational numerics for each phase space point and in most cases it was possible to obtain completely analytic expressions for the integrands of the helicity amplitudes using modest computing resources.

NUMERICAL EVALUATION
The unitarity based method outlined above has been complemented by an approach based on numerical evaluation of Feynman diagrams to determine the coefficients of independent monomial bases. Both of these methods use a momentum twistor [32] parametrisation of the external kinematics to obtain a rational numerical phasespace point. This is extremely important since in order to make use of the finite field reconstruction methods our numerical algorithm must be free of all square roots [39][40][41][42]. The parametrisation in this case was chosen (somewhat arbitrarily) to be, where the columns give the 4-component momentum twistors of the 5 external particles (see, for example, Appendix A of Ref. [15] for more details). These methods have been implemented using a combination of tools including Qgraf [43], Form [44,45], Mathematica and a private implementation of the finite field reconstruction method [31].
We have validated our setup on a number of known cases. Firstly, we have reproduced integrand level expressions for the 'all-plus' helicity sector [15] and against the known integrands in N = 4 Super-Yang-Mills theory [46]. The   To obtain a numerical value for the complete amplitude after integration we perform a sector decomposition of the basis integrals combined with Monte Carlo integration. After applying dimension shifting relations [2,47,48] to rewrite the extra-dimensional ISPs as standard integrals we processed the full set of integrals using both Fiesta [49] and pySecDec [50] packages. This setup was validated with the four-gluon helicity amplitudes and cross-checked against results in the literature [25]. Simple topologies with 2 → 2 kinematics were reduced to the known MIs of Ref. [51] using IBPs from Fire5 [52] and Reduze2 [53] and dimensional recurrence relations from LiteRed [54]. This gave a substantial improvement in the numerical accuracy.
The results for evaluation at a specific phase-space point are given in Tables II and III for the amplitudes with helicities λ i and The leading order amplitudes A LO are the tree-level for the --+++ and -+-++ and rational one-loop amplitudes for the +++++ and -++++. The finite (1-loop) 2 configuration A (2), [2] is presented in Tab. IV. Numerical accuracy is not  The numerical evaluation of finite A (2), [2] (1, 2, 3, 4, 5) helicity amplitudes at the same kinematic point of Tab. II. As only one-loop integrals are required for these amplitudes the integration error is negligible.
an issue here since the integrand level reduction already leads to a basis of one-loop MIs. In addition we find complete agreement with the finite part of the known integrated 'all-plus' amplitude [18].
In cases where the pole structure of the amplitudes is non-trivial we compared with the known universal IR structure [55][56][57][58] including the dependence on d s extracted from the FDH scheme results [59]. The leading pole in 1/ 4 was verified analytically and is therefore quoted exactly in Tabs. II and III. By comparing the agreement in the poles between the (d s −2) 0 and (d s −2) 1 we clearly see the effect of the highest rank tensor integrals which only appear in the latter case. We find convincing agreement between the poles and our amplitudes within the numerical integration error [60]. Since the full amplitude is the sum of all three parts we see in this case that the simple (d s − 2) 0 part dominates and the complete amplitude is evaluated with sub-percent level accuracy. This feature is probably not generic for the whole phase-space however.

CONCLUSIONS
The techniques presented in this letter have allowed the first look at a set of five-point two-loop helicity amplitudes with phenomenological relevance for LHC experiments. We have found that unitarity cutting methods in six dimensions can be combined with finite field reconstruction techniques to compute multi-scale dimensionally regulated two-loop amplitudes in QCD. In many cases it was possible to obtain completely analytic expressions for the integrands of the helicity amplitudes.
While a lot of effort was taken to find manageable expressions, the final integrand form was still extremely large and significantly more challenging than the previously known 'all-plus' helicity configuration. One obvious next step is to include a full set of integration-by-parts identities and reduce the amplitude onto a basis of analytically computed MIs. Promising new approaches that use finite field reconstruction [61] or algebraic geometry analyses [27,[62][63][64][65] could make this possible in the near future. We expect there will be other ways to improve the integrand form by using canonical bases [29] and local integrand representations [66][67][68] though at the present time more work is needed to investigate these approaches.
While there clearly remains a long list of tasks to be completed before predictions of 2 → 3 scattering at NNLO in QCD become a reality, the work presented here is the first example of the evaluation of one of the key ingredients. We hope that the techniques and benchmark results presented here will provide a platform towards this final goal.
We thank Fabrizio Caola, Johannes Henn, Claude Duhr, Donal O'Connell, Nigel Glover, Adriano Lo Presti, Francesco Buciuni, and Johannes Schlenk for many stimulating discussions. We would particularly like to thank Claude Duhr for providing computer readable expressions for the integrals in Ref. [51]. It is also a pleasure to thank the organisers of the AMPDEV2017 program at the MITP, Mainz and also the organisers of the "Automated, Resummed and Effective: Precision Computations for the LHC and Beyond" and "Mathematics and Physics of Scattering Amplitudes" programs at the MI-APP centre, Munich for creating stimulating working environments while this work was on-going. This work is supported by the STFC fellowship ST/L004925/1 and grant ST/M004104/1. The work of T.P. has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No 725110). We are also indebted to Samuel Abreu, Fernando Febres Cordero, Harald Ita, Ben Page and Mao Zeng for pointing out typographical errors in a previous version of this letter.