Exploring the supersymmetric U(1)$_{B-L} \times$ U(1)$_{R}$ model with dark matter, muon $g-2$ and $Z^\prime$ mass limits

We study the low scale predictions of supersymmetric standard model extended by $U(1)_{B-L}\times U(1)_{R}$ symmetry, obtained from $SO(10)$ breaking via a left-right supersymmetric model, imposing universal boundary conditions. Two singlet Higgs fields are responsible for the radiative $U(1)_{B-L}\times U(1)_{R}$ symmetry breaking, and a singlet fermion $S$ is introduced to generate neutrino masses through inverse seesaw mechanism. The lightest neutralino or sneutrino emerge as dark matter candidates, with different low scale implications. We find that the composition of the neutralino LSP changes considerably depending on the neutralino LSP mass, from roughly half $U(1)_R$ bino, half MSSM bino, to singlet higgsino, or completely dominated by MSSM higgsino. The sneutrino LSP is statistically much less likely, and when it occurs it is a 50-50 mixture of right-handed sneutrino and the scalar $\tilde S$. Most of the solutions consistent with the relic density constraint survive the XENON 1T exclusion curve for both LSP cases. We compare the two scenarios and investigate parameter space points and find consistency with the muon anomalous magnetic moment only at the edge of $2\sigma$ deviation from the measured value. However, we find that the sneutrino LSP solutions could be ruled out completely by strict reinforcement of the recent $Z^\prime $ mass bounds. We finally discuss collider prospects for testing the model.


INTRODUCTION
The discovery of the SM-like Higgs boson and nothing else presents a serious challenge to particle phenomenology. On one hand, the SM is incomplete, as it fails to explain properties such as the hierarchy problem, neutrino masses, cosmological inflation and dark matter. On the other hand, a Higgs mass of 125 GeV presents a problem for the Standard Model (SM) (e.g., electroweak vacuum instability), and for most of its extensions. So far, no clear directions for theoretical explorations or experimental solutions are indicated. Constructing and studying models which attempt to solve some of the outstanding problems in SM emerges as a viable alternative. Out of these, supersymmetry presents a partial solution to the hierarchy problem and a clear one for dark matter. However, in its minimal incarnation, the minimal supersymmetric model (MSSM) requires squarks and gluinos in the multi-TeV range to explain such a low Higgs mass, raising a serious challenge for the LHC to find any signals.
This issue may be resolved in models with extended gauge groups. In these models, additional D-term contributions to the Higgs mass matrices weaken considerably MSSM mass limits [1][2][3]. Depending on the models studied, these models can also resolve additional problems of MSSM. For instance, models with left-right symmetry [4] can yield neutrino masses via the seesaw mechanism [5][6][7].
In [8], an extended supersymmetric model based on SU (3) c × SU (2) L × U (1) R × U (1) B−L was proposed. The model can be embedded in SO(10) SUSY-GUT, much like the left-right supersymmetric model, and generate a new seesaw mechanism for neutrino masses. The factor U (1) R can be thought off as remnant of a more complete SU (2) R . Unlike the left-right supersymmetric model, which requires Higgs triplet representations with vacuum expectation values (VEV) v R ∼ 10 15 GeV for obtaining neutrino masses and gauge unification, the symmetry in this model can be broken by singlet Higgs bosons (thought of as remnants of a doublet representation in left-right models), with VEVs in the TeV range, while still allowing for gauge coupling unification. In [8], the smallness of neutrino masses was explained as based on an inverse seesaw mechanism. The general features of the TeV scale soft-supersymmetry breaking parameters were explored in [9], outlining conditions for models with intermediate scales obtained from breaking SO (10). The Higgs sector of the model was further explored, showing that a larger mass than that predicted by MSSM can be obtained. The parameter space was further explored in [10], where benchmarks, branching ratios, as well as lepton violation constraints were analyzed.
In this work, we concentrate on investigating, discriminating, and restricting the parameter space of the model using dark matter studies. We include up-to-date constraints on the spectrum coming from the Higgs signal strengths and mass data, and including LHC restrictions on squark and gluino masses, constraints on flavor parameters from the B sector, as well as recent lower limits on the Z mass. Assuming universal scalar and gaugino masses, we show that the lightest supersymmetric particle (LSP) can be the sneutrino (which is different from the usual in this scenario, being a mixture of the right sneutrino and a gauge singlet fermion introduced to generate the inverse seesaw mechanism); or the lightest neutralino (which is favored to be a mixture of the two U (1) binos). Relic density and indirect dark matter detection severely restrict the parameter space, as indeed does the recent limit on the Z mass [11]. Within the parameter space allowed by dark matter limits, we analyze the consequences on sparticle spectra, the neutral Higgs sector and on the anomalous magnetic moment of the muon, which shows more than 3 σ [12] discrepancy with the SM prediction. Finally we investigate the possibilities of testing the model at the LHC.
Our work is organized as follows. We provide a brief description of the model in Sec. 2, capitalizing on more complete descriptions which have appeared previously. In Sec. 3 we describe in detail the parameters of the model and constraints imposed on them. Dark matter phenomenology is explored in Sec. 4, for both neutralino LSP 4.1 and sneutrino LSP 4.2. We then look at the consequences of our findings and compare the two scenarios in Sec. 5, for the sparticle spectrum, the Higgs sector 5.1 and the anomalous magnetic moment of the muon 5.2, and show that imposing the Z strict mass limitsbasically rules out the sneutrino DM solutions in 5.3. We discuss possibilities for detection in Sec. 6 and conclude in Sec. 7. We leave some relevant formulas for the Appendix.

