Three-loop massive form factors: complete light-fermion and large-Nc corrections for vector, axial-vector, scalar and pseudo-scalar currents

We compute the three-loop QCD corrections to the massive quark form factors with external vector, axial-vector, scalar and pseudo-scalar currents. All corrections with closed loops of massless fermions are included. The non-fermionic part is computed in the large-Nc limit, where only planar Feynman diagrams contribute.


Introduction
Vertex form factors with external fermions play a crucial role in a number of phenomenologically interesting processes. Among them are the massive fermion production in electron positron collisions, and in particular the forward-backward asymmetry, where form factor contributions induced by vector and axial-vector currents are needed. Furthermore, building blocks to the decay rates of scalar and pseudo-scalar Higgs bosons are provided by the corresponding form factors. Last but not least, form factors constitute important toys, which help to investigate the structure of high-order quantum corrections.
In this work we consider vertex form factors where the external current is of vector, axial-vector, scalar or pseudo-scalar type. They are given by where for convenience the heavy quark mass m has been introduced in the scalar and pseudo-scalar currents such that no additional overall (ultraviolet) renormalization constants have to be introduced (as for the vector and axial-vector 1 cases) [1]. If j s or j p are used to compute properties of the Higgs boson there is a one-to-one relation of m to the corresponding Yukawa coupling.
We consider the three-point functions of the currents in Eq. (1) and a quark-anti-quark pair. The corresponding vertex functions can be decomposed into scalar form factors which are defined as with incoming momentum q 1 , outgoing momentum q 2 and q = q 1 − q 2 being the outgoing momentum at j δ . The external quarks are on-shell, i.e., q 2 1 = q 2 2 = m 2 and we have σ µν = i[γ µ , γ ν ]/2. We note that in all cases the colour structure is a simple Kronecker delta in the fundamental colour indices of the external quarks and not written out explicitly.
For later convenience we define the perturbative expansion of the scalar form factors as The two-loop corrections to the vector current contributions F v 1 and F v 2 have been computed for the first time in Ref. [2] (see also Ref. [3] for the fermionic contributions) and have been cross checked by several groups [4][5][6][7][8]. In some cases higher order terms in ǫ have been added. Two-loop axial-vector, scalar and pseudo-scalar contributions have been computed in Refs. [9][10][11] and recently been confirmed in Ref. [7] where O(ǫ) and O(ǫ 2 ) terms have been added. Three-loop corrections are only known for well-defined subsets of the vector form factor: The large-N c limit has been computed in Ref. [5] using the master integrals of [12]. This involves only planar integrals. The complete (planar and non-planar) light-fermion contributions to F v 1 and F v 2 have been obtained in Ref. [8]. In this reference also the results of the relevant master integrals are given. Let us mention that all-order corrections to the massive vector form factor in the large-β 0 limit have been considered in Ref. [13].
For the three-point functions one in general distinguishes singlet and non-singlet contributions. The former includes a closed fermion loop which contains the coupling to the external current. It is connected to the fermions in the final state via gluons as is shown in Fig. 1(a). In case the external current contains γ 5 singlet contributions need special attention since the anti-commuting definition for γ 5 can not be used. Instead prescriptions like the one introduced in Ref. [14] have to be applied.
If the external current does not contain γ 5 the singlet contributions can be treated along the same lines as the non-singlet part. However, in contrast to the latter the singlet contributions have massless cuts which requires modifications of the technique described in [5,8,12] to compute the master integrals. Thus, in this paper we restrict ourselves to non-singlet contributions (cf. Figs. 1(b)-(l)), i.e., the external current couples to the fermions in the final state. At three loops we compute the complete light-fermion contributions and consider the large-N c expansion of the remaining part. At one-and two-loop order all colour factors are computed and agreement with the literature [7] is found.
In the next section we introduce the notation and briefly mention some techniques used for the calculation. Afterwards analytical and numerical results are presented in Sections 3 and 4. We close with a brief summary in Section 5.

