Direct Solution of Integration-by-Parts Systems

Systems of integration-by-parts identities play an important role in simplifying the higher-loop Feynman integrals that arise in quantum field theory. Solving these systems is equivalent to reducing integrals containing numerator products of irreducible invariants to a small set of master integrals. I present a new approach to solving these systems that finds direct reduction equations for numerator terms of a given Feynman integral. As a particular example of its power, I show how to obtain reduction equations for arbitrary powers of irreducible invariants, along with their solutions.


I. INTRODUCTION
The computation and simplification of Feynman integrals play a central role in the evaluation of higher-loop scattering amplitudes, form factors, and correlation functions in quantum field theory. In a frontier calculation, one must often consider a large number of integrals, which are nonetheless related by algebraic identities. Revealing the full set of algebraic relations between integrals reduces the number of integrals which have to be evaluated analytically or numerically. Knowing the full set of algebraic identities is crucial to the unitarity method [1,2] for computing scattering amplitudes beyond one loop [3][4][5][6][7], as well as to computing Feynman integrals using differential equations [8].
The integration-by-parts (IBP) approach within dimensional regularization [9] is currently the method of choice for obtaining such algebraic relations between different Feynman integrals. As applied to integrals beyond two points, the approach generates all possible total derivatives with increasing powers of numerator insertions, generating large systems of equations. One then uses Gaussian elimination, in the careful form introduced by Laporta [10] to solve the system of equations. A number of dedicated automated solvers [11] have been introduced and used over the years, complemented by alternative approaches [12].
Can one reduce the size of the systems, and also find a simpler method to solve them?
The first question was answered affirmatively by Gluza, Kajda, and the author [13], through the introduction of so-called generating vectors. These avoid introducing higher powers of propagators into the system of equations, terms which would later disappear during Gaussian elimination to solve the system. These generating vectors have links to algebraic geometry [14,15], and have seen further development [16] and applications [6,17] recently.
(Competing calculations have made use of a mix of algebraic geometry and more conventional tools [7].) An alternative approach to finding the vectors, less linked to algebraic geometry, may be found in ref. [18].
The goal of this paper is to address the second question, and outline an approach to solving IBP systems directly. As an example of the power of such an approach, I will show how to find closed-form expressions for arbitrary powers of numerator insertions, a question which is largely intractable with current methods.
I focus in this article on planar two-loop integrals, and mostly on the two-loop planar double box with massless external legs. This integral is simple enough to display many formulae explicitly, but nontrivial enough to put the approach to the test. The approach is of course applicable much more generally, to integrals with external or internal masses, and to higher loops as well. In the next section, I review two-loop Feynman integrals, the IBP approach, and generating vectors. In section III, I present a pair of challenges which the new method can address. In section IV, I show how to target simple numerators directly.
Section V is devoted to a basic approach to finding master integrals within the present approach. In section VI, I show how to target numerators with generic powers of irreducible invariants. Section VII discusses higher powers of propagators. In section VIII, I show how to solve the kinds of equations derived in Sect. VI. I present a few concluding remarks in Sect. IX.

