Extended Higgs sector of 2HDM with real singlet facing LHC data

We study the phenomenology of the two Higgs doublet model with a real singlet scalar $S$ (2HDMS). The model predicts three CP-even Higgses $h_{1,2,3}$, one CP-odd $A^0$ and a pair of charged Higgs. We discuss the consistency of the 2HDMS with theoretical as well as with all available experimental data. In contrast with previous studies, we focus on the scenario where $h_2$ is the Standard Model (SM) 125 GeV Higgs, while $h_1$ is lighter than $h_2$ which may open a window for Higgs to Higgs decays. We perform an extensive scan into the parameter space of 2HDMS of type I and explore the effect of the singlet-doublet admixture. We found that a large singlet-doublet admixture is still compatible with the recent Higgs data from LHC. Moreover, we show that $h_1$ could be quasi-fermiophobic and would decay dominantly into two photons. We also study in details the consistency of the non-detected decay of $h_2\to h_1 h_1$ with LHC data followed by $h_1\to \gamma \gamma$ which leads to four photons final state at LHC: $pp\to h_2\to h_1 h_1\to 4\gamma$. Using the results of null searches of multi-photons carried by the ATLAS collaboration, we have found that a large area of the parameter space is still allowed. We also demonstrate that various neutral Higgs of the 2HDMS could have several exotic decays.


Introduction
A Higgs-like particle has been discovered in the first run of Large Hadron Collider (LHC) with 7 and 8 TeV energy in 2012 [1,2] and some of its properties, such as its decay to some Standard Model (SM) particles, have been established. During the second run of LHC with 13 TeV center of mass energy, some Higgs-like decay measurements get improved and new observable such as pp → tth was performed [3,4].
The aforementioned Higgs-like properties established so far will be further improved by the High Luminosity program of the future LHC run (HL-LHC). At the HL-LHC, one can pin down the uncertainties on the Higgs-like couplings to a few percent level in some cases [5] and provide indirect hints to the awaited new physics. Moreover, in the clean environment of the e + e − collider such as ILC and CEPC, which can act like a Higgs factory, one can improve the Higgs-like properties [6,7].
Although all data collected by LHC seems to indicate that the Higgs-like particle is in perfect agreement with the SM predictions, there are theoretical as well as experimental hints that indicate that the SM is only an effective field theory of a more fundamental one. One common feature of those Beyond SM (BSM) theories is an extended Higgs sector with an extra singlet, doublet, and/or triplet.
Most of the higher Higgs representations with an extra doublet/singlet predict in their spectrum extra neutral and/or charged Higgs states. Discovery of another Higgs state or more would be seen as a clear evidence of an extended Higgs sector and a departure from the SM.
Following the discovery of a Higgs-like particle, there has been a large amount of works dedicated to extending the SM Higgs sector by extra Higgs fields. Among the simplest one we mention: A real singlet scalar that has a mixing with the SM Higgs boson [8], the popular Two Higgs Doublet Model (2HDM) with or without natural flavor conservation [8], the popular two Higgs Doublet Model (2HDM) with or without natural flavor conservation [9] and the inert Higgs model that provides dark matter candidate [10,11]. Recently, there have been also phenomenological studies that extend the 2HDM with an additional real gauge-singlet scalar which act as dark matter candidate [12,13]. One can also extend the 2HDM by adding a real scalar singlet with non-vanishing expectation value that can mix with the doublets [14,15], a model which we call 2HDMS. In the two variants of the 2HDMS, the scalar spectrum is richer than the traditional 2HDM, which imply an interesting phenomenology at colliders including but not restricted to scalar-to-scalar decays, exotic decays and fermiophobic scenarios which are precluded in the SM and occurs hardly in the 2HDM. The model can easily accommodate a SM-like Higgs Boson that easily satisfies all the constraints from LHC measurements. Despite the existence of mixing among CP-even mass eigenstates, which would modify the SM-like Higgs couplings to fermions and gauge bosons, constraints from signal strength measurements can be easily satisfied (within the present range of systematic and statistical errors).
In the 2HDMS, the scalar spectrum contains 3 CPeven states h 1,2,3 , one CP-odd A and a pair of charged Higgs H ± . h 1,2,3 are admixtures of doublet and singlet components while A and H ± are pure doublet Higgs. A comprehensive analysis of the 2HDMS has been carried by the authors of Ref. [15] assuming that the lightest CP-even scalar (h 1 ) is the SM-like Higgs boson. In this scenario, and while satisfying all theoretical and experimental constraints, a large singlet-doublet admixture is still compatible with LEP and LHC data [15]. In the present study, we would like to do a comprehensive study for the case where h 2 is the 125 GeV Higgs while h 1 is lighter than h 2 which open a window for Higgs to Higgs decay such as h 2,3 → h 1 h 1 , AZ, W ± H ∓ and also h 3 → h 1 h 2 , h 2 h 2 , which are still compatible with Higgs data.
The paper is organized as follow: section 2 is devoted to the 2HDMS and its parametrization, in section 3 we review the theoretical and experimental constraints that the 2HDMS is subject to. We present our numerical result in section 4 and conclude in section 5. Several technical details are given in the Appendix. 2 The 2HDM with a real singlet field: 2HDMS In this section, we present a review of the 2HDMS. We discuss the scalar potential and derive the spectrum and the parametrization of the model. We also present the Yukawa textures and discuss the natural flavor conservation of the model. Couplings of the Higgs bosons to gauge bosons are also shown and their sum rules are discussed.

