Two-Loop master integrals for heavy-to-light form factors of two different massive fermions

We calculate the full set of the two-loop master integrals for heavy-to-light form factors of two different massive fermions for arbitrary momentum transfer in NNLO QCD or QED corrections. These integrals allow to determine the two-loop QCD or QED corrections to the amplitudes for heavy-to-light form factors of two massive fermions in a full analytical way, without any approximations. The analytical results of the master integrals are derived using the method of differential equations, along with a proper choosing of canonical basis for the master integrals. All the results of master integrals are expressed in terms of Goncharov polylogarithms.


Introduction
The calculation of heavy-to-light form factors have a number of applications in particle physics. For instance, the semileptonic decay of a heavy quark to a light massive quark (t → b + W + (l +ν), b → c + l +ν), the decay of a massive lepton to anther massive lepton (µ → e+ν µ +ν e , τ → µ+ν τ +ν µ ) , and the decay of W bosons into two massive quarks. The NNLO QCD or QED corrections to the decay of a heavy fermion to a light fermion were computed in a number of papers [1][2][3][4][5][6][7][8][9], and the masses of the light fermions were neglected to simplify the calculation. In Ref. [10][11][12], the semileptonic decay b → c + l +ν have been calculated with the mass of charm quark taken into account, an expansion in powers and logarithms of mass ratio mc m b was performed there. It has been shown that in the calculation of similar NNLO corrections to light fermion energy spectrum in heavy fermion decay, the mass of the light fermion can not be neglected since there exist the large logarithm terms ln( m heavy m light ) [13]. These large logarithm terms cancel out in the calculation of the total rate which make the calculation simpler. At order O(α 2 ) or O(α 2 s ), double-logarithmic ln 2 ( m heavy m light ) and single-logarithmic ln( m heavy m light ) enhanced terms will appear, which make it impossible to compute the radiative corrections to quantities such as the light fermion energy spectrum by neglecting the mass of the light fermion from the very beginning. Having the masses as regulators can simplify the treatment of real emission processes, however the computation of virtual corrections is more complicate compared to a purely massless case. The full dependence on the muon and electron masses for electron energy spectrum in muon decay at NNLO QED corrections have first been calculated in a numerical way in [14]. However, the appearance of m 2 light m 2 heavy ≪ 1 will cause numerical instability for the numerical evaluation of multi-dimension integrals, which make the analytical calculations highly desirable. Moreover, the analytic expressions are obviously desirable in order to have control over errors in approximations, especially when the convergence of expansion in powers of mass ratio works slowly.
On the theoretical side, unravelling the mathematical structure of Feynman integrals will be important to handle the complexity of their calculation and may help us to obtain a better understanding of the perturbative quantum field theory. The study of the mathematical properties of Feynman integrals has attracted increasing attention both by the physics and the mathematics communities. Significant progresses were achieved in understanding the analytical computations of multi-loop Feynman integrals in the last years .
The two-loop master integrals of QED electron form factors have been calculated in [37], for on shell electrons of finite squared mass and arbitrary momentum transfer. The calculations of those master integrals have been refined by a suitable choice of basis [23]. The analytical results of two-loop master integrals for form factors of heavy fermion to massless fermion have been obtained in [6, 38? ? ]. All the master integrals for these two processes can be expressed in terms of Harmonic polylogarithms. However, for twoloop vertex integrals with two different type of massive fermions, the master integrals will contain one more scale, to the best of our knowledge, the integrals have not been calculated analytically in the literature. Furthermore, understanding the structure of loop integrals more generally is an interesting and important challenge. In this work, employing the method of differential equations, along with a proper choice of canonical basis, we calculate all the master integrals for two-loop heavy-to-light form factors of two massive fermions, the results are expressed in terms of Goncharov polylogarithms.
The outline of the paper is as follows. In section 2, the kinematics and notations are introduced for the processes we concern. We also present the generic form of the differential equations with respect to the kinematics variables in terms of the derivatives of the external momentum. In section 3, the Goncharov polylogarithms as well as Harmonic polylogarithms are introduced. In section 4, the canonical basis are explicitly presented, followed by the discussion of their solutions. In sections 5, the determination of the boundary conditions, as well as the analytical continuation are explained. Discussions and conclusions are made in section 6. In Appendix A, we present all the rational matrices of the system of differential equations in canonical form. All the analytical results up to weight four from our computation as well as the rational matrices in electronic form are collected in ancillary files accompanying the arXiv version of this publication.

