Roads for Right-handed Neutrino Dark Matter: Fast Expansion, Standard Freeze-out, and Early Matter Domination

Right-handed neutrinos appear in several extensions beyond the Standard Model, specially in connection to neutrino masses. Motivated by this, we present a model of right-handed neutrino dark matter that interacts with Standard Model particles through a new gauge symmetry as well as via mass mixing between the new vector field and the Z boson, and investigate different production mechanisms. We derive the dark matter relic density when the Hubble rate is faster than usual, when dark matter decouples in a matter domination epoch, and when it decouples in a radiation domination regime, which is then followed by a matter domination era. The direct detection rate features a spin-independent but velocity suppressed operators, as well as a spin-dependent operator when the mass mixing is correctly accounted for. We put all these results into perspective with existing flavor physics, atomic parity violation, and collider bounds. Lastly, we outline the region of parameter space in which a weak scale right-handed neutrino dark matter stands as a viable dark matter candidate.


I. INTRODUCTION
The presence of dark matter in our universe is ascertained through a variety of dataset, going from dwarf galaxy scales to cosmological scales.The precision acquired by Cosmic Microwave Background (CMB) probes allowed us to quantify precisely the abundance of dark matter in our Universe.It accounts for nearly 27% of the total energy density [1].Its nature is unknown, though.Massive particles that feature weak interactions with Standard Model (SM) particles dubbed (WIMPs) have stood out from the crowd, driving most experimental efforts for many years [2].Currently, other plausible dark matter candidates such as the axion and axion-like particles, which have been always theoretically compelling, gained lots of interest from the community due to the non-positive signals in WIMP searches [3].However, no solid positive signals have been observed, favoring other dark matter candidates either.The customary assumption that dark matter had a thermal production have given us a misapprehension that a positive signal should appear in the current generation of experiments.There are ways to successfully have a thermal dark matter candidate while yielding no signal at current direct detection experiments.Hence, this sudden disbelief for the WIMP miracle does not seem justified.It has always been worthwhile to explore well-motivated dark matter production mechanisms that lead to attainable detection, specially if the mechanism still leaves imprints at current experiments.
Besides the dark matter siege, neutrino masses stand as a concrete evidence for physics beyond the Standard Model.With the observation of neutrino oscillations, we have concluded that at least two neutrinos are massive [4][5][6][7][8][9][10][11].It should be noted that the sum of the neutrino masses is constrained by CMB observations because the energy density of neutrinos which depends on the sum of neutrino masses affects both the relativistic and non-relativistic energy density relevant to the derivation of the CMB temperature and polarization power spectra [12].In fact, CMB probes give rise to the most stringent bound that reads Σm ν < 0.12 eV [13].Although, the individual masses are unknown.
An elegant solution to theoretically generate neutrino masses is offered by right-handed neutrinos.In particular, if copies of right-handed neutrinos are added to the Standard Model, one can naturally address neutrino masses through a type I seesaw mechanism when both Dirac and Majorana mass terms are included in the Yukawa Lagrangian.The Majorana mass is expected to be very large to play the seesaw and bring the active neutrino masses down below the eV scale.
That said, one may wonder if the mechanisms behind the presence of dark matter in our universe is somehow connected to neutrino masses.This question has driven a multitude of works.In particular, we will be interested in connecting these two new physics landmarks in the context of the type I seesaw mechanism [14][15][16][17][18][19][20][21][22][23][24][25][26][27][28], which features the Lagrangian,

