This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

Freeze-In Dark Matter within the Seesaw mechanism

Michele Lucente michele.lucente@unibo.it Dipartimento di Fisica e Astronomia, Università di Bologna, via Irnerio 46, 40126 Bologna, Italy INFN, Sezione di Bologna, viale Berti Pichat 6/2, 40127, Bologna, Italy Institute for Theoretical Particle Physics and Cosmology (TTK), RWTH Aachen University, D-52056 Aachen, Germany
(September 7, 2025)
Abstract

We show that the minimal Type-I Seesaw mechanism can successfully account for the observed dark matter abundance in the form of a keV sterile neutrino. This population can be produced by the decay of the heavier neutral leptons, with masses above the electroweak mass scale, while they are in thermal equilibrium in the early Universe (freeze-in). Moreover, the implementation of the relevant phenomenological constraints (relic abundance, indirect detection and structure formation) on this model automatically selects a region of the parameter space featuring an approximate lepton number symmetry.

keywords:
Type-I seesaw , dark matter , freeze-in , neutrinos , lepton number symmetry
††journal: Physics Letters B

1 Introduction

The nature of dark matter (DM) and the mechanism at the origin of neutrino masses (and lepton flavour mixing) are among the most pressing questions in particle physics and cosmology. A minimal and natural extension of the Standard Model (SM) of particle physics is the Type-I Seesaw mechanism [1, 2, 3, 4, 5, 6], in which the field content of the SM is extended by including neutrino fields with right-handed chirality (RHN). Being gauge singlets, RHN can have Majorana mass terms whose origin is not related to the SM Higgs mechanism, thus leaving the RHN mass scale as a phenomenological free parameter; for the same reason, the number of RHN fields is not fixed by the requirement of anomaly cancellation, although this is not generally true in the context of gauge extensions of the SM, where the number of RHN has to match the number of SM generations (n=3n=3) to ensure anomaly cancellation (cf. e.g. [7]). We also notice that, in the context of the minimal SM, there exists only one global symmetry which is anomaly free, given by the difference between the baryon number BB and lepton number LL, i.e. B−LB-L.

Sterile neutrinos with masses at the keV scale can be suitable DM candidates (c.f. [8] for a review), subject to a number of phenomenological constraints concerning their current abundance [9], lifetime [10] and impact on structure formation (notably as inferred by observations of the Lyman-α\alpha forest [11, 12, 13, 14, 15, 16]): the first two mainly constrain the mixing θ\theta between active and sterile neutrinos, while the latter strongly depends on the DM production mechanism. These particles are unavoidably produced in the early Universe as long as a non-zero mixing with the active neutrinos is present (Dodelson-Widrow mechanism, DW) [17]. For masses below 28.8 keV (95% C.L.), this production mechanism results in a DM population with a too large free-streaming length (classified as warm dark matter, WDM), incompatible with the observed large scale structure of the Universe [18]. For larger masses, the constraints on θ\theta derived from the non-observation of the sterile neutrino DM decay products [19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40] prevent the DW production to account for the observed DM relic density. The conversion of active into sterile neutrinos can be resonantly enhanced if the primordial plasma features a sizeable lepton asymmetry (Shi-Fuller mechanism, SF) [41]; the effect is twofold: on one hand, smaller values of θ\theta are required to account for the observed DM relic density, thus evading the bounds from indirect detection. On the other hand, the DM free streaming is reduced with respect to DW, thus resulting in a colder DM population and relaxing structure formation bounds. These ideas have been coherently implemented in the so-called ν\nuMSM model [42, 43], which is a phenomenologically motivated implementation of the Type-I Seesaw mechanism with 3 RHN: one of them at the keV scale acting as DM candidate, while the other two have strongly degenerate masses at the GeV scale and, after electroweak symmetry breaking (EWSB), produce a lepton asymmetry which acts as a background for the resonant DM production via SF mechanism. Moreover, the same heavy sterile neutrinos can as well account for the observed Baryon Asymmetry of the Universe via leptogenesis before EWSB [44, 45].

In this Letter, we point out that an alternative implementation of the minimal Type-I Seesaw mechanism exists, which is able to account for the observed DM relic density111Sterile neutrino DM has been proposed as well in models with an extended gauge or scalar sector, see e.g. [46, 47, 48, 39].. The DM candidate is still a keV sterile neutrino, but the other sterile neutrinos have a much larger mass scale, above the electroweak one. They produce the correct DM abundance while in thermal equilibrium by decaying into a keV sterile neutrino and a Higgs boson: this production mechanism has been firstly pointed out in [47] in the framework of a minimal Inverse Seesaw realisation [49].