MODEL DESCRIPTION
In this section, we describe the supersymmetric model under investigation briefly. This model, based on SU (3) c × SU (2) L × U (1) R × U (1) B−L (thereafter referred to as the BLRSSM) was first introduced in [8] and further studied in [9,10,13]. The model emerges from breaking of supersymmetric SO(10) to the SM through the following intermediary steps, The advantages of this model are • It is obtained by breaking of SO(10) through a left-right symmetric model, thus inheriting some of its attractive features [4,14]; • It is able to explain neutrino masses by the inverse seesaw mechanism [8]; • It preserves gauge coupling unification of the MSSM, even when the breaking scale in the last step is of the order of the electroweak scale [9]; • It resolves the MSSM Higgs mass problem by yielding larger Higgs masses through additional D-terms in the soft-breaking potential, without resorting to heavy particles [9]; • It could yield signals differentiating it from MSSM, which may lie in different regions of SUSY parameter space; • It could provide different dark matter candidates and phenomenology, which in turn inform the study of direct and indirect searches.
The particle content of the model contains, in addition to the SM particles: 1. In the fermionic/matter sector, an additional (right-handed) neutrino N c i , required for anomaly cancellation, and an additional singlet fermion S, needed for generating neutrino masses. Both these fermions come in 3 families and are accompanied by their scalar partners; 2. In the bosonic/Higgs sector, two new Higgs fields, X R and X R , remnants of SU (2) R doublets, needed to break U (1) R × U (1) B−L → U (1) Y , and their fermionic partners; 3. In the gauge sector, an additional neutral gauge field, Z , which emerges from the mixing of the neutral gauge fields of SU (2) L , U (1) R and U (1) B−L , (W 0 , B R , B B−L ), and its fermionic partner.
In a sense, the model described here is minimal: however it requires an extra Z 2 matter parity to avoid breaking of R-parity [10].
The superpotential in this model is described by where the first line of Eq.(2.1) contains the usual terms of the MSSM, while the second line includes the additional interactions from the right-handed neutrino N c i and the singlet Higgs fields X R , X R with -1/2 and +1/2 B − L, and +1/2 and -1/2 R charges, respectively. The first term of the second line in superpotential describes the Yukawa interactions between neutrinos, and Y ij ν is the Yukawa coupling associated with these interactions. In a similar manner, Y ij s represents the Yukawa coupling among N c i , X R and S. Moreover, µ R is similar to µ term of the B − L extension of supersymmetric model (BLSSM) and stands for bilinear mixing between X R and X R fields. Note that there is also a µ S term to generate non-zero neutrino masses with inverse seesaw mechanism, and as customary, it is restricted to have a low value, as it cannot give important contributions to any other sector except for neutrinos. Contrary to BLSSM [15][16][17], where neutrinos have Majorana mass terms, N c i fields interact with X R and S through Y ij s N c i X R S term, and lead to SM-singlet pseudo-Dirac mass eigenstates. Besides, the interaction of the SU (2) L singlet Higgs fields X R , S and N c i yield a significant contribution to the masses of the extra Higgs bosons. Implementing the inverse seesaw mechanism into model allows Y ij ν and Y ij s to be at the order of unity. Hence, the contribution from the right-handed neutrino sector to the Higgs boson cannot be neglected and yields a different low scale phenomenology from MSSM and BLSSM with inverse seesaw mechanism [18][19][20].
The soft-breaking Lagrangian terms in the model are which contain triple scalar interactions, scalar masses and masses for the gauginos of all gauge groups, denoted by λ's.
while SU (2) L × U (1) Y is broken further to U (1) EM by the VEVs of the Higgs doublets The spectrum for this model, including particle masses, neutrino seesaw, mixing of gauge bosons and the neutralino sector has been discussed before [13], and we do not repeat it here. In what follows we concentrate on scanning the model parameters first by imposing Higgs sector, particle masses and other low energy restrictions, and then looking for dark matter candidates and resolution of the anomalous magnetic moment of the muon, thus restricting the parameter space to region where these conditions are satisfied.

