Double Higgs boson production at $e^+e^-$ colliders in the two-Higgs-doublet model

TL;DR

This paper investigates how double Higgs production at electron-positron colliders can be significantly enhanced in the two-Higgs-doublet model due to heavy Higgs resonances, especially outside the alignment limit, with implications for collider experiments.

Contribution

It demonstrates the potential for large cross section enhancements in double Higgs production within the 2HDM, linking these effects to Higgs coupling modifications and experimental constraints.

Findings

Cross section can be enhanced by hundreds of percent due to heavy Higgs resonances.

Strong correlation between cross section enhancement and Higgs coupling scaling factors.

Enhancements are significant outside the alignment limit, constrained by current data.

Abstract

We study the double Higgs boson production processes $e^+e^- \to hh f\bar{f}$ ($f\neq t$) with $h$ being the 125 GeV Higgs boson in the two-Higgs-doublet model with a softly-broken $Z_2$ symmetry. The cross section can be significantly enhanced, typically a few hundreds percent, as compared to the standard model prediction due to resonant effects of heavy neutral Higgs bosons, which becomes important in the case without the alignment limit. We find a strong correlation between the enhancement factor of the cross section and the scaling factor of the $hf\bar{f}$ couplings under constraints from perturbative unitarity, vacuum stability and current experimental data at the LHC as well as the electroweak precision data.

Equations

$$\left(\begin{array}{c}{\Phi }_1\\ {\Phi }_2\end{array}\right)=\left(\begin{array}{cc}\cos \beta & -\sin \beta \\ \sin \beta & \cos \beta \end{array}\right)\left(\begin{array}{c}\Phi \\ {\Phi }^{\prime }\end{array}\right),$$

$$\Phi =\left(\begin{array}{c}{G}^{+}\\ \frac{h_1^{\prime }+v+i{G}^0}{\sqrt{2}}\end{array}\right),{\Phi }^{\prime }=\left(\begin{array}{c}{H}^{+}\\ \frac{h_2^{\prime }+iA}{\sqrt{2}}\end{array}\right).$$

$$\left(\begin{array}{c}{h}_1^{\prime }\\ h_2^{\prime }\end{array}\right)=\left(\begin{array}{cc}{c}_{\beta -\alpha }& s_{\beta -\alpha }\\ -s_{\beta -\alpha }& c_{\beta -\alpha }\end{array}\right)\left(\begin{array}{c}H\\ h\end{array}\right),$$

$${\mathcal{L}}_{\text{kin}}=\sum _{i=1,2}|D_{\mu }{\Phi }_i{|}^2=|D_{\mu }\Phi {|}^2+|D_{\mu }{\Phi }^{\prime }{|}^2,$$

$${\mathcal{L}}_{\text{kin}}\supset \left(\frac{2m_W^2}{v}{W}_{\mu }^{+}{W}^{-\mu }+\frac{m_Z^2}{v}{Z}_{\mu }{Z}^{\mu }\right)(s_{\beta -\alpha }h+c_{\beta -\alpha }H).$$

$${\mathcal{L}}_{\text{kin}}$$

$$-ie{H}^{+}{\overleftrightarrow{\partial }}^{\mu }{H}^{-}{A}_{\mu }-i\frac{g\cos 2{\theta }_W}{2\cos {\theta }_W}{H}^{+}{\overleftrightarrow{\partial }}^{\mu }{H}^{-}{Z}_{\mu },$$

$${\mathcal{L}}_Y$$

$$-{\bar{L}}_L^3\frac{\sqrt{2}{m}_{\tau }}{v}\left(\Phi +{\xi }_{\tau }{\Phi }^{\prime }\right){\tau }_R+\text{h.c.},$$

$$\displaystyle\begin{array}[]{ccccc}&\text{Type-I}&\text{Type-II}&\text{Type-X}&\text{Type-Y}\\ (\xi_{b},\xi_{\tau}):&(\cot\beta,\cot\beta)&(-\tan\beta,-\tan\beta)&(\cot\beta,-\tan\beta)&(-\tan\beta,\cot\beta),\end{array}$$

$${\mathcal{L}}_Y$$

$$-\frac{\sqrt{2}}{v}\bar{t}(m_b{P}_R-m_t{P}_L)H^{+}b-\frac{\sqrt{2}}{v}{\bar{\nu }}_{\tau }{m}_{\tau }{P}_R{H}^{+}\tau +\text{h.c.},$$

$$V$$

$$+\frac{{\lambda }_1}{2}({\Phi }_1^{†}{\Phi }_1{)}^2+\frac{{\lambda }_2}{2}({\Phi }_2^{†}{\Phi }_2{)}^2+{\lambda }_3({\Phi }_1^{†}{\Phi }_1)({\Phi }_2^{†}{\Phi }_2)+{\lambda }_4|{\Phi }_1^{†}{\Phi }_2{|}^2+\frac{{\lambda }_5}{2}[({\Phi }_1^{†}{\Phi }_2{)}^2+\text{h.c.}].$$

