and

The method of Mellin-Barnes representation is used to calculate dimensionally regularized massive on-shell double box Feynman diagrams contributing to Bhabha scattering at two loops.


Introduction
The method of Mellin-Barnes representation has turned out to be very successful for the analytical evaluation of two-loop Feynman diagrams with four external lines within dimensional regularization [1]. The problem of the evaluation of massless on-shell double boxes was solved in [2][3][4][5], with multiple subsequent applications to the calculation of two-loop scattering amplitudes in gauge theories [6]. First calculations of massless on-shell fourpoint three-loop Feynman integrals were done in [7].
Complete algorithms for the evaluation of massless on-shell double boxes with one leg off shell were also constructed. Master integrals were calculated using the MB representation (first results in [8]) and the method of differential equations (DE) [9], where a systematic evaluation was described in [10]. The reduction to master integrals was done using a Laporta's idea [11] in [10]. All results are expressed in terms of two-dimensional harmonic polylogarithms [10] which generalize harmonic polylogarithms (HPL) [12]. Let us stress that this very combination of reduction based on [10,11] and DE was successfully applied in numerous calculations, e.g., to various classes of vertex diagrams [13,14]. E-mail address: vsmirnov@mail.desy.de (V.A. Smirnov). 0370-26930370-  20040370- Elsevier B.V. doi:10.10160370- /j.physletb.2004 Open access under CC BY license.
Open access under CC BY license. Another important class of Feynman integrals with one more parameter, with respect to the massless on-shell double boxes, are massive on-shell double boxes relevant to Bhabha scattering at NNLO. Bhabha scattering plays an important role as a luminosity monitor for e + e − colliders. Because of this phenomenological interest, the twoloop corrections to Bhabha scattering have been studied extensively in recent years. In [15], the complete two-loop virtual QED corrections were calculated in the massless fermion approximation. The infrared singularity structure of this result was studied in [16]. The first result for a two-loop master double box with nonzero electron mass was obtained in [17]. Planar two-loop box diagrams with one-loop insertion and without neglecting the fermion masses were calculated by DE in [18]. This work was continued to calculate the subset of all two-loop diagrams involving a closed fermion loop in [19]. The two-loop contribution which factorizes into squares of one-loop diagrams was calculated in [20]. Numerical results for the three master integrals shown in Fig. 1 were given in [21] for Euclidean points. Despite all these efforts, a complete two-loop differential cross section is not available yet [22], mainly because the most complicated master integrals, the two-loop box diagrams, have not all been calculated so far.

