Spin-dependent constraints on blind spots for thermal singlino-higgsino dark matter with(out) light singlets

The LUX experiment has recently set very strong constraints on spin-independent interactions of WIMP with nuclei. These null results can be accommodated in NMSSM provided that the effective spin-independent coupling of the LSP to nucleons is suppressed. We investigate thermal relic abundance of singlino-higgsino LSP in these so-called spin-independent blind spots and derive current constraints and prospects for direct detection of spin-dependent interactions of the LSP with nuclei providing strong constraints on parameter space. We show that if the Higgs boson is the only light scalar the new LUX constraints set a lower bound on the LSP mass of about 300 GeV except for a small range around the half of $Z^0$ boson masses where resonant annihilation via $Z^0$ exchange dominates. XENON1T will probe entire range of LSP masses except for a tiny $Z^0$-resonant region that may be tested by the LZ experiment. These conclusions apply to general singlet-doublet dark matter annihilating dominantly to $t\bar{t}$. Presence of light singlet (pseudo)scalars generically relaxes the constraints because new LSP (resonant and non-resonant) annihilation channels become important. Even away from resonant regions, the lower limit on the LSP mass from LUX is relaxed to about 250 GeV while XENON1T may not be sensitive to the LSP masses above about 400 GeV.


Introduction
Weakly Interacting Massive Particle (WIMP) has been considered as one of the most attractive candidates for dark matter (DM). WIMPs have been intensively searched for in direct detection experiments and the limits on the DM scattering cross-section on nuclei improved by several orders of magnitude over the last decade. Currently the most constraining limits come from the LUX [1] experiment which will be soon improved by the XENON1T [2] experiment which is expected to be superseded by the LZ [3] experiment in near future. In many theories, WIMP interactions with nuclei are mediated by the Higgs and Z 0 bosons. The electroweakstrength couplings of WIMP to the Higgs and Z 0 bosons were excluded already few years ago by XENON100 [4] as noted in Ref. [5]. The recent LUX constraints pushed models to regions of parameter spaces where these couplings are strongly suppressed i.e. to the vicinity of so-called blind spots in direct detection.
Supersymmetric (SUSY) extensions of the Standard Model (SM) generically provide a WIMP candidate in the form of the lightest neutralino which is often the lightest sparticle (LSP). Some earlier studies of neutralino DM include Refs. [6]- [25]. Blind spots for neutralino DM has been identified e.g. in Refs. [5,26,27,28] for Minimal Supersymmetric Standard Model (MSSM) and in Refs. [29,30] for Next-to-Minimal Supersymmetric Standard Model (NMSSM). Several recent papers emphasized a big impact of new LUX constraints on the spin-independent (SI) scattering cross-section on the parameter space of MSSM [31,32] and NMSSM [33,34,35]. A universal conclusion of these papers is that viable points in the parameter space still exist but they reside very close to blind spots for SI scattering cross-section.
In the present article we study implications of the assumption that the SI scattering of the LSP is so small (below the neutrino background) that probably it will never be detected in direct detection of its SI interactions with nuclei. We also demand that the LSP has thermal relic abundance in agreement with the Planck measurement Ωh 2 ≈ 0.12 [36]. These two assumptions lead to interesting predictions for the model parameter space. We focus on NMSSM with singlino-higgsino LSP but many of our conclusions are valid also for more general singlet-doublet DM models (studied e.g. in Refs. [37]- [40]). We investigate how the resulting parameter space can be constrained by direct detection experiments focusing on the spin-dependent (SD) LSP interactions with nuclei. We assess the impact of new LUX results presented at the Moriond 2017 conference [41] on the parameter space, as well as sensitivity of future XENON1T and LZ experiments. Rather than doing huge numerical scans of the NMSSM parameter space we study separately several classes of SI blind spots identified in Ref. [29]. This way it becomes possible to understand which effects have the biggest impact on the constraints. In particular, we emphasize the role of light singlets in relaxing the constraints.
The rest of the article is organized as follows. In section 2 we introduce the model and conventions used. In the next two sections we focus on general NMSSM. In section 3 we discuss the case in which singlet-like states are heavy and SI cross-section is solely determined by the exchange of the Higgs boson. In section 4 we discuss how both the relic density and blind spot condition is affected by the presence of light singlets. In section 5 we analyze the Z 3 invariant NMSSM. We reserve section 6 for our conclusions.

