††institutetext: Institut für Theoretische Physik, Universität Tübingen,
Kepler Center for Astro and Particle Physics,
Auf der Morgenstelle 14, D-72076 Tübingen, Germany

# The massless single off-shell scalar box integral — branch cut structure and
all-order epsilon expansion

Juliane Haug 1 and Fabian Wunder111Corresponding author juliane-clara-
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We investigate the single off-shell scalar box integral with massless internal
lines in dimensional regularization. A special emphasis is given to higher
orders in the dimensional regularization parameter epsilon, its branch cut
structure, and kinematic limits. Common representations of the box integral
introduce superficial branch cuts, which we eliminate to all orders in the
epsilon expansion. In the literature so far a satisfactory solution for this
issue only exists up to finite order in the epsilon expansion. However, for
certain calculations at NNLO in perturbation theory, higher orders in epsilon
of this integral are required. In this paper, we present results to all orders
in epsilon in terms of single-valued polylogarithms and explicitly determine
the real and imaginary part of the box integral in all kinematic regions.

###### Keywords:

Feynman integrals, perturbative QCD, dimensional regularization, epsilon
expansion

††arxiv: 2211.14110

##  1 Introduction

With the upcoming electron-ion collider (EIC) (Accardi:2012qut, ), it is
essential to evaluate the hard-scattering coefficients for semi-inclusive
deep-inelastic scattering (SIDIS) at next-to-next-to-leading order (NNLO) in
QCD (Anderle:2016kwa, ; Abele:2021nyo, ; Borsa:2022vvp, ; Abele:2022wuy, ).
One type of contributions arising at NNLO are the one-loop corrections to the
next-to-leading-order processes of real gluon emission and photon-gluon
fusion. One-loop integrals such as these arise in virtually all higher-order
perturbative calculations and were thoroughly studied in the past 40 years
(tHooft1979, ; vanNeerven1984, ; vanOldenborgh1990, ; Beenakker1990, ;
Denner1991, ; vanOldenborgh1991, ; Hahn1999, ; Fleischer2003, ; Patel2015, ;
Carrazza2016, ; FabiValery, ). Given that tensorial loop integrals can be
systematically reduced to scalar integrals using the Passarino-Veltman
technique (PassarinoVeltman, ), it is sufficient to evaluate scalar integrals.
The case of massless propagators is of particular interest for perturbative
QCD calculations, since light quark masses can be neglected in high-energy
scattering processes and gluons are massless. Loop integrals with vanishing
internal masses generally have soft and collinear divergences, which must be
regularized. The most commonly used method is dimensional regularization,
which simultaneously serves to regularize ultraviolet and infrared divergences
(Bollini:1972ui, ; tHooft:1973mfk, ; Gastmans:1973uv, ; Leibbrandt:1975dj, ).
We will work in $d=4-2\varepsilon$ dimensions.

Both real and imaginary part of a loop integral do in general yield
contributions to the absolute value squared of a matrix element, hence we must
know both. The imaginary part of a loop integral arises in certain kinematic
regions and will be governed by the causal $+\mathrm{i}0$ from the
propagators. For massless two- and three-propagator scalar loop integrals in
general dimensions $d$, it is easy to keep track of this infinitesimal
imaginary part.222By dimensional analysis, massless one-propagator tadpole
integrals must vanish since they depend on no mass scale. Following the
calculation given in ref. (FabiValery, ), one finds that the external momentum
invariants are essentially replaced by $(p_{i}\pm
p_{j})^{2}\rightarrow(p_{i}\pm p_{j})^{2}+\mathrm{i}0$ in the final results.
This simple rule yields no unambiguously well-defined imaginary part for the
four-propagator box integral in general dimension $d$, as we will discuss in
section 2 below.

The scalar box integral with one off-shell external particle (i.e. with non-
vanishing momentum squared) is especially important since it is needed for
higher order (NNLO) perturbative calculations of processes such as SIDIS,
Drell-Yan, or single inclusive annihilation (SIA), which all feature a single
off-shell gauge boson. In ref. (FabriciusSchmitt, ), this integral and its
imaginary part were first given up to finite order $\varepsilon^{0}$ in
dimensional regularization. The result for general dimensions $d$ written in
terms of three Gauss hypergeometric functions was obtained in refs.
(Matsuura1989, ; BernDixonKosower1, ). However, the simple prescription
$(p_{i}\pm p_{j})^{2}\rightarrow(p_{i}\pm p_{j})^{2}+\mathrm{i}0$ given in
ref. (BernDixonKosower1, ) for analytic continuation of the result to all
kinematic regions only applies to the expansion up to finite order. Ref.
(Matsuura1989, , eq. (D.4)), which was used in the calculation of DIS
structure functions at NNLO (Zijlstra1992, ), does not touch on the question
of analytic continuation. In (FabiValery, ), the box integral was expanded to
all orders in $\varepsilon$, however the branch cut structure was also not
discussed. The causal $+\mathrm{i}0$ from the propagators was systematically
kept in ref. (DuplancicNizic, ), but the result was only given up to finite
order in dimensional regularization. To regularize divergences occurring in
kinematic limits, however, we need the full $\varepsilon$-dependence of any
terms divergent in these kinematic limits.

In this paper, we independently recalculate the massless single off-shell
scalar box integral for general dimensions $d$. We explicitly keep the causal
$+\mathrm{i}0$ from the propagators, such that our results are valid for
arbitrary (real) values of the kinematic variables. This allows us to
explicitly determine real and imaginary part to all orders of the
$\varepsilon$-expansion, which goes beyond the finite
$\mathcal{O}(\varepsilon^{0})$ results of refs. (BernDixonKosower1, ;
DuplancicNizic, ) and the Mathematica implementation in Package-X (Patel2015,
). In section 2, we calculate the internal-massless single-off shell scalar
box integral using the well-known Feynman parametrization approach and an
elegant factorization which allows for a trivial reduction to one-dimensional
integrals. As a limiting case our result also includes the internal-massless
box integral with on-shell massless external particles, which is needed for
direct photon production and Compton-like scattering processes. In section 3,
we perform the $\varepsilon$-expansion to all orders and the explicit
determination of real and imaginary parts. The special cases of vanishing
Mandelstam variables, as well as several identities for the Gauss
hypergeometric function, and the definition of single-valued polylogarithms,
are collected in the appendices. Additionally, we provide an example
illustrating the need of higher orders in $\varepsilon$ of the scalar box
integral for phase space regularization in appendix D.

##  2 Calculating the massless single off-shell scalar box integral

Figure 1: The massless scalar box diagram with one off-shell external particle
carrying momentum $q$. It is $p_{1}^{2}=p_{2}^{2}=p_{3}^{2}=0$ and $q^{2}\neq
0$. All external momenta are taken to be incoming. Drawn with TikZ-Feynman
(Ellis2017, ).

The massless single off-shell scalar box integral is defined through the loop
integral

$\displaystyle\mathrm{D}_{0}\,\equiv\,\frac{\mu^{4-d}}{\mathrm{i}\pi^{d/2}}\int\mathrm{d}^{d}l\,\frac{1}{\left[l^{2}+\mathrm{i}0\right]\left[(l+p_{2})^{2}+\mathrm{i}0\right]\left[(l+p_{2}+p_{3})^{2}+\mathrm{i}0\right]\left[(l-p_{1})^{2}+\mathrm{i}0\right]}\,,$
(1)

where the external momenta are labelled as indicated by the Feynman diagram
depicted in figure 1. It is $p_{1}^{2}=p_{2}^{2}=p_{3}^{2}=0$ and
$q^{2}=(p_{1}+p_{2}+p_{3})^{2}\neq 0$. Note that we take all external momenta
to be incoming. We have also explicitly kept the causal $+\mathrm{i}0$ in the
propagators, which is necessary to determine the physical side of the branch
cuts that appear in the final result.333Here, $\mathrm{i}0$, with $0$ taken to
be positive, is an infinitesimal imaginary part indicating on which side of
the branch cut a multivalued function should be evaluated. The prefactor is
chosen as in ref. (Ellis, ).

To evaluate this loop integral, we first combine the propagators into one
generalized propagator of higher power using the well-known Feynman
parametrization (Peskin1995, , eq. (A.39))

$\frac{1}{abcd}\,=\,6\int_{0}^{1}\mathrm{d}x_{1}\,\mathrm{d}x_{2}\,\mathrm{d}x_{3}\,\mathrm{d}x_{4}\,\frac{\delta(1-x_{1}-x_{2}-x_{3}-x_{4})}{\left[ax_{1}\,+\,bx_{2}\,+\,cx_{3}\,+\,dx_{4}\right]^{4}}\,.$
(2)

Next, we complete the square with respect to $l$ and then shift the loop
momentum in the usual manner, which leads to

$\displaystyle\mathrm{D}_{0}\,$
$\displaystyle=\,\mu^{4-d}\int_{0}^{1}\mathrm{d}x_{1}\,\mathrm{d}x_{2}\,\mathrm{d}x_{3}\,\mathrm{d}x_{4}\,\delta(1-x_{1}-x_{2}-x_{3}-x_{4})$
$\displaystyle\hphantom{=}\,\times\,\frac{1}{\mathrm{i}\pi^{d/2}}\int\mathrm{d}^{d}l^{\prime}\frac{6}{\left[l^{\prime
2}+s_{1}x_{1}(x_{2}+x_{3})+s_{2}x_{3}(x_{1}+x_{4})+s_{3}x_{1}x_{3}+\mathrm{i}0\right]^{4}}\,,$
(3)

where we introduced the Mandelstam variables

$\displaystyle s_{1}\,$ $\displaystyle=\,(p_{1}+p_{2})^{2}\,,$ (4)
$\displaystyle s_{2}\,$ $\displaystyle=\,(p_{2}+p_{3})^{2}\,,$ (5)
$\displaystyle s_{3}\,$ $\displaystyle=\,(p_{1}+p_{3})^{2}\,.$ (6)

In the following, we treat all Mandelstam variables as real, but otherwise
arbitrary.

Evaluating the loop integral on the second line of eq. (3) with the help of
(Peskin1995, , eq. (A.44))

$\frac{1}{\mathrm{i}\pi^{d/2}}\int\mathrm{d}^{d}l\,\frac{1}{\left[l^{2}-\Delta+\mathrm{i}0\right]^{\alpha}}\,=\,(-1)^{\alpha}\frac{\Gamma(\alpha-d/2)}{\Gamma(\alpha)}(\Delta-\mathrm{i}0)^{d/2-\alpha}$
(7)

results in

$\displaystyle D_{0}\,$
$\displaystyle=\,\mu^{2\varepsilon}\,\Gamma(2+\varepsilon)\int_{0}^{1}\frac{\mathrm{d}x_{1}\,\mathrm{d}x_{2}\,\mathrm{d}x_{3}\,\mathrm{d}x_{4}\,\delta(1-x_{1}-x_{2}-x_{3}-x_{4})}{\left[-s_{1}x_{1}(x_{2}+x_{3})-s_{2}x_{3}(x_{1}+x_{4})-s_{3}x_{1}x_{3}-\mathrm{i}0\right]^{2+\varepsilon}}\,.$
(8)

Note that this expression is symmetric under exchanging $s_{1}\leftrightarrow
s_{2}$. To reflect this, we write the scalar box integral as
$D_{0}(s_{1},s_{2},q^{2})$, where $q^{2}=s_{1}+s_{2}+s_{3}$ is the squared
momentum of the off-shell external particle. To decouple the Feynman parameter
integrals, which are currently coupled via the delta function, we use the
substitution

$\displaystyle x_{1}\,$ $\displaystyle=\,\eta_{1}\xi_{1}\,,\quad
x_{4}\,=\,\eta_{1}(1-\xi_{1})\,,$ (9) $\displaystyle x_{3}\,$
$\displaystyle=\,\eta_{2}\xi_{2}\,,\quad x_{2}\,=\,\eta_{2}(1-\xi_{2})\,.$
(10)

The Jacobian of this substitution, which was inspired by ref.
(Smirnov_Feynman_Integrals, , page 44), is $\eta_{1}\eta_{2}$. Through the
substitution, the $\eta$\- and $\xi$-integrals factorize,

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)\,$
$\displaystyle=\,\mu^{2\varepsilon}\,\Gamma(2+\varepsilon)\int_{0}^{1}\mathrm{d}\eta_{1}\int_{0}^{1}\mathrm{d}\eta_{2}\;\eta_{1}^{-\varepsilon-1}\,\eta_{2}^{-\varepsilon-1}\,\delta(1-\eta_{1}-\eta_{2})$
$\displaystyle\hphantom{=}\,\times\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[-s_{1}\xi_{1}-s_{2}\xi_{2}-s_{3}\xi_{1}\xi_{2}-\mathrm{i}0\right]^{-\varepsilon-2}$
$\displaystyle=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[-s_{1}\xi_{1}-s_{2}\xi_{2}-s_{3}\xi_{1}\xi_{2}-\mathrm{i}0\right]^{-\varepsilon-2}\,.$
(11)

Here, the delta function was used to evaluate one of the $\eta$-integrals. The
remaining $\eta$-integral thus became a representation of the Beta function,
which was evaluated in terms of Gamma functions,

$B(a,b)\,\equiv\,\int_{0}^{1}\mathrm{d}t\,t^{a-1}(1-t)^{b-1}\,=\,\frac{\Gamma(a)\,\Gamma(b)}{\Gamma(a+b)}\,.$
(12)

For the resulting integral to safely converge, we must take
$\varepsilon<0$.444In the special cases of one or two (or three) vanishing
Mandelstam variables $\varepsilon<-1$ $(\varepsilon<-2)$ is required. We will
first evaluate the integral in this area of convergence and subsequently
analytically continue the end result to larger $\varepsilon$.

So far, our calculation did not differ much from the way the scalar box
integral was calculated before (see e. g. ref. (DuplancicNizic, )) and was
mainly presented for pedagogical reasons and introduction of notation. We will
now deviate from the usual path by factoring the term
$\frac{s_{1}s_{2}}{s_{3}}$ out of the integrand. This is of course not
sensible for ${\frac{s_{1}s_{2}}{s_{3}}=0,\pm\infty}$, that is
$s_{1},s_{2},s_{3}=0$. The special cases of one or several vanishing
Mandelstam variables will be considered in appendix A. Using the identity
(compare (DuplancicNizic, , eq. (19)))

$\displaystyle(a-\mathrm{i}0)^{\alpha}\,=\,(b-\mathrm{i}0)^{\alpha}\left(\frac{a}{b}-\mathrm{i}0\,\mathrm{sgn}\\!\left(b\right)\right)^{\alpha},\quad\text{where
}a\in\mathbb{R}\,,\;b\in\mathbb{R}\\!\setminus\\!\\{0\\}\,,\;\alpha\in\mathbb{C}\,,$
(13)

which is not completely trivial since $\alpha$ is not an integer, we obtain

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$
$\displaystyle=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\left(\frac{s_{1}s_{2}}{s_{3}}-\mathrm{i}0\right)^{\\!-\varepsilon-2}$
$\displaystyle\hphantom{=}\times\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[-\frac{s_{3}}{s_{2}}\xi_{1}-\frac{s_{3}}{s_{1}}\xi_{2}-\frac{s_{3}^{2}}{s_{1}s_{2}}\xi_{1}\xi_{2}-\mathrm{i}0\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}s_{2}}\right)\right]^{-\varepsilon-2}\,.$
(14)

