# 1 Introduction

## 1 Introduction

The massless double box diagram shown in Fig. 1 enters many important physical observables, e.g., amplitudes of the Bhabba scattering at high energies. An experience shows that master diagrams, i.e. with all powers of propagators equal to one, are most complicated for evaluation. In the massless off-shell case, the master double box Feynman integral has been analytically evaluated in  strictly in four dimensions. It is the purpose of the present paper to evaluate it analytically on shell, i.e. for $`p_i^2=0,i=1,2,3,4`$, in the framework of dimensional regularization , with the space-time dimension $`d=42ϵ`$ as a regularization parameter.

To do this, we start, in the next Section, from the alpha-representation of the double box and, after expanding some of the involved functions in Mellin–Barnes (MB) integrals, arrive at a five-fold MB integral representation with gamma functions in the integrand. We use, in Sec. 3, a standard procedure of taking residues and shifting contours to resolve the structure of singularities in the parameter of dimensional regularization, $`ϵ`$. This procedure leads to an appearance of multiple terms where Laurent expansion in $`ϵ`$ becomes possible. The resulting integrals in all the MB parameters but one are evaluated explicitly in gamma functions and their derivatives. In Sec. 4, the last MB integral is evaluated by closing an initial integration contour in the complex plane to the right, with an explicit summation of the corresponding series. A final result is expressed in terms of polylogarithms $`\text{Li}_a\left(t/s\right)`$, up to $`a=4`$, and generalized polylogarithms $`S_{a,b}(t/s)`$, with $`a=1,2`$ and $`b=2`$. Starting from the same one-fold MB integral and closing the contour of integration to the left, we obtain a similar result written through the same class of functions depending on the inverse ratio, $`s/t`$. Furthermore, we obtain, as a by-product, an explicit result for the backward scattering value, i.e. at $`t=s`$, of the double box diagram.

## 2 From momentum space to MB representation

The massless on-shell double box Feynman integral can be written as

$`{\displaystyle \frac{\text{d}^dk\text{d}^dl}{(k^2+2p_1k)(k^22p_2k)k^2(kl)^2(l^2+2p_1l)(l^22p_2l)(lp_1p_3)^2}}`$   
$`{\displaystyle \frac{\left(i\pi ^{d/2}\mathrm{e}^{\gamma _\mathrm{E}ϵ}\right)^2}{(s)^{2+2ϵ}(t)}}K(t/s,ϵ),`$ (1)   
where $`s=(p_1+p_2)^2,t=(p_1+p_3)^2`$, and $`k`$ and $`l`$ are respectively loop momenta of the left and the right box. Usual prescriptions, $`k^2=k^2+i0,s=si0`$, etc are implied. We have pulled out not only standard factors that arise when integrating in the loop momenta but also a factor that makes the resulting function $`K`$ depend on the dimensionless variable, $`x=t/s`$.

The alpha representation of the double box is straightforwardly obtained:

$$K(x,ϵ)=\mathrm{\Gamma }(3+2ϵ)_0^{\mathrm{}}\text{d}\alpha _1\mathrm{}_0^{\mathrm{}}\text{d}\alpha _7\delta \left(\alpha _i1\right)D^{1+3ϵ}\left(A+x\alpha _5\alpha _6\alpha _7\right)^{32ϵ},$$
(2)   
where

$`D`$ $`=`$ $`(\alpha _1+\alpha _2+\alpha _7)(\alpha _3+\alpha _4+\alpha _5)+\alpha _6(\alpha _1+\alpha _2+\alpha _3+\alpha _4+\alpha _5+\alpha _7),`$ (3)   
$`A`$ $`=`$ $`\alpha _1\alpha _2(\alpha _3+\alpha _4+\alpha _5)+\alpha _3\alpha _4(\alpha _1+\alpha _2+\alpha _7)+\alpha _6(\alpha _1+\alpha _3)(\alpha _2+\alpha _4).`$ (4)   
As it is well-known, one can choose a sum of an arbitrary subset of $`\alpha _i,i=1,\mathrm{},7`$ in the argument of the delta function in (2). We choose it as $`\delta \left(_{i6}\alpha _i1\right)`$ and change variables by turning from alpha to Feynman parameters