Model and conventions
In this section we collect formulae useful for the analysis performed in the rest of this work. We adopt conventions as in the previous paper [29] where more details may be found. We start with the most general NMSSM with the superpotential and the soft terms given by where S is a chiral SM-singlet superfield. In the simplest version, known as the scale-invariant or Z 3 -symmetric NMSSM, m 2 3 = m 2 S = ξ S = 0 while f (S) ≡ κS 3 /3. In more general models f (S) ≡ ξ F S + µ S 2 /2 + κS 3 /3.
The mass squared matrix for the neutral CP-even scalar fields, in the basis ĥ ,Ĥ,ŝ related to the interaction basis by a rotation by angle β (see [29] for details), reads: where M 2 hĥ = M 2 Z cos 2 (2β) + λ 2 v 2 h sin 2 (2β) + ∆ĥĥ, ∆ĥ iĥj are radiative corrections, v s , v h sin β and v h cos β are VEVs of the singlet and the two doublets, respectively. The mass eigenstates of M 2 are denoted by h i with h i = h, H, s (h is the 125 GeV scalar discovered by the LHC experiments). These mass eigenstates are expressed in terms of the hatted fields with the help of the diagonalization matrixS: The mass squared matrix for the neutral pseudoscalars, after rotating away the Goldstone boson, has the form where Diagonalizing M 2 p with the matrixP one gets the mass eigenvalues a j = a, A: After decoupling of the gauginos (assumed in this work) the neutralino mass sub-matrix describing the three lightest states takes the form: Trading the model dependent term ∂ 2 S f for one of the eigenvalues, m χ j , of the above neutralino mass matrix we find the following (exact at the tree level) relations for the neutralino diagonalization matrix elements: where N j3 , N j4 and N j5 denote, respectively, the two higgsino and the singlino components of the j-th neutralino mass eigenstate while j = 1, 2, 3 and |m χ 1 | ≤ |m χ 2 | ≤ |m χ 3 |. Using the last two equations and neglecting contributions from decoupled gauginos one can express the composition of the three lighter neutralinos in terms of: λv h , m χ /µ and tan β. Later we will be interested mainly in the LSP corresponding to j = 1, so to simplify the notation we will use m χ ≡ m χ 1 . The physical (positive) LSP mass is given by m LSP ≡ |m χ |.
In the present work we are interested mainly in two properties of the LSP particles: their cross-sections on nucleons and their relic abundance. The spin-dependent LSP-nucleon scattering cross section is dominated by the Z 0 boson exchange and equals σ (N ) where C (p) ≈ 4, C (n) ≈ 3.1 [42]. The spin-independent cross-section for the LSP interacting with the nucleus with the atomic number Z 0 and the mass number A is given by where µ 2 red is the reduced mass of the nucleus and the LSP. When the squarks are heavy, as assumed in the present work, the effective couplings f (N ) (N = p, n) are dominated by the t-channel exchange of the CP-even scalars [43]: Further details may be found in Appendix A. Formulae for the LSP annihilation cross section and its relic density are much more complicated (some of them are collected in Appendices B and C).
The main goal of the present work is to identify regions of the NMSSM parameter space for which the singlino-higgsino LSP particles fulfill three conditions: 1) have very small, below the neutrino background, SI cross-section on nuclei (SI blind spots); 2) have small SD cross-sections to be consistent with present experimental bounds; 3) have relic density close to the experimentally favored value Ωh 2 ≈ 0.12 so the LSP can play the role of a dominant component of DM. 1 Of course, such points in the parameter space must be consistent with other experimental constrains e.g. those derived from the LHC [47] and LEP data [48,49]. In the next sections we discuss solutions fulfilling all the above mentioned conditions, starting with the simplest case of blind spots without interference effects and then for blind spots for which such effects are crucial. We investigate also modifications present in Z 3 invariant NMSSM.

Blind spots without interference effects and relic density
In this section we consider situation when the SI blind spot (BS) takes place without interference effects i.e. when all three contributions to the effective coupling f (N ) in eq. (21) are very small. Two of them are small because the corresponding scalars, s and H, are very heavy while f h is suppressed due to smallness of α hχχ . In the next subsection we discuss the simplest case in which the mixing among the scalars may be neglected. Then the effects of such mixing will be taken into account.  [45,46]) as functions of m LSP and tan β for blind spots with decoupled all Higgs particles except the SM-like one. Dashed and dotted blue lines correspond to eq. (27) and (28), respectively, after substituting (23), for Ωh 2 = 0.12. Red region depicts the points with m χ ± 1 < 103 GeV which are ruled out by LEP [48]. Yellow area is forbidden because of eq. (30), whereas green/cyan one due to LUX/IceCube (IC) limits on spin-dependent LSP interaction with nucleons. Vertical red lines correspond to (from left to right) m LSP = m Z 0 /2, m h /2, m W , m t .