The $\mathrm{i}0$ prescription determines on which side of the branch cut the
respective factors should be evaluated if their argument is negative. The
remaining integral on the second line only depends on the two dimensionless
variables

$x_{1}\,\equiv\,-\frac{s_{3}}{s_{1}}\,,\quad
x_{2}\,\equiv\,-\frac{s_{3}}{s_{2}}\,,$ (15)

which is why factoring the term $\frac{s_{1}s_{2}}{s_{3}}$ out of the
integrand is advantageous. By substituting
$\xi_{1}\rightarrow\zeta_{2}=x_{2}\xi_{1}$ and
$\xi_{2}\rightarrow\zeta_{1}=x_{1}\xi_{2}$, all dependence on the
dimensionless $x_{i}$ is moved into the integration boundaries,

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$
$\displaystyle=\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\,\frac{1}{s_{1}s_{2}}\left(\frac{s_{3}\mu^{2}}{s_{1}s_{2}}+\mathrm{i}0\right)^{\\!\varepsilon}I_{12}(x_{1},x_{2})\,,$
(16) $\displaystyle\text{where}\quad I_{12}(x_{1},x_{2})\,$
$\displaystyle\equiv\,\int_{0}^{x_{1}}\mathrm{d}\zeta_{1}\int_{0}^{x_{2}}\mathrm{d}\zeta_{2}\,\left[1-(1-\zeta_{1})(1-\zeta_{2})-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-2}\,.$
(17)

Note that we inverted the prefactor using

$(a\pm\mathrm{i}0)^{\alpha}\,=\,\left(\frac{1}{a}\mp\mathrm{i}0\right)^{\\!-\alpha},\quad\text{where
}a\in\mathbb{R}\\!\setminus\\!\\{0\\}\,.$ (18)

Additionally, we introduced the abbreviation

$\displaystyle\mathrm{sgn}_{123}\,\equiv\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}s_{2}}\right).$
(19)

And furthermore, we have written the integrand in a form suggesting the
substitution $\zeta_{i}\rightarrow 1-\zeta_{i}$. Subsequently substituting
$\zeta_{2}\rightarrow\zeta_{12}=1-\zeta_{1}\zeta_{2}$ allows us to evaluate
the resulting $\zeta_{12}$ integral,

$\displaystyle
I_{12}(x_{1},x_{2})\,=\,\int_{1}^{1-x_{1}}\mathrm{d}\zeta_{1}\int_{1}^{1-x_{2}}\mathrm{d}\zeta_{2}\left[1-\zeta_{1}\zeta_{2}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-2}$
$\displaystyle=\,-\int_{1}^{1-x_{1}}\frac{\mathrm{d}\zeta_{1}}{\zeta_{1}}\int_{1-\zeta_{1}}^{1-\zeta_{1}(1-x_{2})}\mathrm{d}\zeta_{12}\left[\zeta_{12}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-2}$
$\displaystyle=\,\frac{1}{1+\varepsilon}\int_{1}^{1-x_{1}}\frac{\mathrm{d}\zeta_{1}}{\zeta_{1}}\left\\{\left[1-\zeta_{1}(1-x_{2})-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-\left[1-\zeta_{1}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}\right\\}.$
(20)

While this integrand is finite at $\zeta_{1}=0$, we must regularize the
divergence at $\zeta_{1}=0$ if we split the integral into two separate
integrals. For this, we subtract and add $1$ and subsequently substitute
$\zeta_{1}\rightarrow\zeta=1-\zeta_{1}(1-x_{2})$ in the first integral and
$\zeta_{1}\rightarrow\zeta=1-\zeta_{1}$ in the second integral,

$\displaystyle I_{12}(x_{1},x_{2})\,$
$\displaystyle=\,\frac{1}{1+\varepsilon}\left\\{\int_{1}^{1-x_{1}}\frac{\mathrm{d}\zeta_{1}}{\zeta_{1}}\left(\left[1-\zeta_{1}(1-x_{2})-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-\,1\right)\right.$
$\displaystyle\phantom{=\,\frac{1}{1+\varepsilon}\\{}\left.-\int_{1}^{1-x_{1}}\frac{\mathrm{d}\zeta_{1}}{\zeta_{1}}\left(\left[1-\zeta_{1}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-\,1\right)\right\\}$
$\displaystyle=\,\frac{1}{1+\varepsilon}\left\\{-\int_{x_{2}}^{1-(1-x_{1})(1-x_{2})}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-\,1\right)\right.$
$\displaystyle\phantom{=\,\frac{1}{1+\varepsilon}\\{}\left.+\int_{0}^{x_{1}}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-1\right)\right\\}.$
(21)

Splitting the first of these integrals in two, we find that the scalar box
integral is given by

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)=\,$
$\displaystyle-\frac{1}{\varepsilon}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}\left(\frac{s_{3}\mu^{2}}{s_{1}s_{2}}+\mathrm{i}0\right)^{\varepsilon}$
$\displaystyle\times\left\\{\int_{0}^{x_{1}}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-1\right)\right.$
$\displaystyle\hphantom{\times\\{}\;\;+\int_{0}^{x_{2}}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-1\right)$
$\displaystyle\hphantom{\times\\{}\;\;\left.-\int_{0}^{1-(1-x_{1})(1-x_{2})}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-1\right)\right\\}.$
(22)

Note that we have rewritten the Gamma functions in the prefactor such that
they become unity in the limit $\varepsilon\rightarrow 0$, using
${\Gamma(z+1)=z\,\Gamma(z)}$.

Since the result for the internal-massless single off-shell scalar box
integral is commonly given in terms of Gauss hypergeometric functions, we will
identify them here, even though this will introduce spurious additional branch
cuts. The three integrals appearing in eq. (22) are of the form

$I(\chi)\,\equiv\,\int_{0}^{\chi}\frac{\mathrm{d}\zeta}{1-\zeta}\left(\left[\zeta-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon-1}-1\right).$
(23)

This is trivially zero if $\chi=0$. All three integrals have vanishing $\chi$
for $s_{3}=0$. As mentioned above, this special case of vanishing Mandelstam
variables will be considered separately in appendix A. The upper boundary of
the third integral,

$1-(1-x_{1})(1-x_{2})\,=\,-\frac{s_{3}q^{2}}{s_{1}s_{2}}\,,$ (24)

also vanishes if $q^{2}=0$. We will account for this by setting the third
integral to zero in case $q^{2}=0$. In case $\chi\neq 0$, we can substitute
$\zeta\rightarrow\chi^{-1}\zeta$,

$\displaystyle
I(\chi)\,=\,\int_{0}^{1}\frac{\mathrm{d}\zeta}{1-\chi\zeta}\left(\left[\chi-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}\zeta^{-\varepsilon-1}-\chi\right).$
(25)

As $\zeta>0$, we could factor it out of the square brackets here. To split
this integral into two separate integrals, we must regularize the divergence
that occurs for $\zeta\,=\,1/\chi$ if $\chi>1$. We regulate it by introducing
an infinitesimal imaginary part $\chi\rightarrow\chi+\mathrm{i}\tilde{0}$ in
the denominator, with $\tilde{0}$ taken to be positive.555Here, we write
$+\mathrm{i}\tilde{0}$ as an abbreviation for $\mathrm{i}\tilde{\eta}$ and
taking the limit $\tilde{\eta}\rightarrow+0$. This slightly unusual notation
is chosen in analogy to the causal $+\mathrm{i}0$ originating from the
propagators, where $+\mathrm{i}0$ is a common short hand notation for
$\mathrm{i}\eta$ and taking the limit $\eta\rightarrow 0$. The tilde marks the
new regulator as an analogous but independent limiting process. The
distinction between $\mathrm{i}0$ and $\mathrm{i}\tilde{0}$ allows us to
explicitly check that the end result does not depend on the regulator
introduced to split $I(\chi)$ into two integrals. Note that the end result
will not depend on the regulator $\mathrm{i}\tilde{0}$, since the original
integral $I(\chi)$ was not divergent for $\zeta\,=\,1/\chi$. Proceeding as
described, we obtain

$\displaystyle I(\chi)\,$
$\displaystyle=\left[\chi-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}\int_{0}^{1}\mathrm{d}\zeta\,\zeta^{-\varepsilon-1}\left(1-(\chi+\mathrm{i}\tilde{0})\zeta\right)^{-1}-\,\chi\int_{0}^{1}\mathrm{d}\zeta\,\frac{1}{1-(\chi+\mathrm{i}\tilde{0})\zeta}$
$\displaystyle=\,-\frac{1}{\varepsilon}\left[\chi-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\chi+\mathrm{i}\tilde{0}\right)+\,\ln(1-\chi-\mathrm{i}\tilde{0})\,.$
(26)

Here, the first integral became a representation of the Gauss hypergeometric
function discussed in appendix B (see eq. (77)), while the second integral
could be elementarily evaluated in terms of a logarithm. Both hypergeometric
function and logarithm are evaluated on the branch cut for $\chi>1$. The
regulator $\mathrm{i}\tilde{0}$ tells us which side of these branch cuts to
choose. Also note that the prefactor of the hypergeometric function vanishes
in case $\chi=0$, since we assumed $\varepsilon<0$ for $D_{0}$ to converge.

Applying eq. (26) to the term in curly braces in eq. (22) and using a
different regulator $\mathrm{i}\tilde{0}_{i}$ for each of the three integrals
gives

$\displaystyle-\,\frac{1}{\varepsilon}\left[x_{1}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x_{1}+\mathrm{i}\tilde{0}_{1}\right)$
$\displaystyle-\,\frac{1}{\varepsilon}\left[x_{2}-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x_{2}+\mathrm{i}\tilde{0}_{2}\right)$
$\displaystyle+\,\frac{1}{\varepsilon}\left[1-(1-x_{1})(1-x_{2})-\mathrm{i}0\,\mathrm{sgn}_{123}\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;1-(1-x_{1})(1-x_{2})+\mathrm{i}\tilde{0}_{3}\right)$
$\displaystyle+\,\ln(1-x_{1}-\mathrm{i}\tilde{0}_{1})\,+\,\ln(1-x_{2}-\mathrm{i}\tilde{0}_{2})\,-\,\ln\\!\left((1-x_{1})(1-x_{2})-\mathrm{i}\tilde{0}_{3}\right).$
(27)

Writing the logarithms on the last line of this expression in terms of their
real and imaginary parts, we see that only the imaginary parts yield a
contribution,

$\displaystyle\ln(1-x_{1}-\mathrm{i}\tilde{0}_{1})\,+\,\ln(1-x_{2}-\mathrm{i}\tilde{0}_{2})\,-\,\ln\\!\left((1-x_{1})(1-x_{2})-\mathrm{i}\tilde{0}_{3}\right)$
$\displaystyle=\,\mathrm{i}\pi\left[-\,\mathrm{sgn}\\!\left(\tilde{0}_{1}\right)\Theta(x_{1}-1)\,-\,\mathrm{sgn}\\!\left(\tilde{0}_{2}\right)\,\Theta(x_{2}-1)\right.$
$\displaystyle\left.\hphantom{=\,\mathrm{i}\pi}\,+\,\mathrm{sgn}\\!\left(\tilde{0}_{3}\right)\left\\{\Theta(x_{1}-1)\,\Theta(1-x_{2})\,+\,\Theta(1-x_{1})\,\Theta(x_{2}-1)\right\\}\right].$
(28)

Here, $\Theta(x)$ is the Heaviside step function. Let us now choose the sign
of the regulators $\mathrm{i}\tilde{0}_{i}$ such that the contribution of the
logarithms completely cancels. The kinematic regions in the $x_{1}x_{2}$-plane
in which the respective regulators yield a contribution are depicted in figure
2.

Figure 2: This plot shows in which kinematic regions the respective imaginary
parts $\mathrm{i}\tilde{0}_{i}$ of the three logarithms in eq. (28) yield a
contribution.

For the right-hand side of eq. (28) to vanish, the regulators must satisfy

$\displaystyle\mathrm{sgn}\\!\left(\tilde{0}_{1}\right)\,$
$\displaystyle\stackrel{{\scriptstyle!}}{{=}}\,\mathrm{sgn}\\!\left(\tilde{0}_{3}\right)\quad\text{if
}x_{1}>1\text{ and }x_{2}<1\text{ (yellow-green region)}\,,$ (29)
$\displaystyle\mathrm{sgn}\\!\left(\tilde{0}_{1}\right)\,$
$\displaystyle\stackrel{{\scriptstyle!}}{{=}}\,-\,\mathrm{sgn}\\!\left(\tilde{0}_{2}\right)\quad\text{if
}x_{1}>1\text{ and }x_{2}>1\text{ (yellow-blue region)}\,,$ (30)
$\displaystyle\mathrm{sgn}\\!\left(\tilde{0}_{2}\right)\,$
$\displaystyle\stackrel{{\scriptstyle!}}{{=}}\,\mathrm{sgn}\\!\left(\tilde{0}_{3}\right)\quad\text{if
}x_{1}<1\text{ and }x_{2}>1\text{ (blue-green region)}\,.$ (31)

In case both $x_{1,2}<1$ (white region), the right-hand side of eq. (28)
vanishes independently of the choice of regulators. Note that since these
conditions only determine the relative signs of the regulators, we can choose
the absolute sign of one regulator. Let us select
${\mathrm{i}\tilde{0}_{3}\equiv+\mathrm{i}\tilde{0}}$ everywhere. Then, both
$\mathrm{i}\tilde{0}_{1}$ and $\mathrm{i}\tilde{0}_{2}$ must change their sign
in the yellow-blue region in order to satisfy the conditions. To preserve
symmetry, we let them change signs along the dashed diagonal depicted in
figure 2, such that

$\displaystyle\mathrm{i}\tilde{0}_{1}\,$
$\displaystyle\equiv\,\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(x_{1}-x_{2}\right),$
(32) $\displaystyle\mathrm{i}\tilde{0}_{2}\,$
$\displaystyle\equiv\,\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(x_{2}-x_{1}\right),$
(33) $\displaystyle\mathrm{i}\tilde{0}_{3}\,$
$\displaystyle\equiv\,\mathrm{i}\tilde{0}\,.$ (34)

Combining our results, we find that the internal-massless single off-shell
scalar box integral is given by

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$
$\displaystyle=\,\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}\left(\frac{s_{3}\mu^{2}}{s_{1}s_{2}}+\mathrm{i}0\right)^{\\!\varepsilon}$
$\displaystyle\times\left\\{\left[-\frac{s_{3}}{s_{1}}-\mathrm{i}0\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}s_{2}}\right)\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}}{s_{1}}+\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{2}}-\frac{s_{3}}{s_{1}}\right)\right)\right.$
$\displaystyle\left.\hphantom{\times\\{}\,+\left[-\frac{s_{3}}{s_{2}}-\mathrm{i}0\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}s_{2}}\right)\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}}{s_{2}}+\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}}-\frac{s_{3}}{s_{2}}\right)\right)\right.$
$\displaystyle\left.\hphantom{\times\\{}\,-\left[-\frac{s_{3}q^{2}}{s_{1}s_{2}}-\mathrm{i}0\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}s_{2}}\right)\right]^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}q^{2}}{s_{1}s_{2}}+\mathrm{i}\tilde{0}\right)\right\\}.$
(35)