$$m_{H^{\pm }},m_A,m_H,m_h,M^2,\tan \beta ,v,\alpha ,$$

$${\lambda }_{hhh}$$

$${\lambda }_{Hhh}$$

$$\displaystyle\hskip 28.45274pt+(2m_{h}^{2}+m_{H}^{2}-3M^{2})[2c^{2}_{\beta-\alpha}+s_{\beta-\alpha}c_{\beta-\alpha}(\tan\beta-\cot\beta)]\Big{\}},$$

$${\kappa }_f\equiv \frac{g_{hff}^{\text{THDM}}}{g_{hff}^{\text{SM}}}=s_{\beta -\alpha }+c_{\beta -\alpha }\cot \beta .$$

$$R\equiv \frac{{\sum }_f{\sigma }^{\text{THDM}}(e^{+}{e}^{-}\to f\bar{f}hh)}{{\sum }_f{\sigma }^{\text{SM}}(e^{+}{e}^{-}\to f\bar{f}hh)},$$

Full text

I Introduction

Various signatures of the discovered Higgs boson at the LHC Aad et al. (2012); Chatrchyan et al. (2012) show that its properties such as observed production cross sections and decay branching ratios are consistent with those predicted in the Standard Model (SM) Aad et al. (2016). Although the SM assumes the minimal form of the Higgs sector composed of one isospin scalar doublet, one may ask him- or herself a natural question: whether the discovered Higgs boson comes from just one doublet or not? In fact, the discovered Higgs boson can be regarded as one of Higgs bosons arising from an extended structure of the Higgs sector, and there is no strong reason to restrict the Higgs sector to be minimal. On the other hand, extended Higgs sectors often appear as a low energy effective theory of physics beyond the SM based upon various physics motivations. The important thing is that phenomenological features in the extended Higgs sectors strongly depend on a specific scenario of the underlying theory. Therefore, as a bottom-up approach, it is important to study the structure of the Higgs sector in order to narrow down new physics models.

Among various possibilities of the extended Higgs sector, a two-Higgs-doublet model (THDM) Branco et al. (2012) is one of the simplest but important examples, as it appears in several new physics models. For example, models proposed to solve the gauge hierarchy problem predict THDMs as their low energy effective theories such as the minimal supersymmetric SM Haber and Kane (1985) and composite Higgs models Mrazek et al. (2011); De Curtis et al. (2018a, b). In addition, extra CP violating phases can arise from the multi-doublet structure of the scalar potential, which are needed to realize the successful scenario for the electroweak baryogenesis Turok and Zadrozny (1990); Nelson et al. (1992). Furthermore, the second Higgs doublet is often introduced in models beyond the SM to explain tiny neutrino masses Zee (1980, 1985); Ma (2006); Campos et al. (2017); Camargo et al. (2018) and dark matter in the Universe Barbieri et al. (2006). For these reasons, we consider the THDM as a reference model, in which we impose a softly-broken discrete $Z_{2}$ symmetry to avoid flavor changing neutral currents at tree level Glashow and Weinberg (1977).

There are basically two ways to test the THDM at collider experiments, namely, the direct searches for additional Higgs bosons such as a charged Higgs boson and the indirect searches finding deviations in properties of the discovered Higgs boson ($h$) from the SM prediction. Concerning the former way, if we discover the additional Higgs bosons, it turns out to be a direct evidence for the THDM or at least an extended Higgs sector, but no report has hitherto been provided for the discovery of such new particles at the LHC. Recent studies about the prospects of the direct searches at the high-luminosity LHC (HL-LHC) include Refs. Kling et al. (2018); Adhikary et al. (2018), where possibilities to find new Higgs bosons via cascade decays of heavy scalar resonances are discussed in detail.

Recently the latter way, seeking deviations as the indirect searches, is getting much attention, since $h$ has already been discovered and its properties will be precisely measured in near future. For example, at the HL-LHC, the $h$ couplings to weak gauge bosons ($hVV$, $V=W,Z$) and fermions (e.g., $h\tau\tau$, $hb\bar{b}$) are expected to be measured with a few percent level Dawson et al. (2013), while at the International Linear Collider (ILC), $hVV$ and $hf\bar{f}$ ($f=\tau,~{}c$ and $b$) couplings can be measured to be sub-percent and one percent level, respectively Fujii et al. (2017). In the THDM, these $h$ couplings can deviate from the SM prediction with various patterns depending on a particular scenario as it has been clarified in Ref. Kanemura et al. (2014a) at tree level and in Refs. Arhrib et al. (2004); Kanemura et al. (2014b, 2017, 2015, 2018a) at one-loop level. Therefore, by looking at the possible deviation in the $h$ couplings in the future collider experiments, one can distinguish the scenarios of the THDM, which is still important even if an additional Higgs boson could be discovered via a resonant peak in the direct searches.

