Echoes of 2HDM inflation at the collider experiments

Abstract

We study the correlation between the constraints on general two Higgs doublet model from Higgs inflation and from collider experiments. The parameter space receives meaningful constraints from direct searches at the large hadron collider and from flavor physics if $m_H$, $m_A$, and $m_{H^\pm }$ are in the sub-TeV range, where H, A, and $H^\pm$ are the CP even, CP odd, and charged Higgs bosons, respectively. We find that in the parameter region favored by the Higgs inflation, H, A, and $H^\pm$ are nearly degenerate in mass. We show that such near degeneracy can be probed directly in the upcoming runs of the Large Hadron Collider, while the future lepton colliders such as the International Linear Collider and the future circular collider would provide complementary probes.

1 Introduction

The cosmic inflation [1,2,3] in the early universe is a well established paradigm which can successfully explain the horizon, flatness and exotic-relics problems, and can provide the initial condition for the hot big bang as the reheating process in the early Universe [4]. The slow-roll inflation [5,6,7] can seed the primordial density fluctuations [8, 9] which eventually evolve into large scale structure that we observe today in cosmic microwave background (CMB) anisotropies [10].

Despite of its prevalent success, the underlying mechanism behind the inflationary dynamics still remains unknown. In the simplest inflationary scenario a slowly rolling scalar field (inflaton) can account for the nearly scale-invariant density fluctuation observed in the CMB. In the standard model (SM), the only available scalar field is the Higgs boson, which has a quartic potential. However, it alone, when used in the chaotic inflation, cannot support the observed scalar spectral index and tensor-to-scalar ratio [10].

The Higgs inflation [11,12,13,14,15,16,17,18]; see also [19] is one of the best fit models to the CMB data, and is testable due to its connection to the Higgs physics at the large hadron collider (LHC) and beyond. In the SM Higgs inflation, the Higgs doublet $\varPhi$ is assumed to couple with gravity via Ricci scalar R by $\xi \varPhi ^\dagger \varPhi R$, where $\xi$ is a dimensionless nonminimal coupling of order $10^4$–$10^5$. The successful Higgs inflation requires the stability of Higgs potential up to at least $M_\text{ P }/\xi$. Even if we demand the stability up to $M_\text{ P }$, the required upper bound on the pole mass of the top quark is $m_t^\text{ pole }\lesssim 171.4$ GeV [22], which is perfectly consistent at 1.4 $\sigma$ with the current value $172.4\pm 0.7$ GeV [23].

