Dominant production of heavier Higgs bosons through vector boson fusion in NMSSM

We study the features of the additional Higgs bosons in the Next-to-Minimal Supersymmetric Standard Model where the lightest beyond Standard Model Higgs boson does not dominantly couple to up-type quarks. The new state is dominantly singlet-like while it can also accommodate a small down-type Higgs component. The gluon-gluon fusion cannot be adequate enough for such a Higgs production. We show that the vector-boson fusion may become the leading production mechanism to probe this new scalar at the LHC. Using the existing 13 TeV LHC data for an integrated luminosity $36.1~fb^{-1}$, we show the LHC constraints on the parameter space. Finally, we also study the reach of the planned high luminosity LHC (${\mathcal L}=3~{\rm ab}^{-1}$ at $\sqrt s=$~14 TeV) and the proposed high energy upgrade of the LHC (${\mathcal L}=15~{\rm ab}^{-1}$ at $\sqrt s=$~27 TeV) to probe this singlet-like Higgs scalar.


Introduction
The Higgs boson of the Standard Model (SM) which is considered as one of the main motivations of the Large Hadron Collider (LHC) has already been discovered a few years ago [1,2] with mass m H 125 GeV. Accommodating the observed Higgs scalar has put significant constraints on the allowed parameter space of any existing model of beyond Standard Model (BSM) physics. In the context of the simplest Supersymmetric model (SUSY) -the Minimal Supersymmetric Standard Model (MSSM) [3,4], large stop masses and/or mixing have been required to generate a Higgs mass of 125 GeV (for a review see [5]) which may lead to dangerous charge and color breaking minima [6] and could produce a large fine-tuning on the allowed parameter space [7]. The MSSM has an enlarged Higgs sector with CP-even, CP-odd and charged Higgs bosons [5] which has a particle spectrum and couplings as in a generic Two Higgs Doublet Model of type II [8] but determined by the MSSM parameters [5] 1 .
In the Next-to-Minimal Supersymmetric Standard Model (NMSSM) (for reviews see [12,13]), a SM singlet superfieldŜ is introduced in addition to the MSSM Higgs superfieldsĤ d ,Ĥ u . The solution to the fine-tuning problem can be more easily addressed since the SUSY breaking scale can be made relatively low while satisfying the SM-like Higgs boson mass of 125 GeV. The µ term is generated dynamically by the vacuum expectation value (vev) of the real scalar component of the gauge singlet superfieldŜ. Here one can efficiently enhance the Higgs boson mass either via new tree level contributions [12,13,14] or through singlet-doublet mixings [15]. Thus larget mass and/or t L −t R mixings would not be necessary which can in turn significantly reduce the fine-tuning [16]. NMSSM offers interesting collider phenomenology as lighter SUSY particles are now allowed with a Higgs boson mass of 125 GeV and the fine-tuning criterion [17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35]. Similarly, a singlino-like lightest Supersymmetric particle (LSP) can quantitatively relax the lower bound on SUSY particles through additional SUSY cascades into singlino LSP [36,37,38]. In the scalar sector, a very light Higgs is not yet excluded by the LEP searches [39], provided its couplings to the SM particles are small enough. The extra singlet scalar may lead to new decays of Higgs bosons, particularly into two lighter Higgs scalars which can qualitatively be different than a two Higgs doublet Supersymmetric Standard model [40,41,42] like MSSM.
The leading-order production processes for the Higgs boson production can be divided as: (i) gg → H (ggF), (ii) vector boson fusion (VBF) and associated V + H production (iii) associated tt + H along with a single t + H production (see reviews like [43,44,45]). The common lore is that the production cross section for m H 1 TeV is primarily dominated via the gluon-gluon fusion [46] : gg → H through intermediate quark loops, in particular via the top quark. For m H 100 GeV, the associated-production gg, qq → W +H, Z +H can have significant effects in production. The vectorboson fusion processes W + W − , ZZ → H become important for m H 100 GeV and dominate the production processes for m H 1000 GeV [47]. Other processes like gg, qq → bbH, ttH and in association with a single t ort could also contribute to Higgs productions.
In the case of MSSM, it is known that production of the heavier Higgs scalar has also been dominated by the ggF for small or moderate values of tan β [5,48,49]. At high tan β, in addition to the ggF (where b quark loop now also contributes), the associated Higgs production with the b quark would also significantly contribute, thanks to strong enhancement of the Higgs couplings to the down type fermions. Consequently, the new Higgs can decay into τ τ which has been considered as an important search channel at the LHC. In the case of NMSSM, in the CP-even Higgs sector, in addition to a SM-like Higgs, one may observe a MSSM-like and also a new BSM Higgs scalar. In a specific parameter space, this extra Higgs scalar can be found to be almost independent of H u and determined mainly by S and H d . Due to a very stringent constraint on the CP-even BSM Higgs couplings to the vector bosons, in the phenomenological acceptable region, the new state is found to be mostly singlet like. Then for small tan β, we may observe a few new possibilities for such a CP-even BSM Higgs state H with a mass much smaller than the TeV scale (m H < 1 TeV): • First we review the role of gluon-gluon fusion in the context of production of the scalar H. It depends on the effective Hgg coupling which in turn depends on the masses of the squarks and the heavy fermions in the loops. In the parameter space where squarks are heavier (which has been considered here), the effective Hgg coupling may not receive significant contributions from SUSY scalars in the loops. Similarly, the contributions from the heavy fermion loops in the said coupling are also suppressed as the Higgs state can only have small couplings with the tt and the bb fermions. While the suppression in the Htt coupling is related to an insignificant H u component in H, the smallness of the Hbb coupling can be attributed to the small tan β and large singlet component [13] 2 . Thus in this specific region, the gluon-gluon fusion may not be adequate enough for such a Higgs production. In addition, the associated Higgs production with the b quark would also be insignificant 3 .
• We know that vector-boson fusion processes W + W − , ZZ → H can be unsuppressed if any of the doublet components in the Higgs is non-negligible. Indeed, this can be true in the NMSSM parameter space where ggF may not be sufficient enough to produce a Higgs boson.
Here VBF can be observed as the main production mechanism even for probing a lighter non SM-like Higgs boson.
• Interesting changes can be observed even in the decays of the new Higgs scalar in the NMSSM. H → τ τ and H → tt (when kinematically allowed) would now be extremely suppressed which are promising search channels for MSSM-like heavy Higgs boson. Note that the branching ratios (Brs) into τ τ can be primarily small if the heavier Higgs state is dominantly singlet like. Final states with 2 gauge bosons and/or Higgs→ Higgs decays (when kinematically allowed) may become important to search such a Higgs scalar at the LHC.
The aim of this work would be to study (i) the status of such a Higgs scalar in the light of existing 13 TeV LHC data at 36.1 f b −1 integrated luminosity, (ii) prospects of probing this new scalar in the context of planned high luminosity run of the LHC (L = 3 ab −1 at √ s = 14 TeV) and at the proposed high energy upgrade of the LHC (L = 15 ab −1 at √ s = 27 TeV). We note in passing that there already exists a number of models, where neutral Higgs bosons with suppressed or vanishing couplings to the SM fermions are present at the weak scale. These are commonly known as fermiophobic Higgs bosons and have many generic features which have been considered in the literature [50]. They have already been searched by experimental collaborations [51]. They can be accommodated in type-I two Higgs doublet models [52]. However, our work considers Supersymmetric type-II two Higgs doublet models augmented by a singlet and our conclusion can be generalized to a non Supersymmetric version as well.
The paper is organised as follows. In Sec.2 we primarily discuss NMSSM and its Higgs sectors. In Sec.3 we show parameter space where rate of production of a new singlet dominated Higgs scalar 2 We note in passing that this is hardly possible in the MSSM as the bilinear mixing renders the 2 CP-even Higgses as the admixtures of H d and Hu. Also satisfying Higgs mass at small tan β would need large SUSY scale in MSSM [48] 3 For moderate or relatively larger values of tan β, the situation may be a bit different. Radiative corrections may enhance the non SM-like Higgs couplings to the b quarks which may, in turn, push both ggF and associated Higgs production with the b quark. Similarly, lighter squarks could also boost the ggF contributions to the Higgs production rate at the LHC. This part of the parameter space would not be considered here.
via the VBF channel can be comparable or even dominate over the gluon-gluon fusion process. Sec.4 describes our predictions at the √ s = 14 and √ s = 27 TeV runs of the LHC. Finally we conclude in Sec.5.