In this paper, we focus on the cross section of double Higgs boson production processes at future $e^{+}e^{-}$ colliders as another important observable regarding the indirect search of the THDM. The double Higgs boson production has been discussed in the $pp$ Asakawa et al. (2010); Grober et al. (2017); Basler et al. (2018), $e^{+}e^{-}$ Asakawa et al. (2010); Sonmez (2018) and $\gamma\gamma$ Asakawa et al. (2010) collision in the THDM to extract the Higgs boson self-coupling constant $hhh$, particularly the case with the alignment limit, where all the SM-like Higgs boson couplings become the same as those in the SM prediction. In Ref. Arhrib et al. (2009), the Higgs boson pair production has also been discussed at the LHC in the case with and without the alignment limit. Experimentally, it is a formidable task for the LHC to precisely determine the Higgs trilinear coupling because of the tiny cross section of the double Higgs boson production (see collaboration (2018); Sirunyan et al. (2018) for recent results). This is certainly one of the motivations of construction of new powerful colliders: the HL-LHC will be sensitive only to $\mathcal{O}(100\%)$ for the trilinear coupling while the ILC running at 500 GeV can measure it with an about 25% uncertainty at the 68% confidence level Di Vita et al. (2018).

Our motivation to discuss the double Higgs production is to find the correlation between the deviation in the Higgs boson couplings and the modification of the cross section. These two variables are expected to be strongly correlated with each other, because the deviation in the $h$ couplings appears in the case without the alignment limit at tree level, in which case the additional neutral Higgs bosons can mediate the double Higgs boson production process and can provide sizable enhancement of the cross section. For the LHC, large enhancement of the di-Higgs production cross section due to resonantly produced heavy Higgs bosons in the CP violating THDM has been studied in Ref. Grober et al. (2017). In the present paper, we clarify how large enhancement can be obtained at $e^{+}e^{-}$ colliders in the case without the alignment limit. In the numerical analysis, we focus on a special case of the THDM often referred to as Type-I, because scenarios without the alignment limit are highly constrained by the Higgs boson signal strengths in the other types of Yukawa interactions such as Type-II. Under the current theoretical and experimental constraints on the parameter space, we find a strong correlation between the enhancement of the cross section and the scaling factor of the $hf\bar{f}$ couplings. We note that Higgs boson pair production at $e^{+}e^{-}$ colliders requires the collision energy $\sqrt{s}$ being larger than 250 GeV while the scaling factor can be precisely determined with $\sqrt{s}=250~{}\text{GeV}$. Hence, the Higgs boson couplings will have been known when experiments reach to the collision energy enough for Higgs pair boson production measurement.

This paper is organized as follows. In Sec. II, we define the model and present the relevant Higgs boson interactions. In Sec. III, we discuss general properties of the Higgs boson pair production in $e^{+}e^{-}$ colliders. Then, we explain a mechanism that enhances the cross section in the THDM without the alignment limit and give a rough estimate of the enhancement factor. Detailed numerical analysis on the aforementioned processes at tree level is performed in Sec. IV. Finally, Sec. V is devoted to our conclusions.

II Model

We consider the THDM whose Higgs sector is constructed by two isospin scalar doublet fields $\Phi_{1}$ and $\Phi_{2}$. For simplicity, we consider the CP-conserving case throughout the paper. The vacuum expectation values (VEVs) of the two doublets, i.e., $\sqrt{2}\langle\Phi_{1,2}^{0}\rangle=v_{1,2}$ are parameterized by $v=\sqrt{v_{1}^{2}+v_{2}^{2}}=(\sqrt{2}G_{F})^{-1/2}$ with $G_{F}$ being the Fermi constant and their ratio $\tan\beta=v_{2}/v_{1}$.

It is convenient to introduce the so-called Higgs basis Davidson and Haber (2005); Georgi and Nanopoulos (1979) as follows:

[TABLE]

where

[TABLE]

In this basis, the Nambu-Goldstone bosons $G^{\pm}$ and $G^{0}$, which are absorbed into the longitudinal components of the $W$ and $Z$ bosons, are separated from the physical singly-charged Higgs boson $H^{\pm}$ and the CP-odd Higgs boson $A$, while two CP-even Higgs states $h_{1}^{\prime}$ and $h_{2}^{\prime}$ are generally not mass eigenstates at this stage. By introducing another mixing angle $\alpha$, the mass eigenstates are given by

[TABLE]

where we have abbreviated $\sin X$ and $\cos X$ as $s_{X}$ and $c_{X}$, respectively. We identify $h$ as the discovered Higgs boson with a mass of about 125 GeV.

The kinetic terms for the Higgs doublets are given as

[TABLE]

where $D_{\mu}$ is the covariant derivative for the Higgs doublets. We note that from the first term on the right-hand side, the gauge-gauge-Higgs type interactions are obtained as

[TABLE]

From the second term in the right-hand side of Eq. (4), we obtain the Higgs-Higgs-gauge type interactions as

[TABLE]

with $X\overleftrightarrow{\partial}^{\mu}Y\equiv X(\partial^{\mu}Y)-(\partial^{\mu}X)Y$ and $\theta_{W}$ being the weak mixing angle.