2 The model

We extend the SM field content with the addition of 3 RHN fields NIN_{I}, I=1,2,3I=1,2,3, leading to the most general renormalizable Lagrangian:

ℒ=ℒSM+iNI¯∂̸NI−(12NIc¯MI​JNJ+Fα​IℓLα¯ϕ~NI+h.c.),\mathcal{L}=\mathcal{L}_{\mathrm{SM}}+i\overline{N_{I}}\not{\partial}N_{I}-\left(\frac{1}{2}\overline{N_{I}^{c}}M_{IJ}N_{J}+F_{\alpha I}\overline{\ell_{L}^{\alpha}}\tilde{\phi}N_{I}+h.c.\right), (1)

where ℓLα\ell_{L}^{\alpha} is the left-handed S​U​(2)LSU(2)_{L} lepton doublet of flavour α=e,μ,τ\alpha=e,\mu,\tau, MI​JM_{IJ} is a symmetric matrix of Majorana mass terms, Fα​IF_{\alpha I} are dimensionless Yukawa couplings and Nc=i​γ2​γ0​N¯TN^{c}=i\gamma^{2}\gamma^{0}\overline{N}^{T}, ϕ~=i​σ2​ϕ∗\tilde{\phi}=i\sigma^{2}\phi^{*}. After EWSB, the Yukawa couplings generate a Dirac mass term mD=v​Fm_{D}=vF, where vv is the Higgs vacuum expectation value (vev, taking the value v0=174v_{0}=174 GeV at zero temperature), and Eq. (1) results in a non-vanishing Majorana mass matrix mνm_{\nu} for the light active neutrinos; in the Seesaw limit, |MI​J|≫v​|Fα​I||M_{IJ}|\gg v|F_{\alpha I}|, it reads mν≃−v2​FT​M−1​Fm_{\nu}\simeq-v^{2}F^{T}M^{-1}F, while at leading order in the expansion parameter Θ=v​F/M\Theta=vF/M the heavier mass eigenstates have masses coinciding with the eigenvalues of the matrix MM. Without loss of generality, it is always possible to chose a basis in which the matrix MM is real and diagonal, in which case the mixing elements coupling the active lepton flavour α\alpha to the heavy mass eigenstate II are given by Uα​I∗=mDα​I/MIU_{\alpha I}^{*}=m_{D}^{\alpha I}/M_{I}. In order to ensure the agreement of the model with experimental neutrino mixing parameters, it is convenient to parametrise the Yukawa couplings using the Casas-Ibarra (CI) parametrisation [50]:

mD=−i​UPMNS∗​m^​R​M,m_{D}=-i\ U_{\mathrm{PMNS}}^{*}\sqrt{\hat{m}}\ R\ \sqrt{M}, (2)

where m^\hat{m} is a diagonal matrix containing the mass eigenvalues of the light (mostly active) neutrinos at zero temperature, UPMNSU_{\mathrm{PMNS}} is the unitary PMNS mixing matrix [51] and RR is an orthogonal matrix parametrised by 3 complex angles ωi​j\omega_{ij},

R=V23​V13​V12,R=V_{23}V_{13}V_{12}, (3)

with

V12=(cos⁡ω12sin⁡ω120−sin⁡ω12cos⁡ω120001),V_{12}=\left(\begin{array}[]{ccc}\cos\omega_{12}&\sin\omega_{12}&0\\ -\sin\omega_{12}&\cos\omega_{12}&0\\ 0&0&1\end{array}\right), (4)

and analogous definitions for V23,V13V_{23},V_{13}.

2.1 Model parameters

Low energy neutrino oscillation data fix the value of UPMNSU_{\mathrm{PMNS}} and (squared) mass differences for m^\hat{m} in Eq. (2). In order to be a viable DM candidate, a sterile neutrino must have a mass at the keV scale [8]: we will thus explore realisations with M1M_{1} in the range [1,100][1,100] keV. The heavier mass eigenstates must lie above the electroweak scale for the freeze-in production mechanism proposed in [47] to take place, but not too far from it; for definiteness, we will chose them to be almost degenerate and fix the values M2,3=300M_{2,3}=300 GeV. The large hierarchy of values between M1M_{1} and M2,3M_{2,3} may seem a fine-tuned scenario, however we will show that this can be the natural result of an approximate lepton number symmetry. Most importantly, this symmetry-motivated scenario gets automatically realised once we impose astrophysical and cosmological constraints to the model.