II. INTEGRALS, INTEGRATION-BY-PARTS, AND GENERATING VECTORS
Let us consider a Feynman integral with two or more loops in dimensional regularization, , (2.1) where L is the number of loops, and n d the number of denominator factors. The generic numerator expression is given in terms of dot products of loop momenta i with each other or with the external momenta k i . All integrals with numerators containing dot products of loop momenta with arbitrary external vectors can ultimately be expressed in terms of integrals in eq. (2.1), so they suffice to express the result of any L-loop Feynman diagram, and hence any L-loop amplitude or form factor.
The standard IBP approach proceeds by forming a sufficient number of total derivatives, Denom({ j 5 , k j 6 }) , (2.2) where v µ is taken in turn to be any loop momentum or independent external momentum, in order to close the system of equations. The system will close, as discussed in ref. [10], when one considers polynomials of sufficiently high order, along with all subtopologies where one or more propagators are omitted.
One can instead seek special vectors v µ j such that [13], Given a set of vectors, an infinite tower of IBP equations can be generated by multiplying them by polynomials in Lorentz invariants of the loop momenta, In a certain sense, the use of the vectors v j block-diagonalizes the IBP system. It does not completely solve the system, however, in that we still have to generate multiple equations and solve them together in order to reduce a generic term in an integrand. To discuss the details of reduction, it will be convenient to recall some classes of Lorentz invariants from the literature, and to introduce some additional specializations.
We are interested in 'natural' Lorentz invariants, products of the loop momenta with other loop momenta or the external momenta of the integral. Any Lorentz invariant which can be written as linear combinations of propagator denominators and invariants built out of external momenta is called a 'reducible invariant' or RI. (In the literature these are often called reducible scalar products, however we wish to consider quantities which may not be simply scalar products.) Invariants which can be written purely in terms of propagator denominators, without use of invariants in external momenta, we will denote 'pure reducible invariants' or PRIs. Invariants which cannot be written as a linear combination of propagator denominators and external invariants are called 'irreducible invariants' or IrIs. They first arise at two loops, and play a central role in IBP systems 1 .
The coefficients of terms in the numerator polynomials in eq. (2.5) are all rational functions of = (4 − D)/2 and ratios of external invariants, which we can treat as parameters.
Terms with factors of PRIs reduce to integrals with fewer propagators, that is corresponding to simpler topologies. In this article, I will discuss only the first stage of reduction, and so will set aside such terms. Of course, one can and must deal with the resulting simpler topologies to obtain a complete reduction to a basis of integrals.
Any term in the numerator polynomial in eq. (2.5) that contains a PRI yields nothing interesting for the top-level topology, as it merely cancels against a linear combination of denominators. Accordingly we can take a generic term, without loss of generality, to be a 1 Leaving aside parity-odd terms at one loop.
product of powers of IrIs, in the sense that they have smaller | n| = j n j . Let us call this value the irreducible degree or i -degree for short. (If we need to distinguish between monomials of the same i -degree, we can use any monomial ordering employed in computational algebraic geometry -for example, lexicographic -to determine which is 'simplest'.) A pair of challenges illustrates the power of such an approach. Consider the slashed-box integral, shown in Fig. 1, following the notation of ref. [13]. In this expression, the external momenta k 1···4 are all massless and directed outwards. The first challenge is to simplify, for a generic integer value of n.
Consider also the planar double-box integral, shown in Fig. 2, where the notation again follows ref. [13], and where, The second challenge is to simplify, for a generic integer value of n.

IV. TARGETED REDUCTIONS
Let us begin by studying simple numerators in the double box. There are three linearly independent pairs of IBP-generating vectors when all external legs are massless. All generating vectors can be written as linear combinations, with coefficients taken to be general polynomials in the Lorentz invariants. The first two pairs will suffice for our initial purposes.
The first pair is, where the second index corresponds to the loop-momentum index. (These expressions differ from those in ref. [13], but are equivalent as solutions to eq. (2.3).) In the notation of It will be convenient to introduce a basis of RIs and IrIs and a short-hand notation, for the PRIs, and for the IrIs. With these variables, we can rewrite the first pair of vectors as follows, The second pair is,  where χ 14 = s 14 /s 12 .
If one forms the corresponding differential operators to the vectors, the third vector (given in appendix A) is related to the commutator of the first two 2 , (4.9) Our first task would be to determine the master integrals. We will return to this question in Sect. V; for the moment, let's assume we've already done this. We could in principle do this by generating IBP equations using numerator polynomials of increasing engineering dimension, starting with constants, and solving the equations until the number of independent integrals stabilizes. In the case of the double box, we can choose our masters to be, (4.10) We will consider eq. (2.5) using a variety of polynomials with the two pairs of vectors given in eqs. (4.5) and (4.6). We can expand the integrand in eq. (2.5) and multiply by the denominator to obtain an expression for the polynomial-dependent part of the numerator, (4.11) In this equation, n v is the number of sets or tuples of IBP-generating vectors. The first term in the equation is independent of the derivative of the numerator polynomial, and hence has 2 I thank Harald Ita for pointing this out. a universal structure. We can record the values of the coefficients for the two vector pairs, (4.12) The simplest IBP we can get comes from using the first vector pair (4.5), which allows us to solve for P * * 2,2 [t 14 ] in terms of the masters and integrals of simpler topology. In the present case, the simpler integrals cancel after using their symmetries.
If we look at the next simplest equation, multiplying the first vector pair by t 14 , we obtain, which has two terms of i -degree two, and hence we need a pair of equations to solve for both quadratic powers of IrIs present. If we take a more general polynomial of i -degree one, we find the following IBP, which has all three quadratic terms present. We can remove only one via choices of a 1,2 , if we want to obtain a non-trivial equation.
If we use both vector pairs given in eqs. (4.5) and (4.6), the situation is different. Each can be multiplied by a different polynomial; as they have different engineering dimensions, we must choose the polynomials to have different dimensions too. Taking the simplest possibility, multiplying the first pair by the polynomial in eq. (4.15), and the second pair by a constant expression, we obtain the IBP, We can now isolate each quadratic term separately by choosing a 1,2 and b 1 appropriately; for example, with we obtain an IBP for P * * 2,2 [t 2 14 ] in terms of integrals with simpler (lower i -degree) numerators and reducible integrals, (4.20) Upon substituting eq. (4.13) for the lower i -degree polynomial t 14 , this gives us a direct equation for P * * 2,2 [t 2 14 ], We can view the polynomials (4.15,4.17) with the values of a 1,2 and b 1 given by eq. (4.19) as conjugates to t 2 14 , for the given basis of IBP-generating vectors (4.5,4.6).
Similarly, we can also find direct equations for the other two quadratic IrIs. Taking, Taking , (4.24) in eq. (4.18), we obtain a direct equation for P * * 2,2 [t 2 21 ], We will generalize these choices to higher powers of irreducibles in later sections.