Notations and Kinematics
We consider the process of a heavy fermion decay into a light massive fermion, or the decay of W boson into two massive fermions, with p 2 1 = m 2 1 and p 2 2 = m 2 2 . For the decay of a heavy fermion to a light massive fermion, the squared momentum transfer have the following relations While for the decay of W boson into two different massive fermions, we have We take the dimensionless variables x and y to express the analytical results, they are defined by x , and m 2 = m 1 y. (2.5) In the above equations m 1 is treated as constants, s and m 2 are considered as variable. The derivatives of s and m 2 2 can be written in terms of the derivatives of the external momenta and expressed as (2.6) The corrections to the processes (2.1) and (2.2) could be calculated using Feynman diagram approach. All the amplitudes can be expressed in terms of a set of 40 scalar integrals. The calculation of these scalar integrals always turns out to be the toughest parts in the whole work. We use packages FIRE [39][40][41] to reduce the group of scalar integrals into a minimum set of independent master integrals. FIRE is also adopted in the derivations of differential equations for master integrals.

Goncharov polylogarithms and Harmonic polylogarithms
The Goncharov polylogarithms (GPLs) [42] we use to express the analytical results are defined as follow They can be viewed as a special case belonging to a more general type of functions named Chen-iterated integrals [43]. If all the index a i belong to the set {0, ±1}, the Goncharov polylogarithms turn into the well-known Harmonic polylogarithms (HPLs) [44] H a 1 ,a 2 ,...,an (x) = (−1) k G a 1 ,a 2 ,...,an (x), (3.4) where k equals to the times of element (+1) taken in (a 1 , a 2 , . . . , a n ) . The Goncharov polylogarithms fulfil the following shuffle rules Here, aXb is composed of the shuffle products of lists a and b. It is defined as the set of the lists containing all the elements of a and b, with the ordering of the elements of a and b preserved. Both the GPLs and HPLs can be numerically evaluated within the GINAC implementation [45,46]. A Mathematica package HPL [47,48] is also available to reduce and evaluate the HPLs. Both the GPLs and HPLs can be transformed to the functions of ln, Li n and Li 22 up to weight four, with the algorithms and packages described in [49].

The canonical basis
All the amplitudes of two-loop QCD or QED corrections for the processes we concern can be reduced to a set of 40 master integrals, including two non-planar integrals. The master integrals M i (i = 1 . . . 40) are shown in Fig. 1. The vector of canonical basis F is built up with 40 functions F i (s, m 1 , m 2 , ǫ)(i = 1 . . . 40), defined in terms of the linear combinations of master integrals M i .
F 10 = ǫ 2 s M 10 , (4.10) All the NNLO master integrals for heavy-to-light form factors with different type of fermion masses. The thin lines denote massless propagators; the thick dashed lines denote massive fermion with mass m 1 ; and the thick solid lines denote massive fermion with mass m 2 . A dot on a propagator indicates that the power of the propagator is raised to 2. Two dots represents that the propagator is raised to power 3.
, (4.14) (4.28)  For our choice of master integrals above, the differential equations for F = (F 1 . . . F 34 ) have the following canonical form with A(x, y) = A 1 ln(x) + A 2 ln(x + 1) + A 3 ln(x − 1) + A 4 ln(x + y) + A 5 ln(x − y) +A 6 ln(x y + 1) + A 7 ln(x y − 1) + A 8 ln(y) + A 9 ln(y + 1) + A 10 ln(y − 1) +A 11 ln(x 2 y − 2x + y) + A 12 ln(x 2 − 2y x + 1). (4.42) The notations A i (i = 1 · · · 12) are 40×40 matrices with rational numbers, they are presented in appendix. A. For reader's convenience, the rational matrices in electronic form are also provided in ancillary file accompanying the arXiv version of this paper. We can see in equation (4.42) that in equal masses case with y = 1, the alphabet turns into {x, x+1, x−1}, and the results of all the canonical basis can simply be reexpressed in terms of Harmonic polylogarithms.
The integral M 1 is defined as follow with the measure of the integration is defined as For master integrals without numerators, their definition can readily be read off from Fig.  1, with the normalization defined above. For master integrals with numerator, we first define a series of propagators (4.45) Then, the master integral with numerator (M 34 ) can be expressed as (4.46)