2.2 Phenomenological constraints

2.2.1 Cosmology and astrophysics

There are three broad categories of constraints that apply to the model under discussion. The first requirement is to not overclose the Universe, by producing a larger DM abundance than the observed one, ΩDMobs​h2=0.12\Omega^{\mathrm{obs}}_{\mathrm{DM}}h^{2}=0.12 [9] (where hh is the dimensionless Hubble parameter, H0=100​h​km s−1​Mpc−1H_{0}=100\ h\ \textrm{km s}^{-1}\ \textrm{Mpc}^{-1}). In the absence of a sizeable lepton asymmetry, there are two production mechanisms to consider: the first one is the well known DW, effective at temperatures T≈150T\approx 150 MeV and resulting in a relic density

ΩDMD​W​h2=0.11⋅105​M1keV​∑αCα​(M1)​|(UPMNS∗​m^eV​R)α​1|2,\Omega_{\mathrm{DM}}^{DW}h^{2}=0.11\cdot 10^{5}\frac{M_{1}}{\mathrm{keV}}\sum_{\alpha}C_{\alpha}(M_{1})\left|\left(U_{\mathrm{PMNS}}^{*}\sqrt{\frac{\hat{m}}{\mathrm{eV}}}R\right)_{\alpha 1}\right|^{2}, (5)

where CαC_{\alpha} are dimensionless parameters of order unity, whose exact value depends on the lightest sterile neutrino mass M1M_{1} [52]. It is evident from Eq. (5) that the requirement ΩDMD​W​h2≤0.12\Omega_{\mathrm{DM}}^{DW}h^{2}\leq 0.12 imposes a strong upper bound on the elements of the first column of RR, and/or an upper bound on the light active neutrino masses (see [53] for an extended discussion on this second point). The second production mechanism, i.e. freeze-in production from the decay of heavy neutrinos in thermal equilibrium, has never been included (to the best of our knowledge) in previous studies concerning the minimal Type-I Seesaw mechanism. It is mainly effective at temperatures T≈M2,3T\approx M_{2,3}, but the related couplings are different from zero only after EWSB; its resulting relic abundance can be estimated as,

ΩDMF​I​h2=2.16⋅1022​∑J=2,3|(R†​m^v​R)1​J|2​gJ​(1−mh2MJ2)2​ε​(MJ),\Omega_{\mathrm{DM}}^{FI}h^{2}=2.16\cdot 10^{22}\sum_{J=2,3}\left|\left(R^{\dagger}\frac{\hat{m}}{v}R\right)_{1J}\right|^{2}g_{J}\left(1-\frac{m_{h}^{2}}{M_{J}^{2}}\right)^{2}\varepsilon\left(M_{J}\right), (6)

with mhm_{h} the Higgs boson mass and gJg_{J} the internal degrees of freedom of the heavy neutrinos. The function ε​(MJ)\varepsilon(M_{J}) was defined in [47] (following [54]), and takes into account the evolution of the Higgs vev vv and the suppression of the decay rate in DM due to the electroweak symmetry restoration for decays happening at temperatures larger than the EWSB scale (notice that since mν∝v2m_{\nu}\propto v^{2} the ratio mν/vm_{\nu}/v vanishes for v=0v=0). The factorisation in Eq. (6) neglects the evolution of the Higgs mass with vv for temperatures approximately in the range 160 to 130 GeV, where the sharp electroweak crossover takes place [55], thus slightly underestimating the available phase space in the decay and thus the final DM abundance. For sub-eV active neutrinos, we have (m^/v0)2≲10−23(\hat{m}/v_{0})^{2}\lesssim 10^{-23}, giving the DM abundance in the correct ballpark, provided RR is of order unity. We notice that, although other decay channels mediated by the SM gauge bosons are potentially relevant at zero temperature, thermal effects strongly suppress the production of sterile neutrinos from gauge vector boson decays at electroweak temperatures [56]; we thus neglect these decay channels as they give a negligible contribution to the final abundance.