In order to avoid Higgs boson mediating flavour changing neutral currents at tree level, we impose a discrete $Z_{2}$ symmetry Glashow and Weinberg (1977) into the Higgs sector, where the two doublets are transformed as $\Phi_{1}\to+\Phi_{1}$ and $\Phi_{2}\to-\Phi_{2}$. Under the $Z_{2}$ symmetry, only one of the two Higgs doublets can couple to each up-type, down-type quarks and charged leptons, by which the interaction matrices in the flavour space between neutral Higgs bosons and fermions are diagonalized in the fermion mass eigenbasis. It has been known that there are four independent types of Yukawa interactions so-called Type-I, -II, -X and -Y Aoki et al. (2009) depending on the way to assign the $Z_{2}$ charge for fermions Barger et al. (1990); Grossman (1994).

The Yukawa Lagrangian for the third generation fermions is then given in the Higgs basis by

[TABLE]

where $\tilde{\Phi}=i\tau_{2}\Phi^{*}$ and $\tilde{\Phi}^{\prime}=i\tau_{2}\Phi^{\prime*}$. The factors $\xi_{b}$ and $\xi_{\tau}$ depend on the choice of the types of Yukawa interaction as

[TABLE]

while $\xi_{t}=\cot\beta$ for all the types of Yukawa interaction. The interaction terms of the Higgs bosons and fermions are extracted as

[TABLE]

where $I_{t}\,(I_{b,\tau})=1/2\,(-1/2)$ and $P_{L}\,(P_{R})$ is the projection operator for the left (right) hand chirality.

The Higgs potential is generally written by eight independent real parameters when we include the soft-breaking term of the $Z_{2}$ symmetry:

[TABLE]

As we already mentioned in the above, we assume the CP-invariance of the Higgs sector, so that $m_{3}^{2}$ and $\lambda_{5}$ parameters are taken to real. After imposing the tadpole conditions, i.e., the requirement of vanishing the linear terms of $h$ and $H$ in the Higgs potential, we can eliminate the $m_{1}^{2}$ and $m_{2}^{2}$ parameters. Then, these eight parameters, i.e., six parameters in the potential and two VEVs are expressed as follows:

[TABLE]

where $M^{2}\equiv m_{3}^{2}/(s_{\beta}c_{\beta})$. Among these eight parameters, $v$ and $m_{h}$ are fixed to about 246 GeV and 125 GeV, respectively.

From Eq. (12), we can extract the scalar three-point couplings. In particular, the $hhh$ and $Hhh$ couplings will be important in the later analysis, which are expressed in terms of the parameters shown in (13) as follows:

[TABLE]

where the above quantities are defined by the coefficient in front of the $hhh$ and $Hhh$ vertices in the Lagrangian.

Before closing this section, it would be worth to mention that taking $s_{\beta-\alpha}\to 1$, all the $hVV$ ($V=W,Z$), $hf\bar{f}$ and $hhh$ couplings become the same values as in the SM at tree level. This limit has been known as the alignment limit, where the SM-like Higgs boson $h$ completely comes from the doublet $\Phi$ in the Higgs bases in Eq. (2).

III Double Higgs boson production

In this section, we discuss the general property of double Higgs boson production processes at $e^{+}e^{-}$ colliders, and consider how much the production cross section can be different in the THDM with respect to that in the SM prediction.

The double Higgs boson production $e^{+}e^{-}\to hh+X$ is possible when the collision energy of the electron and positron $\sqrt{s}$ is larger than about 250 GeV. Figs. 1 and 2 show relevant Feynman diagrams, where those in Fig. 1 are common to the SM and THDM, while those in Fig. 2 only appear in the THDM. We here ignore diagrams induced via Yukawa couplings222For the case of $\sqrt{s}$ larger than $2(m_{t}+m_{h})\sim 600$ GeV, the $t\bar{t}hh$ production is allowed and its contribution could be comparable as compared with that shown in Fig. 1. In the numerical analysis of this paper given in the next section, we focus on the collision energy to be 500 GeV, so that this process is not needed to be considered. . In the SM, there are two types of diagrams contributing to the double Higgs boson production, namely, the $s$-channel process $e^{+}e^{-}\to Zhh$ shown as diagrams (a)–(c) in Fig. 1 and the vector boson fusion processes $e^{+}e^{-}\to e^{+}e^{-}hh$ and $e^{+}e^{-}\to\nu_{e}\bar{\nu}_{e}hh$ shown as diagrams (d)–(f) in Fig. 1. The cross section of these processes are calculated as 0.16 (0.12), $2.6\times 10^{-3}$ ($7.2\times 10^{-2}$) and $3.6\times 10^{-4}$ ($9.5\times 10^{-3}$) fb for the $e^{+}e^{-}\to Zhh$, $\nu_{e}\bar{\nu}_{e}hh$ and $e^{+}e^{-}hh$ process at $\sqrt{s}=500\ (1000)$ GeV, respectively. Therefore, the cross section of the double Higgs boson production is mainly determined by the $s$-channel process.