Boundary conditions and analytic continuation
Now, we are ready to perform the calculations of the differential equations. The first step is to specify all the boundary conditions that will completely fix the solutions of the differential equations. The results of basis (F 1 . . . F 5 ) and (F 35 , F 36 ) are already known in the literature [25,[50][51][52]. They can be recalculated with the assistance of Mathematica packages MB [53] and AMBRE [54][55][56]. The integrals (F 6 . . . F 9 ) have been calculated in [50,51], their boundary conditions could be determined by setting y = m 2 m 1 = 1, with the masses of two fermions equaling to each other, and their boundaries are proportional to the results of F 35 .
The integral M 10 does not have singularity at s = 0. Thanks to the normalization factor s that multiplying with M 10 in F 10 , we can readily know that F 10 = 0 at s = 0. Considering the fact that M (11,13,14,16...20,23,24,26,27,29...34,37,39) are regular at s = (m 1 − m 2 ) 2 as well as s = (m 1 +m 2 ) 2 , and their normalization factor to be s , the boundaries of basis F (11,13,14,16...20,23,24,26,27,29...34,37,39) are 0 at s = (m 1 − m 2 ) 2 and s = (m 1 + m 2 ) 2 . The boundary of F 15 at s = (m 1 − m 2 ) 2 and s = (m 1 + m 2 ) 2 can also be determined similarly and then expressed as Considering the fact that M 12 does not have singularity at s = 0, i.e. (x = y) or (x = 1 y ) , the boundary condition of F 12 can be determined from the differential equation of F 12 . To further illustrate it, we consider the differential of F 12 with respect to variable x and find that Both F 11 and F 12 have finite limit at x = y and x = 1 y , this consistency leads to two relations between F 11 and F 12 The boundary of (F 21 , F 22 , F 25 , F 28 , F 38 , F 40 ) can also be determined with the discussion above. By now, all the boundary conditions are fixed.
We then proceed to determinate the analytic continuation of all the master integrals. When considering the processes (2.1,2.2), the variables of the master integrals lie either in Euclidean region or in Minkowski region. In order to obtain the correct numerical results, their analytic continuation should be considered carefully. The proper analytic continuation can be achieved by the replacement s → s+i0 at fixed m 2 1 and m 2 2 . This transfer corresponds to x → x + i0.
After the determination of the boundary conditions and analytical continuation, we can readily obtain the analytical results of all the master integrals using the differential equations we obtain. The results of the master integrals can be written in terms of Goncharov polylogarithms introduced in section 3. The results for F i (i = 1 . . . 40) are calculated up to weight four in this work. All the analytic results up to weight four are collected in the ancillary file "results.m" which is supplied with the paper. Here, for illustration, we present the results for integrals (F 6 , F 10 , F 33 ) up to weight three.
Note that the weight three results of F 33 depend only on variable x, the dependence of F 33 on variable y will start at O(ǫ 4 ).
The calculations are performed with our in house Mathematica code. It is desirable to perform an independent check on the analytical expressions of the master integrals. We check all the analytical results against the numerical results obtained from packages Fiesta [57,58] and SecDec [59,60]. Perfect agreement has been achieved between the analytical and numerical approaches with kinematics in both Euclidean region and Minkowski region. To be more specific, we show the numerical results of the master integrals M 33 (s, m 1 , m 2 ) obtained from packages Fiesta and SecDec and our analytical expressions in two different kinematics, We can see from the equations above that the numerical results obtained from packages Fiesta and SecDec with error estimates are perfectly agree to that precision with the numerical results obtained from our analytical answer. Note that it took Fiesta and SecDec packages each more than ten hours on a workstation with 16 cores (2.3 GHz) processor to reach the above precision for kinematic (s = 5.4, m 1 = 1.0, m 2 = 0.2) in Minkowski region.

Discussions and Conclusions
In summary, applying the method of differential equations, we calculate the full set of the two-loop master integrals for heavy-to-light form factors of two different types of massive fermions and arbitrary momentum transfer in NNLO QCD or QED corrections. With a proper choice of master integrals, it turns out that we can cast all the differential equations into the canonical form, which can straightforwardly be integrated order by order in dimensional regularization parameter ǫ. Under the determination of all boundary conditions, we express the results of all the basis in terms of Goncharov polylogarithms. The integrals allow to determine the two-loop amplitudes for heavy-to-light form factors of two massive fermions in a full analytical way, and thus to compute the NNLO corrections to the decay of heavy massive fermion into light massive fermion, or the decay of W boson into two massive quarks without any approximations. The integrals are also applicable to the calculations of NNLO QCD corrections to the inclusive decay of B meson into charmonium, and two-loop form factors of H + tb coupling vertex.