Without scalar mixing
When the interference effects and scalar mixing may be neglected the SI blind spot condition has the following simple form: which corresponds to vanishing Higgs-LSP-LSP coupling for m s , m H → ∞. In Fig. 1 the dependence of Ωh 2 on m LSP and tan β for such blind spots is shown for some specific values of λ and κ (the latter parameter does not influence the situation as long as the resonance with the lightest pseudoscalar is not considered). Some parts of the (m LSP , tan β) plane are excluded by the upper bounds on σ SD from LUX [41] and IceCube [50] experiments (and also by the LEP data). Particularly important are the new LUX constraints which exclude large part of the parameter space with Ωh 2 = 0.12 ± 0.02. We note that in the allowed part of the parameter space presented in Fig. 1 correct thermal relic abundance is obtained for singlino-dominated LSP. It is not difficult to understand the results shown in Fig. 1 using (approximate) analytic formulae.
For given values of m χ and tan β, the blind spot condition (22) together with eqs. (17) and (18) may be used to obtain the LSP composition. For example, the combination which determines the LSP coupling to the Z 0 boson is given by (see Appendix A) 2 All three components (gauginos are decoupled) may be expressed in terms of the above combination using The above expressions are valid as long as blind spot condition (22) is satisfied. For the LSP masses for which the annihilation cross section is dominated by the s-channel Z 0 exchange this is enough to calculate the LSP relic density to a good accuracy. The approximate formulae are (see the Appendix B for the details): for m LSP of order m Z 0 /2 (and below the W + W − threshold) and above the tt threshold. They reproduce very well the relic density calculated numerically using MicrOMEGAs, as may be seen in Fig. 1. The composition of the LSP is crucial not only for its relic density but also for some of the experimental constraints. One gets the following upper bounds on the combination (23): from the LEP bound for invisible decays of Z 0 [49,18] and   28)). Yellow area is forbidden because of eq. (30). Grey color in the left plot denotes the regions in which h resonance and W + W − /Z 0 Z 0 channels may be important and affect the results. Green/cyan areas correspond to LUX/IceCube limits on spin-dependent LSP interaction with nucleons [41], [50]. Dashed (continuous) green lines in both plots correspond to the precision of the future XENON1T (LZ) experiment of SD direct interaction of the LSP with neutrons [51].
One can find the allowed range of m LSP in the vicinity of m Z 0 /2. The points in the left plot in Fig. 2 at which the blue and red dashed lines enter the green region (excluded by LUX constraints on the SD cross-section) determine the limiting values of m LSP for which the Z 0 resonance may give the correct relic density of the LSP. To be more accurate we used the results for Ωh 2 obtained from MicrOMEGAs and found the limiting LSP masses to be approximately 41 and 46.5 GeV. 3 Substituting the corresponding values of |N 2 13 − N 2 14 | into eq. (23) we find the following linear dependence: The values of tan β necessary to obtain good relic abundance of the LSP become larger when moving to values of m LSP closer to the peak of the resonance (which is slightly below m Z 0 /2). The situation is illustrated in the left plot of Fig. 3. One can see that in the region allowed by LUX large tan β is required unless λ is small (∼ 0.1 in this case). Of course the above limits may become stronger (i.e. for a given λ larger tan β might be required) when SD direct detection experiments (especially based on the interactions with neutrons) gain better precision. For instance, the LZ experiment will be able to explore the entire region considered here -see the left plot in Fig. 2.  (23) whereas dotted lines in the right plot from eq. (28) after substituting (23)) for blind spots with decoupled all Higgs particles except the SM-like one. Left plot: green areas and thick red lines denote the points excluded by LUX and the LEP chargino searches, respectively; thick blue points correspond to eq. 31 and 32. Right plot: green region and thick orange line depict the points excluded by LUX and points for which the stop masses above 5 TeV are necessary to obtain the correct Higgs mass (even when the contribution from the stop mixing is maximized), respectively. The red lines denoted by LP depict regions with a Landau pole below the GUT scale.
In the LSP mass range between the W + W − and tt thresholds the annihilation cross section is dominated by gauge boson (W + W − /Z 0 Z 0 ) final states with the chargino/neutralino exchanged in the t channel. The related couplings are proportional to the higgsino components of the LSP (the gauginos are decoupled) which for the blind spot (22) are related to the LSP-Z 0 coupling by eqs. (24) and (25). The values of N 14 (N 13 is smaller by factor 1/ tan β) necessary to get Ωh 2 ∼ 0.12 lead to too large σ SD and are excluded by both LUX and IceCube data (see Fig. 1). Thus, the LSP masses in the range m W + m LSP m t are excluded. The only way to have correct relic abundance consistent with all experimental constraints is to go to very small values of λ in order to suppress SI cross-section below the LUX constraint also away from the blind spot (22) and increase tan β such that the Higgs mass constraint is fulfilled. 4 More flexibility in the parameter space may appear if some additional particles are exchanged and/or appear in the final state of the LSP annihilation (such situations will be discussed in the next sections).
For m LSP 160 GeV annihilation into tt (via s-channel Z 0 exchange) starts to be kinematically accessible so smaller higgsino component suffices to have large enough annihilation cross-section to fit Ωh 2 ≈ 0.12. In consequence, Ωh 2 ≈ 0.12 is obtained for somewhat larger tan β than between the W + W − and tt thresholds and smaller SD cross-section is predicted. As a result, IceCube [50] constraints are satisfied for m LSP 175 GeV. However, the new LUX constraints exclude m LSP up to about 300 GeV for Ωh 2 = 0.12. This lower bound on m LSP may change by about 50 GeV when the uncertainties in the calculation of the relic abundance are taken into account. It becomes stronger (weaker) for smaller (bigger) values of Ωh 2 . We should also emphasize that the lower bound on the LSP mass from SD constraints are the same for the whole class of singlet-doublet fermion DM as long as it annihilates dominantly to tt. In particular, similar lower bound of 300 GeV on the LSP mass was recently set by LUX on the well-tempered neutralino in MSSM [31]. One can also see in Fig. 2 that the correct relic abundance requires |N 2 13 −N 2 14 | ∼ 0.05 which, depending on tan β, translates to N 2 15 ∼ 0.9−0.95 (see eq. (26)) and such values may be explored by XENON1T. The right panel of Fig. 3 shows values of tan β necessary to get Ωh 2 = 0.12 as a function of m LSP and λ. Contrary to the Z 0 resonance case small values of tan β are preferred and hence moderate or large λ (in order to have big enough Higgs mass at the tree level). However, too small values of tan β lead, for a given big value of λ, to a Landau pole below the GUT scale. Thus, the assumption of perturbativity up to the GUT scale and the requirement Ωh 2 = 0.12 give constraints which result in a λ-dependent upper bound on the mass of the LSP. For example, m LSP 700 GeV for λ = 0.7 (see Fig. 3) and m LSP 800 GeV for λ = 0.6. Let us also note that for large LSP masses coannihilation becomes non-negligible. This effect relaxes the upper bound on tan β and is increasingly important as λ decreases, as can be seen in Fig. 3 from comparison of full result by MicrOMEGAs and the approximated one with only tt included.
Let us comment on two features of the Ωh 2 ≈ 0.12 curves in Fig. 1. First: there are no signs of a resonant annihilation with the h boson exchange in the s-channel. This is simply the consequence of the blind spot condition leading to vanishing (or at least very small) LSP-Higgs coupling. This is characteristic of all blind spots without interference effects. Second: tan β decreases with m χ for all m χ > m Z 0 /2 with the exception of the vicinity of the tt threshold. This is related to the fact that the annihilation cross section is directly (Z 0 in the s channel) or indirectly (the V V final states) connected to the value of (N 2 13 −N 2 14 ) given by eq. (23). The r.h.s. of (23) is a decreasing function of m χ and decreasing function of tan β (for tan β > 1+λv h /m χ ). Thus, in order to keep it approximately unchanged the increase of m χ must be compensated by the decrease of tan β (other parameters determining the annihilation cross section may change this simple relation only close to the tt threshold and below the Z 0 resonance).
Another comment refers to constraints obtained from the indirect detection experiments. The IceCube upper bounds on σ SD change by orders of magnitude depending on what channels dominate the LSP annihilation. This can be already seen in the simple case discussed in this subsection. The lower bound on tan β obtained from the IC data (as a function of m LSP ) visible in Fig. 1 drops substantially above the tt threshold because the tt pairs give softer neutrinos as compared to the W + W − pairs [50]. Moreover, the latest LUX results on σ SD lead to stronger bounds in almost all cases. Only for quite heavy LSP the IC limits are marginally stronger, as may be seen in the right panel of Fig. 2.
To sum up, in this section we identified two crucial mechanisms (Z 0 resonance and annihila-tion into tt) which may give correct relic density and are allowed by the experiments. However, both of them rely on the Z 0 boson exchange in the s channel and therefore are proportional to the LSP-Z 0 coupling, which controls also the SD cross section of the LSP scattering on nucleons. Therefore, the future bounds on such interaction will be crucial in order to constrain the parameter space. In fact, XENON1T is expected to entirely probe regions of the parameter space in which annihilation into tt dominates while LZ will be able to explore the entire region of Z 0 resonance. It is also worth noting that the situation presented in Fig. 1 may change if we consider light pseudoscalar a with mass m a ∼ 2m χ . Such resonance for singlino-dominated LSP (we require κ = 0) is controlled mainly by the mixing in (pseudo)scalar sector and hence may not be so strongly limited by the SD direct detection experiments. For instance, we checked with MicrOMEGAs that for m a in a few hundred GeV range we can easily obtain points in parameter space with correct relic density and σ SD below the future precision of the LZ experiment. In principle, the effect of light pseudoscalar may be also important for 2m χ m h + m a when the LSP starts to annihilate into ha state which, depending on κ, may suppress the tt channel and may weaken the IceCube limits. However, in the case considered in this subsection, the contribution from the ha channel may be important only for large mixing in the pseudoscalar sector. This requires quite large values of A λ which leads to unacceptably small values of the Higgs mass. We will come to these points in the next sections where the annihilation channels involving the singlet-dominated pseudoscalar may play a more important role.