arXiv:2108.11398v2 [hep-ph] 30 Aug 2024
where Φ is the SM Higgs doublet, and M the Majorana mass term for the right-handed neutrinos.Note that this Majorana mass is a 3 × 3 matrix.After spontaneous electroweak symmetry breaking, the first term yields a Dirac mass term for the active neutrinos.The presence of a Majorana mass allow us to employ the so-called type I seesaw mechanism that results into m ν ≃ m T D M −1 m D and m N ∼ M .For Yukawa couplings of order one in Eq. ( 1), this simple and plausible framework requires the right-handed neutrino masses to be very large, above 10 12 GeV.Consequently, this mechanism is hardly testable by existing experiments.It would be interesting if one could successfully explain neutrino masses and dark matter without invoking a very large energy scale suppression.A possible route is to consider light right-handed neutrinos with masses around 1 keV.If one of the three right-handed neutrinos is lighter than the others, then it can be made cosmologically stable as long as the mixing between active and right-handed neutrinos is sufficiently small.This scenario is dubbed sterile neutrino dark matter 1 [30,31].The dark matter relic density is found through a non-resonant production [30], which is based on the active-sterile neutrino oscillations.The approximate relic density is [32][33][34][35][36], where sin2 θ i ∼ a y 2 ab v 2 /M 2 is the active-sterile neutrino mixing and v is the SM Higgs vacuum expectation value.As mentioned earlier, the sterile neutrino is unstable because of its mixing with the active neutrinos.In fact, this mixing leads to the decays N → ννν, and N → νγ decay [37,38].The latter decay has been extensively looked for in X-rays surveys.However, such X-ray searches combined with structure formation requirements rule out this sterile neutrino production [39,40].There are ways to bypass this exclusion by evoking other production mechanisms [41,42].Instead of dwelling on light sterile neutrinos as dark matter, we investigate under which circumstance we can accommodate the type I seesaw mechanism, while hosting a TeV scale right-handed neutrino dark matter.Note, we have renamed sterile neutrino to righthanded neutrino.Theoretically speaking, they are the same field, but when ones refers to sterile neutrino typically small masses and suppressed mixing angles are adopted.When this sterile neutrino is heavy, with masses above 1 MeV, the term right-handed neutrino is used more often 2 .
Anyway, how can we have a weak scale right-handed neutrino dark matter having the presence of the first term in Eq. 1, which is necessary for the type-I seesaw to work but make it unstable?The idea is to invoke a Z 2 symmetry, where only one of the right-handed neutrinos is odd under it [45][46][47][48].In this way, two massive right-handed neutrinos will play a role in the seesaw mechanism, whereas the lightest becomes a plausible dark matter candidate.Notice that in general, one needs to either invoke fine-tuning in the yukawa couplings to make the right-handed neutrino a viable dark matter candidate or symmetries.We focus on the latter.Similar works have been carried out in this direction [49][50][51].Although, our work extends previous studies because we go beyond the standard freeze-out case.Concretely, we investigate the right-handed neutrino dark matter scenario under three different cosmological histories, namely fast expansion, matter-domination during freeze-out, and matter-domination after freeze-out.These cosmological histories change the theoretical predictions for the dark matter relic density compared to the standard freezeout case, and consequently lead to different conclusions.Besides computing the dark matter relic density under different cosmological scenarios, we also compute the dark matter scattering cross-section.We highlight that for right-handed neutrino dark matter, the dark matter scattering cross-section does not fall into the standard classification of spin-independent or spin-dependent dark matter-nucleon scattering.Therefore, we need to compare the resulting scattering rate with the data to derive the correct limits.Lastly, we also include collider bounds on the mediator masses based on the LHC search for heavy dilepton resonances.

II. THE MODEL
2HDMs have been extensively studied in the literature [81][82][83] for naturally keeping the ρ parameter unchanged despite the extended scalar sector.Furthermore, it leads to interesting collider searches that are within reach of the LHC.Although, the canonical version of the 2HDM suffers from flavor changing neutral interactions and does not address neutrino masses.It would be appealing if one could solve both issues via gauge symmetries.This is precisely the idea behind the 2HDM-U(1).An Abelian gauge symmetry is incorporated to the 2HDM in such a way that only one scalar doublet contribute to fermion masses, and three right-handed neutrinos are required to cancel the gauge anomalies and consequently play a role in the type I seesaw mechanism aforementioned.We will choose the new Abelian gauge symmetry to be the Baryon-Lepton (B-L) number, for concreteness, but we highlight that other gauge symmetries are also conceivable.The purpose of our work is to incorporate dark matter without adding new fields and explore the impact of our conclusions under different cosmological histories.Concerning the model itself, which is based on the SU (3) c × SU (2) L × U (1) Y × U (1) B−L , the fermion content features, and, Notice that thus far the only difference to the SM is the presence of right-handed neutrinos.The model has three scalar fields, and where the subscript a = 1, 2, 3 accounts for the three generations.We stress that the presence of the scalar singlet is twofold: (i) it breaks the U(1) gauge symmetry generating mass to the Z ′ gauge boson; (ii) it gives rise to the Majorana mass for the right-handed neutrinos.In the original version of the type I seesaw, we have a Majorana bare mass term for the right-handed neutrino masses, but in our work it is generated through a mechanism of spontaneous symmetry breaking.This fact will be clear below.Fermions get mass through the Yukawa Lagrangian L Y = L Y1 + L Y2 , where, Having a weak scale right-handed neutrino from Eq. ( 6) requires the addition of an extra discrete symmetry.We will invoke a Z 2 symmetry to stabilize, say N 1R , the lightest righthanded neutrino.In this way, the first term in Eq. ( 6) involves only two right-handed neutrinos, whereas the latter remains unaltered as long as y M ij is diagonal.The two right-handed neutrinos in the Dirac mass term act the type I seesaw mechanism leading to two massive active neutrinos and a massless one, in agreement with neutrino oscillation observations [84].This is a well-known fact which can be directly seen using the Casas-Ibarra parametrization where the neutrino masses are directly tied to the Yukawa couplings.As we are basically removing one Yukawa coupling, the first term of Eq.( 6) which is 6x6 matrix, should have at least four non-zero entries.Using the Casas&Ibarra parametrization one can see that the Yukawa matrix will depend on two right-handed neutrinos masses and two active neutrinos masses [85].We emphasize we have assumed that the Majorana mass matrix is purely diagonal.The scalar potential, in agreement with the gauge charges of the doublets, including the singlet scalar, is given by, The neutral components of all scalars acquire vacuum expectation values (vev) as and which give mass to all fermions and gauge bosons.It is important to emphasize that the dark matter candidate N 1R becomes massive through the singlet scalar vev, v s , Concerning the gauge sector, kinetic mixing arises at tree level, with ϵ, the kinetic mixing parameter, θ W , the Weinberg angle.The Lagrangian responsible for gauge bosons' masses is given by where the covariant derivative is, with Q Y and Q X the charges under Yukawa and B − L gauge symmetries, respectively, and g and g ′ the usual gauge couplings associated with SU (2) L and U (1) Y symmetries, respectively, and g X , the B − L coupling.Diagonalizing the matrices for X µ and Z µ bosons, for details we recommend [71], we get the following gauge boson masses where GeV 2 and g Z = g/ cos θ W .The Lagrangian in Eq. ( 11) is also responsible for the interactions between gauge bosons and scalars, for example which are also relevant for the dark matter phenomenology.These couplings can be found in [71].We emphasize that we have taken the kinetic mixing to be zero at tree-level.At one-loop level, it can be safely ignored, as our phenomenology will be driven by the gauge coupling g X .
The dark matter phenomenology of the model is governed by the neutral current involving the Z and Z ′ gauge bosons that reads [49,86], where Q X f = −1 for charged leptons, Q X f = 1/3 for quarks, and with, and, Identifying N 1R as our dark matter candidate, we can straightforwardly carry out the dark matter phenomenology within the standard freeze-out and ΛCDM model.This scenario was investigated somewhere else [71].Our plan is to go beyond and explore other cosmological scenarios.To do so, we start reviewing the important ingredients of early dark energy-like, early radiation, and early matter domination epochs.We have now set the particle physics model, we will review the cosmological background.