SCANNING PROCEDURE AND EXPERIMENTAL CONSTRAINTS
We proceed to analyze the model by scanning the fundamental parameter space of BLRSSM. We use SPheno 3.3.3 package [21,22] obtained from the model implementation in Sarah 4.6.0 [23,24]. This package employs renormalization group equations (RGEs), modified by the inverse seesaw mechanism to evolve Yukawa and gauge couplings from M GUT to the weak scale, where M GUT is determined by the requirement of gauge coupling unification. We do not strictly enforce the solutions to unify at M GUT , since a few percent deviation is allowed due to unknown GUT-scale threshold corrections [25]. M GUT is thus dynamically determined by the requirement of gauge unification, that is g L = g R = g B−L ≈ g 3 , with subindices denoting the gauge couplings associated with SU (2) L , SU (2) R , U (1) B−L and SU (3) C respectively. With boundary conditions determined at M GUT , all the soft supersymmetry breaking (SSB) parameters along with the gauge and Yukawa couplings are evolved to the weak scale.
We performed random scans over the parameter space, as illustrated in Table I, imposing universal boundary conditions for scalar and gaugino masses. We comment briefly first on the parameters chosen, and then on the  m τ 1 105 GeV [33] m b 1 222 GeV [33] m q 1400 GeV [33] m τ 1 > 81 GeV [33] m e 1 > 107 GeV [33] m µ 1 > 94 GeV [33]   constraints included. Here m 0 corresponds the mass terms for all scalars, and M 1/2 represents the mass terms for all gauginos, including the ones associated with the U (1) B−L and U (1) R gauge groups. In setting the ranges for the free parameters, we scan scalar and gaugino SSB mass terms between 0-3 TeV, regions which yield sparticle masses at the low scale, especially the LSP.
Here A 0 is the trilinear scalar interaction coupling coefficient, and we adjusted its range to avoid charge and/or color breaking minima, which translates into |A 0 | 3m 0 [26,27]. Also, tan β is the ratio of vacuum expectation values of the MSSM Higgs doublets v u /v d , while tan β R which denotes the ratio of vacuum expectation values of v X R /v X R , is also free parameter in this model. Practically however, tan β R is required to be close to 1, in order to prevent large D-term contributions to the sfermion masses and to avoid tachyonic solutions. The VEV v R represents the vacuum expectation value which breaks extra U (1) B−L × U (1) R symmetry. Since the breaking scale of the extra symmetry plays a crucial role in determining the Z mass, the gauge boson associated with U (1) B−L × U (1) R symmetry, we scan v R between 6.5 and 20 TeV to obtain Z boson masses consistent with the current experimental bounds.
The parameter µ is the bilinear mixing of the MSSM doublet Higgs fields, while µ R is the bilinear mixing of the SU (2) R remnants Higgs fields, which are singlet under SU (2) L symmetry. The values of µ and µ R can be determined by the radiative electroweak symmetry breaking (REWSB) but their signs cannot; thus, only their signs remain as free parameters. Since the model contributions to muon anomalous magnetic moment are related to the sign of µM 1/2 , we scan over positive µ values, but we accept both negative and positive solutions of µ R , while requiring solutions consistent with experimental predictions, and favoring solutions which improve upon the SM predictions for the muon g − 2 factor. The superpotential of the model also includes a µ S parameter, which yields non-zero neutrino masses via the inverse seesaw mechanism. However, µ S is constrained to be small, so that it cannot effect any supersymmetric particle masses or decays. We also fixed the top quark mass to its central value (m t = 173.3 GeV) [28] in our scan. The Higgs boson mass is very sensitive to the top quark mass, and small changes in its value can shift Higgs boson mass by 1-2 GeV [29,30], although it does not significantly affect sparticle masses [31]. Hence, we scan both diag(Y ij ν ) and diag(Y ij s ) between 0.001-0.99, though the inverse seesaw mechanism prefers values of order 1. In scanning the parameter space, we use the interface which employs Metropolis-Hasting algorithm described in [39]. All collected data points satisfy the requirement of REWSB. After collecting the data, we impose current experimental mass bounds on all the sparticles and SM-like Higgs boson as highlighted in Table II. Although we restrict the SM-like Higgs boson to lie between 122-128 GeV with 3 GeV uncertainty, we also employed HiggsBounds 4.3.1 package [40] to compare our Higgs sector predictions with the experimental cross section limits from the LHC, and we require agreement with Higgs boson decay signal strengths at tree level, h → W W , h → ZZ and h → bb. Thus using the mass-centered χ 2 , and selecting the parametrization for the Higgs mass uncertainty as "box" we employed HiggsSignals 1.4.0 package [41] and bounded the solutions which yield total χ 2 (μ) 3. Another constraint comes from rare B-decay processes, B s → µ + µ − [34], b → sγ [36] and B u → τ ν τ [35]. The B-meson decay into a muon pairs, in particular, constrains the parameter space since there the SM predictions are consistent with the experimental measurements. The supersymmetric contributions are proportional to (tan β) 6 /m 4 Ai and constrained to be small. Hence, m Ai has to be heavy enough (m Ai ∼ TeV) to suppress the supersymmetric contributions for large tan β values. In addition to these limitations, dark matter observations severely restrict the parameter space, requiring the LSP to be stable and electric and color neutral, which excludes a significant portion of parameter space where stau is the LSP. We concentrate on two different data sets, one with the neutralino being the LSP, and one where sneutrino is the LSP, and we shall distinguish these two scenarios throughout our investigations. We employ micrOMEGAs 4.3.1 package [42] and tag the solutions which yield consistent relic density within 20% uncertainty range provided from WMAP data [37,38] as specified in Table II. Apart from relic abundance constraint, we do not impose any restriction from the dark matter experiments. All the experimental restrictions mentioned above are listed in Table II.

DARK MATTER PHENOMENOLOGY
For either neutralino or sneutrino to be viable candidates for dark matter, they must yield the correct level of relic abundance for thermal dark matter production in the early Universe, determined very precisely as the amount of non-baryonic dark matter in the energy-matter of the Universe, Ω DM h 2 = 0.1199 ± 0.0027 [43], with Ω DM being the energy density of the dark matter with respect to the critical energy density of the universe, and h the reduced Hubble parameter.
In addition, as the lack of any dark matter signals in either direct or indirect dark matter detection experiments confront our theoretical expectations, these must satisfy increasingly severe constraints from experiments. The interaction of dark matter with detector nuclear matter can be spin-dependent or spin-independent. The spin-dependent scattering can only happen for odd-numbered nucleons in the nucleus of the detector material, while in spin-independent (scalar) scattering, the coherent scattering of all the nucleons in the nucleus with the DM are added in phase. Consequently, in direct detection experiments, the experimental sensitivity to spin-independent (SI) scattering is much larger than the sensitivity to spin-dependent scattering, and thus we shall concentrate on the former.
We proceed as follows. First, we analyze the consequences of having the lightest neutralino as the dark matter candidate. Using the results in the previous sections, we explore the parameter space of the model which is consistent with this assumption. We follow in the next subsection with the parameter restrictions for sneutrino dark matter.

