Next-to-leading order electroweak corrections to gg → H H and gg → gH in the large-m t limit

We compute two-loop electroweak corrections to Higgs boson pair and Higgs plus jet production, taking into account all sectors of the Standard Model. All diagrams with virtual top quarks are computed in an expansion for large top quark mass up to order 1 /m 8 t or more. We present analytic results for the form factors and discuss the convergence properties. For the process gg → gH we also consider QCD corrections in the large-m t limit.


Introduction
Higgs boson pair production is a crucial process to obtain deeper insight into the symmetry breaking mechanism of the Standard Model (SM).For this reason, it is one of the most important processes studied in detail at the Large Hadron Collider at CERN and similarly at its High Luminosity upgrade which will begin operation within this decade.The main SM production mechanism for Higgs boson pairs is via gluon-gluon fusion and a number of higher-order corrections have been computed, mainly in the context of QCD.As far as electroweak corrections are concerned comparatively very little is known.First steps have been taken in Refs.[1,2].In Ref. [1] the two-loop box diagrams have been considered where a Higgs boson is exchanged between the massive top quarks.It has been shown that a deep expansion in the high-energy limit leads to results for the form factors which are valid in a large part of the phase space.In Ref. [2] top-quark Yukawa corrections have been considered, partly in the infinite top quark mass limit.Electroweak corrections proportional to the Higgs self-couplings have been considered in Ref. [3] using a numerical approach.In the present work we compute the complete NLO electroweak corrections as an expansion in the large top quark mass limit, including sub-leading terms up to 1/m 10 t .The corresponding corrections in the case of QCD have been computed in Refs.[4][5][6].
A similarly important process at the LHC is the production of a Higgs boson in association with a jet.As for Higgs boson pair production the dominant production channel is gluon-gluon fusion, with the partonic process gg → gH.NLO QCD corrections have been considered in a number of works: in the large-m t limit [7], in the high-energy limit [8][9][10] and numerically, including exact dependence on m t [11][12][13].NNLO corrections have even been computed in the infinite top quark mass limit [14][15][16][17][18]. NLO electroweak corrections via massless bottom quark loops have been computed in Ref. [19], and the corrections induced by a trilinear Higgs coupling in the large top mass limit have been recently calculated in Ref. [20].In this work we compute, for the first time, the full NLO electroweak corrections involving virtual top quark loops.We consider all sectors of the Standard Model and perform an expansion for large m t up to order 1/m 8 t .Furthermore, we provide analytic results for the NLO QCD corrections, again expanded up to 1/m 8 t .These expressions will be of interest for cross checks of numerical results and for the construction of approximation formulae involving expansions in different limits.
Calculations in the electroweak sector of the Standard Model are in general much more complicated than in the strong sector since many different mass scales are involved.For the case of QCD corrections it has been shown (see, e.g., Refs.[21][22][23]) that precise approximations can be obtained by combining expansions performed in different regions of the phase space.This motivates developing these techniques beyond QCD to the electroweak sector of the Standard Model.In this work we take a first step in this direction by considering the region in which the top quark mass is larger than all other kinematic invariants.While the radius of convergence of such an expansion is limited only to small values of the centre-of-mass energy, the results will serve as benchmarks for cross checks of other expansions or for numerical results.This paper is organized as follows: in the next section we define the form factors which describe the two processes considered, and the technical details needed for our calculation are presented in Section 3. In particular, we describe the asymptotic expansion and our renormalization procedure.Section 4 contains our results for Higgs boson pair production and Section 5 is dedicated to the electroweak corrections to gg → gH.The QCD corrections to gg → gH are discussed in Section 6.In all cases we study the influence of the higher-order 1/m t terms on the form factors and provide our complete analytic expressions in the ancillary files of this paper [24].A brief summary of our findings is provided in Section 7.
2 Form factors for gg → HH and gg → gH
Here q ij = q i • q j with q 2 1 = q 2 2 = 0 and q 2 3 = q 2 4 = m 2 H . p T is the transverse momentum of the final-state Higgs bosons, given by with the Mandelstam variables Using these definitions we introduce the form factors F 1 and F 2 as where a, b are adjoint colour indices, and α s (µ) is the strong coupling constant evaluated at the renormalization scale µ.We decompose the functions F 1 and F 2 introduced in Eq. ( 5) into "triangle" and "box" form factors. F 1 has contributions with zero, one and two s-channel Higgs boson propagators whereas F 2 only has box contributions.Thus we write In order to obtain this decomposition it is important to re-write the factors of s which occur in the numerators during the calculation using s/(s−m 2 H ) = 1+m 2 H /(s−m 2 H ). Note that at two loops F tri is not the same as the form factor for single Higgs boson production (as is the case for QCD corrections), since loop corrections to the HHH vertex also enter here.
We define the perturbative expansion of the form factors as where α is the fine structure constant and the ellipses indicate higher-order QCD and electroweak corrections.
In Section 4 we discuss the results for the squared matrix element constructed from the form factors F tri , Ftri , F box1 and F box2 .Analytic results for the leading-order form factors (F box1 and F box2 ) are available from [25,26].Two-loop corrections to F (0,1) box1 and F (0,1) box2 originating from the exchange of a virtual Higgs boson have been computed in Ref. [1] in the high-energy limit.
In Fig. 1 we show sample one-and two-loop diagrams contributing to gg → HH.At two-loop order we have: • one-particle irreducible box and triangle diagrams, • one-particle reducible diagrams with a one-loop correction to the HHH vertex or a one-loop self-energy correction to the Higgs propagator of a one-loop gg → H → HH diagram, • one-loop tadpole corrections to one-loop diagrams.
At two-loop order there are also contributions without top quarks which are not suppressed by small Yukawa couplings.In these contributions the gluons couple to light quarks and the connection to the final-state Higgs bosons is mediated via Z bosons.An example is given by diagram (g-1) in Fig. 1 if a light quark runs in the fermion loop.In our expansion these contributions formally contribute to the m 0 t term, however in this work we do not compute such diagrams; they can be computed following the approach of Ref. [19].