Higgs sector in the Next to minimal Supersymmetric Standard Model
In the scale invariant NMSSM, superpotential W can be read as [13]: In the above, the Yukawa couplings h u , h d , h e and the superfields Q, U c R , D c R , L and E c R should be understood as matrices and vectors in family space, respectively. The vacuum expectation value (vev) s of the real scalar component ofŜ generates an effective µ-term which in turn solves the µ-problem of the MSSM.
In the soft-SUSY breaking sector, apart from the L sof t MSSM which contains soft-SUSY breaking terms of the MSSM except the Bµ term, there are new trilinear interactions as well as mass terms involving the singlet field.
Compared to the MSSM, the gauge singlet superfieldŜ augments new degrees of freedom to CPeven and CP-odd Higgs sectors. Hence the spectrum contains (1) 3 CP-even neutral Higgs bosons H i , i = 1, 2, 3, (2) 2 CP-odd neutral Higgs bosons A 1 and A 2 , (3) one charged Higgs boson H ± . At the tree level, one may characterise the Higgs sector by the six parameters, namely, Then after eliminating m 2 H d , m 2 Hu and m 2 S using the minimization equations of the potential, one can read the elements of the 3 × 3 CP-even mass matrix in the basis (H dR , H uR , S R ) [13,28].
In the above, one may use ) and In the neutral CP-even Higgs sector, one of the lighter Higgs scalar H 1 /H 2 can resemble with SM-like Higgs boson whose mass can be expressed as where δ mix [15] represents singlet-doublet mixing. Clearly, for small tan β (≤ 2), there can be an enhancement even at the tree level via ∼ λ 2 sin 2 2β term (see Eq:7) while for larger tan β, δ mix could play the crucial role to obtain SM-like Higgs boson mass ∼ 125 GeV, provided there is a lighter singlet. In the present study, the lightest NMSSM Higgs H 1 becomes SM-like (∼ H 125 ) while the lighter and the heavier BSM Higgs bosons H 2 and H 3 are approximately singlet and MSSM-like respectively. Since the squark masses are assumed to be large in this analysis, the lightest SM-like Higgs boson may receive handful contributions from the radiative corrections. Before we present numerical results, here, for completeness, we list the couplings of the lighter BSM Higgs H 2 with the quarks and the gauge bosons at the tree level [13].
Here S 2,i (i = 1, 2, 3) refers to down, up and singlet components. As already stated, up-type Higgs component in H 2 is considered to be extremely tiny i.e., |S 2,2 | ∼ 0. Thus H 2 couplings with the up-type quarks are vanishingly small. The same couplings with the down-type quarks are small because of singlet dominance (though |S 2,1 | >> |S 2,2 |) and also for small values of tan β. As explained earlier, the couplings of our interest would be H 2 V V and its maximum allowance where V refers to any gauge boson. Now defining reduced couplings C H i V V as the ratio of Higgs couplings to the vector bosons relative to the corresponding couplings of the SM-like Higgs boson, one gets Demanding a SM-like Higgs boson in the spectrum, the said condition leads to a small allowance for a singlet-like or any of the BSM Higgs couplings to the vector bosons. For a quantitative estimate of the allowance, one may find that C H 125 V V ≥ 0.83 [53] at the 3σ CL level which is also quite consistent with Ref: [54]. The future prospects for the measurements of C H 125 V V at the LHC depend upon uncertainty scenarios in particular and have been discussed in detail in Ref: [55,56,57]. Using Ref: [57], one may set a lower bound on C H 125 V V ∼ 0.95 at the 3σ CL level which means that C H BSM V V can accommodate an allowance of ∼ 0.3 (where H BSM refers to any BSM Higgs H 2 or H 3 ). The said limit has been obeyed in our analysis. This, in general, makes VBF processes as somewhat less interesting to probe the properties of any BSM Higgs scalar at the LHC.