Planar double boxes of the first type
The general double box Feynman integral of the first type (see Fig. 1(a)) takes the form B PL,1 a 1 , . . . , a 8 ; s, t, m 2 ; where s = (p 1 + p 2 ) 2 , t = (p 2 + p 3 ) 2 , and k and l are the loop momenta corresponding to the left and the right box. Usual prescriptions k 2 = k 2 + i0, s = s + i0, etc., are implied. For convenience, we consider the factor with (k + p 1 + p 2 + p 3 ) 2 corresponding to the irreducible numerator as an extra propagator but, really, we are interested only in the non-positive integer values of a 8 . A first analytical result for the planar double box of the first type, was obtained in [17] using the following sixfold MB representation: ,4,5,6,7 (a j where a = a 1···8 , a 4567 = a 4 + a 5 + a 6 + a 7 , a 13 = a 1 + a 3 , etc., and integration contours are chosen in the standard way. To evaluate the master double box B PL,1 (1, . . ., 1, 0; s, t, m 2 ; ) the standard procedure of taking residues and shifting contours [2] was used [17], with the goal to obtain a sum of integrals where one may expand the integrands in a Laurent series in = (4 − d)/2 (where d is the space-time dimension within dimensional regularization [1]). Then one can use the first and the second Barnes lemmas and their corollaries to perform some of the MB integrations explicitly. In the last integrations which usually carry a dependence on the external variables, one closes the contour in the complex plane and sums up the corresponding series. (See [24] for details of the method.) In [22] it was announced that the massive on-shell double boxes with two reduced lines have been calculated by DE, with a solution of reduction problem by the method of [10,11]. However, an extension of these results to all master integrals meets problems. Differential equations of third order and higher are encountered there.
On the other hand, the method of MB representation can certainly be applied to arbitrary massive on-shell double boxes. To illustrate this point, we present the evaluation of typical complicated master integrals from this class.
As has been demonstrated in [3], it can be useful to use a double box with a numerator as a master integral in order to avoid coefficients of order 1/ in the course of the tensor reduction when calculating a scattering amplitude. Therefore, we present here results for the planar double box B PL,1 (1, . . . , 1, −1; s, t, m 2 ; ).
Our calculation of this Feynman integral is quite similar to the previous case [17]. It is based on (2). After a preliminary analysis of the integrand one observes that it is reasonable to start the procedure of resolution of the singularities in , by shifting contours and taking residues, from taking care of the gamma function (−1 − z 2 − z 3 − z 4 ). Here one meets the same problem as in the massless case (see the third reference in [3]) connected with spurious singularities in MB integrals. It can also be cured in the same way, by introducing an auxiliary analytic regularization. To do this, we have chosen a 8 = −1 + x. Then the singularities in the MB integrals are first resolved with respect to x and then with respect to when x and tend to zero. With this complication, the key gamma function becomes (x − 1 − z 2 − z 3 − z 4 ). As a result of the procedure of resolution of the singularities in x and , the spurious singularities in x drop out, and we are left only with six MB integrals which can be evaluated by expanding the integrands in and evaluating the corresponding finite MB integrals.
We have obtained the following result: [12]. The two other contributions to the finite part are more cumbersome: − 2m x p 2 xy + m y p 2 xy + p x p 2 xy − 2l y m x p y − 2l y m xy p y − m x m xy p y + 2m x m y p y + 2l y p x p y − 3m xy p x p y − 2m y p x p y + 2l y p xy p y + 3m x p xy p y + p x p xy p y (7) − m x p 2 y + 3m xy p 2 where y = 1/ 1 − 4m 2 /t and S 2,2 is a generalized polylogarithm [25]. The following abbreviations are also used: (3), l z = ln z for z = x, y, 2, p z = ln(1 + z) and m z = ln(1 − z) for z = x, y, xy. The results have been checked numerically in two independent ways: the finite MB integrals, obtained after the resolution of singularities, have been evaluated by numerical integration along the imaginary axis. The parameter space for the Monte Carlo integration was maximally 4-dimensional in this case. In addition, the overall coefficients of the poles have been calculated at Euclidean points by the completely independent method of sector decomposition [21,23].

Planar double boxes of the second type
For the general planar double box of the second type (see Fig. 1 a sixfold MB representation can be derived [26] similarly to the first case: B PL,2 a 1 , . . . , a 8 ; s, t, m 2 ; We have used this representation to calculate the master planar double box of the second type B PL,2 (1, . . . , 1, 0; s, t, m 2 ; ). The resolution of singularities in was performed similar to the previous cases. The number of resulting MB integrals where an expansion in can be performed in the integrand is again equal to six. We have arrived at the following result: B PL,2 1, . . ., 1, 0; s, t, m 2 ; Furthermore, where the functionb 02 (x, y) is given by the braces in Eq. (9) of [17] which defines the contribution b 02 (x, y) to the finite part of the master double box of the first type. Finally, We have controlled this result similarly to the previous case, by numerical evaluation of our finite MB integrals and numerical evaluation by the method of [23].
It is not clear whether the two-parametric integrals present in the finite part can be expressed in terms of HPL or 2dHPL depending on special combinations of s, t and m 2 . If the answer to this question is negative one might think about the introduction of a new class of functions, e.g., a sort of generalized two-dimensional HPL. In this context, let us point out that generalized HPL of various types were introduced in [14]. These new special functions were defined similarly to HPL, with other basic functions, in particular 1/ √ t (t + 4). However, for example, the GHPL and Li 2 (r, θ ) = Re[Li 2 (re iθ )] is the dilogarithm of a complex argument. Still it is not yet clear whether any GHPL can be expressed in terms of polylogarithms and HPL.
From the mathematical point of view it is natural to try to express any new results in terms of known functions, which may depend on the initial variables through some special combinations. Pragmatically, it is sufficient to present new results in terms of some functions which are, probably, new and which can be evaluated with high precision at physical values of the given variables. Anyway, the two-parametrical integrals presented above can be evaluated numerically in a straightforward way at any non-singular point of the variables s, t, m 2 .

Non-planar double boxes
For the general non-planar double box (see Fig. 1(c)), one can derive the following eightfold MB representation [26]: Although the number of integrations is rather high one can proceed also in this case. However, it turns out that the massive non-planar case is rather complicated. We have performed the resolution of singularities for the non-planar master planar double box. Here is our result for its double-pole part in : where z = 1/ 1 − 4m 2 /u. A one-parametric integral which can hardly be expressed in terms of known special functions is present already in the 1/ part. Unfortunately, we are unable to check our preliminary result (for the coefficient of the simple pole and for the finite part) numerically by the method of [23] for the following reason. As in the case of the massless non-planar double boxes (see the second reference in [2]), it is natural to treat the Mandelstam variables as not restricted by the physical condition s + t + u = 4m 2 , because this condition does not simplify the calculation. Then the natural procedure is to perform the calculation for some extension of the given Feynman integral as a function of the variables s, t, u (and, now, m 2 ) to the Euclidean domain, with s, t, u < 0, and perform the analytical continuation to the physical domain in the result. When one starts from the alpha/Feynman parameters of the whole diagram one naturally arrives at an extension of the physical condition just by considering u as an independent variable. Since, in the massless case, a fourfold MB representation was derived starting from the global Feynman parameter representation, the same extension as in [23] was implied there. In the massive case however, if we start from Feynman parameters and try to separate various terms entering functions present in the integrals over Feynman parameters, we do not see any possibility to arrive at a MB representation with the number of integrations less than ten. The eightfold MB representation (15) was derived in another way, by introducing some Feynman parameters upon integrating over the first loop momentum, then integrating over the second-loop momentum and completing the procedure of introducing MB integrations. After this procedure, the variable u is considered as an independent variable. It turns out that this extension to non-physical Euclidean points defined by (15) differs from the version implied within [23], based on the global Feynman parametric representation for the initial diagram. This can be seen by checking (15) in limiting cases where some of the indices a i are zero and where one has explicit analytical results: In some of these limits, it is necessary to use the physical condition s + t + u = 4m 2 to prove agreement with the known result, while the agreement is achieved without this condition when one starts from the global Feynman parametric representation.
Thus, to perform numerical control of our results, a generalization of the method of sector decompositions [23] to points with kinematical invariants of different signs is necessary and will be developed in the near future. Then one will be able to make reliable results for massive non-planar double boxes.