Scalar one-loop four-point Feynman integrals with complex internal masses

Based on the method in Refs.~{\tt [D.~Kreimer, Z.\ Phys.\ C {\bf 54} (1992) 667} and {\tt Int.\ J.\ Mod.\ Phys.\ A {\bf 8} (1993) 1797]}, we present analytic results for scalar one-loop four-point Feynman integrals with complex internal masses. The results are not only valid for complex internal masses, but also for real internal mass cases. Different from the traditional approach proposed by G. 't Hooft and M. Veltman in the paper {\tt[Nucl.\ Phys.\ B {\bf 153} (1979) 365]}, this method can be extended to evaluate tensor integrals directly. Therefore, it may open a new approach to cure the inverse Gram determinant problem analytically. We then implement the results into a computer package which is {\tt ONELOOP4PT.CPP}. In numerical checks, one compares the program to {\tt LoopTools version} $2.12$ in both real and complex mass cases. We find a perfect agreement between the results generated from this work and {\tt LoopTools}.


Introduction
The future colliders, like the Large Hadron Collider (LHC) at high luminosities and the International Linear Collider (ILC) [1][2][3], aim to measure the properties of Higgs boson (or to explore the Higgs sector), of top quark and vector bosons, as well as search for Physics Beyond the Standard Model (BSM). These measurements will be performed at high precision, e.g. the Higgs boson's couplings will be measured at the precision of 1% or better for a statistically significant measurement [3]. In order to match the high precision data in the near future, the higher-order corrections from theoretical calculations are necessary. Therefore, the detailed evaluations for one-loop multi-leg and higher-loop at general scale to the selected scattering cross sections at the colliders are urgent requirements.
In traditional framework, the cross-section of the processes at one-loop corrections will be obtained by integrating over the phase space of squared amplitudes which are decomposed into the tensor integrals. These tensor ones are then reduced into scalar one-loop one-, two-, three-and four-point functions. It is well-known that the traditional tensor reductions may meet the inverse Gram determinant problem [4,5] at several kinematic points in the phase space. Consequently, it leads to numerical instabilities. One may apply the suitable experimental cuts to avoid the problem. However, the situation is completely different when we consider one-loop multi-leg processes, for instance, 2 → 5, 6, etc, where the higher-point functions will be reduced to scalar one-loop one-, two-, three-and four-point functions with In numerical checks, one compares this work to LoopTools version 2.12 [22] in both real and complex mass cases. We find a perfect agreement between the results generated from this work and LoopTools.
The layout of the paper is as follows: In section 2, we present the method for evaluating scalar one-loop four-point functions in detail. In section 3, we will show the numerical checks on the program with LoopTools. Conclusions and plans for future work are presented in section 4. Several useful formulae used in this calculation are shown in the Appendix.

The calculation
In this section, we apply the method which was proposed in Refs. [11][12][13] to calculate scalar one-loop four-point functions. The Feynman integrals for these functions are given by D 0 ≡ D 0 (p 2 1 , p 2 2 , p 2 3 , p 2 4 , s, t, m 2 1 , m 2 2 , m 2 3 , m 2 4 ) := d n l P 1 P 2 P 3 P 4 . (1) Where the inverse Feynman propagators are The term iρ is Feynman's prescription. p i and m i for i = 1, 2, 3, 4 are external momenta and internal masses respectively. The external momenta flow incoming as Fig. 1 and follow momentum conservation law p 1 + p 2 + p 3 + p 4 = 0. The loop momentum is l and n is space-time dimension. In this calculation, we are not going to deal with infrared divergence. Thus, we will work directly in space-time dimension n = 4. In general, D 0 is a function of p 2 1 , p 2 2 , p 2 3 , p 2 4 , s, t, m 2 1 , m 2 2 , m 2 3 , m 2 4 with s = (p 1 + p 2 ) 2 , t = (p 2 + p 3 ) 2 . Fig. 1 The box diagrams.
The Feynman integral is then written explicitly as .
Partitioning the integrand into the form 1 P 1 P 2 P 3 P 4 = 1 P 1 (P 2 − P 1 )(P 3 − P 1 )(P 4 − P 1 ) + 1 P 2 (P 1 − P 2 )(P 3 − P 2 )(P 4 − P 2 ) we then make a shift on P k by l i −→ l i + q ki for k = 1, 2, · · · , 4 and i = 0, 1, · · · , 3, the resulting reads 4/28 We have already introduced the following kinematic variables It is important to note that the kinematic variables a lk , b lk , c lk ∈ R and d lk ∈ C. In the next subsections, we are going to calculate this integral.