Technicalities
The techniques and the setup of the programs, which are used to obtain the results of this paper, are straightforward extensions of the works [5,8] and thus we refrain from repeating in detail the technical descriptions. Note, however, that in contrast to Ref. [5] we do not define a "super family", which includes the eight relevant planar families as sub-cases. Rather, we generated the input files for FIRE [15] from scratch and computed separate tables for each individual family. Let us mention that for the reduction to master integrals and the minimization of the latter it is useful to combine FIRE [15] with LiteRed [16,17], which provides important symmetry information. In fact, for the most complicated integral family the reduction took about a day of CPU time on a computer with 18 cores, even for general gauge parameter.
For the form factors it is useful to introduce the following variable which maps the complex q 2 /m 2 plane into the unit circle, as illustrated in Fig. 2. The low-energy (q 2 → 0), high-energy (q 2 → ∞) and threshold (q 2 → 4m 2 ) limits correspond to x → 1, x → 0 and x → −1, respectively. Furthermore, the interval q 2 < 0 is mapped to x ∈ (0, 1) and q 2 ∈ [0, 4m 2 ] to the upper semi-circle. Note that for x ∈ (0, 1) and x = e iφ with φ ∈ [0, π] the form factors have to be real-valued since the corresponding Feynman diagrams do not have cuts. This is different for the region q 2 > 4m 2 , which corresponds to x ∈ (−1, 0), where the form factors are complex-valued. Note that for negative x we interpret log(x) as log(x + i0) = log(−x) + iπ.
For the threshold limit (q 2 → 4m 2 , x → −1) it is convenient to introduce the velocity of the produced quarks which is related to x via For the analytic three-loop expressions we furthermore define In the practical calculation it is convenient to apply projectors in order to extract the scalar form factors. We refrain to provide them explicitly but refer to Ref. [7] where projectors for the four currents in Eq. (1) can be found.
All one-and two-loop Feynman integrals can be expressed as a linear combination of the 2 + 17 master integrals discussed in Ref. [8]. Note that our two-loop basis is smaller than the one of Ref. [7] where 23 non-singlet master integrals are given. After inserting the ǫ-expanded results for the master integrals into the expressions for the form factors we obtain the one-and two-loop expressions expanded up to order ǫ 4 and ǫ 2 , respectively. Our two-loop results agree with [7]. Let us repeat that we do not consider singlet contributions which occur for the first time at two loops. Note that they vanish for the vector current but give non-vanishing contributions for the other three currents.
At three-loop order we have 89 planar master integrals entering the large-N c expressions and 15 additional master integrals for the complete light-fermion n l and n 2 l contributions, only two of the them are non-planar.
To obtain the renormalized form factors we use the MS scheme for the strong coupling constant and the on-shell scheme for the heavy quark mass and wave function of the external quarks. In all cases the counterterm contributions are simply obtained by rescaling the bare parameters with the corresponding renormalization constants, Z αs , Z OS m and Z OS 2 . The latter is needed to three loops whereas two-loop corrections for Z αs and Z OS m are required. For the scalar and pseudo-scalar form factors also the overall overall factor m has to be renormalized (to three-loop order), which we choose to do in the MS scheme. Note that this is the natural choice if F s or F p are used for Higgs boson production or decay since then m takes over the role of the Yukawa coupling. The MS renormalization constants, of course, only contain pole parts. However, for the on-shell quantities also higher order ǫ coefficients are needed since the one-and two-loop form factors develop 1/ǫ and 1/ǫ 2 poles, respectively. Note that in our case the overall renormalization constants of all currents in Eq. (1) are equal to unity.