Choice of parameters and LHC constraints
We use the code NMSSMTools [58] to compute masses and the couplings for sparticles and Higgses and also branching ratios (BRs) of the Higgs states. In the present context we focus on the small tan β (= 2). This choice effectively minimizes the b-quark loop contributions to Higgs production via ggF processes. For different parameters and soft SUSY breaking terms we make the following choices: • We set squark masses ∼ 2.5 TeV in the present analysis. This choice is motivated by the fact that any SUSY particle has not been yet observed at the LHC. With this choice, the contributions coming from the squarks mediated loop diagrams in the effective Hgg couplings can be neglected 4 .
• Trilinear soft SUSY breaking terms are set at • Masses for all gauginos are set at 2 TeV. Here Higgsino can be the lighter which can in general make Higgsino as the viable dark matter (DM) candidate of the Universe. But in the case of singlino-Higgsino or Higgsino like DM, Direct-detection experiments like LUX [59] can lead to stringent constraints on the parameter space which is particularly severe for a lighter Higgs spectra 5 . Assuming gravitino as the LSP would relax this tension completely. Note that thermally produced gravitinos may form the observed dark matter [62] without conflicting with the particle dark matter search.
• We impose a requirement that Landau pole singularities of the running Yukawa couplings would not occur below 10 TeV [13], which may otherwise lead to conflict with precision electroweak tests of the SM. This typically yields a requirement on the coupling parameter λ (one needs to satisfy λ ≤ 2). The NMSSM with such a large λ (precisely 0.7 < λ < 2) is commonly called the λ−SUSY model [63,64,65,66]. In this case, larger tree level contributions helps one to accommodate a SM-like Higgs boson with mass of 125 GeV more easily. One may also find that the sensitivity of the weak scale to the stop mass scale would be reduced here [63,64,65]. In this analysis we scan λ ∈ 0 − 2 and κ ∈ 0 − 1 and always check the absence of unphysical global minima of the Higgs potential as done in NMSSMTools [13].
• We use the following limit for SUSY Higgs mass m H 1 which is the lightest Higgs particle of the NMSSM spectrum with standard model like couplings.
A 3 GeV theoretical uncertainty around m H 1 125 GeV has been considered due to uncertainties in the computation of loop corrections up to three loops, top quark mass, renormalization scheme and scale dependence etc. [67].
• On the flavor physics side, the constraints from B-physics namely B → X s + γ, B s → µ + µ − have been implemented, as done in NMSSMTools.
• Finally we use 173.1 GeV for the top quark pole mass.
We would now try to see the feasibility of the NMSSM parameter space where the ggF can be observed as a sub-dominant process for the lightest non SM-like Higgs boson production at the LHC. Our aim would be to find the parameter space consistent with the experimental and the theoretical constraints and where H 2 couplings with the up-type quarks are vanishingly small. We first show our results in the λ − κ plane (see Fig.1a) where the input parameters are kept fixed at tanβ = 2, A κ = −700 GeV, A λ = 250 GeV and µ ef f = 600 GeV. While varying the input parameters, constraints coming from different phenomenological observables (like one of the Higgs boson would be SM-like, collider and the flavor physics constraints), as implemented in NMSSMTools have been checked. The gray points in the parameter space in the Fig.1a are in general not compatible by the said constraints, thus fail to qualify as the valid points. The red region has been drawn to represent the allowed spectra where we mainly neglect the landau pole or the dark matter constraint. This region satisfies the aforesaid flavor physics constraints and is consistent with a SM-like Higgs boson mass with 122.1 < m H 1 < 128.1 GeV. In the figure, the green region represents the parameter space where one may obtain |S 2,2 | << |S 2,1 | for the lightest non SM-like Higgs boson H 2 . This region, in general, represents the parameter space where H 2 couplings with the up-type fermions are tiny, in particular, one finds that |S 2,2 | ≤ 0.001. Clearly, the green region mostly falls in the gray part while a small part is consistent with the phenomenological constraints (this part resides on the top of red regions). Such a small value of |S 2,2 | depends on the input parameters that specify the Higgs masses and mixings (see Sec:2). From the figure, one may see that the green strip where one can get a tiny |S 2,2 | corresponds to λ ∼ κ. Though the exact dependence can only be understood numerically, intuitively the linear dependence in the λ−κ plane can be found if one sets M 2 S,32 0 in the tree level mass matrix for the CP-even Higgs scalars in Eq:5. The exact dependence depends on the other input parameters e.g., A λ , tan β and µ ef f . Besides, we also plot contours for H 2 mass. Now qualitatively, a tiny |S 2,2 | refers to suppression in the reduced H 2 gg coupling, noted by C 2 H 2 gg . As already discussed, this is the specific region of our interest where the VBF mode may turn out to be the main production channel. For an estimate of the relative dominance of the VBF mode over the ggF channel in the computation of the production rate of H 2 scalar we also present the ratio of the their reduced couplings respectively, i.e., Notably, one may also view the ratio as ) and σ R (H ggF 2 ) are the reduced production cross-sections in the VBF and the ggF mode respectively 6 . We find that the ratio becomes maximum along the green strip and can be as large as ∼ 1000. Interestingly, this also includes the region which has been satisfying all phenomenological constraints.  fig.(a), but now the parameter space has been studied in the κ − A κ plane for a fixed λ = 0.667. Now, it may also be interesting to know the allowance in A κ to obtain |S 2,2 | << |S 2,1 | or a tiny |S 2,2 | for a given κ (see Fig.1b). To understand the same, we set the input parameters at tan β = 2, µ ef f = 600 GeV, A λ = 250 GeV which are the same as in Fig.1a. With this choice of input parameters, we set λ ∼ 0.66 and vary A κ (-1TeV to 1TeV) to study the parametric dependence in the A κ − κ plane for obtaining |S 2,2 | << |S 2,1 |. Again, the valid parameter space is shown by the red colored regions where one may easily accommodate the SM-like Higgs boson mass constraint. Similarly, green region corresponds to the parameter space where one gets |S 2,2 | << |S 2,1 |. Varying 6 The absolute value of the cross-section may be obtained if we multiply it with the respective production crosssection of the H2 scalar assuming it SM-like. Typically, for a heavy Higgs having SM-like properties with mass mH 2 ∼ O(100) GeV, the ggF dominates the production process over the VBF by an order of magnitude i.e., O(10) [68]. Thus a large value of may clearly define the VBF dominance over the ggF in producing the H2 scalar in the NMSSM.
A κ would effect the singlet component of Higgs mass matrix which in turn can influence mass and couplings of the Higgs bosons. Indeed observing the green points, one finds that |A κ | may need to be enhanced with κ to satisfy |S 2,2 | << |S 2,1 |. Here one may also find that |S 2,2 | ≤ 0.001. Similarly, a large value of may be obtained along the green region. Additionally, we also present the contours of H 2 mass for representative values (=500, and 800). In summary, a tiny |S 2,2 | or a large value of can be observed in the λ − κ or in the A κ − κ plane (see Fig:1) which may correspond to a H 2 scalar with mass m H 2 ∼ 600 GeV.
We would also analyse if such a large value of the can be accessed over a wide range of H 2 mass (see Fig.2a). Here we scan the parameters: 0.1 ≤ λ ≤ 2, 0 ≤ κ ≤ 1, 100 ≤ µ ef f ≤ 2000 (GeV), −1000 ≤ A κ ≤ 1000 (GeV). The other input parameters are tanβ = 2, A λ = 250 GeV. All the points satisfy the Higgs mass constraint and other phenomenological constraints as discussed above. The region includes both large and small values of |S 2,2 |. In the parameter space where |S 2,2 | << |S 2,1 | or where one may obtain a tiny |S 2,2 | a large relative enhancement of the can be found. In practice, as discussed already that |S 2,2 | ∼ 0 refers to a tiny value of C 2 H 2 gg . Keeping this in mind, it may be imperative to estimate the magnitude of the reduced couplings in the VBF and ggF mode which we show in Fig:2b. Here we impose a constraint that As can be seen from Fig:2b for the same set of input parameters with the mass of the Higgs of our interest in Fig:3. This may give an idea whether H 2 scalar of a given mass can be probed at the LHC.
In order to demonstrate the observable effects where VBF contribution to H 2 production can be promising, which can potentially be seen at the high luminosity run or at the upgraded LHC, we present a few benchmark points (BPs) in Table:1. All the points, shown here, are consistent with constraints related to the observables in the Higgs sectors as done in HiggsBounds-4.3.1 [69]. We set tan β = 2 for all these points. In addition, for each BP we show masses for Higgs scalars, compositions of scalars (only for the SM-like and the lighter BSM Higgs), reduced couplings of H 2 with the electro-weak gauge bosons, with the gluons and with the fermions. In all cases, H 2 is dominantly singlet like with vanishingly small H u components. Thus, without any surprise, the reduced coupling of H 2 scalar to the top quark C 2 H 2 tt is extremely suppressed. Similarly, the reduced couplings with the down type fermions are smaller for tan β = 2. As a natural consequence, one would expect larger C 2 H 2 W W/H 2 ZZ in comparison to C 2 H 2 gg which may, in turn, enhance the relative contribution of the VBF processes in H 2 production compared to the ggF channel. Indeed one can see it from Table:2, where for the said BPs (as shown in Table:1), we present effective production cross-section for the lighter BSM Higgs boson H 2 through the ggF, the VBF and in association with the b quark channels. Here σ SM refers to the production cross-section corresponding to a Higgs states with mass m H 2 having SM characteristics [68] at the center of mass energy √ s = 13 TeV. For the ggF, σ SM is obtained assuming NNLO+NNLL QCD accuracy while for the VBF processes only NNLO QCD accuracy is used for the same. Here for the calculation of effective production cross-section we rescale the SM cross-sections with the effective squared couplings.   As H 2 is dominated by singlet component, it can hardly decay into τ leptons which is otherwise a standard search channel for MSSM heavy Higgs searches [49]. Here, we see that H 2 dominantly decays into the gauge bosons. A considerable amount of efforts have already been put forward to discover BSM Higgs scalars through its decay into the vector bosons modes at the LHC [70,71,72,73,74] and non observation of any new physics would eventually lead to the upper bounds (UB) on (σ × Br) BSM Higgs at 95% CL. Using Ref: [72,73] which assumes √ s = 13 TeV and integrated luminosity L = 36.1 f b −1 , we study the status of these BPs focusing on the promising W W and ZZ final states and the results have been displayed in Table:3. Here the production of H 2 in association with the b quark has not been considered where one may find upper limits on (σ × Br) BSM Higgs assuming H 2 decaying into τ τ [75] and bb [76] at √ s = 13 TeV run of the LHC. In Table:3, different upper bounds (σ × Br) exp−U B BSM Higgs for the ggF and the VBF modes have also been shown respectively. Here one assumes W W decaying leptonically to eνµν. Similarly, for ZZ final states, 4l and 2l2ν channels are combined. In the last column the necessary improvements for the future runs of the LHC to probe our BPs in terms of the predicted (σ × Br) N M SSM have been calculated. Clearly, the ZZ final state has been appeared as the best channel for observing the H 2 scalar at the future runs of the LHC. In fact the relatively light Higgses with m H 2 (=250,500 GeV) in BP-A and in BP-B can be probed via ZZ final state with roughly an order of magnitude enhancement in the experimental sensitivity while somewhat less enhancement would be necessary to probe m H 2 = 700 GeV in the case of BP-C. Additionally, in the latter case, one may also observe that the ggF contribution to H 2 production becomes comparable with the VBF production channel though the latter still  Table 1: Input parameters for the three benchmark points along with masses for the SM-like and the lighter BSM Higgs scalars are shown. Different components of H 2 scalar along with the reduced couplings for Higgs production via the ggF , the V BF and associated production with the b quark have been presented. In addition, we also show the relevant branching ratios.   Table:1 through the ggF, the VBF and in association with the b quark have been presented.
Here σ SM refers to the production cross-section corresponding to a Higgs states with mass m H 2 having SM characteristics [68] at center of mass energy √ s = 13 TeV. For the ggF, σ SM is obtained assuming NNLO+NNLL QCD accuracy while for the VBF processes only NNLO QCD accuracy is used for the same.
dominates. This is because for this parameter point C 2 H 2 gg is somewhat larger compared to the previous BPs which can again be attributed to somewhat larger S 2,2 . Thus BP-C serves as a good example where the ggF production mode can compete with the VBF production channel. Here one may observe that a roughly equal amount of enhancement in sensitivity would be required for both production channels to discover the Higgs boson at the future runs of the LHC. In the next section, we would study the prospects of observing this Higgs scalar at the high luminosity run of the LHC at √ s = 14 TeV and at a future pp collider with 27 TeV center of mass energy considering the fact that the rate of production of the said scalar is dominated by the VBF channel.
So far we discuss about a BSM Higgs scalar which is primarily singlet dominated with nonnegligible H d components in it. It is the lightest of the two BSM Higgs scalars in the present context. We will see that the problem of discovering H 2 persists as a motivation for experimental searches. In this context, it may be interesting to know the status of other Higgs states. (i) Heavier CP-even Higgs state: The heavier one is mostly doublet dominated (approximately MSSM-like) whose properties can be seen from the Table:4. The effective production cross-section is completely dominated by the ggH channel through t quark loop and has been calculated at NNLO with SusHi [77]. In the shown BPs, lightest H 3 can be obtained for BP-A, σ ggH 3 N M SSM is approximately 1pb. Afterwards, it also dominantly decays into tt final states and non-negligibly to H 1 H 2 or A 1 Z (when kinematically allowed). In case of tt final states, using the simplified analysis as carried out in ref: [49], it may be possible to discover the heavy MSSM-like Higgs H 3 at √ s = 14 TeV LHC run with 300 f b −1 data. Additionally, heavier Higgs to lighter Higgs decays could be an important search channel in this case. In fact Higgs→ Higgs decays can be an important search channel for the lighter BSM Higgs scalar as well (see Table:1). These can in turn lead to different possibilities in the final states with additional jets in case of the Higgs productions through the VBF process. For the decay modes H 1 H 2 and H 1 H 1 , stringent constraints can be obtained from the resonant SM Higgs pair productions where Higgs can subsequently decay into into bbbb, bbτ τ , bbγγ and bblνlν [78,79] [72,73]. In the last column the necessary improvements for the future runs of LHC to probe the BPs in terms of the predicted (σ × Br) N M SSM through the ggF and the VBF production modes have been presented.
well below the present bounds as obtained by √ s = 13 TeV LHC run with L = 35.9 f b −1 (see e.g., the summary plot in ref: [80]). (ii) CP-odd Higgs state: The lightest CP-odd Higgs is dominantly singlet-like while the heavier one is doublet-like. Their masses have been shown in Table:4. The production of the CP-odd states would be dominated by the gluon-gluon fusion process. Clearly, the coupling of relevance would be A i gg. For the lightest CP-odd Higgs scalar, we find that A 1 gg becomes maximum for BP-A and the reduced coupling becomes C 2 A 1 gg ∼ 0.1. Similarly, for the same BP, one may read C 2 A 2 gg ∼ 0.36. In case of other BPs, one gets somewhat smaller values for C 2 A 2 gg . Thus BP-A may potentially be seen as best example to search for the CP-odd states at the LHC in the present context. After production, for BP-A and BP-B, the CP-odd scalars A 1 and A 2 dominantly decay in to tt with branching ratio more than 90%. In fact the MSSM-like CP-odd Higgs state A 2 may be probed at √ s = 14 TeV LHC run with 300 f b −1 data just like the MSSM-like heavy CP even Higgs state H 3 [49]. In case of BP-C, A 2 can decay in to ZH 2 with 20% BR and A 1 H 1 with 10% BR. However, with m A 2 ∼ 1 TeV, the production rate would be suppressed to search it through the cascades of lighter Higgs scalars.  be the most promising channel) where Z can decay into leptons, Z → l + l − (where l = e, µ). One can characterize the VBF induced processes, pp → jjH 2 (j refers to light jets) by the presence of two energetic jets with a large rapidity gap. For a detailed collider simulation, we generate parton level events using MadGraph5 with PDF choice NN23lo1 [81,82,83,84] and subsequently passed the events into Pythia (v8) [85] for hadronization taking into account the initial state radiation/final state radiation (ISR/FSR) and multiple interactions. We use MadGraph5 for implementing decays of the lightest BSM Higgs scalar. The jets, leptons, photons have been reconstructed using fast detector simulator Delphes-v3.3.3 [86,87,88]. The jets and leptons are reconstructed by anti-k t algorithm [89] implemented in the Fastjet [90,89,91] with a cone of ∆R = 0.4 and minimum transverse momentum of 20 GeV. We generate our new physics signal events associated with H 2 production at the leading order (LO) using the UFO file based on Ref: [92] which considers VBF at NLO-QCD accuracy. Since signal does not include much missing transverse energy, the dominant background appears to be pp → ZZ + 2j (where j stands for light jets) where Z decays leptonically [73]. There may be other sub-dominant background processes like pp → γγ + 2j, pp → γZ + 2j or pp → W Z + 2j . Then we use the following selection cuts [73] to optimize the significance (S = S √ B+S ).

