Multiple polylogarithms with algebraic arguments and the two-loop EW-QCD Drell-Yan master integrals

TL;DR

This paper demonstrates that certain complex Feynman integrals with algebraic singularities can be expressed using multiple polylogarithms, enabling efficient numerical evaluation of two-loop corrections in Drell-Yan processes.

Contribution

It introduces a method to evaluate integrals with unrationalizable roots using multiple polylogarithms, applied to two-loop Drell-Yan master integrals.

Findings

Successfully evaluated two-loop master integrals for Drell-Yan production.

Optimized functional basis for fast numerical evaluation.

Extended the class of integrals expressible with multiple polylogarithms.

Abstract

We consider Feynman integrals with algebraic leading singularities and total differentials in $\epsilon\,\mathrm{d}\ln$ form. We show for the first time that it is possible to evaluate integrals with singularities involving unrationalizable roots in terms of conventional multiple polylogarithms, by either parametric integration or matching the symbol. As our main application, we evaluate the two-loop master integrals relevant to the $\alpha \alpha_s$ corrections to Drell-Yan lepton pair production at hadron colliders. We optimize our functional basis to allow for fast and stable numerical evaluations in the physical region of phase space.

Tables

Table 1. Table 1: Integral families for the master integrals which appear in the differential equations for the two-loop mixed EW-QCD Drell-Yan master integrals with two massive internal lines.

Family $A$	Family $C$	Family $E$
$k_1^2$	$k_1^2$	$k_1^2$
$k_2^2$	$k_2^2$	$k_2^2-m^2$
${\left(k_1-k_2\right)}^2$	${\left(k_1-k_2\right)}^2$	${\left(k_1-k_2\right)}^2$
${\left(k_1-p_1\right)}^2$	${\left(k_1-p_1\right)}^2$	${\left(k_1-p_1\right)}^2$
${\left(k_2-p_1\right)}^2$	${\left(k_2-p_1\right)}^2$	${\left(k_2-p_1\right)}^2$
${\left(k_1-p_1-p_2\right)}^2$	${\left(k_1-p_1-p_2\right)}^2$	${\left(k_1-p_1-p_2\right)}^2$
${\left(k_2-p_1-p_2\right)}^2$	${\left(k_2-p_1-p_2\right)}^2-m^2$	${\left(k_2-p_1-p_2\right)}^2-m^2$
${\left(k_1-p_3\right)}^2$	${\left(k_1-p_3\right)}^2$	${\left(k_1-p_3\right)}^2$
${\left(k_2-p_3\right)}^2$	${\left(k_2-p_3\right)}^2$	${\left(k_2-p_3\right)}^2$

Equations

$$\displaystyle\includegraphics[valign={m},raise={.1 cm},height=28.45274pt,keepaspectratio]{BalgLSunmarked}\big{(}s,t,m^{2}\big{)}\qquad\qquad\includegraphics[valign={m},raise={.05 cm},height=28.45274pt,keepaspectratio]{int32}\big{(}s,t,m^{2}\big{)}$$

$$\includegraphics[valign={m},raise={.4 cm},height=28.45274pt,keepaspectratio]{BalgLS}\left(s,t,m^{2}\right)=-e^{2\gamma_{E}\epsilon}\Gamma(1+2\epsilon)\Bigg{[}\prod_{i=1}^{5}\int_{0}^{\infty}{\rm d}\alpha_{i}\Bigg{]}\delta(1-\alpha_{5})\,\mathcal{U}^{-1+3\epsilon}\mathcal{F}^{-1-2\epsilon},$$

$$\mathcal{U}$$

$$\mathcal{F}$$

$$-m^2\left({\alpha }_1{\alpha }_2{\alpha }_4+{\alpha }_1{\alpha }_2{\alpha }_5+{\alpha }_1{\alpha }_3{\alpha }_4+{\alpha }_1{\alpha }_3{\alpha }_5+{\alpha }_1{\alpha }_4{\alpha }_5+{\alpha }_2{\alpha }_3{\alpha }_4+{\alpha }_2{\alpha }_3{\alpha }_5+{\alpha }_3{\alpha }_4{\alpha }_5\right)$$

$$\displaystyle{\rm L}_{\alpha_{1}}=\Big{\{}\vphantom{\left(1+\alpha_{2}+\alpha_{3}\right)^{2}}1+\alpha_{2}+\alpha_{3},m^{2}-t\,\alpha_{2}\,\alpha_{3},\alpha_{2}\left(1+\alpha_{3}\right)+\left(1+\alpha_{2}+\alpha_{3}\right)\alpha_{4},$$

$$m^2{\alpha }_2+\left(m^2+\left(2m^2-s\right){\alpha }_2\right){\alpha }_4+m^2\left(1+{\alpha }_2+{\alpha }_3\right){\alpha }_4^2,$$

$$\displaystyle\alpha_{2}^{2}\left(m^{2}+t\,\alpha_{3}\left(1+\alpha_{3}\right)\right)+\left(2\,m^{2}-s+t\,\alpha_{3}\right)\left(1+\alpha_{2}+\alpha_{3}\right)\alpha_{2}\,\alpha_{4}+m^{2}\left(1+\alpha_{2}+\alpha_{3}\right)^{2}\alpha_{4}^{2}\Big{\}}.$$

$${\alpha }_4=\frac{x_4{\alpha }_2}{1+{\alpha }_2+{\alpha }_3}$$

$${\alpha }_2^2\left(m^2+t{\alpha }_3\left(1+{\alpha }_3\right)\right)+\left(2m^2-s+t{\alpha }_3\right)\left(1+{\alpha }_2+{\alpha }_3\right){\alpha }_2{\alpha }_4+m^2{\left(1+{\alpha }_2+{\alpha }_3\right)}^2{\alpha }_4^2$$

$$={\alpha }_2^2\left(m^2+t{\alpha }_3\left(1+{\alpha }_3\right)+\left(2m^2-s+t{\alpha }_3\right)x_4+m^2{x}_4^2\right).$$

$${\alpha }_2\left(1+{\alpha }_3\right)+\left(1+{\alpha }_2+{\alpha }_3\right){\alpha }_4$$

$$m^2{\alpha }_2+\left(m^2+\left(2m^2-s\right){\alpha }_2\right){\alpha }_4+m^2\left(1+{\alpha }_2+{\alpha }_3\right){\alpha }_4^2$$