Neutralino Dark Matter
In this subsection, we concentrate on analyzing the consequences on the mass spectrum of the BLRSSM obtained by scanning over the parameter space given in Table I where lightest neutralino ( χ 0 1 ) is always the LSP, and highlight the solutions compatible with the constraints showed in Table II. We start with Fig. 1 which displays the allowed parameter regions, with plots in m 0 − M 1/2 , m 0 − A 0 /m 0 and M 1/2 − tan β planes. Throughout the graphs, all points satisfy REWSB. Blue points satisfy all experimental mass bounds, signal strengths of SM-like Higgs boson and rare B-decay constraints given in Table I. Red points obey the above mentioned constraints, as well as relic density bounds, 0.09 ≤ Ω DM h 2 ≤ 0.14. The m 0 − M 1/2 plane shows that solutions for M 1/2 800 GeV are excluded by the constraints in Table II, and the requirement of consistent relic density (red points) excludes a significant portion of the LHC allowed region (blue points). For  [16], where negative A 0 solutions for m 0 ≥ 1 TeV do not satisfy REWSB, here all LSP constraints can be fulfilled for this portion of parameter space, while only the relic density constraint imposes positivity of A 0 . The M 1/2 − tan β plot indicates that it is possible to find solutions with 0.09 ≤ Ω DM h 2 ≤ 0.14 only for large tan β values, 40 ≤ tan β ≤ 60, although it is easier to satisfy LHC limitations for low tan β values.
In Fig. 2, we show specific results for the determination of sparticle mass spectrum, with plots in The color coding is the same as Fig. 1. Furthermore, the mass region in which the two masses are degenerate is displayed as a solid green line. We note that the LSP neutralino solutions consistent with the relic density bound can be obtained only when 300 GeV  Table II. Red points form a subset of blue, and represent solutions consistent with the relic density constraint.
≤ m χ 0 1 ≤ 800 GeV. As can be seen from m t1 − m χ 0 1 and m b1 − m χ 0 1 planes, we find that stop and sbottom masses have to be at least ∼ 1.5 TeV and 2 TeV respectively to fulfill all the restrictions. Even though it is possible to find light stop solutions (m t1 ≤ 1 TeV) when 340 GeV ≤ m χ 0 1 ≤ 550 GeV, the relic density condition is not satisfied for these solutions. Moreover, unlike the results of BLSSM [16] where the lightest chargino masses are always above 600 GeV, here the m χ ± 1 − m χ 0 1 plot shows that there is a region of parameter space where lightest chargino solutions is nearly degenerate with the lightest neutralino when 300 GeV ≤ m χ 0 1 ≤ 500 GeV. These solutions correspond to the case where the lightest chargino decays into the neutralino LSP and W/W boson ( χ ± 1 → χ 0 1 + W ± (W ± )), and the branching ratio for this channel is almost 1. On the bottom right panel, the m τ1 − m χ 0 1 plane illustrates the stau mass along with the LSP neutralino mass. There is a parameter space around mχ0 ∼ 600 GeV, where stau mass is almost degenerate with the LSP neutralino and becomes the next to lightest supersymmetric particle (NLSP), but also for a region of the parameter space, the stau can be much heavier than the neutralino LSP. The lightest stau NLSP solutions compatible with the relic density constraint occur around 500 GeV. One can choose one of these solutions and study relevant neutralino annihilation processes mediated by a light stau [44].
The bottom plots in Fig. 2, show our results for the sparticle spectrum for the gluino and sneutrinos, with the plots in m q − m g , (where q represents squarks from the first two families), and m ν1 − m χ 0 1 planes. The m q − m g plane shows that squarks masses for the first two families and gluino masses should be heavier than 2 TeV but lighter than 4 TeV (light blue points). Although the relic density condition and the current ATLAS experimental limit [45] strictly constrain the crucial portion of the parameter space, most of the solutions are consistent with this experimental exclusion. Finally, the m ν1 − m χ 0 1 plane reveals that it is hard to find solutions with sneutrino as the supersymmetric NLSP if we require consistency with the relic density bound, and the lightest sneutrino solutions satisfying all bounds can be obtained at around 1 TeV.
Note that the graphs contain also information on the composition of the neutralino LSP. As can be seen from gluino vs squarks panel, light red points, which represent the mixed or higgsino-like neutralino LSP solutions consistent with the relic density bounds, are mostly found under the yellow curve (the excluded region). Light blue points representing mixtures of R-bino and B −L bino (gauginos of U (1) R and U (1) B−L , respectively) neutralino LSPs are mostly located within the 1 sigma error of the yellow line.
To continue the investigation of the neutralino LSP composition, in Fig. 3 we plot the correlation between the neutralino mass and gaugino and higgsino mass ratios with (top left) M 4 /M 1 , (top right) M 1 /µ, (bottom left) M 2 /µ, and (bottom right) µ R -µ, for correct relic density. The color coding is the same as Fig. 1. According to the M 4 /M 1 − m χ 0 1 plane, there must be a clear relation between the B − L bino B and B R masses so that the ratio of B R / B should be at around ∼ 1.8, decreasing slightly when the neutralino LSP mass increases. The next two plots compare the bino-higgsino (top right) and wino-higgsino (bottom left) masses, respectively, by looking at their mass ratio. In the top right plot, almost all solutions satisfying LHC collider bounds, and all solutions satisfying relic density constraints have M 1 /µ 1, that is the bino is lighter than the higgsino mass parameter. The left bottom plane shows that, despite allowing for light higgsinos, the wino is mostly lighter than the higgsino over all the parameter space where relic density bounds are satisfied. The µ R − µ plot (bottom right) shows that solutions prefer positive µ R to the negative ones, and µ R can take values in a large range between 500 GeV-7 TeV while the relic density bound can only be fulfilled with the low µ values. As can be seen from µ R − µ plane, the relic density constraint can be satisfied mostly when µ 0.5 TeV and 0.7 TeV µ 1.5 TeV.
The neutralino LSP content consistent with all constraints (including relic density) is as follows: its mass is constrained as 300 GeV m χ 0 1 500 GeV, and for those parameter points, the neutralino LSP content is a B R -ino, H-ino and B-ino mixture, in this region the wino masses are heavier than the higgsino masses for solutions consistent with the relic density bound. Since M 1 /µ 1, the bino mixes more than the higgsinos to form the LSP neutralino. In the region 500 GeV m χ 0 1 800 GeV, the LPS neutralino is about 60% B R − 40% B admixture, consistent also with the top left plot in Fig. 4. In Fig. 4 we present results specific to dark matter phenomenology, plotting the relic density and spin-independent cross section as a function of the lightest neutralino mass. In addition, we plot the correlation between the lightest pseudoscalar and the third lightest neutral Higgs boson h 3 , to highlight the fact that dark matter annihilation proceeds through these two funnels. We show (top left) Ω DM h 2 − m χ 0 1 , (top right) σ SI nucleon − m χ 0 1 , (bottom left) m A1 − m χ 0 1 , and (bottom right) m h3 − m A1 plots. In the top left and top right plane color coding is indicated in the insert, while for the bottom plots the color coding is the same as Fig. 3. The top left plot confirms our previous results on the content of LSP neutralino between 500 -800 GeV is composed of 60% B R -ino and 40% B-ino, whereas when 300 GeV m χ 0 1 500 GeV, its content is shared among B R -ino, H-ino and B-ino . The top left plot shows the dependence of the relic density, and the right plot shows the dependence of the spin-independent proton and neutron cross section, with neutralino LSP mass. The solid green line represents the current exclusion limit for XENON1T experiment [46]. As can be seen from the graph, most solutions consistent with the relic density constraint can be found below the XENON1T exclusion bound, specifically between 10 −10 pb -10 −11 pb. Hence they can be detected by the next generation DM detectors such as XENONnT [47], LZ and DARWIN [48]. Note that we also have a substantial amount of solutions consistent with the relic density above XENON1T exclusion limit. These solutions correspond to the region where 300 GeV with h 3 , with mass between 1 and 3 TeV. In this energy scale A 1 and h 3 provide the main funnel channels of this model. Apart from these, we have also verified the relation of the relic density with the IceCube confidence level exclusion and the neutrino flux, and all neutralino LSP solutions surviving relic and cross section bounds are within 1% confidence level of the experimental result.

