Three-loop massive tadpoles and polylogarithms through weight six

We evaluate the three-loop massive vacuum bubble diagrams in terms of polylogarithms up to weight six. We also construct the basis of irrational constants being harmonic polylgarithms of arguments ekiπ/3.


Introduction
More than two decades ago, the integration-by-parts relations [1] and asymptotic expansions [2,3] became common in the Feynman diagram calculus. The combination of these methods provides a powerful tool for the evaluation of multiloop diagrams. In particular, massive propagator diagrams through the three-loop order can be reduced with the help of asymptotic-expansion methods to three-loop massive tadpoles, which can be done, e.g., using the FORM [4] package MATAD [5] (see also ref. [6]). 1 There are a lot of physical applications, where the above-mentioned technique was applied. Just to mention but a few examples, it was applied to the evaluation of the three-loop ρ parameter in QCD [7,8] and the electroweak theory [9], the three-loop QCD corrections to heavy-quark production [10], and many other quantities. Integral topologies with all lines massive find applications in calculations of renormalization group functions [11][12][13] at the three-loop order and also at higher orders of the epsilon expansion in four-loop [14] and even five-loop [15] calculations.
In his work [16], Broadhurst noticed that all three-loop single-scale vacuum diagrams at order O((4 − d)/2) in dimensional regularization can be related to the elements of the JHEP08(2017)024 algebra of the sixth root of unity. This observation allowed him to evaluate all the threeloop integrals up to their finite parts in terms of a few constants, being polylogarithms of weight four.
In this paper, we proceed by studying three-loop vacuum integrals with a single mass scale at weights five and six. On the one hand, this is a necessary ingredient in evaluations beyond the three-loop approximation, where the three-loop master integrals have to be expanded to higher powers in d − 4. On the other hand, we would like to test the basis of the algebra of the sixth root of unity through weight six.

Notation
We use dimensional regularization with the dimension of space-time being d = 4 − 2ε in euclidean space. Each loop integration is normalized as follows: where γ = 0.577216 . . . is the Euler-Mascheroni constant. Defining the general three-loop vacuum bubble diagram with six scalar propagators, where it is implied that the masses either take the values m or zero, may be written as In addition to diagrams with six lines, we also have three-loop digrams with five and four lines, 3 Polylogarithms, algebra of the sixth root of unity, and its subalgebras In our study, the key role is played by multiple polylogarithms [17][18][19], defined recursively as repeated integrals, where a 1 , a 2 , . . . , a w and z are complex numbers. The definition in eq. (3.1) is modified in the case of q trailing zero indices in the following way: The integer number w is called the weight of the polylogarithm. The functions G obey the so-called shuffle and stuffle relations. In particular, any product of two G functions with the same argument and weights w 1 and w 2 can be rewritten as a linear combination of G functions of weight w 1 + w 2 . In other words, polylogarithms form a graded algebra.
The algebra of the sixth root of unity A ω is obtained from general polylogarithms by restricting all a j to the seven-letter alphabet {0, ω 0 , ω 1 , ω 2 . . . , ω 5 }, where is a primitive sixth root of unity. At arbitrary argument z, such functions include the so-called inverse-binomial-sums functions [20][21][22][23][24][25] and are related to cyclotomic polylogarithms [26]. At z = 1, they represent a set of irrational 2 constants, which is relevant for the description of some single-scale massive diagrams (in particular, three-loop vacuum bubble and two-loop on-shell self-energy diagrams). The complete basis of the algebra of the sixth root of unity A ω through weight 6 has recently been constructed in ref. [27].
In this work, we construct the basis of the subalgebra of A ω formed by the harmonic polylogarithms [28] H n 1 ...np (z) of arguments z k = ω k . We shall call such an algebra A H(ω k ) . The harmonic polylogarithms are defined similarly to eqs. (3.1)-(3.2), but now the parameters a j can only take the values −1, 0, +1. For historical reasons, there is also a difference in the overall sign. Specifically, the harmonic-polylogarithm functions H n 1 ...np (z) are related to the generalized polylogarithms G via H n 1 n 2 ...nw (z) = (−1) (number of n j = 1) G n 1 n 2 ...nw (z) , (3.4) where n j = −1, 0, +1. Using the scaling properties of the polylogarithms together with the shuffle relations, it is easy to show that any element of the form H n 1 ...np (ω k ) can be rewritten as G ω k 1 ...ω kp (1)
We shall denote the uniform bases of Re{A} H(ω) and Im{A} H(ω) for fixed weight w as Re H w and Im H w , respectively. In the next sections, we apply the constructed bases Re H w and Im H w to the evaluation of the three-loop massive vacuum bubble diagrams.

