WUE-ITP-2005-014 6 0 Automatized analytic continuation of 0 2 Mellin-Barnes integrals n a J 3 2 M. Czakon ∗ 2 v Institut fu¨r Theoretische Physik und Astrophysik, Universita¨t Wu¨rzburg, 0 Am Hubland, D-97074 Wu¨rzburg, Germany 0 2 Department of Field Theory and Particle Physics, Institute of Physics, 1 University of Silesia, Uniwersytecka 4, PL-40007 Katowice, Poland 1 5 0 / h Abstract p - p e I describe a package written in MATHEMATICAthat automatizes typical operations h performed during evaluation of Feynman graphs with Mellin-Barnes (MB) tech- : v niques.Themain procedureallows to analytically continue aMB integral in a given i X parameterwithoutanyinterventionfromtheuserandthustoresolvethesingularity r structure in this parameter. The package can also perform numerical integrations a at specified kinematic points, as long as the integrands have satisfactory conver- gence properties. I demonstrate that, at least in the case of massive graphs in the physical region, the convergence may turn out to be poor, making na¨ive numerical integration of MB integrals unusable. I present possible solutions to this problem, but argue that full automatization in such cases may not be achievable. e-mail: [email protected] ∗ Preprint submitted to Elsevier Science 2 February 2008 PROGRAM SUMMARY Title of program: MB Version: 1.1 Catalogue identifier: Programobtainablefrom:http://theorie.physik.uni-wuerzburg.de/~mczakon Computers: All Operating systems: All Programming language used: MATHEMATICA, Fortran 77 for numerical evalu- ation Memory required to execute with typical data: Sufficient for a typical installa- tion of MATHEMATICA. No. of bytes in distributed program, including test data: 337900 Distribution format: ASCII Libraries used: CUBA [1] for numerical evaluation of multidimensional integrals and CERNlib [2] for the implementation of Γ and ψ functions in Fortran. Keywords: Mellin-Barnes integrals, analytic continuation, numerical evalua- tion, Feynman integrals. Nature of physical problem: Analytic continuation of Mellin-Barnes integrals in a parameter and subsequent numerical evaluation . This is necessary for evaluation of Feynman integrals from Mellin-Barnes representations. Method of solution: Recursive accumulation of residue terms occurring when singularities cross integration contours. Numerical integration of multidimen- sional integrals with the help of the CUBA library. Restrictions on the complexity of the problem: Limited by the size of the avail- able storage space. Typical running time: Depending on the problem. Usually seconds for moder- ate dimensionality integrals. 2 1 Introduction The synergy between experiment and theory in the area of elementary parti- cle physics is constantly driving perturbative calculations to higher and higher orders. This is particularly true close to the beginning of the Large Hadron Collider’s operation. Therefore, recent years have seen the emergence of sev- eral powerful methods of evaluation of subsequent terms of the perturbative expansion. As far as multiloop Feynman integrals are concerned, the method of differential equations Ref. [3,4] and Mellin-Barnes integral representations Ref. [5,6] have proved to be the most successful. Some complicated problems turnedout toeven require amixed approach,asadvocated, forexample, inthe case of Bhabha scattering in Ref. [7]. In parallel to analytical approaches, new numeric techniques have been devised, among which the sector decomposition method Ref. [8] occupies a prominent place. Very recently in Ref. [9], the role ofMellin-Barnesintegralrepresentations assourcesofnumericapproximations in the physical region has also been stressed. In the present work, I will concentrate on Mellin-Barnes integral representa- tions.There aretwo advantagesofthisapproach.First,itallowsforsystematic extraction of singularities. Second, the dimensionality of the representation is not directly connected to the number of lines in the graph and therefore, one often arrives at integrals of low dimensionality even for complicated graphs. The calculation of a Feynman integral proceeds in this method in three steps. At first, one derives a representation, then performs the analytic continuation in ǫ, where d = 4 2ǫ is the dimension of spacetime, and finally evaluates the − resulting integrals. The first step above can be performed in several different ways, aiming at the simplest possible representation. The various possibilities are described in Ref. [10]. A general algorithm is here only interesting in the case of subsequent numeric integration, see Ref. [9]. The third step cannot be generalized apart from numerical integration, even though some classes of problems can be solved algorithmically, e.g. by reduction to nested sums, see Ref. [11]. It is only in the second step, the analytic continuation, that one can provide an algorithmic solution that would be satisfactory for both ana- lytic and numeric evaluation. This solution is provided by the MATHEMATICA package MB introduced in the present work. Numeric evaluation of MB integrals has already been mentioned more than once above. Whether just for testing or for the actual calculation, automatiza- tionofthisstepisofvaluebyitself.ThepackageMBcanperformthenecessary integration by means of FORTRAN, the CUBA library [1] of integration rou- tines, and the CERN library implementation of gamma and psi functions [2]. Since the integrals are infinite range and multidimensional, their feasibility depends strongly on their convergence. In all tested examples, where invari- ants are in the Euclidean range, the behaviour is exponential and therefore 3 poses no problems. In [9], physical kinematics have also been considered, but the presented examples were restricted to massless graphs exclusively. Here, I notice that massive graphs have worse properties. In fact, I give examples of integrals, which are not even absolutely integrable, and the integral is similar to the Fourier transform of the inverse square root. Such cases can still be treated, but some initial analysis is necessary and it is difficult to see how it could beautomatized. Moreover, the techniques will rapidlybecome inefficient for higher dimensional integrals. The paper is organized as follows. In the next section, I define the main con- cepts and present the algorithm for analytic continuation. Subsequently, I describe the package starting with the user interface, low level routines, ex- amples, numerical integration routines and some additional tools. Finally, I briefly summarize and conclude the paper. 2 Analytic continuation of Mellin-Barnes integrals At the core of the Mellin-Barnes method lies the following representation i∞ 1 1 1 Az = dz Γ( z)Γ(ν +z), (1) (A+B)ν Γ(ν)2πi Bν+z − −Zi∞ where the contour is chosen in such a way, that the poles of the Γ function with +z are separated from the poles of the Γ function with z. − This representation can be used in Feynman integral computations in several ways. The easiest is to turn massive propagators into massless and integrate the massless integral, if a formula for general powers of propagators exists. In more complicated cases, one can use some parametric representation of the Feynman integral, which is usually an integral of a product of polynomials raised to some powers, and split the polynomials into pieces that are then integrable by some generalization of the Euler formula 1 Γ(α)Γ(β) dx xα−1(1 x)β−1 = . (2) − Γ(α+β) Z 0 AnextensivediscussionofthemethodswithexamplescanbefoundinRef.[10]. Irrespective of the method, however, the expression for any Feynman integral assumes the form 4 i∞ i∞ 1 Π Γ(A +V +c ǫ) j j j j ... Π dz f(z ,...,z ,s ,...,s ,a ,...,a ,ǫ) , i i 1 n 1 p 1 q (2πi)n Π Γ(B +W +d ǫ) −Zi∞ −Zi∞ k k k k (3) where s are some kinematic parameters and masses; a are the powers of the i i propagators; A , B are linear combinations of the a ; V , W are linear com- i i i i i binations of z ; and c , d are some numbers. The function f is analytic, in i i i practice a product of powers of the s , with exponents being linear combina- i tions of the remaining parameters. Because of the assumptions inherent in Eq. (1), the above equation is well defined and corresponds to the original Feynman integral, if the real parts of all of the Γ functions have positive arguments. If these conditions cannot be satisfied with ǫ = 0, then the integral may develop divergences and analytic continuation to 0 is necessary to make an expansion in ǫ. The purpose of the presented package is to perform the analytic continuation of Eq. (3) in ǫ to some chosen value ǫ . The algorithm requires to generalize 0 Eq. (3) to allow for ψ functions in the fraction, with ψ(z) = dlogΓ(z)/dz and ψ(n)(z) = dnψ(z)/dzn, with the same structure of arguments as those of the Γ functions. 2.1 The algorithm There are two known ways to perform the analytic continuation. The first, introduced in Ref. [5] consists in deforming the integration contours and then shifting them past the poles of the Γ functions, which results in residue in- tegrals. It is not clear how to make this method algorithmic, although some attempts in the specific case of massless on-shell double boxes have been un- dertaken in [12]. The second method, introduced in Ref. [6] assumes fixed contours parallel to the imaginary axis, and the analytic continuation consists in accounting for pole crossings past the contours. As described in Ref. [6], this method is an algorithm. I make one modification with respect to the original, namely I assume that the contours are such that no two contours can be crossed simultaneously. This assumption can always be satisfied by infinitesimal shifts of one of the concerned contours. It should be clear from the above considerations, that the imaginary parts of the involved variables do not play any role. It is therefore assumed that z , a i i and ǫ are real. With z = (z ,...,z ), a = (a ,...,a ), and (I,ǫ ) some MB 1 n 1 q I integral with fixed contours and the value of ǫ fixed at ǫ , the algorithm can I 5 let z, a and ǫ be such that all the arguments of the Γ functions be positive. O (I,ǫ) , C ← { } ← ∅ while O = 6 ∅ do R ← ∅ for all (I,ǫ ) O I ∈ do C C I for←all ǫ′∪{[ǫ},ǫ ], such that there is a F(z,a,ǫ), where Γ(F(z,a,ǫ)) I and FI (∈z,aI,ǫ′0) = n, n N 0 ∈ I − ∈ ∪{ } do if ǫ′ = ǫ or ǫ′ = ǫ I I I 0 then halt contour starts or ends on a pole of Γ(F(z,a,ǫ)) fi let z′ F(z,a,ǫ) let z0′ ∈be such that F(z,a,ǫ)|z′=z0′ = −n, n ∈ N∪{0} let I = 1 i∞ dz′G(z′) 2πi −i∞ sz′ signRof the coefficient of z’ in F(z,a,ǫ) RI′ ←←sRz′ sig(nI(′,Fǫ′()z,a,ǫI)−F(z,a,ǫ0)) Resz′=z0′G(z′) ← ∪{ I } od od O R ← od Fig. 1. Analytic continuation algorithm. be formalized as in Fig. 1. The algorithm has been written for Γ functions, but one should add ψ func- tions, wherever Γ functions occur. Upon termination, the set C contains all the integrals following from the analytic continuation. It should be clear that it is the “if” clause that does not allow for crossings of two different contours at a time. A comment about the choice of the contours is in order. Even though all the choices are equivalent, one would like to have the smallest possible num- ber of contributions. An improvement implemented in the package is to first gather all the residue points, and then try to add additional constraints on the contours such that these residues would not occur. If some subset of these con- straints can be satisfied, then the number of residues will be reduced. This is notanalgorithmthatleadstoanabsoluteminimum ofthenumber ofresidues, it gives, however, at least some reduction of the number of contributions. 6 Finally, one should notice that the technique of Ref. [6] has been similarly formalized in Ref. [9]. 3 The package 3.1 User interface The main routine performing the analytic continuation is MBcontinue[integrand, limit, fixedVarRules, intVarRules , options] { } where the input arguments are integrand: any object accepted by MATHEMATICA. Notice that the singu- • larities are determined by analyzing Γ and ψ functions only. limit: a rule, x -> x0, which specifies at the same time the variable, x, in • which the analytic continuation is performed and the point, x0, which the user wants to reach. fixedVarRules: a list of rules giving the values of the real parts of the • variables, which are not integrated over. In particular, it must contain the starting value of the variable, in which the analytic continuation is per- formed. intVarRules:alistofrulesgivingtherealpartsoftheintegrationvariables. • options: • Level: an integer specifying the level at which the recursive analytic · continuation will be stopped. By default, it is set to infinity. Skeleton: a boolean value. If True, the residues will be identified, but · not calculated. This is achieved by replacing all Γ and ψ functions by a dummy functionMBgam. The purpose ofthis optionisto quickly determine the total number of integrals. By default this option is set to False. Residues: a boolean value. If True, the output will also contain the list · of Residue points besides the actual values of the residues. This is mainly for internal use and is set by default to False. Verbose: a boolean value. If True, the level is printed as well as the · position on the list of the currently continued integral and the residue points together with the signs of the residues. This option is switched on by default. The output is a nested list obtained by replacing, at every level, the integral to be continued by its residues and the original integral at the limit. The elements are 7 MBint[integrand, fixedVarRules, intVarRules ] { } objects,wheretheintegrandcanbeexpandedaroundthelimit,whichisplaced on the fixedVarRules list. If the user specified a finite level, then there might also occur MBitc[integrand, limit, fixedVarRules, intVarRules , Options] { } objects,where“itc”standsfor“integraltocontinue”.Thesearenotyetregular at the limit and require further recursive analytic continuation. Furthermore, if the user set the Residues option to True, there will also be a list of MBres[sign, var, val], objects, which signal that there was a residue taken in the variable, var, at the value, val, with sign. Restricted input checking has been implemented, and as long as the input is syntactically correct, the only error that may occur is (see Section 2.1 for further details) contour starts and/or ends on a pole of Gamma[z] In this case the procedure stops and gives an inequality for an integration variable that is sufficient to remove the problem. The integration contours are found with MBoptimizedRules[integrand, limit, constraints, fixedVars, options] For a description of the integrand and limit see MBcontinue. The remaining input parameters are as follows constraints: a list of additional constraints (inequalities) specified by the • user. This should usually be left empty, but might be used for experimen- tation in order to search for contours that might possibly give less residues. fixedVars: a list of variables, which should be considered fixed during an- • alytic continuation. The integration variables are determined automatically from the arguments of the Γ and ψ functions. options: • Level: specifies the level up to which optimization of the contours will · be performed. This option should only be used for very large calculations. Since in this case, the contours are only partially tested, the user will have to correct them himself, if poles lying on a contour are encountered. In practice, independent, small shifts should be sufficient for this purpose. 8 The output matches precisely the form needed in the input of MBcontinue, i.e. fixedVarRules, intVarRules { } Notice that this procedure not only reduces the number of residues, but also generates such contours that, during analytic continuation, no contours will start or end on a pole. During the determination of the real parts, warning messages are generated. These can be ignored apart from the case when there is a single message no rules could be found to regulate this integral and the output is an empty list. In this case, the integral cannot be regulated andtheuserhastoprovide anotherone,e.g.byintroducing afurtherregulator parameter, for example a propagator power, and performing two subsequent analytic continuations. Once the integrals are determined, they can be either merged, i.e. those that have the same contour will be added by linearity; preselected, i.e. those that would vanish in a given order of expansion in some parameter are rejected; or expanded. These tasks are achieved with the following utilities. MBmerge[integrals] Merges MBint objects on the integrals list by linearity, if they have the same contours. Vanishing integrals are rejected. MBpreselect[integrals, x, x0, n ] { } Rejects those MBint objects on the integrals list that would vanish after expansion in the variable x, around the point x0, up to order n. MBexpand[integrals, norm, x, x0, n ] { } Expands MBint objects on the integrals list around the point x0, in the variable x, up to order n. A normalization factor, norm, is included in every integrand. 9 m p m Fig. 2. The B (s,ms,ms) function, with s = p2 and ms = m2. 0 3.2 Low level routines The routines described in the previous section form the interface. It might happen that the user would like to use the low level routines, which actually perform the calculation. MBresidues[integrand, limit, fixedVarRules, intVarRules , options] { } Performs a single step in the recursive analytic continuation algorithm, i.e. it findsalltheresiduesforagivenintegral,butdoesnotproceedwiththeanalytic continuation of the resulting integrals. All the arguments and options are the same as in MBcontinue, apart from Level, which is in this case meaningless. MBrules[integrand, constraints, fixedVars] Finds the real parts of all the fixed and integration variables, such that the real parts of the arguments of all the Γ and ψ functions be positive. The difference to MBoptimizedRules is that no attempt is made to optimize the number of residues or even check whether the contours will not lead to prob- lems with MBcontinue. To perform these tests, MBoptimizedRules needs the limit of the continuation, which is left unspecified here. This routine is of par- ticular interest, because one may use it to write another contour optimization algorithm. MBrules[integrand, limit, constraints, fixedVars] Same as MBrules, but check the contours, so that a complete analytic contin- uation with MBcontinue can be performed. 3.3 Examples As a first example, I consider the B function with two equal masses, Fig. 2. 0 After introduction of two MB integrations (the integral can be further simpli- fied by the use of the first Barnes lemma, see Section 3.5) and normalization 10