III. NON-STANDARD COSMOLOGY HISTORIES
In this section, we explore how the dark matter relic density may be affected by a non-standard cosmological history, where the right-handed neutrino dark matter candidate decouples from thermal equilibrium in different scenarios.The key aspect of different cosmology histories is the effect on the Hubble expansion rate, because it directly impacts the Boltzmann equations.The Hubble rate evolves differently for matter, radiation and dark energy-like fields.Thus, depending on which component dominates H, distinct solutions to the dark matter relic abundance are found.
The standard cosmology predicts that the Universe is radiation dominated at early times.However, before the Big Bang Nucleosynthesis (BBN), we have freedom to choose different cosmological histories for short period of times.For example, the expansion rate may be governed by matter or a new field, which leads to a faster expansion rate.That said, we tease out the dark matter phenomenology when the dark matter particles decouples during: (i) an expansion rate faster than radiation (Sec.III A); (ii) a radiation-dominated era (ERD) followed by a matter-dominated period (Sec.III B); and (iii) in a matter-dominated phase (EMD) (Sec.III C).We start reviewing the key ingredients of a faster expansion episode and how it changes the Boltzmann equation.

A. Faster Than Usual Early Expansion
Assuming that a new scalar field ϕ has an energy density that grows with the scale factor as, where n encodes the non-standard cosmological evolution.Notice that for n = 0, we recover the radiation energy density ρ R (t) ∝ a(t) −4 that corresponds to the standard case.Then, the scalar field energy density ρ ϕ dominates over radiation at early times and redshifts faster than the usual as the Universe expands.
Let us start by assuming that at some period of the early Universe, the Hubble expansion rate was driven by radiation and ϕ energy densities.Thus, the total energy density can be written as Nevertheless, the radiation energy density ρ R must overcome ρ ϕ at a time before BBN to avoid any tension with observations [87].Then, we define the temperature T r at which ρ R (T r ) = ρ ϕ (T r ).To be consistent with BBN observations, we impose, When the scalar energy density ρ ϕ is negligible, the radiation energy density evolves with temperature as usual with, where g ⋆ (T ) is the effective number of relativistic degrees of freedom at temperature T .
The next step is to express the energy density of ϕ in terms of temperature.Firstly, it is important to remember that the entropy density is given by, where g ⋆ s (T ) is the effective degrees of freedom of the SM entropy density.Using entropy conservation per comoving volume, sa 3 = const., we get, Inserting this result into Eq.( 18), and taking the ratio ρ ϕ (T )/ρ ϕ (T r ), the ϕ energy density can be written as a function of temperature, Finally, using the definition ρ R (T r ) = ρ ϕ (T r ), the total energy density at temperature T > T r becomes, From this equation, we conclude that ϕ dominates the Universe for temperatures T > ∼ T r .Once we know the evolution of the energy density, we can determine the Hubble rate using the Friedmann equation [88], where M P l = (8πG) −1/2 = 2.4 × 10 18 GeV.Assuming that g ⋆ (T ) = g ⋆ is a constant for temperatures T ≫ T r , i.e, for the period in which ρ ϕ completely dominates over ρ R , the Hubble rate becomes [87], with g ⋆ = 106.75accounting for the entire SM degrees of freedom.Knowing the Hubble rate in the ΛCDM radiationdominated cosmology [89], we can rewrite Eq. ( 27) as, where we can explicitly see that for temperatures T > T r non-standard cosmologies (n > 0) make the Universe expand faster than in the standard radiation-dominated epoch.Setting T f to be the freeze-out temperature, having a dark matter particle freezing-out in a fast expanding universe, requires T f > T f rad .Note that from our reasoning for temperatures larger than T r , a non-standard cosmology is at play.With the Eq. ( 27) at hand, we can compute the dark matter relic abundance for the right-handed neutrino N 1R when the universe expands too fast.In the standard radiation-dominated scenario, the Boltzmann equation that describes the evolution of the comoving dark matter number density Y N ≡ n N /s is, with x ≡ m N /T and ⟨σv⟩ being the thermally averaged annihilation cross-section, and H(x) the standard radiationdominated Hubble rate.The comoving equilibrium number density can be written as a function of x for a Maxwell-Boltzmann distribution as follows [87], where g N accounts for the degrees of freedom of the righthanded neutrino dark matter and and then from Eq. ( 30) we find, where ⋆ M P l m N .Approximate analytical solutions are given by [87], Eqs. ( 34) are valid for s-wave annihilation processes.In other words, when ⟨σv⟩ does not depend on temperature.Thus, these approximate analytical solutions can only be extrapolated up to x = x r .In Fig. 1, we present the numerical result for the yield, assuming the dark matter mass m N = 100 GeV and the annihilation cross-section ⟨σv⟩ = 10 −9 GeV −2 , for the following values of T r = 0.1 GeV (brown line), T r = 1 GeV (blue line), T r = 10 GeV (red line), and n = 4.For completeness, we also include the yield for the standard case (black line).
The n > 0 cosmologies can produce a higher comoving dark matter number density at the time of freeze-out x f < x r .However, for s-wave annihilation cross-section and n ≥ 2, the annihilation rate after freeze-out scales as Γ ∝ T 3 , whereas For completeness, we also include the yield for the standard case (black line).Note that for Tr = 10 GeV, we have a superposition with the standard case.The shaded gray region is the approximate range for the correct dark matter relic abundance as measured by Planck [13].It is clear that the longer the period of ϕ dominance, the larger the impact on the yield.FIG. 2. Illustration of the density energy evolution and thermal history for faster early expansion than usual cosmologies.
H ∝ T 2+n/2 .Hence, the annihilation rate redshifts either slower (or equally) than the Hubble rate.Consequently, after freeze-out the dark matter keeps annihilating until the temperature T r , i.e. until radiation starts to govern the universe expansion.In this epoch, H R ∝ T 2 , meaning that the dark matter will stop annihilating, leading to a constant energy density.
In Fig. 2, we briefly illustrate the thermal history for the n ≥ 2 cosmologies.The dark matter decouples at T f whilst the scalar field dominates.The period T r < T < T f defines the relentless phase in which the dark matter tries relentlessly to go back to thermal equilibrium, unsuccessfully [87].For n = 2, Eq. ( 34) translates the relentless effect via the slow logarithmic decrease of the comoving dark matter number density.Although for n > 2, the comoving dark matter number density decreases following a power law.
Finally, solving Eq. ( 33) numerically, we compute the dark matter relic density taking x → ∞ via the following expres-sion, where s 0 stands for the SM entropy density today, ρ 0 the critical energy density, and h is the scale factor for Hubble expansion rate [90].In summary, our dark matter phenomenology is ruled by five free parameters.They encode the interplay between particle physics and cosmology, namely the new gauge coupling, g X , the dark matter mass, m N , the mediator mass, m Z ′ , and the (n, T r ), which are the cosmological input parameters that are related to the energy density of the scalar field that eventually governs the expansion rate of the universe.
We emphasize that the key cosmological input here is the energy density of the scalar field.It suffices to determine (n, T r ) and thus describe the entire background in which the dark matter was thermally produced [87].However, many single scalar field cosmologies can reproduce the behavior in Eq. (25).For n = 2, with positive scalar potential, it is associated to theories of quintessence fluids.In our work, we will assume that the energy density scales as ρ ϕ ∝ a −6 in the kination regime [91][92][93].Also, Chaplying gas theories can reproduce the redshift behavior in Eq. ( 25) and possibly explain dark matter and dark energy [94,95].Although, for n > 2, the theories with negative scalar potentials are needed [87,96].