gg → gH
The amplitude for the process g(q 1 )g(q 2 ) → g(q 3 )H(q 4 ) (8) Figure 1: One-and two-loop Feynman diagrams contributing to gg → HH.Dashed, solid, wavy and curly lines correspond to scalar particles, fermions, electroweak gauge bosons and gluons, respectively.can be decomposed into four physical Lorentz structures [8] 1 The corresponding four form factors F 1 , . . ., F 4 are defined through Figure 2: One-and two-loop Feynman diagrams contributing to gg → gH.Dashed, solid, wavy and curly line correspond to scalar particles, fermions, electroweak gauge bosons and gluons, respectively.Diagrams are also shown which contribute to the NLO QCD corrections.
where c is the adjoint colour index of the final-state gluon, X gggH 0 is given by and the perturbative expansions of the form factors are defined as in (7).The Mandelstam variables are defined as in Eq. ( 4); the only difference with respect to gg → HH is that here q 2 3 = 0 and p 2 T = u t/s.Sample Feynman diagrams for gg → gH are given in Fig. 2. The classification is similar to gg → HH, we again include all one-particle reducible and all tadpole contributions.Note that for the QCD corrections, we also include the one-loop self-energy corrections to the gluon propagators and the one-loop vertex corrections to the triple-gluon vertex of the one-loop diagrams.The corrections to the quartic-gluon vertex do not appear at the two-loop order of this process.
3 Technical setup