Remember that the term in curly braces is a function of the two variables
$x_{1}$ and $x_{2}$ defined in eq. (15). The global prefactor
$(\ldots)^{\varepsilon}$ and the prefactors $[\ldots]^{-\varepsilon}$ of the
hypergeometric functions can be recombined by first inverting the factor in
square brackets using eq. (18) and subsequently applying eq. (13) from right
to left, resulting in

$\displaystyle D_{0}\\!$
$\displaystyle\left(s_{1},s_{2},q^{2}\right)=\,\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}$
$\displaystyle\times\left\\{\left[\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right]^{\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}}{s_{1}}+\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{2}}-\frac{s_{3}}{s_{1}}\right)\right)\right.$
$\displaystyle\left.\hphantom{\times\\{}\,+\left[\frac{\mu^{2}}{-s_{1}-\mathrm{i}0}\right]^{\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}}{s_{2}}+\mathrm{i}\tilde{0}\,\mathrm{sgn}\\!\left(\frac{s_{3}}{s_{1}}-\frac{s_{3}}{s_{2}}\right)\right)\right.$
$\displaystyle\left.\hphantom{\times\\{}\,-\left[\frac{\mu^{2}}{-q^{2}-\mathrm{i}0}\right]^{\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\,-\frac{s_{3}q^{2}}{s_{1}s_{2}}+\mathrm{i}\tilde{0}\right)\right\\}.$
(36)

For convergence, we required $\varepsilon<0$ during the calculation. For the
expansion about the physical four dimensions, the result is analytically
continued to larger $\varepsilon$.

The result was calculated for
${s_{1},s_{2},s_{3}\in\mathbb{R}\\!\setminus\\!\\{0\\}}$. In appendix A, we
will see that it is also valid for vanishing Mandelstam variables in the sense
of appropriate limits. As discussed underneath eq. (24), the term on the last
line is equal to zero for ${q^{2}=0}$, since $\varepsilon<0$. Remember that in
dimensional regularization kinematic limits need to be taken before analytic
continuation in $\varepsilon$ to assure consistent results.

For the kinematic variables, the causal $\mathrm{i}0$ and the regulator
$\mathrm{i}\tilde{0}$ specify on which side of the occurring branch cuts the
respective functions are evaluated. Remember that by construction this result
does not depend on the regulator $\mathrm{i}\tilde{0}$. Hence, the
individually discontinuous imaginary parts of the hypergeometric functions
must cancel each other, and the branch cuts of the hypergeometric functions
are spurious. We will see that this is indeed the case in eq. (55) below.

In the region where $s_{1},s_{2},q^{2}<0$ and all three hypergeometric
functions are away from their branch cut, our result agrees with those found
in (Matsuura1989, , eq. (D.4)) and (FabiValery, , eq. (B.10)). The result in
(BernDixonKosower1, , eq. (4.36)) agrees upon application of
(Abramovitz_Stegun, , eq. (15.3.4)) to the hypergeometric functions.

Expanding the hypergeometric functions in $\varepsilon$ gives a result in
terms of polylogarithms (see eq. (78)). As long as one avoids the branch cuts
of the polylogarithms, the analytic continuation works via $s_{1}\rightarrow
s_{1}+\mathrm{i}0$, $s_{2}\rightarrow s_{2}+\mathrm{i}0$, $q^{2}\rightarrow
q^{2}+\mathrm{i}0$ as given in (BernDixonKosower1, ). However, treating the
branch cuts of the polylogarithms by this prescription leads to inconsistent
results. Remarkably, the prescription works correctly for the result expanded
up to $\mathcal{O}(\varepsilon^{0})$ if the Abel identity
(Lewin_Polylogs_Associated_Functions, ) is applied to the dilogarithms first.
Also the analytic continuation procedure introduced in (Beenakker1990, ) and
put forward in (Ellis, ) only applies to the dilogarithm as it is based on the
Euler identity (Lewin_Polylogs_Associated_Functions, ). It is worth mentioning
that all these prescriptions for analytic continuation make use of the causal
$+\mathrm{i}0$ and apply to each term individually. However from eq. (36) it
is clear that the imaginary parts of the hypergeometric functions respectively
polylogarithms are determined by a distinct regulator $\mathrm{i}\tilde{0}$
and are interdependent.

##  3 Epsilon expansion of the scalar box integral

The scalar box integral given in eq. (35) can be expressed as

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$
$\displaystyle=\,\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}\left|\frac{s_{3}\mu^{2}}{s_{1}s_{2}}\right|^{\\!\varepsilon}$
$\displaystyle\hphantom{=}\times\,\exp\\!\left[\mathrm{i}\pi\varepsilon\,\Theta\\!\left(-\frac{s_{3}}{s_{1}s_{2}}\right)\right]\mathcal{D}_{0}\\!\left(\varepsilon;-\frac{s_{3}}{s_{1}},-\frac{s_{3}}{s_{2}}\right),$
(37)

where $\mathcal{D}_{0}$ abbreviates the term in curly braces. By defining

$\displaystyle\mathcal{F}_{\pm}(\varepsilon;x)\,\equiv\,\left(x-\mathrm{i}0\,\mathrm{sgn}_{123}\right)^{\\!-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x\pm\mathrm{i}\tilde{0}\right),$
(38)

we can write

$\displaystyle\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})=\,\mathcal{F}_{\pm}\\!\left(\varepsilon;x_{1}\right)\\!\,+\,\mathcal{F}_{\mp}\\!\left(\varepsilon;x_{2}\right)\\!\,-\,\mathcal{F}_{+}\\!\left(\varepsilon;1-(1-x_{1})(1-x_{2})\right)\\!\,,$
(39)

where the upper sign choice is valid for $x_{1}\geq x_{2}$ and else the lower
sign choice. To expand the scalar box integral around $\varepsilon=0$, we need
the $\varepsilon$-expansion of
$\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$, which is performed in
section 3.1. Using this key ingredient, we subsequently discuss the expansion
of $\mathcal{D}_{0}\\!\left(\varepsilon;x_{1},x_{2}\right)$ in section 3.2 and
of $D_{0}(s_{1},s_{2},q^{2})$ in section 3.3. Tabulated results for the real
and imaginary parts of
$\exp\\!\left[\mathrm{i}\pi\varepsilon\,\Theta\\!\left(-\frac{s_{3}}{s_{1}s_{2}}\right)\right]\mathcal{D}_{0}\\!\left(\varepsilon;x_{1},x_{2}\right)$
in the various kinematic regions are given in table 1 in section 3.3.

###  3.1 Epsilon expansion of $\mathcal{F}_{\pm}(\varepsilon;x)$

To perform the $\varepsilon$-expansion of
$\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$ , which was defined in eq.
(38), we observe that
$\left(x-\mathrm{i}0\,\mathrm{sgn}_{123}\right)^{-\varepsilon}$ develops an
imaginary part for $x<0$, while the hypergeometric function does so for $x>1$.
Hence, there are 3 regions of interest:

  * 1.

$x\in(-\infty,0)$, where we need to expand
$\left(x-\mathrm{i}0\,\mathrm{sgn}_{123}\right)^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\right)$,

  * 2.

$x\in(0,1)$, where we need to expand
$x^{-\varepsilon}\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\right)$,

  * 3.

$x\in(1,\infty)$, where we need to expand
$x^{-\varepsilon}\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\pm\mathrm{i}\tilde{0}\right)$.

####  3.1.1 Expansion for $x\in(-\infty,0)$

For $x<0$, $\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$ reads

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,$
$\displaystyle=\left(x-\mathrm{i}0\,\mathrm{sgn}_{123}\right)^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\right)$
$\displaystyle=\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}}|x|^{-\varepsilon}\left(1-\sum_{n=1}^{\infty}\varepsilon^{n}\,\mathrm{Li}_{n}\\!\left(x\right)\right)$
$\displaystyle=\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}}\left[|x|^{-\varepsilon}-\sum_{k=0}^{\infty}\frac{(-\varepsilon)^{k}}{k!}\ln^{k}|x|\sum_{n=1}^{\infty}\varepsilon^{n}\,\mathrm{Li}_{n}\\!\left(x\right)\right],$
(40)

where the $\varepsilon$-expansion of the hypergeometric function given in eq.
(78) was used in the first step. Next, we substitute $n\rightarrow N=n+k$ and
then interchange the $N$\- and $k$-sum,

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,$
$\displaystyle=\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}}\left[|x|^{-\varepsilon}-\sum_{N=1}^{\infty}\varepsilon^{N}\sum_{k=0}^{N-1}\frac{(-1)^{k}\ln^{k}|x|}{k!}\,\mathrm{Li}_{N-k}\\!\left(x\right)\right].$
(41)

For $N\geq 2$, we now introduce the single-valued polylogarithms
$\mathcal{L}_{n}\\!\left(x\right)\,$ discussed in appendix C by adding and
subtracting $\frac{(-1)^{N-1}}{N!}\ln^{N-1}|x|\ln|1-x|$ in the sum. This
yields

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,$
$\displaystyle=\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}}\left[\vphantom{\frac{(-1)^{N-1}}{N!}}|x|^{-\varepsilon}\,-\,\varepsilon\,\mathrm{Li}_{1}\\!\left(x\right)\right.$
$\displaystyle\qquad\qquad\quad\left.\,-\,\sum_{N=2}^{\infty}\varepsilon^{N}\left(\mathcal{L}_{N}\\!\left(x\right)\,-\,\frac{(-1)^{N-1}}{N!}\ln^{N-1}|x|\ln|1-x|\right)\right]$
$\displaystyle=\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}}\left[1\,+\,\ln\frac{x}{x-1}\sum_{n=1}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\ln^{n-1}|x|-\sum_{n=2}^{\infty}\varepsilon^{n}\mathcal{L}_{n}\\!\left(x\right)\right].$
(42)

Note that since the single-valued polylogarithms are bounded on $\mathbb{R}$,
any divergences in the kinematic variable $x$ are explicitly contained in the
logarithms.

####  3.1.2 Expansion for $x\in(0,1)$

In the area where $0<x<1$, neither
$(x-\mathrm{i}0\,\mathrm{sgn}_{123})^{-\varepsilon}$ nor
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x\pm\mathrm{i}\tilde{0}\right)$
develop an imaginary part. Following the same steps as before, we find

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,$
$\displaystyle=\,x^{-\varepsilon}\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\right)$
$\displaystyle=\,|x|^{-\varepsilon}\left(1-\sum_{n=1}^{\infty}\varepsilon^{n}\,\mathrm{Li}_{n}\\!\left(x\right)\right)$
$\displaystyle=\,|x|^{-\varepsilon}\,-\,\varepsilon\,\mathrm{Li}_{1}\\!\left(x\right)\,-\,\sum_{n=2}^{\infty}\varepsilon^{n}\left(\mathcal{L}_{n}\\!\left(x\right)\,-\,\frac{(-1)^{n-1}}{n!}\ln^{n-1}|x|\ln|1-x|\right)$
$\displaystyle=1+\ln\frac{x}{1-x}\sum_{n=1}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\ln^{n-1}|x|-\sum_{n=2}^{\infty}\varepsilon^{n}\mathcal{L}_{n}\\!\left(x\right).$
(43)

Again, any divergences in the kinematic variable are contained in the
logarithms since the single-valued polylogarithms are bounded on $\mathbb{R}$.

####  3.1.3 Expansion for $x\in(1,\infty)$

In the area where $x>1$, the hypergeometric function
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\pm\mathrm{i}\tilde{0}\right)$
develops an imaginary part. To make it explicit, we employ the inversion
formula given in eq. (86),

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,$
$\displaystyle=\,x^{-\varepsilon}{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon,x\pm\mathrm{i}\tilde{0}\right)$
$\displaystyle=\,|x|^{-\varepsilon}\left[-{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{x}\right)\,+\,1\,+\,(-x\mp\mathrm{i}\tilde{0})^{\varepsilon}\frac{\pi\varepsilon}{\sin(\pi\varepsilon)}\right].$
(44)

By replacing $\varepsilon\rightarrow-\varepsilon$, we can use the
$\varepsilon$-expansion from eq. (78) for the resulting hypergeometric
function. Additionally writing the last term in terms of its real and
imaginary part, we obtain

$\displaystyle\mathcal{F}_{\pm}\\!\left(x;\varepsilon\right)\,$
$\displaystyle=\,|x|^{-\varepsilon}\,\sum_{n=1}^{\infty}(-\varepsilon)^{n}\,\mathrm{Li}_{n}\\!\left(\frac{1}{x}\right)\,+\,\pi\varepsilon\cot(\pi\varepsilon)\,\mp\,\mathrm{i}\pi\varepsilon\,.$
(45)

The $\varepsilon$-expansion of the cotangent term is given by
(Abramovitz_Stegun, , eqs. (4.3.70) & (23.1.18))

$\displaystyle\pi\varepsilon\cot(\pi\varepsilon)\,=\,1\,-\,2\sum_{n=1}^{\infty}\zeta_{2n}\,\varepsilon^{2n}\,=\,1-\sum_{n=2}^{\infty}\,[1+(-1)^{n}]\,\zeta_{n}\varepsilon^{n}\,,$
(46)

where $\zeta_{n}=\zeta(n)=\sum_{k=1}^{\infty}k^{-n}$ are the integer values of
the Riemann zeta function. Additionally expanding $|x|^{-\varepsilon}$ and
collecting orders of $\varepsilon$ yields

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,=\,$
$\displaystyle
1\,-\,\varepsilon\,\mathrm{Li}_{1}\\!\left(\frac{1}{x}\right)\,\mp\,\mathrm{i}\pi\varepsilon$
$\displaystyle-\,\sum_{n=2}^{\infty}\varepsilon^{n}\left[(-1)^{n-1}\sum_{k=0}^{n-1}\frac{(-1)^{k}}{k!}\,\ln^{k}\left|\frac{1}{x}\right|\,\mathrm{Li}_{n-k}\\!\left(\frac{1}{x}\right)+[1+(-1)^{n}]\zeta_{n}\right].$
(47)

Upon adding and subtracting the term
$\frac{(-1)^{n-1}}{n!}\ln^{n-1}\left|\frac{1}{x}\right|\ln\left|1-\frac{1}{x}\right|$,
we can identify the single-valued polylogarithms
$\mathcal{L}_{n}\\!\left(\frac{1}{x}\right)$, as defined in eq. (94).
Inverting them with the help of eq. (98), we find

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,=\,1\,-\,\varepsilon\,\mathrm{Li}_{n}\\!\left(\frac{1}{x}\right)\,\mp\,\mathrm{i}\pi\varepsilon\,-\,\sum_{n=2}^{\infty}\varepsilon^{n}\left[\mathcal{L}_{n}\\!\left(x\right)-\,\frac{(-1)^{n-1}}{n!}\ln^{n-1}\left|x\right|\ln\left|\frac{x-1}{x}\right|\right].$
(48)