B. Freeze-out During Early Radiation-dominated Era
Followed by a Matter Domination Period In a similar vein, we assume the existence of a scalar field ϕ that dominates the Hubble expansion rate during some period in the early universe before its decay into SM radiation.However, at this time, ϕ behaves as a pressureless fluid 3 .After ϕ decay, the universe is again dominated by radiation.
To describe the evolution of this system, we have to use the Boltzmann equations that couples the time evolution of the ϕ energy density ρ ϕ , the SM entropy density s and the dark matter number density n N [97][98][99][100][101], where E ≃ m 2 + 3T 2 is the averaged energy per dark matter particle 4 .From Eqs. ( 37) and (38), one could see that the injected entropy into SM radiation produced by the decay of ϕ can dilute the thermally produced dark matter components.Actually, it happens for any freeze-out scenario before the end of ϕ decays [102].
In Fig. 3, we briefly recall the early thermal history prior to BBN.The ERD takes place up to the temperature T eq at which the ρ ϕ domination starts.Thereafter, ρ ϕ effectively dominates over ρ R at the temperature T c , and the decay of ϕ ends at temperature T end at which the radiation dominates back again.FIG. 3. Illustration of the density energy evolution and thermal history for a scalar matter-field dominating in the early Universe.Here, the radiation-dominated dark matter freeze-out, discussed in Section III B, takes place in the green region (before T −1 eq ), while the matter-dominated freeze-out case, described in Section III C, happens in the blue region (after T −1 eq ) still before the energy density of the matter field takes a critical value ρ c ϕ at a critical temperature Tc.
The temperature T end is one of the free parameters that specify the cosmological background.It is defined as [98,102] where Γ ϕ stands for the total decay rate of ϕ into SM radiation.BBN restricts T end > ∼ 4 MeV [105][106][107].The other free parameter is, which is the ratio between the energy densities of the scalar field and radiation at T = m N .As we will be varying the dark matter mass, this ratio needs to be recomputed accordingly.Both parameters T end and κ are essentials to determine the dilution resulted from the ϕ decay.Notice that the Hubble rate is generally given by, with the dark matter freeze-out happening at temperature T f during the ERD as usual.
As the dark matter freezes-out at temperature higher than T eq , i.e, T eq ≪ T f , the Hubble in Eq. ( 41) is reduced to Eq. (28), which in terms of the time variable x becomes [102], We emphasize that dark matter decouples much before the decay of ϕ.Hence, the SM entropy is conserved, ds/dt = 0, and Eq. ( 38) turns into, which yields, which is the standard solution for the comoving dark matter number density long after the freeze-out and much before ϕ decays.In this way, the freeze-out temperature depends neither on T end nor on κ [102].As the dark matter abundance is firstly computed in the standard freeze-out case, but it will be changed due to ϕ decay (as showed in Fig. 4), this modification is set by the entropy injection episode that yields a dilution factor defined as, The dilution factor is simply the ratio between the SM entropy density at temperatures immediately after T 2 and before T 1 , the decay of the scalar field.We remind the reader that a similar reasoning is done in textbooks to obtain the temperature difference between photon and neutrinos after e + e − annihilations.Assuming the instantaneous decay approximation, the conservation of energy results into, As ρ ϕ (m) = κρ R (m), and taking T 2 = T end , the dilution factor is, In principle, D can take a wide range of values, but κ varies with the dark matter mass, and T end should be larger than 4 MeV.Moreover, κ should be smaller than one at high temperatures to guarantee that the freeze-out happens during a radiation phase.Anyway, combining Eqs. ( 44) and ( 47), we obtain the final comoving dark matter number density [102], and consequently, the overall dark matter relic density is found, where Ω std N h 2 stands for the dark matter relic abundance computed in the standard radiation freeze-out scenario.Notice that one cannot randomly pick D at will, pick any dilution factor wished to salvage dark matter models from exclusions, because the value chosen for D is still connected to the dark matter mass and properties of the scalar field.Very large dilution factors are not feasible.Within this production mechanism, we vary the free parameters (g X , m N , m Z ′ , T end , D) to outline the region of parameter space consistent with existing bounds.In the next section, we obtain the dark matter relic density for freeze-out occurring whilst ρ ϕ drives the expansion rate instead.FIG. 4. Illustration of the yield versus x (red line) provided by a dark matter freeze-out happening during or before a matter-dominated period in the early universe.For this case, we choose the dark matter mass mN = 100 GeV and the annihilation cross-section ⟨σv⟩ = 10 −9 GeV −2 .Choosing T end = 0.007 and D = 550.For completeness, we also include the yield for the standard case (black line).Again, the shaded gray region is the approximate range for the correct dark matter relic abundance as measured by Planck [13].

