Computer algebra tools for Feynman integrals and related multi-sums
DDESY 18-160, DO-TH 18/21
Computer algebra tools for Feynman integrals andrelated multi-sums ∗ Johannes Blümlein
Deutsches Elektronen–Synchrotron, DESY,Platanenallee 6, D-15738 Zeuthen, GermanyEmail:
Carsten Schnneider † Research Institute for Symbolic Computation (RISC)Johannes Kepler University, Altenbergerstraße 69, A-4040 Linz, AustriaE-mail:
In perturbative calculations, e.g., in the setting of Quantum Chromodynamics (QCD) one aimsat the evaluation of Feynman integrals. Here one is often faced with the problem to simplifymultiple nested integrals or sums to expressions in terms of indefinite nested integrals or sums.Furthermore, one seeks for solutions of coupled systems of linear differential equations, that canbe represented in terms of indefinite nested sums (or integrals). In this article we elaborate themain tools and the corresponding packages, that we have developed and intensively used withinthe last 10 years in the course of our QCD-calculations.
Loops and Legs in Quantum Field Theory (LL2018)29 April 2018 - 04 May 2018St. Goar, Germany ∗ This work was supported in part by the Austrian Science Fund (FWF) grant SFB F50 (F5009-N15) † Speaker. c (cid:13) Copyright owned by the author(s) under the terms of the Creative CommonsAttribution-NonCommercial-NoDerivatives 4.0 International License (CC BY-NC-ND 4.0). https://pos.sissa.it/ a r X i v : . [ c s . S C ] S e p omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider
1. Introduction
We aim at the simplification of loop Feynman integrals to expressions in terms of specialfunctions and constants. More precisely, we consider µ -loop massive Feynman parameter integrals(e.g., µ =
3) emerging in renormalizable Quantum Field Theories, like Quantum Electrodynamicsor Quantum Chromodynamics, see e.g. [24, 70], which can be transformed to a linear combinationof s -fold multiple integrals of the form F ( m , . . . , m r , n , ε ) = (cid:90) · · · (cid:90) f ( m , . . . , m r , n , ε , x , . . . , x s ) dx . . . dx s ; (1.1) m , . . . , m r are the involved masses, the discrete parameter n stands for the Mellin variable, and D : = + ε with ε ∈ R is the analytic continuation of the space-time dimension. Often one canconsider m , . . . , m r and ε as indeterminates (and reinterprets them, only if necessary, as concreteparameters that can be evaluated within certain ranges). In the following we assume that the massesare contained in a field K of characteristic 0. A crucial property is that the integrand f of sucha Feynman integral is hyperexponential in each of the integration variables x i ( ≤ i ≤ s ) andhypergeometric in the discrete parameter n .Given such an integral F ( n , ε ) , one seeks for the first coefficients of their Laurent seriesexpansion F ( n , ε ) = F l ( n ) ε l + F l + ( n ) ε l + + · · · + F r ( n ) ε r + O ( ε r + ) where the coefficients are given in terms of special functions and constants; usually, we have l = − µ for a µ -loop Feynman integral and r = , ,
2. More precisely, we applied the following generalstrategies in the course of our calculations:Feynman integral symbolic integration (cid:15) (cid:15) p F q -/Mellin-Barnestechnologies (cid:35) (cid:35) multi-sum expressions symbolic summation (cid:113) (cid:113) special function expressionsFirst one may apply directly symbolic integration algorithms (see Section 2.1) in combinationwith recurrence solving technologies (see Section 2.2) to simplify these integrals in terms of spe-cial functions. Second one can transform in a preprocessing step these integrals to expressions interms of multiple nested sums. Namely, by applying successively Newton’s binomial theorem andMellin-Barnes decompositions on the integrand, implemented in different packages [33, 36, 65],one can carry out all integrals by introducing Mellin-Barnes integrals. Finally, applying the residuetheorem to these introduced Mellin-Barnes integrals yields multiple sums over hypergeometric ex-pressions. We emphasize that this mechanism has to be applied very carefully in order to arrive at h ( x ) is hyperexponential (resp. hypergeometric) in x if h (cid:48) ( x ) h ( x ) (resp. h ( x + ) h ( x ) ) is a rational function in x . In the following we will suppress the mass dependencies. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider expressions that are suitable for our symbolic summation methods. In general, we will end up at alinear combination of hypergeometric multi-sums of the form S ( n , ε ) = L ( n ) ∑ k = . . . L v ( n , k ,..., k v − ) ∑ k v = f ( n , ε , k , . . . , k v ) (1.2)over K ( n , ε ) . Here the upper bounds L ( n ) , . . . , L v ( n , k , . . . , k v − ) are integer linear (i.e., linearcombinations of the variables over the integers) in the dependent parameters or ∞ , and f is hyper-geometric in n and the summation variables k i . Then given such a sum representation, we are inthe position to apply our symbolic summation tools summarized in Section 2.3.In many applications, one is faced with thousands (even millions) of Feynman integrals thatdescribe an underlying physical problem. To treat them directly with the above methods is usu-ally out of scope. As a preprocessing step, one utilizes integration by parts (IBP) methods [31,42, 46, 47, 67]. They crunch the given expression to a more compact form in terms of only fewintegrals, that have to be treated individually. Similarly, while transforming Feynman integralsto multiple sums, one obtains enormous expressions consisting of up to thousands of multiplesums. To simplify these sums successively using these summation tools is not only problematicbecause of time limitations, but also because of the following intrinsic problem. Often the scat-tered sums themselves cannot be simplified in terms of indefinite nested sums, but only a suitablecombination of them can be simplified in this way. In order to bypass these problems, the package SumProduction [23, 61] built on
Sigma [57, 60] can be utilized. It reduces the sum expres-sions to compact forms where the arising sums are merged appropriately. Afterwards the symbolicsummation can be applied to these expressions within reasonable time and without dealing withsums that cannot be handled within the difference ring setting [49, 63, 64]. Summarizing, the fol-lowing simplification techniques are applied in addition in order to reduce the given problem to areasonable size of integrals or sums. huge expression ofFeynman integrals
IBP methods (cid:1) (cid:1) p F q -/Mellin-Barnestechnologies (cid:28) (cid:28) huge expression in terms offew master integrals huge expressionof (small) multi-sums merging/reduction (cid:15) (cid:15) huge expression in termsof few (large) multi-sumsAn extra advantage of the IBP approach is that most of the produced master integrals are describedas solutions of coupled systems of linear differential equations. Only few integrals (usually of sim-pler type) arise in the inhomogeneous part of these systems which themselves are not representedas solutions of coupled systems. Thus the remaining task is to simplify these few integrals by sym-bolic summation/integration methods and to derive the remaining (and usually more complicated)integrals by solving the provided coupled systems of linear differential equations.2 omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider
In a nutshell, one is faced with the problem to simplify Feynman integrals by means of sym-bolic integration (see Section 2.1 in combination with Section 2.2), symbolic summation (see Sec-tion 2.3) or solving coupled systems of linear differential equations (see Section 2.4).
2. Computer algebra methods
In the following we will present the most relevant computer algebra tools, that were crucial tocarry out the challenging QCD-calculations described in [2–4, 8, 18–20, 22]. For a comprehensivesummary of further available tools we refer to the recent survey article [28].
Concerning symbolic integration tools we heavily utilized the multivariate Almkvist-Zeil-berger algorithm [16, 17], in particular an optimized and improved version [1, 4] for Feynmanintegrals that can tackle the following problem.
Recurrence finding.
Given an integral in the form (1.1) where f is hyperexponential in x i for i = , , . . . , s and hyper-geometric in n . Compute a linear recurrence of the form a ( n , ε ) F ( n , ε ) + a ( n , ε ) F ( n + , ε ) + · · · + a d ( n , ε ) F ( n + d , ε ) = a i ( n , ε ) ∈ K [ ε , n ] where a d ( n , ε ) (cid:54) = K ( n , ε ) . The complexityto solve the underlying systems depends heavily on the involvement of the constructed field K , thenumber of the integration variables x , . . . , x s and on the structure of the integrand f (in particularhow the variables are intertwined). However, for various concrete situations this method worksvery well. Example 1.
If one applies the package
MultiIntegrate [1, 4] to the multi-integral F ( n , ε ) = (cid:90) (cid:90) (cid:90) (cid:90) − u ( x + y − ) n x ε / ( − x ) ε / y ε / ( − y ) ε / ( − u − v ) n (cid:0) − u xx − − v yy − (cid:1) − + / ε u + ε / v + ε / dx dy du dv , (2.2)one can compute a recurrence of the form (2.1) with order d = Remark 1.
For some special cases we could also utilize an extended version [45] of the hyperlog-arithm method [30, 50] in order to tackle massive Feynman integrals.
Suppose we calculated a recurrence of the form (2.1) or more generally of the form a ( n ) F ( n ) + a ( n ) F ( n + ) + · · · + a d ( n ) F ( n + d ) = b ( n ) (2.3)with a i ( n ) ∈ K [ n ] for 0 ≤ i ≤ d where a d ( n ) (cid:54) = b ( n ) ∈ K for n ≥
0. Let δ = max ( { i ∈ N | a d ( i ) = } ∪ {− } ) +
1. Then it follows immediately that there is a unique sequence solution3 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider F ( n ) with n ≥ δ of (2.3) if one fixes the first d initial values F ( δ ) , . . . , F ( δ + d − ) ∈ K .An important task is then to represent, if possible, such a solution in terms of special functions.Within the summation package Sigma [57, 60] difference ring algorithms have been encoded [14,15, 37, 39, 49, 53, 55, 56, 58, 59, 62–64], that find such a representation in terms of indefinitenested sums over hypergeometric products, whenever this is possible.
Definition 1.
Let K be a field of characteristic 0. A product ∏ kj = l f ( j ) , l ∈ N , is called hypergeo-metric in k over K if f ( y ) is an element from the rational function field K ( y ) where the numeratorand denominator of f ( j ) are nonzero for all j ∈ Z with j ≥ l . An expression in terms of indefi-nite nested sums over hypergeometric products in k over K is composed recursively by the threeoperations ( + , − , · ) with • elements from the rational function field K ( k ) , • hypergeometric products in k over K , • and sums of the form ∑ kj = l f ( j ) with l ∈ N where f ( j ) is an expression in terms of indefinitenested sums over hypergeometric products in j over K ; here it is assumed that the evaluationof f ( j ) | j (cid:55)→ λ for all λ ∈ Z with λ ≥ l does not introduce any poles.If K and k are clear from the context we simply say that f ( k ) is an expression in terms of indefinitenested sums (over hypergeometric products) . Example 2.
The class of indefinite nested sums over hypergeometric products in n over K coversas special cases harmonic sums [25, 69] like S , ( n ) = n ∑ i = i i ∑ j = j , generalized harmonic sums [12, 48] like n ∑ k = k k k ∑ i = − i i i ∑ j = j , cyclotomic harmonic sums [11] like n ∑ k = ( + k ) k ∑ j = j j ∑ i = + i or nested binomial sums [10] like n ∑ k = (cid:0) kk (cid:1) k ∑ j = j j ∑ i = (cid:0) ii (cid:1) + i . We can solve the following problem with the summation package
Sigma for the class of indefinitenested sums. 4 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider
Recurrence solving.
Given a i ( n ) ∈ K [ n ] for 0 ≤ i ≤ d where a d ( n ) (cid:54) = b ( n ) that can be calculated by anexpression in terms of indefinite nested sums over hypergeometric products; given δ as aboveand c δ , . . . , c δ + d − ∈ K . Decide constructively if the solution F ( n ) of (2.3) with F ( i ) = c i for δ ≤ i ≤ δ + d − ε . In this case we could treat ε just as an extra parameter which is contained in K and seek for a solution within the class ofindefinite nested sums over hypergeometric products over K . However, in most cases such an all n and ε representation does not exist. A more flexible strategy is to hunt for a solution F ( n , ε ) givenin its ε -expansion F ( n , ε ) = F l ( n ) ε l + F l + ( n ) ε l + + . . . (2.4)More precisely, consider a recurrence of the form a ( n , ε ) F ( n , ε ) + a ( n , ε ) F ( n + , ε ) + · · · + a d ( n , ε ) F ( n + d , ε ) = b ( n , ε ) (2.5)with a i ( n , ε ) ∈ K [ n , ε ] where not all a i ( n , ) are zero and with a right-hand side given in its ε -expansion b ( n , ε ) = ∞ ∑ i = l b i ( n ) ε i (2.6)where b i ( n ) ∈ K for i ≥ l and n ≥
0. Let o be the largest integer such that a o ( n , ) (cid:54) = δ = max ( { i ∈ N | a o ( i , ) = } ∪ {− } ) +
1. Then by a slight variation of [24] it follows that there is aunique Laurent series solution (2.4) for n ≥ δ when fixing F i ( j ) ∈ K for i ≥ l and δ ≤ i ≤ δ + o .In particular, the following holds. If there are two solutions, say F ( n , ε ) with (2.4) and F (cid:48) ( n , ε ) = F (cid:48) l ( n ) ε l + F (cid:48) l + ( n ) ε l + + . . . and if F i ( j ) = F (cid:48) i ( j ) for l ≤ i ≤ r and δ ≤ i ≤ δ + o , then F i ( n ) = F (cid:48) i ( n ) for all l ≤ i ≤ r and n ≥ δ . In short, knowing the first initial values determines uniquely the firstcoefficients of a Laurent series solution. In particular, they can be prolonged stepwise to a fullLaurent series solution by fixing further initial values of higher ε -orders.With this preparation step we can now introduce the following machinery from [24] implementedin Sigma . Recurrence solving for ε -expansions. Given a i ( n , ε ) ∈ K [ n , ε ] for 0 ≤ i ≤ d where not all a i ( n , ) are zero and given b ( n , ε ) with (2.6)where b l ( n ) , . . . , b r ( n ) are represented in terms of indefinite nested sums over hypergeometricproducts; given o , δ as given above, and given c i , j ∈ K for l ≤ i ≤ r and δ ≤ i ≤ δ + o . Decide constructively if the first coefficients F l ( n ) , . . . , F r ( n ) of a Laurent series solution (2.4) of therecurrence (2.5) with F i ( j ) = c i , j can be calculated by expressions in terms indefinite nestedsums over hypergeometric products. Example 3.
We continue Example 1 by taking the recurrence of order 5 (see [4, pp. 49–52]) forthe integral F ( n , ε ) given in (2.2). In addition, we take the five initial values F ( , ε ) , . . . , F ( , ε ) If all a i ( n , ) are zero, we can divide the equation by ε r for some r ∈ N yielding again polynomial coefficients onthe left-hand side where not all are zero when setting ε to 0. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider expanded up to ε − . E.g., we have F ( , ε ) = ε − ε + ε (cid:18) + ζ (cid:19) + O ( ε ) . Then using the above algorithm we can calculate the first three coefficients in F ( n , ε ) = F − ( n ) ε − + F − ( n ) ε − + F − ( n ) ε − + O ( ε ) of the ε -expansion within a few seconds. More precisely, we get F − ( n ) = ( − ) n ( n + )( n + ) + ( n + ) ( n + ) ( n + ) , F − ( n ) = − ( − ) n (cid:0) n + n + n + (cid:1) ( n + ) ( n + ) − (cid:0) n + n + n + (cid:1) ( n + ) ( n + ) , F − ( n ) = ( − ) n (cid:32) (cid:0) n + n + n + n + n + (cid:1) ( n + ) ( n + ) + ζ ( )( n + )( n + ) (cid:33) + (cid:0) n + n + n + n + n + (cid:1) ( n + ) ( n + ) + ( n + ) ζ ( )( n + ) ( n + )+ (cid:18) − ( n + ) ( n + ) + ( − ) n ( n + )( n + ) (cid:19) S ( n )+ (cid:18) ( − ) n ( n + )( n + ) − ( n + ) ( n + ) ( n + ) (cid:19) S − ( n ) ;here S a ( n ) = ∑ nk = k k a are the generalized harmonic numbers and ζ ( z ) = ∑ ∞ k = k z denotes theRiemann zeta function.As illustrated in the previous example, one can combine the Almkvist-Zeilberger method with theabove recurrence solver to obtain a method for the following problem. Simplifying multiple integrals.
Given an integral in the form (1.1) where f is hyperexponential in x i for i = , , . . . , s and hy-pergeometric in n . Compute the first coefficients of its ε -expansion in terms of indefinite nestedsums over hypergeometric products. Given a sum representation of the form (1.2), symbolic summation algorithms in the setting ofdifference rings and fields [14, 15, 37, 39, 49, 53, 55, 56, 58, 59, 62–64] can be utilized to derivean alternative representation in terms of indefinite nested sums over hypergeometric products.
Simplifying multiple sums.
Given a multiple sum of the form (1.2) where f is hypergeometric in n and k i for i = , , . . . , v . Find a simplified version in terms of indefinite nested sums over hypergeometric products. If the integrand has a particularly nice shape (e.g., if it is proper hyperexponential in the x i and proper hypergeo-metric in n ), the multivariate Almkvist Zeilberger method will provide such a recurrence (2.5) with b ( n , ε ) =
0. In short,for such special input the method under consideration turns to a complete decision procedure. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider
The proposed method works from inside to outside of such a multiple sum and transforms ineach step the arising definite summation quantifier to an indefinite nested version by using theMathematica package
Sigma [57, 60]. Namely, suppose that we transformed already a sub-sum of (1.2) to an expression in terms of indefinite nested sums w.r.t. k yielding the expression f ( . . . , m , k ) . Moreover, suppose that we have to deal with an extra summation quantifier, say with F ( . . . , m ) = ∑ L ( ..., m ) k = f ( . . . , m , k ) . If this summation quantifier is the outermost sum in (1.2), then m is precisely the Mellin variable n . Otherwise, m will be the summation variable of the next summa-tion quantifier that is applied to F ( . . . , m ) . Then the following three steps lead often to the desiredsimplification.1. Given the definite sum F ( m ) (we suppress further variables), try to compute a linear recur-rence relation of the form (2.3) ( n replaced by m ) using the creative telescoping paradigm [51]in the setting of difference rings [58, 62–64]. Here a i ( m ) are polynomials in m and b ( m ) isan expression in terms of indefinite nested sums over hypergeometric products.2. Solve afterwards the recurrence in terms of indefinite nested sums over hypergeometric prod-ucts w.r.t. m using the toolbox described in Section 2.2.3. Finally, try to combine the solutions to derive a simpler representation of the original inputsum in terms of indefinite nested sums over hypergeometric products w.r.t. m .Suppose that all three steps can be carried out. If n = m is the Mellin parameter in (1.2) then weare done. Otherwise, F ( . . . , m ) will take over the role of f ( . . . , k ) and we repeat this game to treatthe next summation quantifier in (1.2).In some rare cases, this machinery even works if the input sum depends on ε . However, in mostcases one will arrive at linear recurrences that depend on ε and that does not have sufficiently manysolutions in terms of indefinite nested sums. As a consequence, step 3 from above will fail. Avery successful strategy is based on calculating an ε -expansion of the summand of (1.2) up to thedesired order and applying the summation quantifiers to each of the coefficients (which are freeof ε ). Afterwards, the summation mechanism described above is applied to each of the arisingsummation problems. In many applications, the derived recurrences (now free of ε ) are completelysolvable in terms of indefinite nested sums and all three steps 1–3 from above can be carriedout successfully. The package EvaluateMultiSums [60, 61] based on
Sigma combines allthese steps and variants thereof in a very efficient way in order to carry out such simplificationsautomatically. This will works in general if one finds d linearly independent solutions of the homogeneous version of the recur-rence and one particular solution of the recurrence itself. If infinite sums are involved, properties such as uniformal convergence have to be taken extra into account. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider
Example 4.
Consider the 2-mass 3-loop Feynman diagram , (2.7)that arose in the calculation of the gluonic operator matrix element A ( ) gg , Q [7]. The correspondingintegral is transformed to an expression of size 95MB given as a linear combination of multiplesums of the form (1.2). More precisely, the expression is built by 150 single sums, 1000 doublesums, 12160 triple sums and 1555 quadruple sums. A typical triple sum is n ∑ j = j ∑ i = i ∑ k = ( + ε )( − + n )( − + n ) n π ( − ) − k + ε × − + ε e − εγ η k ×× Γ ( − ε − i + j + k ) Γ ( − − ε ) Γ ( + ε ) Γ ( + n ) Γ ( + ε + i − k ) Γ ( − ε + k ) Γ ( − ε + k ) Γ ( − ε + k ) Γ ( − − ε + k ) Γ ( − − ε ) Γ ( + ε ) Γ ( + i ) Γ ( + k ) Γ ( − i + j ) Γ ( − ε + k ) Γ ( − ε + k ) Γ ( − ε + k ) Γ ( + ε + n ) where η = m m is the quotient of the two arising masses. Applying the above tools for the ε − -contribution yields a simplification in terms of indefinite nested sums in about 30 minutes. Thissuggests that the full ε − contribution of the above diagram will require about 1 year of calculationtime. Interestingly enough, our summation tools failed to produce the ε contribution. More pre-cisely, we obtained recurrences that could not be solved in terms of indefinite nested sums.As a consequence, we used the package SumProduction to crunch the derived expression con-sisting of 14865 sums to an expression consisting only of 8 sums (which are rather large). After-wards, we applied our summation tools to these 8 sums. The calculations can be summarized inthe following table [7]. sum size of sum summand size of time of number of(with ε ) constant term calculation indef. sums n − ∑ i = i − ∑ i = i ∑ i = ∞ ∑ i = n − ∑ i = i − ∑ i = ∞ ∑ i =
232 MB 1646.4 MB 980756 s (11.4 days) 747 n − ∑ i = ∞ ∑ i = ∞ ∑ i = n − ∑ i = i − ∑ i = i ∑ i = i ∑ i = n − ∑ i = i − ∑ i = i ∑ i = n − ∑ i = i ∑ i = n − ∑ i = For instance, consider the quadruple sum of the form ∑ n − i = ∑ i − i = ∑ i i = ∑ ∞ i = . In total the summandrequires 17.7MB of memory. Taking its ε -expansion, the constant term requires 266.3 MB. Ap-plying afterwards the package EvaluateMultiSums to the quadruple sum of the constant term The graph has been drawn using
Axodraw [68]. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider yields within 2.1 day an expression in terms of 1188 indefinite nested sums. In order to treat all 8sums, we needed 3 month of calculation time. In total we obtain an expression for diagram (2.7)that uses 154MB of memory and is composed by 4110 indefinite nested sums. For instance one ofthe most complicated sums is n ∑ h = − h ( − η ) h (cid:18) hh (cid:19) (cid:32) h ∑ i = i ( − η ) − i i ( ii ) (cid:33) (cid:32) h ∑ i = ( − η ) i ( ii ) i (cid:33) h ∑ i = i ( − η ) − i i ∑ j = j ∑ k = ( − η ) k k ji (cid:0) ii (cid:1) . Finally, one can simplify this expression further by using the available difference ring algorithmsof
Sigma . The final result is an expression of only 8.3 MB size in terms of 74 binomial sums thatare all algebraically independent among each other.
Remark 2.
We also started to explore the usage of holonomic summation tools [32, 54] in the set-ting of QCD. We strongly expect that new ideas of [26] will push forward the available summationtechniques.
After applying IBP methods [31, 42, 46, 47, 67] physical quantities are crunched to expressionsin terms of master integrals. In particular, these master integrals are described as solutions ofcoupled systems of linear difference or differential equations. E.g., in the univariate differentialequation case they are of the form D x f ( x , ε ) ... f λ ( x , ε ) = A f ( x , ε ) ... f λ ( x , ε ) + g ( x , ε ) ... g λ ( x , ε ) , (2.8)where A is a λ × λ matrix with entries from K ( x , ε ) and the entries of the right-hand side vec-tor g ( x , ε ) = ( g ( x , ε ) , . . . g λ ( x , ε )) are given as a linear combination of simpler master integrals h ( x , ε ) , . . . , h u ( x , ε ) over K ( x , ε ) . They can be either determined by other coupled systems (byusing the proposed method recursively) or have to be tackled, e.g., by our symbolic summationand integration tools introduced above. In the following we suppose that the unknown functions f i ( x , ε ) have a power series representation f i ( x , ε ) = ∞ ∑ k = F i ( k , ε ) x k , (2.9)and we seek for symbolic representations of the coefficients of the ε -expansions F i ( k , ε ) = ∞ ∑ j = l F i , j ( k ) ε j (2.10)with 1 ≤ i ≤ λ . Likewise we assume that the simpler master integrals h i ( x , ε ) with 1 ≤ i ≤ u arising in g ( x , ε ) have power series representations whose coefficients can be represented in termsof indefinite nested sums. More precisely, their ε -expansions are assembled by coefficients that9 omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider can be expressed in terms of indefinite nested sums. Then one can utilize the following algorithmthat is implemented within the package
SolveCoupledSystem [4, 6] and is based on
Sigma . Solving coupled systems.
Given a system as above (2.8) with initial values of the desired solution and l , r ∈ Z with l ≤ r . Decide constructively if the coefficients F i , j ( k ) for i = , . . . , λ , j = l , . . . , r of the power series F i ( k , ε ) in (2.9) which form the coefficients of the ε -expansions (2.10) can be expressed in termsof indefinite nested sums.1. Uncouple the system by using, e.g., Zürcher’s algorithm [29, 71] implemented in the package OreSys [35]. Usually one gets one scalar linear differential equation in one of the unknownfunctions, say f ( x , ε ) where the right-hand side can be given in terms of a power serieswhose coefficients are given explicitly in terms of indefinite nested sums. In addition, theremaining functions f ( x , ε ) , . . . , f λ ( x , ε ) can be expressed by f ( x , ε ) and the simpler masterintegrals h i ( x , ε ) with 1 ≤ i ≤ u in the following form: f i ( x , ε ) = ∑ r α i , r ( x , ε ) D rx f ( x , ε ) + u ∑ j = ∑ r β i , j , r ( x , ε ) D rx h j ( x , ε ) (2.11)for some rational functions α i , r ( x , ε ) , β i , j , r ( x , ε ) ∈ K ( x , ε ) .2. By holonomic closure properties [41] compute a linear recurrence of F ( k , ε ) from the givenscalar differential equation of f ( x , ε ) . The recurrence is of the form (2.5) where the coeffi-cients in (2.6) can be given explicitly in terms of indefinite nested sums.3. Use the tools from Section 2.2 with the corresponding initial values in order to decide algo-rithmically if the first coefficients F , j ( k ) for j = l , . . . , r can be given in terms of indefinitenested sums. If this is not possible, return that such a representation is not possible.4. Otherwise, plug these coefficients into (2.10) (for i =
1) and afterwards plug it into (2.11).Similarly, proceed to plug the known coefficients of the simpler master integrals h i ( x , ε ) with 1 ≤ i ≤ u into (2.11). Finally, extract the first coefficients of the ε -expansions of F ( k , ε ) , . . . , F λ ( k , ε ) . Remark 3.
Often one has to calculate the ε -expansions for F ( k , ε ) and for the simpler masterintegrals h j ( x , ε ) up to a certain order that is higher than r due to factors ε . The correspondingbounds can be determined after the uncoupling has been carried our. We neglect further details onthese technical aspects and refer to [6].Using this toolbox (and slight variants of it [4, 6, 34]) we performed already rather advanced QCD-calculations [2–4, 8, 18–20, 22]. Remark 4.
In certain instances one can also use ideas from [38, 44] to find solutions in terms ofindefinite nested integrals. In general one might obtain several scalar linear differential equations. Then the described steps in 2–4 below arecarried out for each equation. omputer algebra tools for Feynman integrals and related multi-sums Carsten Schnneider
An interesting feature is that this method can be adapted to calculate a large number mo-ments [27], say µ = SolveCoupledSystem . Large moment method.
Given a system as above (2.8) with initial values of the desired solution, l , r ∈ Z with l ≤ r and µ ∈ N . Calculate efficiently F i , j ( ) , F i , j ( ) , . . . , F , j ( µ ) for i = , . . . , λ , j = l , . . . , r of the powerseries F i ( k , ε ) in (2.9) which form the coefficients of the ε -expansions (2.10).1. Suppose that one has calculated already the first µ moments of the simpler master integralsthat arise in the right-hand side of g i ( x , ε ) up to a certain order. This can be accomplished,e.g., by applying this method again to this simpler integrals or by utilizing the represen-tations in terms of indefinite nested sums coming from our symbolic summation and inte-grations tools described above. As a consequence one can calculate the sequence of values b i ( ) , b i ( ) , . . . , b i ( µ ) of (2.6) in (2.5) up to a certain order, say l ≤ i ≤ r .2. Using the recurrence (2.5) with the moments given on the right-hand side one can calculate inlinear time the moments F , j ( ) , F , j ( ) , . . . , F , j ( µ ) for j = l , . . . , r of the ε -expansion (2.10)( i =
1) of F ( k , ε ) in (2.9) provided that F , j ( ) , F , j ( ) , . . . , F , j ( d ) is given (using our sum-mation tools or using from above or procedures like Mincer [43] or
MATAD [66]).3. Finally, one plugs these values into (2.11) and extracts the values F r , j ( ) , F r , j ( ) , . . . , F r , j ( µ ) for r = , . . . , λ and j = l , . . . , r . Remark 5.
The comments in Remark 3 are also here relevant to calculate the correct momentsup to the desired ε -order r . In addition, it might happen that the number of moments of f ( x , ε ) in calculation step (2) and the moments of the simpler master integrals h j ( x , ε ) must be chosenslightly higher than µ in order to get all moments correctly. This is owing to the fact that extrafactors x might occur (also during the uncoupling process) which introduce negative shifts on thecoefficient level.Recall that the IBP methods usually crunch physical expressions to simper expressions in terms ofthe master integrals under consideration. Assembling all these moments produces large moments ofthe physical quantities. With these moments one is now in the position to guess, e.g., homogeneousrecurrence relations using packages like [40]. Afterwards, one can use our recurrence solvers fromSection 2.2 to derive representations in terms of indefinite nested sums. First non-trivial casestudies have been carried out, like the first re-computation of the 3-loop splitting functions [3] andsimpler parts of the massive 3-loop form factor [9].
3. Conclusion
We presented the central methods for symbolic summation/integration, solving linear recur-rences and solving coupled systems of linear differential equations that have been encoded withinthe Mathematica packages
Sigma , EvaluateMultiSum , MultiIntegrate , SumProduc-tion and
SolveCoupledSystem . Besides these computer algebra tools also special functionalgorithms [10–13, 21, 52, 69] implemented within the package
HarmonicSums [1] are used in11 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider order to speed up the above methods; in particular the limit n → ∞ can be performed for expres-sions in terms of indefinite nested sums in n by utilizing asymptotic expansion algorithms.It will be interesting to see how all these computer algebra and special function algorithms canbe generalized to tackle also more complicated special functions like iterative-non-iterative inte-grals and sums [5] as they appeared in recent calculations. In this regard, we expect that the largemoment machinery (see Section 2.4) will play a decisive role within our future calculations. References [1] J. Ablinger.
Computer algebra algorithms for special functions in particle physics . PhDthesis, J. Kepler University Linz, April 2012.[2] J. Ablinger, A. Behring, J. Blümlein, A. De Freitas, A. von Manteuffel, and C. Schneider.The 3-loop pure singlet heavy flavor contributions to the structure function F ( x , Q ) and theanomalous dimension. Nucl. Phys. B , 890:48–151, 2014.[3] J. Ablinger, A. Behring, J. Blümlein, A. De Freitas, A. von Manteuffel, and C. Schneider. Thethree-loop splitting functions P ( ) qg and P ( , N F ) gg . Nucl. Phys. B , 922:1–40, 2017.[4] J. Ablinger, A. Behring, J. Blümlein, A. D. Freitas, A. von Manteuffel, and C. Schneider.Calculating three loop ladder and V-topologies for massive operator matrix elements by com-puter algebra.
Comput. Phys. Comm. , 202:33–112, 2016.[5] J. Ablinger, A. Behring, J. Blümlein, A. D. Freitas, E. Imamoglu, M. van Hoeij, A. vonManteuffel, C. Raab, C. Radu, and C. Schneider. Iterative and iterative-noniterative integralsolutions in 3-loop massive QCD calculations. In A. Hoang and C. Schneider, editors,
PoS(RADCOR2017) 069 , pages 1–13, 2017.[6] J. Ablinger, J. Blümlein, A. de Freitas, and C. Schneider. A toolbox to solve coupled systemsof differential and difference equations.
PoS (RADCOR2015) 060 , 2015.[7] J. Ablinger, J. Blümlein, A. D. Freitas, A. Goedicke, C. Schneider, and K. Schönwald. Thetwo-mass contribution to the three-loop gluonic operator matrix element A ( ) gg , Q . Nucl. Phys.B , 932:129–240, 2018.[8] J. Ablinger, J. Blümlein, A. D. Freitas, A. Hasselhuhn, A. von Manteuffel, M. Round,C. Schneider, and F. Wißbrock. The transition matrix element A gq ( N ) of the variable fla-vor number scheme at O ( α s ) . Nucl. Phys. B , 882:263–288, 2014.[9] J. Ablinger, J. Blümlein, P. Marquard, N. Rana, and C. Schneider. Heavy quark form factorsat three loops in the planar limit.
Phys. Lett. B , 782:528–532, 2018.[10] J. Ablinger, J. Blümlein, C. Raab, and C. Schneider. Iterated binomial sums and their associ-ated iterated integrals.
J. Math. Phys. , 55(112301):1–57, 2014.[11] J. Ablinger, J. Blümlein, and C. Schneider. Harmonic sums and polylogarithms generated bycyclotomic polynomials.
J. Math. Phys. , 52(10):1–52, 2011.12 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider [12] J. Ablinger, J. Blümlein, and C. Schneider. Analytic and algorithmic aspects of generalizedharmonic sums and polylogarithms.
J. Math. Phys. , 54(8):1–74, 2013.[13] J. Ablinger and C. Schneider. Algebraic independence of sequences generated by (cyclo-tomic) harmonic sums.
Ann. Comb. , 22(2):213–244, 2018.[14] S. Abramov and M. Petkovšek. D’Alembertian solutions of linear differential and differenceequations. In J. von zur Gathen, editor,
Proc. ISSAC’94 , pages 169–174. ACM Press, 1994.[15] S. Abramov and M. Petkovšek. Polynomial ring automorphisms, rational ( w , σ ) -canonicalforms, and the assignment problem. J. Symbolic Comput. , 45(6):684–708, 2010.[16] G. Almkvist and D. Zeilberger. The method of differentiating under the integral sign.
J.Symbolic Comput. , 10(6):571–591, 1990.[17] M. Apagodu and D. Zeilberger. Multi-variable Zeilberger and Almkvist-Zeilberger algo-rithms and the sharpening of Wilf-Zeilberger theory.
Adv. Appl. Math. , 37:139–152, 2006.[18] A. Behring, J. Blümlein, G. Falcioni, A. D. Freitas, A. von Manteuffel, and C. Schneider.The asymptotic 3-loop heavy flavor corrections to the charged current structure functions F W + − W − L ( x , Q ) and F W + − W − ( x , Q ) . Phys. Rev. D , 94(11):1–19, 2016.[19] A. Behring, J. Blümlein, A. D. Freitas, A. Hasselhuhn, A. von Manteuffel, and C. Schneider.The O ( α s ) heavy flavor contributions to the charged current structure function xF ( x , Q ) atlarge momentum transfer. Phys. Rev. D , 92(114005):1–19, 2015.[20] A. Behring, J. Blümlein, A. D. Freitas, A. von Manteuffel, and C. Schneider. The 3-loopnon-singlet heavy flavor contributions to the structure function g ( x , Q ) at large momentumtransfer. Nucl. Phys. B , 897:612–644, 2015.[21] J. Blümlein. Structural relations of harmonic sums and Mellin transforms at weight w=6. InA. Carey, D. Ellwood, S. Paycha, and S. Rosenberg, editors,
Motives, Quantum Field Theory,and Pseudodifferential Operators , volume 12 of
Clay Mathematics Proceedings , pages 167–186, 2010.[22] J. Blümlein, A. Hasselhuhn, S. Klein, and C. Schneider. The O ( α s n f T F C A , F ) contributionsto the gluonic massive operator matrix elements. Nucl. Phys. B , 866:196–211, 2013.[23] J. Blümlein, A. Hasselhuhn, and C. Schneider. Evaluation of multi-sums for large scaleproblems. In
PoS(RADCOR2011)32 , pages 1–9, 2012.[24] J. Blümlein, S. Klein, C. Schneider, and F. Stan. A symbolic summation approach to Feynmanintegral calculus.
J. Symbolic Comput. , 47:1267–1289, 2012.[25] J. Blümlein and S. Kurth. Harmonic sums and Mellin transforms up to two-loop order.
Phys.Rev. , D60, 1999. 13 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider [26] J. Blümlein, M. Round, and C. Schneider. Refined holonomic summation algorithms inparticle physics. In C. Schneider and E. Zima, editors,
Advances in Computer Algebra.WWCA 2016. , volume 226 of
Springer Proceedings in Mathematics & Statistics , pages 51–91. Springer, 2018.[27] J. Blümlein and C. Schneider. The method of arbitrarily large moments to calculate singlescale processes in quantum field theory.
Phys. Lett. , B771:31–36, 2017.[28] J. Blümlein and C. Schneider. Analytic computing methods for precision calculations inquantum field theory.
International Journal of Modern Physics A (IJMPA) , 33(1830015):1–35, 2018.[29] A. Bostan, F. Chyzak, and É. de Panafieu. Complexity estimates for two uncoupling algo-rithms. In
Proc. ISSAC’13 , Boston, June 2013.[30] F. Brown. The massless higher-loop two-point function.
Commun. Math. Phys. , 287:925–958,2009.[31] K. Chetyrkin and F. Tkachov. Integration by parts: The algorithm to calculate beta functionsin 4 loops.
Nucl. Phys. B , 192:159–204, 1981.[32] F. Chyzak. An extension of Zeilberger’s fast algorithm to general holonomic functions.
Dis-crete Math. , 217:115–134, 2000.[33] M. Czakon. Automatized analytic continuation of Mellin-Barnes integrals.
Comput. Phys.Commun. , 175:559–571, 2006.[34] A. De Freitas, J. Blümlein, and C. Schneider. Recent symbolic summation methods to solvecoupled systems of differential and difference equations.
PoS (LL2014) 017 , 2014.[35] S. Gerhold. Uncoupling systems of linear ore operator equations. Master’s thesis, RISC, J.Kepler University Linz, 2 2002.[36] J. Gluza, K. Kajda, and T. Riemann. AMBRE: A Mathematica package for the constructionof Mellin-Barnes representations for Feynman integrals.
Comput. Phys. Commun. , 177:879–893, 2007.[37] P. Hendriks and M. Singer. Solving difference equations in finite terms.
J. Symbolic Comput. ,27(3):239–259, 1999.[38] J. Henn. Multiloop integrals in dimensional regularization made simple.
Phys. Rev. Lett. ,110:251601, 2013.[39] M. Karr. Summation in finite terms.
J. ACM , 28:305–350, 1981.[40] M. Kauers, M. Jaroschek, and F. Johansson. Ore polynomials in Sage. In J. Gutierrez,J. Schicho, and M. Weimann, editors,
Computer Algebra and Polynomials , Lecture Notes inComputer Science, pages 105–125, 2014. 14 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider [41] M. Kauers and P. Paule.
The concrete tetrahedron . Texts and Monographs in SymbolicComputation. SpringerWienNewYork, Vienna, 2011.[42] S. Laporta. High precision calculation of multiloop Feynman integrals by difference equa-tions.
Int. J. Mod. Phys. , A15:5087–5159, 2000.[43] S. Larin, F. Tkachov, and J. Vermaseren. The FORM version of Mincer. Technical ReportNIKHEF-H-91-18, NIKHEF, Netherlands, 1991.[44] R. Lee. Reducing differential equations for multiloop master integrals.
JHEP , 04:108, 2015.[45] A. van Manteuffel, E. Panzer, and R. Schabinger. A quasi-finite basis for multi-loop Feynmanintegrals.
JHEP , 02:120, 2015.[46] A. van Manteuffel and C. Studerus. Reduze 2 - distributed Feynman integral reduction. arXiv:1201.4330 [hep-ph] , 2012.[47] P. Marquard and D. Seidel. The crusher algorithm. unpublished.[48] S. Moch, P. Uwer, and S. Weinzierl. Nested sums, expansion of transcendental functions, andmultiscale multiloop integrals.
J. Math. Phys. , 6:3363–3386, 2002.[49] E. Ocansey and C. Schneider. Representing (q-)hypergeometric products and mixed versionsin difference rings. In C. Schneider and E. Zima, editors,
Advances in Computer Algebra.WWCA 2016. , volume 226 of
Springer Proceedings in Mathematics & Statistics , pages 175–213. Springer, 2018.[50] E. Panzer. Algorithms for the symbolic integration of hyperlogarithms with applications toFeynman integrals.
Comput. Phys. Commun. , 188:148–166, 2015.[51] M. Petkovšek, H. Wilf, and D. Zeilberger. A = B . A. K. Peters, Wellesley, MA, 1996.[52] E. Remiddi and J. Vermaseren. Harmonic polylogarithms. Int. J. Mod. Phys. , A15:725–754,2000.[53] C. Schneider. Symbolic summation in difference fields. Technical Report 01-17, RISC-Linz,J. Kepler University, November 2001. PhD Thesis.[54] C. Schneider. A new Sigma approach to multi-summation.
Adv. in Appl. Math. , 34(4):740–767, 2005.[55] C. Schneider. Product representations in ΠΣ -fields. Ann. Comb. , 9(1):75–99, 2005.[56] C. Schneider. Simplifying sums in ΠΣ -extensions. J. Algebra Appl. , 6(3):415–441, 2007.[57] C. Schneider. Symbolic summation assists combinatorics.
Sém. Lothar. Combin. , 56:1–36,2007. Article B56b.[58] C. Schneider. A refined difference field theory for symbolic summation.
J. Symbolic Comput. ,43(9):611–644, 2008. 15 omputer algebra tools for Feynman integrals and related multi-sums
Carsten Schnneider [59] C. Schneider. A symbolic summation approach to find optimal nested sum representations. InA. Carey, D. Ellwood, S. Paycha, and S. Rosenberg, editors,
Motives, Quantum Field Theory,and Pseudodifferential Operators , volume 12 of
Clay Mathematics Proceedings , pages 285–308. Amer. Math. Soc, 2010.[60] C. Schneider. Simplifying multiple sums in difference fields. In C. Schneider and J. Blümlein,editors,
Computer Algebra in Quantum Field Theory: Integration, Summation and SpecialFunctions , Texts and Monographs in Symbolic Computation, pages 325–360. Springer, 2013.[61] C. Schneider. Modern summation methods for loop integrals in quantum field theory: Thepackages Sigma, EvaluateMultiSums and SumProduction. In
Proc. ACAT 2013 , volume 523of
J. Phys.: Conf. Ser. , pages 1–17, 2014.[62] C. Schneider. Fast algorithms for refined parameterized telescoping in difference fields. InM. W. J. Guitierrez, J. Schicho, editor,
Computer Algebra and Polynomials , number 8942 inLecture Notes in Computer Science (LNCS), pages 157–191. Springer, 2015.[63] C. Schneider. A difference ring theory for symbolic summation.
J. Symb. Comput. , 72:82–127, 2016.[64] C. Schneider. Summation theory II: Characterizations of R ΠΣ -extensions and algorithmicaspects. J. Symb. Comput. , 80(3):616–664, 2017.[65] A. Smirnov and V. Smirnov. On the resolution of singularities of multiple Mellin-Barnesintegrals.
Eur. Phys. J. , C62:445–449, 2009.[66] M. Steinhauser. MATAD: A Program package for the computation of MAssive TADpoles.
Comput. Phys. Commun. , 134:335–364, 2001.[67] C. Studerus. Reduze-Feynman Integral Reduction in C++.
Comput. Phys. Commun. ,181:1293–1300, 2010.[68] J. Vermaseren. Axodraw.
Comput. Phys. Commun. , 83:45–58, 1994.[69] J. Vermaseren. Harmonic sums, Mellin transforms and integrals.
Int. J. Mod. Phys. ,A14:2037–2976, 1999.[70] S. Weinzierl. Feynman graphs. In C. Schneider and J. Blümlein, editors,
Computer Algebrain Quantum Field Theory: Integration, Summation and Special Functions , Texts and Mono-graphs in Symbolic Computation, pages 381–406. Springer, 2013.[71] B. Zürcher.