Evaluation of the three-loop vacuum bubble integrals
Using integration-by-parts relations, it is possible to reduce any three-loop bubble integral with a single scale to a set of twelve three-loop master integrals, two two-loop integrals, and one one-loop bubble. These diagrams are shown in figures 1 and 2. JHEP08(2017)024  It is the goal of this paper to evaluate these master integrals analytically in terms of polylogarithms through weight six.
In the previous section, we discussed the construction of the bases of Re{A} H(ω) and Im{A} H(ω) . We now use these bases to reconstruct the analytic expressions for the ε expansions of these diagrams using the PSLQ algorithm [30]. For that purpose, we first need a precise numerical value of each diagram.
Specifically, for the fully massive diagrams D 6 and E 5 , we make use of the series obtained with the help of the DRA method, based on dimensional recurrence relations and analyticity, presented in ref. [31] and summed with the help of the SummerTime package [32]. Within a few hours, we were able to get 20,000 decimal figures of precision for these diagrams.
The general method of calculation which is used in this work and is applicable to all the considered diagrams consists in writing the systems of differential equations for the integrals and solving them later by the Frobenius method. In the first step, instead of a single-scale diagram with mass m, we introduce a similar diagram with two different masses, m 1 and m 2 . Then, using integration by parts, we can write the system of differential equations in the mass ratio z = m 2 1 /m 2 2 . In general, these equations cannot be solved analytically. We solve them as series of the form n z n c n . The unknown coefficients c n are to be determined by substituting the series in the differential equations.

JHEP08(2017)024
We also need the boundary conditions. The easiest choice in our case is the boundary conditions at z = 0, which correspond to a single-scale bubble integral with a smaller number of massive lines.
Finally, we set z = 1 in the series solution in order to recover the original diagram. The summation of the series is done numerically. In this way, we are able to evaluate integrals to an accuracy of typically 4,000 to 10,000 decimal figures depending on the diagram. For that purpose, we need to sum up to 20,000 terms in the n sum in some cases.
Let us consider as an example the integral D N . There are two massive lines in this diagram. Instead of two equal masses, we introduce now two different masses, one of which we set to unity. Thus, we set in eq. (2.3) m 1 = z, m 5 = 1, and m 2 = m 3 = m 4 = m 6 = 0. With such masses, we have the following set of master integrals, which depend on z: where we use the definition in eq. (2.3). Let us denote, for brevity, the integrals in eq. (4.1) as f 1 , . . . , f 7 in this very order. Then, the functions f 1 (z), . . . , f 7 (z) obey a system of linear differential equations in the variable z, which reads: where f j = df j /dz. To solve the system in eq. (4.2), we substitute the following collective ansatz: The exponent shifts µ k are determined as usual in the Frobenius method from the indicial polynomials. Actually, it is easy to establish that µ j can take the values 0, −ε, −2ε, −3ε. Therefore, we have, for each value of j, four different solutions f j , corresponding to the different values of µ k , and the solution we are looking for is the linear combination

JHEP08(2017)024
with unknown constants C j,k , which should be determined from the boundary conditions at z = 0. Thus, for each value j = 1, . . . , 7, we need four boundary conditions, one for each value of µ k . The boundary conditions correspond to the expansions of the integrals in eq. (4.1) about z = 0. Following the standard rules of the large-mass expansion [2,3], we should take into account the four hard subgraphs {123456}, {23456}, {356} + {245}, and {5}. These four subgraphs provide the four boundary conditions for eq. (4.4).