C. Freeze-out During Early Matter-dominated Era
Pre-BBN thermal history is the same as discussed in the last section.Similarly, we consider an unstable scalar field ϕ in the early universe which behaves as matter, ω ϕ = 0.Moreover, before and after the ρ ϕ dominates the expansion, the Universe is radiation-dominated.Now, we investigate the freeze-out process taking place during the period at which the scalar field ϕ governs the expansion rate.However, freeze-out happens much earlier than the decay of ϕ.Taking as reference Fig. 3, the freeze-out in this case happens between T −1 eq and T −1 c .The Hubble parameter is found in Eq. ( 41).Although, we define a temperature T ⋆ at which the ρ ϕ begins to evolve.It allows us to parameterize the relative fraction of the total energy density at T = T ⋆ in terms of [108], where r ∈ [0, 1].In this way, for T > T ⋆ > m N , we have, where H ⋆ ≡ H(T ⋆ ) and a ⋆ ≡ a(T ⋆ ).The conservation of entropy leads to, Inserting this result into Eq.( 51), the Hubble parameter be-comes approximately, (53) Thus, a matter-domination phase arises immediately at T ⋆ , (See Eq. ( 50)) for r ≪ 1, which leads to H ∝ T 3/2 .For r ≃ 1, the Hubble rate goes as H ∝ T 2 , which is associated to a radiation-domination epoch.Anyway, the Universe is radiation-dominated at T ⋆ , and the matter-domination phase will naturally arise at temperature, at which ρ ϕ accounts for 50% of the total density energy in agreement with [102,108].
The matter field that drives the expansion decays only when it effectively dominates over radiation at temperature, We stress that the freeze-out takes place in a matterdomination period, which occurs much before the scalar field decays.In other words, the freeze-out temperature should lie in the range T c ≪ T f ≪ T eq [102,108].We solve the Boltzmann equation for a matter-dominated universe, Eq. ( 38), with the Hubble given in Eq. ( 53).Taking the limiting case of matter-domination, r ≪ 1, and assuming a s-wave annihilation cross-section, we find the comoving dark matter number density after freeze-out long after the matter-domination epoch [108,109], where x ⋆ = m N /T ⋆ and g ⋆ s ≃ g ⋆ ≃ cte at the time of freezeout.However, it will be diluted due to entropy injection of the scalar field decay.In a similar vein, the dilution factor is the ratio between the entropy density computed at temperatures immediately before and after the scalar field decay, which in terms of T end and T ⋆ leads to, 5 [109], Therefore, the comoving dark matter number density becomes, and consequently, which is the dark matter relic density for a particle that freezes-out during a matter-domination era, which later experiences the decay of the scalar field ϕ.This scenario can also be represented by Fig. 4, where after freeze-out we see a dilution in the dark matter yield.In the next section, we present the constraints for our particle physics model.