Sneutrino Dark Matter
The BLRSSM contains, in addition to the three left sneutrinos, six additional singlet states, three right sneutrinos and three S, the scalar partners of S. The latter two provide candidates for sneutrino dark matter, as they do not suffer from too large an annihilation cross section (thus small relic density) from interacting through Z or W bosons. Sneutrinos thus provide alternative candidates for dark matter in this model, and we analyze their consequences in this subsection. In the left and right panels of Fig. 5 we show the dependence of the relic density Ω DM h 2 as a function of the lightest scalar neutrino mass. The color bars to the right of each plot indicate the right-handed sneutrino and the S content, respectively. As can be seen from the plot, even though it is possible to find sneutrino LSP solutions for almost all values of m ν1 between 0-1400 GeV, requiring consistency with the relic density bound constraints LSP sneutrinos to be between 200-400 GeV. Thus the indication would be that sneutrino LSP case allows lighter LSP masses compared to the neutralino LSP scenario. The right-handed content of the sneutrino LSP solutions changes between 45%-80%, while S composition varies between 20%-52%. Imposing relic density bounds, the mixed sneutrino LSP is about 50-50 % between right-handed and S. Thus the scalar partner of S, introduced for neutrino seesaw, plays a crucial role in the sneutrino LSP composition.
In Fig. 6 we analyze the dependence of the nucleon spin-independent cross section, σ SI p for both the proton (left panel) and neutron (right panel). The color coding is the same as Fig. 1 and also indicated in the legend of the plots. The plots show the relation for the spin independent cross section for proton and neutron respectively. We note that both dark matter constraints (the relic density and σ SI p ) severely restrict the parameter space where the sneutrino is the LSP in this model.