$`\alpha _3=\alpha _{35}\xi _1,\alpha _5=\alpha _{35}(1\xi _1),\alpha _1=\alpha _{17}\xi _3,\alpha _7=\alpha _{17}(1\xi _3),`$   
$`\alpha _{35}=\xi _5\xi _2,\alpha _4=\xi _5(1\xi _2),\alpha _{17}=(1\xi _5)\xi _4,\alpha _2=(1\xi _5)(1\xi _4).`$ (5)   
to obtain the following parametric integral:

$`K(x,ϵ)=\mathrm{\Gamma }(3+2ϵ){\displaystyle _0^{\mathrm{}}}\text{d}\alpha _6{\displaystyle _0^1}\text{d}\xi _1\mathrm{}{\displaystyle _0^1}\text{d}\xi _5\xi _2\xi _4\xi ^2(1\xi )^2`$ (6)   
$`\times (\alpha _6+\xi _5(1\xi _5))^{1+3ϵ}Q^{32ϵ},`$   
where

$`Q`$ $`=`$ $`x\alpha _6(1\xi _1)\xi _2(1\xi _3)\xi _4(1\xi _5)\xi _5`$ (7)   
$`+\xi _5(1\xi _5)[\xi _5\xi _1\xi _2(1\xi _2)+(1\xi _5)\xi _3\xi _4(1\xi _4)]`$   
$`+\alpha _6[\xi _5\xi _1\xi _2+(1\xi _5)\xi _3\xi _4][\xi _5(1\xi _2)+(1\xi _5)(1\xi _4)].`$   
We are now going to apply five times the MB representation

$$\frac{1}{(X+Y)^\nu }=\frac{1}{\mathrm{\Gamma }(\nu )}\frac{1}{2\pi i}_i\mathrm{}^{+i\mathrm{}}\text{d}w\frac{Y^w}{X^{\nu +w}}\mathrm{\Gamma }(\nu +w)\mathrm{\Gamma }(w),$$
(8)   
where the contour of integration is chosen in the standard way: the poles with the $`\mathrm{\Gamma }(\mathrm{}+w)`$-dependence (let us call them infrared (IR) poles) are to the left of the contour and the poles with the $`\mathrm{\Gamma }(\mathrm{}w)`$-dependence (ultraviolet (UV) poles) are to the right of it.

First, we introduce a MB integration, in $`w`$, using decomposition of the function $`Q`$ with $`Y`$ as the first line in (7). We introduce a second MB integral choosing as $`X`$ the term with $`\alpha _6`$ in the rest part of $`Q`$. After that we can take an integral in $`\alpha _6`$ in gamma functions. The next three MB integrations, in $`z_1,z_2`$ and $`z_3`$, are to separate terms in the following three combinations: $`[\xi _5\xi _1\xi _2+(1\xi _5)\xi _3\xi _4]`$, $`[\xi _5(1\xi _2)+(1\xi _5)(1\xi _4)]`$ and $`[\xi _5\xi _1\xi _2(1\xi _2)+(1\xi _5)\xi _3\xi _4(1\xi _4)]`$.

All the integrals in Feynman parameters are then taken explicitly in gamma functions. Finally, we perform the change of variables $`z_2=w_2z_11,z_3=w_3z_11`$ and arrive at the following nice 5-fold MB integral:

$`K(x,ϵ)={\displaystyle \frac{1}{\mathrm{\Gamma }(13ϵ)}}{\displaystyle \frac{1}{(2\pi i)^5}}{\displaystyle \text{d}w\text{d}w_2\text{d}w_3\text{d}z\text{d}z_1x^{w+1}}`$   
$`\times \mathrm{\Gamma }(1+w)^2\mathrm{\Gamma }(w)\mathrm{\Gamma }(w_2)\mathrm{\Gamma }(12ϵww_2)\mathrm{\Gamma }(w_3)\mathrm{\Gamma }(12ϵww_3)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(1w_2+z_1)\mathrm{\Gamma }(1w_3+z_1)\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)\mathrm{\Gamma }(z_1)}{\mathrm{\Gamma }(1+w+w_2+w_3)\mathrm{\Gamma }(14ϵww_2w_3)}}`$   
$`\times \mathrm{\Gamma }(1ϵ+z)\mathrm{\Gamma }(2+2ϵ+w+w_2+zz_1)\mathrm{\Gamma }(2+2ϵ+w+w_3+zz_1)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(23ϵww_2w_3+z_1z)\mathrm{\Gamma }(z_1z)}{\mathrm{\Gamma }(3+2ϵ+w+z)}}.`$ (9)   
One can interchange the order of integration in an arbitrary way. For each order, the rules of dealing with poles are as formulated above. Note that if we have a product $`\mathrm{\Gamma }(a+v)\mathrm{\Gamma }(bv)`$, for some integration variable $`v=w,w_2,w_3,z,z_1`$ with $`a`$ and $`b`$ dependent on other variables, then the integration in $`v`$ produces a singularity of the type $`\mathrm{\Gamma }(a+b)`$.

## 3 Resolving singularities in $`ϵ`$

Since it looks hopeless to evaluate our MB integral for general $`ϵ`$ let us try to obtain a result in expansion in $`ϵ`$ up to the finite part. There is a factor $`1/\mathrm{\Gamma }(13ϵ)`$ proportional to $`ϵ`$ when $`ϵ`$ tends to zero. Representation (9) is therefore effectively 4-fold because to generate a contribution that does not vanish at $`ϵ=0`$ we need to take a residue at least in one of the integration variables. None of the integrations can however immediately produce an explicit $`ϵ`$-pole. Let us first distinguish the following two gamma functions

$$\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)\mathrm{\Gamma }(23ϵww_2w_3+z_1z)$$

that are essential for the appearance of the poles.

We can write down the integral in $`z_1`$ as minus residue at the point $`z_1=ϵ+w+w_2+w_3`$ (where the gamma function $`\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)`$ has its first pole which is UV, with respect to $`z_1`$) plus an integral with the same integrand where this pole is IR. We can similarly write down the integral in $`z`$ as minus residue at the point $`z=23ϵww_2w_3+z_1`$ (where the gamma function $`\mathrm{\Gamma }(23ϵww_2w_3+z_1z)`$ has its first pole which is UV, with respect to $`z`$) plus an integral with the same integrand where this pole is IR.

As a result we decompose integral (9) as $`K=K_{00}+K_{01}+K_{10}+K_{11}`$ where $`K_{11}`$ corresponds to the two residues, $`K_{10}`$ to the residue in $`z`$ and the integral in $`z_1`$ with the opposite nature of the first pole of $`\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)`$, etc. For example, the contribution $`K_{11}`$ is given by the following 3-fold integral:

$`K_{11}(x,ϵ)={\displaystyle \frac{1}{(2\pi i)^3}}{\displaystyle \text{d}w\text{d}w_2\text{d}w_3x^{w+1}\mathrm{\Gamma }(1+w)\mathrm{\Gamma }(w)}`$   
$`\times \mathrm{\Gamma }(w_2)\mathrm{\Gamma }(ϵw_2)\mathrm{\Gamma }(1+ϵ+w+w_2)\mathrm{\Gamma }(12ϵww_2)`$   
$`\times \mathrm{\Gamma }(w_3)\mathrm{\Gamma }(ϵw_3)\mathrm{\Gamma }(1+ϵ+w+w_3)\mathrm{\Gamma }(12ϵww_3)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(2+3ϵ+w+w_2+w_3)\mathrm{\Gamma }(ϵww_2w_3)}{\mathrm{\Gamma }(1+w+w_2+w_3)\mathrm{\Gamma }(14ϵww_2w_3)}}.`$ (10)   
This contribution is in turn decomposed, in a similar way, as $`K_{11}=_{i,j=0,1,2}K_{11ij}`$. Here the value $`i=1`$ of the first index denotes the residue in $`w_2`$ at the point $`w_2=0`$. The value $`i=2`$ denotes the residue in $`w_2`$ at $`w_2=1ϵw`$ of the integrand where the first pole of $`\mathrm{\Gamma }(w_2)`$ is UV rather than IR. Finally, $`i=0`$ means that both above poles are IR. The second index similarly refers to the integral in $`w_3`$.

