The light-fermion contribution to the exact Higgs-gluon form factor in QCD

An analytical expression for the three-loop form factors for $ggH$ and $\gamma\gamma H$ is derived for the contributions which involve massless quark loops. The result is expressed in terms of harmonic polylogarithms. It fully agrees with previously obtained kinematical expansions, and confirms a recent semi-numerical approximation which extends over the full kinematic range.


Introduction
The study of the Higgs boson is one of the most promising ways to search for physics beyond the Standard Model (SM). A necessary precondition for this to be successful is the precise understanding of the relevant SM predictions. One of the most important quantities in this respect is the cross section for Higgs production in gluon fusion. In fact, significant theoretical efforts have been made to pin down its SM value, and to estimate the associated uncertainties (see Ref. [1] for a recent review). One source of uncertainties is the fact that, up to now, QCD corrections to the Higgs cross section beyond next-to-leading order (NLO) are based on the approximation of an infinitely heavy quark mediating the 1 gluon-Higgs coupling. For the top-quark contribution, which by far dominates the total cross section at the Large Hadron Collider (LHC), comparison of this limit to the full result at NLO shows agreement at the sub-% level for a Higgs mass of M H = 125 GeV, providing confidence in using this approximation also at higher orders of perturbation theory [2,3]. In fact, an explicit calculation of sub-leading terms in 1/m t at next-to-nextto-leading order (NNLO), combined with the high-energy limit of the cross section, further justifies this procedure [4][5][6][7]. Nevertheless, the lack of the exact top mass dependence still requires one to associate with it an uncertainty on the total cross section of the order of 1%. It is thus a non-negligible contribution to the overall uncertainty of about 5%, which also includes uncertainties induced by parton density functions (PDFs) and α s , for example (see Refs. [8,9]).
A related uncertainty arises from the bottom-quark induced Higgs-gluon coupling. While suppressed by the bottom Yukawa coupling, its effect on the leading order (LO) cross section is still a reduction by about 6%. Since the numerical value of the bottom-quark mass prohibits the analogous approximation as for the top quark, QCD corrections to the bottom-quark induced ggH amplitude are known only through NLO, without significant progress since their original calculation of more than 25 years ago [3]. Serious attempts to capture the dominant logarithmic contributions of the form ln m b /M H to higher orders in perturbation theory have been presented only recently [10]. Thus, also for this source, the LHC Higgs Cross Section Working Group assigned an uncertainty of roughly another 1% to the total gluon fusion Higgs cross section [8].
The total cross section at NNLO requires the inclusion of three-loop virtual corrections to the ggH amplitude (the "Higgs-gluon form factor"), two-loop corrections to single-real emission, and the one-loop double-real emission contributions which occur for the first time at this order. The real-emission contributions are sufficient if one aims for Higgs boson production at non-zero transverse momenta p ⊥ . In this case, top-mass effects have been addressed by several groups recently [11][12][13][14]. After estimates based on 1/m t expansions of the cross section which indicated a break-down of this approximation for p ⊥ 150 GeV [15,16], it came as a surprise to find the K-factor of the exact calculation to be fairly independent of p ⊥ [12,13]. This provides yet another indication that also for the total cross section, the QCD corrections are well described by their heavy-top limit. Also bottom-quark mass effects have been considered for finite p ⊥ [17,18] at this order of perturbation theory.
Concerning the virtual corrections, it took about ten years before the original numerical two-loop result for the Higgs-gluon form factor of Ref. [3], contributing to the total cross section at NLO, was expressed in closed analytic form using harmonic polylogarithms [19][20][21]. The analytic result for the γγH amplitude had been obtained one year earlier [22].
For the three-loop form factor, only approximate results are available up to now, most notably through expansions in the heavy-quark limit [23][24][25][26]. While this expansion is expected to work very well for on-shell Higgs production mediated by a top-quark loops, it will break down for the bottom-mediated contribution, or in cases where the Higgs is produced as a virtual intermediate particle, for example in off-shell or double-Higgs production. Knowledge of the general dependence of the three-loop Higgs-gluon form factor on the quark/Higgs mass ratio is thus very desirable.
Very recently, the expansions in 1/m t were combined with the leading behavior of the amplitude at the top-threshold, i.e.ŝ ≈ 4m 2 t (see also Ref. [27]) in order to construct Padé approximants for the three-loop ggH amplitude which should be valid-within intrinsic Padé uncertainties-for general Higgs and quark masses [28].
In this paper, we provide an analytic result for a subset of the virtual three-loop corrections, namely those involving light (massless) quark loops in addition to the massive (topor bottom) quark loop. Using integration-by-parts (IBP) identities, we reduce the occurring Feynman integrals to a set of master integrals, which we manage to solve in terms of harmonic polylogarithms. Comparing our result to Ref. [28], we find full agreement for this light-fermion component within the uncertainty estimate of Ref. [28]. As a byproduct of this calculation, we also obtain the three-loop γγH form factor from which one may directly derive the exclusive photonic decay rate of the Higgs boson through NNLO.