COMPARISON OF THE TWO DARK MATTER SCENARIOS
In the previous section, we analyzed DM phenomenology for both neutralino LSP and sneutrino LSP scenarios in BLRSSM. As discussed in detail, BLRSSM provides quite different mass spectrum for two distinct variants of LSP, and these relatively two different mass spectra change the low scale DM phenomenology in important manner. While we found sneutrino LSP scenario to be highly constrained and statistically unlikely, there are a few parameter points that survive universal boundary conditions, so in this section, we compare results for the two different LSP scenarios. In Fig. 7 we plot in the µ − µ R and tan β − M 2 /µ dependence. Dark blue points satisfy the mass bounds and constraints from the rare B-decays for the neutralino LSP solutions. Red points form a subset of dark blue, and represent neutralino LSP solutions which satisfy the relic density constraint. Light blue solutions are consistent with the mass bounds and the constraints from the rare B-decays for sneutrino LSP solutions, while yellow points form a subset of light blue, and represent sneutrino LSP solutions consistent with the relic density constraint.
The µ − µ R plots compare the higgsino sectors of our model. We note that while the neutralino LSP solution can allow values of µ R between 7-9 TeV, sneutrino LSP solutions prefer low µ R values, mainly between 0-4 TeV for positive µ R . Even this range becomes narrow, around 1.5 TeV, for lighter higgsinos. For the sneutrino LSP solutions, µ < 1.5 TeV, and µ R values favor the region between 4-7 TeV. On the right panel, the tan β − M 2 /µ plane shows the relative wino and higgsino mass ranges for the two LSP scenarios.
From the plots, we conclude that for sneutrino LSP, M 2 /µ 1 and the wino is always lighter than the higgsino over all the parameter space. For the neutralino LSP case, the higgsinos can be lighter or heavier than winos. Also, tan β values for sneutrino LSP solutions are found anywhere in the 0-50 range, and solutions consistent with the relic density constraint can be obtained for either M 2 /µ 1 or M 2 /µ 1. Requiring consistency with the relic density bound solutions with M 2 /µ 1 correspond to neutralino LSP, and tan β values lie in the 10-50 range. Requiring compatibility with the relic density bound, further constrains the region M 2 /µ 1 to correspond to B − B R dominated neutralino LSP solution, where tan β should be between 40-60. In general the model clearly favors solutions with neutralino LSP to those with sneutrino LSP.

The neutral Higgs sector
The choice of LSP affects the heavier states in the Higgs Sector of BLRMSS. For both neutralino and sneutrino LSP solutions, the lightest neutral Higgs boson can be lighter than 150 GeV. Fig. 8 shows the results for the values of Higgs masses for both LSP cases with plots for m h2 relative to m h1 (left) and m A1 dependence of tan β (right), where A 1 is the lightest pseudoscalar. The color coding is described in the legend of these planes. The left plot shows that while the two lightest neutral Higgs bosons can be degenerate when the LSP is neutralino, degenerate solutions cannot be obtained for the sneutrino LSP, where the second lightest Higgs boson mass is between 150 -700 GeV. This phenomenon can be explained as due the contributions obtained from different elements of CP-even Higgs mass matrix. When m h2 > 150 GeV, the dominant contribution comes from the m 2 RR element of CP-even Higgs mass matrix, corresponding to singlet Higgs fields associated with U (1) R × U (1) B−L . Thus there h 2 is mostly a singlet Higgs boson. The off-diagonal term m 2 LR which provides essential mixing between the two sectors becomes important when m h2 < 150 GeV. For the sneutrino LSP solutions, the Yukawa coupling Y s is constrained to be small (as the sneutrino LSP mass is generated mostly through this term), unlike when the LSP is the neutralino. The Y s coupling then imposes lighter h 2 masses, mostly generated by the singlet Higgs field X R . The other Higgs bosons can be quite heavy. This is seen also in the right-hand side of Fig. 8, where we plot the dependence of the mass of the lightest pseudoscalar Higgs boson A 1 (degenerate with h 3 ), with tan β. As before, the region in tan β ∼ 40-60 represents the mixed binos neutralino LSP solutions, while for tan β < 40, regions with larger (smaller) A 1 mass correspond to sneutrino (neutralino) LSP. Thus second lightest Higgs boson is a singlet in both scenarios, but, while the sneutrino LSP scenario favors the 150-700 GeV mass range, for the neutralino LSP solutions the second Higgs mass can be much heavier than 700 GeV.

The muon anomalous magnetic moment
The experimental results for the muon anomalous magnetic moment pioneered by the BNL E821 experiment [49,50] have been improved with the updated results from FNAL E989 [51] and J-PARC E34 [52] experiments. However, the SM prediction for the muon anomalous magnetic moment [53], a µ = (g − 2) µ /2, indicates a 3.5σ deviation from the experimental results, The SM prediction is limited in precision by the evaluation of hadronic vacuum polarization contributions. Calculations exist for the lowest contributions, evaluated using perturbative QCD and experimental cross section data involving e + e − annihilation into hadrons. However, the large discrepancy has motivated possible explanations within new physics scenarios.
In MSSM, if one of the smuons and bino or wino soft masses can be sufficiently light, supersymmetry can ameliorate this discrepancy. However, if the model is required to obey universality conditions at M GUT , obtaining the correct Higgs boson mass is the greatest challenge to explaining the muon g − 2 anomaly. We can expect better results from the BLRSSM model since it includes inverse seesaw mechanism and an extra gauge sector. The effect of inverse seesaw mechanism can be read through RGE for the smuons. As can be seen from the last two terms of Eq. A7, the Yukawa coupling Y ν helps decrease the smuon masses at low scales, as compared to models without inverse seesaw. A similar effect can be read through the RGE of µ Eq. A1 and sneutrinos Eq. A6. The presence of another free Yukawa coupling Y s in addition to Y ν contributes to evolving light sneutrino masses to the low scale via RGE as can be seen from the Eq. A6.
Here we investigate the effects on the muon g − 2 anomaly for both sneutrino and neutralino LSP cases. Fig. 9 displays the correlations between muon a µ and the relevant free parameters in m 0 , M 1/2 , tan β and µ. The color coding is the same as Fig. 8. In addition, the shadowed regions show 1, 2 and 3 σ deviations between the calculated contribution to muon g − 2 factor and its experimental value. The top left side plot shows that ∆a µ favors low values for m 0 (light scalar masses). Similarly light gaugino masses (light electroweakinos) are also required to decrease the ∆a µ discrepancy, as seen from the top right handed plot. The need of light scalars and electroweakinos agrees with large tan β values (bottom left panel). Finally, the ∆a µ depends sensitively on the µ parameter, as in MSSM, and here the contribution to the muon g − 2 factor drops sharply for µ > 1.5 TeV. This is due to one loop contributions effects, where, as the µ term increases, the contributions where the higgsinos run in the loop are suppressed, while bino-smuon loop is left as only effective contributing diagram. However, as the bino masses cannot be as low asB R masses, the contribution from this channel is insufficient. And thus, against expectations, the inverse seesaw mechanism cannot sufficiently enhance muon ∆a µ within universality conditions, and the corrections hardly reach the 2σ region.