V. MASTER INTEGRALS
In the previous section, we found equations to directly reduce quadratic target monomials in IrIs, of the form, 0 = P * * 2,2 [target + simpler IrIs + purely reducible] . (5.1) (The 'simpler' term may contain no IrIs at all, but only powers of external invariants.) In these cases, we only needed one solution each for the different polynomials that can appear in eq. (2.5), rather than the most general solution. Proceeding by plugging in simple ansätze for the solution, and solving for the coefficients, is perhaps not the most elegant way to proceed, but it is adequate.
In contrast, in order to determine that the 'target' is a master integral 3 , we need to show that there is no polynomial solution to the requirement that eq. (2.5) give rise to an IBP equation (5.1) or one with the 'simpler' term missing, To do this in generality, it may be possible to use computational algebraic geometry methods for D-modules. I leave an investigation of this possibility to future work. Here, I will limit myself to showing that there is no solution for polynomials up to some degree, and assume that no solution miraculously appears for higher-degree polynomials.
For our purposes here, instead of writing out all terms in the polynomial in eq. (2.5), it will be more convenient to write out all monomials as a vector, multiplying by appropriate powers of the selected external invariant s in order to make the engineering dimensions of all entries uniform; for example, a vector of 'degree 3' for the double box would be, We need one such vector for each tuple of IBP-generating vectors. The entries in different vectors will of course be of different engineering dimensions in order to ensure that the resulting IBP equations will be of homogeneous engineering dimension. Independently substituting each entry of each vector for Poly in eq. (2.5) leads to a big vector of IBP equations. Setting the purely reducible terms in this vector to zero yields a matrix which can be regarded as a linear transformation of a vector of monomials of the appropriate engineering dimension. Each row of the matrix corresponds to an IBP relation; each column, to a different monomial. The number of possible reductions corresponds to the dimension of the range of this matrix, while the number of master integrals is given by the dimension of its kernel. This latter number is the number of redundant candidate IBP relations. A basis for its kernel, simplified using non-trivial IBP reductions, then gives candidates for the master integrals themselves.
As an example, let us derive the master integrals for the double box. Although the third pair of generating vectors turns out not to be needed for reductions of the double box, we don't know that ahead of time, and so I include it here (taking n v = 3 in eq. (4.11)). The corresponding prefactor for the first term of eq. (4.11) is also given in appendix A.
The simplest construction takes a degree-zero vector for the first pair (4.5), and omits the second (4.6) and third (A1) pairs of IBP-generating vectors. This just yields the matrix form of eq. (4.13), where the columns correspond to the monomials t 14 , t 21 , and s 12 , and the corresponding IBP equation is, The matrix M 1 has a kernel space of dimension 2, generated by the two vectors, corresponding to P * * 2,2 [t 14 +t 21 ] and s 12 P * * 2,2 [1], respectively. Using the non-trivial IBP equation and removing overall constant factors, we then obtain P * * 2,2 [t 21 ] and P * * 2,2 [1] as masters. The construction of M 1 does not make use of the second (or third) pair of vectors, so one may worry that it is missing information. We can proceed to polynomials of one higher degree, using a degree-one vector for the first IBP-generating pair, and degree-zero vectors for the second and third pairs of IBP-generating vectors. This yields for the linear transformation, where the columns now correspond to t 2 14 , t 14 t 21 , t 2 21 , s 12 t 14 , s 12 t 21 , and s 2 12 , respectively. This matrix again has a kernel of dimension 2, and gives rise to the same master integrals. Repeating this procedure with polynomials of one higher dimension again yields the same result, too.