With scalar mixing
Next we consider the case when the mixing among scalars is not negligible and affects the blind spot condition (22), which is now of the form: where η, defined by depends on the LSP composition (some formulae expressing η in terms of the model parameters may be found in Appendix A). In equation (33) we introduced also parameter γ describing the mixing of the SM-like Higgs with the singlet scalar. This mixing can be expressed (for m s m h assumed in this section) in terms of the NMSSM parameters as: where ∆ mix ≡ m h −M hh is the shift of the SM-like Higgs mass due to the mixing [52]. For m s > m h this shift is negative so we prefer it has rather small absolute value. The scenario of higgsino-dominated LSP with Ωh 2 ≈ 0.12 is very similar to the analogous case in the MSSM model and requires |µ| ≈ 1.2 TeV. Even the present results from the direct and indirect detection experiments constrain possible singlino admixture in the higgsino-dominated LSP to be at most of order 0.1. So small singlino component leads to negligible changes of µ necessary to get the observed relic density of DM particles. 5 Thus, similarly as before, we focus on SI blind spots with Ωh 2 ≈ 0.12 for singlino-dominated LSP. In this case and for non-negligible κ (more precisely, when |N 13 N 14 | |κ| λ N 2 15 ) the blind spot condition may be approximated by: This condition is a quadratic equation for µ and has solutions only when One can see from (36) that for γκµ > 0 we have always m χ µ > 0 and the LSP has more higgsino fraction than when condition (22) holds. In the opposite case i.e. γκµ < 0 we can have either m χ µ > 0 with slightly smaller higgsino fraction or strongly singlino-dominated LSP with m χ µ < 0. However, for values of |γ| small enough not to induce large negative ∆ mix the higgsino component of the LSP with m χ µ < 0 is too small to obtain Ωh 2 ≈ 0.12. Therefore, from now on we focus on the case m χ µ > 0. Solving eq. (36) for the ratio m χ /µ and substituting the solution into eq. (23) one can find the difference of two higgsino components of the LSP. For small |γ| it can be approximated as (23) with a small correction: As we already mentioned discussing eq. (23), the first term in the r.h.s. of the last equation is a decreasing function of tan β (unless tan β is very close to 1). So, in order to keep the value of |N 2 13 − N 2 14 | necessary to get Ωh 2 ≈ 0.12, the contribution from the second line of (38) may be compensated by increasing (decreasing) the value of tan β when γκµ > 0 (< 0). 6 This effect from the s-h mixing is illustrated in Fig. 4. One can see that indeed, depending on the sign of γκµ, we can have smaller or larger (in comparison to eq. (22)) values of tan β for a given LSP mass, while keeping Ωh 2 ≈ 0.12. In particular, non-negligible Higgs-singlet mixing may relax the upper bound on the LSP mass arising from the perturbativity up to the GUT scale. It is also important to emphasize that the resonant annihilation with s scalar exchanged in the s channel is quite generic for this kind of blind spot (see the right plot in Fig. 4, where we have chosen m s = 600 GeV) because it is easier to have substantial s-h mixing when the  (27) and (36)). We took m χ > 0, µ > 0, m s = 600 GeV and |γ| ≤ 0.025 (corresponding to |∆ mix | 1 GeV), where γ < 0 (> 0) refers to the lower (upper) limit of a given area. The shift of the orange region with respect to the yellow one for LSP masses above about m s is due to the annihilation channels containing s in the final state and co-annihilation effects which become more important for larger LSP masses.
singlet-dominated scalar is not very heavy. Moreover, the presence of resonant annihilation via s exchange can relax the lower limit on the LSP mass from LUX constraints on the SD crosssection, as seen in the right panel of Fig. 4, because in such a case smaller higgsino component is required to obtain correct relic density.

