Singlino-dominated dark matter in $Z_3$-NMSSM

Singlino-dominated dark matter properties are investigated in the $Z_3$ Next-to-Minimal Supersymmetric Standard Model, producing superweak interactions with nucleons involved in dark matter direct detection experiments. Approximate analytical formulas describing the dark matter abundance and cross section in the scattering with nucleons are used to illustrate a dependence on theoretical parameters in neutralino and Higgs sectors. It is shown that the measured abundance requires a sizable singlet--doublet Higgs coupling parameter $\lambda$, while the experimental detection results prefer a small $\lambda$. The parameter space is then surveyed using a nest sampling technique guided by a likelihood function containing various observables in dark matter, Higgs, and B physics, such as the abundance and the scattering cross section. It is demonstrated that dark matter can achieve the correct abundance through $\tilde{\chi}_1^0 \tilde{\chi}_1^0 \to t \bar{t}$ or co-annihilation with higgsinos. The former process provides significantly larger Bayesian evidence than the latter, but this will be examined by the near-future PandaX-4T experiment. If the experiment shows no signs of dark matter, it will become highly disfavored. Furthermore, four cases are summarized to suppress dark matter scattering with nucleons, namely, a small $\lambda$ and three kinds of cancellation between different contributions.


Introduction
The 2012 discovery of the Higgs boson at the Large Hadron Collider (LHC) [1,2] confirmed the correctness of the Higgs mechanism as the origin of the masses of subatomic particles. However, the hierarchy problem caused by sizeable radiative corrections to the Higgs mass term implies there should be new physics between the electroweak scale and the Planck scale. In addition, the Standard Model (SM) cannot explain the existence of dark matter (DM). However, current astronomical observations have confirmed that the universe is composed of 27% DM [3,4]. Therefore, new DM particle candidates and new physics are required beyond the SM. Establishing the structural nature of DM is one of the most fundamental open questions in cosmology and particle physics.
Among multiple proposed theories, the most widely accepted are supersymmetry (SUSY) models with R-parity conservation, i.e., the minimal supersymmetric standard model (MSSM) [5,6] and the next-to-minimal supersymmetric standard model (NMSSM) [7][8][9], which provide elegant solutions to the hierarchy problem by introducing contributions from superpartners to the Higgs mass term. In addition, R-parity conservation ensures that the lightest neutralino is a stable neutral particle if it is the lightest supersymmetric particle (LSP), and it may be an excellent DM candidate. MSSM, the most economic realization of SUSY, exhibits several attractive features, but also includes some challenges (e.g., the "µproblem" [10] and "little hierarchy problem" [11]) that have become exacerbated in recent years by the first run of the LHC experiments. This was particularly true for the uncomfortably large mass of the discovered Higgs boson m h 125 GeV [12][13][14][15][16][17][18][19][20]. Alternatively, the NMSSM solves the µ-problem by adding a singlet chiral superfieldŜ to the MSSM. In this process, the µ parameter is replaced by a dynamic quantity µ ef f = λv s when S develops a vacuum expectation value (VEV) v s , the magnitude of which is naturally at the electroweak scale [8,9]. Furthermore, the SM-like Higgs squared mass can be enhanced by an additional tree-level contribution λ 2 v 2 sin 2 2β and the singlet-doublet Higgs mixing [12,13,15,17], where tan β ≡ v u /v d with v u and v d representing the VEVs of the doublet Higgs fields and v 2 ≡ v 2 u + v 2 d . As a byproduct, the neutralino sector includes a fermionic partner of S (singlino) in addition to neutral electroweak gauginos (bino and wino) and the neutral fermionic partner of the Higgs doublets (higgsinos). Since the coupling of a singlino field with SM particles may be very weak, this study focuses on the case of a singlino-dominated neutralino as a DM candidate [21][22][23][24][25][26]. This scenario is feasible when the Yukawa-like couplings satisfy λ ≥ 2κ and the gauginos are assumed to be heavier than the higgsinos.
• κ ∼ 0.01, λ 0.2, and at least one light singlet-dominated Higgs boson [56]. Here, χ 0 1 annihilates in certain funnel regions and the higgsinos decay in a complex manner to satisfy the LHC constraints.
These conclusions are applicable for |µ ef f | 500 GeV, or equivalently, the fine-tuning criterion ∆m Z 50, where ∆m Z defined in [59] parameterizes the sensitivity of m Z to the SUSY parameters at the weak scale. Given this situation and the fact that the fine-tuning criteria lacks a confirmed scientific basis and may reflect personal prejudice, we update the study by Cao et al. [55] to improve these conclusions. Specifically, we do not require ∆m Z to be less than 50 but impose the condition that µ ef f ≤ 1 TeV. We adopt an advanced MultiNest algorithm [60,61] to perform a sophisticated scan over the parameter space of the Z 3 -invariant NMSSM with a singlino-dominated DM. This algorithm is much more efficient than the other algorithms (e.g., the Markov Chain method [62] adopted in [55]) in providing comprehensive information of the space to reveal the underlying physics, although it usually involves a tremendous amount of calculation. To the best of our knowledge, few researchers have used it to study the NMSSM phenomenology [63][64][65][66][67][68]. We also present a description of DM annihilation and the mechanisms used to suppress the DM nucleon scattering cross section through analytical formulas and numerical analysis. Some of the formulas are new, and some are consistent with excellent related works [69][70][71][72][73]. Evidently, such an analysis is helpful to understand the DM physics.
The remainder of this paper is organized as follows. In Section 2, we briefly introduce the basic properties of the Z 3 -invariant NMSSM, including the Higgs and neutralino sections. We then demonstrate DM annihilation and scattering cross sections for singlinodominatedχ 0 1 with nucleons using analytical formulas. In Section 3, we provide a brief description of our scanning strategy. In Section 4, we investigate predictions for surviving samples and the properties of singlino-dominated DM scenarios to understand their distinctive features. Finally, Section 5 includes a discussion of the results and corresponding conclusions.
2 Next-to-minimal supersymmetric standard model model