VI. HIGHER POWERS OF IRREDUCIBLE INVARIANTS
Let us continue the approach of finding targeted IBP equations with higher powers of the irreducible invariants. We seek equations that directly reduce them to simpler invariants, that is combinations of invariants of lower engineering dimension. We can do this by taking higher-dimension polynomials in our basic equation (2.5). For example, multiply the first vector pair (4.5) by, and the second pair (4.6) by, Feeding it through the differentiation (making use of eq. (4.11)), we then obtain the IBP equation,

(6.5)
We can find an equation that avoids introducing t 21 by starting with a slightly more general polynomial. Multiply the first vector pair (4.5) by, a 1 t n−1 14 + a 2 t n−2 14 t 21 + a 3 t n−2 14 s 12 , 6) and the second pair (4.6) by the same coefficient as above. Taking the same values as in eq. (6.4) along with, , (6.7) we obtain the following simplified equation for P * * 2,2 [t n 14 ], (6.8) or equivalently, + simpler topologies . (6.9) Eq. (6.9) reduces t n 14 in the numerator to two integrals with lower-dimension irreducible numerators (in addition to integrals with simpler topologies). One may wonder whether it is possible to find an equation that has only one integral with lower-dimension irreducibles.
Even with higher-order polynomials, however, this does not seem possible. (Not too surprisingly, using the third vector pair (A1) does not change this conclusion.) What higher-order polynomials do make possible is greater reduction of the degree in t 14 in one reduction step. Multiply the first vector pair by, a 1 t n−1 14 + a 2 t n−2 14 t 21 + a 3 t n−2 14 s 12 + a 4 t n−3 14 t 2 21 + a 5 t n−3 14 t 21 s 12 + a 6 t n−3 14 s 2 12 , (6.10) and the second pair by, Choosing, with a 4 arbitrary, we find the following equation, (6.13) This result could also be obtained by a partial iteration of eq. (6.9), applying it to the t n−1 14 term on its right-hand side.
While we have implicitly taken n to be an integer in the derivations above, there is nothing that requires it to be one. It can be an arbitrary real value; the difference comes in the stopping conditions -a non-integer n would not ultimately reduce to one of the master integrals, but would require new masters, also with fractional powers of t 14 .

VII. HIGHER PROPAGATOR POWERS
The IBP-generating vectors are designed to avoid introducing doubled propagators (or even higher powers) when they are not present initially. They can of course still be used if such higher powers are present at the beginning of a calculation. The generating vectors will ensure that no powers higher than those present originally will be generated by taking derivatives. When doubled propagators are present, the structure of the IBP equations changes; instead of containing just terms with irreducible numerators along with integrals corresponding to simpler topologies, a new kind of term appears, corresponding to the original topology, but with a lower power of the doubled propogators, 0 = P * * 2,2 [target + lower i -degree + lower propagator powers + purely reducible] . (1 + χ 14 ) s 12 + purely reducible .

