A Singlet Extension of the MSSM with a Dark Matter Portal

The minimal extension of the MSSM (NMSSM) has been widely studied in the search for a natural solution to the $\mu$ problem. In this work, we consider a variation of the NMSSM where an extra singlet is added and a Peccei-Quinn symmetry is imposed. We study its neutralino sector and compute the annihilation cross section of the lightest neutralino. We use existent cosmological and collider data to constrain the parameter space and consider the lightest neutralino, which is very light, as a dark matter candidate.


Introduction
The long standing problem of evading the LEP bound [1] on the mass of the lightest (Standard Model (SM)-like) Higgs boson within supersymmetric extensions of the SM (SUSY) has to now be reinterpreted due to the recent discovery of what appears to be a SM-like Higgs boson with a mass of 125 GeV [2,3]. However, the underlying nature of this problem remains the same. In the Minimal SUSY SM (MSSM) the lightest Higgs must lie below the Z mass at tree level. This mass can be pushed up through radiative corrections arising from the third family of quarks and squarks. However, the discovery of a 125 GeV SM-like Higgs boson has placed the MSSM into a region of parameter space where the hierarchy problem, which SUSY is expected to solve, is reintroduced. That is, reaching a value of 125 GeV requires large stop masses, as heavy as 10 TeV, or a tuned value of the stop mixing parameter at the electroweak scale [4,5]. This version of the original hierarchy problem is well known as the "little hierarchy problem" and it is quite generic within the MSSM (see [6][7][8] and references therein).
One popular route that is taken to alleviate this problem is to extend the Higgs sector of the MSSM. This has the effect of generating new quartic terms in the scalar potential [9,10]. These new quartic terms push up the mass of the SM-like Higgs boson at tree level, removing the need for large radiative corrections. In particular, this can be achieved by extending the MSSM sector with SM gauge singlets. The minimal extension of the MSSM, refered to as the next-to-minimal SUSY Standard Model or NMSSM, incorporates a single gauge singlet and was introduced primarily to address the µ-problem of the MSSM (For reviews, see [11][12][13].) This model has had its fair share of success, but it is not clear how one can naturally generate the Higgs mass on the order of 125 GeV without introducing some degree of fine tuning, arguably as large as in the MSSM [14]. Furthermore, its minimal incarnation may introduce tension between the way the hierarchies are stabilized and the generation of domain walls [15]. However, it has been shown that a stable NMSSM without domain walls is possible if certain discrete R-symmetries are imposed [16,17]. Nonetheless, the discovery of a SM-like Higgs boson has led to many studies on the phenomenology of a 125 GeV SM-like Higgs boson within the NMSSM [18][19][20][21][22][23][24][25][26]. In addition, it has been noted that a generalized version of the NMSSM that follows from underlying discrete Rsymmetries can reduce the amount of fine tuning in the scalar sector [27][28][29]. Within this class of models, additional operators are generated when SUSY and the R-symmetry are broken in the hidden sector and the effects are mediated to the observable sector through Planck-suppressed operators. Alternatively, one may generalize the MSSM by introducing effective dimension four and five operators [30][31][32]. These operators can reduce the amount of fine tuning in the scalar sector and be sensitive to new degrees of freedom at the TeV scale.
One particular generalization of the NMSSM, the S-MSSM [33,34], has been implemented to fully address the little hierarchy problem of the MSSM. Within this framework one gives up any attempt at addressing the origin of the µ-term and the absence of a Z 3 symmetry and incorporates a supersymmetric mass for the SM singlet field which is used to stabilize the singlet's vacuum expectation value (vev). A version of this model has been successfully embedded into a model where SUSY breaking is mediated by gauge interactions [35].
It may be possible to argue that a successful natural solution to the little hierarchy problem and the µ-problem exists within one unified model. In the scenario described in [34], it was shown that one can successfully eliminate the µ-term as a phenomenological parameter as long as a small supersymmetric mass for the singlet is incorporated. The model is described by the following superpotential: (1.1) Of course, the above superpotential reintroduces a µ-problem, a µ S -problem. The model we propose in this work goes one step further as it replaces the µ S -parameter by a second SM gauge singlet superfield, µ S → ρN , with ρ a dimensionless parameter. Within this framework the gauge singletN has no direct couplings to the MSSM fields. In this way the S-MSSM gauge singlet,Ŝ, serves as a portal for the singletN . This is possible due to the existence of a (PQ) symmetry 1 . Additionally, we explore the existence of a dark matter candidate within the model presented in this work. Currently, there is plenty of evidence that points towards the existence of dark matter (DM) in our universe [37,38], providing strong evidence for physics beyond the Standard Model. Recent results from Planck [39] suggest a cold dark matter component with a density of Ωh 2 = 0.1199 ± 0.0027. Furthermore, there exist positive signals in direct detection experiments [40,41] that point to a light dark matter candidate with mass at around 10 GeV. Commonly, supersymmetric models contain a light degree of freedom that is cosmologically stable. In fact, the lightest neutralino in our model has a mass between 1 − 15 GeV and annihilates into SM particles, mainly leptons and light quarks, through the exchange of light CP-even and -odd Higgs bosons and neutral gauge bosons. We take a look at the regions of the parameter space that allow a relic density that agrees with observations and consistent with collider searches and constraints arising from Big Bang nucleosynthesis.
The structure of this paper is as follows: In Section 2 we introduce the model and look at the structure of Electroweak Symmetry Breaking (EWSB). In Section 3 we review the constraints arising from colliders that limit our parameter space while in Section 4 we show the main annihilation channels contributing to the density of dark matter in the universe. In Section 5, we offer concluding remarks on the possibility of a light neutralino in singlet extensions of the SM.