IV. CONSTRAINTS
As the model is based on a 2HDM fashion there are several constraints which have been derived specifically to the canonical 2HDM, but some are also applicable to our model with proper adjustments.We will cover one by one below.

A. Collider
The most relevant collider limits stem from Z ′ searches at the LHC.As the Z ′ coupling to SM fermions is not suppressed, such Z ′ boson can be easily produced at the LHC, giving rise to either dijet or dilepton signal events.If the Z ′ coupling to dark matter were different and larger than the one with SM fermions, the invisible decay of the Z ′ boson could be significant to weaken the LHC lower mass bound.However, in our model, the Z ′ field interacts with equal strength to all fermions.Thus, the dilepton and dijet searches at the LHC are not meaningfully impacted by the presence of Z ′ decays into dark matter pairs [110,111].These two datasets have been considered, and the respective bound are displayed in our plots.

B. Flavor Physics
Interestingly, our model can also be constrained by flavor physics observed despite the absence of flavor changing interactions.The charged Higgs boson contributes to the b → sγ at one-loop level via the Cabibbo-Kobayashi-Maskawa matrix [112,113].The limit is reported in terms of tan β and the charged scalar mass.In our model, the mass of the charged Higgs is proportional to the v s , vev of the singlet scalar, and tan β.Taking tan β = 1, we can convert the lower mass limit on the charged Higgs mass into a bound on v s .Knowing that the Z ′ mass is controlled by v s and g X , for a given value of g X we now have a constraint on the Z ′ mass.We exhibit this limit in our plots.

C. Atomic Parity Violation
Atomic parity violation (APV) effects are often overlooked in physics beyond the Standard Model endeavours [56,57,114,115].The Z ′ gauge boson rising from the B-L gauge group has only vectorial interactions with fermions.In principle, that would lead to no APV, but in the presence of mass mixing with the Z boson this is no longer true.Atomic transitions in Cesium has been proved to be a great laboratory to probe Z ′ contributions to APV via mass mixing.Following [74], the contribution to the weak charge of Cesium in our model is found, (60) where δ is the mass mixing parameter, where we have defined tan β d = v1 vs [49].As the v s controls the Z ′ mass, with tan β = 1, we can display this limit in terms of the Z ′ mass for a given g X .We show this constraint in our plots.

D. Direct Detection
Dark matter particles might scatter off of nuclei, leaving a signal at direct detection experiments.We are considering heavy mediators only, thus the momentum transfer in the scattering process is much smaller than the mediator mass.Therefore, the dark matter interaction with quarks can be treated using effective operators.Our model features the effective Lagrangian, where q = u, d.We point out that only the valence quarks contribute due to the vector current present in Eq. (62).
In models where Z − Z ′ mass mixing is absent, Eq. ( 62) represent the only relevant source of direct detection signal.Interestingly, such term does not give rise to the standard spinindependent (SI) or spin-dependent (SD) dark matter-nucleon scattering signal.Hence, we need to map this Lagrangian onto a more general formalism using effective operators [116][117][118] to understand that this interaction gives rise to the operators, where , ⃗ q is the momentum transfer, and µ the reduced mass.
The first term in Eq. ( 63) features an enhancement with the atomic mass but is velocity suppressed, whereas the second is momentum suppressed and spin-dependent.That said, the first term will be dominant but will yield a distinct energy spectrum at the detector.We properly account for behavior using the code described in [119].
However, our model has a second Higgs doublet charged under the new gauge symmetry.Therefore, the Z − Z ′ mass mixing in unavoidable.This mixing induces the operator, where Λ here encodes both the Z mediator and Z ′ mediator via t-channel diagrams.This axial-vector interactions is known to generate a spin-dependent scattering which is found to be, where g A Z u,d (g A Z ′ u,d ) are the axial couplings of the Z (Z ′ ) boson with quarks, while g ZN (g Z ′ N ) are the axial couplings of the Z (Z ′ ) boson with dark matter.The explicit expressions can be extracted directly from Eq. ( 16), and ∆ p u,d,s are form factors accounting for the light quarks contributions to the nucleon spin [116,117].