(7.3)
Reducibility here again means integrals with fewer propagators (simpler topologies), though one of the surviving propagators will still be doubled.
In order to solve for integrals with doubled propagators, we must generalize the polynomials multiplying the generating vectors to rational functions, with a denominator power corresponding to each doubled propagator. We can repeat the analysis of Sect. V to find master integrals in the presence of doubled propagators. Here we must take appropriate additional powers of a propagator multiplying the numerator insertion. In the case of the double box, we find that the structure of the equations changes. Using just the first IBPgenerating pair (4.5), we find two additional masters beyond those given in eq. (4.10); with a polynomial of engineering dimension 2 multiplying the first pair, and constants multiplying the second and third pairs, we find one additional master; and with a polynomial of engineering dimension 4 multiplying the first pair, and polynomials of engineering dimension 2 multiplying the second and third pairs, we find no additional masters beyond eq. (4.10).
This means that all integrals with doubled propagators can be reduced to linear combinations of integrals with lone propagator powers and integrals with simpler topologies. In this case, we also find that the third IBP-generating pair (A1) is no longer redundant, but is in fact required to obtain a sufficient number of equations.

VIII. SOLVING GENERAL POWERS
In Sect. VI, we saw how to obtain a reduction for an arbitrary power of an irreducible invariant, in the form of eq. (6.9). One could imagine reducing a double-box integral with a given high numerator power of the irreducible invariant by repeatedly applying this reduction, until it ultimately terminates (for integer n) when n = 2. We would then be left with integrals which are either masters or directly expressible in terms of masters.
We could also try to solve the recurrence directly. If we define, w n ≡ s −n 12 P * * 2,2 t n 14 , and drop the purely reducible (simpler-topology) terms in eq. (6.8), that equation takes the form, where ' . =' denotes the dropping of simpler topologies.
We will ultimately turn this recurrence into a differential equation, and solve the latter.
Before doing so, however, let us look at a simpler example. Let us study the sunrise integral P 0,0 , shown in Fig. 3, where K 2 = 0. This two-point topology has just one master integral, which I take to be P 0,0 [1], and two irreducible invariants, t 1 = 1 · K and t 2 = 2 · K. It depends only on the kinematic invariant s = K 2 . There are three linearly independent pairs of IBP-generating vectors, (8.4) where r 11 = 2 1 and r 22 = 2 2 as in eq. (4.3), and Let us try to compute P 0,0 [t n 1 ]. The integral is simple enough that we can compute it directly, using the following expression for the one-loop bubble with arbitrary exponents, Performing the 2 integration first, we obtain, Alternatively, we can proceed using the IBP-generating vectors. Take a linear combination of the first and third vector pairs in eq. (8.4), with coefficients, , (8.8) and 1 4(n + 2 − 3 ) , (8.9) respectively. We then find the following equation, Defining y n ≡ s −n P 0,0 t n 1 , (8.11) we have the recurrence relation, We can solve this equation (for example, using Mathematica), obtaining the result, (8.13) in agreement with the explicit computation in eq. (8.7).

B. Differential Equations
It can be difficult to solve the more general recurrence relations such as eq. (8.2) directly (Mathematica, for example, can solve them but provides the solution in an implicit and rather unenlightening form in terms of DifferenceRoot objects). Instead, introduce the generating function, a n x n , (8.14) and derive a differential equation for it. Once one has solved the differential equation, one can obtain the solution for a n by series-expanding the solution. One approach to obtaining a differential equation is to use the RISC-Linz Mathematica package GeneratingFunctions [19]; but one can also proceed in a more pedestrian fashion, as described here.
First recast the recurrence relation Rec so that the indices of a appearing in it are strictly positive for n ≥ 0, and then sum the recurrence (depending on n) into a generating object, Rec n x n . (8.15) Then apply the substitution rule, ∞ n=0 c n a n+r x n → x −r ∞ n=0 c n−r a n x n − x −r r−1 n=0 c n−r a n x n . (8.16) In this rule, c n is a polynomial in n and r ≥ 0; we need consider only linear functions of n (because the single derivative generating the IBP identity can bring down only a single power of an exponent; though factors of n in coefficients could in principle alter this). Finally, using the operator, 8.17) replace ∞ n=0 n p a n x n → D p x f (x) . (8.18) In the recurrences we consider, this will give an inhomogeneous first-order differential equation. It turns out to be easier to solve (using Mathematica, anyway) a higher-order homogeneous equation obtained by further differentiation. The behavior of f (x) as x → 0 provides the additional boundary conditions needed for the higher-order equation. The Gen-eratingFunctions package produces such a higher-order equation directly. In the next two subsections, I give examples of using differential equations for the generating function to solve recurrence relations for general powers of numerator insertions.
The slashed box has one master integral, which we can choose to be P 1,1 [1].