Organizing this expression in terms of logarithms and single-valued
polylogarithms yields

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,=\,1\,+\,\ln\frac{x}{x-1}\sum_{n=1}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\ln^{n-1}|x|\,-\,\sum_{n=2}^{\infty}\varepsilon^{n}\mathcal{L}_{n}\\!\left(x\right)\,\mp\,\mathrm{i}\pi\varepsilon\,.$
(49)

####  3.1.4 Representation of the expansion of
$\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$ on the entire real axis

Combining the results found in eqs. (42), (43), and (49), we can give a
representation for the $\varepsilon$-expansion of
$\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$ valid for the entire
domain of $x\in\mathbb{R}$,

$\displaystyle\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!\,=\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\,\Theta(-x)}\mathfrak{F}(\varepsilon;x)\,\mp\,\mathrm{i}\pi\varepsilon\,\Theta(x-1)\,,$
(50) where
$\displaystyle\mathfrak{F}(\varepsilon;x)\,\equiv\,1\,+\,\ln\\!\left|\frac{x}{x-1}\right|\sum_{n=1}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\ln^{n-1}|x|\,-\,\sum_{n=2}^{\infty}\varepsilon^{n}\mathcal{L}_{n}\\!\left(x\right)\,.$
(51)

This compact expression is well suited to be used for the
$\varepsilon$-expansion of the complete box integral. The imaginary part is
made explicit and all functions of the real variable $x$ are manifestly real.
The complex exponential stems from
$\left(x-\mathrm{i}0\,\mathrm{sgn}_{123}\right)^{-\varepsilon}$, while the
imaginary part $\pi\varepsilon\,\Theta(x-1)$ originates from the
hypergeometric function.

An important feature of the expansion of $\mathfrak{F}$ as given in eq. (51)
is that it is finite in the limit $x\rightarrow\pm\infty$ in every order in
$\varepsilon$. Using the limits of the single-valued polylogarithms given in
appendix C, we obtain

$\displaystyle\lim_{x\rightarrow\infty}\mathfrak{F}(\varepsilon;x)\,$
$\displaystyle=\,1\,-\,2\sum_{n=1}^{\infty}\varepsilon^{2n}\zeta_{2n}\,,$ (52)
$\displaystyle\lim_{x\rightarrow-\infty}\mathfrak{F}(\varepsilon;x)\,$
$\displaystyle=1\,-\,2\sum_{n=1}^{\infty}\varepsilon^{2n}(2^{1-2n}-1)\zeta_{2n}\,.$
(53)

By collecting the $\ln|x|$ terms we can cast eq. (51) into the alternative
form

$\displaystyle\mathfrak{F}(\varepsilon;x)=|x|^{-\varepsilon}+\varepsilon\ln|1-x|-\sum_{n=2}^{\infty}\varepsilon^{n}\left[\frac{(-1)^{n}\ln|1-x|\ln^{n-1}|x|}{n!}+\mathcal{L}_{n}\\!\left(x\right)\right].$
(54)

From this we can read off that $\mathfrak{F}(\varepsilon;x)$ behaves like
$|x|^{-\varepsilon}$ near $x=0$ and it diverges logarithmically like
$\varepsilon\ln|1-x|$ near $x=1$ . In the complete box integral the $\ln|1-x|$
divergences always cancels.

###  3.2 Epsilon expansion of $\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})$

Figure 3: Physical interpretation of the kinematic regions to be considered
for the $\varepsilon$-expansion of $\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})$.
For DIS and Drell-Yan, the Mandelstam variables are $s>0$, $t<0$, and $u<0$.
In the case of SIA, all Mandelstam variables are positive. The virtuality
$q^{2}$ of the off-shell particle is negative for DIS and positive for both
Drell-Yan and SIA

By plugging the $\varepsilon$-expansion of
$\mathcal{F}_{\pm}\\!\left(\varepsilon;x\right)\\!$ from eq. (50) into eq.
(39), it is straightforward to obtain the $\varepsilon$-expansion of
$\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})$. We must consider seven areas
corresponding to different kinematics, which are depicted in figure 3. Since
$\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})$ is symmetric under exchanging
$x_{1}\leftrightarrow x_{2}$, only five of these areas are distinct. There are

  * 1.

DIS kinematics with $x_{1}$ and $x_{2}$ both positive,

  * 2.

DIS kinematics where $x_{1}$ and $x_{2}$ have different signs,

  * 3.

Drell-Yan kinematics with $x_{1}$ and $x_{2}$ both positive,

  * 4.

Drell-Yan kinematics where $x_{1}$ and $x_{2}$ have different signs,

  * 5.

SIA kinematics with $x_{1}$ and $x_{2}$ both negative.

Observing that

$\displaystyle\mathrm{i}\pi\varepsilon$
$\displaystyle\left[-\,\mathrm{sgn}\\!\left(x_{1}-x_{2}\right)\,\Theta(x_{1}-1)\,-\,\mathrm{sgn}\\!\left(x_{2}-x_{1}\right)\,\Theta(x_{2}-1)\,+\,\Theta(-(1-x_{1})(1-x_{2}))\right]=\,0\,,$
(55)

we see that the imaginary parts originating from the hypergeometric functions
cancel in all kinematic regions, as they should by construction. Hence, we
obtain

$\displaystyle\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})\,=\,$
$\displaystyle\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\Theta(-x_{1})}\mathfrak{F}(\varepsilon;x_{1})\,+\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\Theta(-x_{2})}\mathfrak{F}(\varepsilon;x_{2})$
$\displaystyle-\,\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\Theta(-x_{1}-x_{2}+x_{1}x_{2})}\mathfrak{F}(\varepsilon;x_{1}+x_{2}-x_{1}x_{2})\,.$
(56)

The $\varepsilon$-expansion of $\mathcal{D}_{0}$ in the first region, where
none of the three parts develops an imaginary part, is purely real. The
logarithmic poles at $x_{1,2}=1$ cancel between the three
$\mathfrak{F}$-functions, such that

$\displaystyle\mathcal{D}_{0}(\varepsilon;x_{1},x_{2})\,=\,\mathfrak{F}\\!\left(\varepsilon;x_{1}\right)\,+\,\mathfrak{F}\\!\left(\varepsilon;x_{2}\right)\,-\,\mathfrak{F}\\!\left(\varepsilon;x_{1}+x_{2}-x_{1}x_{2}\right)$
$\displaystyle=\,1\,-\,\varepsilon\ln\left|\frac{x_{1}x_{2}}{x_{1}+x_{2}-x_{1}x_{2}}\right|+\sum_{n=2}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\left[\ln\left|\frac{x_{1}}{x_{1}-1}\right|\ln^{n-1}|x_{1}|\right.$
$\displaystyle\phantom{=}\left.+\ln\left|\frac{x_{2}}{x_{2}-1}\right|\ln^{n-1}|x_{2}|-\ln\left|\frac{x_{1}+x_{2}-x_{1}x_{2}}{(1-x_{1})(1-x_{2})}\right|\ln^{n-1}|x_{1}+x_{2}-x_{1}x_{2}|\right]$
$\displaystyle\phantom{=}-\,\sum_{n=2}^{\infty}\varepsilon^{n}\left[\mathcal{L}_{n}\\!\left(x_{1}\right)+\mathcal{L}_{n}\\!\left(x_{2}\right)-\mathcal{L}_{n}\\!\left(x_{1}+x_{2}-x_{1}x_{2}\right)\right].$
(57)

For the order $\varepsilon^{2}$ coefficient, we can use Abel’s formula (115),
in combination with the Euler identity (102) and inverted Landen identity
(114), to reduce the number of appearing $\mathcal{L}_{2}$ to just two,

$\displaystyle\mathcal{L}_{2}\\!\left(x_{1}\right)+\mathcal{L}_{2}\\!\left(x_{2}\right)-\mathcal{L}_{2}\\!\left(x_{1}+x_{2}-x_{1}x_{2}\right)=\mathcal{L}_{2}\\!\left(\frac{x_{2}(x_{1}-1)}{x_{1}}\right)+\mathcal{L}_{2}\\!\left(\frac{x_{1}(x_{2}-1)}{x_{2}}\right)+\zeta_{2}\,,$
(58)

where we assumed $x_{1},\,x_{2},\,(x_{1}+x_{2}-x_{1}x_{2})>0$. For the higher
orders in $\varepsilon$ there presumably exists no representation with fewer
than three (single-valued) polylogarithms, since only for the dilogarithm
there exists a suitable five term relation (see appendix C.3).

###  3.3 Epsilon expansion of $D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$

Combining eqs. (37) and (56), we see that real and imaginary part of
$D_{0}(s_{1},s_{2},q^{2})$ are determined by the phase factors
$\exp\left[{\mathrm{i}\pi\varepsilon\,\Theta\\!\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\right]\\!\cdot\exp\left[{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\Theta(-y_{i})}\right]$,
where ${y_{i}=-\frac{s_{3}s_{i}}{s_{1}s_{2}}}$ with $s_{i}=s_{1},s_{2},q^{2}$.
We can simplify this factor to

$\displaystyle\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\cdot\mathrm{e}^{\mathrm{i}\pi\varepsilon\,\mathrm{sgn}_{123}\Theta(-y_{i})}=\Theta(-s_{i})+\Theta(s_{i})\,\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,.$
(59)

Therefore,
$\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\mathcal{D}_{0}(\varepsilon,x_{1},x_{2})$
admits for the representation

$\displaystyle\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\mathcal{D}_{0}(\varepsilon,x_{1},x_{2})\,=$
$\displaystyle\left(\Theta(-s_{2})+\Theta(s_{2})\,\mathrm{e}^{\mathrm{i}\pi\varepsilon}\right)\mathfrak{F}(\varepsilon;x_{1})$
$\displaystyle+\left(\Theta(-s_{1})+\Theta(s_{1})\,\mathrm{e}^{\mathrm{i}\pi\varepsilon}\right)\mathfrak{F}(\varepsilon;x_{2})$
$\displaystyle-\left(\Theta(-q^{2})+\Theta(q^{2})\,\mathrm{e}^{\mathrm{i}\pi\varepsilon}\right)\mathfrak{F}(\varepsilon;x_{1}+x_{2}-x_{1}x_{2})\,.$
(60)

Hence, the imaginary part is determined by the signs of $s_{1}$, $s_{2}$ and
$q^{2}$. Explicit results for the real and imaginary part in all areas are
given in table 1.

$s_{1}$ | $s_{2}$ | $q^{2}$ | $\mathrm{Re}\\!\left(\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\mathcal{D}_{0}(\varepsilon,x_{1},x_{2})\right)$ | $\mathrm{Im}\\!\left(\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\mathcal{D}_{0}(\varepsilon,x_{1},x_{2})\right)$
---|---|---|---|---
$-$ | $-$ | $-$ | $\mathfrak{F}_{1}+\mathfrak{F}_{2}-\mathfrak{F}_{12}$ | 0
$-$ | $-$ | + | $\mathfrak{F}_{1}+\mathfrak{F}_{2}-\mathfrak{F}_{12}\cos\pi\varepsilon$ | $-\mathfrak{F}_{12}\sin\pi\varepsilon$
$-$ | + | $-$ | $\mathfrak{F}_{2}-\mathfrak{F}_{12}+\mathfrak{F}_{1}\cos\pi\varepsilon$ | $\mathfrak{F}_{1}\sin\pi\varepsilon$
+ | $-$ | $-$ | $\mathfrak{F}_{1}-\mathfrak{F}_{12}+\mathfrak{F}_{2}\cos\pi\varepsilon$ | $\mathfrak{F}_{2}\sin\pi\varepsilon$
$-$ | + | + | $\mathfrak{F}_{2}+(\mathfrak{F}_{1}-\mathfrak{F}_{12})\cos\pi\varepsilon$ | $(\mathfrak{F}_{1}-\mathfrak{F}_{12})\sin\pi\varepsilon$
+ | $-$ | + | $\mathfrak{F}_{1}+(\mathfrak{F}_{2}-\mathfrak{F}_{12})\cos\pi\varepsilon$ | $(\mathfrak{F}_{2}-\mathfrak{F}_{12})\sin\pi\varepsilon$
+ | + | $-$ | $-\mathfrak{F}_{12}+(\mathfrak{F}_{1}+\mathfrak{F}_{2})\cos\pi\varepsilon$ | $(\mathfrak{F}_{1}+\mathfrak{F}_{2})\sin\pi\varepsilon$
+ | + | + | $\left(\mathfrak{F}_{1}+\mathfrak{F}_{2}-\mathfrak{F}_{12}\right)\cos\pi\varepsilon$ | $\left(\mathfrak{F}_{1}+\mathfrak{F}_{2}-\mathfrak{F}_{12}\right)\sin\pi\varepsilon$
Table 1: Real and imaginary part of
$\mathrm{e}^{\mathrm{i}\pi\varepsilon\Theta\left(-\frac{s_{3}}{s_{1}s_{2}}\right)}\mathcal{D}_{0}(\varepsilon,x_{1},x_{2})$
in the areas specified by the signs of $s_{1}$, $s_{2}$ and $q^{2}$. For
compactness we introduce the abbreviations
$\mathfrak{F}_{1}\equiv\mathfrak{F}(\varepsilon;x_{1})$,
$\mathfrak{F}_{2}\equiv\mathfrak{F}(\varepsilon;x_{2})$ and
$\mathfrak{F}_{12}\equiv\mathfrak{F}(\varepsilon;x_{1}+x_{2}-x_{1}x_{2})$.

Note that the superficially divergent terms $\ln(1-x_{1})$ and $\ln(1-x_{2})$
cancel between $\mathfrak{F}(\varepsilon;x_{1,2})$ and
$\mathfrak{F}(\varepsilon;x_{1}+x_{2}-x_{1}x_{2})$ in the regions where
$x_{1,2}=1$.

Combining the results listed in table 1 and eq. (37) we can directly infer
real and imaginary parts of the complete box integral
$D_{0}(s_{1},s_{2},q^{2})$. For completing the $\varepsilon$-expansion we can
employ (Abramovitz_Stegun, , eq. (6.1.33))

$\displaystyle\Gamma(1+\varepsilon)=\exp\\!\left[-\gamma_{\text{E}}\,\varepsilon+\sum_{n=2}^{\infty}\frac{(-1)^{n}\zeta_{n}\,\varepsilon^{n}}{n}\right]$
(61)

for all Gamma functions in the prefactor, where $\gamma_{\text{E}}$ is the
Euler-Mascheroni constant. However, since the same factor of Gamma functions
multiplies all one-loop integrals and a similar factor is present in
corresponding phase-space integrals it is often not required to expand all of
them explicitly. Expanding the term
$\left|\frac{s_{3}\mu^{2}}{s_{1}s_{2}}\right|^{\varepsilon}$ introduces
additional logarithmic terms in the $\varepsilon$-expansion. One might albeit
need to keep this term unexpanded to regularize singular behavior in the
kinematic limits $s_{1},s_{2},s_{3}\rightarrow 0$. This is exemplified in
appendix D.

To sum up our result, we can cast $D_{0}(s_{1},s_{2},q^{2})$ as given in eq.
(35) into the form

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)=\,$
$\displaystyle\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}\left|\frac{s_{3}\mu^{2}}{s_{1}s_{2}}\right|^{\\!\varepsilon}$
$\displaystyle\times\,$
$\displaystyle\left[(\Theta(-s_{2})+\Theta(s_{2})\mathrm{e}^{\mathrm{i}\pi\varepsilon})\,\mathfrak{F}(\varepsilon;x_{1})+(\Theta(-s_{1})+\Theta(s_{1})\mathrm{e}^{\mathrm{i}\pi\varepsilon})\,\mathfrak{F}(\varepsilon;x_{2})\right.$
$\displaystyle\left.-(\Theta(-q^{2})+\Theta(q^{2})\mathrm{e}^{\mathrm{i}\pi\varepsilon})\,\mathfrak{F}(\varepsilon;x_{1}+x_{2}-x_{1}x_{2})\right].$
(62)