Results and discussion
We present the terms of the ε expansions analytically in terms of the bases Re H j and Im H j , j = 1, . . . , 6 in the appendix. In each case, we take the prefactor in such a way that the terms of the expansion are homogeneous in the weights. In some cases, this requires us to evaluate additional integrals (with dots on lines) and to re-express the original integral with the help of integration-by-parts relations. Moreover, we find that, with the suitable choice of prefactors, the elements of the expansion are expressed in each case either through the Re H basis or the Im H basis. This feature is a convenient property which allows us to reduce the length of PSLQ vector. In addition to the analytic expressions, we also give their numerical values accurate to 50 decimal figures.
There is, of course, a certain degree of arbitrariness in the choice of the basis elements. We just use the lexicographical ordering of the three-letter alphabet {−1, 0, 1}. The sets of basis elements and the transformation between different bases (with arguments ω and ω 2 ), as well as all analytic results can be found in the attachment in a Mathematica-readable form. In addition, we give the numerical values of all basis elements both for A H(ω) and A H(ω 2 ) to an accuracy of 20,000 decimal figures.

Integrals evaluated in terms of Γ functions
Four of the integrals in figure 1 can be evaluated in terms of Γ functions, namely

Two-loop integral T 111
The two-loop integral T 111 was considered in ref. [33], where its representation in terms of the hypergeometric function 4 F 3 was given. The construction of its ε expansion was JHEP08(2017)024 discussed in great detail in ref. [34]. There, the expansion to all orders in ε was found in terms of log-sine integrals. We find that T 111 can be written in terms of our bases as 111 are expressed in terms of the homogeneous bases Im H k+2 . They are presented in the appendix.
It should be noted here that, in eq. (5.6), it is necessary to take out the factor Γ 2 (1 + ε) to avoid the mixing of the Im H and Re H bases. This mixing occurs, since ζ k ∈ Re{A} H(ω) , while ζ k / ∈ Im{A} H(ω) and √ 3ζ k / ∈ Im{A} H(ω) . To the accuracy of 50 decimal figures, we have

Diagram BN
The integral BN belongs to the class of the so-called 'QED-type' integrals. These are the integrals with an even number of massive lines at each vertex. They have an especially simple structure and were considered in ref. [35] though weight six. We have The diagram BN 1 was previously considered in refs. [36,37]. There, its explicit representation in terms of the hypergeometric function Q F P of argument 1/4 was obtained. The corresponding ε expansion, through weight-five polylogarithms, was constructed in terms of log-sine integrals.
In order to keep the property of the weight homogeneity and to separate the real and imaginary bases, we write BN 1 in terms of the additional integrals BN 1 , which is BN 1 with additional dots and V 1 ,

Diagram E 3
The integral E 3 was also considered in refs. [36,37], where the same analysis as for BN 1 can be found. Our representation of this integral reads: x + . . . , (5.14) where allD (k) x are expressed in terms of the homogeneous bases Re H k+4 and are explicitly given in the appendix.
It should be noted that the integral D 5 is sometimes replaced by the fully massive three-loop integral with 5 lines. The latter integral was considered in ref. [38], where, again, the hypergeometric representation was presented.
Numerically we have

JHEP08(2017)024 6 Conclusions
In this work, we considered the three-loop massive vacuum bubble diagrams and constructed the pertaining bases of irrational constants through weight six, which are harmonic polylogarithms of argument ω = e iπ/3 . These bases are smaller than the bases of the algebra of the sixth root of unity. Nevertheless, we found by explicit calculation that such reduced bases are large enough to describe all the three-loop single-scale vacuum integrals. We presented the results for all relevant master integrals both numerically and analytically in terms of the introduced constants. Our basis is universal, and its application is not restricted to three-loop tadpoles. As an example, we succeeded in reconstructing 3 all the three-loop integrals contributing to the massive planar form factor through weight six presented in ref. [39]. As expected, these integrals only involve the real parts of the basis Re{A} H(ω) .

A Master integrals in terms of harmonic polylogarithms of argument e iπ 3
In this appendix, we present the first few terms of the ε expansions of the relevant master integrals. For the elements of the bases Re H k and Im H k , we introduce the following short-hand notation: Re