The Higgs potential
The scalar sector of 2HDMS consists of two weak isospin doublets H i (i = 1,2), with hypercharge Y = 1 and a real singlet field with hypercharge Y = 0 which are given by The most general renormalizable scalar potential for the model that respect SU (2) L ⊗ U (1) Y gauge symmetry has the following form: where m 2 11 , m 2 22 and m 2 S are masses terms. By hermiticity of the scalar potential λ 1,2,3,4,6,7,8 are dimensionless real parameters while λ 5 and µ 2 can be complex to allow CP violation in the scalar sector. In the present study, we assume that all scalar parameters are real. Therefore, the only source of CP violation is in the Cabbibo-Kobayashi-Maskawa matrix. We remind here that we allow a dimension 2 term µ 2 which break softly Z 2 symmetry. This discrete Z 2 symmetry is usually imposed in order to avoid Flavor Changing Neutral Current (FCNC) at tree level in the Yukawa sector.
Assuming that spontaneous electroweak symmetry breaking (EWSB) is taking place at some electrically neutral point in the field space, and denoting the corresponding VEVs by The parameters m 2 11 , m 2 22 and m 2 S can be eliminated by the minimization conditions of the potential Eq.(2): where s x , c x , t x stand for sin x, cos x and tan x respectively, and λ 345 = λ 3 + λ 4 + λ 5 and t β = v 2 /v 1 . After the EWSB of SU (2) L ⊗ U (1) Y down to electromagnetic U (1), three of the nine Higgs degrees of freedom corresponding to the Goldstone bosons are absorbed by the longitudinal components of vector boson W ± and Z 0 . The remaining six degrees of freedom should manifest as physical Higgses: three CP-even scalars (h 1 , h 2 , h 3 with m h1 < m h2 < m h3 ), one CP-odd A and a charged Higgs pair H ± .

Higgs masses and mixing angles
The most general form of the squared mass matrix 7 × 7 of the Higgs sector can be recast, using Eqs. (4,5,6), into a block diagonal form of three submatrices: one 3 × 3 matrices denoted in the following by M 2 CPeven for CP-even sector, one 2 × 2 matrix M 2 CP odd for CP-odd sector and one 2 × 2 matrix denoted by M 2 ± for the charged sector. The squared mass matrix for the charged fields φ ± 1,2 is: with λ + 45 = λ 4 + λ 5 . This matrix is diagonalized by the following orthogonal matrix R β , given by : Among the two eigenvalues of M 2 ± , one is zero and corresponds to the charged Goldstone boson G ± while the other one corresponds to the charged Higgs boson H ± and is given by: The charged Higgs H ± and the charged goldstone G ± are orthogonal rotation of the weak eigenstates φ ± 1 , φ ± 2 , The neutral scalar and pseudo-scalar mass matrices are given by: and The physical states h i = {h 1 , h 2 , h 3 } are obtained by an orthogonal transformation h i = R α1,2,3 φ i , (i = 1, 2, s) that diagonalizes the mass matrix M 2 CPeven , with : with s i = sin α i c i = cos α i . Without loss of generality we assume that m h1 < m h2 < m h3 .
From Eq. (12), it is easy to get the two eigenvalues of M 2 CP odd , one is vanishing and corresponds to the neutral Goldstone boson G 0 while the other one corresponds to the pseudo-scalar A: The CP-odd state A and the neutral Goldstone G 0 are obtained by an orthogonal rotation of the weak eigenstates χ 1 , χ 2 :