In the THDM, additional diagrams contribute to the process as seen in Fig. 2, where the extra neutral Higgs boson $H$ or $A$ appears in the diagrams having the same topology as those of (a), (b), (d) and (e) in Fig. 1. In order to estimate the typical size of these contributions, let us focus on the diagram (a’) in Fig. 2 as an example. When the mass of $H$ is taken to between $250~{}\text{GeV}\lesssim m_{H}\lesssim\sqrt{s}-m_{Z}$, $H$ can be on-shell. In this case, the cross section of the diagram (a’) is approximately calculated by the product of the two-body cross section of $e^{+}e^{-}\to ZH$ and the branching ratio of the $H\to hh$ decay assuming $\Gamma_{H}\ll m_{H}$ with $\Gamma_{H}$ being the total width of $H$.333In our scenario, the typical size of $\Gamma_{H}/m_{H}$ is 1% level or smaller. Therefore, the size of this cross section is typically obtained by multiplying the factor of $16\pi^{2}\times c^{2}_{\beta-\alpha}\times\text{BR}(H\to hh)$ with respect to the cross section of the diagram (a) in the SM, where the factor $16\pi^{2}$ appears due to the typical ratio of the two-body and three-body phase-space factors, while $c^{2}_{\beta-\alpha}$ comes from the $HZZ$ coupling normalized to the $hZZ$ one in the SM. For example, when $s_{\beta-\alpha}$ is fixed to be 0.99 (0.995), the above factor becomes 3.14 (1.56) assuming BR$(H\to hh)=1$. The similar enhancement can also be obtained from the diagram (b’) as long as $A$ is produced with on-shell. We thus can expect that the total cross section of the double Higgs boson production can be several times larger than the SM prediction.

The important thing here is that such enhancement of the cross section happens when departure of the alignment limit, i.e., $s_{\beta-\alpha}\neq 1$, is realized, because the both $HZZ$ and $AZh$ couplings are proportional to $c_{\beta-\alpha}$. Therefore, the enhancement is strongly correlated with the deviation in the Higgs boson couplings from the SM prediction. On the other hand, the $hVV$ ($V=W,Z$) and $hf\bar{f}$ couplings are expected to be precisely measured at future $e^{+}e^{-}$ colliders. For example, at the ILC the $hZZ$, $hWW$, $hb\bar{b}$, $h\tau^{+}\tau^{-}$ and $hc\bar{c}$ couplings may be able to be measured with $0.38$, $1.8$, $1.8$, $1.9$ and $2.4$% at 1$\sigma$ level assuming 250 GeV of the collision energy and 2 ab*-1* of the integrated luminosity Fujii et al. (2017). Therefore, if deviations in the Higgs boson couplings are detected in future, we expect the sizable enhancement of the double Higgs boson production as well.

In the next section, we numerically evaluate how large enhancement can be obtained in the scenario without the alignment limit and show the correlation between the double Higgs boson production cross section and the deviation in the Higgs boson couplings.

IV Numerical results

In this section, we compute the cross section of the $e^{+}e^{-}\to f\bar{f}hh$ process, where $f$ is a fermion except for the top quark. These four-body final states are obtained through the $Zhh$ production with the decay of $Z$ into a fermion pair and the vector boson fusion process, see Figs. 1 and 2. We note that when the $Z$ boson from the $Zhh$ process decays into $e^{+}e^{-}$ ($\nu_{e}\bar{\nu}_{e}$), this process interferes with the $e^{+}e^{-}hh$ $(\nu_{e}\bar{\nu}_{e}hh$) final states from the $Z$ ($W$) boson fusion process. We take into account such interference effects in the numerical analysis. We also note that the dependence of the type of Yukawa interactions slightly appears in the cross section through the decay width of $H$ and $A$. In contrast, the type dependence significantly appears in the region of the parameter space allowed from experimental bounds. In particular, except for the Type-I THDM, the scenario with $s_{\beta-\alpha}\neq 1$ which is considered in this section is highly constrained by the Higgs boson signal strengths, see e.g., Blasi et al. (2017). We thus consider the Type-I THDM in what follows444In the Type-I THDM, constraints from flavor experiments such as those from $B\to X_{s}\gamma$ data on the charged Higgs boson mass are also quite milder than those in the other types especially in the Type-II. For example, ${\cal O}(100)$ GeV of the charged Higgs boson mass is allowed by the $B\to X_{s}\gamma$ data when $\tan\beta\gtrsim 2$ Misiak and Steinhauser (2017) in the Type-I THDM, while $m_{H^{\pm}}\lesssim 600$ GeV is excluded with 95% confidence level in the Type-II THDM. .

As seen in (13), there are six free parameters in the THDM, i.e., three masses of extra Higgs bosons, $M^{2}$, $\tan\beta$ and $\alpha$. Instead of inputting $\alpha$, we take $s_{\beta-\alpha}$ and the sign of $c_{\beta-\alpha}$ as inputs. To avoid large contributions to the electroweak oblique $T$ parameter Peskin and Takeuchi (1990, 1992), we take the masses of $H^{\pm}$ and $A$ to be the same, i.e., $m_{A}=m_{H^{\pm}}$ Deshpande and Ma (1978); Sher (1989); Nie and Sher (1999); Kanemura et al. (1999). In this case, the quadratic dependence of the extra Higgs boson masses completely vanishes, and only the small logarithmic mass dependence remains, which is proportional to $c_{\beta-\alpha}^{2}$.