Plugging in the $\varepsilon$-expansion of $\mathfrak{F}(\varepsilon;x)$ from
eq. (51) we have an all-order $\varepsilon$-expansion of
$D_{0}(s_{1},s_{2},q^{2})$ with explicit real and imaginary parts valid in the
entire kinematic domain of $s_{1},s_{2},q^{2}$. Since the term in square
brackets is finite for vanishing $s_{1,2}$, $D_{0}(s_{1},s_{2},q^{2})$ behaves
like $|s_{1,2}|^{-1-\varepsilon}$ in the limits $s_{1,2}\rightarrow 0$. For a
comprehensive discussion of the behavior of $\mathfrak{F}(\varepsilon;x)$ in
different kinematic limiting cases we refer to section 3.1.4.

By bringing the factor
$\left|\frac{s_{3}\mu^{2}}{s_{1}s_{2}}\right|^{\varepsilon}$ into the
brackets, we obtain a form analogous to eq. (36),

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}\right)$
$\displaystyle=\,\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}$
$\displaystyle\times$
$\displaystyle\left[\left(\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right)^{\varepsilon}\,\left|\frac{s_{3}}{s_{1}}\right|^{\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{s_{3}}{s_{1}}\right)+\left(\frac{\mu^{2}}{-s_{1}-\mathrm{i}0}\right)^{\varepsilon}\,\left|\frac{s_{3}}{s_{2}}\right|^{\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{s_{3}}{s_{2}}\right)\right.$
$\displaystyle\left.-\left(\frac{\mu^{2}}{-q^{2}-\mathrm{i}0}\right)^{\varepsilon}\,\left|\frac{s_{3}q^{2}}{s_{1}s_{2}}\right|^{\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{s_{3}q^{2}}{s_{1}s_{2}}\right)\right].$
(63)

Eqs. (36) and (63) differ only by the replacement
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x\pm\mathrm{i}\tilde{0}\right)\rightarrow|x|^{\varepsilon}\,\mathfrak{F}(\varepsilon;x)$.
This can be understood as replacing the hypergeometric function by a suitable
single-valued version. Hence, in contrast to eq. (36), the new eq. (63) does
no longer suffer from spurious branch cuts.

##  4 Conclusion

We have revisited the massless scalar box integral with one external off-shell
particle, keeping the causal $+\mathrm{i}0$ throughout the calculation. The
main new result is extending the explicit determination of real and imaginary
parts beyond finite order in the dimensional regularization parameter epsilon
in all kinematic regions. The epsilon expansion is expressed to all orders in
terms of newly introduced single-valued polylogarithms. This makes our result
free of the spurious branch cuts common to representations using ordinary
polylogarithms. The single-valued polylogarithms are bounded, such that any
divergent behavior in kinematic limits is contained in logarithmic terms.
Since we have determined the box integral in all kinematic regions, our result
is suitable for all kinds of calculations in higher orders of perturbation
theory, where imaginary parts or higher orders in epsilon become important.

The calculation method put forward in this work can be generalized to the
scalar box integral with two non-adjacent end points off the light cone. A
corresponding paper is in preparation (Haug2023, ).

###### Acknowledgements.

We are grateful to Werner Vogelsang for valuable comments regarding this work.
We thank the referee for raising the question if generalizing the result to
more end points off the light cone is possible. This study was supported in
part by Deutsche Forschungsgemeinschaft (DFG) through the Research Unit FOR
2926 (Project No. 409651613). J. H. is grateful to the
Landesgraduiertenförderung Baden-Württemberg for supporting her research.

##  Appendix A Scalar box integral in case of vanishing Mandelstam variables

In this appendix, we evaluate the internal-massless single off-shell scalar
box integral in case one or several Mandelstam variables vanish, since the
result calculated in the main text was obtained for non-vanishing Mandelstam
variables. For this, we go back to $D_{0}$ as given in eq. (11) and explicitly
set the vanishing Mandelstam variables to zero. We will find that all cases
agree with appropriate limits of the general result.

###  A.1 $s_{1}=0$ or $s_{2}=0$

Since the scalar box integral is symmetric under exchanging
$s_{1}\leftrightarrow s_{2}$, we can without loss of generality consider
$s_{1}=0$. In this case, the integrals in eq. (11) factorize,

$\displaystyle D_{0}\\!\left(s_{1}=0,s_{2},q^{2}=s_{2}+s_{3}\right)$
$\displaystyle=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\\!\int_{0}^{1}\\!\\!\mathrm{d}\xi_{1}\\!\int_{0}^{1}\\!\\!\mathrm{d}\xi_{2}\,\left[\xi_{2}\left(-s_{2}-s_{3}\xi_{1}-\mathrm{i}0\right)\right]^{-\varepsilon-2}$
$\displaystyle=\,\frac{2}{\varepsilon}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\mu^{2\varepsilon}\\!\int_{0}^{1}\\!\\!\mathrm{d}\xi_{1}\left[-s_{2}-s_{3}\xi_{1}-\mathrm{i}0\right]^{-\varepsilon-2}.$
(64)

Here the $\xi_{2}$-integration was performed elementarily. Note that
$\varepsilon<-1$ was imposed for convergence. While the $\xi_{1}$-integral can
in principle be evaluated elementarily as well, first substituting
${\xi_{1}\rightarrow\zeta=-s_{2}-s_{3}\xi_{1}}$ is useful for discussing the
additional limits $q^{2},s_{2}=0$. Splitting the resulting integral into two
separate integrals, we obtain

$\displaystyle D_{0}\\!\left(s_{1}=0,s_{2},q^{2}=s_{2}+s_{3}\right)$
$\displaystyle=\,\frac{1}{\varepsilon}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{3}}\,\mu^{2\varepsilon}$
$\displaystyle\hphantom{=}\times\left\\{\int_{0}^{-s_{2}}\mathrm{d}\zeta\left[\zeta-\mathrm{i}0\right]^{-\varepsilon-2}-\int_{0}^{-q^{2}}\mathrm{d}\zeta\left[\zeta-\mathrm{i}0\right]^{-\varepsilon-2}\right\\}$
$\displaystyle=\,\frac{1}{\varepsilon(1+\varepsilon)}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{3}}$
$\displaystyle\hphantom{=}\times\left\\{\frac{1}{s_{2}}\left(\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right)^{\\!\varepsilon}-\frac{1}{q^{2}}\left(\frac{\mu^{2}}{-q^{2}-\mathrm{i}0}\right)^{\\!\varepsilon}\right\\}.$
(65)

Note that the second integral and thus the second term on the last line
vanishes for $q^{2}=s_{2}+s_{3}=0$. The result can be analytically continued
to values of $\varepsilon$ beyond the original area of convergence.

Let us compare this result to the limit $s_{1}\rightarrow 0$ of the general
result as given eq. (36). The limit of the first and last term in the general
result are of the form

$\lim_{s_{1}\rightarrow
0}\;\frac{1}{s_{1}}\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;\frac{a}{s_{1}}\pm\mathrm{i}\tilde{0}\right),\quad\text{where
}a\in\mathbb{R}\\!\setminus\\!\\{0\\}\,.$ (66)

To evaluate this limit, we first apply the inversion formula for the
hypergeometric function derived in eq. (86) in the appendix, which yields

$\displaystyle\text{eq.\,\eqref{eq:first_and_third_limit_s1_to_zero}}\,=\,\lim_{s_{1}\rightarrow
0}\;\frac{1}{s_{1}}\left\\{-{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{s_{1}}{a}\right)\,+\,1\,+\,\left(\frac{-a}{s_{1}}\mp\mathrm{i}\tilde{0}\right)^{\\!\varepsilon}\frac{\pi\varepsilon}{\sin(\pi\varepsilon)}\right\\}.$
(67)

Imposing the condition $\varepsilon<-1$, which is needed for convergence in
case $s_{1}=0$, the second term vanishes. The remaining limit can be evaluated
using L’Hôpital’s rule,

$\displaystyle\text{eq.\,\eqref{eq:first_and_third_limit_s1_to_zero}}\,=\,\lim_{s_{1}\rightarrow
0}\;\left\\{-\,\frac{1}{a}\left.\frac{\mathrm{d}}{\mathrm{d}x}\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;x\right)\right|_{x=\frac{s_{1}}{a}}\right\\}\,=\,-\,\frac{1}{a}\,\frac{\varepsilon}{1+\varepsilon}\,.$
(68)

The derivative of the hypergeometric function was read off the series
expansion given in eq. (76). The second term in the general result also
vanishes for $s_{1}\rightarrow 0$ since $\varepsilon<-1$. Combining
everything, we find that taking the limit $s_{1}\rightarrow 0$ of the general
result given in eq. (36) is equivalent to setting $s_{1}=0$ in eq. (11).

###  A.2 $s_{1}=0$ and $\boldmath{s_{2}=0}$

In this case, we set $s_{2}=0$ in eq. (64) and obtain

$D_{0}\\!\left(s_{1}=0,s_{2}=0,q^{2}=s_{3}\right)\,=\,-\frac{1}{\varepsilon(1+\varepsilon)}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{3}^{2}}\left(\frac{\mu^{2}}{-s_{3}-\mathrm{i}0}\right)^{\\!\varepsilon}.$
(69)

This is identical to the expression obtained by taking the limit
$s_{2}\rightarrow 0$ of eq. (65), i.e. taking the limits
$s_{1},s_{2}\rightarrow 0$ of eq. (36).

###  A.3 $s_{3}=0$

In this case, eq. (11) reads

$\displaystyle D_{0}\\!\left(s_{1}=0,s_{2},q^{2}=s_{2}\right)$
$\displaystyle=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[-s_{1}\xi_{1}-s_{2}\xi_{2}-\mathrm{i}0\right]^{-\varepsilon-2}$
$\displaystyle=\,\frac{1}{\varepsilon}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{2}}\mu^{2\varepsilon}$
$\displaystyle\hphantom{=}\times\int_{0}^{1}\mathrm{d}\xi_{1}\left\\{\left[-s_{1}\xi_{1}-\mathrm{i}0\right]^{-\varepsilon-1}-\left[-s_{1}\xi_{1}-s_{2}-\mathrm{i}0\right]^{-\varepsilon-1}\right\\},$
(70)

where the $\xi_{2}$-integration was performed elementarily. Note that
$\varepsilon<-1$ was imposed for convergence. Next, we substitute
${\xi_{1}\rightarrow\zeta=-s_{1}\xi_{1}}$ in the first term and
${\xi_{1}\rightarrow\zeta=-s_{1}\xi_{1}-s_{2}}$ in the second term, and
subsequently split the latter integral into two separate integrals,

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}=s_{1}+s_{2}\right)$
$\displaystyle=\,\frac{1}{\varepsilon}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}\,\mu^{2\varepsilon}\left\\{-\int_{0}^{-s_{1}}\mathrm{d}\zeta\left[\zeta-\mathrm{i}0\right]^{-\varepsilon-1}\right.$
$\displaystyle\hphantom{=}\left.-\int_{0}^{-s_{2}}\mathrm{d}\zeta\left[\zeta-\mathrm{i}0\right]^{-\varepsilon-1}+\int_{0}^{-q^{2}}\mathrm{d}\zeta\left[\zeta-\mathrm{i}0\right]^{-\varepsilon-1}\right\\}.$
(71)

Note that the last integral vanishes for $q^{2}=0$. Evaluating the remaining
$\zeta$-integrals yields

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}=s_{1}+s_{2}\right)=\,$
$\displaystyle\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{1}s_{2}}$
$\displaystyle\times\left\\{\left(\frac{\mu^{2}}{-s_{1}-\mathrm{i}0}\right)^{\\!\varepsilon}+\left(\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right)^{\\!\varepsilon}-\left(\frac{\mu^{2}}{-q^{2}-\mathrm{i}0}\right)^{\\!\varepsilon}\right\\},$
(72)

which can be analytically continued beyond $\varepsilon<-1$ which was
originally required for convergence. One can easily see that this result is
equivalent to setting $s_{3}=0$ in the general result as given in eq. (36),
since the resulting hypergeometric functions evaluate to
${{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;0\right)=1}$
according to eq. (76).

###  A.4 $s_{1}=s_{3}=0$ or $s_{2}=s_{3}=0$

As above, we can without loss of generality consider $s_{1}=s_{3}=0$. In this
case, eq. (11) yields

$\displaystyle D_{0}\\!\left(s_{1},s_{2},q^{2}=s_{1}+s_{2}\right)$
$\displaystyle=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[(-s_{2}-\mathrm{i}0)\xi_{2}\right]^{-\varepsilon-2}$
$\displaystyle=\,\frac{1}{\varepsilon}\frac{\Gamma(1+\varepsilon)\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{2}^{2}}\left(\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right)^{\\!\varepsilon}\,.$
(73)

As in the previous cases, one has to impose $\varepsilon<-1$ for convergence
before analytically continuing to larger $\varepsilon$.

Let us compare this expression to the corresponding limit of the general
result as given in eq. (36). For $s_{3}=0$, the hypergeometric functions
evaluate to
${{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;0\right)=1}$.
The limit remaining to be taken is

$\displaystyle\lim_{s_{1},s_{3}\rightarrow
0}D_{0}\\!\left(s_{1},s_{2},q^{2}\right)\,$
$\displaystyle=\,\frac{1}{\varepsilon^{2}}\frac{\Gamma(1+\varepsilon)\,\Gamma^{2}(1-\varepsilon)}{\Gamma(1-2\varepsilon)}\,\frac{2}{s_{2}}$
$\displaystyle\times\,\left\\{\lim_{s_{1}\rightarrow
0}\frac{\left[\frac{\mu^{2}}{-s_{2}-\mathrm{i}0}\right]^{\varepsilon}\,-\left[\frac{\mu^{2}}{-s_{1}-s_{2}-\mathrm{i}0}\right]^{\varepsilon}}{s_{1}}\right.\;\left.+\,\lim_{s_{1}\rightarrow
0}\frac{\left[\frac{\mu^{2}}{-s_{1}-\mathrm{i}0}\right]^{\varepsilon}}{s_{1}}\right\\}.$
(74)

Imposing the condition $\varepsilon<-1$, which is needed for convergence in
case $s_{1}=0$, the second term vanishes. Applying L’Hôpital’s rule to the
first limit, we find that the limit ${s_{1},s_{3}\rightarrow 0}$ of the
general result given in eq. (36) is equivalent to setting $s_{1}=s_{3}=0$ in
eq. (11).

###  A.5 $s_{1}=s_{2}=s_{3}=0$

Setting all Mandelstam variables to $0$ in eq. (11) and requiring
$\varepsilon<-2$ for convergence yields