V. RESULTS
In this section, we present the results for each scenario described above.We compute the dark matter relic density for different cosmological setups and compare aforementioned bounds.The results are presented in the m Z ′ vs m N plane.We cover several cosmological scenarios, varying their cosmological parameters.We also fixed the gauge coupling, g X = 1, since the main focus here is to study the impact of different cosmological parameters on this model.In any case, decreasing the gauge coupling g X will provide lower crosssections, leading to higher relic densities.To compensate for this change and to get the right relic density, we have to diminish the Z ′ mass, usually leading to more constrained scenarios 6 .For completeness, we also include the main Feynman diagrams responsible for the dark matter relic density, collider and direct searches, atomic parity violation and flavor physics in Fig. 5.  6 For a detailed study on the impact of the gauge coupling on the results we recommend [74].
We assess the impact of non-standard cosmology by comparing our findings with the one stemming from standard freeze-out.The colorful shaded regions represent the bounds from direct detection (blue region), dijet (red region), dilepton (cyan region), flavor physics (light green), APV (dark green region), and perturbative unitarity (magenta).We overlay the curves that delimit the parameter, yielding the correct relic density within the standard freeze-out (solid black).

A. Faster than usual early expansion
In the standard freeze-out case, the right-handed neutrino dark matter can nicely reproduce the correct relic density and evade existing bounds around the Z ′ resonance, and away from the resonance when m N > m Z ′ .In the context of a faster expansion rate at early times, the (n, T r ) parameters have an enormous impact on the dark matter relic abundance [87,120].Moreover, the dark matter mass also plays a key role.We explore n = 2, 4, 6 cosmologies, for temperatures T r = 0.1, 1, 10 GeV to have an overall idea of how a faster expansion affects the standard relic density.The largest temperature at the end of the relentless phase should be T r = 10 GeV, because we must enforce T f > T r to keep freeze-out occurring during the scalar field domination.That said, we exhibit in Fig. 6 the region of parameter in the m N vs m Z ′ plane that yields the correct relic density.For comparison, we also display with a black curve the region of parameter space that gives the right relic density in the standard freeze-out regime.For all plots we took g X = 1, in the top-left panel we adopted n = 2, in the top-right panel we assumed and n = 4, and in the bottom plot and n = 6.
It is clear that the parameters n and T r can enhance the dark matter relic density.We explored an interplay between n and T r , keeping one fixed and varying the other.In this cosmological background, a faster expansion, which enhances the dark matter abundance, the parameter space that reproduces the correct relic density is shifted to smaller Z ′ masses.One can understand this conclusion, remembering that the Z ′ mass enters the denominator of the thermal annihilation cross-section.Hence, the smaller the Z ′ mass, the smaller the abundance.This shift to lower Z ′ masses appears to counterbalance the fast expansion present during the dark matter freeze-out.As we increase further the expansion rate and take n = 4 (top-right-panel), and n = 6 (bottom-panel), the regions that reproduce the correct relic are also further shifted to smaller Z ′ masses.Furthermore, for each cosmological background n, the larger the ratio x r = m N /T r , the longer the relentless phase.Then, for a fixed T r , the heavier righthanded neutrino dark matter particles, the larger Z ′ masses (smaller cross-sections) to compensate the suppression to the dark matter relic density due to their longer time relentlessly annihilating.
The slope of the curves is mostly governed by n as it is closely related to the equation of state of the field ϕ [87].Notice that this fast expansion history is completely excluded by existing bounds, even in the Z ′ resonance.In the face of this, we will explore the impact of the early matter dominated field on this model.

B. Matter-Field Domination in the Early Universe
In this section, we deal with the non-standard freeze-out scenarios discussed in Sections III B and III C. Here, we use the cosmological parameters (T end , D) and (T ⋆ , ζ) for ERD and EMD, respectively.We firstly analyze the early radiationdominated freeze-out, thereafter, the dark matter decoupling during EMD in which the field ϕ drives the Hubble expansion, but it is still away from its decay.
a. Freeze-out During Early Radiation-dominated Era Followed by a Matter Domination Period In this case, we explore D = 550 and 2750, fixing T end = 0.007 GeV [98].As discussed earlier, κ is not fixed in order to adjust to the dark matter mass scales (See Eq.( 47)).It is common to evoke a constant and ad hoc dilution factor to bring down the dark matter relic density and circumvent existing bounds.Notice that this dilution factor cannot take any value, because it does depend on T end and κ.Naively, one might think that κ and T end are completely independent parameters.As κ is the ratio between the ϕ and radiation densities, thus depends on g * (T end )/g * (T = m N ) [102].Therefore, κ does feature a dependence on T end .As we need to impose T end > 4 MeV due to BBN bounds, clearly the dilution factor D cannot take arbitrarily large values.
In Fig. 7, we see that the contours that delimit the region of parameter space which yield the right dark matter abundance are free from existing bounds even away from the resonance, for m N < m Z ′ .Notice that the larger the entropy after the decay of the scalar field quantified by the dilution factors D = 550 and 2750, the smaller the cross-sections and the lighter the dark matter masses are required to obtain the correct relic density.In other words, there is a shift towards heavier Z ′ bosons and smaller dark matter masses.It is important to emphasize that we are fixing the value of T end = 0.007 and leaving the κ parameter free to get the right value of the dilution factor.If we fix the mass m and the κ parameter, T end will be inversely proportional to the dilution factor D, according to the Eq.(47).Conversely, the region of parameter space that reproduces the correct relic density in the standard freezeout is rather more constrained.Hence, a successful way to salvage WIMP models from restrictive direct detection and collider bounds is to assume a standard freeze-out followed by a short matter domination stage governed by the scalar ϕ, which then decays and injects entropy, shifting the relic density curve to the right side in Fig. 7.
b. Freeze-out During Early Matter-dominated Era In Fig. 8, we show the results for the early matter-dominated freeze-out.We have chosen ζ = 10 −6 and T ⋆ = 10 6 GeV, and ζ = 10 −8 and T ⋆ = 10 8 GeV.These choices are in agreement with the freeze-out during EMD condition, namely, T c ≪ T f ≪ T eq .We took T end ≃ 10 GeV, and r = 0.99.The blue and magenta curves delimit the region of parameter space that yields the correct relic density.It is visible that most region of parameter space is now consistent with existing bounds.
Our results comparing with the standard case (black lines), leads to scenarios completely free from the bounds, for dark matter masses larger than 200 GeV, for both ERD and EMD FIG. 6. Fast Expansion Case: We exhibit the curves that delimit the parameter space that leads to the correct relic density for Tr = 0.1 (brown), 1 GeV (blue) and 10 GeV (red).We also display the region which reproduces the correct relic density within the standard freeze-out for comparison.The colorful shaded regions represent the bounds from direct detection (blue region), dijet (red region), dilepton (cyan region), flavor physics (light green), APV (dark green region), and perturbative unitarity (magenta).cases.On the other side, for the standard case, the only regions that survive the limits are near the Z ′ resonance.Therefore, for our model, the search for dark matter masses around the T eV -scale may hint the presence of non-standard cosmological histories in the early universe.