Electroweak Symmetry Breaking and Scalar Higgs Sector
In this work, we modify the S-MSSM [33,34] by replacing the supersymmetric mass term for the SM gauge singletŜ by an additional SM gauge singlet superfieldN . This new superfield does not couple directly to the fields in the MSSM, but only through the mixing induced by a superpotential coupling betweenN andŜ. Additionally, we impose a Peccei-Quinn (PQ) symmetry where both MSSM Higgs doublets, H u and H d , have charge 1 and the singlets S and N are given charges −2 and 4 respectively. Furthermore, under this symmetry, quarks and leptons have P Q charges of −1/2. Using this framework, the superpotential we consider, defined at some high energy scale, is given by: where the κ term is introduced to give mass to a Nambu-Goldstone boson that arises from the spontaneous breaking of the PQ symmetry due to EWSB. Additionally, SUSY breaking generates the following contributions to the scalar potential . In this analysis, we are interested in the limit where S and N interact weakly, that is ρ 1. Furthermore, we consider only small values for κ such that the PQ symmetry is only slightly broken. Although there is no symmetry that forbids the PQ symmetry breaking operators α 1N 2Ŝ and (α 2 /3)N 3 and λ NNĤu ·Ĥ d , we have set them to zero at the messenger scale. This is a scale-dependent assumption. Nevertheless, this hypothesis is quite stable under renormalization group effects. Indeed, assuming small values for α 1 , α 2 and λ N at the messenger scale, these couplings remain small and well below ρ and κ at the electroweak scale. The running between a messenger scale given by M mess = 10 12 GeV and the weak scale is shown in Figure 1, where we have used the one-loop renormalization group equations for the dimensionless couplings given in Appendix A. Again, it is important to emphasize that a particular high energy choice for the superpotential in Equation (2.1) was made. The structure may be achieved with additional dynamics above the messenger scale.  Furthermore, we note that once the PQ symmetry is broken by the κŜ 3 operator in Equation (2.1), additional contributions to V sof t breaking this symmetry will be generated at one and two loops. These operators are further suppressed by powers of κ, ρ and λ. In addition, our model has an exact Z 3 symmetry as in the NMSSM. This discrete symmetry can lead to the generation of domain walls. However, our framework (with κ → 0) has a U (1) R symmetry, as the one discussed in [17], which has a Z 5 subgroup that can induce a tadpole term in the scalar potential large enough to avoid a domain wall problem without destabilizing the electroweak hierarchy. Nevertheless, our model contains two singlets and a more detailed analysis of the R-symmetries and their discrete subgroups is imperative for an in depth study of Equation (2.1). Suffice it to say that one may completely avoid an explicit PQ breaking term in the superpotential of Equation (2.1) and find an appropriate discrete R-symmetry that can induce tadpole terms capable of stabilizing the hierarchy and avoid the cosmological domain wall problem.
The superpotential in Equation (2.1) together with the soft-breaking terms in Equation (2.2) give rise to the following scalar potential for the neutral components of the two electroweak Higgs doublets and singlet fields: where g and g are the gauge couplings of the SU (2) W and U (1) Y gauge groups respectively. Minimizing the scalar potential with respect to H 0 u , H 0 d , S and N leads to the following four constraints: GeV and tan β = v u /v d . The vacuum expectation values for the two singlets have been obtained in the limits where both κ and ρ are much smaller than one. Equations (2.4) and (2.5) are analogous to the MSSM minimization conditions with an effective µ-parameter given by µ ef f = λv S and an effective B µ term given by In the absence of explicit CP-violating phases in the Higgs sector, the physical spectrum of the model includes a single charged Higgs boson (H ± ), four neutral scalars that we label h N , h S , h 0 , H and three neutral pseudoscalars (A N , A S , A). The states labeled with a subscript will turn out to have a large singlet component. For the state most resembling the usual pseudoscalar Higgs of the MSSM, the mass is given by In the ρ → 0 limit, ρ · v N is largely suppressed and the effective supersymmetric mass for the singlet S, µ S,ef f , is small. This is interesting since the spectrum of H u , H d and S mimics the one studied in [34]. In the analysis, in the limit where µ 2 S , m 2 S λ 2 v 2 , mixing of the singlet into the light MSSM-like scalar vanishes, yet receives an NMSSM-like enhancement [12]: Therefore, in the ρ, κ → 0 limit, we expect a similar state since we can work within the regime where µ 2 S,ef f , m 2 S λ 2 v 2 . However, this is not the case for finite κ, where a "push-up" effect is expected to increase the mass of the SM-like scalar at tree level [24]. This effect is due to the fact that the singlet h S is lighter than the SM-like state and mixing between the two increases the mass of the latter. This phenomenon is evident if we write the upper 3 × 3 mass matrix of the CP-even scalar sector in the basis (H d cos β + H u sin β, H u cos βH d , S) ≡ (h 0 , H, h S ) as in [24,34], but for finite κ and m 2 S and neglect the small corrections proportional to ρ. In this basis, the mass matrix has the following form: where R = 1 v λ(κv S + 1 2 A λ ) sin 2β . Using Equation (2.6) for v S , one can see that the (1, 3) element of the above matrix vanishes for κ, m 2 S → 0. For finite κ and m 2 S , the SM-like Higgs mass at tree-level is given by where δm 2 h 0 ,mix is a function of κ, ρ and m 2 S and parametrizes the contribution from the h 0 − h S mixing. The mixing between the SM-like Higgs, h 0 , and h S will also have an effect on the couplings of the former to SM matter fields. In particular, it will suppress the coupling of h 0 to gauge bosons while it will enhance the coupling of S. Therefore, for finite κ and m 2 S , hiding the light S state from Higgs searches carried out by LEP becomes a strong constraint on the model's parameter space [1,[42][43][44].
In the limit where κ, ρ 1, the mass matrix for the CP-even scalars in the basis (H uR , H dR , S R , N R ) is given by the following terms: In the limit where ρ 0, the mixing between the singlet N and the other three scalars, M 2 i,4 , is largely suppressed for not too large values of the tri-linear coupling A ρ and the vev of the singlet S, which in this model can be adjusted through A λ . The mass of N will then depend mostly on the soft SUSY breaking mass parameter m 2 N . In our framework, we choose then to work in the following limit µ 2 S,ef f , m 2 S λ 2 v 2 while keeping m 2 N as a less constrained free parameter. In this limit, the masses of the CP-even scalars are given by: In order to keep the above equations as clear and simple as possible, we have not incorporated corrections proportional to ρ. However, the calculation of the masses is done exactly in our numerical routines.
As mentioned earlier in this section, in order to generate an spectrum similar to the one studied in [34], it is important to work in the limit where m 2 S λ 2 v 2 . This condition is somewhat unnatural since there exist contributions to the one-loop renormalization group equation for m 2 S that are proportional to A 2 λ [12], that in our framework is large in order to decouple the MSSM-like pseudoscalar, with mass given in Equation (2.8), from the spectrum. One may alleviate this by embedding the model into a SUSY breaking mediation mechanism where the scale of SUSY breaking is not very high.