$$\displaystyle\frac{\alpha_{2}}{1+\alpha_{2}+\alpha_{3}}\Big{(}m^{2}\left(1+\alpha_{2}+\alpha_{3}\right)+\big{(}m^{2}+\big{(}2\,m^{2}-$$

$$\displaystyle{\rm L}_{\alpha_{1}\alpha_{2}}=\Big{\{}1+\alpha_{3},1+\alpha_{3}+x_{4},m^{2}+t\,\alpha_{3}+t\,\alpha_{3}^{2},m^{2}+\left(2\,m^{2}-s\right)x_{4}+m^{2}x_{4}^{2},$$

$$\displaystyle m^{2}\left(1+x_{4}\right)-s+\left(m^{2}(2+x_{4})-s\right)\alpha_{3},m^{2}+\left(2\,m^{2}-s\right)x_{4}+m^{2}x_{4}^{2}+t(1+x_{4})\alpha_{3}+t\,\alpha_{3}^{2}\Big{\}}$$

$${\alpha }_3=\frac{x_3}{1-x_3}.$$

$$s=-\frac{4m^2(v_2-v_1{)}^2}{(1-v_1)(1+v_1)(1-v_2)(1+v_2)}\mathrm{and}t=-\frac{m^2(v_2-v_1{)}^2}{v_1{v}_2}$$

$$m^2+t{\alpha }_3+t{\alpha }_3^2=\frac{m^2\left(v_1-(v_2-v_1){\alpha }_3\right)\left(v_2+(v_2-v_1){\alpha }_3\right)}{v_1{v}_2},$$

$$m^2+\left(2m^2-s\right)x_4+m^2{x}_4^2=$$

$$\frac{m^2\left((1+v_1)(1-v_2)+(1-v_1)(1+v_2)x_4\right)\left((1-v_1)(1+v_2)+(1+v_1)(1-v_2)x_4\right)}{(1-v_1)(1+v_1)(1-v_2)(1+v_2)}.$$

$$\sqrt{1+\frac{2(4m^2-2s-t)}{4m^2-t}{x}_4+x_4^2}.$$

$$1+\frac{2(4m^2-2s-t)}{4m^2-t}{x}_4+x_4^2={\rho }^2.$$

$$\left(x_4^{(0)},{\rho }^{(0)}\right)=\left(-\frac{2(4m^2-2s-t)}{4m^2-t},1\right)$$

$$\rho =y_4\left(x_4-x_4^{(0)}\right)+{\rho }^{(0)}.$$

$$x_4=\frac{2y_4\left(1+\frac{4m^2-2s-t}{4m^2-t}{y}_4\right)}{(1-y_4)(1+y_4)}.$$

$$\left(v_1,v_2,m^2\right)=\left(1/7,1/5,2\right).$$

$$\includegraphics[valign={m},raise={.2 cm},height=22.76228pt,keepaspectratio]{BalgLSunmarked}\left(-1/36,-8/35,2\right)\approx-6.317550089475753330169497\ldots\,+\,\mathcal{O}\left(\epsilon\right),$$

$$\sqrt{p_1(v_1,v_2)p_2(v_1,v_2)p_3(v_1,v_2)p_4(v_1,v_2)},$$

$$p_1(v_1,v_2)$$

$$p_2(v_1,v_2)$$

$$p_3(v_1,v_2)$$

$$p_4(v_1,v_2)$$

$$q(p_1)\bar{q}(p_2)$$

$$s=(p_1+p_2{)}^2,t=(p_1-p_3{)}^2,u=(p_2-p_3{)}^2,p_1^2=p_2^2=p_3^2=p_4^2=0.$$

$$\displaystyle\mathbf{f}_{1}^{\mathrm{\accentset{\rule{2.79996pt}{0.8pt}}{A}:38}}=\includegraphics[valign={m},raise={.3 cm},height=54.2025pt,width=54.2025pt,keepaspectratio]{int1}(t)\quad\mathbf{f}_{2}^{\mathrm{C:97}}=\includegraphics[valign={m},raise={.3 cm},height=54.2025pt,width=54.2025pt,keepaspectratio]{int2}(s)\quad\mathbf{f}_{3}^{\mathrm{C:76}}=\includegraphics[valign={m},raise={.3 cm},height=21.68231pt,width=21.68231pt,keepaspectratio]{int3}$$

Full text

††preprint: MITP/19-043, MSUHEP-19-012

1 Introduction

The Drell-Yan process Drell:1970wh is one of the most important and basic processes measured at the Large Hadron Collider at CERN. It is used for precision measurements of the $W^{\pm}$ mass, the $Z$ mass, and the weak mixing angle, as well as for new physics searches. In order to interpret the increasingly precise data, higher-order corrections must be included in the theory predictions. In fixed-order perturbation theory, the pure Quantum Chromodynamic (QCD) corrections to the cross section are known through to next-to-next-to-leading order Hamberg:1990np ; Harlander:2002wh ; Anastasiou:2003ds ; Melnikov:2006di , the pure Quantum Electrodynamics (QED) corrections to neutral gauge boson production and decay are known at next-to-next-to-leading order Berends:1987ab ; Blumlein:2019srk , whereas the exact electroweak (EW) corrections are known only at next-to-leading order Baur:2001ze ; Dittmaier:2001ay ; Baur:2004ig . Since these corrections turn out to be significant, the inclusion of exact mixed EW-QCD corrections is well-motivated; at the present time, results are known only for the mixed Quantum Electrodynamic (QED)-QCD corrections Kilgore:2011pa and in the pole approximation for the full EW-QCD corrections Dittmaier:2014qza ; Dittmaier:2015rxo .

The two-loop master integrals relevant to the virtual part of the mixed EW-QCD corrections to the Drell-Yan process were studied in Bonciani:2016ypc ; vonManteuffel:2017myy and all of them admit $\epsilon$-decoupled and, subsequently, $\epsilon\,\mathrm{d}\ln$ differential equations Henn:2013pwa . All of the double-box integrals are actually planar, and for most of them, no intrinsically non-rational prefactors appear in the $\epsilon$-decoupled differential equations written with respect to suitable kinematic variables. It has come to light, however, that the two most complicated two-loop integral topologies (right panels of Figures 1 and 2) involve three square roots at once which can not be rationalized simultaneously mainz ; Besier:2019hqd . Here and in the following, we will refer to a set of roots as unrationalizable if no locally invertible rational variable transformation exists, which turns all of them into rational functions at the same time. The presence of non-rational symbol letters Goncharov:2010jf in the $\epsilon\,\mathrm{d}\ln$ differential equations makes the standard approach to the integration of the differential equations in terms of multiple polylogarithms impossible.

For most Feynman integrals considered in the literature so far, kinematic square roots could be rationalized by a suitable choice of variables, e.g. by using a Landau variable for Feynman integrals with a two-mass threshold Barbieri:1972as , by using momentum twistor variables Hodges:2009hk for multi-leg Feynman integrals Gehrmann:2015bfy ; Bourjaily:2018aeq ; Caron-Huot:2018dsv , or by using diophantine equations to construct suitable variables JohannesAmplitudes2015talk . However, in several processes studied in the literature, no simultaneous rationalization of the square roots which appear could be found. In addition to the Drell-Yan process mentioned above, this is an issue for a subset of the planar two-loop master integrals for Bhabha electron-positron scattering Henn:2013woa ; Festi:2019 and for a subset of the planar two-loop master integrals for the next-to-leading order QCD corrections to Higgs plus jet production with full heavy top quark mass dependence Bonciani:2016qxi .

So far, master integrals satisfying $\epsilon\,\mathrm{d}\ln$ differential equations which involve non-rational symbol letters have typically been treated as generic Chen iterated integrals Chen:1977oja , with a focus on the Euclidean region of phase space. Given our experience from other processes Bonciani:2013ywa ; Gehrmann:2015ora , one may expect significantly faster and also more stable numerical evaluations for results directly expressed in terms of suitably chosen multiple polylogarithms. A constructive algorithm to obtain such a solution in the presence of root-valued symbol letters was used in vonManteuffel:2017hms to calculate Feynman integrals. However, the latter case is special since the symbol alphabet is univariate and simple, and a rational parametrization does in fact exist. In other words, it was unclear until now whether it should generically be possible to find linear combinations of multiple polylogarithms which solve $\epsilon\,\mathrm{d}\ln$ differential equations with symbol letters involving unrationalizable algebraic functions.

In this article, we show for the first time that it is indeed possible to integrate Feynman integrals of current phenomenological interest with unrationalizable roots in their symbol letters in terms of multiple polylogarithms, focusing primarily on the most complicated two-loop mixed EW-QCD Drell-Yan master integrals. Our method allows us to derive results written in terms of multiple polylogarithms well-suited for the physical region of phase space. The two-loop mixed EW-QCD Drell-Yan master integrals with two massive internal lines exhibit a rich structure of thresholds and pseudo-thresholds, which means that one must think carefully about the analytic structure of the involved functions. We discuss in detail a procedure to systematically filter the functions in a given region of phase space, such that no explicit Feynman $+i\,0$ prescriptions are required for the kinematic variables and spurious singularities of individual multiple polylogarithms at pseudo-thresholds of the Feynman integrals are avoided.

The outline of the article is as follows. In Section 2 we discuss the linear reducibility of examples which are well-known in the literature to frustrate the standard machinery of Feynman integral calculus. Specifically, we prove that the integral topologies of Figure 2, a five-line master integral for Bhabha scattering and a six-line master integral for Drell-Yan production, are linearly reducible for the first time. Linear reducibility is a technical criterion which is of interest here because linearly reducible Feynman integrals are guaranteed to be integrable in terms of multiple polylogarithms to all orders in $\epsilon$. In Section 3, we define a normal form basis with $\epsilon\,\mathrm{d}\ln$ differential equations for the two-loop master integrals for the mixed EW-QCD corrections to Drell-Yan production with two massive internal lines (see the right panel of Figure 1) and discuss a partial rationalization of the roots appearing in our integral basis definition. In Section 4, we show how to integrate the differential equations directly in terms of multiple polylogarithms even in the presence of root-valued symbol letters. In Section 5, we review the analytic continuation of multiple polylogarithms and outline our procedure to filter multiple polylogarithms with undesirable analytic properties out of our ansätze. In Section 6, we present results for the most complicated two-loop mixed EW-QCD Drell-Yan master integrals. In particular, we highlight the notable analytic features of our solution for the six-line integral from the right panel of Figure 2. We conclude in Section 7 and, for clarity, we give the set of complete $\epsilon\,\mathrm{d}\ln$ differential equations for the two-loop master integrals considered in this paper in Appendix A. In Appendix B we give an example for the construction of algebraic letters in the presence of multiple root-valued leading singularities.

2 Linear reducibility for algebraic symbol letters

In this section, we discuss the direct integration of a five-line master integral for the two-loop QED corrections to massive Bhabha scattering (left panel of Figure 2). This particular five-line integral is of special interest because its symbol letters are not simultaneously rationalizable Festi:2019 . It has been known for quite some time that the planar master integrals for the two-loop QED corrections to Bhabha scattering (see the left panel of Figure 1) satisfy $\epsilon$-decoupled differential equations Henn:2013woa , but, even at leading order in the $\epsilon$ expansion, it is not at all clear that integral (1) below may be expressed as a linear combination of standard multiple polylogarithms LDanilevsky ; Goncharov:1998kja ; Borwein:1999js . In fact, it was suggested in ClaudeAmplitudestalk that elliptic multiple polylogarithms BrownLevin might actually be required for such cases. It is therefore of some importance to demonstrate that one can actually integrate it to all orders in $\epsilon$ in terms of Goncharov or multiple polylogarithms using a HyperInt-like approach Brown:2008um ; Brown:2009ta ; Panzer:2014caa ; Bogner:2014mha ; Panzer:2015ida ; vonManteuffel:2013uoa ; Gehrmann:2013cxs ; vonManteuffel:2013vja ; Gehrmann:2014bfa ; Duhr:2019tlz ; Becchetti:2019tjy .

In fact, this may be readily achieved by studying the polynomials which appear in the polynomial reduction Brown:2008um ; Brown:2009ta ; Panzer:2015ida of the Symanzik polynomials at intermediate stages. The idea is to make changes of variables which factorize totally quadratic polynomials to allow for further integrations without producing square root-valued functions of the remaining Feynman parameters. Ultimately, the goal is to find a complete sequence of integration variables or, in the language introduced in Brown:2009ta and refined in Panzer:2014caa , to show that (1) is a linearly reducible Feynman integral in the chosen variables. Practically speaking, a Feynman integral is linearly reducible if it can be evaluated in terms of Goncharov polylogarithms, starting from the Feynman parametric representation but allowing for arbitrary variable changes along the way. In this case, the Feynman parametric representation is

[TABLE]

where

[TABLE]

are the integral’s Symanzik polynomials. We implement the delta function constraint at the very beginning by setting $\alpha_{5}=1$ inside $\mathcal{U}$ and $\mathcal{F}$. Our normalization for the integral is chosen to facilitate comparison with FIESTA 4 Smirnov:2015mct .

At the outset, the polynomial reduction of $\mathcal{U}$ and $\mathcal{F}$ indicates that just one integration can safely be performed, with respect to either $\alpha_{1}$ or $\alpha_{3}$. The irreducible polynomials in the remaining Feynman parameters which could potentially appear in the result after integrating out $\alpha_{1}$ are Panzer:2014caa :

[TABLE]

The final totally quadratic polynomial in ${\rm L}_{\alpha_{1}}$ is conveniently dealt with by using the first non-trivial variable change from the discussion of the period $P_{7,11}$ in Panzer:2015ida . It is readily apparent from the form of the polynomial that two powers of $\alpha_{2}$ factorize from it once the variable change

[TABLE]

is applied:

[TABLE]

Quite remarkably, the other polynomials in ${\rm L}_{\alpha_{1}}$ which depend on $\alpha_{4}$ still give rise exclusively to irreducible polynomials linear in $\alpha_{2}$ after changing variables:

[TABLE]

From Eqs. (2), (7), and (8), it is clear that we can now safely integrate out $\alpha_{2}$.

Rerunning the polynomial reduction with HyperInt after changing variables shows that the irreducible polynomials

[TABLE]

could appear in the result after integrating out $\alpha_{1}$ followed by $\alpha_{2}$. ${\rm L}_{\alpha_{1}\alpha_{2}}$ is more non-trivial to treat than ${\rm L}_{\alpha_{1}}$ and requires us to think carefully, both about what sort of result we expect to find and how the Mathematica-based private direct integration code written by one of us actually operates in practice. As should by now be clear, our goal is to delay the appearance of square roots until the final variable of integration is moved into the arguments of the Goncharov polylogarithms which depend on it. With the final integration of the most complicated weight-three functions in mind, it is useful to aim for a final integration domain of $[0,1]$. This will naturally produce Goncharov polylogarithms of argument $1$ which contain a complicated square root in the weights. Actually, due to the way that our integration script is written, it is most natural for us to integrate complicated finite Feynman parameter integrals such as the one considered in this section over the unit hypercube. This is because it is easier for our integration code to take a definite integral on $[0,1]$ than on $[0,\infty)$, due to the fact that no complicated argument inversion formulae for non-trivial Goncharov polylogarithms need to be derived if the upper integration endpoint is always set to $1$.

In this spirit, we initiate our treatment of ${\rm L}_{\alpha_{1}\alpha_{2}}$ by trivially mapping one of the remaining integration domains onto $[0,1]$ via the transformation

[TABLE]

It is also useful to note that two of the quadratic polynomials in Eq. (9) depend on just a single Feynman parameter each, a situation which strongly suggests that a change of kinematic variables would be advantageous. Indeed, with Euclidean $s$ and $t$ in mind, we see that the variable changes

[TABLE]

rationalize two square roots which would otherwise arise from the third irreducible polynomial of ${\rm L}_{\alpha_{1}\alpha_{2}}$,

[TABLE]

and the fourth irreducible polynomial of ${\rm L}_{\alpha_{1}\alpha_{2}}$,

[TABLE]

Curiously, transformation (11) was also found to be useful in the context of the multi-loop QED corrections to light-by-light scattering Caron-Huot:2014lda .

At this stage, it turns out that just a single non-trivial algebraic function of $x_{4}$ is generated if one attempts to naively continue the calculation by integrating out $x_{3}$:

[TABLE]

To our knowledge, a situation such as this was first discussed in the Feynman integral literature in the context of the fully-massless (three-loop) $K_{4}$ integral Panzer:2014gra . Using the method of parametrization by lines Panzer:2014gra ; Besier:2018jen , it is completely straightforward to derive a final variable change for $x_{4}$ which rationalizes (14) and produces an integration domain of $[0,1]$. As usual, one begins by identifying a suitable rational point of the algebraic variety

[TABLE]

For our purposes, the point

[TABLE]

is ideal.111We avoid the other obvious rational point, $\left(x_{4}^{(0)},\rho^{(0)}\right)=(0,1)$, because the resulting variable change does not map the $x_{4}$ integration domain onto the unit interval. Eq. (15) defines a hyperbola which can be rationally parametrized by a family of lines passing through (16),

[TABLE]

We combine Eqs. (15) and (17) to determine $x_{4}$ as a function of $y_{4}$:

[TABLE]

As desired, we see from the above that the integration domain for $y_{4}$ is the unit interval.

By combining together the various variable changes given above, we obtain the complete sequence of integration variables $\{\alpha_{1},\alpha_{2},x_{3},y_{4}\}$ and establish that, up to an overall normalization factor, integral (1) may be evaluated as a rational linear combination of Goncharov polylogarithms to all orders in $\epsilon$. As a sanity check, we explicitly evaluate all Feynman parameter integrals analytically using our direct integration code at $\mathcal{O}\left(\epsilon^{0}\right)$ and then evaluate the result obtained numerically to high precision using GiNaC Vollinga:2004sn ; Bauer:2000cp at the randomly chosen Euclidean phase space point

[TABLE]

We find

[TABLE]

which we could rapidly confirm agrees to five significant digits with an independent FIESTA 4 evaluation of the integral.

Let us stress once more that the unrationalizable square root

[TABLE]

where

[TABLE]

appears in our analytic integrations only at the very end, once all Feynman parameters have been integrated out. The most complicated, weight-four functions which appear in the final result have argument unity and weights which are non-trivial functions of the square root above. For the final integration, our integration code actually utilizes the generalized weights of vonManteuffel:2013vja . The idea is that, for the final integration over $y_{4}$, one can employ non-linear integrating factors of the form $\mathrm{d}y_{4}(\partial f(y_{4};v_{1},v_{2})/\partial y_{4})/f(y_{4};v_{1},v_{2})$, for non-linear irreducible polynomials $f(y_{4};v_{1},v_{2})$ in $y_{4}$. This allows for a more concise representation of the final result without giving up access to the well-tested numerical routines provided by GiNaC. If desired, a representation written strictly in terms of Goncharov polylogarithms may also be derived by fully factoring all of the generalized weights.

Of course, the direct integration approach described in this section is not without its limitations. For one thing, it would at first seem quite non-trivial to find an explicit integration order for e.g. the two-loop double box with two massive internal lines from the two-loop mixed EW-QCD corrections to Drell-Yan lepton production (right panel of Figure 1). Actually, with the help of HyperInt, it is comparatively easy to find a complete sequence of integration variables; for a fixed number of loops and legs, one gets the impression that the number of internal lines is a less important factor in determining the difficulty of a linear reducibility analysis than the total number of massive lines. Although we could directly integrate through to weight four the six-line basic scalar integral (right panel of Figure 2) which first inherits unrationalizable symbol letters from its unrationalizable leading singularity, we find the direct integration of the top-level two-massive-line integrals to be very cumbersome at the technical level. Furthermore, the produced Goncharov polylogarithms are very complicated, not minimal, and suboptimal for numerical evaluations with GiNaC.

An alternative approach which can avoid these issues is to apply the differential equation method. We shall see that the differential equation method has the very appealing feature that one can fit an ansatz of multiple polylogarithmic functions appropriate for each physical kinematic region separately. In the following sections, we discuss in some detail how this is achieved for the mixed EW-QCD Drell-Yan master integrals with two massive lines, including, in particular, the six- and seven-line integrals with an unrationalizable square root in their symbol letters.

3 An $\epsilon$ basis for the Drell-Yan master integrals with two massive lines

In this work, we are primarily concerned with neutral-current lepton pair production222The two-loop mixed EW-QCD corrections to the charged-current Drell-Yan process could be accessed by expanding in $1-m_{W}^{2}/m_{Z}^{2}$, as this would allow one to make effective use of the equal-mass integrals relevant to the neutral-current process. in quark-antiquark annihilation,

[TABLE]

where all external particles are taken massless and on their mass shell. The most complicated master integrals for the two-loop mixed EW-QCD corrections to above process are those with two internal lines of mass $m$, where, depending on the parent Feynman diagram, $m$ may refer to either $m_{W}$ or $m_{Z}$. In particular, as mentioned above, it has been known to us for some time that the integral topologies from the right panels of Figures 1 and 2 actually contain master integrals with unrationalizable symbol letters in their $\epsilon$-decoupled differential equations. To index the seventeen master integrals with two massive internal lines, it suffices to consider a single integral family based on the two-loop planar double box with two massive internal lines from Figure 1 (Family $\mathrm{E}$ from Table 1). For completeness, we will also give definitions for those normal form integrals with zero or one massive internal lines which appear on the right-hand sides of the $\epsilon$-decoupled differential equations for the master integrals with two massive internal lines. To index these auxiliary integrals, we reintroduce two additional integral families which were already studied in the physical region by two of us some time ago vonManteuffel:2017myy (Family $\mathrm{A}$ and Family $\mathrm{C}$ from Table 1).

In order to obtain a closed system of differential equations for the masters with two massive internal lines, we need to consider thirty-six integrals in total. Our notation for Feynman integrals in this section is exactly that of vonManteuffel:2017myy (i.e. dots for doubled propagators, heavy lines for massive propagators, numerator insertions written in square brackets, and $\mathrm{\accentset{\rule{2.79996pt}{0.8pt}}{F}}$:$x$ for the crossed version of sector $x$ from family $\mathrm{F}$). For the kinematic invariants we use

[TABLE]

In the following, we keep the dependence on the internal mass parameter $m$ implicit for the sake of brevity and since it is anyway clear from the thick-line notation. We build up our normal form basis out of the following thirty-six Feynman integrals:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

From the general principles discussed in Gehrmann:2014bfa , one can readily cast (25) above into a normal form basis for the integrals of interest.333Due to the presence of unrationalizable square roots, the available public software packages for the construction of a normal form basis of integrals Prausa:2017ltv ; Gituliar:2017vzm ; Meyer:2017joq are not applicable to the problem at hand. Abbreviating the three square roots which appear as

[TABLE]

we find:

[TABLE]

For $\mathbf{m}_{1},\ldots,\mathbf{m}_{36}$, we employ the integration measure

[TABLE]

with $\epsilon=(4-d)/2$, which allows us to consider our integrals to be functions of two dimensionless kinematic variables.444 The normal form basis (3) closely resembles the one given in Bonciani:2016ypc , but we correct typos and employ a different overall normalization. Note that $r_{1}^{2}=s(s-4\,m^{2})$ is positive for both $s>4\,m^{2}$ and $s<0$. However, one must take care when attempting to simplify $r_{1}$ to make contact with the Euclidean-region analysis of Bonciani:2016ypc ; $r_{1}=s\sqrt{1-4\,m^{2}/s}$ for $s>4\,m^{2}$, but $r_{1}=-s\sqrt{1-4\,m^{2}/s}$ for $s<0$. For instance,

[TABLE]

One achieves a rationalization of two of our three roots, $r_{1}$ and $r_{2}$, with the parametrization  Bonciani:2016ypc

[TABLE]

and we hereafter work primarily with the variables $w$ and $z$. It is immediate from (3) that, in the $(w,z)$ representation, the only normal form integral which involves a square root in its definition is $\mathbf{m}_{32}$ (i.e. the integral from the right panel of Figure 2). Indeed, we see from (32) and Appendix A of Bonciani:2016ypc that

[TABLE]

in the region $s>4\,m^{2}$ (see also Section 5.1 for a detailed discussion of the $(w,z)$ parametrization).

In the following, we will replace $r_{1}$, $r_{2}$ and $r_{3}$ according to (33)-(35) in the normal form definitions $\mathbf{m}_{1},\ldots,\mathbf{m}_{36}$ and use these partially rationalized expressions for the entire physical region of phase space. That is, only the definition of $\mathbf{m}_{32}$ involves a root

[TABLE]

in the prefactor, which requires a non-trivial analytic continuation from $s>4\,m^{2}$ to other regions of phase space on its own. Using these definitions, it is straightforward to obtain the differential equations in the $\epsilon$-decoupled form discussed in the next section. We used Reduze 2  vonManteuffel:2012np ; Studerus:2009ye ; Bauer:2000cp ; fermat to compute the integration-by-parts identities required to derive the differential equations.

4 Integrating root-valued symbols in terms of multiple polylogarithms

Our normal form basis for the mixed EW-QCD Drell-Yan integrals with two massive internal lines, (3), is chosen to bring the associated differential equations into an $\epsilon$-decoupled form Kotikov:2010gf ; Henn:2013pwa :

[TABLE]

where $l_{k}$ are the symbol letters, $A^{(k)}$ are matrices of rational numbers, and $i,j=1,\ldots,36$. It was demonstrated already in Bonciani:2016ypc that a $\epsilon\,\mathrm{d}\ln$ form does in fact exist, but their choice of the symbol letters is not optimal for our purposes. For now, we proceed with the understanding that some algebraic symbol letters appear in Eq. (37), but leave their number and precise form to be determined by the analysis below.

Using $\vec{x}=(w,z)$ and $\vec{\mathbf{m}}=(\mathbf{m}_{i})$, $i=1,\ldots,36$, one can give a formal solution of Eq. (37) in terms of Chen iterated integrals as

[TABLE]

where

[TABLE]

$\vec{\mathbf{m}}(\epsilon,\vec{x}_{0})$ are the boundary constants of the master integrals at the point $\vec{x}=\vec{x}_{0}$, $\gamma$ is a piece-wise smooth path connecting $\vec{x}_{0}$ and $\vec{x}$, and the path-ordered exponential in Eq. (38) is defined in the usual way as an infinite series of integral operators acting to the right,

[TABLE]

In (40), the product of two or more $\mathrm{d}\ln$ terms is understood as an instruction to take the corresponding iterated integral of the kernel along the path $\gamma$. As a concrete example, let us consider a straight-line path $\tau$ on the real axis from [math] to $x$. One can write, e.g.:

[TABLE]

Note that, in the above example, one could equally well identify the iterated integral as a Goncharov polylogarithm,

[TABLE]

If all symbol letters are linear, it is always possible to integrate Chen iterated integrals in closed form using $G$ functions. This is no longer true, however, if one encounters non-linear or non-rational letters. In the former case, it is often possible to choose an appropriate integration order for which the non-linear letters appear only in the final integration kernel. At this stage, one may employ the generalized weights used in Section 2 to obtain concise results. If one encounters non-rational symbol letters, the standard integration algorithms for $G$ functions cannot be applied. In many cases, a transformation can be found which simultaneously rationalizes all letters of the alphabet. Nevertheless, it can be proven that the symbol alphabet of the five most complicated integrals from (3), $\{\mathbf{m}_{32},\ldots,\mathbf{m}_{36}\}$, cannot be rationalized mainz ; Besier:2019hqd . Of course, it is a priori not obvious whether the functional basis for $\{\mathbf{m}_{32},\ldots,\mathbf{m}_{36}\}$ consists solely of multiple polylogarithms; at the outset, it is certainly possible that one could have to deal with a more involved space of functions.

Fortunately, our linear reducibility analysis from Section 2 guarantees that standard multiple polylogarithms suffice; what we need is a better way to integrate $\epsilon\,\mathrm{d}\ln$ differential equations. A clear alternative is to proceed by matching the symbol of the Chen iterated integrals to a suitable ansatz built out of logarithms and ${\rm Li}$ functions. For rational alphabets, a method was provided in Duhr:2011zq by Duhr, Gangl, and Rhodes to construct suitable ${\rm Li}$ function arguments, such that the functional basis contains no spurious letters. At weight one, it is clear that logarithms of the letters are admissible functions. The non-trivial step is to find a sufficiently large set of admissible ${\rm Li}_{n}$ arguments; once these arguments have been found, it is straightforward to construct all admissible arguments for the ${\rm Li}$ functions of depth greater than one. By considering the symbol of these functions,

[TABLE]

we see that it is desirable to admit only those function arguments $f$ which have the property that both $f$ and $1-f$ may also be written as a power product of the symbol letters. In practice, one therefore forms power products $f$ out of the letters and tests if $1-f$ factorizes over the alphabet. The symbols of higher-depth Li functions are more complicated and lead to additional constraints. To treat ${\rm Li}_{2,1}$, ${\rm Li}_{3,1}$, and ${\rm Li}_{2,2}$, let $\mathcal{F}$ be the union of the set of admissible ${\rm Li}_{n}$ function arguments and the set $\{1\}$. Then, the symbol dictates that a possible pair of arguments for ${\rm Li}_{n_{1},n_{2}}$, $(f_{i},f_{j})$ such that $f_{i},f_{j}\in\mathcal{F}$, is admissible if $1-f_{i}f_{j}$ factorizes over the alphabet. It has proven useful for practical applications to impose further constraints on the functional basis to ensure real-valuedness and good numerical performance. An implementation of this method written by one of us in Mathematica has been applied successfully to various processes vonManteuffel:2013uoa ; Gehrmann:2013cxs ; Gehrmann:2014bfa ; Gehrmann:2015ora ; vonManteuffel:2017myy ; Becchetti:2019tjy .

In the presence of square roots, we use a heuristic factorization algorithm to detect admissible function arguments. For a given expression $g$ we are interested in factorizations of the form

[TABLE]

with a rational number $c$ and $a_{n}\in\mathbbm{Q}$. It is non-trivial to find such factorizations using standard computer algebra systems due to the presence of the root $r$ in the symbol letters. We observe that the factorization (44) implies

[TABLE]

Replacing the variables by numerical samples allows us to find these relations using heuristic integer relation finders. To find the required factorizations, we employ the Lenstra-Lenstra-Lovász (LLL) algorithm Lenstra1982 implemented in PARI/GP PARI2 for a parallelized C++ code written by one of us.

We would like to stress that the definition of the symbol letters is not unique. One can replace a letter by power products of letters and it is a priori unclear which choice is optimal for practical purposes. For example, we find that out of the 17 letters presented in Bonciani:2016ypc , only 16 combinations are actually required for the integration of the Drell-Yan integrals. Furthermore, since we consider actual derivatives of Feynman integrals, we are not sensitive to numerical letters like 2. In principle, one needs to include also ad-hoc letters like $-1$, $2$, etc. in the construction of function arguments $f$. This is a problem occurring also in the purely rational case and does not seem to be a major obstacle in practice. In the applications we have considered so far, including $-1$ and $2$ was sufficient. For the case of the Drell-Yan integrals, we were able to absorb the letter $2$ by a redefinition of our symbol alphabet, as will be explained below.

In practice, we encounter two main problems specific to the case of non-rational symbol alphabets:

1. (i)

In general, one needs to allow for non-integer powers, e.g. $1/2$, $1/4$, etc. when forming power products $f$. Consequently, one may have to test many more expressions than in the rational case in order to construct enough function arguments to successfully “integrate the symbol.” In fact, without additional constraints, the inherent combinatorial complexity makes this extension of the Duhr-Gangl-Rhodes method too costly for two-loop calculations of current phenomenological interest. 2. (ii)

Factorization over algebraic functions is not unambiguously defined; in principle, it is not clear where to “stop” factorizing. Consider, for instance, the set of letters $\{\sqrt{x},\sqrt{y},x-y\}$ and note that it is a priori unclear whether one ought to factor the third letter further, i.e. $x-y=(\sqrt{x}+\sqrt{y})(\sqrt{x}-\sqrt{y})$.

In the case of the Drell-Yan integrals, both of these problems can be tackled by the following observations: given an alphabet, we introduce the subset of rational letters, $\mathcal{L}_{R}$, and the subset of intrinsically non-rational letters, $\mathcal{L}_{A}$. Let us assume that in $\mathcal{L}_{A}$ we encounter a square root, $r$. We find it natural to take $r$ itself to be an element of $\mathcal{L}_{A}$. For a given algebraic letter with a non-trivial rational part, $l_{a}$, we define the conjugated letter, $\bar{l}_{a}$, by making the replacement $r\rightarrow-r$ in $l_{a}$. For each $l_{a}\in\mathcal{L}_{A}$ we observe that $\bar{l}_{a}$ can be written as a power product of other symbol letters; that is to say, one could exchange any letter for its conjugate without affecting the singularity structure of the alphabet. Furthermore, we observe that the product $l_{a}\bar{l}_{a}$ always factorizes over the rational part of the alphabet, $\mathcal{L}_{R}$. We find these observations to be quite compelling and therefore conjecture that these structural constraints hold also for other symbol alphabets of this type.

Our conjecture is a strong one, since it implies that one can predict the form of the remaining unrationalizable symbol letters which appear in the presence of the square root $r$: knowing only the rational part of the alphabet, $\mathcal{L}_{R}$, and the unrationalizable square root which appears, it allows one to construct the algebraic part of the alphabet, $\mathcal{L}_{A}$, unambiguously. As we shall see, our new insight allows for a drastic simplification of the algebraic part of the symbol alphabet relative to what was presented in Eqs. (5.13)-(5.19) of Bonciani:2016ypc and what we were able to derive ourselves initially using the above-mentioned heuristic factorization code.

We want to show this in some detail for the EW-QCD Drell-Yan master integrals with two massive internal lines. At the outset, after carefully considering various possibilities, we find the rational alphabet

[TABLE]

and the intrinsically algebraic alphabet

[TABLE]

where $r=\sqrt{(1+w^{2}z^{2})(w+z)^{2}+2w\,z(w-z)^{2}+4w\,z^{2}(1+w^{2})}$, using our heuristic factorization code. Note that $r$ already appears as a letter in Eq. (47). However, there are also two rather complicated-looking symbol letters with terms involving $r^{2}$.

As stated above, knowing only $r$ and $\mathcal{L}_{R}$, we can also construct an improved representation of the algebraic part of the alphabet by making an ansatz of the form

[TABLE]

where $q_{a}$ is a rational function in $w$ and $z$, and then requiring that $l_{a}\bar{l}_{a}$ factorize over $\mathcal{L}_{R}$. In practice, it is more convenient to directly make an ansatz of the form

[TABLE]

and then solve for the unknowns $a_{ij}$. The algorithm to construct (simple) algebraic letters then reads:

In order to also reproduce the original letters we would need to allow for a polynomial prefactor in front of the $r^{2}$ term, which we deliberately avoid here to simplify our construction. In practice, we run the algorithm with a degree $n$ at least sufficiently large to be able to express the differential equations in terms of the new alphabet.

By using this procedure up to degree $4$, we find the following algebraic letters for the Drell-Yan integrals:

[TABLE]

The overall factors of $1/2$ in Eq. (4) are judiciously chosen after the fact to prevent the appearance of explicit factors of $2$ in our Li function arguments; as alluded to above, our original construction of the functional basis appended $2$ as an auxiliary letter. The presence of this factor can also be understood in light of the factorization property: one can easily check that all letters in $\tilde{\mathcal{L}}_{A}$ after multiplication with their conjugate, indeed factorize over the rational part of the alphabet without the presence of the factor $2$ in this product. In order to find the correct normalization, one can run the above algorithm first adding the letter $2$ in the rational part of the alphabet.

Using our heuristic approach to algebraic function factorization, we immediately find that the non-trivial elements of $\mathcal{L}_{A}$ can indeed be expressed as power products of letters drawn from the improved alphabet:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

In summary, we define the full symbol alphabet for the two-massive-line Drell-Yan integrals, $\mathcal{L}=\mathcal{L}_{R}\cup\tilde{\mathcal{L}}_{A}$, in terms of positive definite letters for physical phase space points which satisfy $-1<w<0$ and $w<z<-w^{2}$ (this component of $(w,z)$ space corresponds to part of the $s>4\,m^{2}$ region, see Section 5.1 for details). From Eqs. (4) and (4), we have:

[TABLE]

Most importantly, our improved representation of the alphabet effectively solves the problem of finding power products of high degree, since we observe in practice that we do not need to consider square roots of any letter in $\mathcal{L}$. For the integration of the EW-QCD Drell-Yan master integrals with two massive internal lines through to weight four, it was enough to consider power products of symbol letters derived above up to total degree $9$.

The construction above employed a representation with only a single root-valued leading singularity. In Appendix B, we consider a one-loop integral which involves five different root-valued leading singularities. We demonstrate that the algebraic part of the symbol alphabet can be constructed without any reparametrization by generalizing the procedure above to the case of multiple roots.

5 Analytic continuation and optimization of the functional bases

In this section, we review the salient features of the $(w,z)$ representation for the mixed EW-QCD corrections to Drell-Yan production introduced in Section 3 above, as well as subtleties one encounters when analytically continuing multiple polylogarithms. Our primary goal in Section 5.1 is to motivate the analysis of Section 5.2, where we show how we avoid all explicit analytic continuations and $+i\,0$ prescriptions by partitioning the physical phase space and finding several solutions to the differential equations in terms of well-behaved Li functions valid in judiciously chosen regions. However, we also find it useful to review the fundamentals of analytic continuation of Feynman integrals as well as some details relevant to our specific representation.

5.1 Analytic continuation

Feynman’s $+i\,0$ prescription for the propagators determines the value of a given Feynman integral in a specific region of phase space unambiguously. In principle, one could imagine to solve the Feynman integral in each region of phase space separately. Alternatively, one can try to solve the integral in one region and then use the solution to obtain a result for it in a neighboring region by analytic continuation. The latter method is of particular interest for the method of differential equations, since it typically involves regularity conditions in some regions of phase space, and one needs to transport this knowledge to the region of interest.

If we want to continue a specific representation of the solution based on just the solution itself (without reference to the original Feynman integral), we need to make sure Feynman’s $+i\,0$ prescription is maintained by appropriate complex values of the kinematic parameters. It is is essential to observe that the analytic continuation is along a path in complex phase space and that the $+i\,0$ prescriptions must be respected for all points along this path, not just for the start and endpoint. This is in general non-trivial and needs to be checked for the representation at hand.

What is commonly referred to as “analytic continuation” in the physics literature should really be regarded as a two-step procedure in general:

1. (i)

The actual analytic continuation in the mathematical sense: given a solution for one region, derive a solution for a connected region in some representation. 2. (ii)

A possible change of functional representation in the new region such that no explicit $+i\,0$ prescription is necessary.

We will now work out the details for our current application.

From the second Symanzik polynomials of our Drell-Yan master integrals, we see that their Euclidean region is given by $s<0$, $t<0$ and $m^{2}>0$.555For non-planar topologies with four massless external legs, also cuts in $u$ must be taken into account, which may actually prevent the existence of a Euclidean region with $s+t+u=0$ Tausk:1999vh . Considering other regions of phase space, we observe that the Feynman propagator prescription can effectively be implemented by the replacements $s\rightarrow s+i\,0$, $t\rightarrow t+i\,0$, since these are external scales, and $m^{2}\rightarrow m^{2}-i\,0$, since this is an internal scale. The $+i\,0$ prescription is relevant for $s$ ($t$) whenever $s$ ($t$) is positive. For the discussion of $w$ and $z$, we will see that it is sufficient to view $m^{2}$ as being normalized to 1 without any imaginary part.

Let us focus on the $(w,z)$ parametrization of our integrals in the physical region of phase space, $s>0$, $t<0$ and $m^{2}>0$. It was noted in Bonciani:2016ypc that the $(w,z)$ representation has a point of non-analyticity at the physical two-mass threshold $s=4\,m^{2}$ and a rather different character depending on whether $s$ is above or below this threshold. We find it natural to adopt the definitions

[TABLE]

We can see from the above that $w$ and $z$ are real variables which satisfy $-1<w<0$ and $-1<z<-w^{2}$ in the region $s>4\,m^{2}$. With our choices (56) and (57), we find that the Feynman prescription implies $w\rightarrow w+i\,0$ and $z\rightarrow z+i\,0$ in this region. On the other hand, both $w$ and $z$ become pure phases in the region $0<s<4\,m^{2}$. Therefore, it makes sense to explicitly extract their real and imaginary parts,

[TABLE]

to emphasize that, in the below-threshold region, the imaginary parts of $w$ and $z$ are fixed in terms of the real parts of $w$ and $z$. In particular, we can deduce from Eqs. (58) and (59) that $\mathrm{Re}[w]$ and $\mathrm{Re}[z]$ satisfy $-1<\mathrm{Re}[w]<1$ and $-1<\mathrm{Re}[z]<1-2\left(\mathrm{Re}[w]\right)^{2}$ in the below-threshold region. A visualization of the full physical phase space in the $(w,z)$ representation is given in Figure 3.

In the following, we will study how to analytically continue solutions between different regions in $s$. In practice, we work with simple straight-line paths in the complex $(w,z)$ space, checking after the fact that the chosen paths of analytic continuation always preserve the $+i\,0$ prescription for the original kinematic variables. In what follows, we carefully go through typical elementary examples of analytic continuation in order to clearly illustrate the subtleties for our ${\rm Li}$ functions which one must be wary of to avoid introducing errors.

First, let us illustrate the importance of taking into account the complete path of analytic continuation rather than just the start and endpoint. Consider the analytic continuation of \ln\big{(}w^{2}\big{)} along a straight-line path from $w_{i}=w_{(0)}+i\,\delta$ to $w_{e}=-w_{(0)}+i\,\delta$, for real $w_{(0)}$ and $\delta$ such that $0<w_{(0)}<1$ and $0<\delta<w_{(0)}$ (see Figure 4). Such a path allows one to connect solutions for $s<0$ with solutions for $s>4\,m^{2}$. Naively, one might be tempted to erroneously implement the analytic continuation as \ln\big{(}w_{i}^{2}\big{)}=\ln\big{(}w_{(0)}^{2}\big{)}=\ln\big{(}(-w_{(0)})^{2}\big{)}=\ln\big{(}w_{e}^{2}\big{)} and incorrectly conclude that the analytic continuation is trivial along the chosen path. The problem here is that the logarithm is a multi-valued function and one must therefore carefully check whether or not the specified path of analytic continuation forces the polynomial function argument to cross the branch cut of the logarithm on the negative real axis, $(-\infty,0)$. In the absence of any branch cut crossings, knowledge of the endpoint of the path is sufficient. However, in the presence of one or more branch cut crossings, the function leaves its principal Riemann sheet and one must add an appropriate monodromy contribution, taking into account the details of the path of analytic continuation.

Let us spell out in detail how to correctly analytically continue solutions valid for different values of $s$. We parametrize our chosen path in $w$ space as

[TABLE]

with a parameter $v\in[0,1]$. We then have

[TABLE]

for the argument of \ln\big{(}w^{2}\big{)}. As depicted in Figure 5, our path in $w$ space takes $w^{2}(v)$ from just above the point $w_{(0)}^{2}-\delta^{2}$ on the positive real axis to just below the point $w_{(0)}^{2}-\delta^{2}$ and, crucially, it passes through the negative real axis at $w(1/2)=\delta\,i$. The monodromy contribution to \ln\big{(}w^{2}\big{)} in this case is well-known to be simply $2\pi i$, due to the fact that our path induces just one counter-clockwise branch cut crossing. Therefore, the analytically continued function is \ln\big{(}w^{2}\big{)}+2\pi i in a neighborhood of the endpoint of our chosen path. This conclusion may be quickly checked by rewriting the function \ln\big{(}w^{2}\big{)} as $2\ln(w)$ in the Euclidean region where $0<\mathrm{Re}[w]<1$ before carrying out the analytic continuation along our chosen path. When the function is viewed in this alternative way, no branch cut crossing is induced and one can simply rewrite the function to explicitly extract its imaginary part: for $\mathrm{Im}[w]>0$, 2\ln(w)=2\big{(}\ln(-w)+\pi i\big{)}=2\ln(-w)+2\pi i=\ln\big{(}w^{2}\big{)}+2\pi i. In the two-step picture mentioned above we thus have

[TABLE]

Depending on the type of representation we consider, either the analytic continuation (i) or the rewriting of the function to be independent of $+i\,0$ prescriptions (ii) is trivial in this example.

The relevant monodromy contributions become more complicated for higher-weight ${\rm Li}$ functions. Already for ${\rm Li}_{2}$, a new feature emerges: moving across its branch cut on the positive real axis, $(1,\infty)$, onto a Riemann sheet other than the principal one actually exposes the existence of a hidden branch point at 1. To see this, consider Euler’s identity for the dilogarithm of $1-w^{2}$,

[TABLE]

The representation furnished by the right-hand side of (64) and the above discussion of \ln\big{(}w^{2}\big{)} make it clear that {\rm Li}_{2}\big{(}1-w^{2}\big{)} has non-trivial monodromy for the path $w(v)$ defined above in Eq. (60). For the path $w(v)$, we find:

[TABLE]

This implies that {\rm Li}_{2}\big{(}1-w^{2}\big{)} picks up a monodromy contribution of -2\pi i\ln\big{(}1-w^{2}\big{)} due to the branch point at $w=1$, despite the fact the function is continuous at that point on the principal sheet.

When considering general analytic continuations of ${\rm Li}_{n_{1},\ldots,n_{k}}$ functions, one must take into account all function arguments to obtain the monodromy contributions as one moves along the path of analytic continuation. Building on the work of Goncharov Goncharov:2001iea , it was shown in Duhr:2012fh how one can easily compute the monodromy of an arbitrary ${\rm Li}$ function in terms of monodromies of simple logarithms using the coproduct. In particular, the presence of a monodromy contribution can be detected by studying the first entry of the symbol. For our master integrals, the sheer number of distinct branch cut crossings which can arise from the arguments of the various ${\rm Li}$ functions renders the analytic continuation between different regions rather involved in practice. While our final goal is to obtain representations in terms of well-behaved ${\rm Li}$ functions for kinematic regions, we find it convenient to also employ auxiliary representations for the purpose of analytic continuations.

If available, a representation in terms of $G$ functions with the kinematic variables in the argument avoids many of the subtleties involved in the continuation described above. For the integrals with unrationalizable alphabets, however, this is not an option. For such cases, we find it useful to use expansions around regular and singular points Lee:2017qql for the continuation. For this part of the analysis, we do not include the unrationalizable root of $\mathbf{m}_{32}$ in (A) in the definition of our basis in order to work with rational differential equations. Using high-precision numerical evaluations for two expansion points with an overlapping region of convergence and the PSLQ algorithm PSLQ , one can effectively transport analytic integration constants.

What has been discussed so far allows us to construct a functional basis for a given region of phase space, integrate the symbol in terms of these functions, and relate solutions for different regions by analytic continuation in order to fix the integration constants. We will now discuss how to construct domain-restricted but well-behaved ansätze of multiple polylogarithms, which don’t require an explicit $+i\,0$ prescription and perform well numerically.

5.2 Optimizing the bases of multiple polylogarithms

Starting from a Duhr-Gangl-Rhodes basis of multiple polylogarithms for the integrals $\vec{\mathbf{m}}$, we wish to remove ${\rm Li}$ functions which have suboptimal analytic properties for physical kinematics ($s>0$). In fact, one finds that even the purely rational function arguments allowed by the letters (4) are non-trivial to treat systematically and, we will therefore restrict the discussion in this section to this rational subset. The main idea is to make a partition of the physical phase space into regions $\mathcal{D}_{i}$ such that, inside each region, a solution to the differential equations may be constructed out of ${\rm Li}$ functions which never diverge or move off of their principal Riemann sheets for arbitrary phase space trajectories contained in $\mathcal{D}_{i}$. Although our primary goal is to show how to one can avoid supplying explicit $+i\,0$ prescriptions for $w$ and $z$, we also find it convenient to impose further aesthetic criteria on our polylogarithmic bases in order to simplify our results. For example, we find it useful in each ansatz to give precedence to those functions which do not involve the symbol letter $1-w+w^{2}$, since this ensures that the master integrals which do not depend on $1-w+w^{2}$ are manifestly free of $1-w+w^{2}$ at the level of functions.666The presence or absence of the letter $1-w+w^{2}$ for particular master integrals is linked to the presence or absence of a one-mass threshold at $s=m^{2}$.

The first step of our analysis is to study the letters of the rational alphabet, $\mathcal{L}_{R}$, above and below the physical two-mass threshold at $s=4\,m^{2}$. To proceed, we must determine under what conditions the logarithms of the letters either diverge or move off of their principal Riemann sheets. We find that, in practice, it is simplest to use Eqs. (56) and (57) above threshold, Eqs. (58) and (59) below threshold, and the Mathematica function Reduce to work out how various constraints on polynomials of $w$ and $z$ map back to conditions on the original and more familiar kinematic variables, $s$, $t$, and $m^{2}$. We find that Reduce works in an efficient way when the system of inequalities to be reduced is formulated in terms of real-valued parameters and the solution to the system does not involve roots of high-degree polynomials. Fortunately, these assumptions are always satisfied for the ${\rm Li}$ functions which have symbols built out of letters from $\mathcal{L}_{R}$.777 In Section 6, the method described in this section is also used to prepare an ansatz involving ${\rm Li}$ functions which have symbols built out of letters from both $\mathcal{L}_{R}$ and $\mathcal{L}_{A}$. To check whether it is possible to live without a $+i\,0$ prescription, we must first understand for what values of $s$, $t$, and $m^{2}$ both the real and imaginary parts of the letters vanish simultaneously and for what values of $s$, $t$, and $m^{2}$ the imaginary parts of the letters vanish while their real parts happen to be negative:

• •

$\ln(l_{1})=\ln(1-w)$ diverges at the phase space boundary point $s=0$.

• •

$\ln(l_{2})=\ln(-w)$ has no issues for $m^{2}>0$, except at the phase space boundary point $s=0$ where it becomes ill-defined.

• •

$\ln(l_{3})=\ln(1+w)$ diverges at the two-mass threshold $s=4\,m^{2}$.

• •

$\ln(l_{4})=\ln(1-w+w^{2})$ diverges at the one-mass threshold $s=m^{2}$.

• •

$\ln(l_{5})=\ln(1-z)$ diverges at the phase space boundary point $s=-t=2\,m^{2}$.

• •

$\ln(l_{6})=\ln(-z)$ has no issues for $m^{2}>0$, except at the phase space boundary point $s=-t=2\,m^{2}$ where it becomes ill-defined.

• •

$\ln(l_{7})=\ln(1+z)$ diverges at the phase space boundary where $t=0$ and at the two-mass threshold $s=4\,m^{2}$.

• •

$\ln(l_{8})=\ln(1-w\,z)$ diverges at the two-mass threshold $s=4\,m^{2}$.

• •

$\ln(l_{9})=\ln(1+w^{2}\,z)$ diverges at the phase space boundary $s=-t$ for $0\leq s\leq 2\,m^{2}$ and at the two-mass threshold $s=4\,m^{2}$.

• •

$\ln(l_{10})=\ln(-z-w^{2})$ diverges at the phase space boundary $s=-t$ for $s\geq 2\,m^{2}$ and, in particular, at $s=4\,m^{2}$.

• •

$\ln(l_{11})=\ln(z-w)$ diverges when $t=-m^{2}$ for $s\geq m^{2}$ and, in particular, at $s=4\,m^{2}$. $\ln(l_{11})$ is also not analytic at $s=-\frac{4\,m^{2}t}{m^{2}-t}$ for $-m^{2}\leq t\leq 0$.

The next step is to determine from the above data how many regions of the physical phase space it makes sense to consider separately. In a sense, this step is the most non-trivial because different choices besides the one we ultimately make are possible. Our choices are guided by imposing analytic properties of the integrals onto the basis functions themselves. Due to the absence of $u$ dependence in our planar double boxes, we expect them to have a regular limit as $u$ approaches [math]. We see from the above that, depending on whether $s$ is less than or greater than $2\,m^{2}$, regularity as $u\rightarrow 0$ suggests the absence, respectively, of $\ln(l_{9})$ or $\ln(l_{10})$ at weight one. Of course, we wish to achieve that also the higher weight basis functions lack logarithmic singularities as $l_{9}\rightarrow 0$ for $0<s<2\,m^{2}$ and $l_{10}\rightarrow 0$ for $s>2\,m^{2}$. This can be achieved by choosing the point $s=2\,m^{2}$ as a region boundary and then imposing an appropriate first-entry condition Gaiotto:2011dt ; Dixon:2011ng ; CaronHuot:2011ky ; Duhr:2012fh on the symbol to select suitable functions. In other words, we remove from consideration all ${\rm Li}$ functions which have the letter $l_{9}$ in the first entry of their symbols in $0<s<2\,m^{2}$ regions and all functions which have the letter $l_{10}$ in the first entry of their symbols in $s>2\,m^{2}$ regions. All of the other letters exhibit problems only at the physical one- and two-mass thresholds or at other exceptional points on the boundary of the physical phase space like $s=0$. The only logarithm with a spurious singularity in the physical region is $\ln(l_{11})$, which can be dealt with by imposing an additional first-entry condition on $l_{11}$ in all regions.

Ultimately, we find it natural to partition the physical phase space into four regions:

[TABLE]

where, as suggested by Eqs. (58) and (59), we have eliminated $\mathrm{Im}[w]$ and $\mathrm{Im}[z]$ below the two-mass threshold by exploiting the fact that $w$ and $z$ become pure phases:

[TABLE]

Even for the twelve two-mass master integrals which have rational symbols in the $(w,z)$ representation, it is not obvious that it suffices to consider a partition of the physical phase space into just four separate regions. Fortunately, no further subdivisions are necessary888This conclusion remains unchanged if one passes from $\mathcal{L}_{R}$ to the full alphabet. and it is even possible to consistently impose stronger constraints in each region on the subset of ${\rm Li}$ functions which survive our first-entry cuts.

Our basis of ${\rm Li}$ functions in region $\mathcal{D}_{a}$ may be further refined along the lines described in vonManteuffel:2017myy , where the EW-QCD Drell-Yan master integrals with a single massive line were evaluated in the physical region for the first time. Consider, for example, the set of 192 ${\rm Li}_{n}$ arguments consistent with our first-entry conditions in this region:

[TABLE]

Above threshold, the letters are real-valued functions and it turns out that the ${\rm Li}_{n}$ functions themselves may be chosen to be real-valued in $\mathcal{D}_{a}$. By using the Reduce function of Mathematica, we can rapidly check whether any of the above function arguments assume values larger than one for $(w,z)$ in $\mathcal{D}_{a}$. As one might expect, this is a non-trivial constraint to impose on (5.2) and only 148 function arguments survive. We find

[TABLE]

for our final set of preferred ${\rm Li}_{n}$ function arguments above the two-mass threshold. We can proceed analogously to filter the set of ${\rm Li}_{2,1}$, ${\rm Li}_{3,1}$, and ${\rm Li}_{2,2}$ function argument pairs which survive our first-entry cuts. For the sake of brevity, let us simply state that we find 2496 preferred pairs of function arguments $(f_{i},f_{j})$ after using Reduce to throw away all pairs for which $f_{i}>1$ or $f_{i}f_{j}>1$ for $(w,z)$ in $\mathcal{D}_{a}$. As explained in vonManteuffel:2017myy , we can further improve the numerical performance of our functional basis in GiNaC by reordering the lists of function arguments obtained to ensure that those which lead to convergent power series expansions are given precedence.

In region $\mathcal{D}_{b_{1}}$, the first-entry conditions are the same as for $\mathcal{D}_{a}$ and again lead to (5.2). As mentioned above, we can use Reduce to efficiently check which ${\rm Li}_{n}$ basis functions have suitable analytic properties below the two-mass threshold as well; despite the fact that $w$ and $z$ become complex functions of the original kinematic variables in $\mathcal{D}_{b_{1}}$ and the remaining below-threshold regions, we can proceed as before provided that we first eliminate $\mathrm{Im}[w]$ and $\mathrm{Im}[z]$ via Eq. (70) and work only with $\mathrm{Re}[w]$ and $\mathrm{Re}[z]$. It turns out that a surprisingly large number of exceptional phase space points and trajectories cause problems for particular basis functions in the $(w,z)$ representation. Specifically, we find that certain ${\rm Li}_{n}$ arguments assume values on the real axis greater than one

• •

when  $s=3\,m^{2}$  for  $-3\,m^{2}<t<0$.

• •

when  $s=-\frac{4\,m^{2}t}{m^{2}-t}$  for  $-3\,m^{2}<t<-m^{2}$.

• •

when  $s=-\frac{2\,m^{2}(-2\,t-m^{2})}{t}$  for  $-\left(2+\sqrt{2}\right)m^{2}<t<-m^{2}$.

• •

when  $s=\frac{4\,m^{2}+t+\sqrt{-t\,(8\,m^{2}-t)}}{2}$  for  $-\left(2+\sqrt{2}\right)m^{2}<t<0$.

We also observe that a number of function arguments have imaginary parts which vanish identically below threshold. Any such ${\rm Li}_{n}$ function with an argument which can attain values greater than one on $\mathcal{D}_{b_{1}}$ must be discarded because it will be ill-defined unless small imaginary parts are assigned to $w$ and $z$.

In the end, we find the 178 preferred ${\rm Li}_{n}$ function arguments

[TABLE]

for $\mathcal{D}_{b_{1}}$. In an analogous manner, we find 3314 preferred pairs of ${\rm Li}_{2,1}$, ${\rm Li}_{3,1}$, and ${\rm Li}_{2,2}$ function arguments $(f_{i},f_{j})$ after using Reduce to throw away all pairs for which either $\mathrm{Re}[f_{i}]>1$ and $\mathrm{Im}[f_{i}]=0$ or, alternatively, $\mathrm{Re}[f_{i}f_{j}]>1$ and $\mathrm{Im}[f_{i}f_{j}]=0$ for $(w,z)$ in $\mathcal{D}_{b_{1}}$.

For regions $\mathcal{D}_{b_{2}}$ and $\mathcal{D}_{b_{3}}$, the first-entry conditions lead to the preliminary list of 194 ${\rm Li}_{n}$ function arguments

[TABLE]

The rest of the filtering procedure in the remaining regions is closely analogous to what is described above for $\mathcal{D}_{b_{1}}$, but we nevertheless summarize our findings for the benefit of the reader.

In $\mathcal{D}_{b_{2}}$, besides those functions whose arguments have imaginary parts which vanish identically, we find that certain ${\rm Li}_{n}$ arguments assume values on the real axis greater than one

• •

when  $s=-\frac{4\,m^{2}t}{m^{2}-t}$  for  $-m^{2}<t<-\frac{m^{2}}{3}$.

• •

when  $s=-\frac{2\,m^{2}(-2\,t-m^{2})}{t}$  for  $-m^{2}<t<-\frac{2\,m^{2}}{3}$.

• •

when  $s=\frac{4\,m^{2}+t-\sqrt{-t\,(8\,m^{2}-t)}}{2}$  for  $-\frac{m^{2}}{3}<t<0$.

After discarding all of the offending ${\rm Li}_{n}$ function arguments, we find that 184 preferred ones remain:

[TABLE]

In much the same way, we find 3420 preferred pairs of function arguments for ${\rm Li}_{2,1}$, ${\rm Li}_{3,1}$, and ${\rm Li}_{2,2}$ in $\mathcal{D}_{b_{2}}$.

In $\mathcal{D}_{b_{3}}$, besides those functions whose arguments have imaginary parts which vanish identically, we find that certain ${\rm Li}_{n}$ arguments assume values on the real axis greater than one

• •

when  $s=-\frac{4\,m^{2}t}{m^{2}-t}$  for  $-\frac{m^{2}}{3}<t<0$.

• •

when  $s=-\frac{2\,m^{2}(-2\,t-m^{2})}{t}$  for  $-\frac{2\,m^{2}}{3}<t<-\left(2-\sqrt{2}\right)m^{2}$.

• •

when  $s=\frac{4\,m^{2}+t-\sqrt{-t\,(8\,m^{2}-t)}}{2}$  for  $-\left(2-\sqrt{2}\right)m^{2}<t<-\frac{m^{2}}{3}$.

After filtering, we find that 184 preferred ${\rm Li}_{n}$ function arguments remain:

[TABLE]

Note that (5.2) and (5.2) are not identical even though they are of identical length. As before, we find 3420 preferred pairs of function arguments for ${\rm Li}_{2,1}$, ${\rm Li}_{3,1}$, and ${\rm Li}_{2,2}$ in $\mathcal{D}_{b_{3}}$.

Although it is not a priori clear that such restrictive choices for the ansätze of Li functions in each region are actually allowed, it turns out that we can in fact use the sets of function arguments derived from the considerations outlined above to write explicit solutions which satisfy the differential equations for the twelve two-massive-line master integrals with rational symbols in $\mathcal{D}_{a}$, $\mathcal{D}_{b_{1}}$, $\mathcal{D}_{b_{2}}$, and $\mathcal{D}_{b_{3}}$. We have explicitly checked numerically against the $G$ function representation with kinematic variables in the arguments at a large number of physical phase space points that, as desired, our results are correct independent of whether explicit positive imaginary parts are assigned to $w$ and $z$. In the next section, we will discuss the more complicated master integrals which involve the intrinsically algebraic part of the full symbol alphabet, $\mathcal{L}_{A}$, in a similar manner.

6 Weight-four multiple polylogarithms for Drell-Yan master integrals

In this section, we apply the techniques introduced in Sections 4 and 5.2 to explicitly integrate the $\epsilon\,\mathrm{d}\ln$ form for $\{\mathbf{m}_{1},\ldots,\mathbf{m}_{36}\}$, (A), in the physical region above the two-mass threshold ($s>4\,m^{2}$) through to weight four in terms of real-valued Li functions. To our knowledge, this is the first time a complete solution in terms of standard multiple polylogarithms has been found for a multi-loop Feynman integral with an $\epsilon\,\mathrm{d}\ln$ differential equation containing unrationalizable symbol letters. We fix our integration constants using a variety of established techniques, such as direct evaluation in Feynman parameters, exploiting regular kinematic limits of the differential equations, and exploiting kinematic limits where particular integrals are power-suppressed and must therefore vanish (see Section 5.2 of Bonciani:2016ypc ).

Due to the significant length of our results, we refrain from presenting them explicitly in the text. For the benefit of the reader interested in the detailed structure of our solution, we provide the result through to weight four for the six-line normal form integral from the right panel of Figure 2, $\mathbf{m}_{32}$, as an ancillary file included with our arXiv submission. $\mathbf{m}_{32}$ has a number of notable analytic features. We find a basis of multiple polylogarithms comprised of: 14 $\ln$ functions (by construction, $\ln(l_{10})$ and $\ln(l_{11})$ cannot appear), 46 ${\rm Li}_{2}$ functions, 235 ${\rm Li}_{3}$ functions, 342 ${\rm Li}_{4}$ functions, 28 ${\rm Li}_{2,1}$ functions, 742 ${\rm Li}_{3,1}$ functions, and 324 ${\rm Li}_{2,2}$ functions. To give the reader a feeling for the complexity of our result for $\mathbf{m}_{32}$, we present the complete list of 347 ${\rm Li}_{n}$ function arguments which appear in it:

[TABLE]

It is also worth pointing out that our choice of alphabet, Eq. (4), produces remarkably simple integration constants. We find

[TABLE]

at weight three and

[TABLE]

at weight four.

In lieu of our explicit analytic results, it is straightforward to obtain high-precision numerical results for our master integrals using GiNaC. For the randomly chosen physical phase space point

[TABLE]

we find for the most complicated master integrals with two massive internal lines

[TABLE]

Using the finite integral method vonManteuffel:2014qoa ; vonManteuffel:2015gxa ; vonManteuffel:2017myy and SecDec 3 Borowka:2015mxa , we were able to independently check these numerical results to a few decimal digits.

We find that a double precision evaluation of all 36 master integrals ${m_{1},\ldots,m_{36}}$ at the point (80) using GiNaC takes $0.5$ seconds on one Intel E3-1275 CPU core. The performance in this central point in phase space is therefore not much worse than what is encountered for Feynman integrals involving rational alphabets, see e.g. the discussion in Gehrmann:2015ora . Therefore, we expect the functional representation discussed in this paper to be well-suited for direct usage in Monte-Carlo programs.

7 Outlook

In this article, we considered irreducible Feynman integrals satisfying $\epsilon\,\mathrm{d}\ln$ differential equations with non-rational symbol letters and described methods for their evaluation in terms of standard multiple polylogarithms. We show for the first time that a complete solution in terms of standard multiple polylogarithms can be obtained even in the presence of unrationalizable symbol letters. In particular, new techniques for the construction of an ansatz which matches the symbol in a particular region of phase space allowed us to calculate the two-loop master integrals for the mixed EW-QCD corrections to Drell-Yan lepton pair production in the physical region using ${\rm Li}$ functions through to weight four. We discussed in detail how to optimize the functional basis to allow for fast and stable numerical evaluations, systematically avoiding arguments on branch cuts which would require an explicit $+i\,0$ prescription. The presented techniques and results may be applied directly to the calculation of the virtual amplitudes for the full EW-QCD corrections to Drell-Yan production of relative order $\alpha\alpha_{s}$.

As motivation, we also considered a rather non-trivial master integral for massive Bhabha electron-positron scattering and showed that it may be directly integrated from Feynman parameters to all orders in $\epsilon$ in terms of multiple polylogarithms despite the presence of root-valued symbol letters. Besides the applications discussed in this paper, $\epsilon\,\mathrm{d}\ln$ differential equations with root-valued symbol letters also appear in other interesting contexts, such as the next-to-leading order QCD corrections to $H$+jet production with full heavy top quark mass dependence. One may therefore expect that our techniques could be fruitfully applied also to other problems of current phenomenological interest.

Acknowledgements.

We gratefully acknowledge Erik Panzer for inspiring and stimulating discussions over the last several years about various aspects of this work. RMS is also very grateful to Francis Brown, Erik Panzer, and All Souls College of Oxford University for their hospitality at an early stage of this work. We thank Stefan Müller-Stach, Duco van Straten, and Marco Besier for many enlightening discussions about the theory of rational parametrization. We would also like to thank Hubert Spiesberger for his support and for collaborations on related topics. We thank Roberto Bonciani, Stefano Di Vita, Pierpaolo Mastrolia, and Ulrich Schubert for a comparison of results. MH was supported in part by the German Research Foundation (DFG), through the Collaborative Research Center, Project ID 204404729, SFB 1044, and the Cluster of Excellence PRISMA+, Project ID 39083149, EXC 2118/1. AvM was supported in part by the National Science Foundation under Grant No. 1719863. The authors would like to express a special thanks to the Mainz Institute for Theoretical Physics (MITP) of the Cluster of Excellence PRISMA+ for its hospitality and support. Our figures were generated using Jaxodraw Binosi:2003yf , based on AxoDraw Vermaseren:1994je .

Appendix A Differential equations for the Drell-Yan masters with two massive lines

In terms of the symbol letters $\{l_{1},\ldots,l_{16}\}$ defined in (4), the system of differential equations for $\mathbf{m}_{1},\ldots,\mathbf{m}_{36}$ in (3) can be cleanly written in differential form:

[TABLE]

Appendix B A symbol alphabet with several independent roots

In this appendix, we consider the symbol alphabet of the one-loop box integral shown in Figure 6. This integral occurs in processes like $gg\to hh$ or $gg\to ZZ$ with a top quark running in the loop.

We use the kinematical definitions

[TABLE]

and $m_{t}^{2}$ for the internal mass squared.

It is straightforward to write down a normal form basis,

[TABLE]

which fulfils a differential equation of the form

[TABLE]

with $i,j=1,\ldots,8$. It is not straightforward, however, to actually arrive at (96) with a simple form for the letters $l_{k}$ due to the presence of the five different square roots. Moreover, while it is is not difficult to choose a set of independent letters, it is not obvious how to minimize the number of independent letters.

While one could try to rationalize some of the roots with a suitable reparametrization, we want to take a different approach here and just directly work with the original kinematic invariants. Our primary motivation is to extend the construction of symbol letters presented in Section 4 to the case of several root-valued leading singularities.

Employing the integration measure

[TABLE]

the integrals $m_{i}$ depend only on the dimensionless ratios of scales

[TABLE]

In a conventional approach, one might try to directly integrate the terms in the partial differential equations with e.g.  Mathematica to arrive at the form (96). In this approach, we obtained 24 independent letters, many of which are very complicated and involve high powers of the variables.

Instead, the following construction leads to much simpler and fewer independent letters. Starting from the rational part:

[TABLE]

and the five square roots

[TABLE]

we used the algorithm described in Section 4 for each of these square roots independently. Next, we formed pairs and triples of the roots and used the algorithm on their products. We find a total of 22 independent candidates for letters in this way, which are much simpler than the original set of letters and involve only low powers of the variables:

[TABLE]

The factor $1/2$ appears for the same reason as in Section 4 to prevent any factor of $2$ in the factorization of the product of a letter with its conjugate over the rational part. All of these letters are power products of the 24 letters in the differential equation obtained in the first attempt.

To go to even lower degrees, one can do a further factorization of the letters by factorizing the square roots further, i.e. by splitting them in a similar way to:

[TABLE]

and dividing them out in all letters, if possible. Note, that this leads to a representation with a possible sum of different square roots in each letter. If one would like to derive this representation along the lines described in Section 4, one would need to allow for more general structures than

[TABLE]

with rational $q_{a}$.

Interestingly, in our new representation even more letters drop out of the differential equation, and one finds only $19$ letters in the end:

[TABLE]

with the differential equation given by:

[TABLE]