Blind spots and relic density with light singlets
Now we move to the case when the singlet-dominated scalar is lighter than the SM Higgs. Neglecting the effect from the heavy doublet H exchange for the SI cross section (i.e. setting f H to zero) the blind spot condition may be written in the following form [29]: where and parameters measure (in the large tan β limit) the ratio of the couplings, normalized to the SM values, of the h i (= s, h) scalar to the b quarks and to the Z 0 bosons.
In the rest of this section we will consider singlino-like LSP, because the case of higgsinodominated LSP does not differ much from the one described in section 3.2. The blind spot condition (39), analogously to a simpler case (33), may be approximated by a quadratic equation which has solutions only if The last condition may be interpreted as the upper bound on tan β (lower bound on cos β).
It is nontrivial when its r.h.s. is positive i.e. when the second term under the absolute value is negative but bigger than −1. The bound is strongest when that term equals − 1 2 . However, usually the absolute value of that term is smaller then 1 2 because the h-s mixing measured by γ is rather small. So, typically the bound on tan β becomes stronger with increasing LSP mass or increasing |κ| (with other parameters fixed). We focus again on the more interesting case m χ µ > 0 because for m χ µ < 0 the blind spot condition may be satisfied only for the LSP strongly dominated by singlino which typically leads to too large relic density. 7 Since in this section we consider m s < m h , |A s | is typically larger than |γ| (unless c s and/or c h deviate much from 1, which under some conditions may happen which we discuss in more detail later in this section) -in such a case the condition (43) is always fulfilled for γκµ < 0 (see Fig. 5), whereas for γκµ > 0 (see comment in footnote 6 ) there is an upper bound on tan β which gets stronger for larger values of γκ.
In the following discussion we focus on big values of |γ| because they lead to a relatively big positive contribution to the Higgs mass from the Higgs-singlet mixing [52]. In our numerical analysis we take |γ| = 0.4 which corresponds to ∆ mix ∼ 4 GeV. In order to emphasize new features related to the modification of the BS condition we also consider rather large values of |κ| ∼ O(0.1). For such choices of the parameters there are no viable blind spots for γκµ > 0 in accord with the discussion above so we focus on γκµ < 0.
The fact that our blind spot condition now comes from destructive interference between f h and f s amplitudes (rather than vanishing Higgs-LSP coupling) influences strongly the relation between the DM relic density, especially in small LSP mass regime, and other experimental constraints. Now the Higgs-LSP coupling is not negligible so the LSP mass below m h /2 is forbidden, or at least very strongly disfavored, by the existing bounds on the invisible Higgs decays [53]. Thus, resonant annihilation with Z 0 or s (in this section we consider s lighter than h) exchanged in the s channel can not be used to obtain small enough singlino-like LSP relic abundance. As concerning the h resonance: it may be used but only the "half" of it with m LSP m h /2 (this effect is visible in all panels of Fig. 5). However, we found that other experimental constraints, such as the ones from the LHC and/or LUX exclude even this "half" of the h resonance when the mixing parameter |γ| is large.
In general NMSSM the masses of singlet-like scalar s and pseudoscalar a are independent from each other so let us first consider the situation when a is heavy. The case with m a = 1 TeV is presented in the upper left panel of Fig. 5. The contours of Ωh 2 = 0.10, 0.12, 0.14 above the tt threshold are quite similar to the case with heavy singlet. The only difference is that now somewhat larger values of tan β are preferred but even in this case they cannot exceed about 5.
Let us now check what happens when the lighter pseudoscalar is also singlet-dominated (i.e. a 1 = a) and relatively light. The existence of such light pseudoscalar is very important for both the relic abundance of the LSP and the constraints from the IceCube experiment. Let us now discuss these effects in turn.
The DM relic density is influenced by a pseudoscalar in two ways. First is given by possible resonant annihilation with a exchanged in the s channel. This possibility is interesting only for m a m h because for lighter a one still has problems with non-standard Higgs decays constrained by the LHC data (see above) 8 . However, as for any narrow resonance, the DM relic density may be in agreement with observations only for a quite small range of the DM mass (for a given a mass). One can see this in three panels (except the upper left one) of Fig. 5. Second effect is related to new annihilation final states including the singlet-dominated pseudoscalar, namely sa, ha, aa (in addition to similar channels involving only scalars: ss, sh, hh). It is best illustrated in the upper right panel of Fig. 5 for m a = 260 GeV. In this case the sa threshold roughly coincides with the tt threshold. Near this threshold the curves of constant Ωh 2 ≈ 0.12 go up towards bigger values of tan β and leave the region excluded by the LUX data for smaller LSP mass than in the case with heavy singlets. The reason is quite simple. With increasing contribution from the annihilation channels mediated by (non-resonant) pseudoscalar exchange smaller contribution from the channels mediated by Z 0 exchange is enough to get the desired value Ωh 2 ≈ 0.12. Moreover, smaller LSP-Z 0 coupling is obtained for bigger values of tan β so larger values of tan β are preferred than in the case with heavy a. As a result the lower possible LSP mass consistent with the LUX SD limits is almost 100 GeV smaller when a is relatively light. The precise values depend on the relic density and for Ωh 2 = 0.12 is around 200 GeV instead of around 300 GeV as in the case with heavy singlets. The behavior of the Ωh 2 = 0.12 ± 0.02 curves close to and slightly above the tt threshold depends on the parameters. Particularly important is the sign of γ. We see from the right panels of Fig. 5 that even for the same mass of the pseudoscalar, m a = 260 GeV in this case, the plots are very different for different signs of γ. Most differences originate from the fact that γ > 0 implies c h > 1 and c s < 1 while for γ < 0 the inequalities are reversed. There are two important implications of these correlations which we describe in the following.
Firstly, the LHC constraints from the Higgs coupling measurements are stronger for γ > 0 because in such a case the Higgs coupling (normalized to the SM) to bottom quarks is larger than the one to gauge bosons. In consequence, the Higgs branching ratios to gauge bosons is suppressed as compared to the SM. Moreover, non-zero γ results in suppressed Higgs production cross-section so if |γ| is large enough the Higgs signal strengths in gauge boson decay channels is too small to accommodate the LHC Higgs data which agree quite well with the SM prediction. Moreover, a global fit to the current Higgs data shows some suppression of the Higgs coupling to bottom quarks [47] 9 which disfavors c h > 1, hence also large γ > 0. It can be seen from the upper right panel of Fig. 5 that for γ = 0.4 the LHC excludes some of the interesting part of the parameter space which is allowed by LUX due to the LSP annihilations into sa final state. The LHC constraints can be satisfied for small values of γ but this comes at a price of smaller ∆ mix , hence somewhat heavier stops.
Secondly, |A s | is larger for γ < 0 than in the opposite case (see eq. (40)). Moreover, since deviations of c s and c h from 1 grow with tan β, |A s | increases (decreases) with tan β for negative (positive) γ. For γ > 0, this implies that for large enough tan β the r.h.s. of the blind spot condition (42) changes sign and γκµ < 0 is no longer preferred. Equivalently, the upper bound on tan β in eq. (43) gets stronger as tan β grows so it is clear that at some point condition (43) is violated. The appearance of the violation of the blind spot condition is clearly visible at large tan β in the upper panels of Fig. 5 (the black regions). For γ < 0 instead the blind spot condition may be always fulfilled for γκµ < 0 (by taking e.g. appropriate value of µ). 10 However, there is interesting phenomenon that may happen for γ > 0 above the tt threshold which is well visible in the lower right panel of Fig. 5. The values of tan β corresponding to Ωh 2 = 0.12 grow rapidly just above the tt threshold and the there is a gap in the LSP masses for which there are no solutions with SI BS and observed value of Ω. Such solutions appear again for substantially bigger m LSP (above 300 GeV in this case). The reason why Ωh 2 = 0.12 curve is almost vertical near m LSP ∼ m t and the gap appears is related to the fact that m χ /µ varies very slowly with tan β, which results in fairly constant |N 2 13 − N 2 14 | which determines the tt annihilation cross-section, hence also Ωh 2 . The weak dependence of m χ /µ on tan β originates from the fact that for increasing tan β both sides of the blind spot condition (42) grow. The l.h.s. grows because of decreasing sin(2β) while the r.h.s. due to increasing c s . Of course, the fact that these two effects approximately compensate each other relies on specific choice of parameters and does not necessarily hold e.g. for different values of κ.
The presence of light a influences also the IceCube constraints in a way depending on the LSP mass. For m χ (m a + m s )/2 (assuming m a > m s ) the IceCube constraints are very much relaxed and become practically unimportant for the cases discussed in this section. This is so because the additional annihilation channels (into sa, ha, ss, sh) at v = 0 lead to softer neutrinos as compared to otherwise dominant V V channels (or tt channel for even heavier pseudoscalar). The situation is different (and more complicated) for LSP masses between the W + W − and as thresholds. In this region one can have destructive/constructive interference between Z 0 and a-mediated amplitudes 11 for bb annihilation at v = 0 which strengthens/reduces the constraints (IceCube limit on SD cross section are two orders of magnitude stronger for V V than for bb). In our case the effect depends on the sign of κ: for κ < 0 the IceCube limits are strengthen (relaxed) for m χ m a ( m a ) and vice versa for κ > 0 -see the lower panels of Fig. 5.
In the examples considered in this section and shown in Fig. 5 we have chosen the sign of µ to be positive (then m χ is also positive because, as discussed earlier, there are no interesting solutions with m χ µ < 0). The corresponding solutions with negative µ, and also changed signs of other parameters like κ, γ and A κ , are qualitatively quite similar. Of course there are some quantitative changes. Contours in the corresponding plots are slightly shifted towards smaller or bigger (depending on the signs of other parameters) values of tan β. Typically also the regions of unphysical minima are moved towards bigger values of the LSP mass.
As explained in detail in Ref. [29], one can also get vanishing spin-independent cross section when the standard BS condition (22) is fulfilled. Then, the scalar sector has no effect on the blind spot condition. In such a case, we can have as large ∆ mix as is allowed solely by the LEP and LHC constraints, irrespective of the DM sector. Interestingly, the standard BS condition appears also, in a non-trivial way, when |κ| λ is relatively small (we still consider singlinodominated LSP) and both terms in the denominator of (34) are comparable and approximately cancel each other. The blind spot condition (33) requires η to be fixed and not very large. From eq. (54) we see that in such a case a small denominator, for a singlino-dominated LSP, may be compensated by small factor (m χ /µ − sin 2β) which means that the simplest BS condition (22) is approximately fulfilled. In both cases, the analysis performed in subsection 3.1 holds. However, it should be noted that small |κ| weakens the singlet self-interaction ∼ κS 3 and hence may decrease the above-mentioned effects from the a exchange.