Before going to show the numerical results for the cross section, let us summarize the constraints on the parameter space what we take into account in the analysis. We impose the perturbative unitarity bound Kanemura et al. (1993); Akeroyd et al. (2000); Ginzburg and Ivanov (2003); Kanemura and Yagyu (2015) and the vacuum stability bound Deshpande and Ma (1978); Sher (1989); Nie and Sher (1999); Kanemura et al. (1999) as the constraints from theoretical consistency. These bounds restrict the size of scalar quartic couplings in the potential, which can be translated into the bound on the Higgs boson masses and the mixing angles. As the experimental constraints, we take into account the electroweak oblique $S$ and $T$ parameters Tanabashi et al. (2018), direct searches for additional Higgs bosons at the LEP, Tevatron and LHC experiments as well as the compatibility of the signal strengths for the discovered Higgs boson with a mass of 125 GeV. For the direct searches, we use the HiggsBounds-5.3.0beta Bechtle et al. (2014). For the signal strengths of the discovered Higgs boson, we use the combined data from ATLAS and CMS at the LHC Run-I experiments Aad et al. (2016), and require that the prediction of the signal strengths for the $\gamma\gamma$, $ZZ^{*}$, $WW^{*}$ and $\tau\tau$ modes of the $h$ state does not exceed the given width of the error bar. We impose these experimental bounds at the 95% confidence level.

In Fig. 3, we show the region of the parameter space excluded by the theoretical and experimental bounds explained in the above. Here, we take $m_{H}=m_{A}(=m_{H^{\pm}})=300$ GeV and $s_{\beta-\alpha}=0.99$ (left panel) and 0.995 (right panel). The other parameters $M^{2}$ and $\tan\beta$ are scanned in these plots. Instead of showing the value of $\tan\beta$, we introduce the scaling factor $\kappa_{f}$ for the Yukawa coupling $g_{hff}^{\text{THDM}}$ in the THDM normalized to the SM value $g_{hff}^{\text{SM}}$, see Eq. (7):

[TABLE]

We note that due to the choice of the Type-I THDM, the scaling factor $\kappa_{f}$ does not depend on the choice of a fermion $f$, see Sec. II. From Eq. (16), we see that $\kappa_{f}<s_{\beta-\alpha}$ ($\kappa_{f}>s_{\beta-\alpha}$) is obtained by taking the sign of $c_{\beta-\alpha}$ to be negative (positive), and $\kappa_{f}=s_{\beta-\alpha}$ is given in the limit of $\tan\beta\to\infty$.

We see that the constraint from experiments (shown by the magenta shaded region) becomes important in the region with smaller $\tan\beta$ values, i.e., larger values of $|1-\kappa_{f}|$, where the direct search at the LHC, particularly for $H\to ZZ$ dominantly contributes to the exclusion. This can be understood by the fact that the production cross section $gg\to H$ is proportional to $\cot^{2}\beta$ in the limit of $s_{\beta-\alpha}\to 1$, so that the larger $\tan\beta$ case can avoid the constraint due to the smaller cross section. Another thing we can find in this figure is that the theory bound (shown by the blue shaded region555In fact, the constraints from the $S$ and $T$ parameter are also taken into account in the blue shaded region, but the shape of the exclusion is almost determined by the perturbative unitarity and vacuum stability bounds. ) becomes important in the case with larger $\tan\beta$ and/or larger difference between $m_{H}^{2}$ and $M^{2}$. In particular, for a larger $\tan\beta$ case, quite small area with $M^{2}\sim m_{H}^{2}$ is allowed, because some of the scalar quartic couplings become very sensitive to the value of $m_{H}^{2}-M^{2}$. The typical behavior of the constraints does not change so much between the cases with $s_{\beta-\alpha}=0.99$ (left panel) and 0.995 (right panel), but smaller values of $|1-\kappa_{f}|$ are excluded by experimental constraints. This is simply because the value of $|1-\kappa_{f}|$ given by the same value of $\tan\beta$ becomes smaller in the case of $s_{\beta-\alpha}=0.995$ as compared to the case of $s_{\beta-\alpha}=0.99$, as seen in Eq. (16). For reference, we also show the contours of the branching ratio of the $H\to hh$ mode which is important to understand the behavior of the enhancement of the cross section of the double Higgs boson production.