Fundamental NMSSM properties
As the simplest extension of the MSSM, the NMSSM includes one additional gauge singlet Higgs fieldŜ. The associated superpotential can be expressed as follows [8,9]: where W MSSM is the MSSM superpotential without the µ term, λ and κ are dimensionless parameters, andĤ u andĤ d are the common Higgs superfields. It is the most general Rparity-conservation superpotential satisfying a Z 3 discrete symmetry given the considered field content. Assuming CP-conservation, the Higgs sector of the Z 3 -NMSSM is determined by six parameters at the tree-level [8,74]:  In the Z 3 -NMSSM, mixtures of bino (B 0 ), wino (W 0 ), higgsino (H 0 d,u ), and singlino (S 0 ) fields form neutralinos. Assuming a basis of ψ 0 = (−iB, −iW 0 ,H 0 d ,H 0 u ,S) produces the following neutralino mass matrix [8]: where M 1 , M 2 , and µ ef f denote the soft breaking masses of the bino, wino, and higgsinos, respectively. Diagonalizing the mass matrix with a unitary matrix N yields five mass eigenstates (ordered by mass): The lightest neutralino,χ 0 1 , acting as the DM candidate is the focus of this work. In the limit that |M 1 | and |M 2 | are much larger than |µ ef f | and v,χ 0 1 is approximated byχ (2.6) If |κ/λ| < 1, the dominant composition ofχ 0 1 is a singlino. In this case, κ is related to mχ0 1 as follows [69,70]: and the elements of the matrix N exhibit the following relationships [69][70][71][72][73]: Thus, The higgsino and singlino fractions inχ 0 1 can be defined as Z h = N 2 13 + N 2 14 and Z s = N 2 15 , respectively. The ratio of Z h to Z s can then be expressed as follows: This expression implies that a small λ can suppress the higgsino fraction inχ 0 1 . The couplings of DM to scalar Higgs states, the Z boson, and the Goldstone boson G 0 are included in the calculation of DM annihilation. They take the following form [8]: where the coefficients are given by where the approximation in Eq. (2.6) was applied. In addition, the coupling of a CP-even Higgs h i to two quasi-pure singlet CP-odd Higgs a s , L NMSSM C h i asas h i a s a s , is relevant to our study. Its coefficient C h i asas is given as follows [77,78]: In most cases, C hsasas C hasas , since v s v.