The Higgs inflation is also possible in models with additional Higgs doublet. After the discovery of the Higgs boson h of mass 125 GeV [20, 21], it is conceivable that the Higgs field has an extra generation since all the known fermions in the SM has more than one generations. The general two Higgs doublet model (g2HDM) is one of the simplest renormalizable extensions of the SM where the scalar sector ($\varPhi$) is extended by one extra doublet ($\varPhi '$). The g2HDM would share the same virtue of being one of the best fit inflationary models to account for the CMB data if one has sufficiently large nonminimal couplings to $\varPhi$ and/or $\varPhi '$.

In this article we study the possibility of slow-roll inflation with nonminimal Higgs couplings in general two Higgs doublet model¹ and its implications at the collider experiments. In general we have three nonminimal couplings between the Higgs fields and the Ricci scalar in g2HDM. As a first step, we study two different scenarios in this article. In Scenario-I we switch only on the nonminimal coupling of $\varPhi$, while in Scenario-II we switch only on that of $\varPhi '$. In both scenarios we find the parameter space for inflation satisfying all observational constraints from Planck 2018 [10].

Without the presence of discrete symmetry, in g2HDM, at tree level both the scalar doublets couple with both the up- and down-type fermions. After diagonalizing the fermion mass matrices two independent Yukawa couplings $\lambda ^F_{ij}$ and $\rho ^F_{ij}$ emerge, where F denotes leptons (L), up-type quarks (U), and down-type quarks (D): the $\lambda ^F_{ij}$ matrices are real and diagonal and responsible for mass generation of the fermions, while the $\rho ^F_{ij}$ are in general complex and non-diagonal matrices. The parameter space for inflation receives constraints from several direct and indirect searches, in particular from the LHC and Belle experiments. We show that extra Yukawa couplings $\rho ^U_{tt}$ and $\rho ^U_{tc}$ can provide unique test for the parameter space for inflation at the LHC. Discoveries are possible at the LHC or future lepton colliders such as International Linear Collider (ILC) and the future circular collider (FCC-ee), depending on the magnitude of extra Yukawa couplings $\rho ^U_{tt}$ and $\rho ^U_{tc}$. We also show that $B_{s}$ and $B_{d}$ mixing data as well as future measurements of B meson decay observables would provide sensitive probe to the inflationary parameter space.

In the following we outline the g2HDM framework in Sect. 2 followed by formalism for inflation in Sect. 3. The scanning and parameter space for inflation is summarized in Sect. 4, and direct and indirect constraints are discussed in and Sect. 5. We discuss our results with some outlook in Sect. 6.

2 Model framework

Here we outline the framework of the g2HDM following the notation of Refs. [32, 33]. In the Higgs basis, the most general two Higgs doublet potential can be written as [33], see e.g., [34]

$$\begin{aligned}&V(\varPhi ,\varPhi ') = \mu _{11}^2|\varPhi |^2 + \mu _{22}^2|\varPhi '|^2 - \left( \mu _{12}^2\varPhi ^\dagger \varPhi ' + h.c.\right) \nonumber \\&\qquad + \frac{\eta _1}{2}|\varPhi |^4 + \frac{\eta _2}{2}|\varPhi '|^4 + \eta _3|\varPhi |^2|\varPhi '|^2 + \eta _4 |\varPhi ^\dagger \varPhi '|^2 \nonumber \\&\qquad + \left[ \frac{\eta _5}{2}(\varPhi ^\dagger \varPhi ')^2 + \left( \eta _6 |\varPhi |^2 + \eta _7|\varPhi '|^2\right) \varPhi ^\dagger \varPhi ' + h.c.\right] , \end{aligned}$$

(1)

where the vacuum expectation value v arises from the doublet $\varPhi$ via the minimization condition $\mu _{11}^2=-\frac{1}{2}\eta _1 v^2$, while $\left\langle \varPhi \right\rangle =(0,~v/\sqrt{2})^T$, $\left\langle \varPhi '\right\rangle =0$ (hence $\mu _{22}^2 > 0$), and $\eta _i$s are quartic couplings. A second minimization condition, $\mu _{12}^2 = \frac{1}{2}\eta _6 v^2$, removes $\mu _{12}^2$, and the total number of parameters are reduced to nine. For the sake of simplicity, we assumed CP-conserving Higgs sector. The mixing angle $\gamma$ between the CP even scalars h, H satisfy relations:

$$\begin{aligned} \cos \gamma ^2 = \frac{\eta _1 v^2 - m_h^2}{m_H^2-m_h^2},~\quad \sin {2\gamma } = \frac{2\eta _6 v^2}{m_H^2-m_h^2}. \end{aligned}$$

(2)

The alignment limit corresponds to $c_\gamma \rightarrow 0$ with $s_\gamma \rightarrow -1$, where we used shorthand $c_\gamma = \cos \gamma$ and $s_\gamma = \sin \gamma$ . The current LHC data suggests [35, 36] that $c_\gamma$ to be small i.e. the so called approximate alignment [32].

The physical scalar masses can be expressed in terms of the parameters in Eq. (1),

$$\begin{aligned} m_{A}^2&= \frac{1}{2}(\eta _3 + \eta _4 - \eta _5) v^2+ \mu _{22}^2, \end{aligned}$$

(3)

$$\begin{aligned} m_{h,H}^2&= \frac{1}{2}\bigg [m_A^2 + (\eta _1 + \eta _5) v^2\nonumber \\&\quad \mp \sqrt{\left( m_A^2+ (\eta _5 - \eta _1) v^2\right) ^2 + 4 \eta _6^2 v^4}\bigg ], \end{aligned}$$

(4)

$$\begin{aligned} m_{H^\pm }^2&= \frac{1}{2}\eta _3 v^2+ \mu _{22}^2 . \end{aligned}$$

(5)

We now express the quartic couplings $\eta _1$, $\eta _{3-6}$ in terms of [32, 34] $\mu _{22}$, $m_h$, $m_H$, $m_A$, $m_{H^\pm }$, all normalized to v, and the mixing angle ${\gamma }$,

$$\begin{aligned} \eta _1&= \frac{m_h^2 s_\gamma ^2 + m_H^2 c_\gamma ^2}{v^2}, \end{aligned}$$

(6)

$$\begin{aligned} \eta _3&= \frac{2\left( m_{H^\pm }^2 - \mu _{22}^2\right) }{v^2}, \end{aligned}$$

(7)

$$\begin{aligned} \eta _4&= \frac{m_h^2 c_\gamma ^2 + m_H^2 s_\gamma ^2 -2 m_{H^\pm }^2+m_A^2}{v^2}, \end{aligned}$$

(8)

$$\begin{aligned} \eta _5&= \frac{m_H^2 s_\gamma ^2 + m_h^2 c_\gamma ^2 - m_A^2}{v^2}, \end{aligned}$$

(9)

$$\begin{aligned} \eta _6&= \frac{\left( m_h^2 - m_H^2\right) (-s_\gamma )c_\gamma }{v^2}. \end{aligned}$$

(10)

The quartic couplings $\eta _2$ and $\eta _7$ do not enter in the scalar masses, nor in the mixing angle $\gamma$. Therefore in our analysis we take v, $m_h$, and $\gamma$, $m_A$, $m_H$, $m_{H^\pm }$, $\mu _{22}$, $\eta _2$, $\eta _7$ as the nine phenomenological parameters.

The scalars h, H, A and $H^\pm$ couple to fermions by [33, 34]

$$\begin{aligned} {\mathcal {L}}&= -\frac{1}{\sqrt{2}} \sum _{F = U, D, L} {\bar{F}}_{i} \bigg [\left( -\lambda ^F_{ij} s_\gamma + \rho ^F_{ij} c_\gamma \right) h \nonumber \\&\quad +\left( \lambda ^F_{ij} c_\gamma + \rho ^F_{ij} s_\gamma \right) H -i ~\mathrm{sgn}(Q_F) \rho ^F_{ij} A\bigg ] P_R\; F_{j}\nonumber \\&\quad -\bar{U}_i\left[ (V\rho ^D)_{ij} P_R-(\rho ^{U\dagger }V)_{ij} P_L\right] D_j H^+ \nonumber \\&\quad - {\bar{\nu }}_i\rho ^L_{ij} P_R \; L_j H^+ +\mathrm{H.c.}, \end{aligned}$$

(11)

where $P_{L,R}\equiv (1\mp \gamma _5)/2$, $i,j = 1, 2, 3$ are generation indices, V is Cabibbo–Kobayashi–Maskawa matrix and, $U=(u,c,t)$, $D = (d,s,b)$, $L=(e,\mu ,\tau )$ and $\nu =(\nu _e,\nu _\mu ,\nu _\tau )$ are vectors in flavor space. The matrices $\lambda ^F_{ij}\; (=\sqrt{2}m_i^F/v)$ are real and diagonal, whereas $\rho ^F_{ij}$ are in general complex and non-diagonal. In the following we drop superscript F. For simplicity, we assume all $\rho _{ij}$ are real in our analysis. It is likely that $\rho _{ij}$ follow similar flavor organizing principle as in SM i.e. $\rho _{ii}\sim \lambda _i$ with suppressed off-diagonal elements of $\rho _{ij}$ matrices [32]. Therefore $\rho _{tt}\sim \lambda _t$, $\rho _{bb}\sim \lambda _b$ etc., while as we show below the flavor changing neutral Higgs coupling $\rho _{tc}$ could still be large. In the following, for simplicity we assumed $\lambda _t$, $\rho _{tt}$, and $\rho _{tc}$ to be nonzero and set all other $\lambda _i$ and $\rho _{ij}$ couplings to zero; their impact will be discussed in the later part of the paper.

For inflationary dynamics we chose the $m_H$, $m_A$, and $m_{H^\pm }$ between 200 and 800 GeV. This is primarily because of our aim to find signatures at the collider experiments, in particular at the LHC. In general lighter masses are possible. However they will be subjected to severe bounds from flavor physics as well as direct searches. We remark that heavier masses are also possible for inflationary dynamics. The potential for discovery or probing, although, becomes limited for heavier masses due to rapid fall in the parton luminosity. Thus we focus on sub-TeV mass range and restrict ourselves below 800 GeV.² As discussed earlier it is likely that $\rho _{ii}\sim \lambda _i$. However, as we shall see below for the bulk of the 200–800 GeV mass range $\rho _{tt}=\lambda _t$ is excluded by various direct and indirect searches. In particular we set $\rho _{tt} = 0.5$ at low scale. Furthermore we take $\rho _{tc}=0.2$, which is still allowed by current data and can have exquisite signatures at the LHC.

3 Inflationary dynamics

To study the inflationary dynamics we first write down the action in Jordan’s frame:

$$\begin{aligned} S&= \int d^4 x \sqrt{-g}\bigg [-\frac{M_\text{ P}^2}{2}\bigg (1+ 2 \xi _{11} |\varPhi |^2 + 2 \xi _{22} |\varPhi '|^2 \nonumber \\&\quad + 2 \left( \xi _{12} \varPhi ^\dagger \varPhi '+ h.c.\right) \bigg )R -g^{\mu \nu }\big (\partial _\mu \varPhi ^\dagger \partial _\nu \varPhi \nonumber \\&\quad + \partial _\mu \varPhi '^\dagger \partial _\nu \varPhi '\big )-V(\varPhi ,\varPhi ')\bigg ], \end{aligned}$$

(12)

where $\xi _{11}$, $\xi _{22}$, and $\xi _{12}$ are dimensionless nonminimal couplings; $g^{\mu \nu }$ and g are the inverse and determinant of metric, respectively; and $M_\text{ P }$ is the reduced Planck mass ($\approx 2.4\times 10^{18}$ GeV) with $M_\text{ P } =1$. The action in Eq. (12) can be written in the Einstein’s frame as

$$\begin{aligned} S_E&= \int d^4 x \sqrt{-g_E}\bigg [-\frac{R}{2}+\frac{3}{4}\left( \partial _\mu \left( \log F^2\right) \right) ^2\nonumber \\&\quad -\frac{|\partial _\nu \varPhi |^2+ |\partial _\mu \varPhi '|^2}{F^2}-V_E(\varPhi ,\varPhi ')\bigg ]. \end{aligned}$$

(13)

where,

$$\begin{aligned} F^2=1+ 2 \left( \xi _{11} |\varPhi |^2 + \xi _{22} |\varPhi '|^2 + \left( \xi _{12} \varPhi ^\dagger \varPhi '+ h.c.\right) \right) \end{aligned}$$

(14)

and $V_E(\varPhi ,\varPhi ') = V(\varPhi ,\varPhi ')/F^2$.

For inflationary dynamics we choose the Higgs field in the electromagnetic preserving direction:

$$\begin{aligned} \varPhi = \frac{1}{\sqrt{2}} \begin{pmatrix} 0 \\ \rho _1 \\ \end{pmatrix}\quad \text{ and },\quad \varPhi ' = \frac{1}{\sqrt{2}} \rho _2 \begin{pmatrix} 0\\ e^{i\chi }\\ \end{pmatrix}. \end{aligned}$$

(15)

The Einstein action in terms of field $\phi ^I=\{\rho _1,\rho _2,\chi \}$ becomes

$$\begin{aligned} S_E = \int d^4 x \sqrt{-g_E} \bigg [-\frac{R_E}{2}-S_{IJ}g_E^{\mu \nu }\partial _\mu {\phi _I^\dagger } \partial _\nu \phi _J-V_E(\phi ^I)\bigg ], \end{aligned}$$

(16)

where $S_{IJ}= \delta _{IJ}/F + 3k \; F^\dagger _I F_J/(2F^2)$ with $F_I=\partial F/\partial \phi _I$; $k= 1$ and 0 are for metric and Palatini formulations, respectively. The potential $V_E(\phi ^I)$ can be written as

$$\begin{aligned} V_E(\rho _1,\rho _2,\chi )&= \frac{1}{8 \left( 1+ \xi _{11} \rho _1^2 + \xi _{22} \rho _2^2 + 2 \xi _{12}{c_\chi }\rho _1 \rho _2\right) ^2}\nonumber \\&\quad \times \bigg [{\tilde{\eta _1}} \rho _1^4+{\tilde{\eta _2}} \rho _2^4+ 2\rho _1^2\rho _2^2\big ({\tilde{\eta _3}} +\big ({\tilde{\eta _4}} + {c_{2\chi }} \tilde{\eta _5} \big )\big )\nonumber \\&\quad +4{c_\chi }\rho _1 \rho _2 \left( {\tilde{\eta _6}} \rho _1^2+{\tilde{\eta _7}} \rho _2^2\right) \bigg ], \end{aligned}$$

(17)

where $c_\chi = \cos \chi$ and ${c_{2\chi }} = \cos 2\chi$ and we have only taken into account the quartic terms of the Jordan-frame potential V, discarding the quadratic terms, as we are interested in the inflaton dynamics for very large field values. The $\tilde{\eta _i}$s denote the quartic couplings in Eq. (1) at the inflationary scale.

As we will see below, one nonminimal coupling is sufficient to account for all the observational constraints on the Higgs inflation. Therefore in the following we turn only one nonminimal coupling at a time. In particular we primarily focus on the scenarios when either of $\xi _{11}$ and $\xi _{22}$ are nonzero, while $\xi _{12}=0$ throughout, and denote them as Scenario-I and Scenario-II respectively. The impact of nonzero $\xi _{12}$ will be briefly discussed at the latter part of the paper.

3.1 Scenario-I

In the Scenario-I we set $\xi _{22}=0$. Let us perform following field redefinition [26]:

$$\begin{aligned}&\rho = \frac{\rho _2}{\rho _1}~\text{ and }~\varphi = \sqrt{\frac{3}{2}}\log \left( 1+ \xi _{11} \rho _1^2\right) . \end{aligned}$$

(18)

With this field redefinition, the potential in Scenario-I becomes

$$\begin{aligned} V(\rho , \varphi ,\chi )&= \frac{1}{8\xi _{11}^2}\bigg [{\tilde{\eta _1}} +{\tilde{\eta _2}} \rho ^4+ 2 \rho ^2 \left( {\tilde{\eta _3}} + {\tilde{\eta _4}} + {c_{2\chi }} {\tilde{\eta _5}}\right) \nonumber \\&\quad +4{c_\chi }\rho \left( {\tilde{\eta _6}} +{\tilde{\eta _7}} \rho ^2\right) \bigg ] \left( 1- e^{-2\varphi /\sqrt{6}}\right) ^2, \end{aligned}$$

(19)

where $\varphi$ can play the role of inflaton. To find the slow-roll direction the $\varphi$ independent part of Eq. (19)

$$\begin{aligned} V(\rho ,\chi )&\approx \frac{1}{8\xi _{11}^2}\bigg [{\tilde{\eta _1}} +{\tilde{\eta _2}} \rho ^4+ 2 \rho ^2 \left( {\tilde{\eta _3}} + {\tilde{\eta _4 }}+ \left( 2 c_\chi ^2 -1\right) {\tilde{\eta _5}}\right) \nonumber \\&\quad +4{c_\chi }\rho \left( {\tilde{\eta _6}} +{\tilde{\eta _7}} \rho ^2\right) \bigg ] \end{aligned}$$

(20)

has to be minimized with respect to $\rho$ and $c_{\chi }$. It is hard to find analytical minimization for Eq. (20). Instead we minimize Eq. (20) numerically as follows. The potential has a extremum at $(\rho _0,c_{\chi _0})$, which is found by solving $\partial V/\partial \rho =0$ and $\partial V/\partial c_\chi =0$ simultaneously. The extremum is considered a minimum if both the determinant and trace of the covariant matrix $X_{ij}=\partial ^2 V/\partial x_i \partial x_j$ (with $x_{i,j}= \rho ~\text{ and }~ c_\chi$), calculated at the minima $(\rho _0,c_{\chi _0})$, are $>0$. In total there are three cases of minima $(\rho _0 ,c_{\chi _0})$ which we categorize as $c_\chi \ne \pm 1$, $c_\chi = 1$, and $c_\chi = -1$. In general the case $\rho _0 = 0$ and $c_{\chi _0} = 0$ could be a minimum, however the determinant of the covariant matrix $X_{ij}$ in this case is $\propto -\eta _6^2$. As we assume all $\eta _i$s real, the case $\rho _0 = 0$ and $c_{\chi _0} = 0$ cannot be a minimum in our case. The minima for the case $c_\chi = 1$ and $c_\chi = -1$ are found simply setting $c_\chi = \pm 1$ and demanding $\partial V/\partial \rho =0$ with $\partial ^2 V/\partial \rho ^2 > 0|_{\rho = \rho _0}$.

After stabilizing the potential at the minima $(\rho _0,c_{\chi _0})$, the potential for single Higgs inflation becomes

$$\begin{aligned} V\approx \frac{\eta _\mathrm{{eff}}}{8\xi _{11}^2}\left( 1- e^{-2\varphi /\sqrt{6}}\right) ^2, \end{aligned}$$

(21)

where

$$\begin{aligned} \eta _\mathrm{{eff}}&={\tilde{\eta _1}} +{\tilde{\eta _2 }}\rho _0^4+ 2 \rho _0^2 \left( {\tilde{\eta _3}} + {\tilde{\eta _4}} + \left( 2 c_{\chi _0}^2 -1\right) {\tilde{\eta _5}}\right) \nonumber \\&\quad +4 c_{\chi _0} \rho _0\left( {\tilde{\eta _6 }}+{\tilde{\eta _7}} \rho _0^2\right) \end{aligned}$$

(22)

is required to be positive to have a positive potential energy $V_0$ during inflation.

3.2 Scenario II

In Scenario II, the potential of Eq. (17) after field redefinition becomes

$$\begin{aligned} V(\rho , \varphi ,\chi )&= \frac{1}{8 \left( \rho ^2\xi _{22}\right) ^2}\bigg [{\tilde{\eta _1}} +{\tilde{\eta _2}} \rho ^4+ 2 \rho ^2 \left( {\tilde{\eta _3}} + {\tilde{\eta _4}} + {c_{2\chi }} {\tilde{\eta _5}}\right) \nonumber \\&\quad +4{c_\chi }\rho \left( {\tilde{\eta _6}} +{\tilde{\eta _7}} \rho ^2\right) \bigg ] \left( 1- e^{-2\varphi /\sqrt{6}}\right) ^2, \end{aligned}$$

(23)

where

$$\begin{aligned} \rho = \frac{\rho _2}{\rho _1}\quad \text{ and }\quad \varphi = \sqrt{\frac{3}{2}}\log \left( 1+ \xi _{22} \rho _2^2\right) . \end{aligned}$$

(24)

As in previous section, we minimize the $\varphi$-independent part of potential (23) numerically. Again, there exists three sets of minima: $c_\chi \ne \pm 1$, $c_\chi = 1$ and $c_\chi = -1$. After stabilizing the potential at the minima $(\rho _0,c_{\chi _0})$, the potential for single Higgs inflation becomes

$$\begin{aligned} V\approx \frac{\eta _\mathrm{{eff}}}{8\xi _{22}^2}\left( 1- e^{-2\varphi /\sqrt{6}}\right) ^2, \end{aligned}$$

(25)

where $\eta _\mathrm{{eff}}$ is written as

$$\begin{aligned} \eta _\mathrm{{eff}}&=\frac{1}{\rho _0^4}\bigg [{\tilde{\eta _1 }}+{\tilde{\eta _2}} \rho _0^4+ 2 \rho _0^2 \left( {\tilde{\eta _3}} +{\tilde{ \eta _4}} + \left( 2 c_{\chi _0}^2 -1\right) {\tilde{\eta _5}}\right) \nonumber \\&\quad +4c_{\chi _0} \rho _0\left( {\tilde{\eta _6}} +{\tilde{\eta _7}} \rho _0^2\right) \bigg ], \end{aligned}$$

(26)

calculated at the minimum $(\rho _0,c_{\chi _0})$. As before this is required to be positive for the positive potential energy during inflation.

3.3 Kinetic mixing

If there exist kinetic mixing, the heavy state needs to be integrated out during inflation to get an effective theory [37,38,39,40] such that $\varphi$-independent parts of Eq. (19) or Eq. (23) would induce the slow-roll inflation for the light state ($\simeq \varphi$) while the mass of the heavy state is exponentially suppressed. Let us elaborate on this.

The kinetic terms of the Lagrangian can be written as:

$$\begin{aligned} {\mathcal {L}}_\mathrm{{kin}}&\approx - \frac{1}{2}\left( 1+\frac{\rho ^2+1}{6 \left( \xi _{11}+\xi _{22} \rho ^2\right) }\right) (\partial _\mu \varphi )^2-\nonumber \\&\quad \times \frac{\left( \xi _{11}-\xi _{22}\right) \rho }{\sqrt{6}\left( \xi _{11}+\xi _{22} \rho ^2\right) {}^2} (\partial ^\mu \varphi )(\partial _\mu \rho )- \frac{\xi _{11}^2+\xi _{22}^2 \rho ^2}{2\left( \xi _{11}+\xi _{22} \rho ^2\right) {}^3} (\partial _\mu \rho )^2\nonumber \\&\quad -\frac{1}{2} \frac{\rho ^2}{\xi _{22}\rho ^2 + \xi _{11}}\left( 1- e^{-2\varphi /\sqrt{6}}\right) (\partial _\mu \chi )^2 \end{aligned}$$

(27)

It is clear from Eq. (27) that the kinetic terms are not canonically normalized, i.e., there exist kinetic mixing between $\varphi$ and $\rho$. To find canonically normalized kinetic terms, we closely follow the prescription laid out in Ref. [26]. For finite value of the Higgs ratio $\rho$, we consider a perturbation around the minimum $\rho _0$ as $\rho = \rho _0 + {{\bar{\rho }}}$. The kinetic terms of $\varphi$ and ${{\bar{\rho }}}$ can be rewritten as $-1/2~K_{ij} \partial ^\mu \phi _i \partial _\mu \phi _j$ ($\phi = {{\bar{\rho }}},\varphi$), where $K_{\varphi \varphi }=\alpha$, $K_{{{\bar{\rho }}} {{\bar{\rho }}}}=\beta$ and $K_{\varphi {{\bar{\rho }}}}={\gamma '}$. The potentials of Eqs. (21) and (25) can be expanded around the minima as $V\approx V_0 + A {{\bar{\rho }}}^2$ where, $A = A_\mathrm{{I}}/\xi _{11}^2$ in Scenario-I, whereas $A = A_\mathrm{{II}}/\xi _{22}^2$ in Scenario-II. The quantity $A_\mathrm{{I}}$ and $A_\mathrm{{II}}$ are:

$$\begin{aligned} A_\mathrm{{I}}&=\frac{1}{4} \left( 2 {\tilde{\eta _5}} c_{\chi _0}^2+6 {\tilde{\eta _7}} \rho _0c_{\chi _0} +{\tilde{\eta _3}}+{\tilde{\eta _4}}-{\tilde{\eta _5}}+3 {\tilde{\eta _2}} \rho _0^2\right) ,\end{aligned}$$

(28)

$$\begin{aligned} A_\mathrm{{II}}&= \frac{1}{\left( 4 \rho _0^6\right) } \bigg (5 {\tilde{\eta _1}} + 2 {\tilde{\eta _7}} \rho _0^3 c_{\chi _0}+3 \rho _0 \bigg (4 {\tilde{\eta _6}} c_{\chi _0}\nonumber \\&\quad +\rho _0 \left( 2 {\tilde{\eta _5}} c_{\chi _0}^2+{\tilde{\eta _3}}+{\tilde{\eta _4}}-{\tilde{\eta _5}}\right) \bigg )\bigg ). \end{aligned}$$

(29)

Both $A_\mathrm{{I}}$ and $A_\mathrm{{II}}$ are required to be positive. This is an additional requirement in addition to the conditions for potential minimization as described earlier. We can now diagonalize the kinetic terms via the following transformation:

$$\begin{aligned} \varphi '&= \cos \theta ~\varphi + \sin \theta \;{{\bar{\rho }}} \end{aligned}$$

(30)

$$\begin{aligned} {{\bar{\rho }}}'&= -\sin \theta ~\varphi + \cos \theta \;{{\bar{\rho }}}, \end{aligned}$$

(31)

where

$$\begin{aligned} \theta = \mathrm {ArcTan}\bigg [ \frac{2 {{\gamma '}}}{\alpha -\beta +\sqrt{\alpha ^2-2 \alpha \beta +\beta ^2+4 {{\gamma '}}^2}}\bigg ]. \end{aligned}$$

(32)

The eigenvalues of the kinetic terms can be identified as

$$\begin{aligned} \lambda _\pm =\frac{\alpha +\beta \pm \sqrt{\alpha ^2-2 \alpha \beta +\beta ^2+4 {{\gamma '}}^2}}{2}, \end{aligned}$$

(33)

while the potential can be re-expressed in terms of the new variables as

$$\begin{aligned} V\approx V_0 + A \left( \sin \theta \varphi ' + \cos \theta {{\bar{\rho }}}'\right) ^2. \end{aligned}$$

(34)

By further field redefinition ${{\tilde{\varphi }}} = \sqrt{\lambda _+} \varphi '$ and ${{\tilde{\rho }}} = \sqrt{\lambda _-} {{\bar{\rho }}}'$, the kinetic terms become canonically normalized. The different elements of the mass matrix for ${{\tilde{\varphi }}}$ and ${{\tilde{\rho }}}$ are

$$\begin{aligned} m^2_{{{\tilde{\varphi }}} {{\tilde{\varphi }}}}&= \frac{A \sin ^2\theta }{\lambda _+}, \end{aligned}$$

(35)

$$\begin{aligned} m^2_{{{\tilde{\varphi }}} {{\tilde{\rho }}}}&= m^2_{{{\tilde{\rho }}}{{\tilde{\varphi }}} } = \frac{A \sin \theta \cos \theta }{\sqrt{\lambda _+ \lambda _-}},\end{aligned}$$

(36)

$$\begin{aligned} m^2_{{{\tilde{\rho }}}{{\tilde{\rho }}}}&= \frac{A \cos ^2\theta }{\lambda _-}. \end{aligned}$$

(37)

After diagonalizing the mass matrix we get two eigenvalues $m^2_\mathrm{{light}} = 0$ and $m^2_\mathrm{{heavy}} = A \left( \sin ^2\theta /\lambda _+ + \cos ^2\theta \lambda _-\right)$. The $\varphi$-dependent part of Eq. (21) or Eq. (25) induces slow-roll inflation for the massless mode $m_\mathrm{{light}}$, while the $m_\mathrm{{heavy}}$ mode is exponentially suppressed. In Scenario-I (Scenario-II) for the large value of $\xi _{11}$ ($\xi _{22}$), the mass of the heavy state becomes $m^2_\mathrm{{heavy}}\sim A_\mathrm{{I}}/\xi _{11}$ ($\sim (A_\mathrm{{II}}\ \rho _0^4)/\xi _{22}$). This is much larger than the Hubble parameter ${\mathcal {H}}^2\sim \eta _\mathrm{{eff}}/\xi _{11}^2~(\sim \eta _\mathrm{{eff}}/\xi _{22}^2)$, and heavy states can be integrated out. To find the parameter space for inflation, along with all aforementioned conditions, for both the scenarios additionally we also demanded $m^2_\mathrm{{heavy}}>{\mathcal {H}}^2$ in our numerical analysis.

4 Parameter space for inflation

4.1 Inflationary observables

Let us spell out our notation for basic quantities. The dimensionless slow-roll parameters which measures the slope and curvature are defined as $\epsilon _\varphi = (1/2)\left( V_{,\varphi }/V\right) ^2$ and $\eta _\varphi =V_{,\varphi \varphi }/V$ where $V_{,\varphi }=\partial V/\partial \varphi$ and $V_{,\varphi \varphi }=\partial ^2 V/\partial \varphi ^2$. The quantities $n_s=1 + 2 \eta _\varphi -6 \epsilon _\varphi$ and $n_t=- 2 \epsilon _\varphi$ are the scalar and tensor spectral indices, respectively, while $A_s = \frac{V}{24 \pi ^2 \epsilon _\varphi }$ and $A_t = \frac{2 V}{3 \pi ^2}$ are the scalar and tensor amplitudes, respectively. To first order approximation $r_\varphi = A_t/A_s = 16 \epsilon _\varphi = - 8 \eta _\varphi$.

4.2 Observational constraints on inflation

For consistent inflationary model the observational constraints from Planck 2018 results are [10]

$$\begin{aligned} A_s^*&= (2.099\pm 0.014)\times 10^{-9} \quad {68\% }~\mathrm {CL},\end{aligned}$$

(38)

$$\begin{aligned} n_s^*&= 0.9649\pm 0.0042 \quad {68\% }~\mathrm {CL},\end{aligned}$$

(39)

$$\begin{aligned} r_{\varphi ^*}&< 0.056 \quad {95\% }~\mathrm {CL}, \end{aligned}$$

(40)

where $A_s^*$, $n_s^*$, and $r_{\varphi ^*}$ are the scalar amplitude, the scalar spectral index, and the tensor-to-scalar ratio, respectively, evaluated at $\varphi =\varphi ^*$. The value of $\varphi ^*$ is obtained by solving the number of e-foldings N

$$\begin{aligned} N \approx \int _{\varphi _\mathrm{{end}}}^{\varphi ^*} d\varphi \frac{V}{V_{,\varphi }}, \end{aligned}$$

(41)

where $\varphi ^*$ correspond to the value of inflaton field when number of e-foldings $N=60$, and $\varphi _\mathrm{{end}}$ denotes the end of slow-roll approximation defined as $\epsilon _\varphi (\varphi _\mathrm{{end}}) :=1$. If we approximate that $\varphi _\text{ end }=0$, from Eq. (41) one finds

$$\begin{aligned} N = \frac{3}{4}\left[ e^{\sqrt{\frac{2}{3}}\varphi ^*} - \sqrt{\frac{2}{3}}\varphi ^* -1 \right] , \end{aligned}$$

(42)

while $A_s^*$, $n_s^*$ and $r_{\varphi ^*}$ are

$$\begin{aligned} A_s^*&= \frac{\eta _\mathrm{{eff}} \sinh ^4\left( \frac{ \varphi ^*}{\sqrt{6}}\right) }{16 \pi ^2 \xi ^2},\end{aligned}$$

(43)

$$\begin{aligned} n_s^*&= \frac{1}{3} \left[ 4 \coth \left( \frac{\varphi ^*}{\sqrt{6}}\right) -4\text {csch}^2\left( \frac{\varphi ^*}{\sqrt{6}}\right) -1\right] , \end{aligned}$$

(44)

$$\begin{aligned} r_{\varphi ^*}&= \frac{64}{3 \left( e^{\sqrt{\frac{2}{3}} \varphi ^*}-1\right) ^2}, \end{aligned}$$

(45)

with $\xi = \xi _{11}$ or $\xi _{22}$.

Solving Eq. (42) for $N=60$ we find $\varphi ^* \approx 5.45$. Correspondingly, $n_s^*\approx 0.9675$ and $r_{\varphi ^*}\approx 3.03\times 10^{-3}$, which are within the limits obtained by Planck 2018 [10], as can be seen from Eqs. (39) and (40). Moreover, scalar amplitude $A_s^*$ of Eq. (43) needs to satisfy the constraint as in Eq. (38).

4.3 Scanning and parameter space

At the low scale ($\mu = m_W$), the dynamical parameters in Eq. (1) need to satisfy the unitarity, perturbativity, positivity constraints, for which we utilized 2HDMC [41]. To save computation time we generated the input parameters $\gamma$, $m_A$, $m_H$, $m_{H^\pm }$, $\mu _{22}$, $\eta _2$, $\eta _7$ randomly in the ranges: $c_\gamma =[0,0.05]$, $m_H=[200,800]$ GeV, $m_A=[200,800]$ GeV, $m_{H^\pm }=[200,800]$ GeV, $\mu _{22}=[0,1000]$ GeV, $\eta _2=[0,1]$ and $\eta _7=[-1,1]$, with $v = 246$ GeV, and $m_h= 125$ GeV. We call them parameter points and fed into 2HDMC for scanning in the Higgs basis. The input parameters in 2HDMC [41] are $\varLambda _{1-7}$ and $m_{H^\pm }$ in the Higgs basis with v being an implicit parameter, and we identify $\varLambda _{1-7}$ with $\eta _{1-7}$. To match the convention of 2HDMC, we take $-\pi /2\le \gamma \le \pi /2$. For more details on the convention and parameter counting we redirect readers to Refs. [33, 42]. One has to also consider oblique T parameter [43] constraint, which restricts the hierarchical structures of the scalar masses [44, 45], and therefore $\eta _i$s. We utilize the expression given in Ref. [45]. The parameter points that passed unitarity, perturbativity and positivity conditions from 2HDMC are further needed to satisfy the T parameter constraint within the $2\sigma$ error [46, 47].

We shall see shortly $\eta _\mathrm{{eff}}\gtrsim 1$ is favored by inflationary constraints, which implies that $\xi _{11}$ and $\xi _{22}$ should be ${\mathcal {O}}(10^4)$ to generate the observed spectrum of CMB density perturbations. On the other hand, unitarity is broken at momentum scales $\mu \gtrsim M_\text{ P }/\xi _{11}$ ($M_\text{ P }/\xi _{22}$) for a scattering around the electroweak vacuum in Scenario-I (Scenario-II), and one might expect that higher dimensional operators are suppressed only by $M_\text{ P }/\xi _{11}$ ($M_\text{ P }/\xi _{22}$) rather than by $M_\text{ P }$. We may either assume that the coefficients of higher dimensional operators have extra suppressions or introduce additional scalars at the inflationary scale to restore the unitarity as discussed in Refs. [24,25,26]. The RGE above unitarity scale depends on the ultra-violet (UV) completion of the model [48], and we only perform the RGE computation of the parameters of g2HDM up to the unitarity scale. As for the unitarity scale, we take $y \approx 26$, corresponding to the scale $\mu =1.6\times 10^{13}$ GeV with $y=\ln (\mu / m_W)$ such that unitarity is maintained for the ballpark nonminimal couplings ${\mathcal {O}}(10^4)$ and $\eta _\mathrm{{eff}}\sim 1$. Later we also call this scale the high scale. The high scale parameters are denoted with tilde in order to differentiate them from the corresponding low scale parameters in Eqs. (1) and (11).

For the RGE of the parameters in Eq. (1) as well as $\rho ^F$ and $\lambda ^F$ in the Eq. (11), from low scale $y=0$ to high scale $y= 26$, we utilized the $\beta _x$ functions ($\beta _x= \partial x /\partial y$) for g2HDM given in Ref. (The beta functions for $\eta _i$s, $\lambda _i$ and $\rho _{ij}$ are taken from [49]). The parameter points that survive the low scale constraints from unitarity, perturbativity, positivity, and T parameter are entered in the RG equations. At the high scale we check perturbativity (i.e. couplings are being within $[-\sqrt{4\pi },\sqrt{4\pi }]$) for ${\tilde{\lambda _i}}$s, ${\tilde{\rho _{ij}}}$s, and $|{\tilde{\eta _i}}|$ as well as positivity ${{\tilde{\eta }}_{1,2}} > 0$. We found that parameter points with $|\eta _i| > 1$ at the low energy get generally excluded after imposing perturbativity and positivity criteria at the high scale. Therefore, with limited computational facility to save time while generating parameters at low scale, we more conservatively demanded $|\eta _{i}|\le 1$.

Fig. 1

The scanned points plotted for Scenario-I (Scenario-II) in $\eta _\mathrm{{eff}}$–$\xi _{11}$ ($\eta _\mathrm{{eff}}$–$\xi _{22}$) plane in the left panel (right panel) respectively. The purple and orange scanned points are respectively for $c_{\chi _0}=1$ and $c_{\chi _0}=-1$. See text for further details

Fig. 2

The corresponding values of the minima of the scanned points in Fig. 1 in $\rho _0$ vs $c_{\chi _0}$ plane for the Scenario-I (left panel) and Scenario-II (right panel)

Fig. 3

The scanned points corresponding to Fig. 1 that are traced back and plotted in the $m_A$–$m_H$ and $m_{H^\pm }$–$m_H$ planes for $y=0$ in Scenario-I

Fig. 4

Same figure as Fig. 3 but for Scenario-II

At the high scale, for each parameter point, we require to satisfy all the necessary conditions as described in the previous section, such as $\eta _\mathrm{{eff}}$ and the potential energy $V_0$ being positive at the potential minima. Finally, the points are needed to satisfy inflationary constraints of Eqs. (38)–(40) from Planck 2018 [10]. The parameter points that passed all of the above mentioned required conditions as well as Eqs. (38)–(40) are termed as “scanned points”. The scanned points are plotted in Fig. 1 for Scenarios-I and -II in the $\xi _{11}$ vs $\eta _\mathrm{{eff}}$ (left panel) and $\xi _{22}$ vs $\eta _\mathrm{{eff}}$ (right panel) planes, respectively. Their corresponding minima $(\rho _0,c_{\chi _0})$ are plotted in Fig. 2. In both figures purple scanned points correspond to the minima $c_{\chi _0}=1$, while orange points are for $c_{\chi _0}=-1$. The scanned points in Fig. 1 are traced back to low scale (i.e. to $y=0$), and plotted in the $m_A$–$m_H$ and $m_{H^\pm }$–$m_H$ planes in Figs. 3 and  4.

In Scenario-I, we find that $0.5\lesssim \eta _\mathrm{{eff}}\lesssim 3.5$ with $10^3\lesssim \xi _{11}\lesssim 5\times 10^4$ as can be seen from the left panel of Fig. 1. The scanned points are mostly concentrated around $\xi _{11}\sim 10^4$. This can be understood easily from Eq. (38). For $\varphi ^*\approx 5.45$ and $0.5\lesssim \eta _\mathrm{{eff}}\lesssim 3.5$ one finds $\xi _{11}$ to be ${\mathcal {O}}(10^4)$. A similar pattern is found also for the Scenario-II. The corresponding minima $(\rho _0,c_{\chi _0})$ in the Scenario-I (Scenario-II) found to be in the range $0\lesssim |\rho _0|\lesssim 1$ ($1\lesssim |\rho _0|\lesssim 100$) while $c_{\chi _0}$ is either 1 or $-1$, as can be seen in Fig. 2. Note that, there exist no minima for $-1< c_{\chi _0} < 1$ in both the scenarios. In most cases $\rho _0$ is found to be complex for $-1< c_{\chi _0} < 1$. Indeed, there exist some real $\rho$ and $c_{\chi }$ that solve $\partial V/\partial \rho =0$ and $\partial V/\partial c_\chi =0$ simultaneously, however, the determinant and/or the trace of the covariant matrix $X_{ij}$ are found to be not positive in such cases.

Let us take a closer look at Figs. 3 and 4. We find that at the low scale the parameter space for inflation requires $m_H$, $m_A$ and $m_{H^\pm }$ to be nearly degenerate with $m_h= 125$ GeV. This behavior can be traced back to our choices of parameters at the low scale. As for the inflation, one requires perturbativity and positivity at the high scale. At the low scale while scanning we demanded all $|\eta _i| \le 1$. This is driven by the fact that the parameter points with $|\eta _i| > 1$ at the low scale tend to become nonperturbative at the high scale. Such choices severely restrict mass splittings between $m_H$, $m_A$ and $m_{H^\pm }$ which are primarily determined by the magnitudes of $\eta _i$. This can be understood easily from Eqs. (3)–(5) and Eqs. (6)–(10). With common $\mu _{22}^2$ terms in $m_H$, $m_A$ and $m_{H^\pm }$, the mass splittings are restricted because we require all $|\eta _i| < 1$. Thus we conclude that parameter space favored by inflation requires $m_H$, $m_A$ and $m_{H^\pm }$ to be nearly degenerate, as reflected in Figs. 3 and  4. In what follows we shall show that these mass ranges receive meaningful direct and indirect constraints and may have exquisite signatures at the LHC, ILC, FCC-ee, etc. In addition, we show that such near degeneracy can be directly probed at the LHC.

5 Direct and indirect searches

5.1 Constraints

Having already found the parameter space for inflation we now turn our attention to the constraints on the parameter space. The couplings $\rho _{tt}$ and $\rho _{tc}$ receive several direct and indirect search limits, particularly in the sub-TeV mass ranges of $m_A$, $m_H$ and $m_{H^\pm }$. We now summarize these constraints in detail. In particular we will show that the parameter space chosen for scanning in Sect. 4 is allowed by current data.

First we focus on the h boson coupling measurements by ATLAS [51] and CMS [52]. The results are provided as the ratios of the observed and SM productions and decay rates of h, called signal strengths. For nonvanishing $c_\gamma$ the extra Yukawa $\rho _{ij}$s modify the couplings h to fermions (see Eq. (11)). Therefore, the parameter space for inflation would receive meaningful constraint from such measurements. The ATLAS results [51] are based on Run-2 ($\sqrt{s}= 13$ TeV) 80 $\hbox {fb}^{-1}$ data, while CMS [52] utilized only up to 2016 Run-2 data (35.9 $\hbox {fb}^{-1}$ ). Both the collaborations measured signal strengths $\mu _i^f$ and corresponding errors to different production and decay chains $i\rightarrow h \rightarrow f$. The signal strengths $\mu _i^f$ are defined as [51, 52]:

$$\begin{aligned} \mu _i^f = \frac{\sigma _i{\mathcal {B}}^f}{(\sigma _i)_{\text {SM}}({\mathcal {B}}^f)_{\text {SM}}} = \mu _i \mu ^f, \end{aligned}$$

(46)

where $\sigma _i$ and ${\mathcal {B}}^f$ are the production cross sections of $i\rightarrow h$ and the branching ratio for $h\rightarrow f$ respectively. ATLAS and CMS considered gluon-fusion ($gg\text{ F }$), vector-boson-fusion ($\text{ VBF }$), Zh, Wh, $t{{\bar{t}}}h$ production processes (denoted by index i) and the $\gamma \gamma$, ZZ, WW, $\tau \tau$, bb, and $\mu \mu$ decay modes (by f). For simplicity we utilized the leading order (LO) $\mu _i^f$ and followed Refs. [53,54,55,56] for their explicit expressions. In our analysis, we focus particularly on the $gg\text{ F }$ and the $\text{ VBF }$ production modes because they put the most stringent constraints. In the $gg\text{ F }$ category we find that the most relevant signal strengths are $\mu _{gg\text{ F }}^{WW}$, $\mu _{gg\text{ F }}^{\gamma \gamma }$, $\mu _{gg\text{ F }}^{ZZ}$ and $\mu _{gg\text{ F }}^{\tau \tau }$, whereas in the $\text{ VBF }$ category $\mu _{\text{ VBF }}^{\gamma \gamma }$, $\mu _{\text{ VBF }}^{WW}$ and $\mu _{\text{ VBF }}^{\tau \tau }$. Further, we have also considered the Run-2 flagship observations of top Yukawa (htt) [57, 58] and bottom Yukawa (hbb) [59, 60] by ATLAS and CMS. We call all these measurements together “Higgs signal strength measurements”. Under the assumptions on couplings in Sect. 4, the flavor conserving couplings $\rho _{tt}$ would receive meaningful constraints for $c_\gamma \ne 0$. Allowing $2\sigma$ errors on each signal strength measurements we find that the $|\rho _{tt}|=0.5$ is still allowed by Higgs signal strength measurements for $c_\gamma = 0.05$. While finding the upper limit, we assumed $m_{H^\pm } = 200$ GeV, which enters in the $h\gamma \gamma$ couplings only from one loop level: The constraints have very mild dependence on $m_{H^\pm }$ and the results remain practically the same for the entire $m_{H^\pm } \in [200,800]$ GeV range.

For nonzero $\rho _{tt}$ the charged Higgs and W bosons loop with t quark modifies the $B_q$-${\overline{B}}_q$ ($q=d,s$) mixing amplitudes $M^q_{12}$. The constraint is stringent specially for the sub-TeV $m_{H^\pm }$. Recasting the type-II 2HDM expression of $B_q$-${\overline{B}}_q$ mixing amplitude [61], Ref. [62] found that $M^q_{12}$ can be written as

$$\begin{aligned} \frac{M^q_{12}}{M^{q\;\mathrm {SM}}_{12}} = 1+ \frac{I_{WH}(y_W, y_H, x) + I_{HH}(y_H)}{I_{WW}(y_W)}, \end{aligned}$$

(47)

where $x = m_{H^\pm }^2/m_W^2$, $y_i = m_t^2/m_i^2$ ($i = W, H^\pm$), and $m_t$ and $m_W$ are top quark and W boson masses. The respective expressions for $I_{WW}$, $I_{WH}$, and $I_{HH}$ are given as [62]

$$\begin{aligned} I_{WW}&=1+\frac{9(1-y_W)-6}{(1-y_W)^2}-\frac{6 \ln y_W}{y_W}\left( \frac{y_W}{1-y_W}\right) ^3,\end{aligned}$$

(48)

$$\begin{aligned} I_{WH}&\simeq \left( \frac{y_H|\rho _{tt}|^2}{\lambda _t^2}\right) \bigg [ \frac{(2x-8) \ln y_H}{(1-x) (1-y_H)^2} \nonumber \\&\quad +\frac{6x \ln y_W}{(1-x)(1-y_W)^2}- \frac{8- 2 y_W}{(1-y_W)(1-y_H)}\bigg ],\end{aligned}$$

(49)

$$\begin{aligned} I_{HH}&\simeq \frac{|\rho _{tt}|^4}{\lambda _t^4} \bigg [ \frac{ 1+ y_H}{(1- y_H)^2} + \frac{2 y_H \ln y_H}{(1- y_H)^3}\bigg ]y_H . \end{aligned}$$

(50)

For the quantity $C_{B_q} e^{2i \phi _{B_q}}:=M^q_{12}/M^{q\;\mathrm {SM}}$, one simply has $C_{B_q}=M^q_{12}/M^{q\;\mathrm {SM}}$ for real $\rho _{ij}$ couplings. The summer 2018 results of UTfit [63] found $C_{B_d}\in 1.05\pm 0.11$, $\phi _{B_d}\in -2.0\pm 1.8~~[\text{ in }~^{\circ } ]$, $C_{B_s}\in 1.110\pm 0.090$, and $\phi _{B_s}\in 0.42\pm 0.89~~[\text{ in }~^{\circ } ]$. Allowing $2\sigma$ uncertainties on the $C_{B_d}$ and $C_{B_s}$, the parameter space excluded by $B_{d,s}$-${\overline{B}}_{d,s}$ mixings are shown by the purple shaded region in $|\rho _{tt}|$–$m_{H^\pm }$ plane in Fig. 5. Note that here we have overlaid the excluded regions by $B_{d}$ and $B_s$ mixing in purple color and denote them together as $B_q$ mixings in Fig. 5.

Fig. 5

The purple and green shaded regions are excluded by $B_{q}$ mixings and $pp\rightarrow {\bar{t}} H^+$ search [66] respectively

Nonvanishing $\rho _{tt}$ can induce $V_{tb}$ enhanced $bg\rightarrow {\bar{t}} H^+$ and $gg\rightarrow {\bar{t}} b H^+$ processes (charge conjugate processes are implied). The processes $pp\rightarrow {\bar{t}} (b) H^+$ followed by $H^+\rightarrow t {\bar{b}}$ are the conventional search program for the $H^\pm$ and covered extensively by ATLAS [64] and CMS [65, 66]. The ATLAS search [64] provides model independent 95% CL upper limit on cross section times branching ratio ($\sigma (pp\rightarrow {\bar{t}} b H^+)\times {\mathcal {B}}(H^+ \rightarrow t {\bar{b}})$) based on its $\sqrt{s}=13$ TeV 36 $\hbox {fb}^{-1}$ dataset for $m_{H^\pm }=200$ GeV–2 TeV. Likewise CMS also set 95% CL upper limits on $\sigma (pp\rightarrow {\bar{t}} H^+)\times {\mathcal {B}}(H^+ \rightarrow t {\bar{b}})$, based on $\sqrt{s}=13$ TeV 35.9 $\hbox {fb}^{-1}$ dataset for $m_{H^\pm }=200$ GeV and 3 TeV in the semileptonic t decay [65], and on combination of semileptonic and all-hadronic final states [66]. We first extract these $\sigma \times {\mathcal {B}}$ upper limits (To obtain the 95% CL $\sigma \times {\mathcal {B}}$ upper limit we digitized the figures of Refs. [61,62,63]. A similar digitization strategy was followed in Ref. [67]) from Refs. [64,65,66] in the mass range $m_{H^\pm }=200$–800 GeV. In order to estimate the constraints, we determine the cross sections $\sigma (pp\rightarrow {\bar{t}} (b) H^+)\times {\mathcal {B}}(H^+ \rightarrow t {\bar{b}})$ at LO for reference $|\rho _{tt}|=1$ value for the $m_{H^\pm }=200$–800 GeV via Monte Carlo event generator MadGraph5_aMC@NLO [68] with NN23LO1 PDF set [69]. To obtain the respective 95% CL upper limits on $|\rho _{tt}|$ these cross sections are then rescaled by $|\rho _{tt}|^2$ simply assuming ${\mathcal {B}}(H^+ \rightarrow t {\bar{b}}) \approx 100\%$. We find the upper limits from ATLAS search [64] are in general much weaker than that of CMS searches [65, 66]. The upper limits from the CMS semileptonic final state [65] are mildly weaker compared to those from the combined semileptonic and all-hadronic final states [66]. Hence in Fig. 5 we only provide the regions excluded by the CMS search of Ref. [66], which is shown in green shaded region. While finding the excluded regions we assumed $\rho _{ij}=0$ except for $\rho _{tt}$ for the sake of simplicity. In general nonzero $\rho _{ij}$ couplings would turn on other decay modes of $H^+$ leading to even weaker upper limits on $\rho _{tt}$. E.g., the $\rho _{tc}$ coupling induces $V_{tb}$ proportional $H^+\rightarrow c {\bar{b}}$ decay. For $\rho _{tc}=0.2$ such additional decay mode can suppress the ${\mathcal {B}}(H^+ \rightarrow t {\bar{b}})$ by $20-30\%$ for $m_{H^\pm }=200$–800 GeV. While finding these upper limits, we have implemented the effective model in FeynRules [70].

Fig. 6

Regions excluded from $gg\rightarrow H/A \rightarrow t {\bar{t}}$ searches by ATLAS [71] and CMS [72] are shown in red and blue shaded regions respectively. The cyan shaded regions are excluded by CMS [73] search for heavy H/A production in association with t quarks with $H/A \rightarrow t {\bar{t}}$ decays

The search for heavy Higgs via $gg\rightarrow H/A \rightarrow t {\bar{t}}$ by ATLAS [71] and CMS [72] would be relevant to constrain $\rho _{tt}$ for $m_A/m_H > 2 m_t$. The ATLAS [71] search set exclusion limits on $\tan \beta$ vs $m_A~(\text{ or }~m_H)$ in type-II 2HDM framework starting from $m_A$ and $m_H=500$ GeV for two different mass hierarchies: $m_A = m_H$ and mass-decoupled $m_A$ and $m_H$. The search is based on $\sqrt{s}=8$ TeV (Run-1) 20.3 $\hbox {fb}^{-1}$ data. The CMS has performed similar search [72] but with Run-2 35.9 $\hbox {fb}^{-1}$ data, and provided 95% CL upper limit on coupling modifier (see Ref. [72] for definition) for $m_A$ ($m_H$) from 400–750 GeV based on different values of decay width to mass ratios $\varGamma _A/m_A$ ($\varGamma _H/m_H$) assuming $m_H$ ($m_A$) is decoupled. After reinterpreting ATLAS results for $m_A=m_H$, which are provided only for three benchmark points 500, 550, and 600 GeV [71], we find red shaded exclusion region in Fig. 6. Note that we utilized ATLAS $m_A=m_H$ result (and not the mass-decoupled $m_A$ and $m_H$ scenario) primarily because most scanned points in Figs. 3 and  4 resemble roughly $m_A\approx m_H$ pattern. We remark the actual constraints would be mildly weaker, depending on the value of $|m_A- m_H|$ for the respective scanned points. The limits for the mass-decoupled scenario are much weaker and not shown in Fig. 6. The CMS [72] provides limits only for mass-decoupled scenario which is shown in blue shaded regions in Fig. 6. The limits are weaker than those from ATLAS even though the latter used only Run-1 data. It is reasonable to assume the constraints could be stronger if CMS [72] provided results for $m_A=m_H$ scenario. A CMS analysis with mass degeneracy with full Run-2 dataset is welcome.

Moreover, $\rho _{tt}$ would also receive constraint from CMS search for SM four-top production [73] with 13 TeV 137 $\hbox {fb}^{-1}$ dataset. Apart from measuring SM four-top production, the search also set 95% CL upper limits on $\sigma (pp \rightarrow t {\bar{t}} A/t {\bar{t}} H)\times {\mathcal {B}}(A/H \rightarrow t {\bar{t}})$: $350~\text {GeV}\le m_{A/H} \le 650~\text {GeV}$. The search also included subdominant contributions from $\sigma (pp\rightarrow t W A/H, t q A/H)$ with $A/H\rightarrow t {\bar{t}}$. To find the constraint on $\rho _{tt}$ we generate these cross processes at LO by MadGraph5_aMC@NLO for a reference value of $|\rho _{tt}|$ setting all other $\rho _{ij} = 0$, and finally rescale simply by $|\rho _{tt}|^2$ assuming ${\mathcal {B}}(A/H \rightarrow t {\bar{t}})\approx 100\%$. We find that the constraints from $\sigma (pp \rightarrow t {\bar{t}} A)\times {\mathcal {B}}(A \rightarrow t {\bar{t}})$ are mildly stronger than that of $\sigma (pp \rightarrow t {\bar{t}} H)\times {\mathcal {B}}(H \rightarrow t {\bar{t}})$. The regions excluded by the former process is shown in cyan shaded regions in Fig. 6. We stress that for simplicity we assumed ${\mathcal {B}}(A/H \rightarrow t {\bar{t}})\approx 100\%$. As we chose $\rho _{tc}=0.2$, which will induce $A/H \rightarrow t {\bar{c}}$ decays, ${\mathcal {B}}(A/H \rightarrow t {\bar{t}})$ would be suppressed, hence the limits will be weaker than the shaded regions in Fig. 6. Note that as in before, while setting upper limits, CMS [73] assumed A (or H) is decoupled from H (or A), which is not the case for the scanned points in Figs. 3 and 4. We remark that the actual limit could possibly be stronger.

The coupling $\rho _{tc}$ receives constraints from ${\mathcal {B}}(t\rightarrow c h)$ measurement. For nonzero $c_\gamma$, $\rho _{tc}$ can induce flavor changing neutral current (FCNC) coupling htc (see Eq. (11)) which can induce $t\rightarrow ch$ decay. Both ATLAS and CMS have searched for the $t\rightarrow ch$ decay and provided 95% CL upper limits on ${\mathcal {B}}(t\rightarrow c h)$. The ATLAS upper limit is ${\mathcal {B}}(t\rightarrow c h) < 1.1\times 10^{-3}$ [74], while the CMS one is weaker ${\mathcal {B}}(t\rightarrow c h) < 4.7 \times 10^{-3}$ [75]. Both ATLAS and CMS results are based on 13 TeV $\sim 36$ $\hbox {fb}^{-1}$ dataset. For $c_\gamma =0.05$, which is the largest value considered while scanning, $|\rho _{tc}| \lesssim 1.8$ is excluded at $95\%$ CL. The constraint is weaker for smaller $c_\gamma$.

The constraints on $\rho _{tc}$ from ${\mathcal {B}}(t\rightarrow ch)$ measurement is rather weak. However, it has been found [76, 77] that $\rho _{tc}$ receives stringent constraint from the CMS search for SM four-top production [73] (based on 13 TeV 137 $\hbox {fb}^{-1}$ dataset), even when $c_\gamma$ is small. The search provides observed and expected number of events for different signal regions depending on the number of charged leptons and b-tagged jets with at least two same-sign leptons as baseline selection criteria [73]. It has been shown [77] that the CRW [73], i.e. the Control Region for $t{\bar{t}}W$ background, defined to contain two same-sign leptons and two to five jets with two of them b-tagged (see Ref. [73] for details), is the most relevant one to constrain $\rho _{tc}$. The Ref. [73] reported 338 events observed in CRW whereas the total events expected (denoted as SM expected events) is $335 \pm 18$ [73]. Induced by $\rho _{tc}$ coupling, the processes $cg\rightarrow t H/tA \rightarrow t t {\bar{c}}$ (charge conjugate processes always implied) with the semileptonically decaying same-sign top quarks have similar event topologies and contribute abundantly to the CRW. However, there is a subtlety. If the masses and widths of A and H are degenerate the $cg\rightarrow t H \rightarrow t t {\bar{c}}$ and $cg\rightarrow tA \rightarrow t t {\bar{c}}$ contributions interfere destructively, leading to exact cancellation between the amplitudes [76, 77]. The cancellation weakens if the mass splitting $|m_{H}-m_{A}|$ is large or widths of H and A become nondegenerate. For the scanned points in Figs. 3 and 4 we find that $|m_{H}-m_{A}|$ are small and widths are nearly degenerate. To understand how strong the constraint is we choose two representative $m_H$ and $m_A$ values: $m_H=200$ GeV and $m_A = 220$ GeV (i.e. $|m_{H}-m_{A}|\approx 20$ GeV) and $m_H=200$ GeV and $m_A = 250$ GeV (i.e. $|m_{H}-m_{A}|\approx 50$ GeV).

We first estimate $cg\rightarrow t H/tA \rightarrow t t {\bar{c}}$ contributions for a reference $\rho _{tc} = 1$ assuming ${\mathcal {B}}(H/A\rightarrow t {\bar{c}})=100\%$. Following the same event selection criteria described for CRW analysis [73], we rescale these contributions by $|\rho _{tc}|^2$ and demand that the sum of the events form the $cg\rightarrow t H/tA \rightarrow t t {\bar{c}}$ contributions and the SM expected events in CRW to agree with the number of the observed events within $2\sigma$ error bars for the SM expectation. We find that $\rho _{tc}\gtrsim 0.5$ is excluded at $2\sigma$ for the scenario $|m_{H}-m_{A}|\approx 20$ GeV, whereas $\rho _{tc}\gtrsim 0.4$ for $|m_{H}-m_{A}|\approx 50$ GeV. Due to smaller mass splitting, and therefore larger cancellation between the amplitudes, the constraint is weaker for the $|m_{H}-m_{A}|\approx 20$ GeV case compared to $|m_{H}-m_{A}|\approx 50$ GeV case. Here we simply assumed Gaussian (for more precise estimation of exclusion limits using likelihood function with Poisson counting, see e.g., [78]) behavior for the uncertainty of the SM expected events. Note that nonzero $\rho _{tc}$ will also induce $cc\rightarrow t t$ via t-channel H/A exchange, which we also included in our analysis. The events are generated at LO by MadGraph5_aMC@NLO interfaced with PYTHIA 6.4 [79] for showering and hadronization, and then fed to Delphes 3.4.2 [80] for fast detector simulation with CMS based detector card. For matrix element and parton shower merging we adopted MLM scheme [81].

For the heavier $m_H$ and $m_A$, we find that the constraints on $\rho _{tc}$ from CRW becomes weaker. This is simply because $cg\rightarrow t H/tA \rightarrow t t {\bar{c}}$ cross sections drops rapidly due to fall in the parton lumniosity. In finding the constraint we assumed ${\mathcal {B}}(H/A\rightarrow t {\bar{c}})\approx 100\%$. However, this assumption is too strong given $c_\gamma =0.05$. For nonzero $c_\gamma$ one has $A \rightarrow Zh$ decay for $m_A > m_Z +m_h$ (or $H\rightarrow h h$ decay for $m_H > 2m_h$), which will weaken the constraint further. In addition as we assumed $\rho _{tt}=0.5$. For scanned points in Figs. 3 and 4 where $m_H/m_A>2m_t$ the ${\mathcal {B}}(H/A\rightarrow t {\bar{c}})$ will be diluted further by large ${\mathcal {B}}(H/A\rightarrow t {\bar{t}})$.

In this regard we also note that ATLAS has also performed similar search [82] however we find the limits are weaker due to difference in event topologies and selection cuts. In addition, ATLAS has performed search [83] for R parity violating supersymmetry with similar event topologies. The selection cuts, however, are still too strong to give meaningful constraints on $\rho _{tc}$. Furthermore, $B_{s,d}$ mixing and ${\mathcal {B}}(B\rightarrow X_s\gamma )$, where $\rho _{tc}$ enters via charm loop through $H^\pm$ coupling [84], can still constrain $\rho _{tc}$. A reinterpretation of the result from Ref. [84] finds $|\rho _{tc}|\gtrsim 1$ is excluded from $B_s$ mixing, for the ballpark mass range of $m_{H^\pm }$ considered in our analysis. The constraints are weaker than those from the CRW region.

Before closing we remark that $\rho _{tt} \sim 0.5$ and $\rho _{tc} \sim 0.2$ are still allowed by the current direct and indirect searches for all scanned points in Figs. 3 and 4. So far for simplicity we set all $\rho _{ij}$ to zero in the previous section, however, there exist searches that can also constrain the parameter space if some of them are nonzero. E.g., the most stringent constraint on $\rho _{bb}$ arises from [96] CMS search for heavy H/A production in association with at least one b-jet and decaying into $b{\bar{b}}$ pair for $m_H/m_A$ 300–1300 GeV [85]. Following the same procedure as in before and utilizing $\sigma (pp\rightarrow b A/H +X)\cdot {\mathcal {B}}(A/H\rightarrow b {\bar{b}})$ in we find that $|\rho _{bb}|\sim 0.2$ is still allowed at 95% CL for all scanned points with $m_H/m_A$ $>300$ GeV. ATLAS preformed a similar search [86] but the limits are somewhat weaker. The CMS search for light resonances decaying into $b{\bar{b}}$ [87] provides limits covering also $m_H/m_A = 200$ GeV, however, the constraint are weaker than Ref. [85] for all scanned points. This illustrates that the current exclusion limits are much weaker than our working assumption $\rho _{bb}\sim \lambda _b$. Same is also true for $\rho _{\tau \tau }$ i.e., all scanned points are allowed if $\rho _{\tau \tau }\sim \lambda _\tau$. Moreover, nonvanishing $c_\gamma$ may induce $H\rightarrow Z Z$, $H\rightarrow W^+ W^-$, $H\rightarrow \gamma \gamma$, $A\rightarrow \gamma \gamma$ etc., however, we have checked such decays are doubly suppressed via $c_\gamma \in [0,0.05]$ and large ${\mathcal {B}}(H/A \rightarrow t {\bar{t}})$ and ${\mathcal {B}}(H/A \rightarrow t {\bar{c}} + {\bar{t}} c)$. In general, we assumed off diagonal $\rho _{ij}$s to be much smaller compared to the diagonal elements in the corresponding $\rho$ matrices, however, $\rho _{tu}$ could still be large, with ${\mathcal {O}}(0.1-0.2)$ is still allowed for $m_A/m_H \gtrsim 200$ GeV [88]. Furthermore, if both $\rho _{tu}$ and $\rho _{\tau \tau }$ are nonzero $B\rightarrow \tau \nu$ decay could provide sensitive probe which could be measured by the Belle-II experiment [89]. We leave out a detailed analysis turning on all $\rho _{ij}$s simultaneously for future. We conclude that there exist sufficient room for discovery in near future while non-observation may lead to more stringent constraints on the parameter space.

5.2 Probing near mass degeneracy at the LHC

In this section we discuss how to probe the near degeneracy of $m_H$, $m_A$ and $m_{H^\pm }$ favored by inflation at the LHC. As discussed earlier, there exist exact cancellation between the $cg\rightarrow t H \rightarrow t t {\bar{c}}$ and $cg\rightarrow tA \rightarrow t t {\bar{c}}$ amplitudes if masses and widths are degenerate [76, 77]. The cancellation reduces if the mass splittings are larger, as can be seen from previous section. For the allowed values of $\rho _{tc}$ and $\rho _{tt}$ discussed above, the decay widths of H and A are also nearly degenerate. Therefore cancellation could be significant for the scanned points in Figs. 3 and 4. With semileptonically decaying same-sign top signature, $cg\rightarrow tH/tA \rightarrow t t {\bar{c}}$ can be discovered at the LHC, even with full Run-2 dataset, unless there exist such cancellation [76, 77].

Note that such cancellation does not exist between $cg\rightarrow tH \rightarrow t t {\bar{t}}$ and $cg\rightarrow tA \rightarrow t t {\bar{t}}$ processes [76] if $m_H$ and, $m_A$ are above $2m_t$ threshold. Induced by $\rho _{tc}$ and $\rho _{tt}$ couplings, the processes $cg\rightarrow tH/tA \rightarrow t t {\bar{t}}$ can be discovered in the Run-3 of LHC if $m_A$ and, $m_H$ are in the sub-TeV range [76]. In general, it is expected [76] that $cg\rightarrow tH/tA \rightarrow t t {\bar{c}}$ (same-sign top signature) would emerge earlier than the $cg\rightarrow tH/tA \rightarrow t t {\bar{t}}$ (triple-top signature). For sizable $\rho _{tc}$ and $\rho _{tt}$ one may also have $cg\rightarrow bH^+\rightarrow b t {\bar{b}}$ [90] process which can also be discovered at the LHC as early as in the Run-3. Hence, we remark that vanishing or small same-sign top and, sizable triple-top and $cg\rightarrow bH^+\rightarrow b t {\bar{b}}$ signatures at the LHC would provide smoking gun signatures for the inflation in g2HDM. We leave out a detailed study regarding the discovery potential of these processes in the context of inflation for future.

6 Discussion and summary

We have investigated inflation in g2HDM in the light of constraints arising from collider experiments. We have primarily focused on the two benchmark scenarios. In Scenario-I we assumed nonminimal coupling $\xi _{11}$ to be nonvanishing while in Scenario-II we assumed $\xi _{22}$ nonzero. In both cases the parameter space favored by inflation require the nonminimal coupling ${\mathcal {O}}(10^3{-}10^4)$. We find that parameter space preferred by inflation requires $m_H$, $m_A$ and $m_{H^\pm }$ to be nearly degenerate.

While finding the available parameter space we turned on only one nonminimal coupling at a time. This is primarily driven by the fact that one nonminimal coupling is sufficient to account for all the constraints from Planck data 2018. Throughout we set $\xi _{12}=0$ in our analysis. We find that a similar parameter space for $\xi _{12}$ can be found. We leave out a detailed analysis where all three nonminimal couplings are nonzero for future.

There exist several direct and indirect constraints for the parameter space. The most stringent constraints on the additional Yukawa couplings $\rho _{tt}$ arise from h boson coupling measurements by ATLAS [51] and CMS [52] as well as from heavy Higgs searches such as $bg\rightarrow {\bar{t}} H^+$ [66], $gg\rightarrow A/H \rightarrow t{\bar{t}}$ [71, 72], and $gg\rightarrow t {\bar{t}} A/H \rightarrow t {\bar{t}} t{\bar{t}}$ [73]. The most stringent indirect constraints arise from $B_{d,s}$ meson mixings. We found that $\rho _{tt}\approx 0.5$ is allowed by current data for $m_H$, $m_A$, and $m_{H^\pm }$ for 200–800 GeV. On the other hand the most stringent constraint on $\rho _{tc}$ arise from the control region of $t{\bar{t}}W$ background of CMS search for SM four-top production [73]. We find that $\rho _{tc}\sim 0.2$ are well allowed by current data.

The near degeneracy of $m_H$ and $m_A$, as preferred by inflation, would lead to small same-sign top $cg\rightarrow t H/tA\rightarrow tt {\bar{c}}$ signature, while triple-top $cg\rightarrow t H/tA\rightarrow tt {\bar{t}}$ cross sections could be large. One expects same-sign top to emerge earlier than triple-top, that is unless $m_H$ and $m_A$ are degenerate or nearly degenerate [76]. One may also have $cg\rightarrow bH^+\rightarrow b t {\bar{b}}$ signature which could be discovered as early as in the Run-3 of LHC. Together they will provide unique probes for the inflation in g2HDM at the LHC if $m_H$, $m_A$, and $m_{H^\pm }$ are sub-TeV. Future lepton colliders such as ILC and FCC-ee might also provide sensitive probes to the parameter space. E.g., if $c_\gamma$ is nonzero one may have $e^+ e^-\rightarrow Z^* \rightarrow A h$, followed by $A\rightarrow t {\bar{t}}$ (or $A\rightarrow t {\bar{c}}$) with $h \rightarrow b {\bar{b}}$. This would be studied elsewhere. The future updates of $B_{d,s}$ mixing or, ${\mathcal {B}}(B\rightarrow X_s \gamma )$ of Belle-II [92] could also relevant.

In our analysis we have assumed $\rho _{ij}$ and $\lambda _i$s to be real for simplicity. In general $\rho ^F_{ij}$, $\mu _{12}^2$, $\lambda _5$, $\lambda _6$ and $\lambda _7$ could be complex. We however briefly remark that such complex couplings receive stringent constraints from electron, neutron and mercury electric dipole moment (EDM) measurements [93, 94]. In this regard, asymmetry of CP asymmetry ($\varDelta {\mathcal {A}}_{\text {CP}}$) of charged and neutral $B\rightarrow X_s \gamma$ decays could be relevant [94, 96] even though the observable has associated hadronic uncertainties. The future Belle-II measurement of $\varDelta {\mathcal {A}}_{\text {CP}}$ [92] could reduce the available parameter space for imaginary $\rho _{tt}$ [94,95,96]. Moreover, we set all $\lambda _{i}$s and $\rho _{ij}$s to zero except for $\lambda _t$, $\rho _{tt}$ and $\rho _{tc}$ and, assumed $\rho _{ii}$ could be $\sim \lambda _i$ with suppressed off diagonal $\rho _{ij}$s. If such coupling structure is realized in nature, we find that couplings other than $\lambda _t$, $\rho _{tt}$ and $\rho _{tc}$ have inconsequential effects in inflationary dynamics.

In summary, we have analyzed the possibility of Higgs inflation in general two Higgs doublet model. We find that parameter space for inflation favors nearly degenerate additional scalars. The sub-TeV parameter space receives meaningful constraints from direct and indirect searches. We also find that parameter space required for inflation could be discovered in the future runs of LHC as well as the planned ILC, FCC-ee, etc., while indirect evidences may emerge in flavor factories such as Belle-II. A discovery would not only confirm beyond Standard Model physics, but may also provide unique insight on the mechanism behind inflation in the early Universe.