Asymptotic expansion of the two-loop amplitudes
For the generation of the gg → HH and gg → gH diagrams and the corresponding amplitudes we use qgraf [27].As input we use the Lagrangian file of the full Standard Model shipped with tapir [28], which is derived from the Feynman rules of UFO [29].tapir translates the qgraf output to FORM [30] notation and generates further auxiliary files which are useful for the manipulation of the amplitudes.The large-m t expansion is realized with the help of exp [31,32] which generates the corresponding subdiagrams and maps them to various integral families. 2e apply the large-m t limit as where no additional hierarchy is assumed among the scales on the right-hand side.This leads to the following integral families: • one-and two-loop one-scale vacuum integrals, • one-loop massless triangle integrals where two external lines are massless, • massive vertex integrals where for one external leg we have (q 1 + q 2 ) 2 = s and for the other two legs we have q 2 3 = q 2 4 = m 2 H (for gg → HH) or q 2 3 = 0 and q 2 4 = m 2 H (for gg → gH), • for the QCD corrections to gg → gH we also need massless one-loop box families with one external mass q 2 4 = m 2 H ; explicit analytic results can be found in Ref. [34].
Our FORM-based setup automatically performs a reduction of arbitrary members of each family to master integrals, which are well known in the literature (see, e.g., Refs.[35,36]).
The tadpole integrals are computed by MATAD [37] and the remaining integral families are reduced by IBP reduction rules derived by LiteRed [38] which have been implemented in FORM.Furthermore all of our reduction routines can deal with tensor integrals, avoiding the need to construct additional projection operators.In Fig. 3 we show how the various integral families appear due to the asymptotic expansion in the large-m t limit.In the Feynman gauge we have performed an expansion of the form factors up to order 1/m 10 t (1/m 8 t ) for gg → HH (gg → gH).In order to check our calculation, we also introduce general gauge parameters ξ Z , ξ W and ξ γ for the Z and W bosons and the photon.From the technical point of view ξ γ does not introduce any additional complexity since no new mass scale is introduced.It drops out after summing all bare two-loop diagrams.This is not the case for ξ Z and ξ W since they appear in combination with gauge boson masses in the gauge boson and Goldstone propagators.Furthermore, ξ Z and ξ W only drop out after renormalization.For this check we assume and perform an expansion which includes terms up to order 1/m To check the cancellation of ξ Z and ξ W we have to consider the combination of the bare two-loop diagrams and the counterterm contribution from the wave function of the external Higgs boson (see also below), which also depends on ξ W and ξ Z . 3It is a welcome and non-trivial check of our calculation that up to this expansion depth, ξ W and ξ Z drop out of the gg → HH and gg → gH amplitudes.

Renormalization
In the following we concentrate on the electroweak sector; for the discussion of the renormalization and the treatment of the infra-red divergences which occur for the NLO QCD corrections to gg → gH we refer to Section 6.
For the renormalization we follow the standard procedure as outlined, e.g., in Refs.[39,40].We express our one-loop amplitudes for the form factors in terms of the parameters where e = √ 4πα, and introduce one-loop on-shell counterterms (see, e.g.Eqs. ( 143), ( 153) and (421) of Ref. [40]).Furthermore, we have to renormalize the wave function of the external Higgs boson, which we also perform in the on-shell scheme (see Eq. (144) of Ref. [40]).
We consistently include tadpole contributions in all parts of our calculation (in the twoloop gg → HH and gg → gH amplitudes, and the gauge boson and fermion two-point functions needed for the counterterms).This guarantees that the top quark mass counterterm is gauge-parameter independent.This prescription is equivalent to the so-called Fleischer-Jegerlehner tadpole scheme [41]. 4or the numerical evaluation of the form factors we transform our results into the so-called G µ scheme where the Fermi constant G F and the gauge boson masses m Z and m W are the input parameters, and the fine structure constant α and the weak mixing angle θ W are derived quantities.(see, e.g., Section 5.1.1 of Ref. [40]).In this scheme it is convenient to express the final result in terms of the variable Although we have computed the exact top quark mass dependence of all counterterm contributions it is convenient to expand them in 1/m t and combine the individual terms with the expanded bare two-loop amplitude.We do not expand the (finite) quantity ∆r, which performs the transformation from the α to the G µ scheme, in the large-m t limit but retain its exact dependence on m t .
Note that the NLO electroweak corrections do not produce infra-red divergences.Thus, already after renormalization we obtain the finite results for the form factors.This is not the case for the NLO QCD corrections to gg → gH; the infra-red subtraction necessary to produce a finite result is discussed in Section 6.Let us also mention that our NLO electroweak form factors do not have an explicit dependence on the renormalization scale since all parameters are renormalized in the on-shell scheme.
4 Results for gg → HH

