Complete reduction of two-loop five-light-parton scattering amplitudes

We reduce all the most complicated Feynman integrals in two-loop five-light-parton scattering amplitudes to master integrals, while other integrals can be reduced even more easily. Our results are expressed as a system of linear relations in block-triangle form, which are very efficient for numerical calculation. Based on this, complete next-to-next-to-leading order QCD calculation for three jets, photons, or hadrons production at hadron colliders becomes possible. In order to find out the block-triangle relations, we develop a new method which is efficient and general. The method may provide a practical solution for the bottleneck problem of reducing multiloop multiscale integrals.

We reduce all the most complicated Feynman integrals in two-loop five-light-parton scattering amplitudes to master integrals, while other integrals can be reduced even more easily. Our results are expressed as a system of linear relations in block-triangle form, which are very efficient for numerical calculation. Based on this, complete next-to-next-to-leading order QCD calculation for three jets, photons, or hadrons production at hadron colliders becomes possible. In order to find out the blocktriangle relations, we develop a new method which is efficient and general. The method may provide a practical solution for the bottleneck problem of reducing multiloop multiscale integrals.
Introduction. -Due to the good performance of the Large Hadron Collider (LHC), we now enter the era of precision physics. Some of the most important observables are three light particles or jets production cross sections [1][2][3], which can both test the strong interaction at high energy and determine the QCD coupling constant. On the theoretical side, predictions with compatible precision are needed, which demands the perturbative QCD calculation up to next-to-next-to-leading order (NNLO). Although great progresses have been made in the past few years , a complete NNLO result is still unavailable. The main obstacle right now is the calculation of two-loop amplitudes.
To evaluate a two-loop five-light-parton scattering amplitude, one usually first generates integrand, then reduces all Feynman integrals to linear combinations of relatively simpler master integrals (MIs), and finally calculates these MIs. Because integrands can be obtained either by unitarity method [4][5][6][7][8][9] or by traditional Feynman diagram method and MIs have been calculated analytically [20][21][22][23][24], the bottleneck is the reduction of Feynman integrals. For example, non-planar contribution of two-loop three photons production at the LHC cannot be calculated for lack of reduction for nonplanar integrals [19].
Reduction is usually achieved by integration-by-parts (IBP) identities combined with Laporta's algorithm [25][26][27][28][29][30][31][32][33][34][35]. Even though many new ideas have been proposed to improve the IBP reduction [36][37][38][39][40][41][42][43][44][45][46][47][48][49] in recent years, the problem of reducing multiloop multiscale integrals has not been fully resolved yet. The difficulty is twofold. On the one hand, due to the number of scales, the explicit solution of IBP system is usually too huge in size to be used for numerical calculation, besides it is very hard to obtain [47][48][49][50][51]. On the other hand, although solving IBP system numerically in a single run is tolerable, one usually needs to solve it for a huge number of times for the purpose of either phase space integration or fitting analytical expressions, which is very time-consuming and resource-consuming. For example, to reconstruct the fully analytical two-loop five-gluon all-plus helicity am-plitude [17], one needs to run numerical IBP for about half a million times 1 . If one uses the same method to reconstruct analytical one-minus or maximal-helicityviolation amplitude, much more times of numerical IBP running may be needed, which is hard to achieve.
In Ref. [52], we pointed out that the difficulty of reduction can be overcome if a system of block-triangle relations are found, which has small expression size and can be solved numerically with very high efficiency. Using our proposed series representation of Feynman integrals as input [52,53], we constructed an algorithm to search for block-triangle relations and obtained some preliminary results in Ref. [52].
In this Letter, by further developing and improving the method in Ref. [52], we successfully find out blocktriangle relations to reduce two-loop five-light-parton scattering amplitudes. As expected, the relations are only 148MB in size, and can be numerically solved hundreds of times faster than other methods. Our work makes complete NNLO QCD calculation for three jets, photons, or hadrons production at the LHC a possibility. As our method is efficient and general, it can be straightforwardly applied to any other process, and thus provides a practical solution for the bottleneck problem of Feynman integrals reduction.
Feynman integrals in two-loop five-light-parton scattering amplitudes. -To obtain the badly needed reduction of Feynman integrals in two-loop five-lightparton scattering amplitudes, we only need to consider integrals originated from the four topologies shown in Fig. 1. All other Feynman integrals are one-loop-like which can be dealt with much more easily.
Let us take the most complicated one, topology (a) in Fig. 1, as an example to explain what kind of Feynman integrals do we need to reduce. There are five external momenta p 1 , · · · , p 5 flowing into the diagram, satisfying on-shell conditions p 2 i = 0 (i = 1, . . . , 5) and momentum conservation 5 i=1 p i = 0. As a result, this problem contains five independent mass scales, which can be chosen as s = {s 1 , s 2 , s 3 , s 4 , s 5 } with s i ≡ 2p i · p i+1 and p 6 ≡ p 1 . A complete set of denominators can be chosen as where the first eight are inverse propagators, and the last three, which are called irreducible scalar products, are introduced to make the set complete.
The family of integrals defined by topology (a) can be expressed as where the indexes ν 1 , · · · , ν 8 are integers, ν 9 , ν 10 and ν 11 are nonpositive integers. Two integrals in this family are called in the same sector if positions of their positive indexes are the same. The degree of an integral is defined by the opposite value of the summation of all its negative indexes. Finally, we call a degree-m integral is m n -type if it has n positive indexes and all these positive indexes are 1. For example, I {1,1,1,1,1,1,1,1,−4,0,−1} is a degree-5 integral in the top sector, and it is 5 8 -type. For convenience, we define operatorsm ± (for nonnegative integer m), which generate a set of integrals in the same sector or its subsectors when acting on an integral. For any integral I ν ,0 ± I ν = I ν , m + 1 ± I ν = m ±1± I ν ,1 − I ν generates a set of integrals with one index decreased by 1, and1 + I ν generates a set of integrals with one nonzero index increased by 1. We also definê m as a list of operators, Note that the definition of these operators are a little bit different from those in Ref. [52].
As is well-known, the most complicated integrals to reduce are top-sector integrals with the highest degree. By studying the two-loop five-gluon scattering amplitude diagram by diagram, we find the highest degree is 5 for top-sector integrals belonging to family (a), which are 5 8 -type. Therefore, we define an integral set In the following, we take the reduction of S (a) as an example to explain our two-step search strategy. At the first step, we want to set up a system of reduction relations that can numerically express all integrals in terms of MIs. It is fine if the system is less efficient in numerical calculation, because it only serves as generating input data for searching procedure at the second step. At the second step, the aim is to find a system of block-triangle reduction relations, which needs to be very efficient for numerical use.
Reduction at the step one.
-We want to set up a set of relations, with which we can express all integrals in S (a) in terms of MIs for any given phase space point (rational numbers for both s and ), with coefficients calculated in the finite field of a 63-bit prime number. Although IBP method [25][26][27][28][29][30][31][32][33][34][35] can do this job, we would like to explain in the following that our method proposed in Ref. [52] may provide a better choice.
For each given integral I ν , called a seed, there are 12 IBP relations among the integral set Besides, there are additional 6 relations due to Lorentz invariance [54], which can be interpreted as linear combinations of IBP relations from other seeds [55]. The above IBP relations can also be found out easily using the method proposed in Ref. [52]. To this end, we introduce a parameter η for all integrals in G IBP ν , and then search relations among them using input from series representation [52,53]. Up to d max = 1, where d max is half of the maximal value of mass dimension for coefficients of relations, we can find at least 12 relations; while up to d max = 2 we can find at least 12+6 relations. Because these relations are analytical in η, we can take η → 0 directly and recover the aforementioned 12 + 6 IBP relations.
The advantage of our method in Ref. [52] is that it has freedom to search relations among any set of integrals. As the simplest generalization of G IBP ν , we can define an integral set and search relations among them. Up to d max = 2, there are typically 2 more relations besides 12+6 IBP relations for each seed. With more relations in hand, it is possible to select better relations to achieve a more efficient reduction. For example, our relations from all 4 8 -type seeds can already reduce 15 out of all 5 8 -type integrals to integrals with lower degree (these relations are available at [56]). IBP relations from these seeds cannot do this job because 5 8 -type integrals do not show up. One can certainly explore other integral sets for each seed to find even better reduction efficiency. We did not do that because efficiency of either the IBP set (4) or the generalized set (5) is sufficient for us to deal with the problem in this work.
With integral sets in hand, we generate a system of linear equations from all seeds belonging to m n -type with 3 ≤ n ≤ 8 and 0 ≤ m ≤ 5, and use the package FiniteFlow [35] to trim the system by removing redundant relations and solve the trimmed system numerically, which expresses all integrals in S (a) as linear combinations of 108 MIs (after exploring symmetries between MIs).
Reduction at the step two. -In this step, we want to search for a system of block-triangle linear relations to reduce integrals in S (a) . Our search strategy is similar to that in Ref. [52], but using input information from the step one.
We first explain how to search linear relations among a given integral set G = {I 1 , . . . , I N } ⊆ S (a) . In general, a linear relation among G can be expressed as where Q i ( , s ) can be decomposed as where max is the maximal power of allowed to appear in the relation, Ω di = { λ ∈ N r | λ 1 + · · · + λ r = d i }, d i is half of the mass dimension of Q i which can be fixed by d max ≡ max{d 1 , · · · , d N }, andQ κλ1...λ5 i are unknown rational numbers to be determined. It is crucial to point out that, for given max and d max , the number of unknowns are finite. Therefore, they can be determined by finite number of constraints.
To obtain constraints over the unknowns for given max and d max , we reduce integrals in Eq. (6) to the 108 MIs numerically at a randomly chosen phase space point through the step one. Because MIs are independent of each others, their coefficients in Eq. (6) must vanish, which results in (at most) 108 linear constraints. Repeat this for sufficient number (several thousand in this work) of phase space points, we can always generate enough constraints to determine all the unknowns. As the above unknowns are actually evaluated in the finite field of a given prime number, we need to repeat the procedure for several different prime numbers (at most 15 in this work) and use the Chinese remainder theorem to reconstruct the true results.
If we want to reduce a set of integrals G 1 to another set of simpler 2 integrals G 2 (the reducibility can be tested numerically easily), we set G = G 1 ∪ G 2 and search relations among G with different values of max and d max . For the purpose of the current work, we find it is sufficient to fix max = 3. In order to find out simple relations, we follow the algorithm proposed in Ref. [52] by starting the search procedure with d max = 0 and increasing d max by 1 each time, until enough relations are obtained to reduce G 1 to G 2 . Now we use the above method to search for blocktriangle reduction relations for S (a) . We begin with the reduction of the most complicated 5 8 -type integrals. To this end, we search relations in the set G = S (a) with G 1 chosen as all 21 5 8 -type integrals. We indeed find out 21 independent relations, which can reduce all 5 8type integrals to simpler integrals. The most complicated relation corresponds to d max = 7, which means that the coefficients of 5 8 -type integrals are degree-2 polynomials in s. We then reduce 4 8 -type integrals, which can be realized by setting G =4 I {1,1,1,1,1,1,1,1,0,0,0} with G 1 chosen as all 15 4 8 -type integrals. To reduce the rest of top-sector integrals, we set G =3 I {1,1,1,1,1,1,1,1,0,0,0} with G 1 chosen as 11 top-sector integrals that are not MIs.
Similarly, we can reduce all integrals in subsectors. For example, by setting G =4 I {0,1,1,1,1,1,1,1,0,0,0} with G 1 chosen as all 35 4 7 -type integrals in G, we realize the reduction of 4 7 -type integrals in this sector. Based on the above scheme, we obtain 3801 reduction relations. By introducing additional 5 symmetry relations among MIs [48], we have 3806 relations in total that can express all the 3914 integrals in S (a) as linear combinations of 108 MIs. In Fig. 2, we show a matrix density plot of these relations with black points representing nonzero elements in the matrix. In this plot, each line represents a relation, and each column corresponds to an integral. These relations are ordered from top sector to subsectors, and integrals are ordered from the most complicated one to the simplest one, with MIs put at the end of the system. It can be found that the matrix is in an exact block-triangle form with the largest size of the blocks being 35 × 35. We note that there is a way to further reduce the block size. For example, by setting G =3 I {0,1,1,1,1,1,1,1,−1,0,0} , we can generate a smaller size of block to reduce part of 4 7 -type integrals.
Reduction relations for other families have similar properties. As shown in Tab. I, the file size of all relations is about 148 MB. To obtain them, it costs about 200 CPU-core hours in the search procedure, in addition to hundreds of CPU-core hours to generate input data. Analytic expression of all these relations are available at [56].
Checking our results and comparison with other methods. -Our final reduction relations have been verified numerically by an independent code FIRE6 [29] for randomly chosen phase space points, and the results agree with each other.
For each given phase space point, numerically solving our reduction relations of the four families totally costs 0.4 second using one CPU. The time spent can be divided into two parts: assignment of these relations for the given phase space point, which is proportional to the file size of these relations; and solving the linear system, which depends on both the number of relations and how these relations are coupled with each other. Because our system has a nice block-triangle form, the time spent in the latter part is shorter than the former. Therefore, the efficiency of numerical calculation of our reduction relations can be simply estimated by the file size.
Comparing with explicit solutions, the file size of our reduction relations is much smaller. The file size of explicit solutions of top-sector integrals with degree up to 4 in family (a), 26 integrals in total, is about 2GB [48]; that of top-sector integrals with degree up to 4 in family (b), 32 integrals in total, is about 0.8GB [47]; and that of all integrals in family (c) is in excess of 20GB for compressed format [49]. It can be expected that our relations should be hundreds of times smaller than the complete explicit solution in file size, which results in at least hundreds of times faster in numerical calculation, even if there is no problem for memory to store the huge expression of explicit solutions.
We note that the file size of trimmed IBP relations to reduce all integrals considered in this work is a few GB, which is also much larger than that of our reduction relations. The reason is that, although each IBP relation is simpler than ours, the IBP system involves hundreds of times more equations. Furthermore, the time spent for numerical IBP is dominated by the latter part because IBP relations are coupled in a complicated way. As a result, numerical IBP should be much more inefficient than our method. Through our test, numerical IBP via FiniteFlow [35] combined with LiteRed [34] costs about 2 minutes for each phase space point, which is slower than our method by hundreds of times.
The above comparison reveals the advantage of our method. Numerical evaluation of explicit solutions spends too much time on assignment; while numerical IBP spends too much time to solve linear equations. Our method has improved both parts, and therefore it is much more efficient. Similar to numerical evaluation over field of a prime number, our reduction relations should also be much more efficient for numerical evaluation with floating numbers, which enables phase space integration to obtain physical cross sections.
Summary and outlook. -In this Letter, we achieve the reduction of a set of integrals which covers all the most complicated integrals in two-loop fivelight-parton scattering amplitudes. Our results are expressed as a system of linear relations in block-triangle form, which are very efficient for numerical calculation. The remained integrals involved in amplitudes can be easily reduced using the same method upon demanding. Therefore, complete reduction of integrals in two-loop five-light-parton scattering amplitudes, which challenges all other methods, is available now. As MIs are already known [20][21][22][23][24], our results make the complete calculation of two-loop five-light-parton scattering amplitudes, and thus complete NNLO calculation of three light particles or jets production at the LHC on the horizon.
To obtain the block-triangle relations, we develop and improve our previous method [52]. As our newly developed method is general and efficient, other more complicated problems, like two-loop integrals for tt + jet, ttH, or 4jets hadproduction, are also within reach. Our work opens the door for complete NNLO QCD calculation for three or more particles production at the LHC.
In the current application of our method, most CPU time is cost to generate input data. Although the time spent is tolerable for the current problem, improvement may be needed for more complicated applications. There are different possible ways. Within the method of [52], one needs to explor better integral sets. Another possible choice is to use trimmed IBP systems obtained by solving syzygy equations [44][45][46][47][48]. We leave this for future study.