Dark matter relic density
In the NMSSM, the abundance of the singlino-dominated DM candidate (χ 0 1 ) tends to be unacceptably large, due to small coupling effects with SM particles. However, such a candidate can still achieve the measured abundance [3,4] by a specific mechanism, e.g., via s-channel exchanges of gauge and Higgs bosons and t-channel exchanges of electroweakinos and sfermions. It can also be achieved through co-annihilation with heavier states, such as sleptons, the next-to-lightest neutralino (χ 0 2 ), or the lightest chargino (χ ± 1 ). The thermal abundance of the DM at the freeze-out temperature T F = mχ0 1 /x F is given as follows [73]: with a thermally-averaged annihilation cross section σv x F ≡ a + 3b x F . Dominant contributions to σv x F in acquiring the measured abundance are discussed in Section 2.2 of Baum et al. [73]. The following conclusions were presented: 1)χ 0 1χ 0 1 → tt is usually the most crucial channel for the abundance when mχ0 1 > m t . It proceeds through the s-channel exchanges of Higgs and Z bosons. Since the top quark is massive, the contribution from the Z boson's longitudinal polarization to σv x F is important in this process when the Higgs mediators are far off-shell. σv x F can then be approximately expressed as follows: The measured abundance Ωh 2 ∼ 0.12 [3,4] is achieved through the coupling of DM pairs to the Goldstone boson |Cχ0 is then given as follows [73,79]: where the s-and t-channel contributions are approximated as in the case of |A t | |A s | [73]. Hence, it is evident that once the involved particle masses are fixed, the density is primarily determined by the parameter κ. However, because λ > 2|κ| to ensure a singlino-dominatedχ 0 1 , the measured abundance can set a lower bound on λ.
To date, the sensitivities of the XENON-1T experiments have reached the precision of 10 −47 cm 2 for the SI cross section [50] and 10 −42 cm 2 for the SD cross section [51]. They have strongly restricted the λ 0.3 case (see the discussion about DM-nucleon scattering), which is preferred by the annihilationsχ 0 1χ 0 1 → tt, h s a s to account for the measured abundance. In addition, the LHC searches for new particles [30][31][32][33][34][35][36][37][38][39][40][41][42][43][44][45][46][47][48][49], and its precise measurement of the discovered scalar's properties [75,76] have limited the theoretical prediction on light sparticles and Higgs bosons. This has a significant impact on the DM annihilation channels since they are usually accompanied with light particles to account for the abundance (see the discussion in [55]). Thus, the scenario preferred by the scan of Ref. [55] involves a singlino-dominated DM with λ 0.05. In this case, an effective mechanism to obtain the measured abundance includes co-annihilation with higgsinos. The corresponding reaction isχ iχj → XX , in which XX denotes SM particles andχ iχj may be an LSP-NLSP or NLSP-NLSP annihilation state (e.g.,χ 0 1χ 0 . This mechanism is distinct in that the effective annihilation rate at a temperature T is very sensitive to theχ 0 1 -higgsino mass splitting [79,80], and even for a small λ and κ, it can still explain the measured abundance 3 .