Z mass constraints
To highlight the differences between the two scenarios, we kept the model as general as possible and did not impose Z mass bounds so far. We include here an investigation of implications of the constraints imposed on the Z mass by a recent new study at ATLAS [11], requiring an increase in the lower bound for the BLRSSM model to M Z > 3.9 (3.6) TeV in the ee(µµ) channels. One must be careful when applying these bounds. First, the experiment assumes nonsupersymmetric models, and thus a case where Z does not decay to supersymmetric particles, which will modify its total decay width and thus branching ratios. Second, the parameter choice and unification scale is different from ours: the choice depends on symmetry breaking scales and assumed multiplet composition of the GUT parent. With this note of caution, we explore the parameter space here.
First, we show some of the decay rates of the Z boson in BLRSSM. Fig. 10 displays some of the important decay channels of Z where BR(Z → ll) = BR(Z → ee) + BR(Z → µµ), BR(Z → l l), BR(Z → qq) and BR(Z → χ χ), all plots as a function of m Z . Throughout, all points are consistent with LHC, B-physics bounds, HiggsBounds and HiggsSignals. Dark blue points show neutralino LSP solutions whereas light blue ones stand for sneutrino LSP solutions.
The top left panel in Fig. 10 exhibits the branching ratio of Z into lepton pairs while the top right panel shows the branching for the supersymmetric partners in the same channel. As can be seen from top left plane, the branching ratio of Z into leptons changes between 25% − 37% while its decays into their supersymmetric partners, sleptons, are low, in the range of 0% and 8%. It is interesting to note that these models, unlike E 6 -derived models containing an extra U (1) gauge group, are not likely to be leptophobic as the branching ratio into leptons is significant throughout the parameter space investigated. The bottom panels of Fig. 10 show the branching ratio into quarks (left) and into neutralinos and/or charginos (right). As usual, the largest branching ratio obtained is hadronic (40%-62%), which,  Fig. 8. In addition, the shadowed regions show 1 σ, 2 σ and 3 σ differences between the theoretical contribution to muon g − 2 factor and its experimental value. though significant, is not as large as for U (1) models [54], which will likely adversely affect the Z production cross section. The decay into two charginos or neutralinos occurs above their mass threshold and is very small throughout the whole parameter space (0%-13%). So it appears that the decay of the Z boson is fairly consistent with a non-supersymmetric scenario. Based on this, we shall investigate the effects of setting the mass lower bound to be m Z > 3.5 TeV throughout our analyses.
With these constraints, we revisit the plots for the spin independent cross section for proton and neutron respectively. While in the Fig. 6, we considered m Z ≥ 2.5 TeV, and the spin-independent proton (or neutron) cross sections for sneutrino LSP solutions were satisfied with XENON1T experimental exclusion limit, imposing the the new Z mass limit excludes most of the parameter space for sneutrino LSP solutions, as shown in Fig. 11. Specifically, of about 10 6 scanned parameter points only 18 solutions compatible with the relic density bound are found, and only 10 of them can survive XENON1T experimental exclusion limit. Imposing Z mass constraints, the sneutrino LSP scenario thus emerges as extremely constrained and, realistically, ruled out.