Neutralino Sector
The neutral gauginos of this model (B,W 0 ) mix with the two neutral higssinos and the two singlinos to form the neutralino mass eigenstates due to the electroweak symmetry breaking and the Yukawa couplings. Using the basis the mass terms in the Lagrangian for the neutralino sector are given by (2.17) The corresponding mass eigenstates are given by where the unitary mixing matrix, N i j , diagonalizes Equation (2.17), and where the eigenmasses have been labeled in ascending order. The leading contributions to the masses of the two lightest neutralinos, χ 0 1,2 , are given by Within our framework, χ 0 1 is mostly singlino and couples weakly to the other particles in the spectrum. However, as it will be shown in Section 3, it could have a small but significant bino and higgsino components. This will play an important role in the cosmological evolution of the energy density of this stable particle, since it makes the self-annihilation effective enough to avoid overabundance of the relics. On the other hand, the mass of the next-to-lightest neutralino within the parameter space considered in Section 3, will be at least twice as massive as the lightest supersymmetric particle (LSP); and co-annihilations between χ 0 1 and the heavier neutralino will not relevant for the calculation of its relic abundance. Instead, the relic density will be determined by the annihilation cross-section of the LSP, as explained in Section 4. The couplings of the lightest neutralino to the CP-odd and CP-even Higss scalars, that will be used in the relic abundance calculation, are given by where i = 1, 2, 3, j = 1, 2, 3, 4, N ≡ N −1 and R A ,R H diagonalize the CP-odd and CP-even mass matrices defined in Equations (2.13) and (2.15).