Analytic results
The analytic results for the form factors are expressed in terms of Goncharov polylogarithms (GPLs) [18] with letters −1, 0, +1 and r 1 . They are quite long and we refrain from presenting them in the paper. Rather we collect all relevant expressions in a computerreadable format; the corresponding file can be downloaded from [19]. To fix the notation we provide one-loop results for the six scalar form factors introduced in Eq. (2) up to the constant term in ǫ. For µ 2 = m 2 they are given by where the two GPLs can be written as In Eq. (8) we have the colour factor C F = (N 2 c − 1)/(2N c ). At two-loop order one has C 2 F , C F C A , C F T F n l and C F T F n h where C A = N c , n l counts the massless quark loops and n h = 1 the quark loops with mass m. At three-loop order there are the colour factors The latter is obtained by the replacements C F → N c /2 and C A → N c in the non-n l terms and taking only the leading contribution for large N c . This limit removes in particular all (pure) n h terms.
In the following three subsections we discuss the analytic structure of the form factors in three important kinematical regions where the external momentum q 2 is either small, large or close to the threshold for producing the heavy quarks on-shell. In these limits the expressions are compact and analytic results can be reproduced in this paper. We obtain the expansions of the full result by expanding the GPLs in the respective region. We restrict ourselves to the choice µ 2 = m 2 and refer to [19] for the general results. Subsection 3.4 contains a brief discussion on the infrared structure of the form factors and mentions several checks on our calculation.

Static limit
After expanding the GPLs in the low-energy limit we obtain the expansion of the form factor up to order (1 − x) 6 . In the following we present the results for the first two terms for the axial-vector, scalar and pseudo-scalar cases. The results for the vector currents can be found in Refs. [8,12] (see Sections 4.2.1 and 5.2, respectively). Since the interval q 2 /m 2 ∈ [0, 4] is mapped to the upper semi-unit circle in the complex x plane we use x = e iφ and parametrize our results as a function of φ which is real-valued. Our results read F a,(1) 1 In all cases the limit q 2 → 0 exists (i.e., there are no logarithmic terms in q 2 ) and all form factors become infrared finte. The infrared divergences are present starting from the φ 2 term. Note that for the vector case all quantum corrections vanish for q 2 = 0 and we have F v 1 (0) = 1 whereas for all other form form factors this is not the case. The static form factors for the vector and axial-vector case have been discussed in Ref. [20] up to two-loop order and the physical interpretations have nicely been summarized. Once the complete non-fermionic pieces and all singlet contributions are available the analysis of Ref. [20] can be extended to three loops. Three-loop corrections to F v 2 (0) have been computed in [21].

High-energy limit
In the high-energy limit, i.e. for x → 0 it is convenient to introduce where we have computed seven expansion terms, i.e., up to O(x 6 ), for all six scalar form factors. Note that the leading terms are identical both for F v 1 and F a 1 and for F s and F p since in this limit the quark masses in the numerator can be neglected and γ 5 is anticommuted through an even number of γ matrices to one end of the fermion string. As a consequence we also have that f a,(n,0) 2 = 0 since f v,(n,0) 2 = 0. We illustrate the structure of the analytic expressions by showing the terms of order x 0 and x 1 which are given by f a,(1,1) 2,lar Results for f v,(n,k) i,lar can be found in Eq. (13) of Ref. [8] and Eq. (16) of Ref. [12].
In general one has two powers of l x for each loop-order. At one and two loops we indeed observe l 2 x and l 4 x terms, respectively. However, for the shown colour structures we have at three-loop order at most l 5 x terms in the above expressions. Note that for the vector form factor one has l 6 x terms in the N 3 c term of f v,(3,0) 1,lar [22,23]. For a dedicated analysis of the leading logarithmically enhanced terms in power-suppressed contributions we refer to [24][25][26].