Now we present numerical results on the cross section. For the numerical computation of the cross section, we used the public version of GRACE Yuasa et al. (2000); Fujimoto et al. (2003); Collaboration with some modifications, and all the calculations were performed at tree level with $\sqrt{s}=500\ \text{GeV}$. Let us first illustrate the invariant mass distribution of the di-Higgs system $m_{hh}$ in Fig. 4. Here, we take $m_{H}=300$, 350 and 400 GeV as examples, and the other parameters are specified as described in the caption. It is clearly seen that the sharp peaks appear at around $m_{H}$, because of the on-shell $H$ mediation with the $H\to hh$ decay shown in the diagram ($a^{\prime}$) of Fig. 2. The height of the peak is getting lower as the mass of $H$ increases. Due to this resonant effect, the cross section of the double Higgs boson production is sizably enhanced as we will see in what follows.

Next, we consider the ratio of the di-Higgs production cross section in the Type-I THDM to that in the SM in order to see how large enhancement can be obtained:

[TABLE]

where the summation for $f$ is done over all the fermions except for $f=t$. In the following calculation, we scan $\tan\beta$ and $M^{2}$ with the ranges of $1\leq\tan\beta\leq 30$ and $0\leq M^{2}\leq(500~{}\text{GeV})^{2}$. The scatter plot in Fig. 5 shows the ratio $R$ as functions of the scaling factor $\kappa_{f}$ with using sampling points passing all the constraints, namely, points in the white region of Fig. 3. Red and blue points represent for $s_{\beta-\alpha}=0.99$ and 0.995, respectively. As it is expected in Sec. III, a clear correlation between $R$ and $\kappa_{f}$ is seen from this plot, and a considerable enhancement of the cross section is observed due to the on-shell mediation of the extra neutral Higgs bosons $H$ and $A$. See also the contours of the branching ratio of $H\to hh$ shown in Fig. 3 to figure out the behavior of the correlation. It is also seen that a larger value of $R$ is obtained in the case of $s_{\beta-\alpha}=0.99$ in comparison with the case of $s_{\beta-\alpha}=0.995$, because the $HZZ$ and $AZh$ couplings are proportional to $c_{\beta-\alpha}$. The strong enhancement seen in the large $\tan\beta$ region (i.e., $\kappa_{f}\sim s_{\beta-\alpha}\sim 1$) can be traced back to a large value of the $\lambda_{Hhh}$ coupling, see Eq. (15). This gives a considerable deviation $R>1$ even at $\kappa_{f}=1$. We find that the value of $R$ can be maximally around 5 (3) for $s_{\beta-\alpha}=0.99$ (0.995).

We also perform the similar calculation shown in Figs. 7 and 7, but for the case with $m_{H}=300$ GeV and $m_{A}=m_{H^{\pm}}=500$ GeV. The region allowed by the constraints is almost the same as the previous plots in Fig. 3. On the other hand, by looking at Fig. 7, we find that the shape of the scattering points are shifted to below due to the $A$ mediation being off-shell, and the maximally allowed value of $R$ becomes around 4 and 2.5 for $s_{\beta-\alpha}=0.99$ and 0.995, respectively.

In Fig. 9, we scan the value of $s_{\beta-\alpha}$ with $0.99\leq s_{\beta-\alpha}\leq 1$ under the constraints described in the above. The left (right) panel shows the case with $m_{A}=m_{H^{\pm}}=300$ (500) GeV, while $m_{H}$ is fixed to be 300 GeV in the both panels. The value of $R$ is distributed in the range between 1 and 4.5 (4) for the case with $m_{A}=m_{H^{\pm}}=300$ (500) GeV. In Fig. 9, we also show the dependence on the scaling factor of $hVV$ ($V=W,Z$) couplings denoted as $\kappa_{V}\equiv g_{hVV}^{\text{THDM}}/g_{hVV}^{\text{SM}}=s_{\beta-\alpha}$. It is clarified that the value of $R$ approaches to unity when $\kappa_{V}\to 1$, because $H$ and $A$ do not appear in the diagram in this limit and all the couplings of $h$ become the same value as those in the SM. However, once $\kappa_{V}<1$ is taken, a significant enhancement is driven by the appearance of the $H$ and $A$ mediations in the diagram depending on the value of $\kappa_{f}$, see also Figs. 5 and 7.