Z 3 -symmetric NMSSM
All the analysis presented till this point apply to a general NMSSM. In this section we focus on the most widely studied version of NMSSM with Z 3 symmetry. In this model there is no dimensionful parameters in the superpotential: while the soft SUSY breaking Lagrangian is given by (2) with m 2 3 = m 2 S = ξ S = 0. This model has five free parameters less than general NMSSM which implies that some physical parameters important for dark matter sector are correlated. The main features of Z 3 -symmetric NMSSM relevant for phenomenology of neutralino dark matter are summarized below: • sgn(m χ µ) = sgn(κ) • LSP dominated by singlino implies • Neither singlet-like scalar nor singlet-like pseudoscalar can be decoupled due to the following tree-level relation (for singlino LSP after taking into account the leading contributions from the mixing with both scalars coming from the weak doublets,ĥ andĤ): Masses of both singlet-dominated scalar and pseudoscalar are at most of order m LSP .
• Phenomenologically viable (small) Higgs-singlet mixing leads to the following tree-level relation: where the last approximation is applicable for large tan β and/or singlino-like LSP and forbids resonant LSP annihilation via heavy Higgs exchange. Such resonance may be present only for λ 1 since only in such a case significant deviation from relation (47) is possible. Important constraints on dark matter sector of Z 3 -symmetric NMSSM follow from relation (46) which we discuss in more detail in the following subsections.