Analytic results
It is instructive to begin by discussing the leading contributions in the large-m t expansion, of order m 4 t and m 2 t , which are present in F (0,1) tri and F (0,1) box1 .Our results for the two-loop form factors read For reference, we also provide the large-m t limit of the leading-order form factors which are given by Results for F tri and F box1 have also been presented in Ref. [2], in which leading m 2 t corrections to the ggH and ggHH vertices at two-loop order are taken into account using an effective-theory approach, while one-particle reducible diagrams have been computed with full m t dependence.Furthermore, all one-and two-particle reducible diagrams involving Yukawa couplings have been considered.After extracting the m 4 t and m 2 t terms we find agreement with our results.To make this comparison it is important to consider subleading terms in the expansion of the LO form factors which are factored out in Ref. [2] and contain exact m t dependence.
Using the asymptotic expansion described in Section 3.1 we have obtained expansion terms up to order 1/m 10 t .Up to order 1/m 4 t we have performed the calculation for general gauge parameters and we have verified that they drop out from the renormalized results.The higher-order 1/m t terms have been computed only in the Feynman gauge.The analytic expressions for the form factors can be obtained from [24].
In our analytic expressions we observe poles of the form 1/(s−4m 2 H ) k where k > 0 is larger for the higher-order 1/m t terms.The origin of these terms are massive one-loop triangle (co-)subgraphs, such as the one on the first row of Fig. 3 with external squared momenta s, m 2 H and m 2 H .The expansion of the subgraph leads to numerators in the triangle diagram and the 1/(s − 4m 2 H ) terms result from the subsequent reduction to master integrals.We note that the poles are spurious; for each 1/m t term the limit s → 4m 2 H exists.We also point out that the m 0 t term presented here is not complete, since it should also receive contributions from diagrams without top quarks, for e.g., the first diagram in Fig. 3 where the top quarks are replaced by light quarks.We do not compute such diagrams in this paper.They can be computed following the approach of, e.g., Ref. [19] where similar contributions to gg → gH have been considered, or with the help of expansions as proposed, e.g. in Ref. [1].

Numeric results
For the numerical evaluation of our form factors we adopt the G µ scheme and use the following input values Furthermore, we express the form factors in terms of s and p T and introduce the parameter In the following we choose ρ p T = 0.1 and discuss results for the squared matrix element For the numerical evaluation of the massive two-and three-point functions we use the program Package-X [43].
For reference, in Fig. 4  In Fig. 5 we show the NLO quantity Ũ(0,1) ggHH as a function of √ s.The curves include increasing expansion depths starting from the leading term proportional to m 4 t (which originates from F (0,1) tri ) up to 1/m 10 t .For the √ s axis we choose values from the Higgs pair production threshold at 2m H = 250 GeV up to √ s = 380 GeV.Note that convergence of the expansion is not expected beyond the top quark pair production threshold at 2m t = 344 GeV.Below this value we observe, at first sight, a reasonable convergence.Below √ s ≈ 300 GeV a significant shift is obtained from the constant contribution proportional to m 0 t and higher order 1/m t terms are small up to 1/m 8 t .However, the 1/m 10 t contribution again provides a sizeable shift, which is clearly visible on the right panel which shows the ratio with respect to the m 0 t contribution.This behaviour is due to diagrams with a closed quark loop which contains both top and bottom quarks, see, e.g., the second diagram in Fig. 3.Such diagrams contain cuts through a top quark and W boson and thus the large-m t expansion is expected to break down above √ s = m t + m W ≈ 250 GeV.Diagrams with such a cut contribute to both F 1 and F 2 .To demonstrate this, in Fig. 6 we show the results for Ũ(0,1) ggHH where we set all diagrams containing a bottom quark to zero in the finite parts. 5We indeed observe that after removing these contributions the large-m t expansion converges as expected up to the threshold at √ s = 2m t .We note that the two-loop diagrams have further cuts where no top quark is involved at √ s = 2m W , 2m Z , 2m H .In our approach all of these are taken into account exactly, so they do not affect the convergence of the large-m t expansion.
In view of the above discussion the validity of the leading m t terms (see Section 4.1 and Ref. [2]), and indeed of the deeper large-m t expansion, for a description of the electroweak corrections to gg → HH is questionable.More insight will be provided in a future publication which considers the small-t expansion of these diagrams in the style of Ref. [23].We have compared our one-loop form factors to Ref. [20] and find agreement up to 1/m 14 t .We also compare with the subset of NLO contributions induced by the trilinear Higgs boson coupling considered in Ref. [20], by extracting the corresponding pieces from our bare two-loop form factors.We have compared up to 1/m 2 t and find agreement.Our result provides solid predictions for the energy range m H ≤ √ s ≲ 300 GeV and will thus serve as an important cross check for future (analytic) calculations in different kinematic limits or of numerical evaluations.
6 NLO QCD corrections to gg → gH in the large-m t limit A finite expression for the NLO virtual QCD corrections to gg → gH is obtained after introducing counterterms for the ultra-violet poles and subtracting the infra-red divergences.We first renormalize the strong coupling constant in the MS scheme with six active flavours.The top quark mass and gluon wave functions are renormalized in the on-shell scheme. 6Afterwards we express the form factors in terms of α s (µ), with five active flavours.Finite form factors are then obtained via the subtraction (i = 1, 2, 3, 4) where F (1,0) i,ren are the ultra-violet renormalized form factors.The quantity I (1) g on the right-hand side of Eq. ( 22) is given by [44] with β 0 = 11C A /12 − T F n l /3, where T F = 1/2, C A = n c and n l is the number of massless quarks.
For illustration we present the one-and two-loop expressions for the form factors F (0) 1 and F (1,0) 1,fin to the expansion order 1/m 2 t and m 0 t , respectively.Deeper expansions can be found in the supplementary material [24] of this paper.At one-loop order we have and the two-loop expression is given by For the construction of the squared matrix element the infra-red subtracted form factors (Eq.22) have been used.The panel on the right shows the result normalized to the m 0 t expansion term.
where n c = 3 and Li 2 is the dilogarithm.
In Fig. 9 we show the NLO QCD corrections to U gggH for ρ p T = 0.1 as a function of √ s.We observe a rapid convergence, even beyond the top quark threshold (although the expansion is not expected to produce the correct result in this region).In fact, only the 1/m 2 t terms lead to a shift of a few percent; the higher-order expansion terms are much smaller.This behaviour can be explained by the dominance of the diagrams involving ggH triangle contributions and the suppression of the box-type Feynman diagrams.