Finally, we show the mass dependence of the value of $R$ in Fig. 10 with a fixed value of $\tan\beta$ to be 2 (top), 3 (middle) and 5 (bottom). In these plots, we take the degenerate mass of the extra Higgs bosons, i.e., $m_{H}=m_{A}=m_{H^{\pm}}$, and varying the mass range to be from 200 GeV to 500 GeV. Again, all the points are passed all the constraints explained in the above. For the low $\tan\beta$ case such as $\tan\beta=2$ shown in the top panels, the constraint from direct searches excludes the lower mass region, so that the points appear only in the larger mass region. For the larger $\tan\beta$ case particularly with $c_{\beta-\alpha}>0$ (left panels), the lower mass region is allowed, and it is clearly seen that the value of $R$ drastically increases at the mass of extra Higgs boson just above 215 and 250 GeV, because of the on-shell $A\to Zh$ and $H\to hh$ decays open. We also find that for the larger mass region, e.g., $m_{H}\gtrsim 400$ GeV the value of $R$ becomes smaller than 1, in which the $H$ and $A$ appearing in the $Zhh$ production are getting off-shell666In fact, for the case with $m_{H}\sim 400$ GeV, $H$ can still be on-shell as it is produced with the $Z$ boson at the collision energy of 500 GeV. However, we have checked that the cross section of the $e^{+}e^{-}\to ZH\to Zhh$ with $m_{H}=400$ GeV is smaller than the other subdominant contributions such as the $e^{+}e^{-}\to\nu_{e}\bar{\nu}_{e}hh$ process, because of the phase space suppression. Therefore, we cannot obtain the large enhancement of the cross section at $m_{H}\sim 400$ GeV or larger. , so that these contributions become unimportant. The value $R<1$ is simply explained by the fact that the contribution containing $h^{*}\to hh$ as the diagram (a) in Fig. 1 becomes smaller than the SM prediction because of the smaller $\lambda_{hhh}$ coupling with respect to the SM value, see Eq. (14). We note that in these plots, the value of $\tan\beta$ is fixed, so that the value of $\kappa_{f}$ is determined to be (1.06, 1.04, 1.02; 0.92, 0.94, 0.96) for $s_{\beta-\alpha}=0.99$ and (1.04, 1.03, 1.01; 0.95, 0.96, 0.98) for $s_{\beta-\alpha}=0.995$ where the first (second) three values are the case for $c_{\beta-\alpha}>0~{}(c_{\beta-\alpha}<0)$ with $\tan\beta=2$, $3$ and 5.

Before closing this section, we would like to give a brief comment on the difference of the prediction from the other models, e.g., models with an isospin singlet scalar field, see for instance Lewis and Sullivan (2017). In these models, both $\kappa_{V}\,(V=W,Z)$ and $\kappa_{f}$ are given by the same factor like $\cos\theta$, where $\theta$ is a mixing angle between 2 CP-even Higgs bosons. Therefore, the case ($\kappa_{f}\sim 1$, but $R>1$) what we found in the above does not happen, and so our analysis would give an important result to pin down the possible models.

In addition, we would also like to comment on the possibility of the direct detection of $H$ at future collider experiments such as the HL-LHC. In the parameter space what we considered in this paper, the branching ratio of $H$ is approximately given by $\text{BR}(H\to hh)+\text{BR}(H\to WW)+\text{BR}(H\to ZZ)\sim 1$, where the small portion of the branching ratio is filled by the fermionic decay modes777For $m_{H}\gtrsim 350$ GeV, the branching ratio of the $H\to t\bar{t}$ mode can be 10% level for a smaller value of $\tan\beta$. . Thus, $H$ can be discovered via the $H\to WW/ZZ$ mode. In fact at the HL-LHC, there is a study for direct searches of a neutral scalar boson $X$ using the $X\to ZZ$ mode, and the 95% confidence level upper limit on the cross section ($pp\to X$) times the branching ratio $(X\to ZZ)$ is taken in Ref. Cepeda et al. (2019). This bound can be translated into the constraint on the parameter space in our model, by which smaller values of $\tan\beta$ and/or $m_{H}$ can be excluded. Combining such constraint into our analysis, more restricted results would be obtained.

V Conclusions

We have discussed the correlation between the scaling factor of the Yukawa coupling for the SM like Higgs boson $\kappa_{f}$ and the ratio of the cross section for the $e^{+}e^{-}\to hhf\bar{f}$ ($f\neq t$) process normalized to the SM prediction in the Type-I THDM. We particularly concentrate on the case without the alignment limit, in which resonant effects of the extra neutral Higgs bosons $H$ and $A$ provide a sizable enhancement of the cross section and the value of $\kappa_{f}$ is different from unity at tree level. Under the constraints from perturbative unitarity, vacuum stability, electroweak oblique parameters, direct searches for heavy Higgs bosons at collider experiments and compatibility of the signal strengths of the discovered Higgs boson, we have found that the considerable enhancement of the cross section, typically a few times larger than the SM prediction, can be obtained depending on the value of $\kappa_{f}$ and the masses of the extra Higgs bosons. The value of $\kappa_{f}$ is expected to be precisely measured at future collider experiments such as the high-luminosity LHC and the ILC, typically with a few percent and one percent level, respectively. Therefore, if some deviations in the Higgs boson couplings are found at future colliders, we expect the sizable enhancement of the double Higgs boson production and can extract information of the mass of the extra neutral Higgs boson and/or dimensionful model parameter $M^{2}$ in the Higgs potential.

Although expected not to give so drastic change, radiative corrections to the cross-section enhancement studied in this work may be investigated by using, for example, H-COUP Kanemura et al. (2018b) to incorporate one-loop electroweak-corrected vertices or GRACE-loop Belanger et al. (2006); Fujimoto et al. (2007) for full one-loop computation. We leave it for future works.

Acknowledgements.

All the authors are grateful for fruitful discussions with Masaaki Kuroda, Yoshimasa Kurihara, Kiyoshi Kato, Masato Jimbo, Tadashi Ishikawa, Junpei Fujimoto, Yusaku Kouda and Ryotaro Nara.