The second constraint comes from requiring the DM stability over cosmological timescales: keV sterile neutrinos can indeed radiatively decay into a photon and a light active neutrino [10], with a lifetime proportional to the active-sterile mixing angle, that one can write as

θ12=∑α|(UPMNS∗​m^​R)α​1|2​M1−1.\theta_{1}^{2}=\sum_{\alpha}\left|\left(U_{\mathrm{PMNS}}^{*}\sqrt{\hat{m}}\ R\right)_{\alpha 1}\right|^{2}M_{1}^{-1}. (7)

Even if the sterile neutrino lifetime exceeds the age of the Universe, a fraction of DM can decay at present time, resulting in a monochromatic gamma-ray line with energy M1/2M_{1}/2 (i.e. X-rays). The non-observation of this signal results in an upper bound on θ12\theta_{1}^{2} as a function of the DM mass M1M_{1} [19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40].

Finally, the DM momentum distribution affects the structure formation of the Universe: sterile neutrinos produced via the DW mechanism are classified as WDM, meaning they have a sizeable free streaming length that erases substructures at galactic scales. A WDM dominated Universe is at odd with observation, while mixed scenarios where WDM constitutes a fraction fWDMf_{\mathrm{WDM}} of the current DM relic density are potentially allowed, depending on the value of fWDMf_{\mathrm{WDM}} and on the free-streaming length of the WDM component. In this study we adapt the 95% C.L. results from [18], whose bounds can be translated for the case of sterile neutrinos produced via DW. Sterile neutrinos produced via freeze-in have a colder spectrum, compatible with the structure formation constraints for masses approximately above few keV [57, 47, 8] (see also [58, 59, 60]). In the present Letter we consider DM masses as light as 1 keV, keeping in mind that the validation of solutions with masses smaller than approximately 10 keV requires a dedicated study, that goes beyond the scope of the present work. Finally, a lower bound on fermionic DM mass resulting from the phase space density of DM in dwarf spheroidal galaxies [11] excludes (disfavours) values below 400 eV (2 keV) (see [8] and references therein).

2.2.2 Laboratory

We require the model in Eq. (1) to correctly reproduce neutrino oscillation data [61] at the 3σ\sigma accuracy level. In addition, we impose constraints on the two heavy neutral leptons (HNL) from the non-observation of further beyond the SM signals: these include upper bounds on the active-sterile mixing from direct searches at accelerators [62, 63, 64, 65], deviations from unitarity of the leptonic mixing matrix [66] and neutrinoless double beta decay (0​ν​β​β0\nu\beta\beta) [67].

2.2.3 Theory

The CI parametrisation in Eq. (2) ensures that neutrino oscillation parameters are correctly reproduced. However, since the Yukawa couplings grow exponentially with the imaginary parts of the angles ωi​j\omega_{ij}, these cannot be arbitrarily large, either because the theory becomes non-perturbative or because the Seesaw hypothesis |v​F|≪|M||vF|\ll|M| gets violated, invalidating the CI derivation itself. In our study, we impose an upper bound of 4​π4\pi on each Yukawa coupling |Fα​I||F_{\alpha I}|, and require in addition that ΓI<MI/2\Gamma_{I}<M_{I}/2 for each state [68, 69, 70, 71, 72, 73]. Moreover, for each realisation of the Dirac matrix (2) we explicitly diagonalise the full Lagrangian (1), excluding solutions not reproducing neutrino oscillation data at the 3​σ3\sigma level. Finally, we require the largest Yukawa coupling for each heavy neutrino to be larger than 2⋅10−7\sqrt{2}\cdot 10^{-7} [74], so as to ensure that the state thermalises in the early Universe (a necessary condition for freeze-in production).

Refer to caption
Refer to caption
Refer to caption
Figure 1: Parameter space of solutions for the model with the choice of parameters as described in the main text, for different DM masses. The black line (grey region) reproduces (exceeds) the observed relic density. Cyan and yellow regions are excluded by Lyman-α\alpha and X-ray observations, respectively, while in the blue region the model fails in reproducing neutrino oscillation data or is excluded by HNL searches. The black-dashed line shows the expected solution by considering the DW production mechanism only.

3 Methodology