$D_{0}\\!\left(s_{1}=0,s_{2}=0,q^{2}=0\right)\,=\,\mu^{2\varepsilon}\,\frac{\Gamma(2+\varepsilon)\Gamma^{2}(-\varepsilon)}{\Gamma(-2\varepsilon)}\int_{0}^{1}\mathrm{d}\xi_{1}\int_{0}^{1}\mathrm{d}\xi_{2}\,\left[0\right]^{-\varepsilon-2}\,=\,0\,.$
(75)

This integral also has to vanish since the integrand has a mass dimension of
$-\varepsilon-2$, but does not depend on any mass scale. Taking the limit
${s_{1,2,3}\rightarrow 0}$ of the general result given in eq. (36) leads to
the same conclusion.

##  Appendix B Gauss hypergeometric function
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$

For $|z|<1$ the Gauss hypergeometric function is defined by the series
(Abramovitz_Stegun, , eq. (15.1.1))

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(a,b,c;z\right)=\sum_{n=0}^{\infty}\frac{(a)_{n}(b)_{n}}{(c)_{n}}\,\frac{z^{n}}{n!}\,,$
(76)

where $(a)_{n}=\frac{\Gamma(a+n)}{\Gamma(a)}$ is the Pochhammer symbol. It is
analytically continued by the integral representation (Abramovitz_Stegun, ,
eq. (15.3.1))

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(a,b,c;z\right)\,\equiv\,$
$\displaystyle\frac{\Gamma(c)}{\Gamma(b)\,\Gamma(c-b)}\int_{0}^{1}\mathrm{d}t\,t^{b-1}\,(1-t)^{c-b-1}\,(1-zt)^{-a}\,,$
(77)
$\displaystyle\text{where}\;\mathrm{Re}(c)>\mathrm{Re}(b)>0\,,\;z\in\mathbb{C}\\!\setminus\\!\mathbb{R}_{\geq
1}\,.$

Note that the integral may be ill-defined for $z>1$, since in this case the
term ${(1-zt)}$ changes sign within the range of integration. This change of
sign can, in general, give rise to a discontinuity or branch cut. In the
following, we discuss selected properties (including the branch cut structure)
of the specific hypergeometric function
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$ which
appears in the internal-massless single off-shell scalar box integral.

###  B.1 Epsilon expansion of
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$

The $\varepsilon$-expansion of this hypergeometric function was derived in
ref. (FabiValery, , eq. (B.28)),

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)\,=\,1\,-\,\sum_{n=1}^{\infty}\varepsilon^{n}\,\mathrm{Li}_{n}\\!\left(z\right).$
(78)

Since the polylogarithms $\mathrm{Li}_{n}\\!\left(z\right)$ have a branch cut
along the positive real axis greater than 1
(Lewin_Polylogs_Associated_Functions, ), the hypergeometric function also has
a branch cut there. The branch cut structure will be discussed in section B.3.

###  B.2 Inversion relation for
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$

An equation relating the hypergeometric function
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$ to
another hypergeometric function with argument $\frac{1}{z}$ is useful because
it allows us to transform the argument away from the branch cut. To find this
inversion relation, we write the hypergeometric function in terms of its
$\varepsilon$-expansion and invert the appearing polylogarithms using the
inversion relation for the polylogarithm given in ref.
(Lewin_Polylogs_Associated_Functions, , eq. (7.20)),

$\displaystyle\mathrm{Li}_{n}\\!\left(z\right)\,+\,(-1)^{n}\,\mathrm{Li}_{n}\\!\left(\frac{1}{z}\right)\,=\,-\frac{\ln^{n}(-z)}{n!}\,+\,2\sum_{k=1}^{\lfloor
n/2\rfloor}\frac{\ln^{n-2k}(-z)}{(n-2k)!}\,\mathrm{Li}_{2k}\\!\left(-1\right)\,,$
(79)

which holds for $z\in\mathbb{C}\setminus[0,\infty)$.666For reference, the
corresponding equation on the branch cut, i.e. for $x\in(0,\infty)$, is
$\displaystyle\mathrm{Li}_{n}\\!\left(x\pm\mathrm{i}0\right)\,+\,(-1)^{n}\,\mathrm{Li}_{n}\\!\left(\frac{1}{x\pm\mathrm{i}0}\right)\,=\,-\frac{\ln^{n}(x)}{n!}\,+\,2\sum_{k=1}^{\lfloor
n/2\rfloor}\zeta_{2k}\frac{\ln^{n-2k}(x)}{(n-2k)!}\,\pm\,\mathrm{i}\pi\,\frac{\ln^{n-1}(x)}{(n-1)!}\,.$
Plugging the second term on the left into eq. (78) leads to the closely
related $\varepsilon$-expansion of
${{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{z}\right)$,
while plugging the second term on the right into eq. (78) leads to the series
expansion of an exponential function, such that

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)+\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{z}\right)$
$\displaystyle=$
$\displaystyle\,1\,+\,(-z)^{\varepsilon}-\,2\sum_{n=1}^{\infty}\varepsilon^{n}\sum_{k=1}^{\lfloor
n/2\rfloor}\frac{\ln^{n-2k}(-z)}{(n-2k)!}\mathrm{Li}_{2k}\\!\left(-1\right)\,.$
(80)

To simplify the term on the second line of this expression, we interchange the
two sums and subsequently shift the summation index to $N=n-2k$, such that the
$N$-sum becomes an exponential function,

$\sum_{n=1}^{\infty}\varepsilon^{n}\sum_{k=1}^{\lfloor
n/2\rfloor}\frac{\ln^{n-2k}(-z)}{(n-2k)!}\mathrm{Li}_{2k}\\!\left(-1\right)\,=\,\sum_{k=1}^{\infty}\mathrm{Li}_{2k}\\!\left(-1\right)\varepsilon^{2k}\underbrace{\sum_{N=0}^{\infty}\frac{\left(\varepsilon\ln(-z)\right)^{N}}{N!}}_{=(-z)^{\varepsilon}}\,.$
(81)

Thus, we obtain

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)+\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{z}\right)\,$
$\displaystyle=\,1\,+\,(-z)^{\varepsilon}\left[1\,-\,2\,\sum_{k=1}^{\infty}\mathrm{Li}_{2k}\\!\left(-1\right)\varepsilon^{2k}\right].$
(82)

The remaining polylogarithm constants can be written as
(Lewin_Polylogs_Associated_Functions, , eq. (7.29))

$\displaystyle\mathrm{Li}_{2k}\\!\left(-1\right)\,=\,-(-1)^{k-1}\left(2^{2k-1}-1\right)\mathrm{B}_{2k}\frac{\pi^{2k}}{(2k)!}\,,$
(83)

where $\mathrm{B}_{2k}$ are the Bernoulli numbers, which are given by the
coefficients in the following series expansion777Note that this modern
definition of the Bernoulli numbers differs slightly from the older definition
given in ref. (Lewin_Polylogs_Associated_Functions, , eq. (7.23)), the
relation cited in eq. (83) was adjusted accordingly.

$\displaystyle\frac{t}{e^{t}-1}\,=\,\sum_{n=0}^{\infty}\mathrm{B}_{n}\frac{t^{n}}{n!}\,.$
(84)

Using this, we find that the term in square brackets is the series
representation of
$\pi\varepsilon\csc(\pi\varepsilon)=\pi\varepsilon/\sin(\pi\varepsilon)$
(Abramovitz_Stegun, , eq. (4.3.68)),

$\displaystyle
1\,-\,2\,\sum_{r=1}^{\infty}\mathrm{Li}_{2r}\\!\left(-1\right)\varepsilon^{2r}\,$
$\displaystyle=\,1\,+\,\sum_{r=1}^{\infty}\frac{(-1)^{r-1}\,2\left(2^{2r-1}-1\right)\mathrm{B}_{2r}}{(2r)!}(\pi\varepsilon)^{2r}\,=\,\frac{\pi\varepsilon}{\sin(\pi\varepsilon)}\,.$
(85)

Hence the inversion formula for the hypergeometric function is

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)+\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{z}\right)=\,1\,+\,(-z)^{\varepsilon}\,\frac{\pi\varepsilon}{\sin(\pi\varepsilon)}\,.$
(86)

###  B.3 Branch cut structure of
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$

The hypergeometric function
${{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;z\right)$ has a
branch cut along the positive real axis starting at $1$. This is explicitly
seen through the inversion formula derived in eq. (86). If we take
$x\in\mathbb{R}_{>1}$ and $\mathrm{i}0$ as an infinitesimal imaginary part,
then we have

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x\pm\mathrm{i}0\right)\,=\,-\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{x}\right)\,+\,1\,+\,(-x\mp\mathrm{i}0)^{\varepsilon}\frac{\pi\varepsilon}{\sin(\pi\varepsilon)}\,.$
(87)

Note that the imaginary part could be dropped in the hypergeometric function
on the right-hand side since the function is well-defined for arguments
$1/x<1$. Explicitly writing the factor $(-x\mp\mathrm{i}0)^{\varepsilon}$ in
terms of its real and imaginary part, we obtain

$\displaystyle{{}_{2}}\mathrm{F}_{1}\\!\left(1,-\varepsilon,1-\varepsilon;x\pm\mathrm{i}0\right)\,=\,-\,{{}_{2}}\mathrm{F}_{1}\\!\left(1,\varepsilon,1+\varepsilon;\frac{1}{x}\right)\,+\,1\,+\,\pi\varepsilon\cot(\pi\varepsilon)\,x^{\varepsilon}\,\mp\mathrm{i}\pi\varepsilon\,x^{\varepsilon}\,.$
(88)

##  Appendix C Single-valued polylogarithms

###  C.1 Single-valued polylogarithms in the literature

There are several definitions of so called single-valued polylogarithms in the
literature. Their common feature is that they do not have branch cuts opposed
to regular polylogarithms. Depending on the intended use, different versions
prove themselves to be most suitable (Zhao2016, ).

Bloch and Wigner devised an important version of a single-valued dilogarithm
(Bloch2011, )

$\displaystyle
D_{2}(z)\,=\,\mathrm{Im}\left[\mathrm{Li}_{2}\\!\left(x\right)\right]+\mathrm{arg}(1-z)\ln|z|\,,$
(89)

defined on $\mathbb{C}\setminus\\{0,1\\}$. This was generalized to
polylogarithms by Ramakrishnan (Ramakrishnan1986, ). With a slight
modification for odd $n$ due to Wojtkowiak (Wojtkowiak1989, ) and the notation
$\mathfrak{R}_{n}=\mathrm{Re}$ / $\mathrm{Im}$ if $n$ is odd / even, his
definition reads

$\displaystyle\mathrm{D}_{n}\\!\left(z\right)=\mathfrak{R}_{n}\left[\sum_{k=0}^{n-1}\frac{(-1)^{k}\ln^{k}|z|}{k!}\mathrm{Li}_{n-k}\\!\left(z\right)-\frac{(-1)^{n}\ln^{n}|z|}{2n!}\right].$
(90)

Similarly Zagier proposed (Zagier1991, , eq. (33))

$\displaystyle\mathrm{P}_{n}\\!\left(z\right)=\mathfrak{R}_{n}\left[\sum_{k=0}^{n-1}\frac{2^{k}\mathrm{B}_{k}}{k!}\ln^{k}|z|\,\mathrm{Li}_{n-k}\\!\left(z\right)-\frac{2^{n-1}\mathrm{B}_{n}}{n!}\ln^{n}|z|\right].$
(91)

Both are single-valued real-analytic functions from
$\mathbb{C}\setminus\\{0,1\\}$ to $\mathbb{R}$. Zagier’s version satisfies the
functional equations of the polylogarithm in a clean fashion without product
terms. However, both Ramakrishnan’s and Zagier’s functions are identically
zero on $\mathbb{R}$ for even $n$. Hence, they are not suitable for our
purpose.

Also in an attempt to get rid of product terms in the functional equations,
Lewin defined his version of polylogarithms by (Lewin1991, , eq. (3.19))

$\displaystyle\mathrm{L}_{n}\\!\left(x\right)=\sum_{k=0}^{n-1}\frac{(-1)^{k}}{k!}\,\ln^{k}|x|\,\mathrm{Li}_{n-k}\\!\left(x\right)+\frac{(-1)^{n-1}}{n!}\ln^{n-1}|x|\ln(1-x)$
(92)

for $x<1$, as a generalization of Rogers’ dilogarithm
$\mathrm{L}_{2}\\!\left(x\right)$ (Rogers1907, , eq. (1)). They are related to
Kummer’s functions (Kummer1840, , p. 330)
${\Lambda_{n}(x)=\int_{0}^{x}\mathrm{d}t\,\frac{\ln^{n-1}|t|}{1+t}}$ by
(Lewin1991, , eq. (3.18))

$\displaystyle\mathrm{L}_{n}\\!\left(x\right)=\frac{(-1)^{n}}{(n-1)!}\,\Lambda_{n}(-x)+\frac{(-1)^{n-1}}{n!}\ln^{n-1}|x|\ln(1-x)\,.$
(93)

Kummer’s functions are real-valued for real arguments, which makes them
single-valued polylogarithms predating Bloch and Wigner’s construction by well
over a century. We observe from eq. (93) that
$\mathrm{L}_{n}\\!\left(x\right)$ can be continued to a single-valued function
for $x>1$ by replacing $\ln(1-x)\rightarrow\ln|1-x|$ in definition (92).

###  C.2 Definition and basic properties of continuous single-valued
polylogarithms

We define for $n\in\mathbb{N}_{\geq 2}$ and real $x$ a version of single-
valued polylogarithms by

$\displaystyle\mathcal{L}_{n}\\!\left(x\right)\,\equiv\,\sum_{k=0}^{n-1}\frac{(-1)^{k}}{k!}\,\ln^{k}|x|\,\mathrm{Li}_{n-k}\\!\left(x\right)+\frac{(-1)^{n-1}}{n!}\ln^{n-1}|x|\ln|1-x|\,.$
(94)

For $x>1$ the polylogarithms are taken on their principal branch, which is
inherited from the principal branch of the logarithm, i.e.
${\mathrm{Li}_{n}\\!\left(x\right)\equiv\mathrm{Li}_{n}\\!\left(x-\mathrm{i}0\right)}$.
As a result we have a function of a real variable $x$ with the following
properties:

  * (a)

$\mathcal{L}_{n}\\!\left(x\right)$ is single-valued, i.e. there is no branch
cut and no imaginary part, in contrast to $\mathrm{Li}_{n}\\!\left(x\right)$.

  * (b)

$\mathcal{L}_{n}\\!\left(x\right)$ is continuous for all $x\in\mathbb{R}$.

  * (c)

$\mathcal{L}_{n}\\!\left(x\right)$ is bounded on $\mathbb{R}$, in contrast to
$\mathrm{Li}_{n}\\!\left(x\right)$.

  * (d)

$\mathcal{L}_{n}\\!\left(x\right)$ satisfies clean versions of the functional
equations of $\mathrm{Li}_{n}\\!\left(x\right)$, i.e. without product terms,
as for Lewin’s and Zagier’s functions.