LHC phenomenology at
• We demand the presence of two same-flavour, opposite-sign leptons (electron or muon) with p e T or p µ T > 15 GeV. This criterion is somewhat different compared to what has been assumed in Ref: [73] where the highest-p T lepton in the quadruplet must satisfy p T > 20 GeV, and the second (third) lepton must satisfy p T > 15 GeV (10 GeV). Similarly, the pseudorapidity range of electron or muon has been considered to lie |η| < 2.47 (2.7).
• We require the two leading forward jets to satisfy p T ≥ 30 GeV and |η| < 4.5 . Moreover, these two jets should lie in opposite hemispheres with η j 1 × η j 2 ≤ 0.
• While computing the statistical significance, we impose additional cuts : (i) For BP-A (m H 2 = 250 GeV), the four lepton invariant mass (m 4l ) should lie between 240 − 260 GeV and (ii) for BP-C (m H 2 = 700 GeV), m 4l should be above 600 GeV.
For background simulation, we generate pp → ZZ in association with 0-2jets at the matrix element level using MadGraph5. Here we also implement the MLM matching technique to avoid double counting. The production may arise from quark-antiquark annihilation, gluon-initiated production associated with jets and from EW vector-boson scattering. The last process plays a more important role in the VBF-enriched category which makes signal discrimination very difficult from the backgrounds. Considering the fact that effective production cross-section of the signal process pp → H 2 jj → ZZjj → l + l − l + l − jj is not quite large for singlet like Higgs, we find that probing lighter BSM Higgs with mass m H 2 = 250 GeV (BP-A) is more difficult even at the upgraded LHC. In order to have better significance we devise a small window of 20 GeV around the Higgs boson mass in the first case. For quantitative estimates, we compute the number of signal and background events for the lighter and heavier BSM Higgs scalar for the BP-A and BP-C. As shown in Table:3, for these BPs, the necessary improvements which would be required to probe them at the future runs of LHC would be smaller compared to the BP-B. This is related to smaller C 2 H 2 ZZ which is in turn related to smaller S 2,1 in case of BP-B. The production cross-section for the H 2 scalar has been evaluated [92] at the NLO-QCD level which resembles quite well with [93] while for the background productions pp → ZZjj → l + l − l + l − jj (l = e, µ), we use LO estimates multiplied by K-factor. At the leading order our computation produces σ(pp → ZZjj) = 17.4 × 10 3 fb at the √ s = 14 TeV and 44.3 × 10 3 fb for the √ s = 27 TeV runs of the LHC. In Ref: [94,95], the total ZZ cross-section predictions is computed to be 18 pb with MATRIX [96] at NNLO in QCD for √ s = 14 TeV which is quite consistent with our result. We use a conservative choice for K-factor (∼ 1.0) which is consistent with the limit [97] if one considers the scale uncertainty. We note that for the 27 TeV run of the LHC, we assume the detector acceptance and efficiency are the same as the 14 TeV run. In Table 5 and 6, we show the effective production cross-section for signal (σ ef f = σ(pp → H 2 jj → ZZjj → 4ljj) × ef f iciency) and for the backgrounds (σ ef f B = σ(pp → ZZjj → 4ljj) × ef f iciency) and the number of typical signal and background events obtained for 3 and 15 ab −1 integrated luminosity (L) for √ s = 14 and 27 TeV runs of the LHC respectively. One may see that the effective production cross-section for signal events σ ef f is not quite large, though the situation somewhat improves for higher cms energy. This is in general responsible for having a small number of signal events at the LHC. Hence, though the Higgs production through the VBF mode appears to be an interesting possibility, but it really lacks to provide with the desired sensitivity at the high-luminosity run or even at the upgraded LHC. This is especially true when the lightest non SM-like Higgs scalar is light. Situation improves for a heavier Higgs boson, thanks to the rapid fall of the four lepton backgrounds with its invariant mass.
√   Table 6: Expected number of Signal and background events that qualify the selection criteria, mentioned in the text for the BPs A and C. We take L = 3&15 ab −1 integrated luminosity for the √ s = 14 and 27 TeV runs of the LHC respectively.

Conclusion
In this work, we show that the rate of non SM-like Higgs productions via the gluon-gluon fusion process could be insignificant in some of parts of the NMSSM parameter space. This can be attributed to the smallness of up-type Higgs component in the Higgs scalar and also to the small values of tan β. In this case, the vector-boson fusion may play the leading role in the BSM Higgs searches if the new state is not completely dominated by singlet component. We present three benchmark points where the mass of the non SM-like Higgs scalar may vary from 250 GeV to 700 GeV. Then we study the role of existing LHC constraints on the said BPs in the light of the 13 TeV data. We also discuss the prospects of observing such a state at the high-luminosity L = 3 ab −1 run of the LHC at cms energy √ s = 14 TeV and at the proposed √ s = 27 TeV run of the LHC with L = 15 ab −1 . The dominant background comes from pp → ZZ + 2jets in which VBF may also play a significant role which makes it somewhat difficult to distinguish the backgrounds from the signal, especially for lighter BSM Higgs scalar in BP-A (m H 2 = 250 GeV). For heavier BSM Higgs, we can use a stronger cut on the 4 lepton invariant mass, m 4l and in turn this can yield better significance. Finally, in view of this work, it may be interesting to study the prospects of such a BSM Higgs scalar at the Linear Colliders.