We chose to fix the neutrino oscillation parameters (mass differences, mixing angles and Dirac CP-violating phase) to their best-fit values in [61], and we set the CP-violating Majorana phases for the light neutrino states to zero, as their impact on the addressed dynamics is expected to be subdominant. To enhance the predictivity of the model, we set to zero the lighter neutrino mass as well as the complex angle ω13\omega_{13}, and assume that ω12\omega_{12} and ω23\omega_{23} are purely real and imaginary numbers, respectively. Together with fixing the heavy neutrino masses to 300 GeV, this leaves only 3 real free parameters: M1M_{1}, Reω12\real\omega_{12} and Imω23\imaginary\omega_{23}. This choice allows us to gain useful insights on the role played by the different parameters in the model.

4 Results

We report in Fig. 1 the parameter space of solutions for a Normal Ordering (NO) of light active neutrino masses, and three choices of DM mass: one where Lyman-α\alpha bounds are dominant, M1=2M_{1}=2 keV, one that reproduces the (not yet confirmed) monochromatic gamma-ray line at 3.55 keV [75, 76] (M1=7.1M_{1}=7.1 keV), and finally one where X-ray bounds are dominant, M1=30M_{1}=30 keV. In the plots, the black lines represent the choices of parameters that reproduce the observed DM relic abundance, ΩDM​h2=0.12\Omega_{\mathrm{DM}}h^{2}=0.12, while the grey region above is excluded by the overproduction of DM. The cyan and yellow regions are excluded by Lyman-α\alpha and X-ray observations, respectively, while in the blue region (on the top) the model does not reproduce neutrino oscillation data or is excluded by the HNL searches previously described: we find that the most constraining bounds in the present scenario are the ones imposed by the HNL contributions to charged lepton flavour violating processes, parametrised as deviations from unitarity of the leptonic mixing matrix (cf. e.g. [77, 66] and references therein). The black-dashed line represents the expected solution if the freeze-in production is neglected, i.e. if all the relic abundance is due to the DW production: as it is well known, this solution is at odds with astrophysical observations. We also notice that the freeze-in production is less effective for decreasing values of Imω23\imaginary\omega_{23}. The results for Inverted Ordering (IO) of light active neutrino masses are qualitatively similar, but this scenario corresponds to larger Yukawa couplings for fixed CI parameters, cf. Eq. (2). This implies that the upper bound on Imω23\imaginary\omega_{23} from HNL constraints is stronger, Imω23<9.87\imaginary\omega_{23}<9.87 (while Imω23<10.56\imaginary\omega_{23}<10.56 for NO).

By increasing the value of the DM mass, X-ray constraints on θ12\theta^{2}_{1} become increasingly stronger (given that the sterile neutrino decay width is proportional to the fifth power of its mass): we find that in the considered scenario this closes the parameter space of solutions for masses M1M_{1} above 60 keV (49 keV) for NO (IO).

Refer to caption
Refer to caption
Refer to caption

Refer to caption

Figure 2: Allowed and excluded regions for the complex angles ωi​j\omega_{ij} in the CI parametrisation, with the other parameters chosen as described in the main text: green points are allowed by experimental constraints and theoretical considerations, while other colours denote exclusion by at least one of the criteria reported in the Figure. Finally, black points reproduce the observed DM relic density at the 3​σ3\sigma precision.

4.1 Phenomenologically emergent approximate lepton number symmetry

The parameter space of solutions discussed in the previous section may seem strongly fine tuned, featuring a large hierarchy of values between the DM mass and the heavier HNLs, as well as between Imω23\imaginary\omega_{23} and the other parameters of the RR matrix. However, this is exactly the limit in which global lepton number conservation is recovered in the CI parametrisation: in the limit M1→0M_{1}\rightarrow 0, M2→M3M_{2}\rightarrow M_{3}, ω12,ω13→0\omega_{12},\omega_{13}\rightarrow 0 and ω23→±i​∞\omega_{23}\rightarrow\pm i\infty, the new fields arrange in a decoupled massless state N1N_{1} and a Dirac fermion given by a linear combination of the fields N2N_{2} and N3N_{3}, thus resulting in a theory with lepton number conserving interactions only. In the exactly symmetric limit light active neutrinos are massless, m^=0\hat{m}=0, thus the symmetry must be broken at some level. This also provides a possible explanation for the lightness of active neutrinos when compared to the other SM fermions [78, 79, 80, 81, 82, 83]. Concerning the specific values of MiM_{i}, the symmetry only imposes relations among them, while the absolute scales M1M_{1} and M2,3M_{2,3} must derive from the (unknown) mechanism at the origin of these Majorana mass terms. Without assuming a specific mechanism, the restrictions on these parameters are mostly phenomenological in nature. In particular, to be a viable DM candidate M1M_{1} must lie around the keV scale (smaller values are excluded by structure formation, larger ones give a too short DM lifetime), while the scale of M2,3M_{2,3} is determined by the assumed DM production mechanism (much smaller values imply too small Yukawa couplings to produce a sizeable DM population, while for much larger values the parent particles decay at temperatures above the electroweak phase transition, when the active-sterile mixings are vanishing).