Yukawa texture
There are different types of Higgs couplings to fermions. If we do like in the SM and allow both Higgs fields to couple to all fermions through the following lagrangian: where Q 0 L is the weak isospin quark doublet, U 0 R and D 0 R are the weak isospin quark singlets and η U,0 1,2 , η D,0 1,2 are matrices in flavor space, then the above lagrangian will generate Flavor Changing Neutral Currents (FCNC) at the tree level which can invalidate some low energy ob-servables in B, D and K physics. In order to avoid such FCNC, it is customary to invoke a Z 2 symmetry that forbids FCNC couplings at the tree level [16]. Depending on the Z 2 assignment, we have four type of models [17]. In the present study, we focus only on type-I where all fermions couple only to one of the two Higgs doublets. In this case: Neutral Higgs couplings to a pair of fermions are: where f designate any type of fermions.

Higgs couplings to gauge bosons and sum rules
We present shortly here the Higgs couplings to gauge bosons and discuss the sum rules required by unitarity [18,19] that are subject to. The normalized couplings of neutral Higgs to a pair of gauge bosons V = Z, W are given by: which satisfy the following sum rule: This sum rule imply that each coupling g hiV V is requested to satisfy: | g hiV V |≤ 1.
For the couplings between two Higgs bosons and one gauge boson, we can distinguish two cases, a neutral case which corresponds to h i AZ vertex and charged case associated with h i H ∓ W ± vertex. From the kinetic terms of the Higgs fields, one can derive the various trilinear couplings among neutral, charged Higgses and gauge bosons. In units of (p hi − p A ) µ for neutral Higgs, respectively in units of λ c = ∓ g 2 (p hi − p H ± ) µ for charged Higgs, we have: where V=Z and S = A for the neutral case and V = W ± and S = H ∓ for charged one.
In the 2HDMS, one can easily derive the following sum rules: Where R i3 is the singlet component of the Higgs h i . An other type of sum rule which relate h i f f and h i V V can be derived from the Feynman rules and is given by [20]: where g hiV V and g hif f are the normalized couplings of h i to gauge bosons and fermions. From above, it follows that: and h i AZ must be very suppressed, and this will present a real challenge for the production and detection of such Higgs bosons.
then both singlet component R i3 as well as g hiSV couplings must vanish. This scenario could happen only when h i have no singlet component.
This would imply from Eq.(29) that the reduced coupling to fermions must satisfy g hif f = 1.

Theoretical and experimental constraints
The Two Higgs Doublets Model plus a Singlet possesses a large freedom in the scalar sector, coming from the large number of free parameters of the scalar potential. In order to obtain a viable model, many theoretical constraints have to be imposed on the scalar potential like perturbative unitarity, vaccum stability and electroweak precision observables. In what follows, we will describe briefly these constraints.

Boundedness from below (BFB) of the potential
In order to ensure a stable vacuum, the scalar potential has to be bounded from below in any directions in the field space as the field strength becomes extremely large. At large field values, the scalar potential is fully dominated by quartic couplings whose the BFB will depend only. At large field strength, the potential defined by Eq.(2) is generically dominated by the quartic terms: The study of V (4) (H 1 , H 2 , S) will thus be sufficient to obtain the main constraints. The full BFB constraints reads as for λ 7 > 0 and λ 8 > 0.
If λ 7 or λ 8 < 0, we have to satisfy two additional constraints: Full technical details on the proof of these constraints can be found in Appendix.(A).

Perturbative unitarity
To constrain further the scalar potential parameters of the 2HDMS one can ask that tree-level perturbative unitarity is preserved for a variety of scattering processes: gauge boson-gauge boson scattering scalar-scalar scattering, and also scalar-gauge boson scattering. Moreover, according to the equivalence theorem which states that at high energy limit √ s the amplitudes of a scattering process involving longitudinally polarized gauge bosons V are asymptotically equal, up to correction of the order m V / √ s, to the corresponding scalar amplitudes in which longitudinally polarized gauge bosons are replaced by their corresponding Goldstone bosons. We conclude that perturbative unitarity constraints can be implemented by considering pure scalar-scalar scattering only.
In order to derive the perturbative unitarity constraints on the scalar parameters of 2HDMS we follow refs [21,22]. According to [21,22], one computes the scattering amplitude in the weak eigenstate basis where the quartic couplings are less involving (does not involve mixing angles α i and β). The important point is that the amplitude expressed in the mass eigenstate fields can be transformed into the amplitude for the non-physical fields by making a unitary transformation. The eigenvalues for the scattering amplitude should be unchanged under such a unitary transformation.
In the Appendix.(B) we present the technical details of the different 2 → 2 scattering amplitudes. The explicit forms of the eigenvalues at tree level are given by: Other eigenvalues are coming from the cubic polynomial equation associated to the submatrix M 2 corresponds to scattering with one of the following initial and final states: Moreover, we also force the potential to be perturbative by imposing that all quartic couplings of the scalar potential satisfy |λ i | ≤ 8π (i = 1, ..., 8).