Heavy singlet scalar
Let us first discuss the case of heavy singlet scalar in which only the Higgs exchange is relevant in the SI scattering amplitude and the SI blind spot has the standard form (22). In this case |γ| must be close to zero to avoid large negative correction to the Higgs mass and eq. (46) implies m LSP > m s . This is demonstrated in Fig. 6 where it is clearly seen that for m LSP m s there are no solutions (due to a tachyonic pseudoscalar). 12 We also note that eq. (46) implies that resonant LSP annihilation via singlet-like scalar or pseudoscalar is typically not possible in this case. 13 On the other hand, eq. (46) implies that the LSP annihilation channel into sa via a exchange is almost always open (for small mixing γ and m s > m h this channel is kinematically forbidden only in a small region of the parameter space for which m a ≈ 3m s ). This allows for smaller annihilation rate into tt, hence also for smaller higgsino component of the LSP and larger tan β. In consequence, larger LSP masses consistent with Ωh 2 = 0.12 and perturbativity up to the GUT scale are possible than in the case with both singlets decoupled (compare Fig. 6 to Fig. 3). For the same reason large enough LSP masses are beyond the reach of XENON1T, as seen from Fig. 6.

Light singlet scalar
The situation significantly changes when singlet-like scalar is light, especially if the Higgssinglet mixing is not small (which enhances the Higgs mass if m s < m h ). This is because the blind spot condition changes to eq. (42). Moreover, for light singlet the loop corrections to condition (46) can no longer be neglected which under some circumstances allows for resonant LSP annihilation via the s-channel exchange of a.
In the Z 3 -symmetric NMSSM the singlet-dominated pseudoscalar a plays quite important role for the relic density of the singlino-dominated LSP. First we check if and when the s-channel exchange of a may dominate the LSP annihilation cross section and lead to the observed relic density. Of course, this may happen if we are quite close to the resonance, i.e. when m a ≈ 2m LSP . It occurs that it is not so easy to fulfill this requirement in the Z 3 -symmetric model. This is related to the condition (46) which, for m a ≈ 2m LSP and after taking into account the loop corrections in eqs. (6) and (14), may be rewritten in the form The l.h.s. of the above expression is positive so this condition can not be fulfilled without the loop contributions. The last equation may be treated as a condition for the size of the loop corrections necessary to have resonant annihilation of the LSP mediated by the pseudoscalar a.
In order to understand qualitatively the impact of condition (48) on our analysis it is enough to consider the following simple situation: We assume that the scalar mixing γ is negligible and the BS is approximated by (22). On the r.h.s. of eq. (48) we take into account only the first term of the loop correction ∆ŝŝ [56] ∆ŝŝ ≈ 1 2π 2 λ 2 µ 2 ln (the second term is subdominant because m 2 LSP µ 2 for a singlino-dominated LSP and κ 2 < 1 4 λ 2 due to condition (45)). In such approximation and for tan β 1 condition (48) simplifies to This should be compared to the general case when such curves are almost vertical (narrow ranges of the LSP mass but wide ranges of tan β) -see Fig. 5. This difference comes from the fact that in the general model there are more parameters and eq. (46) is not fulfilled. Fig. 7 shows that there are two situations for which BS and correct value of DM relic density are still compatible with the latest bound on DM SD cross-section. One is the above discussed case of resonant annihilation with the light pseudoscalar exchanged in the s channel. The second one occurs for smaller tan β but bigger m LSP and corresponds to annihilation via non-resonant exchange of particles in the s channel. Usually the main contribution to the annihilation crosssection in such a case comes from the exchange of Z 0 boson decaying into tt final state. This process allows to avoid the LUX bounds on σ SD for m LSP above about 300 GeV but is not sufficient to push σ SD below sensitivity of XENON1T, as discussed in section 3. The situation changes when new final state channels, especially as, open. Then not only the present bounds on σ SD may be easily fulfilled but some parts of the parameter space are beyond the XEXON1T reach. We see from Fig. 7 that for light singlets the lower limit on the LSP mass from LUX may be relaxed to about 250 GeV. The effect of annihilation into light singlets is even more important for even heavier LSP so XENON1T may not be sensitive to LSP masses above about 400 GeV. In both cases discussed above the allowed values of tan β are correlated with the LSP mass. The exact form of such correlation depends on theŝ-ĥ mixing parameter γ. Quite generally values of tan β decrease with γ. This is illustrated in Fig. 8 where the bands of allowed tan β as functions of γ are shown for a few values of the LSP mass. This correlation between tan β and γ can be easily understood from eqs. (40)- (42). The first factor on the r.h.s of (42) grows in the first approximation like −γ. This can not be compensated by decreasing κ because in the Z 3 -symmetric NMSSM κ is fixed by the LSP mass. The BS condition (42) with increasing r.h.s. may be fulfilled by decreasing the absolute value of the negative contribution to its l.h.s. i.e. by increasing tan β.