Threshold expansion
In the threshold region we can use our expressions for the form factors to obtain physical results for decay rates and productions cross sections since the corresponding real radiation is suppressed by a relative order β 3 . In the expressions we present in this subsection, the factor m in the definition of the scalar and pseudo-scalar currents, see Eq.
(1), has been transformed from the MS to the on-shell scheme which allows for a more straightforward comparison with results present in the literature.
In Refs. [5,8] the vector form factors F v 1 and F v 2 have been used to obtain results for the cross section σ(e + e − → QQ) in the limit of small quark velocities. In principle these results can be extended in order to incorporate the Z-boson contribution with vector and axial-vector couplings. We prefer to represent the results in a slightly different way, namely as the decay rate of a (hypotetical) boson with either vector, axial-vector, scalar or pseudo-scalar couplings which is related to the form factors via R δ are defined such that the exact result at tree level reads Note that the R δ enter physical quantities as building blocks. For example, we have where H and A are scalar and pseudo-scalar Higgs bosons with masses M H and M A , respectively, and the ellipses indicate quantum corrections from real radiation. In Eq. (21) G F is Fermi's constant, M Q is the heavy quark on-shell mass, σ 0 = 4πα 2 Q 2 Q /(3s), α is the fine structure constant and Q Q is the electric charge of the quark Q. The decay rates in (21) can be obtained from an effective Lagrangian of the form where v = 1/ √ 2G F is the vacuum expectation value.
We cast R δ in the form with K v = 3/2, K a = K s = K p = 1, n v = n p = 1 and n a = n s = 3. Then the threshold expansion of ∆ (i),δ starts at 1/β i and the leading term of the real radiation contribution to ∆ (i),δ is of order β 3−i . For the three-loop fermionic results we can provide results including β 0 . In the following we only show results for δ = a, s and p since the expressions for the vector current can be found in Refs. [5,8]. The one-and two-loop results for ∆ (i),δ are given by For the C 2 F /β 2 and the {C A C F , C F T F n l }/β terms of ∆ (2),δ we find agreement with Ref. [27]. At three-loop order we have − 320 9 π 2 log 2 (2β) + 128 9 π 4 log(2β) − 512 9 π 2 log 2 (2β) + 128 9 π 4 log(2β) + 8384 27 π 2 log(2β) At one-and two-loop order the leading terms are proportional to C F /β and (C F /β) 2 , respectively. At three-loop order we observe 1/β 3 terms only in the axial-vector and scalar case but not for the vector and pseudo-scalar currents. Our findings are in agreement with considerations in the non-relativistic limit which can be used to predict the leading terms of order (α s /β) n . In fact, in this limit ∆ δ can be written as a combination of the Sommerfeld factor and a factor taking into account P -wave scattering [28] (see also [27]): After the second equality sign y = C F α s π/β has been used and P δ is zero for S-wave and unity for P -wave processes, i.e., we have P v = P p = 0 and P a = P s = 1. The ellipses in the above equations represent subleading terms.