We report in Fig. 2 the results of a general scan on the parameter space of the model, performed for NO, fixing the DM mass M1=5M_{1}=5 keV and randomly varying the parameters Reωi​j\real\omega_{ij} and Imωi​j\imaginary\omega_{ij} in the range [10−10,2​π][10^{-10},2\pi] and [10−10,20][10^{-10},20] respectively, drawing from a log-flat distribution. We also allow the physical Majorana phase in the active neutrino sector to vary in the [0,2​π]\left[0,2\pi\right] range. We conclude that astrophysical constraints set an upper bound |ω1​i|≲6⋅10−3\left|\omega_{1i}\right|\lesssim 6\cdot 10^{-3} (i=2,3i=2,3) (clearly, values of Reω1​i\real\omega_{1i} above (2​π−10−3)(2\pi-10^{-3}) are as well allowed), while no bound is identified for Reω23\real\omega_{23} and the upper bound on Imω23\imaginary\omega_{23} (Imω23≲11\imaginary\omega_{23}\lesssim 11) is solely determined by the requirement of not overproducing the DM relic density. We thus conclude that experimental constraints force the parameters ω12\omega_{12} and ω13\omega_{13} to be much smaller than unity. In addition, the requirement of reproducing the observed DM relic abundance via freeze-in decay requires HNLs heavier than the electroweak scale, thus introducing a hierarchy of masses, as well as Imω23≈10\imaginary\omega_{23}\approx 10, recovering the approximate symmetric limit. Concerning the mass degeneracy between M2,3M_{2,3}, this suppresses the 0​ν​β​β0\nu\beta\beta decay rate222We thank Jordy de Vries for having encouraged us to explicitly verify that very large mass splittings in this region of the parameter space can result in 0​ν​β​β0\nu\beta\beta rates exceeding current experimental bounds. which could be sizeable for large values of Imω23\imaginary\omega_{23}, although mild mass splittings of 𝒪​(10)\mathcal{O}(10) GeV are nevertheless allowed; relaxing this condition relaxes the symmetry, thus resulting in a less justified framework from the model building point of view, while opening at the same time to the possibility of large loop corrections to the active neutrino masses [84]. The approximate lepton number symmetry is known to play a role in justifying the structure of the ν\nuMSM [80], as well as in other low-scale leptogenesis models [85, 86, 87].

5 Conclusion

We have shown that, in the minimal SM extended with 3 RHN, the freeze-in production of a sterile neutrino DM can account for the observed relic abundance, while complying with existing constraints. This production channel has been, to the best of our knowledge, neglected in previous studies, and effectively opens a new parameter space of solutions for the sterile neutrino DM hypothesis. This space of solutions is characterised by a hierarchy between different parameters of the model (sterile neutrino masses as well as complex angles in the CI parametrisation) that points towards an approximate lepton number symmetry, which could be embedded in a larger B−LB-L framework. When presenting solutions, we have worked under a simplifying assumption in order to reduce the number of free parameters and gain a deeper insight on the model; the actual parameter space of solutions is thus larger than the one presented here, and we postpone a dedicated statistically comprehensive study to a future work.

Acknowledgments

I am deeply grateful to Asmaa Abada, Giorgio Arcadi, Gioacchino Piazza and Salvador Rosauro-Alcaraz for insightful discussions on the impact of thermal effects on the production mechanisms. I thank Marco Drewes for useful discussions on the connection between the CI parametrisation and global lepton number symmetry, as well as Asmaa Abada and Felix Kahlhoefer for useful comments on the manuscript. This work is funded by the European Union under the Horizon Europe’s Marie Sklodowska-Curie project 101068791 - NuBridge, and has been supported by the Alexander von Humboldt Foundation.