DM-nucleon cross sections
Serving for a WIMP,χ 0 1 might be detected by measuring the recoil of a nucleus after an elastic scattering ofχ 0 1 on a nucleus taking place. In the non-relativistic limit, only two different kinds of interactions between a neutralino and a nucleon need to be considered [81]: the spin-dependent interaction (SD) where the WIMP couples to the spin of the nucleus, and the spin-independent interaction (SI) where the WIMP couples to the mass of the nucleus.
When mq 2 TeV, only the t-channel Z exchange diagram contributes significantly to the spin-dependent (SD) scattering cross section at the tree level, which is approximated is less important than theχ 0 1 -mediated contribution for two reasons. One is that, if |κ| is comparable to λ, |Cχ0 , the former contribution is relatively suppressed by the propagator. 3 Note that the co-annihilation mechanism applies under the premise thatχ 0 1 and the higgsinos remained in thermal equilibrium in the early universe [79,80]. In the Z3-NMSSM, many processes, such asχ 0 1χ 0 1 ↔χiχj,χ 0 1 X ↔χiX , andχi ↔χ 0 1 XX , could keepχ 0 1 in chemical equilibrium withχi, and the conversion rates of some of them might be enhanced if the mediator were around its mass-shell. We add that maintaining the thermal equilibrium does not necessarily require the involved couplings to be moderately large. For example, the equilibrium condition was discussed in Eqs. (4.3) and (4.4) of Ref. [80] in the framework of the DM model ST11. It was found that the involved coupling may be as low as 10 −4 to maintain the equilibrium. by [82,83] σ SD with N = p(n) denoting protons (neutrons) and C p 2.9 × 10 −41 cm 2 (C n 2.3 × 10 −41 cm 2 ) [71,72]. From the expression for Cχ0 4 . Furthermore, in the co-annihilation case, the degeneracy of mχ0 1 and µ ef f leads to a minuscule denominator in Eq. (2.12), which requires a small value of λv/µ ef f to satisfy the DM-DD experimental constraints.
In contrast, the spin-independent (SI) scattering cross section in the heavy squark limit is dominated by a t-channel exchange of CP-even Higgs bosons h i [84][85][86][87] and can be expressed as [88] σ SĨ represent the normalized light quark contribution to the nucleon mass, and f influences other heavy quark mass fractions in nucleons [86,87]. In this study, the default settings for f N q were used in the micrOMEGAs package [85], and they predict Currently, non-SM doublet Higgs bosons H are preferred to be heavier than several hundreds of GeV in LHC experiments. In this case, the contribution from H to the SI cross section will be suppressed by (a q ) 2 H ∝ 1/m 4 H , and it is much smaller than that from h for a not exceedingly large tan β. As such, the primary contribution to SI scattering comes from the t-channel exchange of SM-like Higgs bosons h and the singlet Higgs boson h s . The latter contribution may be crucial when h s is much lighter than h. Since the non-SM doublet components of h and h s are approximately zero, (2. 26) In these expressions, if only the contribution from a pure SM Higgs state is considered, A can be simply expressed as (2.27) It is immediately evident that σ SĨ χ 0 1 −N will vanish for mχ0 1 /µ ef f = sin 2β, which corresponds to a blind spot (BS) condition in [71,73,90].
The above analytical formulas for the SD and SI scattering cross sections suggest that σ SD ∝ λv /µ ef f ) 2 and (by contrast) σ SI depends on λ, µ ef f , and mχ0 1 in a complex way. In general, a larger λ, a smaller µ ef f , and mχ0 1 /µ ef f → 1 will increase σ SI and thus strengthen the DM-DD experimental constraints.

Model scans and constraints
The NMSSMTools-5.4.1 package [77,91] was used to produce samples of singlino-dominated DM scenarios in the Z 3 -NMSSM and to model the corresponding features in detail 4 . A sophisticated scan was first performed over the following ranges in the parameter space: 0 < λ ≤ 0.7, |κ| ≤ 0.7, 1 ≤ tan β ≤ 60, 100 GeV ≤ µ ef f ≤ 1000 GeV, |A κ | < 1 TeV, 0 < A λ ≤ 5 TeV, |A t | ≤ 5 TeV, |M 1 | ≤ 500 GeV, (3.1) in which all parameters were defined at the scale Q = 1 TeV. Upper bounds of 0.7 were imposed on λ and |κ| to maintain a perturbable theory up to the grand unification scale. A lower bound of 100 GeV was placed on µ ef f by the LEP search for electroweakinos [92], and an upper bound of 1000 GeV for µ ef f is large enough to allow us to consider various possibilities (see the discussion presented below). In addition, noting that the LHC search for SUSY prefers massive charged sparticles, the following assumptions were made concerning unimportant SUSY parameters. The electroweak gaugino masses were set to M 2 = 2 TeV, and the gluino masses were set to M 3 = 5 TeV. Soft SUSY-breaking parameters in the squarks sector were fixed at 2 TeV, excluding trilinear couplings A t = A b used as free parameters to adjust the Higgs mass spectrum to coincide with relevant experimental measurements at the LHC. In addition, all slepton soft parameters were set to 2 TeV, as we did not want to explain the muon g-2 anomaly. We also required λ ≥ 2|κ| in the scan to achieve a singlino-dominatedχ 0 1 .

