Three-loop vacuum integrals with arbitrary masses

Three-loop vacuum integrals are an important building block for the calculation of a wide range of three-loop corrections. Until now, analytical results for integrals with only one and two independent mass scales are known, but in the electroweak Standard Model and many extensions thereof, one often encounters more mass scales of comparable magnitude. For this reason, a numerical approach for the evaluation of three-loop vacuum integrals with arbitrary mass pattern is proposed here. Concretely, one can identify a basic set of three master integral topologies. With the help of dispersion relations, each of these can be transformed into one-dimensional or, for the most complicated case, two-dimensional integrals in terms of elementary functions, which are suitable for efficient numerical integration.


Introduction
The need for higher-order radiative corrections is growing more and more important due to the increasing precision of measurements at the LHC and planned future colliders. The anticipated precision of future experiments will require the evaluation of three-loop corrections with arbitrary masses, e.g. for precision electroweak quantities (for a recent review see ref. [1]) or a detailed understanding of the Higgs potential and its stability [2]. In this article, the calculation of general three-loop vacuum integrals is considered, i.e. integrals with vanishing external momentum and arbitrary propagator masses. Such integrals may arise in low-energy observables, in the coefficients of low-momentum expansions (see e.g. ref. [3]) or as building blocks in more general three-loop calculations. At the two-loop level, analytical formulae for general vacuum integrals have been known for some time [3][4][5]. When expanding in powers of ǫ = (4−D)/2 within dimensional JHEP11(2016)145 regularization, they can be written in terms of polylogarithms. At the three-loop level, results for vacuum integrals are only available for one [6-11, 13, 14] and two [13][14][15][16][17][18][19] 1 independent mass scales. The derivation of analytical results for the class of two-scale three-loop vacuum integrals requires the introduction of harmonic polylogarithms [21], and some cases are only known numerically [13,14].
In light of these facts, a numerical approach to three-loop vacuum integrals with general mass pattern appears most promising. In ref. [22][23][24] a numerical technique for the calculation of the four-propagator topology has been presented. In the following, a method for the evaluation of all relevant master integrals is proposed, which is based on dispersion relations. This technique has been previously used for the numerical evaluation of two-loop self-energy and vertex integrals [25][26][27]. For the master integrals considered in this paper, the dispersion relation approach leads to simple numerical integrals for their finite part. For two of the three master integral topologies, one can obtain one-dimensional numerical integral representations in terms of elementary functions. For the most complicated case, the six-propagator master integral, one may construct a two-dimensional integral in terms of elementary functions. 2 Note that in some applications it may be necessary to evaluate the master integrals to higher orders in ǫ. This happens when a master integral is multiplied by a coefficient that has poles in 1/ǫ. The method described in this paper, in its present form, is not suitable for such situations. This article begins by defining the set of three-loop vacuum master integrals in section 2. Each master integral topology is discussed in turn in sections 3-5. Several special cases, which require a modification of the integral representations, are treated separately in sections 3 and 4. For the most complicated master integral, which is the subject of section 5, no such special case has been identified so far. The paper finishes with some comments on the implementation of the numerical integrations in section 6 before concluding in section 7. Some useful formulae are collected in the appendix.

Definition of basic integrals
After trivial cancellations of numerator and denominator terms, a general scalar three-loop vacuum integral may be written in the form small set of master integrals with the help of integration-by-parts identities [30]. 3 In most cases, not involving any special mass patterns, one can choose the following basis of three master integrals, see figure 1, besides integrals that factorize into products of one-and two-loop contributions. Two other simple integrals that are often encountered (see figure 2) can be reduced to these three with the help of integration-by-parts identities:  where B 0 is the standard scalar one-loop self-energy function (see appendix C), and B 0,m 1 (p 2 , m 2 1 , m 2 2 ) = ∂ ∂m 2 1 B 0 (p 2 , m 2 1 , m 2 2 ). The discontinuities of these two functions are denoted by ∆B 0 and ∆B 0,m 1 , respectively. In D = 4 dimensions, they are given by where λ(x, y, z) = x 2 + y 2 + z 2 − 2(xy + yz + xz) , (3.5) and Θ(t) is the Heaviside step function.