In particular, we have

$`K_{1111}(x,ϵ)=\mathrm{\Gamma }(ϵ)^2{\displaystyle \frac{1}{2\pi i}}{\displaystyle \text{d}wx^{w+1}}`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(1+ϵ+w)^2\mathrm{\Gamma }(2+3ϵ+w)\mathrm{\Gamma }(12ϵw)^2\mathrm{\Gamma }(ϵw)\mathrm{\Gamma }(w)}{\mathrm{\Gamma }(14ϵw)}},`$ (11)   
$`K_{1112}(x,ϵ)=K_{1121}(x,ϵ)={\displaystyle \frac{\mathrm{\Gamma }(1+2ϵ)\mathrm{\Gamma }(ϵ)}{\mathrm{\Gamma }(3ϵ)}}{\displaystyle \frac{1}{2\pi i}}{\displaystyle \text{d}wx^{w+1}}`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(1+w)^2\mathrm{\Gamma }(1+ϵ+w)^2\mathrm{\Gamma }(12ϵw)\mathrm{\Gamma }(ϵw)\mathrm{\Gamma }(w)}{\mathrm{\Gamma }(2+ϵ+w)}},`$ (12)   
$`K_{1122}(x,ϵ)=\mathrm{\Gamma }(ϵ)^2{\displaystyle \frac{1}{2\pi i}}{\displaystyle \text{d}wx^{w+1}}`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(1+w)^3\mathrm{\Gamma }(1+ϵ+w)^2\mathrm{\Gamma }(ϵw)^2\mathrm{\Gamma }(ϵw)\mathrm{\Gamma }(w)}{\mathrm{\Gamma }(2+ϵ+w)\mathrm{\Gamma }(12ϵ+w)\mathrm{\Gamma }(12ϵw)}}.`$ (13)   
The next contribution is

$`K_{10}(x,ϵ)={\displaystyle \frac{1}{\mathrm{\Gamma }(13ϵ)}}{\displaystyle \frac{1}{(2\pi i)^4}}{\displaystyle \text{d}w\text{d}w_2\text{d}w_3\text{d}z_1x^{w+1}\mathrm{\Gamma }(1+w)^2\mathrm{\Gamma }(w)}`$   
$`\times \mathrm{\Gamma }(w_2)\mathrm{\Gamma }(ϵw_2)\mathrm{\Gamma }(12ϵww_2)\mathrm{\Gamma }(w_3)\mathrm{\Gamma }(ϵw_3)\mathrm{\Gamma }(12ϵww_3)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(2+3ϵ+w+w_2+w_3)\mathrm{\Gamma }(1w_2+z_1)\mathrm{\Gamma }(1w_3+z_1)}{\mathrm{\Gamma }(1+w+w_2+w_3)\mathrm{\Gamma }(14ϵww_2w_3)}}`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(14ϵww_2w_3+z_1)\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)\mathrm{\Gamma }(z_1)}{\mathrm{\Gamma }(1ϵw_2w_3+z_1)}},`$ (14)   
where the first pole of $`\mathrm{\Gamma }(ϵ+w+w_2+w_3z_1)`$ is IR, with respect to $`z_1`$, rather than UV. We further decompose this contribution by changing the nature of the first pole of $`\mathrm{\Gamma }(14ϵww_2w_3+z_1)`$ in $`z_1`$. We obtain $`K_{10}=K_{100}+K_{101}`$, where the new index 1 corresponds to the residue and has the form