C. The relentless phase versus the entropy injection
In this section, we put all results into perspective and discuss the key findings concerning the overall dark matter density in the presence of non-standard cosmologies.In Fig. 9, we present the results for the scenarios in which a quintessence fluid (red curve) with n = 2 and T r = 10 GeV dominates the expansion rate (Faster Expansion), the standard early radiation domination occurs (ERD), but it is followed by the decay It is remarkable that for a quintessence-dominated freezeout, the contour arises in the left-hand side of the standard freeze-out.This reflects the fact that the larger Hubble rate, the larger the annihilation cross-section to reproduce the correct relic density.
Instead, for freeze-out happening during periods in which the Hubble rate is equivalent to or slower than radiation, the contour appear in the right-hand side, which is less constrained by data.Hence, the presence of a scalar field that drive the Hubble rate during or after the dark matter freezeout represents an important direct to be explored by the next generation of experiments.

VI. CONCLUSIONS
In this work, we investigated the dark matter phenomenology of a weak scale right-handed neutrino dark matter under three different cosmological settings.We concretely incorporated this right-handed neutrino into a Two Higgs Doublet Model featuring an Abelian gauge symmetry.The model itself is well-motivated for being able to account for neutrino masses and absence of flavor changing interactions in the scalar sector.
The dark matter phenomenology is driven by the presence of a Z ′ field, as well as by the Z − Z ′ mass mixing that arises because a Higgs doublet is involved in the spontaneous symmetry breaking of the new Abelian gauge symmetry.Concerning non-standard cosmology, we explored the case where the universe expands faster than usual, and the scenarios where the dark matter freeze-out takes place during a matter domination epoch, as well as during a usual radiation domination, but it is then followed by a matter domination phase.
For the faster expanding case, we found a very constrained scenario, since larger cross-sections are necessary to reproduce the right relic density, whereas when a matter field dominates the expansion rate in the early universe either during or after the dark matter freeze-out the relic density curves are shift away from existing bounds.After putting our results into perspective with flavor bounds, direct detection, collider and atomic parity violation, we solidly conclude that weak scale right-handed neutrino stands for a plausible dark matter candidate in the presence of an early matter domination epoch.

ΩNh² ~0. 11 ⋯ 4 YeqFIG. 1 .
FIG. 1.We exhibit of the yield versus x for a freeze-out happening during the fast expansion regime.For this case, we choose the dark matter mass mN = 100 GeV and the annihilation cross-section ⟨σv⟩ = 10 −9 GeV −2 , for the following values of Tr = 0.1 (brown line), Tr = 1 GeV (blue line), Tr = 10 GeV (red line), and n = 4.For completeness, we also include the yield for the standard case (black line).Note that for Tr = 10 GeV, we have a superposition with the standard case.The shaded gray region is the approximate range for the correct dark matter relic abundance as measured by Planck[13].It is clear that the longer the period of ϕ dominance, the larger the impact on the yield.

FIG. 5 .
FIG. 5.Here we have the Feynman diagrams related to the dark matter phenomenology of this work.The diagram (a) gives the dark matter relic abundance, (b) corresponds to collider searches, (c) is associated to the atomic parity violation bounds, (d) to direct detection searches, and (e) & (f) comes from flavor physics.

FIG. 7 .
FIG. 7. The correct dark matter relic abundance for dark matter freeze-out happening during the early radiation-dominated epoch.For gX = 1, the cosmological parameters are fixed to T end = 0.007 GeV and D = 550 and 2750 (blue and red curves, respectively).The constraints are the same as in Fig.6.