Γ cusp and checks
An interesting feature of the (renormalized) form factors is the presence of infrared poles that can be described by a universal function, the cusp anomalous dimension Γ cusp [29][30][31], which is independent of the external current. This means that we can write where F is any of our (six) scalar form factors and F f is the corresponding (ultraviolet and intrared) finite version. The factor Z, which is defined in the MS scheme and thus only contains poles in ǫ, absorbs the infrared divergences and F f is finite. The single ǫ pole of the n-loop corrections of Z contains the n-loop expression of Γ cusp . Using the notation where α (n l ) s is the strong coupling constant with decoupled heavy quark, we have We have used Eq. (28) for the four form factors F v 1 , F a 1 , F s and F p , and have determined the corresponding Z factor to three loops. This requires to use decoupling relations for α s up to two loops (including higher order ǫ terms [21]) since the calculation of the form factors described above has been peformed in the full theory with α s ≡ α which constitutes a strong cross check of our calculations. Furthermore, our results agree with Refs. [29,[31][32][33] where dedicated calculations of Γ cusp have been performed.
Besides the correct infrared structure there are several other checks which our analytic expressions fulfill: • As mentioned in Section 3.2 we observe that the results for the vector (scalar) and axial-vector (pseudo-scalar) form factors agree in the high-energy limit.
• We have furthermore performed the calculation for general gauge parameter which drops out for the renormalized form factors. Note that the cancellation is non-trivial and only occurs in the proper interplay between bare three-loop expression and wave function and quark mass counterterm contributions.
• We have also numerically cross checked all three-loop master integrals up to the finite term in ǫ using the program FIESTA [34].
• As a further check, our results fulfill the axial Ward identity which is given by Using Eq. (2) this leads to which is satisfied by our explicit results up to three loops after transforming the mass parameter m in Eq. (1) into the on-shell scheme.

Numerical results
In this section we evaluate both the exact result and the approximations in the various limits numerically. For illustration we plot for each of the six form factors the ǫ 0 term both as a function of x ∈ [−1, 1] and φ ∈ [0, π]. This means we cover the whole real q 2 axis. In the plots we restrict ourselves to the real part of the form factors. Using the results from [19] it is straightforward to obtain plots for the imaginary parts, as well.
At one-and two-loop order we plot the complete (non-singlet) results and set N c = 3, n l = 5 and n h = 1. These values are also used at three loops where the sum of the complete n l and the N 3 c terms are shown. For the numerical evaluation of the GPLs we use the program ginac [35,36] which is straightforward for real values of x. GPLs with complex arguments (in our case for x = e iφ with φ ∈ [0, π]) are evaluated with the help of transformation rules given in Ref. [36]. Some of the GPLs involving r 1 require extraordinary long run times. In some cases the results are even unstable. For this reason we generate in a first step for each GPL, which is present in our anayltic result, a data At three-loop order we add the complete light-fermion part for n l = 5 and the N 3 c contribution. Short-(blue), medium-(red) and long-(green) dashed lines correspond to the low-energy, high-energy and threshold approximation, respectively base for φ ∈ [0, π] and construct an interpolation function. Afterwards the numerical evaluation of the form factors is fast and stable.
The approximations shown in the plots contain terms up to order x 6 and (1 − x) 6 in the high-and low-energy expansion, respectively. At threshold we only include terms up  Figure 4: Same as Fig. 3 but for the the scalar and pseudo-scalar currents.
to order β 3 although higher order terms are available [19]. However, for the N 3 c term we observe a bad convegence behaviour which is the reason that we drop β 4 and higher terms.
In the high-energy limit the form factors exhibit logarithmic singularities. Thus, for the plots in the range x ∈ [−1, 1] we subtract the leading high-energy behaviour, i.e., all terms which are not power-suppressed by x, in order to ensure a smooth behaviour for x → 0. Furthermore, we multiply by (1+x) 4 to ensure that at threshold (i.e. for x = −1) the one-, two-and three-loop expressions become zero. This leads to a numerical enhancement for x = 1, however, also in this limit finite results are obtained. Thus, the function we use for the plots readsF Our results for the six scalar form factors are shown in Figs. 3 and 4. The exact result is shown as solid (black) curve and the approximations are plotted as dashed lines. Note that in all cases the whole range x ∈ [−1, 1] can be covered by the approximations, i.e., for each x-value there is at least one dashed curve on top of the (black) solid line.
It is interesting to mention that some of the plots show a peculiar bump-like structure at or close to x = 0. This is not a numerical artifact but probably related to the way the high-energy limit is subtracted. Note that the exact results and the high-energy expansions (red dashed curve), which is simple to evaluate numerically, perfectly agree with each other. which corresponds to q 2 /m 2 0.25 and q 2 /m 2 2.8, respectively.

Conclusions
We have considered the vertex form factors induced by vector, axial-vector, scalar and pseudo-scalar heavy quark currents, which play an important role both in the Standard Model but also in extensions. The form factors are parametrized by six scalar functions, which we have computed up to three-loop order. Our results are expressed in terms of GPLs with letters {−1, 0, 1, r 1 = e iπ/3 } and the argument x defined via the relation q 2 /m 2 = −(1 − x) 2 /x. After expanding the GPLs for small and large q 2 and around the threshold given by q 2 = 4m 2 we obtain compact and easy to evaluate expansions in the corresponding kinematical regions. We have discussed the convergence properties by comparing to the exact expressions. On the way to our three-loop result we have obtained the two-loop form factors including order ǫ 2 terms. This work extends the considerations of Refs. [5] and [8] to axial-vector, scalar and pseudo-scalar currents. Obvious next steps towards the full result are singlet contributions and the subset of Feynman diagrams containing closed massive fermion loops. However, it can be expected that even these sub-classes show a more involved mathematical structure and it is likely that not all pieces of the final result can be expressed in terms of GPLs.