$`K_{101}(x,ϵ)={\displaystyle \frac{1}{(2\pi i)^3}}{\displaystyle \text{d}w\text{d}w_2\text{d}w_3x^{w+1}\frac{\mathrm{\Gamma }(1+w)^2\mathrm{\Gamma }(w)}{\mathrm{\Gamma }(2+3ϵ+w)}}`$   
$`\times \mathrm{\Gamma }(w_2)\mathrm{\Gamma }(ϵw_2)\mathrm{\Gamma }(2+4ϵ+w+w_2)\mathrm{\Gamma }(12ϵww_2)\mathrm{\Gamma }(w_3)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(ϵw_3)\mathrm{\Gamma }(2+4ϵ+w+w_3)\mathrm{\Gamma }(12ϵww_3)\mathrm{\Gamma }(2+3ϵ+w+w_2+w_3)}{\mathrm{\Gamma }(1+w+w_2+w_3)}}.`$ (15)   
Each of the contributions $`K_{100}`$ and $`K_{101}`$ is then decomposed using the change of the nature of poles $`w_2=0`$ and $`w_3=0`$. We obtain $`K_{10j}=K_{10j00}+K_{10j01}+K_{10j10}+K_{10j11}`$, for $`j=0,1`$. Here the value $`i=1`$ of the last index denotes the residue in $`w_3`$ at the point $`w_3=0`$ and the $`i=0`$ an integral where the first pole of $`\mathrm{\Gamma }(w_3)`$ is considered UV. The second index from the end similarly refers to $`\mathrm{\Gamma }(w_2)`$. For example,

$$K_{10111}(x,ϵ)=\mathrm{\Gamma }(ϵ)^2\frac{1}{2\pi i}\text{d}wx^{w+1}\mathrm{\Gamma }(2+4ϵ+w)^2\mathrm{\Gamma }(1+w)\mathrm{\Gamma }(12ϵw)^2\mathrm{\Gamma }(w).$$
(16)   
Then we have