Linearization and l 0 -integration
To integrate over l 0 , we first linearize l 0 by applying a transform as then the Jacobian of this transformation is |J| = ∂(l 0 , l 1 , l 2 , l ⊥ ) ∂(z, y, x, t) = 1. We arrive at with the integrand is corresponding to . (18) In this integrand, we have already used new variables The integrand now depends linearly on x. Thus, the x-integration will be taken easily by applying the residue theorem. For this purpose, one should first analyze the x-poles of this integrand. These poles are We realize that Im which depends on the sign of z. The location of x 0 in the x-complex plane will be determined by the sign of z, see Fig. 2 for more detail. We therefore should rewrite D 0 as following: with The contour integration in the x-plane.
By closing the integration contour over upper x plane when z > 0 and vise versa, as shown in Fig. 2, one then applies the residue theorem, the resulting then reads and The f lk and f − lk functions indicate the residue contributions from the x-poles in Eq. (20). These functions are defined as 6/28 The new kinematic variables introduced in this step are listed as following The delta function is defined It is important to note that A mlk , B mlk ∈ R, and C mlk ∈ C. Combining with the definition of f lk and f − lk in Eq. (27), we verify easily that

The y-integration
We are now going to evaluate the three-fold integrals which arrived the previous subsection, see Eqs. (25,26). In order to work out the y-integration, one has to linearize y by using the Euler shift t → t + y. However, we realize that the terms proportional to t 2 and y 2 in the integrand have the same sign. One can first make t 2 and y 2 having opposite sign by applying a complex rotation in the t-plane. The integrand written in terms of t have two poles which are Because of Eq. (32), we find that t 1,2 are located in the first and the third quarters of t-complex plane, as described in Fig. 3. As a matter of this fact, one should choose the integration contour on the fourth quarter of the t-complex plane, as Fig. 3. There are no residua of t-poles contribute to the t-integration contour. We, therefore, derive the following relation We are now applying the relation in Eq. (34). One then makes the rotation like t −→ it. The resulting reads After obtaining the opposite sign of t 2 and y 2 in these integrands, we proceed the linearization of y by performing the above Euler transformation, or t → t + y. We arrive at The y-integration will be taken by using the residue theorem. The locations of y-poles are more complicated than in the case of x-integration. According to Eqs. (37,38), combining with Eq. (32), it is easy to check that the imaginary parts of the y-poles in these integrands are which depend on the sign of t − blk AClk z. Similar to the x-integration, we should cut the integration of t into two segments t α lk z and t α lk z with α lk = blk AClk . The imaginary parts of the remaining poles which are the roots of the following equation become more complicated. They then contribute to the residua of the taken integrations. Closing the contour on the upper plane of y if t α lk z and vice versa, one takes into account the residua of the y-poles in Eq. (40). Finally, we arrive at with 9/28 and the integrand We have already introduced following kinematic variables: The g mlk and g − mlk functions will indicate the locations of y-poles in Eq. (40) which contributed to the integrations. They are defined as and We now make a shift t −→ t ′ = t − α lk z. The Jacobian of this shift is 1 and the t-integrals change the border to [0, ±∞]. The resulting reads with the new notations We note that D mlk ∈ R and F nmlk ∈ C.
With the definitions of f kl , f − kl in Eq. (27) and g mlk , g − mlk in Eq. (50), one again confirms that We have just arrived at the two-fold integrations. In the next subsections, we will present the approach to calculate these integrals (51, 52, 53, 54) in detail.

The t-integration
To linearize t, we perform a shift The Jacobian of this shift is To remove the quadratic term of t, we have to choose β mlk as the roots of the equation or these roots are written explicitly as With this, the integrand written in terms of t depends linearly on t which is We will choose ϕ mlk as the root of the equation Z mlk = 0, or their solutions are given by We note that the final result of D 0 is independent of the parameters ϕ mlk and β mlk . Without the loss of generality, we choose ϕ mlk = ϕ (1) mlk and β mlk = β (1) mlk in the following calculation. In this case, we have The relations between D mlk with the external momenta are shown in Table 1.
Table 1 D mlk are written in terms of external momenta.
To follow the calculation easily, we would like to omit the index for the kinematic variables which appear in Eq. (62) in the remaining text of this paper. 12 Fig. 4 The integration region.
2.3.1. In the case of D mlk < 0. In this case, β mlk 0 and ϕ mlk 0. The integration region now looks as Fig. 4.
To integrate over t, one first splits the integrations written in terms of t as follows We next rewrite the t-integrand in the form of We are going to apply the formula (138) for calculating the t-integrations. In order to use the formula (138), the integrand I nmlk (t, z) must have no poles in real t-axes. However, in the case of real masses, F nmlk are real, and so I nmlk (t, z) may have poles in negative real t-axes. To treat this problem, we make F nmlk → F nmlk + iρ ′ with ρ ′ 0. The final result is obtained by taking ρ ′ → 0. Using the master integral (138), one gets σz −∞ dt I nmlk (t, z) = 1 β(Q + P z) 13/28 and ∞ σz dt I nmlk (t, z) = 1 β(Q + P z) Where σ will be −ϕ mlk or − 1 βmlk . With the help of (73, 74), one obtains Here the G(z) function is given With the help of iρ ′ , all the logarithmic functions which appear in the z-integrands now are well-defined in z-complex plane. Summing up the above terms, one obtains