Conclusions
In this work we consider the gluon-fusion induced processes gg → HH and gg → gH and compute complete NLO electroweak corrections in the large top quark mass limit and present results for the form factors up to order 1/m 10 t and 1/m 8 t , respectively.We discuss the renormalization procedure in detail and compute all counterterm contributions without assuming any mass hierarchy.Thus, this part can also be applied to expansions in other kinematic limits or an exact (numerical) calculation.
Partial electroweak results for gg → HH are already available in the literature [1,2]; in this work we provide sub-leading terms in the large-m t expansion.
For gg → HH the expansion in 1/m t does not show a convergent behaviour in the physical region 2m H ≲ √ s ≲ 2m t .We have demonstrated that this is due to diagrams involving a cut through a W boson and a top quark.If these diagrams are omitted, we observe reasonable convergence below √ s ≈ 330 GeV.Despite the limited applicability of the large-m t expansion we believe that our results serve as reference for future expansions in other kinematic regions or exact (numerical) calculations.Despite the convergence issues, if we assume that the order of magnitude is at least correct, in the large-m t region the electroweak contribution provides a correction of a few tens of percent with respect to the leading order.
For the NLO electroweak corrections to gg → gH we observe very good convergence below the top quark threshold.In particular, for √ s < 300 GeV we can provide precise predictions on the basis of an expansion which includes corrections up to 1/m 8 t .In this region the electroweak corrections are small, below the percent level with respect to the leading order.
We also provide NLO QCD corrections for the four form factors needed for gg → gH up to 1/m 8 t .Here a rapid convergence is also observed up to the top quark threshold.

Figure 3 :
Figure 3: Asymptotic expansion of two sample Feynman diagrams.The subgraphs left of the stars have to be expanded in the small quantities, i.e., masses, external momenta or loop momenta of the co-subgraphs, which are to the right of the stars.
Results are shown up to order 1/m 10 t .The panel on the right shows the result normalized to the m 0 t expansion term.
we show the LO contribution to ŨggHH as a function of √ s.Below the top quark threshold the expansion converges well, however it converges more slowly as √ s gets closer to 2m t .W , see text for details.The panel on the right shows the result normalized to the m 0 t expansion term.
ggHH without contributions involving a cut at √ s = m t + m Right: Ratio with respect to the m 0 t expansion term.The various colours correspond to the inclusion of different expansion terms.