Calculation
The amplitude for the processes ggH and γγH can be parameterized with the momenta q 1,2 of the two external vector bosons as where we already implied Bose symmetry.
Here and in what follows, a and b denote color indices of the adjoint representation, while µ and ν are d-dimensional Lorentz indices. Because of the trivial color structure in eq. (2.1a) which can be projected out using (δ ab /N A )M ab;µν ggH , where N A is the number of gauge generators, we ignore the color structure in the following and focus only on the Lorentz structure of the amplitudes. For both amplitudes the Ward identity 1 q 1µ j,ν (q 2 )M µν = 0 (2.2) yields a constraint In case of photons in the external state the even stronger Ward identity q 1µ M µν γγH = 0 leads in addition to a vanishing form factor A γγH .
For physical quantities the only contribution stems from the form factors C. Therefore, the physical part of the amplitudes can be written as Since gluons and photons do not directly couple to the Higgs boson, the Feynman diagrams contributing to the ggH and γγH amplitudes always involve at least one closed massive quark loop if higher orders in the electroweak coupling are neglected. In this paper, we address the calculation of the component of the form factors C which, in addition to this massive quark loop, involve a closed loop of a light quark (assumed massless here). Since the corresponding Yukawa coupling vanishes, the Higgs boson will still only couple to the diagram via the massive quark loop.
In section 2.1 we describe the toolchain used to express the contribution from light quarks to the form factors C in terms of master integrals. In section 2.2 the method for the calculation of the master integrals is explained.

Toolchain
For the calculation of the light-quark contribution to the ggH and γγH form factors it is required to evaluate the Feynman diagrams in fig. 2.1. These Feynman diagrams are generated in a first step using the tool qgraf [29]. After the insertion of Feynman rules in R ξ -gauge with the help of q2e [30,31] the diagrams are mapped to a set of seven topologies via exp [30,31].
A custom code for the computer algebra system FORM [32] was written in order to further process the output of exp. We use the projector to project out the form factor C which already implies the validity of the Ward identity in eq. (2.3). Moreover, we use to project out all the form factors in eq. (2.1) in order to check our calculational setup by explicitly verifying the validity of the Ward identities, see eq. (2.3) and below. The color factor of each diagram is determined via the FORM package color [33].
After projecting out the form factors, the results can be expressed in terms of scalar Feynman integrals, which are subsequently reduced to 45 master integrals using integrationby-parts identities [34,35] and the Laporta algorithm [36], implemented in the computer program Kira 2 [37,38].
After the reduction to master integrals the dependence on the gauge parameter ξ drops out and the validity of eq. (2.3) as well as A γγH = 0 is confirmed.

Calculation of master integrals
A very successful technique for the evaluation of two-scale Feynman integrals is based on the method of differential equations [39][40][41]. The solution of the resulting coupled system of differential equations simplifies significantly if it can be written in the canonical form proposed in Ref. [42], where the right hand side of the system is proportional to = (4 − d)/2. An algorithm to compute a basis transformation to such a canonical form was presented by Lee in Ref. [43]. We utilize its implementation in the computer program epsilon [44] in order to evaluate the relevant master integrals.
The class of transformations Lee's algorithm is able to find is restricted to be rational in the kinematic variable. Hence, a proper choice for the kinematic variable is inevitable to obtain a canonical form. For our purposes, an appropriate variable is where τ = m 2 H /(4m 2 q ) with the mass m H of the Higgs boson and the mass m q of the massive quark.
Ordering the master integrals by the number of lines in their topology yields a blocktriangular structure of the system of differential equations. For most applications it is sufficient to transform only the diagonal blocks into the previously described canonical form. The differential equations for master integrals of a certain block f (x, ) can then be written as where g(x, ) consists of already solved master integrals of a lower topology, and M (x) is fuchsian, i.e. it possesses only simple poles in x. For the master integrals entering the light-quark contributions, these poles lie at x = −1, 0, 1. The homogeneous part of (2.8) can be solved using an evolution operator U (x, x 0 ; ) which fulfills via iterated integrations in terms of multiple polylogarithms. This evolution operator allows expressing the full solution as where f (x 0 , ) are the boundary conditions of the master integrals at x = x 0 . In order to simplify the integral in (2.10) we chose to transform the off-diagonal blocks B(x, ) of the system of differential equations into fuchsian form via epsilon. Doing that ensures f (x, ) to be a linear combination of multiple polylogarithms without rational function prefactors in case g(x, ) is of this form as well.
The boundary conditions f (x 0 , ) are calculated as an asymptotic expansion around x 0 = 1, which corresponds to a limit where the quark mass m q is large compared to the Higgs mass m H . For this purpose, we expand the master integrals by subgraphs [45][46][47][48] as it is implemented in the computer program exp [30,31].
Via the method described in this section we were able to solve not only the 45 master integrals relevant for the light-quark contributions, but in total 202 master integrals for the full amplitudes (including n 0 l terms). The master integrals entering the lightquark contributions were cross checked against numerical results obtained by the package FIESTA [49].

Results
In this section, we define the parameterization of our results, which we have evaluated for a general gauge group with fundamental and adjoint quadratic Casimir eigenvalues C F and C A , and fundamental trace normalization T F . For QCD, it is C F = 4/3, C A = 3, and T F = 1/2. The actual analytic expressions are deferred to the appendix for the sake of readability of the main text. In this section, we restrict ourselves to a numerical presentation of the results.