15/28
In the next steps, we are now going to calculate the z-integrals. We realize that and Where we use the notations ++, +−, · · · , −− which are corresponding to the appearance of f lk g mlk , f lk g − mlk , · · · , f − lk g − mlk in the mentioned formulae. From Eqs. (81, 82), we also confirm that From the formula (71), we realize that P mlk is a nonzero real in the case D mlk < 0. We rewrite G(z) as follows with Defining the arguments of logarithmic functions are with 16/28 In order to perform the z-integrals, we decompose logarithmic functions in the z-integrands, as Eq. (131). In special, one has ln S(σ, z) By determining that Im ± S(σ, z) P z + Q is independent of σ and using the formulae (92, 93), D 0 can be presented in the form We first emphasize that ln(P z + Q) might have poles in negative real-axes in the real mass cases. However, in these cases f lk = f − lk = 1 and g mlk = g − mlk = 1. As a result, one checks that (f − lk g − mlk − f lk g mlk ) ln(P z + Q) = 0. Thus we don't need to make Q −→ Q + iρ ′ as F nmlk 17/28 case. Secondly, the z-integrals now are splitted into three basic integrals which are These integrals can be calculated in concrete in the Appendix.

2.3.2.
In the case of 0 < D mlk Amlk Bmlk − α lk 2 and Amlk Bmlk − α lk < 0. In this case, we have It is easy to check that to D mlk < 0 case. As a result, the analytical calculation of D 0 in this case is same to the case of D mlk < 0. Bmlk − α lk > 0. In this case, β mlk and ϕ mlk are positive and we confirm that Therefore, the integration region now looks like Fig. 6. Fig. 6 The integration region.
Applying the same procedure, D 0 reads Where the new kinematic variables introduced in this formula are 19/28 and The last integrals written in terms of z will be evaluated by means of the basic integrals which are presented in Appendix.

In the case of
. From Eq. (80), each term relating to the following integral will be presented as Where C σ are coefficients in front of the mentioned integrals. We also apply the same trick for each term relates to We have already added to D 0 the extra terms which the sum of them is up to zero. These extra terms will contribute to the residue of z-poles when β mlk , ϕ mlk become complex.
Where new kinematic variables are introduced Instead of linearizing t, we are going to calculate the z-integrals directly. The resulting reads The H(t) is defined as Where W (1,2) nmlk are given Finally, one gets 21/28 It is easy to confirm that Moreover, one also verifies that Im t 2 + L mlk t − m 2 k + iρ > 0 and Im K mlk t + M mlk > 0.
Noting that X(t) is the second order polynomial written in terms of t in the argument of logarithmic functions. In detail, it is with Because Im(X(t)) and Im(K mlk t + M mlk ) have the same sign, the logarithmic functions are decomposed as follows ln ± X(t) and ln[X(t)] = ln(t − X (1) Finally, we arrive at Where the new kinematic variables are given The remaining integrals (written in terms of t) will be integrated by using the basic integrals which are devoted in the Appendix.
In the next step, we will extend this work for evaluating tensor one-loop four-point functions with complex internal masses. In the parallel and orthogonal space [11,12,23], a tensor one-loop N -point integral with rank M can be decomposed as The traditional tensor reduction for one-loop integrals has been proposed by Passarino and Veltman [4], later developed by Denner et al [5]. In these schemes, the form factors will be obtained by contracting the Minkowski metric (g µν ) and external momenta into the tensor integrals. At this stage, we have to solve a system of linear equations where the Gram determinants appear in the denominator. If the Gram determinants will vanish or become very small, the reduction method will break or spoil numerical stability (this 23/28 Let us consider f (z) is a rational function with lim z→∞ f (z) = 0 and no poles on the negative real z-axes, z k are poles of f (z). By closing the integration contour as Fig. 7, one then can derive the following relation
We consider three basic integrals in the following paragraphs.
(1) Basic integral I: The basics integral I is defined as 26/28 with x, y ∈ C.
The integral in right-hand side of this equation is nothing but it will be reduced to basic integral I.