JHEP11(2016)145
Inserting the dispersion relation eq. (3.1) into the three-loop integral U 4 , one obtains This integral is divergent, and thus a numerical integration of (3.7) in D = 4 dimension is not possible. Instead one may consider the sum The U 4 master integrals with one or two non-zero masses can be calculated analytically, with results collected in appendix B. The first four terms on the right-hand side of eq. (3.8) precisely reproduce the divergencies of the general U 4 integral on the left-hand side, such that the remainder U 4,sub is finite and can be integrated numerically. It is given by which has been written in a way that makes its finiteness manifest. The divergent part of the A 0 function in eq. (3.9) integrates to zero and thus can be ignored.

Special case: m 1 = 0
A special treatment is required for the case m 1 = 0, when U 4 (m 2 1 , m 2 2 , m 2 3 , m 2 4 ) develops an infrared divergence. This singularity occurs in the limit when the first loop momentum goes on-shell, q 2 1 → 0. In this limit U 4 (0, m 2 2 , m 2 3 , m 2 4 ) factorizes into the product B 0 (0, 0, 0) T 3 (m 2 2 , m 2 3 , m 2 4 ). Note that B 0 (0, 0, 0) = 0 in dimensional regularization, but it has been kept here for illustration. By alternatively considering a mass regulator for the infrared divergence, one arrives at the identity where T 3 is the basic scalar two-loop vacuum integral (see appendix C), and δ is an infinitesimal mass parameter. In principle, one could directly evaluate U 4 (δ 2 , m 2 2 , m 2 3 , m 2 4 ) for a small numerical value of δ according to the prescription in the previous subsection, although one has to pay the price of having numerical cancellations between the two log δ terms from the last two terms in eq. (3.11).

JHEP11(2016)145
Alternatively, it is possible to extract the log δ dependence explicitly from U 4 (δ 2 , m 2 2 , m 2 3 , m 2 4 ). For this purpose, let us consider the following small δ expansions: where f (s) is some arbitrary well-behaved function that does not depend on δ. For the remaining terms in the integrand, involving B 0 and ∆B 0 functions, one can simply set δ to zero. One then obtains

General case
The master integral U 5 can be addressed with a dispersion relation similar to eq. (3.1). As in the previous section, one first needs to subtract the divergencies to arrive at a finite integral suitable for numerical evaluation. For this purpose, it is useful to consider the following relation, which has been derived from integration-by-parts identities: where F [. . .] is a linear combination of U 4 functions and products of A 0 and T 3 functions, whose explicit form is too lengthy to be included here. It is provided as an ancillary file in Mathematica format with this article. Furthermore, Following the approach of section 3.1, M (2, 1, 1, 2, 1, 0) can be expressed in terms of a dispersion integral. The relevant dispersion relation reads   The formula (4.2) is not valid for the case m 2 1 = m 2 2 + m 2 5 . 4 As long as m 5 > 0, this problem can be avoided by simply flipping m 1 ↔ m 2 . A special treatment is needed, 4 The similar case m 2 4 = m 2 3 + m 2 5 can be mapped to the former by flipping the two sub-bubbles of U5, i.e. {m1 ↔ m3, m2 ↔ m4}.

JHEP11(2016)145
however, for the case m 5 = 0 and m 1 = m 2 . In this case, instead of eq. (4.2) one can use the modified formula (2, 1, 1, 2, 1, 0) . (4.6) The full expression for the terms involving A 0 , T 3 and U 4 functions is again provided in the ancillary file. The finite remainder M (2, 1, 1, 2, 1, 0) can be evaluated as above. While the integral U 5 is infrared finite for m 1 = m 2 = 0, the application of the procedure in section 4.1 would lead to spurious divergencies in several terms in eq. (4.2) that only cancel in the sum. Instead, the following integration-by-parts relation proves more useful for this situation:

U 6
Without loss of generality, it is assumed in the following that m 6 ≥ m i for i = 1, . . . , 5. Following ref. [9], one can write the U 6 master integral in terms of where ∆T 5 is the discontinuity of the two-loop five-propagator self-energy master integral shown in figure 4. Evidently, eq. (5.1) is divergent, but one can arrive at a finite integral by considering the difference between the general U 6 integral and a simpler U 6 with only one independent mass scale: The special case of U 6 with identical masses in all six propagators is known analytically [9], where ζ(x) is the Riemann zeta function, and Cl 2 is the second Clausen function. The discontinuity ∆T 5 can be written as a sum over the possible two-and three-particle cuts [27,32], where ∆C 0 (s, p 2 2 , p 2 In general, eq. (5.5) also contains a contribution from an anomalous threshold [27], but for the mass ordering m 6 ≥ m i it is never encountered. Combining eqs. (5.2)-(5.11) one arrives at a two-dimensional integral representation for U 6 in terms of elementary functions, such a logarithms and square roots.
Alternatively, the two-particle cut contributions (2a) and (2b) can also be written directly in terms of the tree-point function C 0 , Explicit formulae for the C 0 function in terms of polylogarithms can be found in ref. [34][35][36]. Therefore, by using eq. (5.12) one arrives at one-dimensional numerical integrals for the (2a) and (2b) terms, but at the cost of requiring non-elementary functions in the integrand. The numerical dispersion relation based on eq. (5.2) proved to be applicable to a large variety of different mass configurations, without the need to consider special cases as in the previous two sections.

Numerical integration and checks
The dispersion relation techniques introduced in the previous sections lead to one-and two-dimensional integrals, which can be efficiently evaluated with numerical methods. In many cases, one encounters integrals of the form For s ′ > s 0 , these can be split into a residuum contribution and principal value integral, resulting in In this expression, the integrand of the first remaining integral is now regular at the point s = s ′ . If f (s) is real, the residuum contribution can be dropped since only the real part is needed for the evaluation of the three-loop vacuum integrals.
Results have been tested both within Mathematica 10 [37] and in C++ using the Gauss-Kronrod algorithm from the Quadpack library [38]. The following checks were carried out for the U 4 and U 5 master integrals: at least ten digits agreement were obtained JHEP11(2016)145 when comparing with the numbers in ref. [9], the integrals eq. (24) of ref. [10], and the results for the J 10a and J (3) 10b integrals in ref. [13,14]. 5 These correspond to U 4 , U 5 and U 6 integrals with one or two different mass scales. Similarly good agreement was obtained for the J 8b integrals in ref. [13,14], which are not among the set of master integrals in figure 1, but which can easily be derived using eq. (2.5).
A public computer code, currently under development, will be presented in a future publication.

Conclusions
A general three-loop vacuum integral can be reduced to one of the three master integral topologies shown in figure 1. The master integrals can be evaluated analytically for all cases with one independent mass scale and for several known cases with two independent mass scales. For general mass patterns, however, numerical integration techniques need to be employed.
Using a method based on dispersion relations, the four-propagator master integral U 4 and the five-propagator master integral U 5 can be expressed as one-dimensional numerical integrals in terms of elementary functions, such a logarithms and square roots. Similarly, the six-propagator master U 6 can be representated by a two-dimensional numerical integral in terms of elementary functions. These numerical integrals can be efficiently evaluated numerically, yielding results with at least ten digits precision for the U 4 and U 5 functions and eights digit precision for the U 6 function.
To ensure that the numerical integrals are UV-finite, the divergent pieces of the U 4 , U 5 and U 6 must be subtracted beforehand. This can be achieved by subtracting suitable linear combinations of special cases of these integrals, which can be evaluated analytically. A public computer code that carries out these subtractions and the numerical integrations is in development and will be presented in a future publication.
The technique described in this article has been tailored for the evaluation of the UVdivergent and finite parts of the master integrals. In principle, the same dispersion relation approach could be used for the calculation of higher order terms in ǫ, but its concrete realization would require a significant amount of additional work.

A Divergent parts of master integrals
This section presents analytic results for the divergent part of the master integrals in figure 1. where "cycl ijk " refers to cyclic permutations of {m i , m j , m k }, and Furthermore [9,33], where ζ(x) is the Riemann zeta function. For completeness, the divergent terms of the four-propagator integral M (1, 1, 1, 1, 0, 0) are also quoted here:

B Analytic results for some master integrals with one and two massive propagators
Various special cases of the U 4 and U 5 integrals with one and two massive propagators can be computed analytically to all orders in ǫ using Mellin-Barnes representations. The ǫ expansions of the hypergeometric 2 F 1 functions below were performed with the help of the HypExp package [39], which internally utilizes the HPL package [40].  Finally, the two-loop vacuum integral is given by [3] T 3 (m 2 1 , m 2 2 , m 2 x log x + y log y + x log 2 x + y log 2 y − (1 − x − y) log x log y