On the Resolution of Singularities of Multiple Mellin–Barnes Integrals A.V. Smirnov1 Scientific Research Computing Center of Moscow State University and 9 0 V.A. Smirnov2 0 2 Nuclear Physics Institute of Moscow State University n a J Abstract 4 ] h One of the two existing strategies of resolving singularities of multifold p - Mellin–Barnes integrals in the dimensional regularization parameter, or p e a parameter of the analytic regularization, is formulated in a modified h form. The corresponding algorithm is implemented as a Mathematica [ code MBresolve.m 1 v 6 8 3 0 . 1 0 9 0 : v i X r a 1E-mail: [email protected] 2E-mail: [email protected] 1 Introduction The method of Mellin–Barnes (MB) representation is one of the most powerful meth- ods to evaluate multiloop Feynman integrals. It is based on a very simple formula, 1 1 1 +i∞ Yz = dzΓ(λ+z)Γ(−z) , (1) (X +Y)λ Γ(λ)2πi Xλ+z Z−i∞ which is called MB representation and is applied to replace a sum of two terms raised to some power by their products in some powers. The contour of integration should be chosen so that the poles with a Γ(...+z) dependence are to the left of the contour and the poles with a Γ(...−z) dependence are to the right of it. The simplest way to apply (1) is to represent massive propagators as continuous superpositions of massless ones, at the cost of introducing additionalintegrations over MB variables. More often, MB integrals are introduced at the level of Feynman or alpha parameters. Anyway, (1) is used in an appropriate way, with the goal to obtain integrals over loop momenta or Feynman/alpha parameters which can be taken in terms of gamma functions. As a result one obtains a multiple MB integral with an integrand involving gamma functions and powers of kinematic invariants. ItisimportanttoderiveaMBrepresentationforgeneralvaluesofthepowersofthe propagators. First, this provides unambiguous prescriptions for a choice of contours of integration. The resulting rule is similar to the rule for (1): when we understand a multiple MB integral iteratively and consider an integration over a variable z it is i implied that the poles with a Γ(...+z ) dependence are to the left of the contour and i the poles with a Γ(...−z ) dependence are to the right of it. Second, such a general i MB representation provides crucial checks and can be applied at any chosen values of the powers of the propagators. When evaluating a Feynman integral for specific powers of the propagators (in- dices), one starts from the general MB representation and obtains an integral of the form 1 +i∞ +i∞ iΓ ai +biε+ j cijzj n ... xdk dz , (2) (2πi)n Z−i∞ Z−i∞ QiΓ(cid:16)a′i +b′iε+Pj c′ijzj(cid:17) Yk k Yl=1 l (cid:16) (cid:17) Q P ′ whereε = (4−d)/2isthedimensionalregularizationparameter,a ,...,c arerational i ij numbers, x are ratios of kinematic invariants and/or masses, and their exponents, k d , are linear combinations of ε and z-variables. Typically, c = ±1. k ij Although results for dimensionally regularized Feynman integrals are practically needed in a Laurent expansion in ε, one could try to evaluate a given MB integral (2) for general ε. Then the ε-expansion is performed in the result. However, for sufficientlycomplicatedmultiloopFeynmanintegrals,theevaluationatgeneralεturns out to be impossible, and one proceeds with an ε-expansion. 1 One of the advantages of this method is that the singularity structure in ε can be resolved in a simple way. This procedure basically consists of taking residues and shifting contours, with the goal to obtain a sum of integrals where one can expand integrands in Laurent series in ε. To do this one can apply two strategies formulated in [1] and [2] which we are going to call Strategy A and Strategy B, respectively. Strategy A is described and illustrated in numerous examples in Chapter 4 of [3]. It was applied, e.g., in [4]. According to strategy A, one performs an analysis of the integrand to reveal how poles in ε arise. The guiding principle is that the product Γ(a + z)Γ(b − z), where a and b can depend on the rest of the integration variables, generates, due to the integration over z, the singularity of the type Γ(a+ b). Indeed, if we shift an initial contour of integration over z across the point z = −a we obtain an integral over a new contour which is not singular at a + b = 0, while the corresponding residue involves an explicit factor Γ(a+b). This observation shows that any contour of one of the following integrations over the rest of the MB variables should be chosen according to this dependence, Γ(a+b). Hence one thinks of integrationsinvarious ordersandthen identifies some ‘key’ gammafunctions which arecrucialforthegenerationofpolesinε. Thenonetakesresiduesandshiftscontours, starting from the first poles of these key gamma functions. The same analysis and procedure is applied to the contributions of the residues. Within Strategy B, one chooses an initial value of ε and values of the real parts of the integration variables, z ,... in such a way that the real parts of all the arguments i of the gamma functions in the numerator are positive and one can integrate over straight lines. Then one tends ε to zero and whenever the real part of the argument of some gamma function vanishes one crosses this pole and adds a corresponding residue which has one integration less and is treated as the initial integral within the same procedure. StrategyBisalgorithmicinitscharacter and, indeed, twoalgorithmicdescriptions were formulated in [5, 6]. A public code called MB.m was presented in [6]. Strategy B was successfully applied, e.g., in [7, 8] and in many other papers. The purpose of this letter is to formulate Strategy A in a slightly modified form and to present the corresponding algorithm implemented in Mathematica. 2 The modified Strategy A To present a modified Strategy A let us explicitly formulate what was implied in the initial Strategy A. When we take care of one of the key gamma functions we shift a contour and take a residue. Let Γ(A ) with A = a + b ε + c z be one of the i i i i j ij j key gamma functions in (2). Without loss of generality, we can consider ε real. Then P changing the nature of the first pole of this gamma function means changing the rule for an admissible contour, i.e. that, instead of the condition ReA > 0 when crossing i the real axis in the process of the integration, we have the condition −1 < ReA < 0. i Let us denote this transition by replacing Γ(A ) by Γ(1)(A ). The initial rule for the i i 2 contour can be changed again and then we have the condition −n < ReA < −n+1 i for n = 2,3,... with the notation Γ(n)(A ). i Within Strategy B one has straight contours in the beginning. Rather, in the modified Strategy A, we will be oriented at straight contours in the end. Apparently, it is desirable to achieve a minimal number of terms after the resolution of the singu- larities in ε. To do this, let us try to search for contours which are going to have in the end of this procedure and for which the gamma functions in the numerator are changed, in the above sense, in a minimal way. To formalize this requirement, let us introduce the function σ(x) = [(1 − x) ] + where [...] is the integer part of a number and x = x for x > 0 and 0, otherwise. In + other words, if −n < x < −n+1 then σ(x) = n for n > 0 and σ(x) = 0 for n ≤ 0. So, let us set ε = 0 and search for contours, i.e. Rez , for which the sum i σ(ReA | ) ≡ σ a + c Rez i ε=0 i ij j ! i i j X X X is minimal. After such a choice is done we identify gamma functions which should be changed, in the above sense, in order to arrive at a final integral where a Laurent expansion in ε is possible. In fact, this step replaces the first step in the primary Strategy A where one identified such key gamma functions after the analysis characterized above. Then the second step in Strategy A is the same as in the old version: we take care of the distinguished gamma functions, i.e. take a residue and replace Γ by Γ(1)(A ) i (and, possibly, Γ(1)(A ) by Γ(2)(A ) etc.) We proceed iteratively, as in the previous i i strategy: every residue is considered from scratch, i.e. treated in the same way as the initial MB integral. Let us emphasize that although this strategy aims to minimize the number of resulting terms, we cannot exclude that there is another way of resolving the singu- larities in ε that is the best one in this sense. For example, it can happen, in rather complicated examples, that different orders of changing the key gamma functions lead to different numbers of resulting terms. Still we believe that such a difference is negligible and that the Strategy A provides a resolution of the singularities in ε at least very close to the theoretically best one. The difference of the new and the old Strategies A is minor. In the examples of [4] one can see that resulting contours were straight indeed. This difference can still be seen in the following simple example of the integral 1 +i∞ Γ(1+ε+z)Γ(−1/2+ε+z)Γ(3/2−ε−z)Γ(−z) dz 2πi Z−i∞ which, of course, can be evaluated at general ε by the first Barnes lemma. Within the ‘old’ Strategy A, one can observe that there is no gluing of poles so that one can expand the integrand in ε. However, a resulting contour cannot be chosen as 3 a straight line. Rather, within the ‘new’ Strategy A, we have to choose a gamma function to be modified and, as a result, we obtain a residue and an integral over a straight line which both can be expanded in ε. 3 Implementation in Mathematica The code should be loaded together with the package MB [6] by Czakon. It uses the MBresidue routine and some other private functions from that package; after obtaining the list of integrals one can continue the contour optimization and the numerical evaluation with the standard routines of MB. Formally the code can be described the following way: the initial function results in contour prescriptions (the arguments of Γ-functions have to be positive in the end). Now MBresolve[function,cont pr] 1. find a point with the minimal number of contour prescriptions broken. 2. If this number is equal to zero, then return function as an answer 3. l = list of broken contour prescriptions 4. For each {x,n} in l do 5. Try to evaluate t1 = MBresolve[function,cont pr′] where cont pr′ are obtained by replacing {x,n} with {x,n+1} 6. On exception move to the next For cycle 7. represent x as cα+r, where c 6= 0, α is one of the integration variables and r does not depend on α 8. Return t1 - Sign[c]MBresolve[MBresidue[function,{α,−r/c}],cont pr′′], where cont pr′′ are inherited from cont pr by removing {x,n} and performing the substitution α → −r/c in all other functions 9. End 10.Throw exception The function should be a function of integration variables such that the poles are only due to the Gamma and PolyGamma factors (with linear arguments). The contour prescriptions cont pr are a list of pairs x,n where x is a linear function of integration variables and n is an non-negative integer. If n is equal to zero, such a term means that x has to be positive; if n is positive, then x has to be greater that −n and smaller than −n+1. If the cont pr parameter is missing, it is created automatically by considering all arguments of Gamma functions and their derivatives and assuming them all to be positive. An uncaught exception at the top level appears if there is a degenerate case and one is required to introduce extra regularization parameters. Hence, instead of using MBoptimizedRules[] or MBrules[] and, subsequently, MBcontinue[], within MB.m [6], one can now apply MBresolve[]. 4 The corresponding Mathematica code MBresolve.m is public and can be found at http://projects.hepforge.org/mbtools/ together with other tools for evaluating MB integrals. Here is an example of a tenfold MB representation derived loop by loop for the four-loop ladder massless on-shell diagram with p2 = 0, i = 1,2,3,4, where s = i (p +p )2 and t = (p +p )2. 1 2 1 3 In[1]: << MB.m; << MBresolve.m; F = -(((-s)^(-5 - 4*ep - z7)*(-t)^z7*Gamma[1 + z1]* Gamma[-1 - ep - z1 - z2]*Gamma[-z2]*Gamma[-1 - ep - z1 - z3]* Gamma[-z3]*Gamma[1 + z1 + z2 + z3]*Gamma[2 + ep + z1 + z2 + z3]* Gamma[z10 - z4]*Gamma[-z1 + z4]*Gamma[-ep + z1 + z2 - z4 - z5]* Gamma[-z5]*Gamma[-ep + z1 + z3 - z4 - z6]*Gamma[-z6]* Gamma[1 + z4 + z5 + z6]* Gamma[1 + ep - z1 - z2 - z3 + z4 + z5 + z6]*Gamma[-z7]* Gamma[1 + z7]*Gamma[-z10 + z7]*Gamma[-ep - z10 + z4 + z5 - z8]* Gamma[-z8]*Gamma[-ep + z10 - z7 + z8]* Gamma[-ep - z10 + z4 + z6 - z9]* Gamma[1 + ep - z10 + z7 - z8 - z9]*Gamma[-z9]* Gamma[-ep + z10 - z7 + z9]*Gamma[1 + z10 + z8 + z9]* Gamma[1 + ep + z10 - z4 - z5 - z6 + z8 + z9])/(Gamma[-2*ep]* Gamma[1 - z2]*Gamma[1 - z3]*Gamma[1 - 2*ep + z1 + z2 + z3]* Gamma[1 - z5]*Gamma[1 - z6]*Gamma[1 - 2*ep + z4 + z5 + z6]* Gamma[1 - z8]*Gamma[1 - z9]*Gamma[1 - 2*ep + z10 + z8 + z9])); Fcont = MBresolve[F, ep]; Length[Fcont] MB 1.2 by Michal Czakon improvements by Alexander Smirnov more info in hep-ph/0511200 last modified 2 Jan 09 MBresolve 1.0 by Alexander Smirnov more info in arXiv:0901.0386 last modified 4 Jan 09 CREATING RESIDUES LIST..........653.5156 seconds EVALUATING RESIDUES..........15.2969 seconds Out[1]: 656 5 4 Discussion and perspectives Let us remind that, sometimes, a given MB representation can be ill-defined even for a well-defined Feynman integral, in the sense that a gluing of poles is present, i.e. one can distinguish a subproduct of the gamma functions such that the sum of there arguments, at general ε, is equal to a non-positive integer number. To cure such a MB representation, one can introduce an auxiliary analytic regularization into some index, i.e. a → a +y and then analytically continue the given MB integral, first, in y i i to the point y = 0 and then ε to the point ε = 0. Let us stress that in such situation one starts, within the new Strategy A, with setting y = 0,ε = 0 and then searches for appropriate contours for which the gamma functions in the numerators are changed in the minimal way, similarly to the case without such analytic regularization. Let us emphasize that both Strategy A and Strategy B are based on the fact that singularities in ε can be generated, because of gluing of poles of different nature, by the integration at compact regions. At least for planar diagrams, for which one can applythecodecalledAMBRE[9]basedontheloop-by-loopstrategy, thislookstobethe only source of the singularities. However, for MB representations derived within the loop-by-loop strategy for nonplanar diagrams, there can be another source of poles.3 This feature can be exemplified by the massless two-loop nonplanar diagram with two external legs on-shell.4 In the corresponding MB representation derived within the loop-by-loop strategy one meets, in particular, the following onefold MB integral 1 +i∞ Γ(1+2ε+z)Γ(−z) e−iπzdz . 2πi 1+ε+z Z−i∞ There is no gluing of poles so that a pole in ε cannot be generated by the integration over finite regions. Still a pole is generated and this can be seen by an explicit evalu- ation of this integral by closing the integration contour to the right and summing up the resulting series. This can be seen also by analyzing the asymptotic behaviour of the integrand at infinity. Setting z = x+iy and using formulae of the asymptotic be- haviour of gamma functions at large arguments in the complex plane, we can observe that the leading asymptotic behaviour when y → +∞ is 1/y1−2ε which explains the appearance of the pole. Let us still mention that, for this concrete diagram, there is a better way to obtain a ‘good’ MB representation — see, e.g., Chapter 4 of [3]. For simple MB integrals, say, up to six-fold, the codes MB.m and MBresolve.m work more or less in the same way. For higher-fold MB integrals, e.g. ten-fold ones, MBresolve.m produces at first fewer integrals after the resolution of singularities in ε. However, the command MBmerge of MB.m combines many terms, and the number of resulting integrals becomes smaller that produced by MBresolve.m. We should bear in mind, though, that in the current version of MB.m, this command combines 3See also a similar discussion in [8]. 4It was first calculated in an expansion in ε up to ε0 in [10]. Now, results of expansion up to ε2 [11], up to ε4 [12] and even a result for general ε [12] are available. 6 all integrals which have the same integration contours, while it is more natural to consider separately integrals with different patterns of gamma functions, both from analytical and numerical points of view. Presumably, the MBmerge command can be improved in this respect. Acknowledgments. We are grateful to M. Czakon and D. Kosower for helpful discussions and to M. Czakon for careful reading of a draft version of the paper. Many thanks to N. Glover for kind hospitality during our visit to the Institute of Particle Physics Phenomenology in Durham, where a part of this work was done. The work was supported by the Russian Foundation for Basic Research through grant 08-02-01451. References [1] V.A. Smirnov, Analytical result for dimensionally regularized massless on-shell double box, Phys. Lett. B 460 (1999) 397 [arXiv:hep-ph/9905323]. [2] J.B. Tausk, Non-planar massless two-loop Feynman diagrams with four on-shell legs, Phys. Lett. B 469 (1999) 225 [arXiv:hep-ph/9909506]. [3] V.A. Smirnov, Evaluating Feynman Integrals, Springer Tracts Mod. Phys. 211 (2004) 1; V.A. Smirnov, Feynman integral calculus, (Berlin, Germany, Springer, 2006). [4] V.A. Smirnov and O.L. Veretin, Analytical results for dimensionally regularized massless on-shell double boxes with arbitrary indices and numerators, Nucl. Phys. B 566 (2000) 469 [arXiv:hep-ph/9907385]; V.A. Smirnov, Analytical Result for Dimensionally Regularized Massless Mas- ter Double Box with One Leg off Shell, Phys. Lett. B 491 (2000) 130 [arXiv:hep-ph/0007032]; Analytical result for dimensionally regularized mass- less master non-planar double box with one leg off shell, Phys. Lett. B 500, (2001) 330 [arXiv:hep-ph/0011056]; The leading power Regge asymptotic be- haviour of dimensionally regularized massless on-shell planar triple box, Phys. Lett. B 547 (2002) 239 [arXiv:hep-ph/0209193]; Analytical result for dimension- ally regularized massless on-shell planar triple box, Phys. Lett. B 567 (2003) 193 [arXiv:hep-ph/0305142]; Analytical result for dimensionally regularized massive on-shellplanardouble box, Phys. Lett.B 524(2002)129[arXiv:hep-ph/0111160]; Evaluating multiloop Feynman integrals by Mellin–Barnes representation, Nucl. Phys. Proc. Suppl. 135 (2004) 252 [arXiv:hep-ph/0406052]; G. Heinrich and V.A. Smirnov, Analytical evaluation of dimensionally reg- ularized massive on-shell double boxes, Phys. Lett. B 598 (2004) 55 [arXiv:hep-ph/0406053]; Z. Bern, L.J. Dixon, and V.A. Smirnov, Iteration of planar amplitudes in max- imally supersymmetric Yang–Mills theory at three loops and beyond, Phys. Rev. 7 D 72 (2005) 085001 [arXiv:hep-th/0505205]; A.G. Grozin, A.V. Smirnov and V.A. Smirnov, Decoupling of heavy quarks in HQET, JHEP 0611 (2006) 022 [arXiv:hep-ph/0609280]; B. Jantzen and V.A. Smirnov, The two-loop vector form factor in the Sudakov limit, Eur. Phys. J. C 47 (2006) 671 [arXiv:hep-ph/0603133]; J.M. Drummond, J. Henn, V.A. Smirnov and E. Sokatchev, Magic identities for conformal four-point integrals, JHEP 0701 (2007) 064 [arXiv:hep-th/0607160]. [5] C. Anastasiou and A. Daleo, Numerical evaluation of loop integrals, JHEP 0610 (2006) 031 [arXiv:hep-ph/0511176]. [6] M. Czakon, Automatized analytic continuation of Mellin–Barnes integrals, Com- put. Phys. Commun. 175 (2006) 559 [arXiv:hep-ph/0511200]. [7] E.W.N. Glover and M.E. Tejeda-Yeomans, Progress towards 2 to 2 scattering at two-loops, Nucl. Phys. Proc. Suppl. 89 (2000) 196 [arXiv:hep-ph/0010031]; C. Anastasiou, J.B. Tausk and M.E. Tejeda-Yeomans, The on-shell massless pla- nar double box diagram with an irreducible numerator, Nucl. Phys. Proc. Suppl. 89 (2000) 262 [arXiv:hep-ph/0005328]; C. Anastasiou, T. Gehrmann, C. Oleari, E. Remiddi and J.B. Tausk, The tensor reduction and master integrals of the two-loop massless crossed box with light-like legs, Nucl. Phys. B 580 (2000) 577 [arXiv:hep-ph/0003261]; M. Czakon, J. Gluza and T. Riemann, On the massive two-loop corrections to Bhabha scattering, Acta Phys. Polon. B 36 (2005) 3319 [arXiv:hep-ph/0511187]; T. Becher and M. Neubert, Toward a NNLO calculation of the B¯ → X γ decay s rate with a cut on photon energy. II: Two-loop result for the jet function, Phys. Lett. B 637 (2006) 251 [arXiv:hep-ph/0603140]; Z. Bern, M. Czakon, D.A. Kosower, R. Roiban and V.A. Smirnov, Two-loop iteration of five-point N = 4 super-Yang–Mills amplitudes, Phys. Rev. Lett. 97 (2006) 181601 [arXiv:hep-th/0604074]; Z. Bern, M. Czakon, L.J. Dixon, D.A. Kosower and V.A. Smirnov, The four-loop planar amplitude and cusp anomalous dimension in maximally supersymmetric Yang–Mills theory, Phys. Rev. D 75 (2007) 085010 [arXiv:hep-th/0610248]; T. Gehrmann, G. Heinrich, T. Huber and C. Studerus, Master integrals for massless three-loop form factors: One-loop and two-loop insertions, Phys. Lett. B 640 (2006) 252 [arXiv:hep-ph/0607185]; G. Heinrich, T. Huber and D. Maitre, Master integrals for fermionic contributions to massless three-loop form factors, Phys. Lett. B 662 (2008) 344 [arXiv:0711.3590 [hep-ph]]; A.V. Smirnov, V.A. Smirnov and M. Steinhauser, Applying Mellin–Barnes tech- nique and Gr¨obner bases to the three-loop static potential, PoS RADCOR2007 (2007) 024 [arXiv:0805.1871 [hep-ph]]; Evaluating the three-loop static quark po- tential, Nucl. Phys. Proc. Suppl. 183 (2008) 308 [arXiv:0807.0365 [hep-ph]]; M.B. Green, J.G. Russo and P. Vanhove, Modular properties of two-loop max- 8 imal supergravity and connections with string theory, JHEP 0807 (2008) 126 [arXiv:0807.0389 [hep-th]]; M. Czakon, A. Mitov and S. Moch, Heavy-quark production in massless quark scattering at two loops in QCD, Phys. Lett. B 651 (2007) 147 [arXiv:0705.1975 [hep-ph]]; M. Czakon, Tops from light quarks: full mass dependence at two-loops in QCD, Phys. Lett. B 664 (2008) 307 [arXiv:0803.1400 [hep-ph]]; M.Beneke, T. Huber and X.Q.Li, Two-loop QCD correction to differential semi- leptonic b → u decays in the shape-function region, [arXiv:0810.1230 [hep-ph]]; A. Pak and A. Czarnecki, Heavy-to-heavy quark decays at NNLO, [arXiv:0808.3509 [hep-ph]]; I. Bierenbaum, J. Blumlein, S. Klein and C. Schneider, Two-loop massive opera- tor matrix elements for unpolarized heavy flavor production to O(ǫ), Nucl. Phys. B 803 (2008) 1 [arXiv:0803.0273 [hep-ph]]; G.SomogyiandZ.Trocsanyi, A subtraction scheme for computing QCD jet cross sections at NNLO: integrating the subtraction terms I, JHEP 0808 (2008) 042 [arXiv:0807.0509 [hep-ph]]. [8] M. Czakon, A. Mitov and S. Moch, Heavy-quark production in gluon fusion at two loops in QCD, Nucl. Phys. B 798 (2008) 210 [arXiv:0707.4139 [hep-ph]]. [9] J. Gluza, K. Kajda and T. Riemann, AMBRE - a Mathematica package for the construction of Mellin–Barnes representations for Feynman integrals, Comput. Phys. Commun. 177 (2007) 879 [arXiv:0704.2423 [hep-ph]]. [10] R.J. Gonsalves, Dimensionally regularized two loop on-shell quark form-factor, Phys. Rev. D 28 (1983) 1542. [11] S. Moch, J.A.M. Vermaseren and A. Vogt, The quark form factor at higher orders, JHEP 0508 (2005) 049 [arXiv:hep-ph/0507039]. [12] T. Gehrmann, T. Huber and D. Maitre, Two-loop quark and gluon form factors in dimensional regularisation, Phys. Lett. B 622 (2005) 295 [arXiv:hep-ph/0507061]. 9