These properties make the function $\mathcal{L}_{n}\\!\left(x\right)$ well-
suited for various applications, where one is interested in having a well-
behaved version of polylogarithms with real arguments. Plots of the first few
single-valued polylogarithms are shown in figure 4.

Figure 4: Plots of the single-valued polylogarithms defined in eq. (94), on
the left for even $n$, on the right for odd $n$.

The qualitative behaviour of the functions is different for even and odd $n$.
Here, $\zeta_{n}=\zeta(n)=\sum_{k=1}^{\infty}k^{-n}$ are the integer values of
the Riemann zeta function.

  * For even $n$: $\mathcal{L}_{n}\\!\left(x\right)$ is monotonously increasing with $\mathcal{L}_{n}\\!\left(-\infty\right)=2(2^{1-n}-1)\zeta_{n}$ and $\mathcal{L}_{n}\\!\left(\infty\right)=2\zeta_{n}$. Special values are $\mathcal{L}_{n}\\!\left(-1\right)=(2^{1-n}-1)\zeta_{n}$, $\mathcal{L}_{n}\\!\left(0\right)=0$, $\mathcal{L}_{n}\\!\left(1\right)=\zeta_{n}$. These functions are discontinuous at infinity, meaning $\mathcal{L}_{n}\\!\left(\infty\right)\neq\mathcal{L}_{n}\\!\left(-\infty\right)$. Therefore functional equations mapping infinity to $0$ or $1$ will be discontinuous at those points.

  * For odd $n$: It is $\mathcal{L}_{n}\\!\left(-\infty\right)=0$, $\mathcal{L}_{n}\\!\left(x\right)$ is decreasing in $(-\infty,-1)$ to a minimum value of $\mathcal{L}_{n}\\!\left(-1\right)=(2^{1-n}-1)\zeta_{n}$, increases in $(-1,1)$, with $\mathcal{L}_{n}\\!\left(0\right)=0$, to a maximum of $\mathcal{L}_{n}\\!\left(1\right)=\zeta_{n}$ and monotonously falls in $(1,\infty)$ to $\mathcal{L}_{n}\\!\left(\infty\right)=0$. Since ${\mathcal{L}_{n}\\!\left(\infty\right)=\mathcal{L}_{n}\\!\left(-\infty\right)}$ these functions are continuous at infinity. Therefore functional equations mapping infinity to $0$ or $1$ are continuous.

For both odd and even $n$, $\mathcal{L}_{n}\\!\left(x\right)$ is real analytic
on $\mathbb{R}\backslash\\{0,1\\}$. At $x=0$ its derivative diverges
logarithmically, at $x=1$ it is $(n-2)$-times differentiable, the $(n-1)$th
derivative diverges logarithmically at this point.

###  C.3 Functional equations of continuous single-valued polylogarithms

The functional equations satisfied by single-valued polylogarithms are
inherited from the corresponding polylogarithms
(Lewin_Polylogs_Associated_Functions, ). Single-valued polylogarithms obey
these equations in a cleaner fashion, where all product terms vanish. There
are equations satisfied for all $n$ and additional identities for small values
of $n$. We explicitly give the identities for di- and trilogarithm, identities
for higher order polylogarithms can be found in
(Lewin_Polylogs_Associated_Functions, ; Lewin1991, ).

####  C.3.1 Equations of all single-valued polylogarithms

For all $n$ there is an inversion relation and a duplication formula.

  * (i)

Inversion relation:

$\displaystyle\mathcal{L}_{n}\\!\left(x\right)+(-1)^{n}\mathcal{L}_{n}\\!\left(\frac{1}{x}\right)=\left\\{\begin{array}[]{cl}0&\text{
for odd }n\,,\\\ 2\zeta_{n}&\text{ for even }n\text{ and }x>0\,,\\\
2(2^{1-n}-1)\zeta_{n}&\text{ for even }n\text{ and }x<0\,.\end{array}\right.$
(98)

Alternatively we can write this in the suggestive form

$\displaystyle\mathcal{L}_{n}\\!\left(x\right)+(-1)^{n}\mathcal{L}_{n}\\!\left(\frac{1}{x}\right)=\mathcal{L}_{n}\\!\left(\mathrm{sgn}\\!\left(x\right)\infty\right).$
(99)

  * (ii)

Duplication formula:

$\displaystyle\mathcal{L}_{n}\\!\left(x\right)+\mathcal{L}_{n}\\!\left(-x\right)-2^{1-n}\mathcal{L}_{n}\\!\left(x^{2}\right)=0\,.$
(100)

####  C.3.2 Equations of the single-valued dilogarithm

For the single-valued dilogarithm, defined as

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)=\mathrm{Li}_{2}\\!\left(x\right)+\ln|x|\ln(1-x)-\frac{1}{2}\ln|x|\ln|1-x|\,,$
(101)

there is the richest structure. We have identities connecting the arguments
$x$, $1-x$, $\frac{x}{x-1}$, $\frac{1}{x}$, $\frac{1}{1-x}$ and
$1-\frac{1}{x}$. Additionally there is the famous 2-variable 5-term relation
discovered independently by Spence (Spence1809, , p.9) and Abel (Abel1881, ,
XIV.(9)). All discontinuities in these identities are induced by the mapping
of the discontinuity at infinity.

  * (i)

Euler identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)+\mathcal{L}_{2}\\!\left(1-x\right)=\zeta_{2}\,.$
(102)

  * (ii)

Landen identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)+\mathcal{L}_{2}\\!\left(\frac{x}{x-1}\right)=\left\\{\begin{array}[]{cl}0&\text{for
}x<1\,,\\\ 3\zeta_{2}&\text{for }x>1\,.\end{array}\right.$ (105)

  * (iii)

Inversion identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)+\mathcal{L}_{2}\\!\left(\frac{1}{x}\right)=\left\\{\begin{array}[]{cl}-\zeta_{2}&\text{for
}x<0\,,\\\ 2\zeta_{2}&\text{for }x>0\,.\end{array}\right.$ (108)

  * (iv)

Inverted Euler identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)-\mathcal{L}_{2}\\!\left(\frac{1}{1-x}\right)=\left\\{\begin{array}[]{cl}-\zeta_{2}&\text{for
}x<1\,,\\\ 2\zeta_{2}&\text{for }x>1\,.\end{array}\right.$ (111)

  * (v)

Inverted Landen identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)-\mathcal{L}_{2}\\!\left(1-\frac{1}{x}\right)=\left\\{\begin{array}[]{cl}-2\zeta_{2}&\text{for
}x<0\,,\\\ \zeta_{2}&\text{for }x>0\,.\end{array}\right.$ (114)

  * (vi)

Abel identity:

$\displaystyle\mathcal{L}_{2}\\!\left(x\right)+\mathcal{L}_{2}\\!\left(y\right)+\mathcal{L}_{2}\\!\left(\frac{1-x}{1-xy}\right)+\mathcal{L}_{2}\\!\left(1-xy\right)+\mathcal{L}_{2}\\!\left(\frac{1-y}{1-xy}\right)=\pm
3\zeta_{2}\,,$ (115)

where the sign on the right is negative if $1-xy<0$ and $x,y<0$, else
positive.

####  C.3.3 Equations of the single-valued trilogarithm

Besides the inversion relation, the single-valued trilogarithm, defined as

$\displaystyle\mathcal{L}_{3}\\!\left(x\right)=\mathrm{Li}_{3}\\!\left(x\right)-\ln|x|\,\mathrm{Li}_{2}\\!\left(x\right)-\frac{1}{2}\ln^{2}|x|\ln(1-x)+\frac{1}{6}\ln^{2}|x|\ln|1-x|\,,$
(116)

satisfies a single-variable identity connecting three terms and a three-
variable 22 term relation discovered by Goncharov (Goncharov1995, ). From the
latter there follow a plethora of two variable identities as special cases.
Since the single-valued trilogarithm is continuous at infinity, these
identities are continuous.

  * (i)

Inversion relation:

$\displaystyle\mathcal{L}_{3}\\!\left(x\right)-\mathcal{L}_{3}\\!\left(\frac{1}{x}\right)=0\,.$
(117)

  * (ii)

Landen identity:

$\displaystyle\mathcal{L}_{3}\\!\left(x\right)+\mathcal{L}_{3}\\!\left(1-x\right)+\mathcal{L}_{3}\\!\left(\frac{x}{x-1}\right)=\zeta_{3}\,.$
(118)

  * (iii)

Goncharov identity:

$\displaystyle\mathcal{L}_{3}\\!\left(xyz\right)+$
$\displaystyle\bigoplus_{\mathrm{Cyc}(x,y,z)}\left[\mathcal{L}_{3}\\!\left(x\right)-\mathcal{L}_{3}\\!\left(xy\right)+\mathcal{L}_{3}\\!\left(\frac{x(1-y)}{x-1}\right)+\mathcal{L}_{3}\\!\left(\frac{y(1-x)}{y-1}\right)\right.$
$\displaystyle+$
$\displaystyle\left.\mathcal{L}_{3}\\!\left(\frac{1-x}{1-xyz}\right)+\mathcal{L}_{3}\\!\left(\frac{xy(1-z)}{1-xyz}\right)-\mathcal{L}_{3}\\!\left(\frac{x(1-y)(1-z)}{(x-1)(1-xyz)}\right)\right]=3\zeta_{3}\,,$
(119)

where cyclic summation for a function $f(x,y,z)$ is defined as

$\displaystyle\bigoplus_{\mathrm{Cyc}(x,y,z)}f(x,y,z)=f(x,y,z)+f(y,z,x)+f(z,x,y)\,.$
(120)

##  Appendix D Example: Regularizing kinematic limits for Drell-Yan
differential in rapidity

In this appendix, we illustrate the regularization of kinematic divergences
occurring in the box integral using the Drell-Yan process differential in
momentum transfer $Q^{2}$ and rapidity $y$ at
$\mathcal{O}(\alpha\alpha_{s}^{2})$ as an example. This process was first
calculated to NNLO in (Anastasiou2003, ) using the reversed unitarity method.
Instead, for the purpose of our example, we consider the NNLO calculation
using the traditional treatment of the phase space integral as was done for
the NLO calculation in (Altarelli1979, ). We will explicitly see that orders
beyond $\varepsilon^{0}$ of $D_{0}$ are required and the full
$\varepsilon$-dependence of any terms divergent in the kinematic limit is
needed here. To illustrate these points it is sufficient to look at the
contribution to $\frac{\mathrm{d}\sigma}{\mathrm{d}Q^{2}\mathrm{d}y}$ from the
quark-antiquark process. At NNLO the box integral appears through the
interference term of one-loop and tree level diagrams for
$q\bar{q}\rightarrow\gamma^{\ast}g$, which are depicted in figure 5.

Figure 5: NNLO contribution to $q\bar{q}\rightarrow\gamma^{\ast}$ from
interference of box and tree-level amplitudes. There are three different boxes
and two different tree-level diagrams for this process.

The two particle phase space in $d=4-2\varepsilon$ dimensions is given by (see
eqs. (79) and (80) in (Altarelli1979, ))

$\displaystyle\mathrm{d}\mathrm{PS}_{2}=\frac{1}{8\pi}\left(\frac{4\pi}{Q^{2}}\right)^{\varepsilon}\frac{z^{\varepsilon}(1-z)^{1-2\varepsilon}}{\Gamma(1-\varepsilon)}\int_{0}^{1}\mathrm{d}y\,y^{-\varepsilon}(1-y)^{-\varepsilon}\,,$
(121)

with

$\displaystyle
s=\frac{Q^{2}}{z}\,,\quad\,t=\frac{-Q^{2}}{z}(1-z)(1-y)\,,\quad\,u=-\frac{Q^{2}}{z}(1-z)y\,.$
(122)

For simplicity, we only perform the regularization for $y\rightarrow 1$ in the
following. Treatment of $y\rightarrow 0$ and $z\rightarrow 1$ would be
analogous. The relevant part of the phase space in the limit $y\rightarrow 1$
is $(1-y)^{-\varepsilon}$. Calculating the interference terms of the type
depicted in figure 5, we find that the scalar box integral contributes through

  * (i)

$(1-y)^{-\varepsilon}D_{0}\\!\left(\frac{Q^{2}}{z},-\frac{Q^{2}}{z}(1-z)(1-y),Q^{2}\right)\,f_{1}(y,\varepsilon)$
,

  * (ii)

$(1-y)^{-\varepsilon}D_{0}\\!\left(\frac{Q^{2}}{z},-\frac{Q^{2}}{z}(1-z)y,Q^{2}\right)\,f_{2}(y,\varepsilon)$
,

  * (iii)

$(1-y)^{-1-\varepsilon}D_{0}\\!\left(\frac{Q^{2}}{z},-\frac{Q^{2}}{z}(1-z)y,Q^{2}\right)\,f_{3}(y,\varepsilon)$
,

  * (iv)

$(1-y)^{-\varepsilon}D_{0}\\!\left(-\frac{Q^{2}}{z}(1-z)y,-\frac{Q^{2}}{z}(1-z)(1-y),Q^{2}\right)\,f_{4}(y,\varepsilon)$
,

where the $f_{i}$ are functions finite in the limit $y\rightarrow 1$ order by
order in $\varepsilon$. The behavior of the scalar box functions $D_{0}$ for
$y\rightarrow 1$ can be obtained from eq. (62). We find

$\displaystyle D_{0}(s,t,Q^{2})$
$\displaystyle\sim\frac{1}{st}\left|\frac{u\mu^{2}}{st}\right|^{\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;-\frac{u}{s}\right)+\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{u}{t}\right)-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{uQ^{2}}{st}\right)\right]$
$\displaystyle\sim(1-y)^{-1-\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;(1-z)y\right)+\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{y}{1-y}\right)-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{yz}{1-y}\right)\right],$
(123)

where $\sim$ means that we dropped factors non-singular for $y\rightarrow 1$.
In the limit $y\rightarrow 1$ the argument of the first $\mathfrak{F}$ stays
finite, for the latter two the argument tends to $-\infty$. According to eq.
(53), $\mathfrak{F}$ is finite order by order in $\varepsilon$ in this case.
Hence, we can write

$\displaystyle D_{0}(s,t,Q^{2})=(1-y)^{-1-\varepsilon}A_{st}(y)\,,$ (124)

where $A_{st}(y)$ collects all factors which stay finite in the limit
$y\rightarrow 1$. Similarly,

$\displaystyle D_{0}(t,u,Q^{2})$
$\displaystyle\sim\frac{1}{tu}\left|\frac{s\mu^{2}}{tu}\right|^{\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;-\frac{s}{t}\right)+\mathfrak{F}\\!\left(\varepsilon;-\frac{s}{u}\right)-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{sQ^{2}}{tu}\right)\right]$
$\displaystyle\sim(1-y)^{-1-\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;\frac{1}{(1-z)(1-y)}\right)+\,\mathfrak{F}\\!\left(\varepsilon;-\frac{1}{(1-z)y}\right)\right.$
$\displaystyle\hphantom{\sim(1-y)^{-1-\varepsilon}}\,\,\left.-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{1}{(1-z)^{2}y(1-y)}\right)\right],$
(125)

where again the bracket is finite for $y\rightarrow 1$, since $\mathfrak{F}$
stays finite if its argument tends to $\pm\infty$ (compare eqs. (52) and
(53)). Hence, analogous to eq. (124)