Electroweak precision test observables (EWPT)
The oblique parameters S, T, and U are known to provide an indirect probe of new physics BSM for theories that process SU (2) × U (1) symmetry [23]. These parameters quantify deviations from the SM in terms of radiative corrections to the W, Z and the photon self-energies. In the framework of 2HDMS, the Higgs doublet couples to the W and Z gauge bosons via the covariant derivative. Due to singlet and doublet admixtures in the scalar sector, the singlet field will also couple to the gauge bosons W and Z. Therefore, both neutral Higgs h i , A and charged Higgs will contribute to S and T parameters which are very well constrained by electroweak precision test observables. These EWPT constraints will be converted to limit on the mixing angles and/or masses splitting among the 2HDMS spectrum. The extra-contribution to S, T and U parameters for 2HDMS are given in Appendix. (C).
In order to study the correlation between S and T, we perform the χ 2 test over the allowed parameter space of 2HDMS. Our χ 2 S,T is defined as: where S and T are the computed quantities within 2HDMS framework [24,25,26].Ŝ andT are the measured values of S and T,σ 1,2 are their one-sigma errors and ρ their correlation [27], S = 0.04 ± 0.11, T = 0.09 ± 0.14, ρ S,T = 0.92 (37) It is worth noting here that we have checked the limits on the oblique parameters with the THDM, in this sense, our results match exactly to those outlined in [28,29].
In addition, we have indirect experimental constraints from B physics observables on the contribution of the 2HDMS such as tan β and m H ± . In the 2HDMS, the charged Higgs coupling to fermions is not at all affected by the singlet component of the additional Higgs. Therefore, constraints from B → X s γ and B q mixing will be the same as for the usual 2HDM model. We remind the reader that the recent experimental results presented by the Heavy Flavor Averaging Group (HFAG) [30] of B(B → X s γ) have changed in a significant way the bounds on the charged Higgs boson mass. For instance, in the Type-II 2HDMS, the measurement of the BR(B → X s γ) constrains charged Higgs mass to be larger than about 570 GeV [31], while in Type-I 2HDMS one can still obtain a charged Higgs boson with a mass as low as 100 − 200 GeV provided that tan β ≥ 2.

Constraints from Higgs data
Both ATLAS and CMS experiments of the LHC Run1 with 7 and 8 TeV and Run2 with 13 TeV confirmed the discovery of a Higgs-like particle with a mass around 125 GeV. Both groups performed several measurements on the Higgs-like particle couplings to the SM particles. Recently, both ATLAS and CMS Collaborations have announced the observation of Higgs bosons produced together with a top-quark pair [3,4]. All these measurements seem to be in perfect agreement with SM predictions.
In the case of 2HDMS, all tree level Higgs couplings to fermions and gauge bosons are modified with the mixing parameters α i and β. The loop mediated processes such as gg → h i , h i → γγ and h i → γZ will be affected both by the mixing angles as well as by the additional charged Higgs H ± loops which depend on the triple scalar coupling To study the effects of ATLAS and CMS measurements on 2HDMS, we take into account experimental data from the observed cross section times branching ratio divided by the corresponding SM predictions for the various channels, i.e. the signal strengths of the Higgs boson defined by: where i stand for different modes of Higgs production. The dominant mechanisms of Higgs production are gluon fusion (ggF), followed by vector boson fusion (VBF), Higgsstrahlung (Vh) and associated production with top-quark pairs (tth). All these various signal strength channels are included in our analysis through the public code HiggsBounds and HiggsSignals which also include previous LEP and Tevatron experimental searches.
As said previously, in our analysis, we will assume that h 2 is the 125 GeV Higgs-like particle discovered while h 1 would be lighter than h 2 . Therefore, once the decay channels h 2 → h 1 h 1 and/or h 2 → AA are open, the subsequent decays of h 1 /A into fermions, photons or gluons, will lead either to invisible or undetected h 2 decays that can be constrained by using global analysis to the present AT-LAS and CMS data to Higgs couplings.
We stress here that there is also searches for nondetected decays of the SM Higgs boson both by ATLAS and CMS. CMS look for the following SM Higgs production channels: gluon fusion, vector boson fusion, and Higgsstrahlung process pp → V H (V=W or Z) with subsequent invisible Higgs decays. Upper limits are placed on Br(H → invisible), as a function of the assumed production cross sections. The combination of all the above channels, assuming SM production, yields to an upper limit of 0.24 on the BR(H → invisible) at the 95% confidence level [32]. ATLAS collaboration performs a search for an invisible decay of the Higgs through pp → ZH process with a leptonic subsequent decay of the Z [33]. Their limit is slightly weaker than CMS results.
In our study, we will use the fact that the total branching fraction of the SM-like Higgs boson into undetected BSM decay modes is constrained, as mentioned, by BR(H → invisible) ≤ 0.24 where BR(H → invisible) designate BR(h 2 → h 1 h 1 ) or the sum of BR(h 2 → h 1 h 1 ) and BR(h 2 → AA), BR(h 2 → Z * A) and BR(h 2 → W * H ± ) if the later is open.