19.5
Specifically, the MultiNest algorithm [60,61] with flat distributions for all the parameters in Eq. (3.1) and nlive = 20000 were adopted during the scan to ensure that the conclusions were as complete as possible, and more than 20 thousands CPU hours were spent on the calculations 5 . Several constraints were imposed by constructing the following corresponding likelihood function to guide the process: where L m h and L h,coupling are likelihood functions for the experimentally measured SM-like Higgs boson mass and couplings, respectively. The computation of m h included leading electroweak corrections, two loop terms, and propagator corrections, as in Ref. [93]. Its experimental central value was taken as m h = 125.09 GeV [20], and a total experimental and theoretical uncertainty of 3 GeV was assumed. L h,coupling works in a seven-parameter κ-framework with related experimental measurements, such as the central values and uncertainties of the Higgs couplings and their correlation coefficients, taken from the AT-LAS analysis, using 80fb −1 data collected during the LHC Run-II [75]. Some knowledge about probability and statistics was used in constructing L h,coupling (see the introduction in Ref. [92]). L B is the likelihood function for the measurement of the branching ratio for the decays B → X s γ and B s → µ + µ − . These ratios were calculated by the formulae in Refs. [94,95], and their experimental values were taken from Ref. [92]. L Ωh 2 and L DD are likelihood functions for the measured abundance from the WMAP/Planck experiments [3,4] and the detection of both spin-independent (SI) and spin-dependent (SD) DM-nucleon scattering in the XENON-1T experiment [50,51]. Relevant quantities were calculated using the micrOMEGAs package [78,[96][97][98][99]. In addition, L EW denotes a likelihood function for precision electroweak observables of i (i=1,2,3) [100][101][102] or, equivalently, S, T , and U parameters [103,104] calculated using the formulas from Cao and Yang [105] and fitted to corresponding measurements by the procedure presented in Ref. [106]. Each of these likelihood functions was assumed to follow a Gaussian distribution, with explicit representations provided in Cao et al. [107]. The acquired samples were further refined using the following criteria: the SM-like Higgs mass was within the range of 122-128 GeV, the observed DM relic abundance was within ±10% of the measured central value Ωh 2 = 0.1187 [4], the upper bound was a 90% confidence level of the XENON-1T experimental results of the SI cross section [50] and the SD cross section [51], and all other constraints implemented in NMSSMTools, including various B-physics observables in corresponding experimentally allowed ranges, were at the 2σ level. We also required the samples to satisfy χ 2 EW ≤ 7.8 and χ 2 h,coupling ≤ 14.1, which corresponded to 95% confidence level exclusion limits for three and seven degrees of freedom, 5 The multi-nest sampling algorithm explores a high-dimensional parameter space by determining the iso-likelihood contour in each iteration with nlive active points (the integer nlive is an input parameter of the algorithm, and it usually takes a value larger than 1000. In general, the larger value nlive adopts, the more accurate the scan result becomes.) [60,61]. It is good at dealing with the case in which the samples' posterior distribution is multi-modal or degenerate, which is frequently encountered in new physics studies. In contrast, the Markov Chain method [62] is highly inefficient for such a situation, and thus, it usually provides incomplete information about the distribution. respectively. Since a more significant deviation of the electroweak precision observables (the Higgs couplings) from their measured values would enhance χ 2 EW (χ 2 h,coupling ), this requirement delineated the experimentally allowed range of these observables to further limit the Z 3 -NMSSM.
Constraints were also implemented from the LHC search for electroweakinos using the SModelS-1.2.3 code [108][109][110]. The final states listed in Table 1 from all the electroweakino pair production processes were considered during this process 6 . It is worth noting that a large portion of the samples satisfying these constraints were characterized by a small mass splitting betweenχ 0 1 andχ 0 2 . However, the latest ATLAS analysis of the compressed mass spectra, acquired by searching soft di-lepton signals [111], was not included in the SModelS-1.2.3 code. As such, the analysis constraints are validated in the Appendix, and they were applied to each sample by elaborate Monte Carlo simulations. We verified that they were very effective in excluding the co-annihilation case. We illustrated this point in our recent publication [88].
We add that we did not consider constraints from indirect DM detection experiments (i.e., the Fermi-LAT search for DM annihilation from dwarf spheroidal galaxies) as they become loose for mχ0 1 > 100 GeV [54]. In addition, in the co-annihilation case encountered in this work, the annihilation rate of singlino-dominated DM in present day is very small, which weakens the constraints.

Numerical Results
A closer analysis suggested that the eventual surviving parameter space could be classified into the following three types: • Type-I samples: the lightest CP-even Higgs boson h 1 as the SM-like Higgs boson, 0.4 λ 0.7, 0.13 κ 0.23, 1.5 tan β 6, 450 GeV µ ef f 720 GeV, and the Bayesian evidence ln Z 1 = −24.2. Due to the limited capability of the Markov Chain algorithm, this type of samples were neglected in [55]. It is evident that the three types of parameter spaces are extremely narrow. The Bayesian evidence Z i (i = 1, 2, 3) produced the Jeffreys' scales [112,113] δ 12 ≡ ln Z 1 −ln Z 2 = 3.3 and δ 13 ≡ ln Z 1 − ln Z 3 = 2.8. These results suggest that the considered experiments slightly preferred Type-I samples to Type-II and -III samples. However, as discussed below, Type-I samples will be examined by the near-future PandaX-4T experiment [114]. Furthermore, Type-II and -III samples share the following features in the parameter space: both λ and |κ| are small, 2|κ|/λ 1, and the higgsinos are moderately light. These facts lead to the similarities of the Type-II and -III samples in many aspects of DM physics (see the discussions below for detailed similarities and differences). Nevertheless, they are still distinguished from each other in Higgs physics [13][14][15][16][17][18].
In the following, we investigate the characteristics of the singlino-dominated DM based on numerical results. In Figure 1, Type-I, -II, and -III samples are projected on the sin 2β − mχ0 1 /µ ef f and λv/µ ef f − Z h planes. In Figures 2, 3, and 4, Type-I, -II, and -III samples are projected on the σ SĨ χ 0 1 −p − mχ0 1 and σ SD χ 0 1 −n − mχ0 1 planes, respectively, with different colors indicating the value of λv/µ ef f . In Table 2, a selection of benchmark points are shown to further clarify features in each of the three scenarios.
The following points about Type-I samples can be determined from Figures 1 and 2 and the point P 1 in Table 2.
• Type-I samples are characterized by a relatively large λv/µ ef f ranging from 0.14 to 0.21 (see right panel of Figure 1). This will increase the higgsino composition inχ 0 1 through Eq. (2.11) and theχ 0 1χ 0 1 Z coupling through Eq. (2.12). As indicated in the right panel of Figure 2, the SD scattering rates are thus larger than 7 × 10 −42 cm 2 , which are exceedingly close to the near-future PandaX-4T exclusion limit.
• As shown in the left panel of Figure 1, mχ0 1 /µ ef f and sin 2β are correlated by mχ0 1 /µ ef f sin 2β. In this case, Cχ0 1χ 0 1 H SM in Eq. (2.14) is suppressed by mχ0 1 /µ ef f and sin 2β cancellation, favored by the stringent bound of the XENON-1T experiment on the SI cross section. We further explore its implication by focusing on the point P 1 in Table 2, which predicts the following four terms in Eq. (2.14): • Based on mχ0 1 /µ, λv/µ, sin 2β, and Z h in Figure 1, it is usually predicted that Cχ0  Table 2, we found Cχ0 1 → h s a s contributed to the total rate by 15%. We add that our estimation roughly agrees with the results in Table 2, calculated by the micrOMEGAs package.
Next we consider Type-II samples. Figures 1 and 3 demonstrate following features: Figure 3. Same plots as in Fig.2, but for the results of the Type-II samples.
• As indicated in the left panel of Figure 1, mχ0 1 /µ ef f and sin 2β are of opposite sign, and thus, their contributions to Cχ0 We add that the SI cross sections in Figure 3 are larger than 3 × 10 −49 cm 2 because smaller cross sections require smaller values of λ, which are not readily available since Bayesian evidence is suppressed significantly. We also add that the characteristics of σ SĨ  Table 2. The four terms in Eq. (2.14) and the two contributions in Eq. (2.25) were as follows:  • Since 2|κ|/λ 1, the singlino-dominatedχ 0 1 co-annihilated with the higgsino-dominated neutralinos and charginos to provide the measured abundance. In addition, because λ and κ must be small to suppress Cχ0 Finally, we investigate the Type-III samples. Figure 1 indicates that some samples correspond to the same parameter space as the Type-II samples, and thus, they predict similar DM physics. This conclusion was verified by studying some points in the two scenarios, e.g., points P2 and P4 in Table 2. In the following, we only concentrate on the samples with same sign mχ0 1 /µ and sin 2β. The following results were obtained from Figures 1 and 4.
• Compared with the Type-II samples, λv/µ ef f may have a greater value of up to 0.12. The reason is there are accidental cancellations contributing to σ SĨ χ 0 1 −p , similar to the Type-I samples. This will relax the XENON-1T constraint. We show this characteristic by considering the point P 3 in Table 2 and finding The cancellation also explains why σ SĨ χ 0 • Due to a relatively large λv/µ ef f and |mχ0 1 /µ ef f | → 1 in Eq. (2.12), σ SD χ 0 1 −n may be on the border of 10 −41 cm 2 (see right panel of Figure 4). Furthermore, similar to the Type-II samples,χ 0 1 co-annihilated with the higgsino-dominated neutralinos and charginos to achieve the measured abundance.
We add that the point P 3 predicts a significantly larger λ than points P 2 and P 4 . Consequently, this will be readily tested in the near-future PandaX-4T experiments for both SI and SD scattering.

Conclusion
This work provided updates to previous studies of singlino-dominated DM mainly in four aspects: (1) We no longer require the fine-tuning measurement ∆m Z as a model selection criterion, since it may only reflect personal prejudice.
(2) We adopted an advanced MultiNest algorithm to perform a sophisticated scan over the Z 3 -NMSSM parameter space to ensure that the obtained conclusions were as complete as possible.
(3) We utilized the latest experimental results to restrict singlino-dominated DM scenarios, including the XENON-1T search for both SI and SD DM-nucleon scattering since 2018, ATLAS analyses of sparticle signals with 139-fb −1 data, and measurements of the Higgs couplings with 80-fb −1 data.
(4) We provided simplified analytical formulas for both the DM annihilation cross sections and the SI and SD cross sections of DM-nucleon scattering, and we numerically scrutinized each contribution in these formulas.
As a result, a new singlino-dominated DM scenario (Type-I samples) was found. More model information, such as its Bayesian evidence, was obtained, and the current and future statuses of the scenarios were presented. More importantly, this study provided clear insight into singlino-dominated DM scenarios and explained why they had been tightly limited in the Z 3 -NMSSM. Specifically, this study indicated that the surviving samples can be categorized into three types: • For Type I samples, 0.4 λ 0.7, 0.13 κ 0.23, 1.5 tan β 6, 450 GeV µ ef f 720 GeV, and the annihilationχ 0 1χ 0 1 → tt is primarily responsible for the DM abundance.
The Bayesian evidence (Z) for the three sample types showed that the experiments slightly preferred Type-I samples to Type-II and -III. However, Type-I samples will be examined in the near-future PandaX-4T experiment. They will become highly disfavored if the experiment shows no signs of DM. It should be emphasized that DM-annihilation by a singlet scalar or pseudo-scalar funnel, a Z-boson funnel, and the SM-like Higgs funnel presented by Abdallah et al. [56] were not observed in this study, due to the small Bayesian evidence. In addition, both analytical formulas and numerical results were used to summarize the theory's four primary cases in significantly suppressing the SI scattering cross section for DM-nucleons. This included 1) a small λv/µ ef f , 2) cancellation between m |tildeχ 0 1 /µ ef f and sin 2β, 3) cancellation between two contributions from the SM doublet component (H SM ) and singlet component (H s ) within theχ 0 1χ 0 1 h coupling strength itself, and 4) cancellation between two contributions from h and h s .
In summary, the interaction of the singlino-dominated DM with nucleons in the Z 3 -NMSSM has been tightly restricted in current DM-DD experiments, while the measured abundance favors its involvement in weak interactions. Given the theory's natural preference for electroweak symmetry breaking, it has become increasingly difficult to represent these two seemingly paradoxical features using neutralino DM, due to the limited theoretical structure. Thus, the benefits of singlino-dominated DM are waning unless one extends the Z 3 -invariant theory. Recent studies on the MSSM and NMSSM imply that DM candidates should be gauge singlet fields or that singlet components should at least be naturally dominant over other components [58]. This requirement motivates us to extend the Z 3 -NMSSM to a general NMSSM to increase the Bayesian evidence of the scenarios significantly [88]. It also motivates us to augment the Z 3 -NMSSM with a seesaw mechanism to generate neutrino masses and select the lightest sneutrino as a DM candidate [58,107,115,116]  In our simulation, we generated 100,000 events for the production process. The results, shown in Table 3, indicate that we can reproduce the ATLAS analysis at the 20% level for most cases. Table 3. Cut flow for the analysis in [111]. We considered the point m(˜ ,χ 0 1 ) = (150, 140) GeV in the calculations.