COLLIDER SIGNALS
Lastly, we would like to analyze the production and decays for this scenario at the LHC. We choose benchmarks from the parameter scan results which satisfy all experimental bounds, including the relic density constraint and XENON1T exclusion limits, and favor light neutralino LSP solutions as the only ones surviving all constraints. We proceed by exporting the BLRSSM to the UFO format [55] and use MG5 aMC@NLO framework version 2.5.5 [56] to simulate hard-scattering LHC collisions and evaluate the cross sections for various signals. For the calculation of cross sections, we select four benchmarks with different features, which could showcase different features of the model  11: Dependence of the spin independent cross section for the proton σ SI p (left) and neutron σ SI n (right) as a function on the sneutrino LSP mass m ν 1 , for , m Z ≥ 3.5 TeV. All points are consistent with REWSB and sneutrino LSP. The color coding in each plane is the same as Fig. 1. for detection at the LHC.
The first benchmark, benchmark 1 has H R -like neutralino LSP. (Even though parameter scans allow Higgsino-like and higgsino-binos mixed LSP neutralino solutions between 300-500 GeV, no benchmark in this range can be found as these states are completely excluded by the XENON1T exclusion limit.) We thus select benchmarks with mixed B R − B content. For benchmarks 2-3, BR( χ 0 2 → χ 0 1 h 1 ) and BR( χ ± 1 → χ 0 1 W ± ) are almost unity. Sparticle masses are similar in both cases, with the exception of the lightest charging, which is heavier for benchmark 3. Also, for benchmark 3, BR( τ 1 → τ 1 χ 0 1 ) ∼ 1 while this is much smaller for benchmark 2. Benchmark 4 is selected for light stau masses, leading to increased stau-stau production cross sections. Note that benchmarks satisfy all the constraints, including the Icecube22 exclusion. Our results are shown in Table Sec. (6).
Even though LSP neutralino mass is quite light (67 GeV) for benchmark 1, we find that both chargino-chargino and neutralino-chargino production cross sections are quite low, due to the fact that the neutrino is mostly higgsino. For the other benchmarks, with neutralino contents of mixed binos, the second lightest neutralino and chargino masses are degenerate. We estimated the cross sections for chargino/neutralinos and stau production as being the most promising. The highest cross-section values for chargino-chargino production and chargino-neutralino production are obtained for benchmark 2 whose neutralino and chargino masses are 470 GeV and 767 GeV, respectively. As can be seen from the Sec. (6), chargino-chargino production and neutralino-chargino production cross sections are 4.623 fb and 2.249 fb, respectively. The cross-section values decrease in benchmark 3 (with respect to benchmark 2) when neutralino and chargino masses are 506 GeV and 954 GeV (versus 470 and 767 GeV), respectively. Finally, the last benchmark is selected to enhance σ(pp → τ 1 τ 1 ) where each stau can decay into a tau and a LSP neutralino ( τ 1 → τ 1 χ 0 1 ). Note that here the stau is NLSP, and that the branching ratio of stau into a tau and a LSP neutralino is 1. The relevant cross-section is 0.5713 fb, a factor of 10 larger than for benchmarks 2 and 3, but still too small to be observed at the LHC. For all benchmarks, Z masses are above 4 TeV, consistent with the latest ATLAS result. Note that gluino masses are about 2.5 TeV for benchmarks 2, 3 and 4, making gluino results testable at the HL-LHC or by the next generation colliders [57,58].
Including all the constraints, we conclude that production of supersymmetric particles in BLRSSM fall below detector sensitivity. Especially because the final signals will have even lower production cross sections, as they will be suppressed by branching ratios of chargino/neutralinos to missing energy + leptons. A way to improve our results is to relax some or most universality constraints, and looking for effective cuts which would enhance the signal over the background. We shall return to this in a future work.

SUMMARY AND CONCLUSION
We analyzed the predictions of the mass spectrum in the BLRSSM framework with universal boundary condition, highlighting the solutions consistent with the DM restrictions (relic density and spin independent cross sections with nucleons) for both neutralino and sneutrino LSP scenarios. We found that the stop and sbottom masses are between 2-3 TeV, and the chargino can be degenerate with the LSP neutralino between 300-500 GeV. In addition, the relic density constraint can be satisfied for masses in the range 300 m χ 0 1 800 GeV. H R dominated LSP neutralino solution can be obtained below 300 GeV, however these solutions are ruled out by the XENON1T spin independent cross section exclusion curve. When all DM constraints are taken into account, the model favours LSP neutralinos with masses between 500 m χ 0 1 800 GeV, bino-dominated, and with composition 60% B R -40%B. We also that, when the LSP is neutralino, A 1 and h 3 are funnel channels for pair-producing them.
In addition, the model allows in principle a sneutrino LSP where its content can be either right-handed dominated or mixed, ν R and S, with masses between 250-1300 GeV. In this sense, sneutrino LSP solutions can be lighter than the neutralino LSP ones. Purely right-handed dominated sneutrino LSP solutions have difficulty to satisfy the relic density constraint and only mixed ones survive. In addition most of the sneutrino LSP solutions are consistent with the XENON 1T spin independent cross section exclusion curve. However, strict imposition of the Z mass bounds basically rule out the sneutrino solutions, while not having any effect on the neutralino LSP parameter space. This is one of the most important predictions of the model.
The parameter spaces corresponding to neutralino and sneutrino are quite different. If allowed, sneutrino LSP solutions favor low singlet higgsino mass parameter, µ R , and the second lightest neutral Higgs boson a singlet, while neutralino LSP favor larger µ R parameters. Sneutrino LSP solutions are spread out over the whole range of tan β, while neutralino solutions are restricted in the 40 tan β 60. Neutralino LSP solutions allow for degenerate masses of the two lightest neutral Higgs bosons, while the sneutrino LSP, although favoring a light m h2 , does not. The anomaly in the anomalous magnetic moment of the muon favors neutralino LSP contributions, where for a large range of scalar masses, and a more restricted one for gauginos and higgsinos, the corrections are within 2σ of the experimental result, while sneutrino LSP solutions can at best produce results within 3σ of the desired values.  We analyzed collider signatures of this proposed scenario, including all constraints, and they are not promising. The largest cross sections are obtained for chargino/neutralino production, and they are at most of O(4) fb, without including cascade decays into leptons which would reduce them further. In the future, collider signals could be enhanced by relaxing some of the severe constraints on the model, such as the universality conditions, and finding suitable cuts to enhance signal versus background. This may extend the parameter space, allowing neutrino LSP back into the consideration. Work in these directions is underway.