(8.32)
Here, F 1 is the first Appell function. The first term is not well-defined for < 0 as x → 0, so c 1 must vanish. The second constant of integration, c 2 , is fixed by the requirement that The last constant, c 3 , is fixed in terms ofw 1 via f (0) =w 1 ; butw 1 in turn is not independent, because for n = −1, the recursion (8.29) becomes a two-term relation, (8.34) We ultimately find that c 3 = 0. The solution with the desired boundary behavior is thus, We can then extract the n-th term of this function to obtain an expression forw n . After a bit of algebra and simplification, we find the following expression, (1 + χ 14 ) n 4 Γ(n 1 − n 4 + 2 ) Γ(n 4 + ) (n 1 − n 4 − 1)! n 4 ! . (8.36) By a bit of guesswork, we can find a more compact expression 4 , This form is manifestly a rational function of χ 14 and , as the hypergeometric function terminates for integer n.
Either form meets the challenge posed in eq. (3.3), up to terms arising from simpler topologies.

D. The Double-Box Integral
Let us return to the double-box integral, and the recurrence relation given in eq. (8.2) (with the definition of w n in eq. (8.1)). Rewriting the equation to make all indices positive for n > 0, we obtain, Again using the approach described in Sect. VIII B, we can obtain the corresponding firstorder differential equation, 8.39) As in the case of the slashed box in the previous subsection, we can differentiate twice with respect to x to obtain a more convenient third-order equation, We can solve this latter equation to obtain, . Once again, the first term is not well-defined for < −1/2 as x → 0, and so c 1 must vanish. The second constant of integration c 2 is fixed by the requirement that f (0) = w 0 to be, The third constant c 3 is fixed by the requirement that f (0) = w 1 , The solution with desired boundary behavior is then, .
It isn't obvious how to write down an analog of eq. (8.37), an expression which is given purely in terms of hypergeometric functions and yet is manifestly rational in χ 14 and .

IX. CONCLUSIONS
Finding linear relations between Feynman integrals plays a key role in higher-loop calculations in quantum field theory. Integration by parts has become the method of choice for finding such relations, but the conventional approach to using them leads to equations involving many unwanted integrals with doubled propagators. Moreover, the standard method for solving them requires cumbersome handling of large systems of equations. The first issue can be addressed using the generating-vector approach first introduced in ref. [13]. In this paper, I presented an approach to simplify the second issue. It eliminates the need to handle large systems of equations by allowing one to target desired numerator terms, and derive direct reduction equations for them. A specific numerator can be isolated by choosing appropriate polynomial prefactors for each of the generating-vector tuples for the integral topology under study. One can do this for specific terms, as in the examples of eqs. (4.20-4.25). One can also do this for general powers of irreducible invariants, something not possible in the conventional approach. I gave examples in eqs. (6.3) and (6.5). As an example of the power of the new approach, I showed how to obtain closed-form reductions to master integrals for such arbitrary powers, in eqs. (8.36) and (8.45). It is also possible to find master integrals within the new approach, as seen in Sect. V, though the strategy outlined there can undoubtedly be improved with more insight from algebraic geometry. The generalization of generic-power equations to multiple irreducible invariants, not discussed explicitly in the present paper, is straightforward. Solving the corresponding differential equations, as in Sect. VIII is less straightforward, as one-variable differential equations are replaced by systems of partial differential equations, but should be possible using appropriately designed series Ansätze. Even without explicit solutions to generic powers, the approach described in this paper will greatly simplify integral reductions to masters, and should make possible new calculations at the high-loop frontier in a variety of quantum field theories.