A planar four-loop form factor and cusp anomalous dimension in QCD

We compute the fermionic contribution to the photon-quark form factor to four-loop order in QCD in the planar limit in analytic form. From the divergent part of the latter the cusp and collinear anomalous dimensions are extracted. Results are also presented for the finite contribution. We briefly describe our method to compute all planar master integrals at four-loop order.


Introduction
One of the important tasks of modern high-energy particle physics is the development of new methods to compute quantum corrections to physical cross sections. This is particularly important in the context of Quantum Chromodynamics (QCD) where higher order corrections often have a significant numerical impact. In this article we provide the first next-to-next-to-next-to-next-to-leading order (N 4 LO) contribution to a three-point function within QCD. To be precise, we consider the photon-quark form factor, which is a building block for N 4 LO cross sections. Namely, it is a gauge-invariant part of virtual forth-order corrections for the process e + e − → 2 jets, or for Drell-Yan production at hadron colliders.
Denoting the photon-quark vertex function by Γ µ q the scalar form factor is obtained via where D = 4 − 2ǫ is the space-time dimension, q = p 1 + p 2 and p 1 (p 2 ) is the incoming (anti-)quark momentum. We consider the large-N c expansion of F q (q 2 ). As a consequence we only have to consider the contributions of planar Feynman diagrams. Results for F q can be used to probe the infrared structure of gauge theories. Form factors encapsulate universal infrared contributions coming from soft exchanges between two partons. The general form of the latter is known [1][2][3][4][5][6] and depends on cusp and collinear anomalous dimensions.
Two-loop corrections to F q have been computed more than 25 years ago [7][8][9][10]. The first three-loop result has been presented in ref. [11] and has later been confirmed in ref. [12]. Analytic results for the three-loop form factor integrals were presented in ref. [13]. In ref. [14], the results of ref. [13] have been used to compute F q at three loops up to order ǫ 2 , i.e., transcendental weight eight, as a preparation for the four-loop calculation. In this paper we compute the fermionic corrections to F q in the large-N c limit, to the four-loop order. Sample Feynman diagrams which have to be computed for this purpose are shown in figure 1.

JHEP05(2016)066
Over the last decades, powerful methods for determining loop integrands based on generalized unitarity have become common. However, form factors are simple enough that a direct Feynman diagram approach for determining the loop integrand is perfectly possible. The expression for the integrals contributing to the form factors is then reduced to a set of so-called master integrals, exploiting integration-by-parts identities [15]. This requires rather involved computer algebra, and can be achieved using the latest version of the program FIRE [16][17][18].
This leaves the evaluation of the master integrals as the main technical difficulty. In a previous paper [19], three of the present authors proposed a new technique for computing such integrals. Massless form-factor integrals have a trivial scale dependence, so the powerful method of differential equations [20][21][22][23][24] cannot be used directly. Rather, one first introduces an auxiliary parameter (corresponding to a second off-shell external leg), in which differential equations are set up. The main idea of [19] is that the boundary value of the differential equations can be fixed trivially from a value of the new parameter that corresponds to propagator-type integrals. This boundary value is then related to the original problem via the differential equations. This last step is especially easy in the canonical form [24] of the differential equations. It was suggested in that paper that in order to reach the canonical form it is helpful to select basis integrals that have constant leading singularities [25]. The latter are essentially multidimensional residues of the loop integrand and can be computed algorithmically. This connection makes it easy to reach the canonical form of the differential equations, as was demonstrated in many recent papers.
We classified all massless planar four-loop form-factor integrals and determined the corresponding master integrals. We found a total of 99 master integrals. We then computed them as described in the previous two paragraphs. A very welcome by-product of the approach of [24] is that the results are typically expressed in terms of uniform weight functions. Examples of uniform weight form factor integrals were previously considered in JHEP05(2016)066 refs. [13,26]. Here we systematically found a uniform weight basis for all planar integrals, and expanded them to weight eight. While this is the weight needed for typical four-loop computations, it is also possible to expand our result to higher weight.
Other attempts to calculate similar form factors or master integrals were reported on in refs. [27][28][29]. The evaluation of the master integrals in refs. [27,28] was performed only by numerical methods while ref. [29] presents results only for some individual integrals in an analytical form.
The remainder of the paper is structured as follows: in the next section we briefly outline our calculation and present results for the form factor and for the cusp and collinear anomalous dimensions. Sections 3 and 4 are dedicated to the classification and evaluation of the master integrals. Our conclusions are contained in section 5.