Appendix A Explicit numerical example

We report hereafter the numerical values for an example point that reproduces the correct DM abundance via Eq.s (5, 6).

The input values in the CI parametrisation are

mνlightest\displaystyle m_{\nu}^{\textrm{\tiny{lightest}}} =\displaystyle= 0,\displaystyle 0, (8)
M1\displaystyle M_{1} =\displaystyle= 5​keV,\displaystyle 5\ \textrm{keV}, (9)
M2\displaystyle M_{2} =\displaystyle= 300​GeV,\displaystyle 300\ \textrm{GeV}, (10)
M3\displaystyle M_{3} =\displaystyle= 300​GeV,\displaystyle 300\ \textrm{GeV}, (11)
ω12\displaystyle\omega_{12} =\displaystyle= 5.184×10−10+i​ 3.867×10−8,\displaystyle 5.184\times 10^{-10}+i\ 3.867\times 10^{-8}, (12)
ω13\displaystyle\omega_{13} =\displaystyle= 9.678×10−7+i​ 6.241×10−8,\displaystyle 9.678\times 10^{-7}+i\ 6.241\times 10^{-8}, (13)
ω23\displaystyle\omega_{23} =\displaystyle= 2.439×10−5+i​ 8.398,\displaystyle 2.439\times 10^{-5}+i\ 8.398, (14)
α1\displaystyle\alpha_{1} =\displaystyle= 3.763,\displaystyle 3.763, (15)
α2\displaystyle\alpha_{2} =\displaystyle= 2.021,\displaystyle 2.021, (16)

while the light neutrino mass differences and PMNS parameters are fixed to the best fit values from [61].

With this input the following mass matrix results

ℳ=(0mDmDTMM)​GeV,\mathcal{M}=\left(\begin{array}[]{cc}0&m_{D}\\ m_{D}^{T}&M_{M}\end{array}\right)\ \textrm{GeV}, (17)

with

mD=(𝐦𝐃𝟏𝐦𝐃𝟐𝐦𝐃𝟑),m_{D}=\left(\begin{array}[]{ccc}\mathbf{m_{D}^{1}}&\mathbf{m_{D}^{2}}&\mathbf{m_{D}^{3}}\end{array}\right), (18)

and the individual vectors are given by

𝐦𝐃𝟏\displaystyle\mathbf{m_{D}^{1}} =\displaystyle= (−1.233×10−12+i​ 3.558×10−122.401×10−11+i​ 2.399×10−111.673×10−11+i​ 5.020×10−12),\displaystyle\left(\begin{array}[]{c}-1.233\times 10^{-12}+i\ 3.558\times 10^{-12}\\ 2.401\times 10^{-11}+i\ 2.399\times 10^{-11}\\ 1.673\times 10^{-11}+i\ 5.020\times 10^{-12}\\ \end{array}\right), (22)
𝐦𝐃𝟐\displaystyle\mathbf{m_{D}^{2}} =\displaystyle= (−2.786×10−2−i​ 7.779×10−3−1.726×10−1+i​ 1.954×10−1−3.060×10−2+i​ 1.307×10−1),\displaystyle\left(\begin{array}[]{c}-2.786\times 10^{-2}-i\ 7.779\times 10^{-3}\\ -1.726\times 10^{-1}+i\ 1.954\times 10^{-1}\\ -3.060\times 10^{-2}+i\ 1.307\times 10^{-1}\\ \end{array}\right), (26)
𝐦𝐃𝟑\displaystyle\mathbf{m_{D}^{3}} =\displaystyle= (7.779×10−3−i​ 2.786×10−2−1.954×10−1−i​ 1.726×10−1−1.307×10−1−i​ 3.060×10−2),\displaystyle\left(\begin{array}[]{c}7.779\times 10^{-3}-i\ 2.786\times 10^{-2}\\ -1.954\times 10^{-1}-i\ 1.726\times 10^{-1}\\ -1.307\times 10^{-1}-i\ 3.060\times 10^{-2}\\ \end{array}\right), (30)
MM\displaystyle M_{M} =\displaystyle= (5.×10−60.0.0.3.×1020.0.0.3.×102).\displaystyle\left(\begin{array}[]{ccc}5.\times 10^{-6}&0.&0.\\ 0.&3.\times 10^{2}&0.\\ 0.&0.&3.\times 10^{2}\\ \end{array}\right). (34)