$`K_{01}(x,ϵ)={\displaystyle \frac{1}{\mathrm{\Gamma }(13ϵ)}}{\displaystyle \frac{1}{(2\pi i)^4}}{\displaystyle \text{d}w\text{d}w_2\text{d}w_3\text{d}zx^{w+1}\mathrm{\Gamma }(1+w)^2\mathrm{\Gamma }(w)}`$   
$`\times \mathrm{\Gamma }(w_2)\mathrm{\Gamma }(1+ϵ+w+w_2)\mathrm{\Gamma }(12ϵww_2)\mathrm{\Gamma }(w_3)\mathrm{\Gamma }(1+ϵ+w+w_3)`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(12ϵww_3)\mathrm{\Gamma }(ϵww_2w_3)\mathrm{\Gamma }(1ϵ+z)}{\mathrm{\Gamma }(1+w+w_2+w_3)\mathrm{\Gamma }(14ϵww_2w_3)}}`$   
$`\times {\displaystyle \frac{\mathrm{\Gamma }(2+ϵw_2+z)\mathrm{\Gamma }(2+ϵw_3+z)\mathrm{\Gamma }(ϵ+w+w_2+w_3z)\mathrm{\Gamma }(22ϵz)}{\mathrm{\Gamma }(3+2ϵ+w+z)}},`$ (17)   
where the first pole of $`\mathrm{\Gamma }(22ϵz)`$ is IR, with respect to $`z`$, rather than UV. Using the change of variables $`w_212ϵww_2,w_312ϵww_3,z_1ϵww_2w_3+z`$ in $`K_{10}`$ we see that $`K_{01}K_{10}`$. Finally, the contribution $`K_{00}`$ is similarly decomposed: $`K_{00}=K_{0000}+K_{0001}+K_{0010}+K_{0011}`$.

Now we observe that, in each of the obtained contributions, the only additional (with respect to explicit gamma functions depending on $`ϵ`$) source of the poles in $`ϵ`$ is the last integration, in $`w`$, where the first (UV) pole of the gamma function $`\mathrm{\Gamma }(12ϵw)`$ glues with an IR pole of $`\mathrm{\Gamma }(1+w)`$ or $`\mathrm{\Gamma }(1+ϵ+w)`$ when $`ϵ0`$ — see such examples in (1113) and (16). Therefore we further decompose each of the contributions into two pieces: minus residue at the point $`w=12ϵ`$ plus an integral where we can integrate in the region $`1<\mathrm{Re}w<0`$. In each of these pieces, we now can expand an integrand in a Laurent series in $`ϵ`$ up to the finite part. In particular, no poles in $`ϵ`$ arise in $`K_{0000}`$ so that it it zero at $`ϵ=0`$ because of the overall factor $`1/\mathrm{\Gamma }(13ϵ)`$.

We collect separately the pieces from these last residues and from the last integration at $`1<\mathrm{Re}w<0`$. The first collection gives the leading order term in the expansion of the double box in the limit $`t/s0`$ while the second collection involves the rest of the terms of this expansion. A remarkable fact is that, in all these multiple contributions, the integrations in $`w_2,w_3,z,z_1`$ can be performed analytically, with the help of the first and the second Barnes lemmas

$`{\displaystyle \frac{1}{2\pi i}}{\displaystyle _i\mathrm{}^{+i\mathrm{}}}\text{d}w\mathrm{\Gamma }(\lambda _1+w)\mathrm{\Gamma }(\lambda _2+w)\mathrm{\Gamma }(\lambda _3w)\mathrm{\Gamma }(\lambda _4w)`$ (18)   
$`={\displaystyle \frac{\mathrm{\Gamma }(\lambda _1+\lambda _3)\mathrm{\Gamma }(\lambda _1+\lambda _4)\mathrm{\Gamma }(\lambda _2+\lambda _3)\mathrm{\Gamma }(\lambda _2+\lambda _4)}{\mathrm{\Gamma }(\lambda _1+\lambda _2+\lambda _3+\lambda _4)}},`$   
$`{\displaystyle \frac{1}{2\pi i}}{\displaystyle _i\mathrm{}^{+i\mathrm{}}}\text{d}w{\displaystyle \frac{\mathrm{\Gamma }(\lambda _1+w)\mathrm{\Gamma }(\lambda _2+w)\mathrm{\Gamma }(\lambda _3+w)\mathrm{\Gamma }(\lambda _4w)\mathrm{\Gamma }(\lambda _5w)}{\mathrm{\Gamma }(\lambda _1+\lambda _2+\lambda _3+\lambda _4+\lambda _5+w)}}`$ (19)   
$`={\displaystyle \frac{\mathrm{\Gamma }(\lambda _1+\lambda _4)\mathrm{\Gamma }(\lambda _2+\lambda _4)\mathrm{\Gamma }(\lambda _3+\lambda _4)\mathrm{\Gamma }(\lambda _1+\lambda _5)\mathrm{\Gamma }(\lambda _2+\lambda _5)\mathrm{\Gamma }(\lambda _3+\lambda _5)}{\mathrm{\Gamma }(\lambda _1+\lambda _2+\lambda _4+\lambda _5)\mathrm{\Gamma }(\lambda _1+\lambda _3+\lambda _4+\lambda _5)\mathrm{\Gamma }(\lambda _2+\lambda _3+\lambda _4+\lambda _5)}}`$   
and their corollaries. These are two typical examples of such corollaries:

$`{\displaystyle \frac{1}{2\pi i}}{\displaystyle _i\mathrm{}^{+i\mathrm{}}}\text{d}w{\displaystyle \frac{\mathrm{\Gamma }(\lambda _1+w)\mathrm{\Gamma }(\lambda _2+w)^2\mathrm{\Gamma }(\lambda _2w)\mathrm{\Gamma }(\lambda _3w)}{\mathrm{\Gamma }(\lambda _1+\lambda _2+\lambda _3+w)}}`$ (20)   
$`={\displaystyle \frac{\mathrm{\Gamma }(\lambda _1\lambda _2)\mathrm{\Gamma }(\lambda _2+\lambda _3)\left[\psi ^{}\left(\lambda _1+\lambda _3\right)\psi ^{}\left(\lambda _2+\lambda _3\right)\right]}{\mathrm{\Gamma }(\lambda _1+\lambda _3)}},`$   
where the pole $`w=\lambda _2`$ is considered IR while other poles are treated in the standard way, and

$$\frac{1}{2\pi i}_{1/2i\mathrm{}}^{1/2+i\mathrm{}}\text{d}w\mathrm{\Gamma }(1+w)\mathrm{\Gamma }(w)\mathrm{\Gamma }(w)\mathrm{\Gamma }(1w)\psi (1+w)^2=\frac{\gamma _\mathrm{E}^2\pi ^2}{3}+6\gamma _\mathrm{E}\zeta (3)+\frac{\pi ^4}{45}.$$
(21)   
Here $`\gamma _\mathrm{E}`$ is the Euler constant, $`\psi (z)`$ the logarithmical derivative of the gamma function, and $`\zeta (z)`$ the Riemann zeta function.

After taking these integrations and summing up the resulting contributions into the two above collections we obtain the following result

$`K(x,ϵ)`$ $`=`$ $`K_{0t}(x,ϵ)+K_{1t}(x,ϵ)+o(ϵ),`$ (22)   
$`K_{0t}(x,ϵ)`$ $`=`$ $`{\displaystyle \frac{4}{ϵ^4}}+{\displaystyle \frac{5\mathrm{ln}x}{ϵ^3}}\left(2\mathrm{ln}^2x{\displaystyle \frac{5}{2}}\pi ^2\right){\displaystyle \frac{1}{ϵ^2}}`$ (23)   
$`\left({\displaystyle \frac{2}{3}}\mathrm{ln}^3x+{\displaystyle \frac{11}{2}}\pi ^2\mathrm{ln}x{\displaystyle \frac{65}{3}}\zeta (3)\right){\displaystyle \frac{1}{ϵ}}+{\displaystyle \frac{4}{3}}\mathrm{ln}^4x+6\pi ^2\mathrm{ln}^2x{\displaystyle \frac{88}{3}}\zeta (3)\mathrm{ln}x+{\displaystyle \frac{29}{30}}\pi ^4,`$   
$`K_{1t}(x,ϵ)`$ $`=`$ $`{\displaystyle \frac{2}{\pi i}}{\displaystyle \frac{\text{d}wx^{w+1}}{1+w}\mathrm{\Gamma }(1+w)^3\mathrm{\Gamma }(w)^3}`$ (24)   
$`\times \left[{\displaystyle \frac{1}{ϵ}}{\displaystyle \frac{5}{1+w}}+3\psi (1+w)4\psi (w)\gamma _\mathrm{E}\right].`$   
Let us stop for a moment and observe that this result provides, in a very easy way, not only numerical evaluation of the double box diagram for general values of $`s`$ and $`t`$ but also asymptotic expansions in the limits $`t/s0`$ and $`s/t0`$ which are obtained by taking series of residues respectively to the right or to the left.

## 4 Evaluating the last MB integral

The last MB integration, in (24), is performed analytically by taking the sum of the residues at the points $`w=0,1,2,\mathrm{}`$ and summing up the resulting series. In this last step, we use, in particular, summation formulae derived in . Here is the final result:

$`K_{1t}(x,ϵ)`$ $`=`$ $`\left[2\text{Li}_3\left(x\right)2\mathrm{ln}x\text{Li}_2\left(x\right)\left(\mathrm{ln}^2x+\pi ^2\right)\mathrm{ln}(1+x)\right]{\displaystyle \frac{2}{ϵ}}`$ (25)   
$`4\left(S_{2,2}(x)\mathrm{ln}xS_{1,2}(x)\right)+44\text{Li}_4\left(x\right)4\left(\mathrm{ln}(1+x)+6\mathrm{ln}x\right)\text{Li}_3\left(x\right)`$   
$`+2\left(\mathrm{ln}^2x+2\mathrm{ln}x\mathrm{ln}(1+x)+{\displaystyle \frac{10}{3}}\pi ^2\right)\text{Li}_2\left(x\right)`$   
$`+\left(\mathrm{ln}^2x+\pi ^2\right)\mathrm{ln}^2(1+x){\displaystyle \frac{2}{3}}\left(4\mathrm{ln}^3x+5\pi ^2\mathrm{ln}x6\zeta (3)\right)\mathrm{ln}(1+x),`$   
where $`\text{Li}_a\left(z\right)`$ is the polylogarithm  and

$$S_{a,b}(z)=\frac{(1)^{a+b1}}{(a1)!b!}_0^1\frac{\mathrm{ln}^{a1}(t)\mathrm{ln}^b(1zt)}{t}\text{d}t,$$
(26)   
a generalized polylogarithm introduced in . Note that any (generalized) polylogarithms involved can be expanded in a Taylor series at $`x=0`$ with the radius of convergence equal to one.

We can similarly close the integration contour to the left and obtain a result in a form of functions depending on the inverse ratio, $`y=1/x`$:

$`K(x,ϵ)`$ $`=`$ $`K_{0s}(x,ϵ)+K_{1s}(x,ϵ)+o(ϵ),`$ (27)   
$`K_{0s}(1/y,ϵ)`$ $`=`$ $`{\displaystyle \frac{4}{ϵ^4}}{\displaystyle \frac{5\mathrm{ln}y}{ϵ^3}}\left(2\mathrm{ln}^2y{\displaystyle \frac{5}{2}}\pi ^2\right){\displaystyle \frac{1}{ϵ^2}}`$ (28)   
$`+\left({\displaystyle \frac{7}{2}}\pi ^2\mathrm{ln}y+{\displaystyle \frac{65}{3}}\zeta (3)\right){\displaystyle \frac{1}{ϵ}}+{\displaystyle \frac{1}{3}}\pi ^2\mathrm{ln}^2y+{\displaystyle \frac{76}{3}}\zeta (3)\mathrm{ln}y{\displaystyle \frac{83}{90}}\pi ^4,`$   
$`K_{1s}(1/y,ϵ)`$ $`=`$ $`\left[2\text{Li}_3\left(y\right)2\mathrm{ln}y\text{Li}_2\left(y\right)\left(\mathrm{ln}^2y+\pi ^2\right)\mathrm{ln}(1+y)\right]{\displaystyle \frac{2}{ϵ}}`$ (29)   
$`4\left(S_{2,2}(y)\mathrm{ln}yS_{1,2}(y)\right)36\text{Li}_4\left(y\right)4\left(\mathrm{ln}(1+y)5\mathrm{ln}y\right)\text{Li}_3\left(y\right)`$   
$`2\left(\mathrm{ln}^2y2\mathrm{ln}y\mathrm{ln}(1+y)+{\displaystyle \frac{10}{3}}\pi ^2\right)\text{Li}_2\left(y\right)`$   
$`+\left(\mathrm{ln}^2y+\pi ^2\right)\mathrm{ln}^2(1+y)+2\left(\mathrm{ln}^3y+{\displaystyle \frac{2}{3}}\pi ^2\mathrm{ln}y+2\zeta (3)\right)\mathrm{ln}(1+y).`$   
As a by-product, we obtain an explicit result for the backward scattering value of (1), i.e. at $`t=s`$,

$$\frac{\left(i\pi ^{d/2}\right)^2\mathrm{e}^{2\gamma _\mathrm{E}ϵ}}{(s)^{3+2ϵ}}\left[\frac{4}{ϵ^4}\frac{9\pi ^2}{2ϵ^2}\frac{53\zeta (3)}{3ϵ}+\frac{22\pi ^4}{9}\pi i\left(\frac{5}{ϵ^3}\frac{25\pi ^2}{6ϵ}\frac{148\zeta (3)}{3}\right)\right].$$
(30)   
The presented algorithm is applicable to massless on-shell box Feynman integrals with any integer powers of propagators.

Acknowledgments. I am grateful to M. Beneke and the CERN Theory Group for kind hospitality during my visit to CERN in April–May 1999 where this work was completed. Thanks to M.B. for careful reading the manuscript and to O.L. Veretin for involving me into this problem and useful discussions. This work was supported by Volkswagen Foundation, contract No. I/73611, and by the Russian Foundation for Basic Research, project 98–02–16981.