Calculation and results
We generate the Feynman amplitudes with the help of qgraf [30] and transform the output to FORM [31, 32] notation using q2e and exp [33,34]. For the reduction to master integrals we use the program FIRE [16][17][18] which we apply in combination with LiteRed [35,36]. Relations between primary master integrals occurring in the reduction tables are revealed with the help of tsort, which is part of the latest FIRE version [18], and based on ideas presented in ref. [17]. This leads to 78 master integrals needed for the fermionic part. More generally, we find that a total of 99 master integrals are sufficient for arbitrary planar integrals. They are all computed as described in sections 3 and 4.
In our calculation we allow for a generic QCD gauge parameter ξ and expand the Feynman diagrams around ξ = 0, which corresponds to Feynman gauge, up to linear order. We checked that ξ drops out before inserting explicit results for the master integrals.
In the following we present results for the form factor F q and the related anomalous dimensions. F q is conveniently shown in term of the bare strong coupling constant. In that case the perturbative expansion of F q can be cast in the form Analytic results for F (n) q , with n ≤ 3, expanded in ǫ up to transcendental weight eight can be found in ref. [14]. We refrain from repeating them here.
The main result of this letter is the fermionic contribution to F (4) q in the large-N c limit. It is given by where the ellipses stand for n f -independent contributions. The cusp and collinear anomalous dimension is conveniently extracted from log(F q ) (after renormalization of α s ). The pole part of the latter has the generic structure (see, e.g., refs. [12,37]) where µ 2 = −q 2 has been chosen and the coefficients of the β function are given by The coefficients of the cusp and collinear anomalous dimensions are defined through From eq. (2.3) it is evident that γ cusp can be extracted from the coefficient of the quadratic, and γ q from the first-order pole in ǫ. In the large-N c limit we obtain for γ cusp where the ellipses in γ 3 cusp indicate non-n f terms which are not yet known. For γ q we have For the finite part of log(F q ) we obtain The expressions in eqs. (2.6) and (2.7) up to three-loop order confirm the results in the literature [11,12,[37][38][39][40][41] and the N 3 c n 3 f term of γ 3 cusp agrees with the result of refs. [42,43]. All other terms in the four-loop results γ 3 cusp and γ 3 q and the finite part in eq. (2.8) are new.

Integrals with constant leading singularities
Our calculation involves planar four-loop form-factor integrals. We classified all such integrals and performed an integral reduction, resulting in 99 master integrals. Before discussing their evaluation, we devote this section to our basis choice for these integrals.

Leading singularities and d-log forms
In recent years it has become standard to use a basis, whenever possible, of integrals having constant leading singularities. Leading singularities [25] are essentially defined as multidimensional residues of the Feyman loop integrand.

JHEP05(2016)066
The usefulness of integrals with constant leading singularities was first noticed in the context of maximally supersymmetric gauge theory, where the answer appears to be naturally written in terms of them. Building on experience with such integrals in the literature, their systematic use was advocated in ref. [44]. A particular highlight is an all-n expression, where n is the number of external legs for the integrand of two-loop maximally helicity violating amplitudes in N = 4 super Yang-Mills theory. In fact, it turns out that the appearance of integrals with simple leading singularities in this theory is very natural, as can be seen in the twistor approach of [45], or when expressing leading singularities as certain Grassmannian integrals [46]. Although more examples are known in the planar case, the concept of constant leading singularities also carries over to the non-planar sector, see [26,47,48] for examples.
The use of such integrals is not limited to supersymmetric amplitudes, as was pointed out in ref. [24]. Since then, they were applied to countless calculations of scattering amplitudes required for phenomenology, see, e.g., ref. [49]. Of course, more integrals are needed in QCD compared to supersymmetric theories. In this context, it is perhaps interesting to point out that many of the additional integrals needed for QCD can be thought of as integrals being defined in a dimension shifted by two units. As is well known, integrals in D ± 2 and D dimensions are related. The picture that emerges is that one should not only classify integrals having constant leading singularities in four dimensions, but in all integer (in particular even) dimensions, and then relate them to the four-dimensional case.
Let us give some one-loop examples of such integrals. We define the triangle integral near four dimensions and the propagator-type integral near two dimensions, where p 2 1 = p 2 2 = 0. In the following we will consider leading singularities at ǫ = 0. It is convenient to change variables. For the bubble, we set k µ = αp µ 1 + βp µ 2 , which leads to where the proportionality sign means that the equation holds up to kinematic-independent factors. While there are various locations of the leading singularities, we can see that all poles are kinematic-independent. A similar analysis was done for the triangle integral, see [48]. We mention that one can rewrite the integrands (algebraically) in a form where this property is manifest, namely, dα dβ (α + 1)αβ(β − 1) Figure 2. Twelve-propagator form-factor integral that has unit leading singularities. The numerator (k 2 4 ) 2 normalization factor is implied.
Here k ± denotes the two solutions to the maximal cut condition, (k ± + p 1 ) 2 = 0, (k ± − p 2 ) 2 = 0, which are given by k + = −p 1 + p 2 and k − = 0. Equation (3.4) implies that there exist variables in which the integrand is just d log x 1 d log x 2 , with unit normalization. More formulas of this type, called d-log forms, can be found in refs. [47,48].
Following these ideas, we wrote down a basis of integrals with constant leading singularities for planar four-loop form-factor integrals. The whole basis will be presented elsewhere. Here we give one example, for the twelve-propagator integral that was needed in the n f -calculation. See figure 2 and the first diagram of figure 1 for a representative Feynman diagram. We choose as basis element The normalization factors were chosen for later convenience. We first would like to illustrate that this integral has indeed constant leading singularities. While this can be done algorithmically, it is instructive to rewrite the integrand in a form where this is obvious, namely in terms of d-log forms of the type discussed above. A very useful feature is that this analysis can be done loop by loop, which allows one to recycle formulas. This is very similar to an analysis via cuts, although here we do not assume that any loop momenta are on-shell. First, we note that the box subintegrals with three off-shell legs, i.e. the ones depending on loop momenta k 1 and k 3 (see figure 2), can be written in a d-log form. For the subsequent calculation, only the normalization factor of these subintegrals is relevant. The latter can be obtained from any of their leading singularities. For example, for the box integral on the left, we have the following integrand (3.6)

JHEP05(2016)066
After taking a multi-dimensional residue one obtains either zero, or a term proportional to For the box on the right a similar expression is obtained. Next, we consider the k 2 integration. Taking into account the factors obtained from the k 1 and k 3 integrals, we arrive at a generalized box integral . (3.8) Again, it can be shown that this has a d-log form, with the normalization factor 1/(k 2 4 ) 3 . We refer the interested reader to [44,48] for more detailed examples of leading singularity and d-log calculations. We now see that the numerator in eq. (3.5) cancels the excessive factors of k 2 4 . Indeed, putting everything together, we see that the remaining k 4 integral is exactly of the form of the one-loop triangle integral of eq. (3.1). In summary, this proves that (3.5) has a d-log representation with unit normalization.
We would like to emphasize again that the classification of integrals having constant leading singularities can be done algorithmically. Let us expand on this point. First of all, for a given propagator structure, one makes an ansatz for all possible numerator terms allowed by power counting (or subject to other criteria). It is convenient to parametrize the loop momentum in such a way that the integration parameters are scalars. We illustrated this in the case of the bubble integral, cf. eq. (3.3). Next, one evaluates all leading singularities of this ansatz (i.e., one computes the residues at all poles of the integrand). Requiring that the residues be kinematic-independent yields a system of equations, which is then solved. It is important to realize that this analysis only depends on the integrand at hand, and can be done before attempting to compute the integral.

Transcendental weight properties
One nice feature of integrals with constant leading singularities is that, conjecturally, they evaluate to so-called pure functions, i.e. iterated integrals of uniform weight.
For iterated integrals, such as multiple polylogarithms, the weight is defined as the number of integrations, e.g. one for logarithms, n for classical polylogarithms Li n , etc. Similarly, transcendental constants such as zeta values, ζ n , have weight n. Finally, when considering Laurent expansions in the dimensional regularization parameter ǫ, one can assign weight −1 to ǫ. This is natural since 1/ǫ would be represented by logarithm in a cutoff regularization.
With these definitions, we see that the triangle integral of eq. (3.1) has uniform weight 2, More generally, L-loop integrals with constant leading singularities in 4 dimensions are expected to evaluate to weight 2L functions.

JHEP05(2016)066
Beyond maximally supersymmetric Yang-Mills theory, also functions of weight smaller than 2L are needed. Perhaps the best way to understand the additional integrals is to consider them in different dimensions. Take as an example the 2 − 2ǫ dimensional bubble integral of eq. (3.2). In fact, in this simple example, the bubble and triangle integrals are related by an integration-by-parts relation (and, dimensional shift relation [50,51]), which implies where the integrals are defined in eqs. (3.1) and (3.2). From this formula we see that its weight is shifted by one compared to the triangle. It evaluates to a uniform weight-one function.
More generally, at higher loops one can also generate integrals of various weights, in particular by writing subintegrals formally in different dimensions. For example, all the uniform weight integrals presented in ref. [24] and elsewhere can be understood in this way. See also the lecture notes [52] for more details.
Returning to our form-factor integrals, we can verify the uniform weight conjecture for the most complicated twelve-propagator integral of eq. (3.5). It turns out that the uniform weight property can be best understood systematically using differential equations that are discussed in the next section. Here we anticipate a result of that calculation, = ζ 8 + ζ 5,3 = 1.041785 . . . and ζ 5,3 is a multiple zeta value [53]. Reinstating the 1/ǫ 8 from eq. (3.5), one sees that this is a uniform weight eight integral, as expected from a four-loop integral. As an independent check of eq. (3.11), we derived a Mellin-Barnes representation (see Chapter 5 of [54] for a review) for this integral, which we used to verify the first three terms in the ǫ expansion analytically.

Differential equation bootstrap for single-scale integrals
In this section we discuss the first analytic computation of all planar four-loop on-shell formfactor integrals which are defined in the kinematic regime p 2 1 = p 2 2 = 0, with q 2 ≡ p 2 3 = (p 1 + p 2 ) 2 . Following the strategy of [19] we introduce an auxiliary parameter by taking a second external leg off-shell, i.e. p 2 2 = 0, and x = p 2 2 /p 2 3 , and derive differential equations with respect to x. The main idea can be explained via figure 3. It turns out that the singular points of the differential equations are x = 0, 1, ∞. The point x = 0 corresponds to the original on-shell form-factor integrals, figure 3(b). Similarly, for x = ∞ we have p 2 3 = 0, which again leads to form-factor integrals, as for x = 0, and thus this case does JHEP05(2016)066 not have to be considered separately. On the other hand, the boundary value at x = 1 corresponds to propagator-type integrals, see figure 3(c). They can be determined easily: in most cases, the boundary value is zero due to kinematical factors. Otherwise one can use results for propagator type integrals available in the literature, see, in particular, four-loop analytic results in [55][56][57]. This information is then transported back via the differential equation to x = 0, see figure 3(a). Let us now see how this works in a bit more detail. A pedagogical example is given in [52].
In reference [24], a canonical form of differential equations for Feynman integrals was suggested. Conjecturally, this form can be reached whenever the master integrals can be chosen to have the property that their leading singularities are constant, as explained in section 3. This reduces the problem of finding a canonical basis for the differential equations to a simple classification of integrals having this property. The latter can be done algorithmically.
For the planar form factor with p 2 2 = 0 and p 2 3 = 0 we find a total of 504 master integrals (some of them related by symmetry). After choosing a canonical basis f , we found the following system of differential equations, where a and b are some constant (i.e. x-and ǫ-independent) 504 × 504 matrices. The special features of this form are the manifest Fuchsian property of the singularities, i.e. only single poles in x = 0, 1, ∞ are present on the right-hand side of eq. (4.1), and the fact that the right-hand side is proportional to ǫ. The latter property can be achieved for iterated integrals. Here, it implies that the solution, to any order in ǫ, can be written in terms of iterated integrals over the kernels dx/x and dx/(x − 1), i.e. in terms of harmonic polylogarithms [58]. The former property is true for any Feynman integral. Making it manifest allows us to describe the boundary behavior in a simple way, namely f (x, ǫ)

JHEP05(2016)066
where a and b are the matrices from eq. (4.1) and the coefficients matrices p k (ǫ) and q k (ǫ) in the expansion can be obtained recursively [59]. We fix the boundary value at x = 1 by demanding regularity of the integrals in this limit and using explicit results for some propagator type integrals. This determines f 1 (ǫ).
We then use the differential equation (4.1) to transport this boundary value back to x = 0. (In mathematical language, we construct the Drinfeld associator, perturbatively in ǫ.) This allows us to determine f 0 (ǫ). Finally, unlike the x → 1 limit, the x → 0 limit is singular, in the sense that the matrix exponential x ǫa contains several terms x ǫα , with α = 0. These non-zero values of α correspond to contributions of various regions [60][61][62] to the asymptotic expansion in the given limit. The on-shell integrals we would like to compute correspond to the so-called "hard" region with α = 0.
In order to determine to the on-shell integrals, we reduce the basis f for on-shell kinematics, expressing it in terms of 99 on-shell master integrals. We then match the expression so obtained to the hard region at x = 0. We find that this determines all the 99 integrals (naturally, some of the 504 equations are redundant). In order to carry out these algebraic manipulations, we found the Mathematica implementation HPL.m [53] useful.
In summary, we analytically computed all planar form-factor integrals with two offshell legs (504 master integrals), and with one off-shell leg (99 master integrals). The result for the most complicated on-shell form-factor integral with twelve propagators that appeared in the n f -piece of our calculation was given in eq. (3.11) as an example. The full analytic results for all planar master integrals will be given elsewhere.

Conclusions
In this paper we report the calculation of all 504 master integrals which are needed for a generic planar massless form factor with two off-shell legs. They are obtained by a proper choice of basis integrals, together with boundary conditions where the form factor degenerates to a two-point function. From the generic basis we derive analytic results for the 99 master integrals that are needed for the planar on-shell form factor. 78 out of the 99 master integrals are needed for the fermionic part of the planar photon-quark form factor. The latter is considered in section 2 of this paper, where analytic results up to four loops are presented for the cusp and collinear anomalous dimension and the finite part of F q .
A natural extension of this work is to apply the planar master integrals we computed to evaluate the non-fermionic planar contribution, where the integral reduction is more complicated. Furthermore, we expect that the methods discussed in this paper can also be applied to non-planar form factor integrals.