Conclusions
Motivated by the recent strong LUX constraints we investigated consequences of the assumption that the spin-independent cross-section of singlino-higgsino LSP scattering off nuclei is below the irreducible neutrino background. We determined constraints on the NMSSM parameter space assuming that the LSP is a thermal relic with the abundance consistent with Planck observations and studied how present and future constraints on spin-dependent scattering cross-section may probe blind spots in spin-independent direct detection.
In the case when all scalars except for the 125 GeV Higgs boson are heavy the new LUX constraints exclude the singlino-higgsino masses below about 300 GeV unless the LSP mass is very close to the half of the Z 0 boson mass (between about 41 and 46 GeV). In the allowed region LSP dominantly annihilates to tt and tan β must be below about 3.5 (assuming perturbative values of λ up to the GUT scale) with the upper bound being stronger for smaller λ and heavier LSP. There is also an upper bound of about 700 GeV assuming perturbativity up to the GUT scale. We found that XENON1T has sensitivity to exclude the entire region of dark matter annihilating dominantly to tt. This conclusion apply to general models of singlet-doublet dark matter. On the other hand, the LSP resonantly annihilating via Z 0 boson exchange is possible only for large tan β unless λ is very small e.g. for λ > 0.5, tan β 20. Only small range of LSP masses around the resonance of about 2 GeV is beyond the XENON1T reach while LZ is expected to probe Z 0 resonance completely. In all of the above cases the LSP is dominated by singlino. Current and future constraints can be avoided also for very pure higgsino with mass in the vicinity of 1.1 TeV.
The situation significantly changes when singlet-like (pseudo)scalars are light. Firstly, the presence of light CP-even singlet scalar modifies the condition for spin-independent blind spot when its mixing with other Higgs bosons is non-negligible. Depending on the sign of the mixing angle between the singlet and the 125 GeV Higgs preferred values of tan β may be either increased or decreased, as compared to the case with heavy singlet. Interestingly, tan β is increased when the Higgs coupling to bottom quarks is smaller than that to gauge bosons which is somewhat favored by the LHC Higgs coupling measurements.
Secondly, the presence of light singlets opens new annihilation channels for the LSP. As a result, correct relic abundance requires smaller higgsino component of the LSP which relaxes spin-dependent constraints. We found that resonant annihilation via exchange of singlet pseudoscalar is possible even in the Z 3 -invariant NMSSM. Interestingly, even far away from the resonant region the lower limit on the mass of LSP annihilating mainly to tt may be relaxed to 250 GeV. For larger LSP masses sa may become dominant annihilation channel and the LSP masses above 400 GeV may be beyond the reach of XENON1T.
The couplings of the i-th scalar to the LSP and to the nucleon, appearing in the formula (19) for the SI cross-section, after decoupling the gauginos are approximated, respectively, by The LSP couplings to pseudoscalars, important for the relic abundance calculation, are approximated by citereviewEllwanger whereP ij are elements of the matrix diagonalizing the pseudoscalar mass matrix defined in eq. (15). Parameter η defined by eq. (34) and convenient for the discussion of SI blind, using eqs. (17) and (18) With the help of eqs. (17) and (18), the combination of the LSP components crucial for σ SD , (N 2 13 − N 2 14 ), may be written as: We can see immediately that the cross-section disappears in the limit of tan β = 1 or a pure singlino/higgsino LSP. The ratio of the higgsino to the singlino components of the LSP may be calculating from eqs. (17) and (18): Using this relation we may rewrite formula (55) in the form:

B The LSP (co)annihilation channels
In this Appendix we will use the following expansion of σv around v = 0: Then, the relic density may be written as [57]: where x f ≈ 25.

B.1 Resonance with the Z 0 boson (unitary gauge)
Let us consider the LSP annihilation into the SM fermions (except the t quark 14 ) via Z 0 exchange in s channel. The expansion coefficient a and b in eq. (58) are equal to: where g ≡ (g 2 1 + g 2 2 )/2, c F = 1 for leptons and 3 for quarks, whereas β F = 2|q F | sin θ 2 W . The 0 index in b 0 parameter means that we put fermion masses to 0 (which is a very good approximation for m χ ∼ m Z 0 /2; of course a 0 = 0). The sum over the SM fermions (except the t quark) in (61) equal ∼ 14.6. It is worth noting that b 0 ∼ m 2 χ and a ∼ m 2 F (1 − 4m 2 χ /m 2 Z 0 ) which means that b a ∼ 0 (in contrary to naive expectation). Moreover, the terms proportional to higher powers of v 2 in σv (for m χ m F ) are suppressed with respect to bv 2 term in geometric way by v 2 /4. Therefore we can approximate σv ≈ b 0 v 2 and hence expressed the relic density in the form of eq. (59). We will however improve slightly this approach (see Appendix C) and write our formula in the following form: where the term proportional tov 2 ≈ 0.5 2 stems from the fact that the dark matter particles posses some thermal energy during the freeze-out. Eq. (62) reproduces very well the results obtained from MicrOMEGAs far from the resonance (see eg. Fig. 1), however very close to the resonance, especially for m χ m Z 0 /2, the difference may be sizable (Fig. 3).

B.2 Annihilation into tt via Z 0
In this case the dominant contribution also comes from Z 0 exchange in s channel but in contrary to the previous paragraph m χ ∼ m F (= m t ). Therefore the statement that b a is now longer true. It becomes clear when we write down the expression for a and b terms in the limit m χ m Z 0 : One can see that for m χ ∼ m t both terms are comparable whereas for larger m χ we have a/b ≈ 4 and eq. (63) suffices (as we would expect, the terms proportional to higher powers of v 2 are suppressed for m χ m t , m Z 0 as v 2 /4). Similarly to eq. (62) we can find the expression for Ωh 2 . Combining (59) with (63) and (64) we get: The above equation works well for m χ 175 GeV (see Fig. 1), however for m χ ≈ m t we have to be more careful because the expansion in v 2 breaks down. One can see that for m χ m t , m Z 0 the square bracket in (65) equals roughly 1 and Ωh 2 depends on |N 2 13 − N 2 14 | only. Similarly to the case of the resonance with Z 0 , the crucial experimental bounds comes from SD direct detection (see right plot in Fig. 2).
It is worth pointing out that both a and b coefficients in (63) and (64) come purely from −p µ p ν /m 2 Z 0 term in Z 0 propagator. It was noticed long time ago [58,59] that taking into account this term is also crucial for DM annihilation in galactic halos (v 2 ≈ 10 −6 ) for m χ ∼ m Z 0 /2. This is because the a coefficient in (60) vanishes which causes large dip in the annihilation cross section.
C Improved formula for Ωh 2 near a resonance The method described below may be found in [60]. Let us consider a general expression for σv for scalar dark matter (with mass m) annihilating via s channel exchange of a particle with mass M and total decay width Γ: For simplicity we assume we assume α = const which is generally not the case, however we are mainly focused on the effect on Ωh 2 coming from the denominator in (66). Using dimensionless quantities δ ≡ 4m 2 /M 2 − 1, γ ≡ Γ/M and considering non-relativistic approximation s = 4m 2 /(1 − v 2 /4) ≈ 4m 2 (1 + v 2 /4) we get: σv = α/M 4 (δ + v 2 /4) 2 + γ 2 . (67) Let us now define Y (x) ≡ n s , where x = T /m, and write For x d = 25 we havev ≈ 0.5. The above method effectively includes the fact that the dark matter particles posses some thermal energy during their freeze-out. Other cases of f (v) can be also easily analyzed and compared with numerical results.