$\displaystyle D_{0}(t,u,Q^{2})=(1-y)^{-1-\varepsilon}A_{tu}(y)\,.$ (126)

For $D(s,u,Q^{2})$ we obtain

$\displaystyle D(s,u,Q^{2})$
$\displaystyle\sim\frac{1}{su}\left|\frac{t\mu^{2}}{su}\right|^{\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;-\frac{t}{s}\right)+\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{t}{u}\right)-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{tQ^{2}}{su}\right)\right]$
$\displaystyle\sim(1-y)^{\varepsilon}\left[\mathfrak{F}\\!\left(\varepsilon;(1-z)(1-y)\right)+\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{1-y}{y}\right)\right.$
$\displaystyle\hphantom{\sim(1-y)^{-1-\varepsilon}}\,\,\left.-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\mathfrak{F}\\!\left(\varepsilon;-\frac{(1-y)z}{y}\right)\right].$
(127)

In the limit $y\rightarrow 1$ the arguments of all three $\mathfrak{F}$
functions become zero. From eq. (54) we infer that for
$\tilde{\mathfrak{F}}(\varepsilon;x)=\mathfrak{F}(\varepsilon;x)-|x|^{\varepsilon}$
the limit $\lim_{x\rightarrow 0}\frac{\tilde{\mathfrak{F}}(\varepsilon;x)}{x}$
is finite order by order in $\varepsilon$. Therefore, we can write

$\displaystyle
D(s,u,Q^{2})\sim(1-z)^{-\varepsilon}+\mathrm{e}^{\mathrm{i}\pi\varepsilon}y^{\varepsilon}(1-z^{-\varepsilon})$
$\displaystyle\quad+(1-y)^{\varepsilon}\left[\tilde{\mathfrak{F}}\\!\left(\varepsilon;(1-z)(1-y)\right)+\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\tilde{\mathfrak{F}}\\!\left(\varepsilon;-\frac{1-y}{y}\right)-\mathrm{e}^{\mathrm{i}\pi\varepsilon}\,\tilde{\mathfrak{F}}\\!\left(\varepsilon;-\frac{(1-y)z}{y}\right)\right].$
(128)

Note that the term in square brackets stays finite in the limit $y\rightarrow
1$ even when divided by $(1-y)$ and hence

$\displaystyle D(s,u,Q^{2})$
$\displaystyle=A_{su}(y)+(1-y)^{1+\varepsilon}B_{su}(y)\,,$ (129)

with both $A_{su}(y)$ and $B_{su}(y)$ finite in the limit $y\rightarrow 1$
order by order in $\varepsilon$.

To convert the singular behavior for $y\rightarrow 1$ into poles in
$\varepsilon$, which is required to combine real and virtual corrections
ensuring cancellation of infrared poles, we use the expansion (see e.g. eq.
(110) in (Altarelli1979, ) or eq. (69) in (Campbell2017, ))

$\displaystyle(1-y)^{-1-k\varepsilon}=-\frac{1}{k\varepsilon}\,\delta(1-y)+\sum_{n=0}^{\infty}\frac{(-k\varepsilon)^{n}}{n!}\left[\frac{\ln^{n}(1-y)}{1-y}\right]_{+}\,,$
(130)

where the plus distribution is defined by

$\displaystyle\int_{0}^{1}\mathrm{d}y\,f(y)\,[g(y)]_{+}\,\equiv\int_{0}^{1}\mathrm{d}y\,[f(y)-f(1)]\,g(y)\,.$
(131)

Using eqs. (124), (126), and (129), we find for the four types of
contributions

  * (i)

$(1-y)^{-\varepsilon}D_{0}(s,t,Q^{2})=\left(-\frac{1}{2\varepsilon}\,\delta(1-y)+\sum_{n=0}^{\infty}\frac{(-2\varepsilon)^{n}}{n!}\left[\frac{\ln^{n}(1-y)}{1-y}\right]_{+}\right)A_{st}(y)\,f_{1}(y)$
,

  * (ii)

$(1-y)^{-\varepsilon}D_{0}(s,u,Q^{2})=\left(\sum_{n=0}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\ln^{n}(1-y)\right)(A_{su}(y)+(1-y)B_{su}(y))\,f_{2}(y)$
,

  * (iii)

$(1-y)^{-1-\varepsilon}D_{0}(s,t,Q^{2})=\left(-\frac{1}{\varepsilon}\,\delta(1-y)+\sum_{n=0}^{\infty}\frac{(-\varepsilon)^{n}}{n!}\left[\frac{\ln^{n}(1-y)}{1-y}\right]_{+}\right)A_{su}(y)\,f_{3}(y)$
$(1-y)^{-1-\varepsilon}D_{0}(s,t,Q^{2})=$$\frac{1}{1}$$+B_{su}(y)\,f_{3}(y)$ ,

  * (iv)

$(1-y)^{-\varepsilon}D_{0}(t,u,Q^{2})=\left(-\frac{1}{2\varepsilon}\,\delta(1-y)+\sum_{n=0}^{\infty}\frac{(-2\varepsilon)^{n}}{n!}\left[\frac{\ln^{n}(1-y)}{1-y}\right]_{+}\right)A_{tu}(y)\,f_{4}(y)$
.

We see that the order $\varepsilon^{1}$ of the box integral contributes to the
order $\varepsilon^{0}$ of the full result, since the
$\frac{1}{\varepsilon}$-poles generated by eq. (130) multiply the order
$\varepsilon^{1}$ of the functions $A_{st}(y)$, $A_{su}(y)$, and $A_{tu}(y)$.
Hence, it is required to perform the $\varepsilon$-expansion of the box
integral to one order higher than needed if no regularization would be
required. In the case of multiple regularizations, e.g. for $y,z\rightarrow 1$
as in eq. (111) of (Altarelli1979, ), multiple poles multiply the box
integral, requiring knowledge of the box integral to even higher orders in
$\varepsilon$. Furthermore it was crucial to sum up any terms of the box
integral singular in the kinematic limit to all orders, i.e. summing up terms
of the form $\varepsilon^{n}\ln^{n}(1-y)$ to $(1-y)^{-\varepsilon}$. This
changed $\frac{1}{\varepsilon}$ to $\frac{1}{2\varepsilon}$ in (i) and (iv).

## References

  * (1) A. Accardi et al., _Electron Ion Collider: The Next QCD Frontier: Understanding the glue that binds us all_ , _Eur. Phys. J. A_ 52 (2016) 268 [1212.1701].
  * (2) D. Anderle, D. de Florian, and Y. Rotstein Habarnau, _Towards semi-inclusive deep inelastic scattering at next-to-next-to-leading order_ , _Phys. Rev. D_ 95 (2017) 034027 [1612.01293].
  * (3) M. Abele, D. de Florian, and W. Vogelsang, _Approximate NNLO QCD corrections to semi-inclusive DIS_ , _Phys. Rev. D_ 104 (2021) 094046 [2109.00847].
  * (4) I. Borsa, D. de Florian, R. Sassot, M. Stratmann, and W. Vogelsang, _Towards a Global QCD Analysis of Fragmentation Functions at Next-to-Next-to-Leading Order Accuracy_ , _Phys. Rev. Lett._ 129 (2022) 012002 [2202.05060].
  * (5) M. Abele, D. de Florian, and W. Vogelsang, _Threshold resummation at $\text{N}^{3}$LL accuracy and approximate $\text{N}^{3}$LO corrections to semi-inclusive DIS_, _Phys. Rev. D_ 106 (2022) 014015 [2203.07928].
  * (6) G. ’t Hooft and M. Veltman, _Scalar one-loop integrals_ , _Nucl. Phys. B_ 153 (1979) 365.
  * (7) W. L. Van Neerven, and J. A. M. Vermaseren, _Large loop integrals_ , _Phys. Lett. B_ 137 (1984) 241.
  * (8) G. J. Van Oldenborgh and J. A. M. Vermaseren, _New algorithms for one-loop integrals_ , _Z. Phys. C_ 46 (1990) 425.
  * (9) W. Beenakker and A. Denner, _Infrared divergent scalar box integrals with applications in the electroweak standard model_ , _Nucl. Phys. B_ 338 (1990) 349.
  * (10) A. Denner and U. Nierste, and R. Scharf, _A compact expression for the scalar one-loop four-point function_ , _Nucl. Phys. B_ 367 (1991) 637.
  * (11) G. J. Van Oldenborgh, _FF – a package to evaluate one-loop Feynman diagrams_ , _Comp. Phys. Comm._ 66 (1991) 1.
  * (12) T. Hahn and M. Perez-Victoria, _Automated one-loop calculations in four and D dimensions_ , _Comp. Phys. Comm._ 118 (1999) 153 [arXiv:hep-ph/9807565].
  * (13) J. Fleischer, F. Jegerlehner, and O. V. Tarasov, _A new hypergeometric representation of one-loop scalar integrals in d dimensions_ , _Nucl. Phys. B_ 672 (2003) 303 [arXiv:hep-ph/0307113].
  * (14) H. H. Patel, _Package-X: A Mathematica package for the analytic calculation of one-loop integrals_ , _Comp. Phys. Comm._ 197 (2015) 276 [arXiv:1503.01469].
  * (15) S. Carrazza, R. K. Ellis, and G. Zanderighi, _QCDLoop: a comprehensive framework for one-loop scalar integrals_ , _Comp. Phys. Comm._ 209 (2016) 134 [arXiv:1605.03181].
  * (16) V. E. Lyubovitskij, F. Wunder and A. S. Zhevlakov, _New ideas for handling of loop and angular integrals in D-dimensions in QCD_ , _JHEP_ 06 (2021) 066 [2102.08943].
  * (17) G. Passarino and M. Veltman, _One-loop Corrections for $e^{+}e^{-}$ Annihilation into $\mu^{+}\mu^{-}$ in the Weinberg Model_, _Nucl. Phys. B_ 160 (1979) 151.
  * (18) C. G. Bollini and J. J. Giambiagi, _Dimensional Renormalization: The Number of Dimensions as a Regularizing Parameter_ , _Nuovo Cim. B_ 12 (1972) 20.
  * (19) G. ’t Hooft, _Dimensional regularization and the renormalization group_ , _Nucl. Phys. B_ 61 (1973) 455.
  * (20) R. Gastmans and R. Meuldermans, _Dimensional regularization of the infrared problem_ , _Nucl. Phys. B_ 63 (1973) 277.
  * (21) G. Leibbrandt, _Introduction to the Technique of Dimensional Regularization_ , _Rev. Mod. Phys._ 47 (1975) 849.
  * (22) K. Fabricius and I. Schmitt, _Calculation of Dimensionally Regularized Box Graphs in the Zero Mass Case_ , _Z. Phys. C_ 3 (1979) 51.
  * (23) T. Matsuura, S. C. van der Marck, and W. L. Van Neerven, _The calculation of the second order soft and virtual contributions to the Drell-Yan cross section_ , _Nucl. Phys. B_ 319 (1989) 570.
  * (24) Z. Bern, L. Dixon, and D. A. Kosower, _Dimensionally-regulated pentagon integrals_ , _Nucl. Phys. B_ 412 (1994) 751.
  * (25) E. B. Zijlstra and W. L. Van Neerven, _Order- $\alpha_{s}^{2}$ QCD corrections to the deep inelastic proton structure functions $F_{2}$ and $F_{L}$_, _Nucl. Phys. B_ 383 (1992) 525.
  * (26) G. Duplančić and B. Nižić, _Dimensionally-regulated one-loop scalar integrals with massless internal lines_ , _Eur. Phys. J. C_ 20 (2001) 357.
  * (27) J. P. Ellis, _TikZ-Feynman: Feynman diagrams with TikZ_ , _Comp. Phys. Comm._ 210 (2017) 103 [arXiv:1601.05437].
  * (28) R. K. Ellis and G. Zanderighi, _Scalar one-loop integrals for QCD_ , _JHEP_ 02 (2008) 002 [0712.1851].
  * (29) M. E. Peskin and D. V. Schroeder, _An Introduction to Quantum Field Theory_ , 1995.
  * (30) V. A. Smirnov, _Analytic Tools for Feynman Integrals_ , Springer (2012).
  * (31) M. Abramowitz and I. A. Stegun, _Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables_ , National Bureau of Standards (1972).
  * (32) L. Lewin, _Polylogarithms and Associated Functions_ , Elsevier (1981).
  * (33) J. Haug and F. Wunder, _The massless non-adjacent double off-shell scalar box integral – branch cut structure and all-order epsilon expansion_ , _JHEP_ 05 (2023) 059 [2302.01956].
  * (34) J. Zhao, _Multiple zeta functions, multiple polylogarithms and their special values_ , vol. 12, World Scientific (2016).
  * (35) S. J. Bloch, _Higher regulators, algebraic $K$-theory, and zeta functions of elliptic curves_, vol. 11, American Mathematical Soc. (2011).
  * (36) D. Ramakrishnan, _Analogs of the Bloch-Wigner function for higher polylogarithms_ , _Contemp. Math_ 55 (1986) 371.
  * (37) Z. Wojtkowiak, _A construction of analogs of the Bloch-Wigner function_ , _Mathematica Scandinavica_ 65 (1989) 140.
  * (38) D. Zagier, _Polylogarithms, Dedekind zeta functions, and the algebraic K-theory of fields_ , in _Arithmetic Algebraic Geometry_ , G. van der Geer, F. Oort, J. Steenbrink, ed., no. 89 in Progress in Mathematics, p. 391, Springer (1991), DOI.
  * (39) L. Lewin, _Structural properties of polylogarithms_ , no. 37 in Mathematical Surveys and Monographs, American Mathematical Soc. (1991).
  * (40) L. J. Rogers, _On Function Sum Theorems Connected with the Series $\sum_{n=1}^{\infty}\frac{x^{n}}{n^{2}}$_, _Proceedings of the London Mathematical Society_ 2 (1907) 169.
  * (41) E. E. Kummer, _Ueber die Transcendenten, welche aus wiederholten Integrationen rationaler Formeln entstehen._ , _Journal für die reine und angewandte Mathematik_ Band 21 (1840) 74.
  * (42) W. Spence, _An Essay on the Theory of the Various Orders of Logarithmic Transcendents_ , John Murray and Archibald Constable and Company (1809).
  * (43) N. H. Abel, _Œuvres complètes_ , vol. 2, Grøndahl (1881).
  * (44) A. B. Goncharov, _Geometry of configurations, polylogarithms and motivic cohomology_ , _Advances in Mathematics_ 114.2 (1995) 197.
  * (45) C. Anastasiou, L. Dixon, and K. Melnikov, F. Petriello, _Dilepton Rapidity Distribution in the Drell-Yan Process at Next-to-Next-to-Leading Order in QCD_ , _Phys. Rev. Lett._ 91 (2003) 182002 [0306192].
  * (46) G. Altarelli, R. K. Ellis, and G. Martinelli, _Large perturbative corrections to the Drell-Yan process in QCD_ , _Nuclear Physics B_ 157 (1979) 461.
  * (47) J. M. Campbell, R. K. Ellis, R. Mondini, and C. Williams, _The NNLO QCD soft function for 1-jettiness_ , _Eur. Phys. J. C_ 78 (2018) 234 [1711.09984].