Results for C γγH
The form factor C γγH is presented as a perturbative series in the strong coupling constant α s , renormalized in n l -flavor QCD in the MS scheme: where v denotes the vacuum expectation value and α the electromagnetic coupling constant. In order to fix the notation we provide the one-loop result as with the electric charge Q q of the massive quark and H 0,0 = ln 2 (x)/2. The three-loop result can be parameterized via where Q j are the electric charges of the n l light quarks. The term C

Results for C ggH
The form factor C ggH can also be written as a series in α s , We provide again the one-loop result to fix the notation, In contrast to the C γγH form factor, the purely virtual ggH result is not finite after the ultraviolet renormalization procedure. This is due to infrared divergences which cancel against real corrections, or are absorbed into PDFs. In Ref. [50] it is shown that the structure of these infrared divergences is universal and can be subtracted: The factors I (1) g and I (2) g are given by [50,51] and The finite three-loop terms can be parameterized as whereC (2,0) ggH denotes the contribution of Feynman diagrams without a light quark loop. Since at the three-loop level there are no diagrams with more than two closed quark loops, theC (2,2) ggH contribution originates solely from the subtraction terms of eq. (3.6b). We further decomposeC (2,1) ggH into its contributions to different color factors, Explicit results for C   [28]. The latter is associated with a systematic uncertainty due to the approximation procedure, which is dominated at large τ by the absence of any input from this kinematical region into the Padé approximants. Our results indeed confirm the associated uncertainty estimate for the light-quark terms up to rather large values of τ (corresponding to large Higgs masses/virtualities or small quark masses).

Conclusions
The Higgs-gluon form factor is an essential component for the theoretical description of Higgs physics at hadron colliders. It enters the total cross section for single and double Higgs production at the LHC, for example. In QCD, it involves a massive quark loop   The blue (dashed) line shows the Padé approximation found in Ref. [28]; the associated uncertainty estimate is indicated by the blue shaded band. The lower panel shows the difference between the Padé approximation and our result.
which mediates the Higgs-gluon/photon coupling. The fact that, until recently, the threeloop form factor has been known only in the limit of a very heavy mediating quark has restricted its applicability to top-quark mediated on-shell (or not-too off-shell) single production of the SM-like Higgs boson. Even there, the lack of an exact result implied non-negligible uncertainties. In cases such as double-Higgs production, off-shell or bottomquark mediated single-Higgs production, or the production of heavy beyond-the-SM (BSM) Higgs bosons, the expansion fails and one had to resort to the NLO result.
In this paper, we provided an analytic three-loop result for the component of the gluon-Higgs form factor which, in addition to the massive quark loop, involves a closed massless quark loop. We showed that it can be computed in closed form for a general quark mass, and presented it in terms of harmonic polylogarithms. Comparison of this result to a recent semi-numerical evaluation of the full Higgs-gluon form factor shows very good agreement at the level of the estimated numerical uncertainties. As a byproduct of our calculation, we also presented the amplitude for the decay rate of the Higgs boson to photons at NNLO for general quark (and Higgs) masses.
A large portion of our technical setup is applicable also to the full form factor. However, the calculations are much more expensive than for the light-quark terms considered here. Moreover, one encounters elliptic integrals which cannot be expressed in terms of harmonic polylogarithms. A fully analytical result for the three-loop form factor therefore requires further efforts.

A Results for C γγH
In this appendix, we provide explicit formulas for the two-loop result and the newly computed light-quark contributions to the three-loop result of the γγH form factor, cf. eqs. (3.1,3.3). We use H a ≡ H( a; x) to denote harmonic polylogarithms [54,55], ζ n ≡ ∞ j=1 j −n for Riemann's zeta function, and the short-hand notation L µ ≡ ln(µ 2 /m 2 q ) with the renormalization scale µ.
The two-loop result has been know for about 15 years [22]. In our notation, it reads γγH ) is the result for a quark mass renormalized in the on-shell (MS) scheme. At three loops, the non-singlet contribution, with the quark mass renormalized in the on-shell scheme, is given by The singlet contribution is The results with an MS renormalized quark mass can be written as where ∆C (non-sing) γγH Note that the singlet component appears for the first time at the three-loop order and is therefore renormalization scheme independent.

B Results for C ggH
In this appendix, we provide explicit formulas for the two-loop and the newly computed light-quark contribution to the three-loop result of the ggH form factor, cf. eqs. (3.4,3.10).
The notation is the same as in Appendix A.
Again, the two-loop result has been known for about 15 years [19][20][21]. In our notation, it reads where again the symbols with (without) a bar on top denote the results for an MS (onshell) renormalized quark mass.

C Ancillary File
The ancillary file ggh-aah-nl.m contains the main results of this paper in an electronic form readable by Mathematica [56]. The following table describes its notation: The OS (MSbar) versions correspond to a quark mass renormalized in the on-shell (MS) scheme. The harmonic polylogarithms are in a format compatible with the Mathematica package HPL.m [55].