Parameters scan
The scalar potential Eq.(2) have 15 independent parameters: four masses, 8 quartic couplings λ 1,...,8 and 3 vaccum expectation values. Three masses can be eliminated by the use of the 3 minimization conditions Eq. (6). Moreover, after electroweak symmetry breaking, from the kinetic terms of the Higgs doublets, the W and Z gauge bosons acquire masses which are given by m 2 W = 1 2 g 2 (v 2 1 + v 2 2 ) and m 2 Z = 1 2 (g 2 + g 2 )(v 2 1 + v 2 2 ), where g and g are the SU L (2) and U Y (1) gauge couplings. The combination v 2 1 + v 2 2 is thus fixed by the electroweak scale through the well known relation v 2 = v 2 1 + v 2 2 = (2 √ 2G F ) −1 , and we are left with 11 free parameters. By simple algebraic calculations, from the mass matrix relations, one can express all the quartic couplings λ i as a function of the physical masses, µ 2 , tan β and the mixing angles α i 1 . One can then take the following set of free independent parameters: Note that the usual 2HDM is recovered from 2HDMS by taking the following limits: In the previous studies on the 2HDMS [34,35], they concentrate on h 1 being the 125 GeV SM-like Higgs while in our analysis we will study the consistency of having the second Higgs h 2 as the 125 GeV SM-like Higgs. This scenario open more decay channels for h 2 such as h 2 → h 1 h 1 and probably h 2 → AA, AA * , Z * A, W * H ± .
In order to display the allowed regions for the parameters, we have considered both of the exclusions from both HiggsBounds-5.1beta and HiggsSignals-2.1.0beta to compute the value of χ 2 min considering the combination of 8 TeV and 13 TeV Higgs signal strength data from run-I and run-II. Thus, we show the best fit at 95.5% C.L, which corresponds to 2.3 ≤ ∆χ 2 (χ 2 − χ 2 min ) ≤ 5.99. In the left panel of Figure.(1), we show the allowed parameter for λ 7,8 for various values of λ 6 taking into account perturbative unitarity and BFB constraints. We note first that there is a complete overlap between the three colors for positive λ 7,8 . One can see that for negative λ 7,8 , the theoretical constraints restrict their strength for small λ 6 . The restriction is relaxed for large λ 6 ≈ 4π.
In the right panel of Figure.(1), we present the correlation between S and T, after taking into account the theoretical constraints and the exclusion from Higgs Bounds at 95% C.L. The red, green and black points represent the points with ∆χ 2 value within the 1σ, 2σ and 3σ interval.
In order to be consistent with the EW precision measurements, such light h 1 and A are naturally accompanied also by a light charged Higgs which is consistent with B → X s γ constraint in 2HDMS of type I.
We study first the decay of h 1 into SM particles. As one can see from the couplings of h 1 to fermions given in Eq. (19), h 1 could be fermiophobic if R 12 ∝ cos α 2 sin α 1 vanish. This scenario might happens if we take α 1 ≈ 0 and/or α 2 ≈ π/2 which is possible since both α 1 and α 2 are free parameters in this model.
In Figure.(2), we show the branching ratio of h 1 → W + W − as a function of h 1 → γγ with Br(h 1 → bb + τ + τ − ) represented on the horizontal axis on the left panel while on the right panel we show m h1 . Since m h1 ≤ 125 GeV, h 1 → W + W − will proceed with one or both W being off-shell. We first mention that the singlet component of h 1 does not exceed 50% in our case which makes h 1 dominated by doublet components.
As one can see, in most of the cases h 1 would decay significantly into a bottom pair unless α 1 vanish which is the fermiophobic limit for h 1 . In this case, it is clear that h 1 → γγ reach its maximum value when Br(h 1 → bb) and Br(h 1 → τ + τ − ) are very suppressed. When h 1 is fermiophobic, In what follow discuss only h 1 → W W since h 1 → ZZ is smaller. In fact, h 1 → W * W * which is open for m h1 < m W is very suppressed due to phase space while for m h1 ≥ m W , h 1 → W W * is open and could strongly compete with h 1 → γγ. This is shown in the left and right panel of Figure.(2) where we can see Br(h 1 → W W * ) as a function of Br(h 1 → γγ). Close to the fermiophobic limit where h 1 → bb and h 1 → τ + τ − are suppressed, if the mass of h 1 is larger than the W boson mass then h 1 → W W * can dominate h 1 → γγ in some cases.
One could have also the following scenario: both h 1 → f f , h 1 → V V * and h 1 → γγ, γZ are rather small while the branching ratio of Br(h 1 → AZ * ) + Br(h 1 → H ± W * ∓ ) becomes significant which can be understood from the sum rules eqs. (28). Due to the smallness of h 1 V V coupling and R 13 component, the sum rule eqs. (28) imply that h 1 AZ and h 1 W ± H ∓ could becomes significant. This is illustrated in the left panel of Figure.(2) with the blue dots in the left-down corner where both Br(h 1 → W W * ) ≈ Br(h 1 → γγ) ≈ 10 −4 and also Br(h 1 → τ + τ − + bb) are rather small.
The above configuration is illustrated clearly in Figure.(3) where we show the correlation between Br(h 1 → τ + τ − +bb) and Br(h 1 → γγ +ZZ +W + W − ) as a function of Br(h 1 → H ± W ∓ + AZ) on the horizontal axis.
It can be seen that when the fermionic (τ + τ − , bb) and bosonic (γγ, W + W − , ZZ) decays of h 1 are suppressed, the higgs to higgs decays h 1 → H ± W ∓ and/or h 1 → ZA becomes significant. In the right panel of Figure.(3) it can be read that h 1 is most of the case having large doublet component.
In Figure.(4) we illustrate κ h1 f as a function of κ h1 V with R 1i (i = 1, 2, 3) component of h 1 on the horizontal axis. From this plot, one can read that the doublet component is rather large in most of the cases leaving only small singlet component which is less than 50% . One can also learn that when κ h1 f and κ h1 V are suppressed, the doublet component is very large. Which means that h 1 is mainly coming from the doublet components. According to the sum rule Eq.(23), | κ h1 V |≤ 1 is requested which is consistent with Figure.(4). On the other hand, in large area of parameter space | κ f (h 1 ) |≤ 1 while one can have a scenario with | κ f (h 1 ) |≥ 1. This happens when h 1 f f ∝ R 12 is rather large (small R 11 ). On the left panel of Figure.(4) we show the sensitivity to tan β where we can see a linear correlation between κ h1 f and κ h1 V at large tan β. Let us mention that in this scenario with suppressed h 1 f f and h 1 V V couplings, h 1 can not be produced in the usual channel such as gluon fusions, vector boson fu- Fig. 3 We now discuss the decay of the SM-like h 2 . We first show the consistency of h 2 → V V and h 2 → f f with LHC data. For this purpose we illustrate in Figure.(5) the correlation between κ h2 V and κ h2 f as a function of R 2 2i . According to the sum rule Eq.(23), κ h2 V < 1, and this is clearly illustrated in the plot. One can see from the plot that when κ h2 V ≈ 1 we have also κ h2 f ≈ 1, this is a consequence of the sum rule Eq.(29). However, the suppression of κ h2 V could be of the order of 12% and it could happens both for κ h2 f < 1 or κ h2 f > 1. Note that the suppression of both κ h2 f and κ h2 V takes place when the singlet component of h 2 is rather large R 2 23 > 0.1. One can see that κ h2 f could reach a values less than 0.8 for R 2 23 ≈ 0.25. It is also clear from the plot that one can have an enhancement of κ h2 f in the range of [1.05 − 1.15] for small singlet component of h 2 (R 2 23 ≈ 0.1) and moderate tan β.
In Figure.(6) we show the correlation between κ h2 gg and κ h2 γγ on the left panel and the correlation between κ h2 γγ and κ h2 γZ on the right panel. The SM value is indicated as a black box. One can see that both in κ h2 gg , κ h2 γγ and κ h2 γZ deviation from SM value can reach 15%. Note that both in κ h2 γγ and κ h2 γZ one can see that in most of the case we have a suppression of the rate compared to SM. The figure shows also that for h 2 with relatively large singlet component we have suppression of κ h2 gg , κ h2 γγ and κ h2 γZ rate. We also stress that most of the cases κ h2 γZ < κ h2 γγ As we have seen in Figure.(5), decays of h 2 into SM particles such as W W , ZZ, bb and τ + τ − are consistent with LHC measurements with deviations from SM predictions that could goes up to 10-15%. However, these deviations are mainly due to experimental uncertainties on all the LHC measurements which could be larger than 10% in some channels. Therefore, taking into account these uncertainties, there is still a room for the non-detected SM higgs decays such as h 2 → h 1 h 1 , AA, AA * , AZ * , H ± W * . In our scan we assume that m A > 62.5 GeV, therefore h 2 → AA will not be open and h 2 → AA * is rather suppressed. We are only left with h 2 → h 1 h 1 , AZ * , H ± W * channels. As explained above, all these additional decays of the SM Higgs should not exceed 24% . We show in Figure. Note that the couplings h 1 AZ and h 1 W ∓ H ± are exactly the same (see Eq. (28)). Therefore, if m A ≈ m H± then Br(h 2 → Z * A) and Br(h 2 → W * H ± ) are of the same order. The total amount for Br(h 2 → h 1 h 1 ) + Br(h 2 → Z * A) + Br(h 2 → W * H ± ) should not exceed 24% as requested from the non-detected decay of the SM Higgs, and this is rather clear from Figure.(7). The plots also display the correlation between Br(h 2 → h 1 h 1 ) and Br(h 2 → Z * A)+Br(h 2 → W * ∓ H ± ). When Br(h 2 → h 1 h 1 ) is maximized, Br(h 2 → Z * A) + Br(h 2 → W * ∓ H ± ) is minimal and vice verse. One can have also a configuration where both Br(h 2 → h 1 h 1 ) and Br(h 2 → Z * A) + Br(h 2 → W * ∓ H ± ) are of the same size. In the case where both A and H ± are heavier than 125 GeV, only h2 → h 1 h 1 would contribute to the non-detected decay of h 2 .
On the middle panel of Figure.(7) it is clear that when the reduced coupling of h 2 h 1 h 1 is large, the branching ratio Br(h 2 → h 1 h 1 ) is substantial which would provide an important production channel for h 1 from h 2 decay: gg → h 2 → h 1 h 1 which could compete with the other production channels such as pp → W h 1 and/or pp → {h 1 A, h 1 H ± }. On the right panel of Figure.(7) we illustrate the correlation between Br(h 2 → Z * A) and κ h2 V as a function of κ h2AZ . As one can see from the plot, and according to the sum rule Eq.(28), when h 2 V V is full strength, then Br(h 2 → Z * A) is suppressed.    We have seen previously that Br(h 2 → h 1 h 1 ) could be significant and can reach 20% in some case. In the case where h 1 is dominated by the singlet component, it is well known that it is hard to produce it through the con-ventional channels such as ggF, VBF ect. Therefore, the process gg → h 2 followed by the decay h 2 → h 1 h 1 could be an important process for the production of h 1 . In the case where h 1 is dominated by the singlet component, its decay to SM particle would be suppressed. In such case, it may be possible that h 1 would decay to a pair of photons which could proceed through charged Higgs loops. Therefore, the process gg → h 2 → h 1 h 1 could lead to a spectacular 4 photons final states. In Figure.(8)(left) we illustrate the branching fraction Br(h 2 → h 1 h 1 ) × Br(h 1 → γγ) 2 as a function of m h1 . As can be seen, such branching fraction could reach 10% in some cases. On the right panel of Figure.(8), we show the production cross section for We note that for very small singlet component R 13 ≈ 0 where h 1 is fully dominated by the doublet components, one could have sizeable Br(h 2 → h 1 h 1 ) as it has been noticed in the usual 2HDM [36,37].
Recently, ATLAS publish their results for the search of new phenomena in events with at least three photons [38] based on 8 TeV CM energy with 20.3 fb −1 . This search was used to put constraint on an N-MSSM scenario which leads to four photons final states gg → H → a 1 a 1 → 4γ where a light pseudo-scalar, if dominated by singlet component, can decay fully into two photon with 100% branching ratio. Following this work, it has been demonstrated in [37] that the kinematic distributions for qq → H → a 1 a 1 → 4γ and qq → H → h 1 h 1 → 4γ with h 1 being CP-even are identical. Ref [37] also provide a projection for 14 TeV CM energy. Therefore results from [38] can be applied to our four photons final states. In Figure.(9), we present our predictions for pp → h 2 → h 1 h 1 → 4γ both for 8 TeV and 14 TeV together with the 8 TeV exclusion from ATLAS analysis. ATLAS projection for 14 TeV is also shown in the lower band. It is clear that some benchmark points are already excluded by the 8TeV data and the 14 TeV projection. However, several benchmarks are still alive.
We now discuss h 3 decays. We show in Figure. One can see that before the opening of the tt threshold, h 3 → W W could be the dominant decay mode of h 3 with a branching which can reach up to 80%, while h 3 → ZZ goes up to 20% and in such case h 3 → h 1 h 1 is suppressed. After opening of tt threshold, h 3 → tt can be slightly larger than 10% for large m h3 mass.
We now discuss higgs to higgs decays, such as h 3 → h 1 h 1 , h 1 h 2 , h 2 h 2 and h 3 → ZA, W ± H ∓ . In Figure.(11) (upper plot) we illustrate the branching ratio of h 3 → h 1 h 1 (left), h 3 → h 2 h 2 (middle) and their correlation (right). From the left panel, one can see that Br(h 3 → h 1 h 1 ) can be substantial and becomes the dominant decay mode, while from the middle panel it is clear that Br(h 3 → h 2 h 2 ) can reach 30% as a maximal value. In the case where Br(h 3 → h 1 h 1 ) is the dominant decay, then one can have a new production mechanism for h 1 , namely: pp → h 3 → h 1 h 1 . This production channel might be useful for the case where h 1 has large singlet component in which case it will be challenging to produce it in the conventional channels. In the lower panels of Figure.(11) we display the correlation between Br( It is clear that one can have a scenario where both Br(h 3 → h 2 h 2 ) and Br(h 3 → AA) are rather large with branching fractions of the order 40%. It is also clear that when Br(h 3 → h 1 h 1 ) and Br(h 3 → h 2 h 2 ) are suppressed then Br(h 3 → W W, ZZ) would become slightly large.
In Figure.(12) we show the branching fractions of h 3 → AZ and h 3 → H ± W ∓ versus R 2 33 . As one can see from the plots and according to the sum-rule Eq.(28), when h 3 is dominated by the singlet component R 33 ≈ 1, then both h 3 ZA and h 3 W ∓ H ± couplings are suppressed resulting in a small branching ratio for both channels. For R 33 away from 1, the branching fraction h 3 → AZ and h 3 → H ± W ∓ can be in the range of 10-40% in some cases.
Given that Br(h 3 → h 1 h 1 ) can be sizeable, one can look to the amount of production cross section that one can get from h 3 production followed by h 3 decay into a pair of h 1 . We illustrate in Figure.
One can see that the production rate of h 1 is large especially in the mass range 200 GeV m h3 250 GeV when the decay channel h 3 → h 1 h 1 is kinematically allowed and both h 3 → h 2 h 2 and h 3 → tt are closed. The same behavior is observed in the middle panel, where the magnitude in the cross section is larger in the mass range 250 GeV m h3 350 GeV when the decay channel h 3 → h 2 h 2 is opened and h 3 → tt mode is closed.

C Oblique Parameters
In order to study the contribution of the oblique parameters in the framework of 2HDMS, we have used the general expressions presented in [24,25,26] for the SU (2) × U (1) electroweak model with an arbitrary number of scalar SU (2) doublets, with hypercharge ± 1 2 , and an arbitrary number of scalar singlets.
In this study, we have four real fields that are related to the physical fields h 1,2,3 and A through, where R ij are defined in terms of the mixing angle β and the elements of the rotation matrix given by eq, R ij , as follows, We recall that the charged sector is the same as 2HDM, It contains only a pair of charged scalars H ± . As a result the charged field is related to the physical charged scalar field through the unit matrix.
Therefore, the oblique parameters S, T and U in the 2HDMS are given by :