Di-photon amplitudes in three-loop Quantum Chromodynamics

We consider the three-loop scattering amplitudes for the production of a pair of photons in quark-antiquark annihilation in Quantum Chromodynamics (QCD). We use suitably defined projectors to efficiently calculate all helicity amplitudes. We obtain relatively compact analytic results, that we write in terms of harmonic polylogarithms or, alternatively, multiple polylogarithms of up to depth three. This is the first calculation of a three-loop four-point scattering amplitude in full QCD.

Multiloop scattering amplitudes in Quantum Chromodynamics (QCD) are a crucial ingredient for the precision physics program carried out at particle colliders, such as the Large Hadron Collider (LHC) at CERN. When combined with the corresponding real radiation, they make it possible to compute high-precision predictions for a multitude of key processes. In the recent past, such calculations played a fundamental role for Higgs characterisation [1], for the extraction of fundamental parameters of the SM [2,3] and for the study of electroweak bosons [4][5][6] and the top quark [7], to mention a few examples. These investigations allow for a deep scrutiny of the Standard Model (SM), which is essential for establishing its validity and for revealing possible tensions pointing towards new physics.
Together with their importance for particle collider phenomenology, scattering amplitudes in QCD also provide us with an important laboratory to investigate formal properties of the perturbative expansion of realistic Quantum Field Theories (QFTs). With the increase of the perturbative order, one faces the computation of increasingly complicated multiloop Feynman diagrams which can depend on multiple scales according to the number of particles involved in the collision. This translates into increasingly involved analytic properties of the special functions required to express the corresponding scattering amplitudes in closed analytic form. The calculation of multi-loop scattering amplitudes have contributed to reveal intricate and fascinating mathematical structures, whose further investigation has the potential to improve our understanding of perturbative QFT.
In parallel to the effort to overcome the barrier of twoloop scattering amplitudes for 2 → 3 processes, similar work is required to push the calculation of 2 → 2 QCD processes up to three loops. Indeed, while the first 2 → 2 three-loop results have appeared in N = 4 Super Yang-Mills [41] and in N = 8 Super Gravity [42], three loop calculations in realistic QFTs have only been performed for simple 2 → 1 processes [43][44][45]. Only recently some of the ingredients for obtaining three-loop 2 → 2 amplitudes in massless QCD have started to appear [46,47].
In this letter, we document the calculation of the threeloop QCD corrections to the scattering amplitude for di-photon production in quark-antiquark annihilations. This is the first three-loop amplitude for a 2 → 2 scattering process in full-colour QCD. The qq → γγ process is arguably a very natural place to start the investigation of three-loop corrections to 2 → 2 scattering in a realistic QFT. Indeed, these amplitudes involve only massless external and virtual particles and, moreover, have a much simpler colour structure than the corresponding amplitudes for the production of coloured partons. The relevant one-loop and two-loop amplitudes have been known for a long time [10,12], followed by a multitude of phenomenological studies up to NNLO in QCD [48][49][50]. However, we would like to stress that despite the aforementioned simplifications, this process still contains all of the analytic complexity of a generic massless 2 → 2 scattering process. Hence, we expect that the results obtained here can be extended to compute all three-loop scattering amplitudes for the production of two massless partons in hadronic collisions.
We now describe our calculation. We consider the process and obtain its physical equivalent qq → γγ via crossing symmetry p 3,4 → −p 3,4 . We parametrise the kinematics of Eq. (1) by defining the usual Mandelstam invariants We also introduce the dimensionless variable x = −t/s, arXiv:2011.13946v1 [hep-ph] 27 Nov 2020 such that in the physical scattering region, one has We use the helicity of the incoming quark λ q and of the incoming photons λ 3,4 to label the scattering amplitude of the process. We denote the scattering amplitude between well-defined helicity states by A λqλ3λ4 . Using parity, charge-conjugation and symmetry relations, it is easy to see that there are only two independent helicity configurations [10,12]. In what follows, for clarity we will compute the over-complete set of four helicity configurations {λ q λ 3 λ 4 } = {L − −}, {L − +}, {L + −}, {L + +}, which allow us to obtain the remaining ones for righthanded quarks by a simple charge-conjugation transformation, as it will be discussed below.
In order to compute the helicity amplitudes, we regulate both infrared and ultraviolet divergences using dimensional regularisation, i.e. we work in d = 4 − 2 dimensions. By choosing a gauge such that it is easy to see that, at any order in QCD perturbation theory, Lorentz covariance dictates that the amplitude can be parametrised as where u andū are the spinors for the incoming quark and antiquark, respectively, and the five Lorentz tensors Γ µν i are defined as The functions F i (s, t) are scalar form factors which only depend on the Mandelstam invariants and on the number of dimensions d. Since the colour structure is straightforward, we keep colour indices implicit here. For ease of typing, we define the five independent structures and from now on, with a slight abuse of language, we refer to the T i as the independent tensors for the problem at hand. At first sight, it may seem puzzling that we find five independent tensor structure when we have only four helicity amplitudes (not considering charge conjugation). This mismatch is however easy to explain: the decomposition Eq. (5) is valid for arbitrary dimension d. For fourdimensional external states, the five tensor structures T i are no longer independent, and four of them are enough to span the whole space [51,52]. Since eventually we are only interested in the d → 4 limit, it is convenient to reorganise the tensors T i and choose for T 5 a linear combination which is identically zero when four-dimensional external states are considered. This can be achieved by choosing We then write the scattering amplitude as where the form factors F i are suitable linear combinations of the original F i , whose explicit form will be irrelevant in the following. All the non-trivial information for the process of Eq. (1) is encoded in the form factors F i . We stress once more that, while all five form factors are in general required for the result in Conventional Dimensional Regularisation (CDR), only the first four are enough to compute the helicity amplitudes in tHV, see the supplement material for more details. Eq. (9) can be inverted to select individual form factors. This is done by introducing suitable projectors operators defined such that 1 where we use d dimensional polarization sums. Since the five T i form a basis for our space, we can write the projectors as A straightforward calculation, reported in the supplement material, shows that the matrix c ik is blockdiagonal As a consequence, the fifth tensor T 5 decouples from the other four. This, combined with the fact that it always evaluates to zero for four-dimensional external states, shows that the helicity amplitudes receive contributions only from the first four form factors in Eq. (9). In other words, there is a one-to-one correspondence between helicity amplitudes and form factors in the 't Hooft-Veltman scheme (tHV) [53]. This statement is quite general and it is always possible to construct a basis where the decoupling of tensors leading to vanishing results in d = 4 is manifest. Details of how to achieve this in general will be presented elsewhere [52]. In order to compute the form factors, we generate the relevant tree-level, one-, two-and three-loop Feynman diagrams with QGRAF [54] and apply on each of them the the projectors defined in Eqs (11) and (12). We perform all required colour and Dirac algebra algebra with FORM [55]. There are 2 diagrams at tree level, 10 at one loop, 143 at two loops and 2922 at three loops. At a given number of loops, each Feynman diagram can be identified with a graph, which is specified by the permutation of the external legs and by the number and topology of the internal lines. In particular, at three loops, where the number of Feynman diagrams becomes rather large, it is convenient to sum together those diagrams that are associated to the same graph and manipulate them together to avoid as much as possible to perform the same operation multiple times. Performing these steps allows us to express the three-loop form factors in terms of a large number of scalar integrals of the form 2 with γ E ≈ 0.5772 the Euler constant. The propagators D j can be drawn from three different families of integrals fam = {PL, NPL1, NPL2} and their crossings. The definition of the three families is provided in electronic files attached to this letter. A key features of integrals of the form Eq. (13) is that they are not all independent. Instead, through a by-now standard use of integration-byparts identities, most of them can be expressed in terms of a much smaller set of so-called master integrals [56,57]. This can be done in principle algorithmically for any process [58] for integrals of the form Eq. (13). Despite being well-understood in principle [58], this reduction step is computationally highly non-trivial for the problem at hand, as it requires the symbolic solution of large systems of linear equations. We have achieved it through a combined use of the public code Reduze 2 [59,60], and an in-house implementation, Finred, employing finite field sampling [61][62][63][64] and syzygy techniques [65][66][67][68][69][70]. In this way, the three loop form factors F j (s, t) can be expressed as linear combinations of a total of 486 (132) independent master integrals, including (excluding) integrals related by crossings of the external legs.
A basis of master integrals, excluding integrals obtained from them by crossings of the external legs, has been computed for the first time in Ref. [47] using the method of differential equations [71][72][73][74][75] augmented by the choice of a canonical basis [76,77]. There it was shown that all these integrals can be expressed in terms of a particularly simple and well understood class of functions, harmonic polylogarithms (HPLs) [78][79][80][81][82]. The present calculation requires results through to transcendental weight six. 3 We have re-derived the differential equations fulfilled by the master integrals including all their non-trivial crossings required for the calculation of the actual scattering amplitudes. We have verified that the results of Ref. [47] do fulfil the differential equations. As a check, we also have recomputed various boundary constants for the integrals of Ref. [47] independently. For manipulating the harmonic polylogarithms, we used inhouse routines and PolyLogTools [83].
Extra care has to be put in computing the boundary conditions and in analytically continuing the integrals to the physical scattering region in Eq. (3). Indeed, it is well known that scattering amplitudes are multivalued complex functions of the external kinematics and that, in order to obtain physical results, one must define the integrals on the correct Riemann sheet whenever one or more of the Mandelstam invariants are positive. The standard approach consists, whenever possible, in computing the relevant scattering amplitude for non-physical values of the kinematics where all Mandelstam variables are negative and a real result is expected. If this is possible, it is then often simple to analytically continue the amplitudes to physical kinematics using the Feynman prescription s i → s i + i0 for each Mandelstam invariant s i which crosses a branch cut. Unfortunately, this approach fails for 2 → 2 massless scattering where the momentum conservation relation s + t + u = 0 prevents us from finding such a non-physical region. As a consequence, the scattering amplitudes must be computed directly in a region where at least one Mandelstam invariant is positive.
In the present case, extra care has to be taken since the boundary conditions provided in Ref. [47] are computed on the unphysical Riemann sheet reached by giving a small and negative imaginary part to the Mandelstam invariant s → s − i0. Of course, it is expected that these results should be related to the ones for s → s + i0 by complex conjugation. We have continued the results of Ref. [47] back to the physical Riemann sheet and checked that this holds. We then used them to obtain analytic expressions for all required crossings (see [84] for details). We used these integrals to compute the three loop form factors F The results computed with the procedure discussed so far contain both ultraviolet and infrared singularities. Up to three loops, they can be written as follows where e = √ 4πα is the electric charge, e q is the charge of the incoming quark in units of e, α s,b is the bare strong coupling and δ kl carries the colour indices of the two incoming quarks. We remove ultraviolet singularities by expressing our result in terms of the MS renormalised coupling α s (µ): The relation between renormalised and bare coupling is given by with S = (4π) − e −γ E and The explicit form of the β-function coefficient β 0,1 is reported in the supplemental material. For definiteness, we will present our results for µ 2 = s. It is straightforward to obtain results at any other scale using renormalisationgroup arguments. The renormalised form factors F (k) i of Eq. (15) still contain infrared singularities. Their form is universal and can be expressed in terms of the lower-order scattering amplitudes as follows [85,86] In these equations, the I j are universal factors that only depend on the center of mass energy of the coloured partons s, on and on the QCD Casimirs C F = 4/3, C A = 3, as well as on the number of light fermions n f . We spell them out explicitly in the supplemental material. For our discussion, it is only important to note that the function I i contains infrared poles up to order 2 i , i.e. I 1 ∼ 1/ 2 , (18) are finite in the → 0 limit, and contain all the non-trivial physics information for the process Eq. (1).
As we have already mentioned, it is straightforward to obtain the helicity amplitudes from our form factors by evaluating the tensor structures T i for well-defined helicity states, see Eq. (8). We stress once more that for any helicity configuration one has T 5,λqλ3λ4 = 0. We write for left-handed spinorsū L (p 2 ) = 2| and u L (p 1 ) = |1] and for the the photon j of momentum p j µ j,− (q j ) = q j |γ µ |j] √ 2 q j j , With these definitions we find [23] α(x) , A L−+ = 2 24 [13] 23 [24] β(x) , 24 [32] γ(x) , A L++ = 2 34 2 31 [23] δ(x) . with Note that the remaining amplitudes for right-handed quarks can be obtained by a charge-conjugation transformation as follows where λ * i indicates the opposite helicity of λ i . Symmetry under exchange of the two photons requires (22) and at tree-level we find α = δ = 0 and β = γ = 1.
Using the relations Eqs (19,20), we obtain the threeloop renormalised finite remainders α (3,fin) , β (3,fin) , γ (3,fin) , δ (3,fin) , defined in analogy to Eqs (15,18). This represents the first three-loop calculation of a full QCD amplitude. Our result can be expressed in terms of harmonic polylogarithms of weights 0 and 1, or, alternatively, in terms of a more compact functional basis consisting of 14 classical polylogarithms and the 9 functions which are two and three-fold nested sums in the conventions of [80]. 4 In the latter representation, our three-loop amplitudes can be evaluated in fractions of a second in one phase space point. Although the result is relatively compact, it is too long to be presented here. We attach it in electronic format to the arXiv submission of this letter. In Fig. 1 we show a numerical evaluation of our analytic result for all independent finite remainder functions through to three loops. While only the real parts are shown here, we note that the imaginary parts of the loop corrections can actually exceed their real parts substantially in the case of β.
Before concluding, we list some of the checks that we have performed to verify the correctness of our calculation. We have performed the tree-level, one-and twoloop computation from scratch using our setup, and obtained agreement for the finite reminders at O( 0 ) with results in the literature [10]. From Eq. (18) it follows, that to extract F (3,fin) i , one also requires the one and two loop results expanded up to order 4 and 2 , respectively. To check their correctness, we have abelianised the result of [46] and checked against the abelian part of our results up to weight six. As we have mentioned, Bose symmetry relates amplitudes of different helicities. We have verified that our helicity amplitude coefficients up to three loops display the Bose symmetry properties in Eq. (22).
A non-trivial aspect of our calculation is the analytic continuation of all the required integrals to the physical region. To check our procedure, we verified that our solutions are consistent with the boundary values for the master integrals that can be inferred by imposing regularity conditions on their differential equations. This allowed us to relate complicated four-point integrals to simpler two-and three-point ones (see [87] for a compilation), whose analytic continuation is straightforward. Moreover, we used Reduze 2 to find non-trivial symmetry relations among the master integrals and their crossings and verified that they are all satisfied by our analytic results. In addition, we checked the finite part of the two-loop integrals against the literature and some of the simple three-loop integrals against secdec [88] . Finally, the most powerful check of the correctness of our result is that the remainders F (3,fin) i are in fact finite in the → 0 limit. This also tests in a non-trivial way our analytic continuation procedure, as it links three-loop integrals with lower loop ones.
In conclusion, in this letter we have reported the first calculation of a three-loop four-point scattering amplitude in full QCD. We have obtained our results by defining a minimal set of projectors that allowed us to extract all the information required to reconstruct the helicity amplitudes. Our result can be expressed in terms of classical polylogarithms plus a handful of multiple polylogarithms, and it is very compact for a QCD amplitude of this type. The methods we employed for this calculation are generic and can be used to compute the three-loop helicity amplitudes for 2 → 2 scattering of massless partons in QCD. We leave this for future investigations. + 11π 4 90 while the one for the quark anomalous dimension reads