Upon diagonalisation of ℳ\mathcal{M}, the following values for the mass eigenstates are obtained

(08.6139397443927×10−35.0169629398723×10−25.0000000002942×1033.0000028935617×10113.0000028935621×1011)​eV,\left(\begin{array}[]{c}0\\ 8.6139397443927\times 10^{-3}\\ 5.0169629398723\times 10^{-2}\\ 5.0000000002942\times 10^{3}\\ 3.0000028935617\times 10^{11}\\ 3.0000028935621\times 10^{11}\\ \end{array}\right)\ \textrm{eV}, (35)

with the PMNS matrix (i.e. the 3×33\times 3 upper-left block of the full mixing matrix 𝒰\mathcal{U})

UPMNS=(0.830.17−0.52​i−0.11−0.1​i−0.27+0.03​i0.20−0.57​i0.40+0.63​i0.49+0.02​i−0.16+0.56​i0.34+0.55​i),U_{\textrm{PMNS}}=\left(\begin{array}[]{ccc}0.83&0.17-0.52i&-0.11-0.1i\\ -0.27+0.03i&0.20-0.57i&0.40+0.63i\\ 0.49+0.02i&-0.16+0.56i&0.34+0.55i\\ \end{array}\right), (36)

and the following values for the active-sterile mixings (i.e. the elements in the upper-right block of 𝒰\mathcal{U})

𝒰α​𝟒\displaystyle\mathbf{\mathcal{U}_{\alpha 4}} =\displaystyle= (−2.466×10−7−i​ 7.115×10−74.801×10−6−i​ 4.797×10−63.346×10−6−i​ 1.004×10−6),\displaystyle\left(\begin{array}[]{c}-2.466\times 10^{-7}-i\ 7.115\times 10^{-7}\\ 4.801\times 10^{-6}-i\ 4.797\times 10^{-6}\\ 3.346\times 10^{-6}-i\ 1.004\times 10^{-6}\\ \end{array}\right), (40)
𝒰α​𝟓\displaystyle\mathbf{\mathcal{U}_{\alpha 5}} =\displaystyle= (−9.286×10−5+i​ 2.593×10−5−5.753×10−4−i​ 6.512×10−4−1.020×10−4−i​ 4.355×10−4),\displaystyle\left(\begin{array}[]{c}-9.286\times 10^{-5}+i\ 2.593\times 10^{-5}\\ -5.753\times 10^{-4}-i\ 6.512\times 10^{-4}\\ -1.020\times 10^{-4}-i\ 4.355\times 10^{-4}\\ \end{array}\right), (44)
𝒰α​𝟔\displaystyle\mathbf{\mathcal{U}_{\alpha 6}} =\displaystyle= (2.593×10−5+i​ 9.286×10−5−6.512×10−4+i​ 5.753×10−4−4.355×10−4+i​ 1.020×10−4).\displaystyle\left(\begin{array}[]{c}2.593\times 10^{-5}+i\ 9.286\times 10^{-5}\\ -6.512\times 10^{-4}+i\ 5.753\times 10^{-4}\\ -4.355\times 10^{-4}+i\ 1.020\times 10^{-4}\\ \end{array}\right). (48)

Following Eq.s (5, 6), this point gives a total DM abundance of

Ω​h2=0.1173,\Omega h^{2}=0.1173, (49)

composed by ΩFI​h2=0.1081\Omega^{\textrm{\tiny FI}}h^{2}=0.1081 from freeze-in decay, and ΩDW​h2=0.0092\Omega^{\textrm{\tiny DW}}h^{2}=0.0092 from Dodelson-Widrow production.

One can notice in above equations the approximate lepton-number conserving structure, with M1≪M2,3M_{1}\ll M_{2,3}, M2∼M3M_{2}\sim M_{3}, mD1≪mD2,3{m_{D}^{1}}\ll{m_{D}^{2,3}} and mD3≃i​mD2m_{D}^{3}\simeq i\ m_{D}^{2}.

We stress that when numerically diagonalising such kind of solutions it is important to use a sufficiently large precision, since in a symmetry-protected scenario the sub-eV neutrino masses can result from cancellations between terms at much larger energy scales (such cancellations are indeed imposed by the symmetry itself).

References