LEP Constraints
One important constraint on the parameter space is due to the LEP bound on the chargino mass, m χ + > 104 GeV. This bound translates into a bound on µ ef f given by |µ ef f | > 104 GeV. Using Equation (2.6), this can be re-casted into a bound on A λ given by A λ > GeV. For tan β = 2 and m 2 S λ 2 v 2 , A λ is bounded from below by 260 GeV. However, constraints on the singlet-like scalar fields yield a finite value for m 2 S and the lower bound for A λ lies slightly above 260 GeV.
Constraints on light scalars also limit the parameter space of this model. In particular, searches by LEP [1,[42][43][44] place strong upper bounds on the two main scalar production mechanisms: e + e − → HZ and e + e − → HA, where H and A denote any of the CP-even or -odd scalars respectively. In the HZ channel these constraints assume that each scalar decays to bb or τ + τ − with a branching fraction equal to one. In general, these bounds will soften since the scalars in our framework can also decay to lighter scalars with a significant branching fraction.
The possibility of Higgs cascade decays has also been searched for at LEP [43,44]. They place strong bounds on two channels: (1) Associated Higgs production with a Z, e + e − → ZH i , H i → A j A j and (2) Scalar-pseudoscalar pair production, e + e − → H i A j . In (1), the analysis assumes a BR(H j , A j → bb) = 1 and BR(H i → H j H j , A j A j ) = 1. In (2), five different final states were analyzed: In our analysis, we calculate the normalized cross section for scalar-pseudoscalar pair production which is given by whereλ is a kinematic factor given bȳ and s is the center of mass energy squared. We multiply the normalized cross-section, σ H i A i /σ SM HZ , by the appropriate branching fractions in the decay chain. Furthermore, we implement the constraint found in the channel e + e − → ZH that is independent of the H decay mode [42].

Meson Decays
A pseudoscalar, with a mass in the range between 1 and 40 GeV, has a coupling to fermions that is highly constrained by meson decays and collider data. The couplings can be extracted from the following Lagrangian: where denote the couplings of the pseudoscalar mass eigenstates, A i , to up-type and down-type quarks respectively and R A ij is the unitary matrix that diagonalizes the CP-odd mass matrix introduced in Equation (2.15). For masses below the upsilon threshold of ∼ 9.46 GeV, an analysis by [45] found that Υ → γA i imposes that C A i dd < 0.5 for tan β ∼ 1. Above this mass threshold the same analysis found the strongest constraint on the pseudoscalar mass coming from the process e + e − → bbA i → bbbb measured by DELPHI [44], setting the following limit, C A i dd < O (10). Additional constraints on light pseudoscalars arise form rare B and K decays such as: B → K + invisible, K → π + invisible, B → Ke + e − , K → πe + e − and K → π + + X as well as the muon g − 2 . The analysis in [46] on an NMSSM light pseudoscalar concludes that pseudoscalar masses of m A i < 2m µ are excluded unless the coupling C A i dd lies below 10 −4 . Within our framework, the parameter space consistent with bounds on light scalars and supersymmetric particles yields pseudoscalar masses above the 2m µ threshold.

LHC Constraints
The discovery of a SM-like Higgs boson with mass around 126 GeV provides a new set of constraints that must be addressed in order for the known production cross sections and decay rates to be in agreement with those measured at the LHC [2,3]. The authors in [47] have proposed a method of calculating the total width of a SM-like Higgs boson using data from the LHC and the Tevatron as well as the properties of the SM-like Higgs boson as a benchmark. Furthermore, they provide a method for estimating the branching fraction of the SM-like Higgs boson to dark matter. In addition, the authors in [48] have carried out a global fit to the data and found a total width of a Higgs relative to the SM prediction given by Γ tot /Γ SM tot ∈ [0.5, 2] and an invisible branching fraction of roughly 38% at 95% CL. These results were obtained by varying the Higgs couplings to SM particles independently of each other. More conservative results were obtained by setting the couplings of the Higgs to SM particles to their SM values. They find a Γ tot /Γ SM tot ∈ [1, 1.25] and Br(h 0 → inv) ≤ 19% at 95% CL. In our analysis we calculate the total width of the SM-like Higgs boson, since this gets contributions from light singlet-like scalars and pseudoscalars as well as the light singlet-like neutralinos, and look for deviation from the SM value of Γ SM tot = 4.1 MeV [49]. We require that 0.5 ≤ Γ tot /Γ SM tot ≤ 2 and a Higgs invisible branching fraction of Br(h 0 → inv) 40%.

Γ inv
Z and Neutralino sector The neutralino sector of this model contains two states with a large singlet component, however, the next-to-lightest neutralino may also have a significant amount of Higgsino component. If this neutralino is lighter than m Z /2, Z decays to a pair of next-to-lightest neutralinos could violate bounds on the invisible decay width of the Z. The decay of the Z into a pair of neutralinos is given by: where N n,3 and N n,4 are the down-and up-type Higgsino components of the n th neutralino mass eigenstate respectively as described in Section 2. However, the bound on the invisible Z decay width, ∆Γ inv Z < 2.3 MeV [50], sets more stringent constraints on the next-tolightest neutralino since it has a larger higgsino component. Furthermore, we find that the next-to-lightest neutralino has a mass below 90 GeV and thus the production process e + e − → χ 0 1 χ 0 2 was kinematically accessible at LEP 2. The strongest bound was found by the OPAL collaboration [51]. Since we are considering a lightest neutralino with a mass below 20 GeV, the cross section for the process e + e − → χ 0 1 χ 0 2 is bounded from above by 0.05 pb. To calculate the cross section we follow the analysis in [52] where s = 209. 2 GeV 2 is the center of mass energy at LEP 2 and N ij is the matrix that diagonalizes the neutralino weak eigenstates introduced in Equation (2.19).

Parameter scan considerations
We analyze the parameter space of this model necessary to generate a SM-like Higgs boson with a mass of 126 GeV and light singlet-like states that are consistent with Higgs searches carried out by LEP [1,[42][43][44]. The SM-like Higgs mass is given by: h 0 ,loop parametrizes the leading radiative corrections to the SM-like Higgs mass from third generation of quarks/squarks. This correction is given by wherem t = m t /(1 + 4α s /3π), m t is the pole mass of the top quark, g s is the strong coupling constant, Mt is the geometric mean of the two top squark mass eigenvalues and X t parametrizes the mixing between top squarks: (3.10) In order to maximize this value at tree level, we consider large values of λ. However, we insist that λ remains perturbative at all scales up to the grand unification scale M GU T = 2 × 10 16 GeV. This places an upper bound on λ which peaks for values of tan β between 2 and 3 as in the models described in [12,33,34]. Our analysis is carried out with tan β = 2. Our calculations of the Higgs masses are done using a full one-loop effective potential. Furthermore, in order to maximize the SM-like Higgs mass we use a large MSSM-like pseudoscalar mass, m A . This has the effect of decoupling one of the Higgs doublets from the Higgs sector. In the analysis, we use the four minimization conditions introduced in Equations (2.4)-(2.7) and solve for m 2 Hu , m 2 H d , v S and v N . The remaining parameters of the model are varied as in Table 1. In the scan we fix the mass of the Bino at half the Wino mass, and we set the gluino mass at 3.0 TeV. Based on the constraints introduced in the previous sections we focus on a subset of the parameter scan introduced in Table 1  • κ is scanned in order to minimize the invisible branching fraction contribution to the total width of the SM-like Higgs boson. The coupling of neutralinos to the SM-like Higgs boson is given by: where R H ij is a unitary matrix that diagonalizes the CP-even mass matrix in Equation (2.13) and N is the inverse of N which was introduced in Equation (2.18) and diagonalizes the neutralino sector. Additionally, κ sets the mass of the next-to-lightest neutralino, which in our model sits well above the lightest neutralino mass.
• ρ is scanned in order to generate a lightest neutralino with a mass below ∼ 15 GeV.
We also use a small value of ρ such that the lightest scalar/pseudoscalar in the spectrum have very little mixing with the MSSM-like scalar/pseudoscalar states. The value of κ is scanned between −0.1 and −0.01 and ρ between 0.01 and 0.05. We run our numerical routines considering two values of the Wino mass, [500, 1500] GeV. In Figure 2 we show the different components of the lightest neutralino (left figure) and the next-to-lightest neutralino (right figure) for M 2 = 500 GeV. Both figures are consistent with a Higgs mass of roughly 126 GeV, the invisible Z width, and LEP bounds on charginos. The contribution to the invisible decay width of the SM-like Higgs will arise mainly from the h 0 → χ 0 2 χ 0 2 and h 0 → χ 0 1 χ 0 2 decay channels. This is due to the fact that the next-tolightest neutralino has a large amount of mixing with the Higgsinos and a fine cancellation between the parameters in the model, λ, κ and ρ, is needed. Furthermore, if we compare Figure 2 which corresponds to M 2 = 500 GeV with Figure 3 which corresponds to M 2 = 1.5 TeV, the next-to-lightest neutralino in the former has a larger component along the Wino and Bino directions. Therefore, one will expect that for M 2 = 500 GeV, the values of κ and ρ are more fine tuned for the model to satisfy the constraints from the invisible width of the Higgs as well all other LEP constraints introduced in the previous sections. In fact, this can be seen in Figures 4 and 5 which correspond to contours of the lightest pseudoscalar and lightest neutralino masses in the ρ − κ plane. Both figures were obtained with M 2 = 1.5 TeV. Within the plot on the left, only the SM-like Higgs mass constraint, the invisible Z width, and the chargino mass bound were taken into consideration. The plot on the right was obtained after applying the entire set of constraints. It is evident from the figures that the range of κ becomes narrow as ρ changes. This is due to the fact that some cancellations have to happen between λ, κ and ρ in order for Br(h 0 → inv) 40%. However, we still manage to get a large enough range of χ 0 1 masses for a wide enough range of κ and ρ parameter values. This is also true for the mass of the lightest pseudoscalar, and as we will show in the next section, the annihilation χ 0 1 χ 0 1 → A N →ll,qq can be efficient enough to generate the right density of dark matter. The situation is a bit more constrained for M 2 = 500 GeV. In this case, the next-to-lightest neutralino has a larger wino and bino component and a finer cancellation among parameters is necessary to satisfy the constraint on the invisible decay width of h 0 . Within this benchmark scenario, after all constraints have been applied, the lightest neutralino has a mass of 8 GeV and the allowed values for κ and ρ are −0.05 and 0.023 respectively. The corresponding value of the lightest pseudoscalar mass is 30 GeV.  In the following section we study the cosmological abundance of a light neutralino with mass below 15 GeV that annihilates into SM particles to produce the observed density of dark matter.

A Dark Matter Candidate
The LSP of SUSY models with exact R−parity is known to be a good candidate for cold dark matter [53]. In general, the LSP is a weekly interacting massive particle (WIMP) and, depending on the specifics of the model, it could be a neutralino, the gravitino, an sneutrino or an axino, among others. Particularly, in the context of the NMSSM, the LSP is commonly the lightest neutralino, which has a large fraction of singlino [54][55][56][57][58]. This favors a light DM candidate, with mass below 20 GeV, in the PQ limit or when a continuous R-symmetry is imposed. In such class of models, the relic density is obtained through annihilation into a light scalar or a pseudoscalar Higgs boson [52,[59][60][61][62].
The abundance of thermal relics, X, in the universe is determined by their selfannihilation in relation to the expansion rate of the universe. In the early universe, these particles are abundant and are in thermal equilibrium with the rest of degrees of freedom. When the expansion of the universe dominates over the annihilation rate, and the universe cools down to a temperature below m X , the interaction among DM particles is less efficient and their density "freezes out". The evolution of the comoving particle density is given by the Boltzmann equation [63] dn where H is the Hubble rate and σ XX v is the thermally averaged annihilation cross section. The freeze-out temperature, T F O , at which the particles depart from equilibrium, can be found by solving numerically equation (4.1). This is, approximately, where g * is the number of relativistic degrees of freedom at the freeze-out temperature. Subsequently, the present day relic abundance is given by Here, h is the Hubble parameter in units of 100 km s −1 Mpc −1 . It is convenient to express this relic abundance in terms of the Taylor expansion of the cross section σ XX v ≈ a+bv 2 , .

(4.4)
We now apply this analysis to our model by considering the lightest neutralino, χ 0 1 , as the DM particle. In order to find its abundance, we calculate the annihilation cross section in the same fashion as it was done in [64,65], where the neutralino relic density was computed for the MSSM. As we have already mentioned, since the next-to-lightest neutralino is much heavier than the LSP, we do not include co-annihilations in our calculations.
For convenience, the function w(s) is defined where s is the Mandelstam variable. For annihilations into a two-body final state χχ → f 1 f 2 , w(s) is given by where θ(x) is the Heaviside function, c is a color factor (3 a quark-antiquak final state, 1 otherwise), andw with (4.8) Once σ(s) is obtained, the thermally averaged cross section can be computed using where K 1,2 are modified Bessel functions.
The LSP in the model presented in this work is mostly singlino and very light. This implies that the kinematically allowed annihilation processes are those where the final states are light MSSM fermions. Thus, in the final state, we consider u, d, c, s, b quark-antiquark pairs and lepton ¯ -pairs. The important processes involved in the calculation of the χ 0 1 relic abundance are s-channel annihilations through a Higgs-like scalar (h i and A i ) or a Z boson. In the case of a CP-even scalar, h i , exchange, the contribution is given bỹ where the couplings C f f j S are obtained by inserting the mixing matrix in Equation (2.18) in the Lagrangian. The values of C χχ j S are given in Equation (2.22). On the other hand, the (CP-odd) pseudo-scalar A i exchange yields the s−wave contributionw And, finally, the Z exchange contribution is given bỹ These results altogether give us the cross section that determines the density of χ 0 1 as expressed in Equations (4.5) and (4.6), wherẽ We explore what values of the parameters in our model yield a relic abundance that is consistent with the measured DM density, which corresponds to a thermally averaged annihilation cross section σv ∼ 3 × 10 −26 cm 3 /s. To do so, we scan the parameter space over the ranges presented in Table 2 and impose the invisible Z decay constraints. Also, we require a realistic Higgs mass, m h ≈ 126 GeV, and that m χ + > 104 GeV to be consistent with collider results.
Our findings show that a lightest neutralino with mass between 4 GeV and 9 GeV yields the appropriate relic density, as shown in Figure 6(a). for this mass range, there is a significant component of χ 0 1 alongH u , as depicted in Figure 6(b). Additionally, the annihilation is dominated by the (s-wave) interchange of a light CP-odd scalar, given in Equation (4.11), and as depicted in Figure 7(b), while the contribution from the CP-even scalar mediated annihilation is p-wave suppressed. Therefore, the presence of light pseudoscalars in this model aides in making the annihilation of the relics efficient, avoiding an overabundance of the DM particles. The DM mass obtained in this analysis is significantly smaller than most MSSM neutralino-like proposals, and it is also consistent with the studies of light DM in the NMSSM [52,[59][60][61][62]. The final products from the annihilation of this light neutralino are¯ pairs or light qq pairs. This is shown in Figure 7(a), where the dominant process is that of annihilation into a pair of quarks, specially for neutralino masses close to 10 GeV. In particular, the dominant process for m χ 0 1 4 GeV yields abb pair final state, whereas below this mass the most relevant products aredd,ūu andss pairs.
Let us now take a look at the effects of the wino mass, M 2 , on the allowed neutralino mass values consistent with a thermalized cross section of 3 × 10 −26 cm 3 /s. For M 2 =  which in our case yields a value of the order of 10 −48 − 10 −46 cm 2 for m χ 0 1 ≈ 10 GeV. This is below the range of cross sections that the current direct detection experiments are able to measure, σ SI ∼ 10 −46 cm 2 [66]. However, this cross section sensitivity might be achieved in future detectors. Notice that cross sections of this magnitude are just below the values for which the irreducible neutrino background would affect the discrimination capabilities of the detector [67]. Despite the fact that in our model the DM particle has a mass that is close to the mass hinted by signals detected recently in CDMS [41], where the DM-nucleon scattering cross section in their detected events is about 10 −41 cm 2 , which is far from the expected cross section in our model. Recently, the LUX experiment has released the results of their first WIMP search [68]. They find an upper bound for the annihilation cross section of 7.6 × 10 −46 cm 2 for a WIMP mass of 30 GeV. The limits corresponding to the range of masses considered in this work are between 10 −44 cm 2 and 10 −45 cm 2 .

Conclusions
Extensions of the MSSM have been extensively used in the literature to solve the µ and little hierarchy problems. In this article, we have explored the Higgs and neutralino sectors for an extension of the MSSM, in which those problems are easily addressed. We performed a scan of the parameter space and found the regions that are consistent with collider constraints and a Higgs mass around 126 GeV. In the Higgs sector, we have found two singlet-like scalars that are allowed by present constraints. In the neutralino sector, we have investigated the existence of a light dark matter candidate and its annihilation cross section. In fact, the dark matter particle is "mostly" the fermionic partner of a singlet scalar that does not couple directly to the ordinary matter, but only through a small coupling to the usual singlet present in the NMSSM.
This relic particle turns out to have a mass in the range 8 GeV < m χ < 15 GeV, which is considerably lighter than candidates for dark matter in the MSSM. Its interaction is also remarkably weak, more than the expected interaction in the usual WIMP scenarios. However, the presence of the new singlet-like scalars, and specially the lightest pseudoscalar, favors the annihilation process, and the right relic abundance can be obtained for a wide region of the parameter space. This provides an example of a scenario where the dark matter is somewhat hidden, with the singlet field S field acting as a portal to the MSSM matter content. Along these lines, we found that the cross section sensitivities of the current direct detection experiments are just above the estimated scattering cross section of this dark matter particle with the nucleons, which makes the detection of this type of relic unachievable at present, but it could be tested in future experiments. Finally, this model has been studied at the phenomenological level; it would be interesting to explore the completion at high energies such as embedding this construction in a gauge mediated SUSY breaking scenario, similar to that presented in [35] for the S-MSSM.

A Renormalization Group Equations
In this appendix we give the renormalization group equations to one loop order using the conventions found in [13]: β y ijk = dy ijk dt = γ i n y njk + γ j n y ink + γ k n y ijn .