# A seesaw model for large neutrino masses in concordance with cosmology

Miguel Escudero<EMAIL_ADDRESS>Theoretical Physics Department, CERN,
1211 Geneva 23, Switzerland  Thomas Schwetz<EMAIL_ADDRESS>Institut für
Astroteilchenphysik, Karlsruher Institut für Technologie (KIT), Hermann-von-
Helmholtz-Platz 1, 76344 Eggenstein-Leopoldshafen, Germany  Jorge Terol-Calvo
<EMAIL_ADDRESS>Instituto de Astrofísica de Canarias, C/ Vía Láctea, s/n E38205
- La Laguna, Tenerife, Spain  Universidad de La Laguna, Dpto. Astrofísica,
E38206 - La Laguna, Tenerife, Spain

(August 27, 2024)

###### Abstract

Cosmological constraints on the sum of the neutrino masses can be relaxed if
the number density of active neutrinos is reduced compared to the standard
scenario, while at the same time keeping the effective number of neutrino
species $N_{\rm eff}\approx 3$ by introducing a new component of dark
radiation. We discuss a UV complete model to realise this idea, which
simultaneously provides neutrino masses via the seesaw mechanism. It is based
on a $U(1)$ symmetry in the dark sector, which can be either gauged or global.
In addition to heavy seesaw neutrinos, we need to introduce $\mathcal{O}(10)$
generations of massless sterile neutrinos providing the dark radiation. Then
we can accommodate active neutrino masses with $\sum m_{\nu}\sim 1$ eV, in the
sensitivity range of the KATRIN experiment. We discuss the phenomenology of
the model and identify the allowed parameter space. We argue that the gauged
version of the model is preferred, and in this case the typical energy scale
of the model is in the 10 MeV to few GeV range.

††preprint: CERN-TH-2022-180

###### Contents

  1. 1 Introduction
  2. 2 Review of the mechanism
  3. 3 A Seesaw Model for large neutrino masses and dark radiation
    1. 3.1 Neutrino mixing
    2. 3.2 Global $U(1)_{X}$
    3. 3.3 Gauged $U(1)_{X}$
    4. 3.4 Parameter summary
  4. 4 Viable parameter space of the model
  5. 5 Discussion
    1. 5.1 Heavy right-handed neutrinos, $\Delta N_{\rm eff}$, and leptogenesis
    2. 5.2 Active neutrino decay
  6. 6 Conclusions
  7. A The Mechanism after recombination
  8. B Entropy production and the value of $N_{\rm eff}$
  9. C Decay Rates

##  1 Introduction

Neutrinos play a crucial role during various stages of the cosmological
evolution, see e.g., Dolgov (2002); Lesgourgues and Pastor (2006) for reviews.
In particular, non-zero neutrino masses will influence cosmological structure
formation. Indeed, within the standard $\Lambda$CDM model, cosmology provides
a tight bound on the sum of neutrino masses. From combined CMB and BAO
observations, the Planck collaboration obtains Aghanim  _et al._ (2020):

$\sum m_{\nu}\equiv\sum_{i=1}^{3}m_{i}<0.12\,{\rm eV}\,(95\%\,{\rm CL})\ ,$
(1.1)

where $m_{i}$ are the masses of the three neutrino mass states. Depending on
the specific cosmological data used, even stronger limits can be obtained, see
e.g. Di Valentino  _et al._ (2021). With data from ongoing/upcoming large-
scale structure surveys by DESI Aghamousa  _et al._ (2016) and Euclid Amendola
_et al._ (2018), sensitivities to $\sum m_{\nu}$ of $0.02\,\mathrm{eV}$ could
be achieved in the very near future, see e.g. Brinckmann  _et al._ (2019);
Chudaykin and Ivanov (2019). From neutrino oscillation data, a minimal value
of $\sum m_{\nu}\approx 0.06\,\mathrm{eV}$ for the normal neutrino mass
ordering and $0.1\,\mathrm{eV}$ for inverted ordering is predicted for $m_{\rm
lightest}=0$. Hence, we can expect a positive detection of non-zero neutrino
mass from cosmology in $\Lambda$CDM in the next 5–7 years.

On the other hand, the search for the absolute neutrino mass scale is also one
of the top priorities in terrestrial experimental physics. The best limit on
the kinematical mass relevant for beta decay comes currently from the KATRIN
experiment Aker  _et al._ (2019, 2022) and reads

$m_{\beta}\equiv\left(\sum_{i=1}^{3}|U_{ei}|^{2}m_{i}^{2}\right)^{1/2}<0.8\,{\rm
eV}\,(90\%\,{\rm CL})\,,$ (1.2)

where $U_{ei}$ is the mixing matrix element of the mass state $i$ with the
electron neutrino. The final sensitivity goal of KATRIN is 0.2 eV for
$m_{\beta}$. For these mass ranges, neutrinos are quasi-degenerate and the
quoted current limit and sensitivity correspond to $\sum m_{\nu}\approx 2.4$
eV and 0.6 eV, respectively, already excluded by the cosmological bound from
eq. (1.1).

If neutrinos are Majorana particles, they will induce neutrinoless double-beta
decay, see Ref. Agostini  _et al._ (2022) for a review. In the absence of
cancellations due to other exotic physics, the corresponding decay rate is
proportional to an effective Majorana mass:

$m_{\beta\beta}=\left|\sum_{i=1}^{3}m_{i}U_{ei}^{2}\right|\ .$ (1.3)

The current constraint from the KamLAND-Zen experiment Gando  _et al._ (2016)
reads

$m_{\beta\beta}<0.061-0.165\,{\rm eV}\quad(90\%\,{\rm CL})\ ,$ (1.4)

where the indicated range corresponds to the uncertainty from nuclear matrix
elements. Comparable limits are obtained in Agostini  _et al._ (2020); Adams
_et al._ (2021); Anton  _et al._ (2019), and there is strong ongoing
experimental effort to reach sensitivities in the range $m_{\beta\beta}\approx
0.01-0.02$ eV Giuliani  _et al._ (2019); Agostini  _et al._ (2022). Depending
the the assumed nuclear matrix element and the unknown complex Majorana phases
in $U_{ei}$, the bound in eq. (1.4) corresponds to values of $\sum
m_{\nu}\gtrsim 0.6$ eV, largely excluded by the cosmological bound eq. (1.1).

Hence, if the $\Lambda$CDM bound is taken at face value, terrestrial neutrino
mass experiments are expected to obtain only upper limits for the neutrino
mass, and detection prospects for neutrino-less double beta decay are
challenging (depending on the neutrino mass ordering, see Agostini  _et al._
(2022); Ettengruber  _et al._ (2022)). In view of this situation, several non-
standard scenarios have been considered in the literature, to relax
cosmological neutrino mass bounds and make cosmology consistent with “large”
(i.e., observable) neutrino masses. Such scenarios include neutrino decays
Chacko  _et al._ (2020); Escudero  _et al._ (2020); Escudero and Fairbairn
(2019); Barenboim  _et al._ (2021); Chacko  _et al._ (2021); Abellán  _et al._
(2021); Chen  _et al._ (2022), neutrinos with a time-varying mass Dvali and
Funcke (2016); Dvali  _et al._ (2021); Lorenz  _et al._ (2019, 2021); Esteban
and Salvado (2021), neutrinos with a temperature much lower than the thermal
one supplemented with dark radiation Farzan and Hannestad (2016); Renk  _et
al._ (2021), and neutrinos with a momentum distribution function that deviates
from the canonical Fermi-Dirac distribution Cuoco  _et al._ (2005); Oldengott
_et al._ (2019); Alvey  _et al._ (2022a). For a discussion of phenomenological
consequences of some of these scenarios, see Alvey  _et al._ (2022b).

In this paper we focus on a mechanism similar to the one introduced by Farzan
and Hannestad in Farzan and Hannestad (2016), where the number density of
massive (active) neutrinos is reduced while simultaneously populating a new
component of dark radiation, see also Renk  _et al._ (2021). The goal of our
work is to realise this mechanism within a UV complete model which
simultaneously provides a model for neutrino masses. The model is based on a
variant of the type-I seesaw mechanism and is sometimes called Minimal
Extended Seesaw Chun  _et al._ (1995); Barry  _et al._ (2011); Zhang (2012);
Heeck and Zhang (2013); Ballett  _et al._ (2019); Escudero  _et al._ (2020).

The outline of our paper is as follows. In section 2 we review the mechanism
of Ref. Farzan and Hannestad (2016) and highlight some differences in our
realisation compared to the original paper. In section 3 we introduce the
model and the relevant parameters, where we discuss two versions, based either
on a global or a local $U(1)$ symmetry. The phenomenology is worked out in
section 4, where we study various cosmological and astrophysical constraints
and derive the available parameter space of the model. In section 5 we provide
further discussions of the viable parameter space and mention various
phenomenological and cosmological consequences and predictions. We draw our
conclusions in section 6. The interested reader can find supplementary
material in the appendices: in App. A we comment on the possibility that the
mechanism is active after recombination; in App. B we consider in detail the
effective number of neutrino species within this mechanism; in App. C we
collect a few useful formulae.

##  2 Review of the mechanism

Cosmological observations are not directly sensitive to the neutrino mass, but
rather to the energy density in neutrinos and its evolution,
$\Omega_{\nu}(z)$. As usual, $\Omega_{\nu}$ denotes the energy density
relative to the critical density. Once a cosmological model is specified there
is a direct connection between $\Omega_{\nu}(z)$ and $\sum m_{\nu}$ and a
bound on the neutrino mass can be placed. The energy density in neutrinos is
$\rho_{\nu}=\sum\left<E_{\nu}\right>n_{\nu}$, where $n_{\nu}$ denotes number
density and the sum is over the three neutrino mass states. The sensitivity to
the energy density in neutrinos in cosmology appears in two regimes, when they
are relativistic and when they are non-relativistic. When neutrinos are
ultrarelativistic $\left<E_{\nu}\right>\simeq\left<p_{\nu}\right>$, and the
contribution of neutrinos to the energy density is parameterised by the number
of effective ultrarelativistic neutrino species:

$\displaystyle N_{\rm
eff}\equiv\frac{8}{7}\left[\frac{11}{4}\right]^{4/3}\left[\frac{\rho_{\rm
rad}-\rho_{\gamma}}{\rho_{\gamma}}\right]\,,$ (2.1)

where $\rho_{\rm rad}$ is the energy density in relativistic species. We can
clearly see that when only neutrinos and photons are present $N_{\rm
eff}\propto\left<p_{\nu}\right>n_{\nu}$ Lesgourgues and Pastor (2006). This
quantity, at the time of recombination is measured to be Aghanim  _et al._
(2020):

$\displaystyle N_{\rm eff}=2.99\pm 0.17\,,$ (2.2)

which is in good agreement with the Standard Model prediction of $N_{\rm
eff}^{\rm SM}=3.044(1)$ Escudero Abenza (2020); Akita and Yamaguchi (2020);
Froustey  _et al._ (2020); Bennett  _et al._ (2021).

When neutrinos become non-relativistic, which happens when
$\left<p_{\nu}\right>\simeq m_{\nu}$, neutrinos start contributing to the
expansion rate as dark matter and their energy density is given by
$\rho_{\nu}=\sum m_{\nu}n_{\nu}$. It has been explicitly shown that current
cosmological observations are insensitive to the exact distribution function
of neutrinos Alvey  _et al._ (2022a). This means that CMB observations can
only place a bound on the energy density in non-relativistic neutrinos which
Planck CMB data constrains to be Aghanim  _et al._ (2020); Alvey  _et al._
(2022a):

$\displaystyle\Omega_{\nu}h^{2}<1.3\times 10^{-3}\quad[{\rm 95\%\,CL}]\,.$
(2.3)

Since $\Omega_{\nu}h^{2}\equiv{\sum m_{\nu}n_{\nu}^{0}}h^{2}/{\rho_{\rm
crit}}$ this bound can be seen as a bound on the product of the neutrino mass
and the neutrino number density today:

$\displaystyle\sum m_{\nu}\times\left[\frac{n_{\nu}^{0}}{56\,{\rm
cm^{-3}}}\right]<0.12\,{\rm eV}\quad[{\rm 95\%\,CL}]\,,$ (2.4)

where $n_{\nu}^{0}$ refers to the background number density of neutrinos today
per helicity state, which in the Standard Model is $n_{\nu}^{0}\simeq 56\,{\rm
cm^{-3}}$ Lesgourgues and Pastor (2006).

Eq. (2.4) highlights a way to relax the cosmological neutrino mass bound.
Since what is constrained is a product of number density and mass, reducing
the number density of neutrinos would relax the neutrino mass bound
accordingly. Importantly, since $N_{\rm
eff}\propto\left<p_{\nu}\right>n_{\nu}$, if one reduces the number density of
neutrinos $N_{\rm eff}$ will decrease, but from eq. (2.2) we see that $N_{\rm
eff}$ measurements are compatible with the Standard Model prediction. This
means that if one wants to reduce the neutrino number density before
recombination one should also add new light or massless species beyond the
Standard Model to compensate for the decrease of $N_{\rm eff}$ due to the
decrease of $n_{\nu}$. This was precisely the idea of Farzan and Hannestad in
Farzan and Hannestad (2016). For this mechanism to work, both the reduction of
the neutrino number density and the addition of new massless dark radiation
should happen before recombination111In Appendix A we study the possibility of
actually realizing the mechanism after recombination. We show that while it is
in principle possible the regions of parameter space is significantly more
restricted than if the mechanism operates before recombination.. In addition,
this should certainly happen after proton to neutron conversions have frozen
out in the early Universe (around $T_{\gamma}\sim 0.7\,{\rm MeV}$), because
otherwise the successful predictions of Big Bang Nucleosynthesis (BBN) will be
spoiled. Nevertheless, since CMB observations are only sensitive to the
Universe’s evolution at $z\lesssim 2\times 10^{5}$, see e.g. Allahverdi  _et
al._ (2020), or equivalently $T_{\gamma}\lesssim 10\,{\rm eV}$, there is
plenty of time for this to happen.

Farzan and Hannestad Farzan and Hannestad (2016) pointed out a way to achieve
the two requirements outlined above: have a large number, $N_{\chi}$, of
massless particles that thermalise with neutrinos after BBN but before
recombination, at $10\,{\rm eV}\lesssim T_{\gamma}\lesssim 100\,{\rm keV}$.
Since after neutrino decoupling at $T_{\gamma}\lesssim 2\,{\rm MeV}$ neutrinos
do not interact with the Standard Model plasma, neutrinos cannot be produced
anymore and therefore the production of new particles will be at the expense
of neutrinos. In this case, the number of effective relativistic neutrino
species in the early Universe is almost unchanged from its SM value $N_{\rm
eff}\simeq 3$222In Farzan and Hannestad (2016) it was mentioned that $N_{\rm
eff}$ does not change in this mechanism. However, the production of particles
out of equilibrium always leads to some entropy generation which does indeed
make $N_{\rm eff}$ slightly larger than 3.044. The small difference is however
negligible for practical purposes, see Appendix B for more details., but the
number density decreases and the current cosmological neutrino mass bound
becomes:

$\displaystyle\sum m_{\nu}<0.12\,{\rm eV}\,(1+g_{\chi}N_{\chi}/6)\quad[{\rm
95\%\,CL}]\,.$ (2.5)

Here $g_{\chi}$ corresponds to the number of internal degrees of freedom of
the massless BSM particle $\chi$ per species and $N_{\chi}$ is the number of
species. Fig. 1 explicitly shows the number of new massless species needed to
relax the cosmological neutrino mass bound as a function of the true neutrino
mass for $g_{\chi}=4$, as this is the case for the most relevant model of the
two we will present later. We see that for example, for the case of $\sum
m_{\nu}=0.6\,{\rm eV}$ (which is the sensitivity limit of KATRIN),
$N_{\chi}\sim 6$ would be needed to avoid the current Planck bound.

Figure 1: Number of massless fermions $\chi$ with $g_{\chi}=4$ degrees of
freedom needed to make the standard cosmological bound (shown on the
horizontal axis) consistent with different values of the sum of the neutrino
masses $\sum m_{\nu}$. The vertical dashed lines indicate the current
cosmological bound from Planck+BAO data, eq. (1.1), and the prospect for
future cosmological observations (0.02 eV). Figure 2: Illustration of the
mechanism of Farzan and Hannestad Farzan and Hannestad (2016) to reduce the
neutrino number density between BBN and recombination. We show the relative
number densities of active neutrinos (red), $N_{\chi}=10$ generations of
massless sterile fermions (blue), and the mediator boson $X$ with mass
$m_{X}=1$ keV (purple). For reference we show relevant events taking place in
the early Universe, see e.g. Pospelov and Pradler (2010), as well as the
region of temperatures at which neutrinos (or other relativistic species)
should be freestreaming Taule  _et al._ (2022).

An important question in this mechanism is how can neutrinos thermalise with a
large number of new massless BSM species between BBN and recombination. Ref.
Farzan and Hannestad (2016) considers resonantly enhanced scattering between
neutrinos and these new species, mediated by a new boson $X$ with a mass
$10\,{\rm eV}\lesssim m_{X}\lesssim 100\,{\rm keV}$ via the process
$\bar{\nu}\nu\to(X)\to\chi\chi$ where $\chi$ here represents one of the
massless states. This requirement can actually be relaxed and what is really
minimally required is that the new boson thermalises with neutrinos and that
it interacts efficiently with a large number of massless species beyond the
Standard Model (SM). Thus, the two requirements for this to work are:

$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!{\rm
1)}\quad\quad\left<\Gamma(\bar{\nu}\nu\to X)\right>>H\,,$ (2.6)

and

$\displaystyle{\rm 2)}\quad\quad\left<\Gamma(X\to\chi_{i}+{\rm
anything})\right>>H\,,$ (2.7)

both for $10\,{\rm eV}\lesssim T\lesssim 100\,{\rm keV}$. To illustrate the
mechanism and its main ingredients we show in fig. 2 the evolution of neutrino
and dark-sector particle densities as a function of photon temperature. For
the parameters chosen in the plot, the bound on the sum of neutrino masses can
be relaxed to 0.9 eV.

##  3 A Seesaw Model for large neutrino masses and dark radiation

In this section we discuss a specific model realisation of the mechanism
described in the previous section, which in addition provides a framework to
generate neutrino masses, following closely the discussion of Ref. Escudero
_et al._ (2020), section 4. The beyond-SM ingredients of the model are:

  * •

three fermion singlets $N_{R}$ (“right-handed neutrinos”) which play the usual
role to generate active neutrino masses as in the type-I seesaw,

  * •

a new abelian symmetry $U(1)_{X}$ which can be either global or local,

  * •

a scalar $\Phi$ with $U(1)_{X}$ charge $+1$, and

  * •

a set of $N_{\chi}$ fermions $\chi$ with $U(1)_{X}$ charge $-1$.

With these assignments we can write the following BSM terms in the Lagrangian:

$-\mathcal{L}=\overline{N_{R}}\,Y_{\nu}\,\ell_{L}\,\widetilde{H}^{\dagger}+\frac{1}{2}\,\overline{N_{R}}\,M_{R}\,N_{R}^{c}+\overline{N_{R}}Y_{\Phi}\,\chi_{L}\,\Phi+\,\text{h.c.}\,.$
(3.1)

Here $H$ and $\ell_{L}$ are the SM Higgs and lepton doublets, respectively,
and $\widetilde{H}=i\,\tau_{2}\,H^{*}$, $M_{R}$ is the $3\times 3$ Majorana
mass matrix for $N_{R}$, and $Y_{\nu}$ and $Y_{\Phi}$ are $3\times 3$ and
$3\times N_{\chi}$ Yukawa matrices, respectively. As we are interested in
“large” neutrino masses, possibly in the quasi-degenerate regime, we need 3
right-handed neutrinos $N_{R}$333We note that the mixing pattern of very
degenerate neutrinos is particularly sensitive to radiative corrections Georgi
and Glashow (2000); Ellis and Lola (1999); Casas  _et al._ (1999). In specific
flavor models this poses constraints on the scale of the origin of neutrino
masses, see e.g. for some constructions de Medeiros Varzielas  _et al._
(2009); Boudjemaa and King (2009). . Here and in the following we keep
$SU(2)_{L}$ and flavour indices contractions implicit. The scalar potential is

$\displaystyle V$
$\displaystyle=\mu^{2}_{H}H^{\dagger}H+\lambda_{H}\left(H^{\dagger}H\right)^{2}$
$\displaystyle+\mu_{\Phi}^{2}|\Phi|^{2}+\lambda_{\Phi}|\Phi|^{4}+\lambda_{H\Phi}|\Phi|^{2}H^{\dagger}H\,,$
(3.2)

with $\mu^{2}$ and $\mu_{\Phi}^{2}$ parameters with dimensions of [mass]2 and
$\lambda_{H},\lambda_{\Phi},\lambda_{H\Phi}$ dimensionless. We assume
$\lambda_{H\Phi}=0$, i.e., no mixing between the two scalar fields. With this
assumption we avoid that $\Phi$ gets thermalised in the early Universe due to
its interactions with the SM Higgs. Electroweak symmetry breaking takes place
in the usual way, with

$\langle H\rangle=\frac{1}{\sqrt{2}}\,\left(\begin{array}[]{c}0\\\ v_{\rm
EW}\end{array}\right)\,,$ (3.3)

with $v_{\rm EW}\simeq 246$ GeV denoting the SM Higgs vacuum expectation value
(VEV). The breaking of the $U(1)_{X}$ takes place when $\Phi$ develops a VEV

$\langle\Phi\rangle=\dfrac{v_{\Phi}}{\sqrt{2}}\,,$ (3.4)

with $v_{\Phi}^{2}=-\mu_{\Phi}^{2}/\lambda_{\Phi}$.

###  3.1 Neutrino mixing

After symmetry breaking, several terms in the Yukawa Lagrangian in eq. (3.1)
induce mixing in the neutral lepton sector. In the basis
$n=\left(\nu_{L}^{c},N_{R},\chi_{L}^{c}\right)$, the fermion mass terms can be
written as

$-\mathcal{L}_{m}=\frac{1}{2}\overline{n^{c}}\,\mathcal{M}_{n}\,n+\,\text{h.c.}\,,$
(3.5)

with the $(6+N_{\chi})\times(6+N_{\chi})$ mass matrix given by

$\mathcal{M}_{n}=\begin{pmatrix}0&m_{D}&0\\\ m_{D}^{T}&M_{R}&\Lambda\\\
0&\Lambda^{T}&0\end{pmatrix},$ (3.6)

where $m_{D}=\frac{v_{\rm EW}}{\sqrt{2}}\,Y_{\nu}$ and
$\Lambda=\frac{v_{\Phi}}{\sqrt{2}}\,Y_{\Phi}$. We assume the following
hierarchy between the entries of the mass matrix:

$\Lambda\ll m_{D}\ll M_{R}\,,$ (3.7)

where these relations are understood for the typical scales relevant for the
matrices.

The block-diagonalisation of the mass matrix leads to the masses of the 3
active neutrinos, the 3 heavy neutrinos and the $N_{\chi}$ massless sterile
neutrinos

$\mathcal{M}_{n}^{D}=\begin{pmatrix}m_{\rm active}&0&0\\\ 0&m_{\rm heavy}&0\\\
0&0&m_{\rm sterile}\end{pmatrix},$ (3.8)

with

$\displaystyle m_{\rm active}$ $\displaystyle\approx
m_{D}\,M_{R}^{-1}\,m_{D}^{T}+\Lambda\,\Lambda^{T}\,M_{R}^{-1}\approx
m_{D}\,M_{R}^{-1}\,m_{D}^{T},$ $\displaystyle m_{\rm heavy}$
$\displaystyle\approx
M_{R}+m_{D}\,M_{R}^{-1}\,m_{D}^{T}+\Lambda\,\Lambda^{T}\,M_{R}^{-1}\approx
M_{R},$ $\displaystyle m_{\rm sterile}$ $\displaystyle=0,$ (3.9)

where $m_{\rm active}=U_{\nu}^{*}\,\widehat{m}_{\nu}\,U_{\nu}^{\dagger}$.
Adopting the diagonal mass basis for charged lepton, $U_{\nu}$ is the PMNS
mixing matrix, given in terms of $3$ mixing angles and $3$ CP-violating phases
(including Majorana phases), while $\widehat{m}_{\nu}$ is a diagonal matrix
containing the physical neutrino mass eigenvalues $m_{i}$. There are
$N_{\chi}$ states which are exactly massless at tree level, due to the rank of
the matrix (3.6). Loop contributions to $m_{\rm sterile}$ are small enough to
consider the $N_{\chi}$ states effectively massless Escudero  _et al._ (2020).

The mass basis is obtained by rotating the fields with the unitary matrix $W$
which induces a mixing between the different states:

$\begin{pmatrix}\tilde{\nu}\\\ \tilde{N}\\\
\tilde{\chi}\end{pmatrix}=W^{\dagger}\begin{pmatrix}\nu_{L}^{c}\\\ N_{R}\\\
\chi_{L}^{c}\end{pmatrix}\,,$ (3.10)

where we have introduced the notation $\tilde{\nu},\,\tilde{N},\,\tilde{\chi}$
to denote the active neutrinos, heavy neutrinos, massless sterile neutrino in
the mass basis, respectively. Following e.g., Grimus and Lavoura (2000) one
can find the mixing matrix at leading order, taking into account the hierarchy
in eq. (3.7):

$\small
W=\begin{pmatrix}1&m_{D}^{*}\,(M_{R}^{-1})^{\dagger}&-(m_{D}^{-1})^{T}\,\Lambda\\\
-M_{R}^{-1}\,m_{D}^{T}&1&0\\\
\Lambda^{\dagger}\,(m_{D}^{-1})^{*}&0&1\end{pmatrix}\begin{pmatrix}U_{\nu}&0&0\\\
0&1&0\\\ 0&0&1\end{pmatrix}.$ (3.11)

Without loss of generality, we have adopted a basis where the right-handed
neutrino mass matrix $M_{R}$ is diagonal.

In order to simplify the discussion, we will adopt below the one-flavour
approximation for the active and heavy neutrinos and introduce mixing angles

$\displaystyle\theta_{\nu
N}=\frac{m_{D}}{M_{R}}\,,\quad\theta_{\nu\chi}=\frac{\Lambda}{m_{D}}\,,$
(3.12)

describing the mixing between active neutrinos and the heavy and massless
states, respectively. With our assumption eq. (3.7), both angles are small. We
need to keep $N_{\chi}$ flavors of massless sterile states and
$\theta_{\nu\chi}$ represents the mixing between each of them and the active
neutrinos. Finally, using the seesaw relation
$m_{\nu}=m_{D}^{2}/M_{R}=\theta_{\nu N}^{2}M_{R}$ we will eliminate $m_{D}$
(or $\theta_{\nu N}$) and $\Lambda$ and consider $m_{\nu}$, $M_{R}$ and
$\theta_{\nu\chi}$ as independent parameters.

In the following we discuss the relevant interaction terms and distinguish the
particularities of the global and gauged versions of the model.

###  3.2 Global $U(1)_{X}$

Let us decompose the complex scalar $\Phi$ into two real fields as
$\Phi=\dfrac{1}{\sqrt{2}}\left(v_{\Phi}+\rho+i\phi\right)$, where we take
$v_{\Phi}$ real without loss of generality. The real part $\rho$ has a mass
$m_{\rho}$ of order $|\mu_{\Phi}|$, while $\phi$ corresponds to the Goldstone
boson. We assume that in addition to the spontaneously breaking of the
$U(1)_{X}$ global symmetry also explicit breaking terms are present, e.g.
arising from higher-dimensional terms of the scalar potential, inducing a mass
term for the imaginary part $\phi$. Hence, the pseudo-Goldstone mass
$m_{\phi}$ is an additional independent parameter in the global version of the
model.

The relevant processes for our mechanism are $X\leftrightarrow\nu\,\nu$ and
$X\leftrightarrow\nu\,\chi$, where for the global case $X$ can be the scalar
$\rho$ or the pseudoscalar $\phi$. These interactions arise from the third
term in eq. (3.1) through the mixing of the neutral particles
$\nu_{L},\,N_{R}$ and $\chi_{L}$. In the mass basis and after Spontaneous
Symmetry Breaking (SSB) we have for the interaction of the scalars with two
active neutrinos

$\displaystyle\overline{N_{R}}Y_{\Phi}\,\chi_{L}\,\Phi+\text{h.c.}\supset-\overline{\tilde{\nu}}\,\lambda_{\rho/\phi}^{\nu\nu}\dfrac{1}{\sqrt{2}}(\rho-i\gamma_{5}\phi)\,\tilde{\nu}^{c}+\text{h.c.}\,$
(3.13)

with the coupling

$\displaystyle\lambda_{\rho/\phi}^{\nu\nu}$
$\displaystyle=\dfrac{1}{\sqrt{2}}U_{\nu}^{\dagger}\,m_{D}^{*}\,\left(M_{R}^{-1}\right)^{\dagger}\,Y_{\Phi}\,\Lambda^{T}m_{D}^{-1}\,U_{\nu}^{*}$
$\displaystyle=\frac{1}{v_{\Phi}}\,\widehat{m}_{\nu}\,U_{\nu}^{T}\,\left(m_{D}^{-1}\right)^{\dagger}\,\Lambda\,\Lambda^{T}m_{D}^{-1}\,U_{\nu}^{*}$
$\displaystyle\to\frac{m_{\nu}}{v_{\Phi}}\,\theta_{\nu\chi}^{2}\,.$ (3.14)

We have used that
$\widehat{m}_{\nu}=U_{\nu}^{T}\,m_{D}\,M_{R}^{-1}\,m_{D}^{T}\,U_{\nu}$ and
$\Lambda=Y_{\Phi}\,v_{\Phi}/\sqrt{2}$. The last line holds in the one-flavour
approximation. An analogous calculation for the coupling
$\rho/\phi\leftrightarrow\nu\,\chi$ leads to (again in the one-flavour
approximation)
$\lambda_{\rho/\phi}^{\nu\chi}=\theta_{\nu\chi}m_{\nu}/v_{\Phi}$, whereas
there is no direct vertex for $\rho/\phi\leftrightarrow\chi\,\chi$, as at
leading order there is no mixing between $N_{R}$ and $\chi_{L}$. Summarizing,
in the global realisation of the model we have:

$\displaystyle\lambda_{\rho/\phi}^{\nu\nu}$
$\displaystyle=\frac{m_{\nu}}{v_{\Phi}}\,\theta_{\nu\chi}^{2}\,,$ (3.15a)
$\displaystyle\lambda_{\rho/\phi}^{\nu\chi}$
$\displaystyle=\frac{m_{\nu}}{v_{\Phi}}\,\theta_{\nu\chi}\,,$ (3.15b)
$\displaystyle\lambda_{\rho/\phi}^{\chi\chi}$ $\displaystyle=0\,.$ (3.15c)

We observe that the couplings $\lambda_{\rho/\phi}^{\nu\chi}$ and
$\lambda_{\rho/\phi}^{\nu\nu}$ are suppressed by the ratio $m_{\nu}/v_{\Phi}$
as well as one or two powers of the active-sterile mixing $\theta_{\nu\chi}$,
respectively.

###  3.3 Gauged $U(1)_{X}$

Let us now consider the case of $U(1)_{X}$ being a gauge symmetry. In this
case, besides the particle content of the global case above we have to add
another set of $N_{\chi}$ fermions charged under the new symmetry with the
same charge but opposite chirality, $\chi_{R}$, in order to cancel the gauge
anomaly introduced by the $\chi_{L}$. With the introduction of $\chi_{R}$, two
new terms arise in the Lagrangian, an interaction with the singlet $N_{R}$ and
the scalar $\Phi$ analogous to the last term of eq. (3.1), and a vector-like
mass term for the $\chi$ field. These terms would change the picture of the
neutrino mass generation. These terms can be forbidden by postulating a
discrete $\mathbb{Z}_{2}$ symmetry under which all particles are even except
$\chi_{R}$. Hence, in the following we will assume that these terms are
absent.

The breaking of the $U(1)_{X}$ will give mass to the associated gauge boson,
$Z^{\prime}$, and the would-be Goldstone boson $\phi$ becomes the longitudinal
polarisation of the $Z^{\prime}$. We assume no tree-level mixing between the
SM $U(1)_{\rm em}$ and the $U(1)_{X}$ gauge fields; there will still be loop
contributions Gherghetta  _et al._ (2019); Bauer and Foldenauer (2022) though
negligible for phenomenology in our case as there is no particle content
charged simultaneously under $U(1)_{Y}$ or $SU(2)_{L}$ and $U(1)_{X}$.

The interactions of fermions with the new gauge boson are described by

$\mathcal{L}=\sum_{f=\chi_{L},\chi_{R}}\,Q_{f}\,g_{X}\,Z^{\prime}_{\mu}\,\bar{f}\,\gamma^{\mu}\,f,$
(3.16)

where $g_{X}$ is the $U(1)_{X}$ gauge coupling and $f$ are the $N_{\chi}$
fermions charged under the $U(1)_{X}$ with charges
$Q_{\chi_{L}}=Q_{\chi_{R}}=-1$. For the interaction of active neutrinos with
the $Z^{\prime}$ we have

$\displaystyle
g_{X}\,Z^{\prime}_{\mu}\,\overline{\chi_{L}}\,\gamma^{\mu}\,\,\chi_{L}+\text{h.c.}\supset$
(3.17)
$\displaystyle\dfrac{m_{Z^{\prime}}}{v_{\Phi}}\,Z^{\prime}_{\mu}\,\overline{\tilde{\nu}^{c}}\,U_{\nu}^{T}\,\left(m_{D}^{-1}\right)^{\dagger}\,\Lambda^{*}\gamma^{\mu}\,\Lambda^{T}\left(m_{D}^{-1}\right)\,U_{\nu}^{*}\,\tilde{\nu}^{c}+\text{h.c.},$

where we have used $m_{Z^{\prime}}=g_{X}\,v_{\Phi}$. In analogy to the global
case we introduce couplings for the
$Z^{\prime}\leftrightarrow\nu\,\nu/\nu\,\chi/\chi\,\chi$ interactions in the
one-flavour approximation:

$\displaystyle\lambda_{Z^{\prime}}^{\nu\nu}$
$\displaystyle=\frac{m_{Z^{\prime}}}{v_{\Phi}}\,\theta_{\nu\chi}^{2}\,,$
(3.18a) $\displaystyle\lambda_{Z^{\prime}}^{\nu\chi}$
$\displaystyle=\frac{m_{Z^{\prime}}}{v_{\Phi}}\,\theta_{\nu\chi}\,,$ (3.18b)
$\displaystyle\lambda_{Z^{\prime}}^{\chi\chi}$
$\displaystyle=\frac{m_{Z^{\prime}}}{v_{\Phi}}\,.$ (3.18c)

In comparison to eq. (3.15), here also a $Z^{\prime}\leftrightarrow\chi\,\chi$
vertex is present. Note that these interactions cannot be induced by the
$\chi_{R}$ gauge interaction, since there is no mixing between these particles
and the active neutrinos under our assumption of a $\mathbb{Z}_{2}$ symmetry.
However, $\chi_{R}$ will contribute to the mechanism described in sec. 2 as
relativistic species in thermal equilibrium with the neutrinos through the
process $\nu\,\nu\leftrightarrow
Z^{\prime}\leftrightarrow\overline{\chi_{R}}\,\chi_{R}$. That means that in
this case $g_{\chi}=4$ in eq. (2.5), and therefore requiring a smaller
$N_{\chi}$ compared to the global case. Note that while the interaction
$\phi/\rho\to\nu\,\chi$ is crucial for the mechanism in the global case to
create the new light species, in the gauge case the $Z^{\prime}\to\chi\,\chi$
interaction will be more efficient to do so than the analogous with
$Z^{\prime}$ due to the suppression with $\theta_{\nu\chi}$, see eqs. (3.18).

The real scalar $\rho$ is also present in the gauged version and will induce
scalar-mediated interactions with $\lambda_{\rho}^{\nu\nu}$,
$\lambda_{\rho}^{\nu\chi}$ according to eqs. (3.15a), (3.15b). In fact, it
thermalises, however, for values in the parameter space that are ruled out by
other constraints on the interaction $\nu\,\nu\leftrightarrow Z^{\prime}$ as
we will discuss later.

###  3.4 Parameter summary

Let us now summarise the discussion above and list the independent parameters
of the model, adopting the one-flavour approximation. As discussed in sec.
3.1, the neutrino sector can be parameterised by the three parameters

$\displaystyle m_{\nu},\,M_{R},\,\theta_{\nu\chi}\,,$ (3.19)

being the active neutrino mass, the heavy right-handed neutrino mass, and the
mixing between active and massless sterile neutrinos, respectively. For
numerical estimates we fix $m_{\nu}=0.2$ eV, at the final KATRIN sensitivity.
As we will see below the value of $M_{R}$ is irrelevant for the phenomenology
of the mechanism itself. However, the value of $M_{R}$ can be constrained by
perturbativity requirements and we will be able to make predictions for the
required range of $M_{R}$, which can in turn lead to cosmological consequences
in the very early Universe. Regarding $\theta_{\nu\chi}$, we work within the
one-flavor approximation and thus consider that the mixing between all the
$\chi_{L}$ states and $\nu_{a}$ is comparable in size. As we will see, this
parameter is constrained mainly by BBN observations and we will place
constraints that take into account the number of $\chi_{L}$ states,
$N_{\chi}$.

In the scalar and interaction sector we take as independent parameters

$\displaystyle
v_{\Phi},\,m_{\rho},\,m_{X}=\left\\{\begin{array}[]{l}m_{\phi}\quad\text{(global)}\\\
m_{Z^{\prime}}\quad\\!\text{(gauge)}\end{array}\right.\,.$ (3.22)

We use the VEV and the real-scalar mass $m_{\rho}$ as independent parameters
of the scalar potential. In the global $U(1)_{X}$ case, the pseudo-Goldstone
boson mass $m_{\phi}$ is another independent parameter, whereas in the gauge
version, we trade the gauge coupling constant $g_{X}$ by taking the VEV and
the $Z^{\prime}$ mass as independent.

In summary, the most relevant parameters to be determined by phenomenology are
$\theta_{\nu\chi}$, $v_{\Phi}$ and $m_{X}$, where $X=\phi/\rho$ ($Z^{\prime}$)
in the global (gauge) symmetry case. As outlined in sections 3.2 and 3.3, all
the relevant interaction rates can be expressed in terms of these few
parameters. In what follows we study the available parameter space in terms of
these parameters.

##  4 Viable parameter space of the model

In this section we consider all relevant phenomenological and cosmological
constraints in order to narrow down the parameter space to a region in which
the mechanism can be fully realised and the bounds on the neutrino mass
significantly relaxed. Figs. 3 and 4 show the results for the global and gauge
cases, respectively, in terms of the mass of the new boson, $m_{X}$, and the
vacuum expectation value of the $\Phi$ field, $v_{\Phi}$, for two
representative values of the mixing between active and sterile states
$\theta_{\nu\chi}=10^{-3}$ and $10^{-4}$. Below we elaborate upon all of these
constraints and considerations.

Figure 3: Regions of the parameter space of the global $U(1)_{X}$ model
excluded by several cosmological bounds for a value of the mixing between
active and massless sterile neutrinos, $\theta_{\nu\chi}=10^{-3}$ (left) and
$10^{-4}$ (right). The white region is allowed. Vertical dotted black lines
correspond to the maximum $M_{R}$ value in GeV given by the requirement of
perturbativity for $Y_{\Phi}$, see eq. (5.3). The purple line indicates the
region where $m_{\phi}>v_{\Phi}$, where the explicit breaking of the
$U(1)_{X}$ symmetry by the scalar mass would dominate over the spontaneous
breaking. The vertical green line highlights parameter space excluded from
neutrino freestreaming in the specific case $m_{\phi}=0$.

Figure 4: Regions of the parameter space of the gauge $U(1)_{X}$ model
excluded by several cosmological bounds for a value of the mixing between
active and massless sterile neutrinos, $\theta_{\nu\chi}=10^{-3}$ (left) and
$10^{-4}$ (right). The white region is allowed. Dotted black lines correspond
to the maximum $M_{R}$ value in GeV given by the requirement of perturbativity
for $Y_{\Phi}$, see eq. (5.3). The grey dotted lines indicate regions of
constant value of the gauge coupling constant $g_{X}=m_{Z^{\prime}}/v_{\Phi}$.
We also indicate the region where standard thermal leptogenesis can work
(purple shading).

$\bullet$$\bullet$$\bullet$ Thermalisation: As discussed in sec. 2 the
mechanism requires the $X$ boson to thermalise with neutrinos and $\chi$
states after proton-neutron freeze-out and before recombination at $10\,{\rm
eV}\lesssim T_{\gamma}\lesssim 0.7\,{\rm MeV}$. Calculating the thermally
averaged decay rate and comparing it to the expansion rate of the Universe
allows us to narrow down the range of couplings for which this happens (see
appendix C for formulae). In our models, the rate $\Gamma(X\to\nu\chi)$ is
always larger than $\Gamma(X\to\nu{\nu})$ due to the mixing angle suppression,
see eqs. (3.15) and (3.18). However, we assume that initially no $\chi$s are
present in the plasma, and therefore the thermalisation requirement applies to
the process $X\leftrightarrow\nu{\nu}$, controlled by the coupling
$\lambda_{X}^{\nu{\nu}}$. Since the thermally average decay rate peaks at
$T\sim m_{X}/3$ we demand $\left<\Gamma(X\leftrightarrow\nu\nu)\right>\gtrsim
H(T=m_{X}/3)$, which approximately implies

$\displaystyle\lambda_{X}^{\nu{\nu}}\gtrsim 4\times
10^{-12}\,\sqrt{\frac{m_{X}}{{\rm keV}}}\,.$ (4.1)

In figs. 3 and 4 we highlight in red the region of parameter space for $m_{X}$
and $v_{\Phi}$ for which the $X$ boson will not reach thermal equilibrium with
active neutrinos in the early Universe. The thermalisation requirement leads
to an upper bound on the scalar VEV, where the dependence on $m_{X}$ follows
from combining eq. (4.1) with either (3.15a) or (3.18a).

$\bullet$$\bullet$$\bullet$ BBN Constraints on $X$–$\nu$ interactions: The new
boson $X$ cannot be in thermal equilibrium at temperatures above
$T_{\gamma}>0.7\,{\rm MeV}$ for two reasons: 1) It would reduce the number
density of neutrinos which in turn participate in $p\leftrightarrow n$
conversions, and 2) it would contribute by itself and with the new $\chi$
states to the expansion rate at the time of BBN since neutrino decoupling
happened at a very similar temperature $T\sim 2\,{\rm MeV}$. In order to
ensure that the $X$ boson does not spoil the success of BBN we can place a
bound again on the interaction rate between the $X$ boson and active
neutrinos:

$\displaystyle\left<\Gamma(\nu{\nu}\to X)\right>\lesssim H(T=0.7\,{\rm
MeV})\,,$ (4.2)

which in terms of the relevant coupling approximately reads:

$\displaystyle\lambda_{X}^{\nu{\nu}}\lesssim 10^{-7}\,\frac{{\rm
keV}}{m_{X}}\,.$ (4.3)

This constraint is shown in blue in figs. 3 and 4. In passing, we note that
the coupling $\lambda_{X}^{\nu{\nu}}$ is also bounded by possible $\nu\nu\to
XX$ processes occurring before neutrino decoupling and leading to a large
$\Delta N_{\rm eff}$. For this not to happen, then
$\lambda_{X}^{\nu{\nu}}\lesssim 10^{-5}$ Escudero and Fairbairn (2019).
However, we never probe such large couplings in the parameter space of
interest in our study.

$\bullet$$\bullet$$\bullet$ CMB Constraints on $X$–$\nu$ interactions: The
interaction between $X$ particles and neutrinos and sterile massless states
can leave an imprint on CMB observations if it occurs sufficiently close to
recombination as this would alter neutrino freestreaming and distort the CMB
power spectra. A recent model-independent analysis of Planck legacy data has
shown that provided that the $\nu$–$X$ interactions are not efficient at
$z<10^{5}$ there are no CMB constraints Taule  _et al._ (2022). We will use
this as a constraint on the parameter space, requiring that

$\displaystyle\left<\Gamma({\nu}\nu\to
X)\right><H\,\,\,\text{at}\,\,\,z<10^{5}\,.$ (4.4)

In addition, since in our scenario a large fraction of the energy density in
relativistic particles is in the form of massless sterile states $\chi$ we
will also require

$\displaystyle\left<\Gamma(X\leftrightarrow\chi+\chi/{\nu})\right><H\,\,\,\text{at}\,\,\,z<10^{5}\,.$
(4.5)

The combination of these bounds is shown in light green in figs. 3 and 4 and
essentially implies $m_{X}\gtrsim 1$ keV.

Furthermore, there are constraints coming from neutrino freestreaming from
possible $2\to 2$ processes. In particular, one needs to make sure that
$\chi\chi\leftrightarrow\chi\chi$ interactions as mediated by a $Z^{\prime}$
or $\phi$ are not efficient at $z<10^{5}$. For the scalar case this is
automatically fulfilled because $\lambda_{\phi}^{\chi\chi}=0$. However, for
the gauge case it is not since
$\lambda_{Z^{\prime}}^{\chi\chi}=m_{Z^{\prime}}/v_{\Phi}=g_{X}$. By enforcing:

$\displaystyle\left<\Gamma(\chi\chi\leftrightarrow\chi\chi)\right><H\,\,\,\text{at}\,\,\,z<10^{5}\,,$
(4.6)

we find a constraint on $v_{\Phi}$:

$\displaystyle v_{\Phi}>3\times 10^{4}\,{\rm keV}\,,$ (4.7)

for the gauge case. This bound is shown with a darker green colour in fig. 4
and restricts or even close the allowed parameter space for smaller values of
$\theta_{\nu\chi}$.

For the scalar case we still should consider the process
$\chi\nu\leftrightarrow\chi\nu$. By demanding the same requirement on this
process we find a small region of parameter space to be excluded which is
highlighted in green in the left corner of the available parameter space and
that follows a diagonal shape in fig. 3.

$\bullet$$\bullet$$\bullet$ Astrophysical considerations: Since our $X$ boson
interacts with neutrinos it can be subject to constraints from astrophysical
considerations, in particular from core collapse supernova. Supernova cores
have temperatures $T\sim 30\,{\rm MeV}$ and release almost all their binding
energy in the form of neutrinos on a timescale of $t\sim\mathcal{O}(10)\,{\rm
s}$. The observed neutrino spectrum of SN1987A is in broad agreement with that
expected from standard core collapse supernova simulations, see e.g. for a
review Janka (2012). In this context, there are two bounds one can place.
Firstly, the $X$ particle should not be copiously produced and escape on a
timescale shorter than $t\sim\mathcal{O}(10)\,{\rm s}$, otherwise the
supernova will cool much faster than what has been observed in SN1987A Raffelt
(1996). Requiring that the luminosity in $X$ states is smaller than the one
from neutrinos in the standard scenario rules out couplings in the range
Fiorillo  _et al._ (2022):

$\displaystyle 4\times 10^{-6}\,\frac{\rm
keV}{m_{X}}\lesssim\,\lambda_{X}^{\nu\nu}\,\lesssim 10^{-4}\,\frac{\rm
keV}{m_{X}}\,,$ (4.8)

for ${\rm keV}$ scale bosons. Secondly, it has been very recently pointed out
Fiorillo  _et al._ (2022) that even if the luminosity of $X$ particles emitted
by the supernova is substantially smaller than that of active neutrinos there
could still be constraints due to the lack of high energy events in the time
window where the SN1987A signal was observed. For $m_{X}<{\rm MeV}$, the
reported exclusion range corresponds to

$\displaystyle 3\times 10^{-7}\,\frac{\rm
keV}{m_{X}}\lesssim\,\lambda_{X}^{\nu\nu}\,\lesssim 10^{-4}\,\frac{\rm
keV}{m_{X}}\,.$ (4.9)

As such, we find that these constraints are weaker than the one we impose from
BBN consistency in our scenario.

$\bullet$$\bullet$$\bullet$ BBN Constraints on the mixing between active and
sterile neutrinos: The massless sterile neutrinos we consider are subject to
BBN constraints on their own because they mix with active neutrinos and
therefore can be produced via collisions and oscillations in the early
Universe. In addition, since these states are lighter than active neutrinos
they feature an enhanced resonant production Barbieri and Dolgov (1991). The
production rate for these sterile neutrinos peaks at $T\simeq 10\,{\rm
MeV}({|\Delta m^{2}|/0.1{\rm eV^{2}}})^{1/6}$ Dolgov (2002), which is well
above neutrino decoupling and BBN. The main effect of these additional states
will then be to contribute to the energy density in the Universe both at the
time of BBN and recombination. Here we take the production rate of sterile
neutrinos from Abazajian (2006) and integrate it up to the time of neutrino
decoupling, $T_{\nu}^{\rm dec}\simeq 2\,{\rm MeV}$, in order to obtain this
contribution to the number of ultrarelativistic neutrino species in the early
Universe:

$\displaystyle\left.\Delta N_{\rm eff}\right|_{\chi}\simeq
0.014\sum_{\chi=1}^{N_{\chi}}\frac{|\theta_{e\chi}|^{2}+0.8(|\theta_{\mu\chi}|^{2}+|\theta_{\tau\chi}|^{2})}{10^{-6}}\left(\frac{m_{\nu}}{0.1\,{\rm
eV}}\right)\,.$ (4.10)

Note that for simplicity to obtain this expression we have neglected
destruction of sterile neutrinos in the collision rates. By solving the
relevant Boltzmann equations we have explicitly checked that this is a good
approximation provided that $\Delta N_{\rm eff}\lesssim 0.3$ for a given new
species.

Assuming that the mixing is similar for each species and applying a bound of
$\Delta N_{\rm eff}<0.3$ which is representative of both Planck data Aghanim
_et al._ (2020) and global BBN analyses Pisanti  _et al._ (2021), we can find
a bound on $\theta_{\nu\chi}$ and $N_{\chi}$ which reads:

$\displaystyle|\theta_{\nu\chi}|\lesssim
10^{-3}\sqrt{\frac{10}{N_{\chi}}}\,\sqrt{\frac{0.2\,{\rm eV}}{m_{\nu}}}\,,$
(4.11)

where here $m_{\nu}$ refers to the mass of an individual and almost degenerate
active neutrino. In terms of the parameters of interest in our study this
means that:

$\displaystyle\theta_{\nu\chi}=\frac{\Lambda}{m_{D}}\lesssim
10^{-4}-10^{-3}\,,$ (4.12)

where the ranges are taken by varying $0.1\,{\rm eV}<m_{\nu}<1\,{\rm eV}$ and
$N_{\chi}\lesssim 50$ as relevant for a range of scenarios as seen in fig. 1.
This explains why in Figures 3 and 4 we take as benchmarks
$|\theta_{\nu\chi}|=10^{-3}$ and $|\theta_{\nu\chi}|=10^{-4}$. Choosing
smaller values of $\theta_{\nu\chi}$ would move the allowed regions to smaller
values of $v_{\Phi}$, as a seen from eqs. (3.15), (3.18). This would lead to
non-perturbative gauge couplings or similar inconsistencies in the global case
(see discussion below). For the gauged version, for mixing angles
significantly below $10^{-4}$ the allowed region above the free streaming
bound on $v_{\Phi}$, eq. (4.7), would disappear. Therefore, the preferred
parameter region is close to the upper bounds for $\theta_{\nu\chi}$ discussed
above.

With active neutrinos close to the eV scale and massless sterile neutrinos, we
obtain a mass-squared difference $\Delta m^{2}\sim 1$ eV2, potentially
relevant for short-baseline oscillation experiments Böser  _et al._ (2020);
Dentler  _et al._ (2018). However, mixing angles in the range indicated in eq.
(4.12) are too small to be tested in oscillation experiments.

##  5 Discussion

Let us summarise the main results from the various constraints discussed in
the previous section, referring to figs. 3 and 4. We find a closed region of
parameter space for the mediator mass $m_{X}$ and scalar VEV $v_{\Phi}$, where
the mechanism can work. The mediator mass is restricted roughly to the range 1
keV $\lesssim m_{X}\lesssim 1$ MeV, both in the global and gauged version of
the model. For the global case, the scalar VEV is roughly in the range 1 eV
$\lesssim v_{\Phi}\lesssim 10$ keV for $\theta_{\nu\chi}=10^{-3}$, and roughly
one order of magnitude smaller for $\theta_{\nu\chi}=10^{-4}$. In the gauged
version, we obtain larger scales for $v_{\Phi}$, ranging roughly from 10 MeV
to few GeV for the case $\theta_{\nu\chi}=10^{-3}$, but only viable for 20–200
MeV for $\theta_{\nu\chi}=10^{-4}$ as a result of the $v_{\Phi}$ dependence of
the neutrino freestreaming constraint on $\chi$ self-interactions.

The different scales for the VEV in the global and gauged versions, as well as
the dependence on the mixing angle follow from the parametric dependence of
the coupling constants shown in eqs. (3.15) and (3.18), being proportional to
$m_{\nu}/v_{\Phi}$ ($m_{Z^{\prime}}/v_{\Phi}$) in the global (gauged) version.
For the global case, a shift in the neutrino mass would have a similar effect
as changing $\theta_{\nu\chi}$, with lower masses moving the excluded regions
to smaller $v_{\Phi}$. This follows directly from eqs. (3.15).

Let us now discuss in some more detail the parameter region for the global
$U(1)_{X}$ version. In the scalar sector, the parameters of the model are the
scalar VEV $v_{\Phi}$, the mass of the real scalar, $m_{\rho}$, and the mass
of the pseudo-Goldstone boson $m_{\phi}$, which we consider as an independent
parameter responsible for the explicit breaking of the $U(1)_{X}$ symmetry.
The diagonal line in fig. 3 indicates the condition $m_{X}=v_{\Phi}$.
Considering the case $X=\phi$, we see that in large regions of the parameter
space (for $\theta_{\nu\chi}=10^{-4}$ even all of the viable parameter space),
we have $m_{\phi}>v_{\Phi}$. This means actually that the explicit symmetry
breaking happens at scales higher than the spontaneous breaking, which in some
sense contradicts the notion of the global symmetry. We consider this
configuration as theoretically inconsistent, and we remain only with the small
triangle to the right of the diagonal line in fig. 3 (left), which however,
still implies a mass for the pseudo-Goldstone not too far from the VEV, and
also requires “large” mixing angles $\theta_{\nu\chi}$, saturating the bound
from BBN.

Alternatively, we can consider the situation $X=\rho$, i.e., the real part of
the complex scalar plays the role of the mediator particle. In this case we
can avoid explicit symmetry breaking at all and keep $\phi$ to be a massless
Goldstone. Note that $\rho$ and $\phi$ have the same couplings to the
fermions, eqs. (3.15), and therefore, the phenomenological discussion from
sec. 4 applies equally to $\rho$ and $\phi$. For the mass of $\rho$ we have
$m_{\rho}=\sqrt{2\lambda_{\Phi}}v_{\Phi}$. Therefore, the perturbativity
requirement for the quartic coupling, $\lambda_{\Phi}<\sqrt{4\pi}$, implies
that $m_{\rho}$ cannot be much larger than the VEV: $m_{\rho}\lesssim
2.7v_{\Phi}$. Hence, we see that again we are restricted to the small region
to the right of the diagonal line in fig. 3. In addition, in this case there
is an additional constraint related to the presence of the massless pseudo-
scalar $\phi$, which now can lead to 2-to-2 processes such as
$\chi\chi\to\phi\phi$ and $\nu\nu\to\phi\phi$, which can suppress neutrino
freestreaming. For this not to happen, these processes cannot be in thermal
equilibrium at $z<10^{5}$ Taule  _et al._ (2022). Explicit calculations show
that the coupling mediating this process then should be $\lambda\lesssim
7\times 10^{-7}$ Forastieri  _et al._ (2019). In particular, this bound will
apply to $\lambda_{\phi}^{\nu\chi}$, see eq. (3.15b), and can be interpreted
as a lower bound on the scalar VEV that reads:

$\displaystyle v_{\Phi}>0.3\,{\rm
keV}\,\frac{|\theta_{\nu\chi}|}{10^{-3}}\left(\frac{m_{\nu}}{0.2\,{\rm
eV}}\right)\,.$ (5.1)

This constraint is shown as a vertical green line in fig. 3 and we see that
this bound, together with the perturbativity requirement for the quartic
coupling, puts also severe restrictions to the parameter space in the case of
global $U(1)_{X}$ with no explicit symmetry breaking ($m_{\phi}=0$).

We conclude that the global symmetry version of the model is severely
constrained by perturbativity and theoretical consistency arguments. This is
not the case for the gauged version: we indicate in fig. 4 by diagonal lines
the values of the gauge coupling $g_{X}=m_{Z^{\prime}}/v_{\Phi}$, and we can
see that the allowed region does not have any perturbativity problem since the
couplings are sufficiently small (unless $\theta_{\nu\chi}$ is not chosen much
smaller than the values adopted in the figure, see also the discussion in sec.
4).

In the following subsections we discuss further cosmological consequences of
the model.

###  5.1 Heavy right-handed neutrinos, $\Delta N_{\rm eff}$, and leptogenesis

The mass of the heavy neutrinos $M_{R}$ is essentially unconstrained, as for
given $m_{\nu},v_{\Phi}$ and $\theta_{\nu\chi}$ we still can choose the Yukawa
couplings $Y_{\Phi}$. However, we can derive an upper bound on $M_{R}$ by
requiring perturbativity for $Y_{\Phi}$, i.e. that $Y_{\Phi}$ has to be
smaller than $\sqrt{4\pi}$. In the one flavour approximation this reads

$Y_{\Phi}=\dfrac{\sqrt{m_{\nu}M_{R}}}{v_{\Phi}}\theta_{\nu\chi}\leq\sqrt{4\pi}\,,$
(5.2)

which gives an upper bound on $M_{R}$:

$M_{R}\leq 4\pi\dfrac{v_{\Phi}^{2}}{m_{\nu}\,\theta_{\nu\chi}^{2}}\,.$ (5.3)

This bound is shown in figs. 3 and 4 as vertical dotted lines. We see that in
the global case the heavy neutrinos are relatively light, with masses below
the TeV scale in a large fraction of the parameter space. However, the mixing
between active and heavy neutrinos given by eq. (3.12) is too small to be
tested in collider searches, see e.g. Abada  _et al._ (2022). In the gauge
case we have upper bounds on $M_{R}$ in the range between $10^{6}\,{\rm GeV}$
and $10^{17}\,{\rm GeV}$. In general the heavy neutrinos are not going to
interfere with the mechanism altering the neutrino number density before
recombination, but nevertheless we can wonder about other cosmological
implications they can lead to.

Right handed neutrinos that give masses to the active neutrinos are expected
to thermalise in the early Universe with the SM plasma via inverse decays of
the form $H\,L\to N$, see e.g. Buchmuller  _et al._ (2005). That means that if
the reheating temperature of the Universe was $T_{\rm RH}>M_{R}$ then we
expect these sterile neutrinos to have cosmological implications. Perhaps the
most widely studied consequence of them is baryogenesis via leptogenesis, see
e.g. Davidson  _et al._ (2008). In addition, in our scenario these states
interact with a large number of very light species, which can lead to a
primordial $\Delta N_{\rm eff}$ that is effective both at the time of BBN and
recombination and would therefore exclude the model. The key parameter
controlling the contribution to $\Delta N_{\rm eff}$ and also the possible
effect on the generation of a lepton asymmetry is the branching ratio into the
dark sector:

$\displaystyle{\rm
BR}(N\to\phi\,\chi)=\dfrac{2\,N_{\chi}\,\theta_{\nu\chi}^{2}\,v_{\rm
EW}^{2}}{2\,N_{\chi}\,\theta_{\nu\chi}^{2}\,v_{\rm EW}^{2}+v_{\Phi}^{2}}\,.$
(5.4)

The interaction rate of heavy sterile neutrinos with the SM plasma is directly
determined by the Yukawa coupling $Y_{\nu}$ which in the type-I seesaw
mechanism is intimately related to the neutrino mass, $m_{\nu}\simeq
Y_{\nu}^{2}v_{\rm EW}^{2}/2M_{R}$. We thus expect the sterile neutrinos to
reach thermal equilibrium in the early Universe for $K_{\rm SM}\gtrsim 1$,
where

$K_{\rm SM}\equiv\frac{\Gamma(N\to HL)}{H(T\simeq M_{R})}\simeq
200\,\left(\frac{m_{\nu}}{0.2\,{\rm eV}}\right)\,.$ (5.5)

Clearly, if these states have efficient interactions with $\chi_{L}$ and
$\phi$ then they will generate a large number density of massless states that
will contribute to $\Delta N_{\rm eff}$. In addition, these new interactions,
if efficient, can make the distribution function of heavy sterile neutrinos
very close to thermal and therefore suppress the generation of a primordial
lepton asymmetry.

$\Delta N_{\rm eff}$ can be written as

$\displaystyle\Delta N_{\rm eff}$
$\displaystyle\equiv\frac{8}{7}\left(\frac{11}{4}\right)^{4/3}\left[\frac{\rho_{\rm
DS}}{\rho_{\gamma}}\right]\,,$ (5.6)

where $\rho_{\rm DS}$ refers to the energy density in the dark sector, i.e.
$\chi_{L}$ and $\phi$. This contribution is bounded from above because the
maximum temperature the $\phi-\chi$ population could have obtained is the SM
one, in which case we can use entropy conservation to find:

$\displaystyle\Delta N_{\rm
eff}\leq\frac{8}{7}\left(\frac{11}{4}\right)^{4/3}\left[\frac{g_{\chi}\frac{7}{8}N_{\chi}+g_{\phi}}{2}\right]\left[\frac{g_{\star
S}^{\rm SM}(T_{\rm today})}{g_{\star S}^{\rm SM}(T\simeq
M_{R})}\right]^{4/3}\,.$

Using $g_{\star S}^{\rm SM}(T_{\rm today})\simeq 3.9$ and $g_{\star S}^{\rm
SM}(T\simeq M_{R})\simeq 105$, this yields a maximum contribution to $\Delta
N_{\rm eff}$ of

$\displaystyle\Delta N_{\rm eff}|_{\rm max}\simeq
0.027\left[1+\frac{7}{8}g_{\chi}N_{\chi}\right]\,.$ (5.7)

For instance, for typical values we have $\Delta N_{\rm eff}|_{\rm
max}(N_{\chi}=10,g_{\chi}=2)\simeq 0.5$, which is already excluded.

By solving explicitly the Boltzmann equations for the number density of heavy
sterile neutrinos in the early Universe allowing them to decay into these new
very light species we find that in order to satisfy the CMB and BBN bound
$\Delta N_{\rm eff}\lesssim 0.3$ we require:

$\displaystyle K_{\rm SM}\times{\rm BR}(N\to\phi\,\chi)\lesssim 1\,,$ (5.8)

where this expression is rather insensitive to $N_{\chi}$ provided that
$N_{\chi}\gtrsim 5$.444This simply means that the rate of interactions of
sterile neutrinos with the dark sector should be slower than the expansion
rate at $T\simeq M_{R}$: $\left<{\Gamma(N\to\phi\,\chi)}\right>\lesssim
H(T\simeq M_{R})$. From eq. (5.5) we obtain $10<K_{\rm SM}<800$, where the
limits come from taking $\sqrt{|\Delta m_{\rm sol}^{2}|}<m_{\nu}<0.8\,{\rm
eV}$. Putting this together, we see that the branching ratio of sterile
neutrinos is bounded to be

$\displaystyle{\rm BR}(N\to\phi\,\chi)\lesssim(10^{-3}-0.1)\,.$ (5.9)

Using eq. (5.4), we find that for the global version of the model the full
viable parameter space is incompatible with this condition. Therefore in the
global case we need to require $T_{\rm RH}<M_{R}$, such that the heavy
neutrinos cannot thermalise. Since in the global case the right handed
neutrinos could have masses $M_{R}<100\,{\rm GeV}$ this means that one would
need to invoke a mechanism that explains the baryon asymmetry of the Universe
without the help of sphalerons processes. This argument again makes the global
version of the model less attractive compared to the gauged case, in addition
to the discussion above.

For the gauge case there is some parameter space fulfilling eq. (5.9) which is
displayed as a vertical band in purple labelled “Leptogenesis” in fig. 4. For
these values of $v_{\Phi}$ the contribution to $\Delta N_{\rm eff}$ remains
small and we can be confident that the usual thermal leptogenesis mechanism
will not be distorted because the interactions with the new dark sector will
always be slower than the expansion rate of the Universe and therefore will
not impact the potential generation of a primordial lepton asymmetry. Of
course, for the standard thermal leptogenesis mechanism to be successful the
sterile neutrino mass cannot be arbitrarily low and should be $M_{R}\gtrsim
10^{8}\,{\rm GeV}$ Davidson and Ibarra (2002). For the rest of the parameter
space (to the left of the purple band), in order to be in agreement with BBN
and CMB bounds on $\Delta N_{\rm eff}$, we need to require $T_{\rm RH}<M_{R}$
and the usual thermal leptogenesis will not be operative. However, this
parameter space is compatible with having $T_{\rm RH}$ large enough to allow
for sphaleron processes.

###  5.2 Active neutrino decay

The existence of interacting fermions lighter than active neutrinos allow
neutrinos to decay. Furthermore, if active neutrinos were to decay on a
cosmological timescale this will diminishing their density today, see e.g.
Escudero  _et al._ (2020). Two body decays
$\nu\to\chi\,\,\phi/\rho/Z^{\prime}$ are forbidden in our scenario because the
$\phi/\rho/Z^{\prime}$ bosons are heavier than active neutrinos. As such, the
only possible decay channels at tree-level are $\nu\to 3\chi$ and
$\nu_{i}\to\nu_{j}\chi\chi$. The first cannot happen in the global case
because $\lambda^{\chi\chi}_{\phi\rho}=0$, see eq. (3.15). The rate for the
second is very small specially for the region of degenerate active neutrino
masses we are interested in.

However, the $\nu\to 3\chi$ decay is possible for the gauge mediated case
through the gauge interactions
$\lambda_{Z^{\prime}}^{\nu\chi}\overline{\nu}\gamma^{\mu}\chi$ and
$g_{X}\,Z^{\prime}_{\mu}\,\overline{\chi}\,\gamma^{\mu}\,\,\chi$. Assuming the
mixing between active and sterile neutrinos is similar for all the $N_{\chi}$
sterile species, we get $2\,N_{\chi}^{2}$ possible decay channels and
therefore the total active neutrino decay rate will be $\Gamma_{\rm
tot}^{\nu}=2\,N_{\chi}^{2}\,\Gamma_{\nu\to 3\chi}$. Taking the decay rate
calculated in Escudero  _et al._ (2020) we can compute the lifetime of the
active neutrinos with the parameters of our model in terms of the age of the
universe, $t_{U}\simeq 13.8$ Gyr, and the number of new sterile species, which
reads:

$\displaystyle\tau_{\nu}=4\times
10^{3}\,t_{U}\,\frac{\left(\dfrac{v_{\Phi}}{1\,\rm
MeV}\right)^{4}}{\left(\dfrac{N_{\chi}}{6}\right)^{2}\,\left(\dfrac{m_{\nu}}{0.2\,\rm
eV}\right)^{5}\left(\dfrac{|\theta_{\nu\chi}|}{10^{-3}}\right)^{2}}\,.$ (5.10)

Given this equation it can be seen from fig. 4 that neutrinos do not decay on
cosmological timescales in the regions of parameter space in which the
mechanism is viable.

##  6 Conclusions

We have discussed a UV complete model to realise a mechanism to relax
cosmological neutrino mass bounds to the level that the neutrino mass becomes
observable in terrestrial experiments such as KATRIN and/or neutrinoless
double-beta decay. The mechanism is based on an idea put forward by Farzan and
Hannestad in Ref. Farzan and Hannestad (2016). The main ingredient is the
depletion of the cosmological abundance of active neutrinos by populating
instead a sufficient number of massless species, see eq. (2.5) and figs. 1 and
2.

Our model realisation is based on a version of the seesaw mechanism sometimes
called “minimal extended seesaw” with the following main features: we have 3
heavy right-handed neutrinos, responsible for generating active neutrino
masses via the type-I seesaw. In addition we introduce a large number of
additional sterile neutrino species, which remain massless and provide the
dark radiation. In order to achieve the required relaxation of the neutrino
mass bound to $\simeq 1$ eV for the sum of neutrino masses, we need
$\mathcal{O}(10)$ species of massless sterile neutrinos. To realise the
conversion of active to sterile neutrinos we introduce a $U(1)_{X}$ symmetry,
which can be either global or gauged and is spontaneously broken by a complex
scalar field. The (pseudo-)scalar or vector bosons act as mediators between
active neutrinos and the dark sector.

The depletion of active neutrinos should happen between BBN and recombination.
This restricts the mass of the mediator particle to the range between 1 keV
and 1 MeV. We provide a detailed discussion of the relevant phenomenology and
cosmological constraints in sec. 4 and find that the mechanism can work only
in a closed region of parameter space. The version with the global symmetry is
severely constrained by theoretical consistency arguments and perturbativity
requirements. Furthermore, it is incompatible with thermal leptogenesis and
requires a reheating temperature below the electro-weak scale. In contrast,
the gauged version of the model is safely in the perturbative regime, is
partially compatible with thermal leptogenesis and allows for reheating
temperatures above the electro-weak scale. In the gauged version of the model,
the VEV breaking the $U(1)_{X}$ symmetry can be between 10 MeV and a few GeV
when the mixing $|\theta_{\nu\chi}|=10^{-3}$ while for
$|\theta_{\nu\chi}|=10^{-4}$ only a narrow window of VEVs is allowed $20\,{\rm
MeV}\lesssim v_{\Phi}\lesssim 100\,{\rm MeV}$.

In conclusion, we have presented a relatively simple extension of the Standard
Model, explaining active neutrino masses close to the eV scale, featuring a
non-standard neutrino cosmology. The main signature of the model is the
possibility of observing a non-zero neutrino mass in KATRIN and/or a signal in
neutrinoless double-beta decay, which would be excluded within the standard
$\Lambda$CDM cosmology. The model requires the presence of $\mathcal{O}(10)$
species of massless sterile neutrinos, whose theoretical motivation remains to
be identified, see however Arkani-Hamed  _et al._ (2016); Buchmuller  _et al._
(2007); Ellis and Lebedev (2007) for extensions of the Standard Model with
large number of sterile neutrinos or BSM species.

Looking forward, in this scenario the number density of active neutrinos
forming the cosmic neutrino background (CNB) is smaller than in $\Lambda$CDM.
However, these neutrinos have a smaller temperature and therefore can cluster
much more efficiently in the solar system. This in turn enhances the direct
detection prospects of the CNB within this cosmology Alvey  _et al._ (2022b).
In addition, while in this work we have focused on relaxing the current
cosmological neutrino mass bound, the essence of the model remains relatively
unchanged and will be relevant for the interpretation of upcoming cosmological
neutrino mass measurements from DESI/Euclid, see Fig 1.

Acknowledgements. This project has received support from the European Union’s
Horizon 2020 research and innovation programme under the Marie Sklodowska-
Curie grant agreement No 860881-HIDDeN. The work of JTC is supported by the
Ministerio de Ciencia e Innovación under FPI contract PRE2019-089992 of the
SEV-2015-0548 grant. JTC is very grateful for the hospitality of the
Theoretical Astroparticle Physics group of KIT during the months spent there
when most of this project was carried out and to CERN-TH where this project
was finalised. JTC would also like to thank Andrés Castillo for useful
discussion. ME and TSM thank the Instituto de Fisica Teorica (IFT UAM-CSIC) in
Madrid for support via the Centro de Excelencia Severo Ochoa Program under
Grant CEX2020-001007-S, during the Extended Workshop “Neutrino Theories”,
where part of this work was developed.

##  Appendix A The Mechanism after recombination

In the main text we have considered the possibility of decreasing the number
density of neutrinos before recombination in order to relax the cosmological
neutrino mass bound. However, a priori it is not really needed that the
depletion of the neutrino number density happens before recombination, it can
also happen after. In this appendix we comment on this possibility.

First, let us note that the mean momentum of a neutrino at the time of
recombination is $\left<p_{\nu}\right>\simeq 3T_{\nu}\simeq 0.6\,{\rm eV}$. As
such, the primary CMB spectra can only be affected by neutrinos that are
$m_{\nu}\gtrsim 0.6\,{\rm eV}$. Indeed, the WMAP bound on the neutrino mass
within $\Lambda$CDM was $\sum m_{\nu}<1.3\,{\rm eV}$ at 95% CL Hinshaw  _et
al._ (2013), very close to the estimate given above. Nevertheless, the CMB is
indirectly sensitive to less massive neutrinos via their effect on the growth
of the large scale structures in the Universe. This sensitivity steams from
the fact that such large scale structures can lens the CMB photons emitted at
the surface of last scattering and detected today, see e.g. Lattanzi and
Gerbino (2018). The energy density in non-relativistic neutrinos suppresses
this growth and CMB observations are therefore sensitive to significantly
lighter neutrinos, e.g. $\sum m_{\nu}<0.24\,{\rm eV}$ at 95% CL is the bound
from Planck CMB observations alone Aghanim  _et al._ (2020). Importantly, the
main contribution to the lensing of CMB photons happens only at moderately low
redshifts, $z\lesssim 6$, and this means that provided that $m_{\nu}<0.5\,{\rm
eV}$ Planck CMB observations are really sensitive to $\Omega_{\nu}h^{2}$ at
$z\lesssim 6$. Similar considerations apply for data from Baryon Acoustic
Oscillations. Put together, this means that in principle the mechanism
described in the main text could work equally well provided that the neutrino
number density is reduced at $z\gtrsim 10$.

To concretely see how this could work we first note that neutrinos become non-
relativistic at $z_{\rm nr}\simeq 200\,m_{\nu}/0.1\,{\rm eV}$. In the scenario
we considered in the main text we need a particle $X$ that can decay to
neutrinos and therefore should have a mass $m_{X}>2\,m_{\nu}$. That means that
the mass of this particle would need to be in a rather narrow range $0.1\,{\rm
eV}\lesssim m_{X}\lesssim 1\,{\rm eV}$. In addition, this particle should
interact efficiently with neutrinos at $z\lesssim 10^{3}$ but not at $z\gtrsim
10^{3}$. This is simply because otherwise this new state will suppress
neutrino freestreaming and distort the CMB Taule  _et al._ (2022). Thus, we
can appreciate the three requirements are that the interactions are not
effective until $z\lesssim 10^{3}$, that the mediator is in the right mass
window $0.1\,{\rm eV}\lesssim m_{X}\lesssim 1\,{\rm eV}$, and that the
interactions occur rapidly as compared with the timescale for expansion. We
believe that these three conditions can be met although we think that this
will only occur in a very limited region of parameter space for $m_{\nu}$ and
$m_{X}$ and for the interactions between $X$ and neutrinos. To confidently
study the viable parameter space in this case a solution of the relevant
Boltzmann equations would be needed which is beyond the scope of this study.

Finally, one can consider the possibility of depleting the neutrino number
density by neutrino annihilations into massless states,
$\nu\bar{\nu}\to\chi\bar{\chi}$. This in fact was considered a long time ago
Beacom  _et al._ (2004), but was shown indeed not to be possible because of
the suppression of neutrino freestreaming from these interactions Hannestad
(2005).

##  Appendix B Entropy production and the value of $N_{\rm eff}$

The thermodynamic processes discussed in the main text that depletes the
number density of neutrinos in the early Universe leads to entropy production.
This in particular means that we do not expect $N_{\rm eff}$ to be exactly
3.044. However, for large $N_{\chi}$ we expect the difference not to be
significant. The reason that $N_{\rm eff}\neq 3.044$ is because the production
and decay of a massive particle, $X$, always leads to some extra energy
density production stored in its decay products. Of course, the more $\chi$
states there are the less significant will be this release of energy because
the fractional number density of $X$ states will be small. Since we are
working in a regime where the relevant processes are very efficient with
respect to the characteristic expansion timescale $t_{U}\sim 1/H$, we can use
equilibrium thermodynamics to elucidate this matter.

We should split the discussion in two events that we consider happen
instantaneously. First, we consider that at $T\gg m_{X}$ the $X$ bosons and
the $\chi$ fermions thermalise instantaneously generating a coupled system
with temperature $T_{\rm eq}$ and chemical potential $\mu_{i}^{\rm eq}$. Since
all the particles are massless energy conservation implies:

$\displaystyle\\!\\!\\!\\!\\!\rho_{\nu}(T_{\nu},0)=\rho_{\nu}(T_{\rm
eq},\mu_{\nu}^{\rm eq})+\rho_{\chi}(T_{\rm eq},\mu_{\chi}^{\rm
eq})+\rho_{X}(T_{\rm eq},\mu_{X}^{\rm eq})\,.$ (B.1)

In addition, since we are considering processes $X\leftrightarrow\nu{\nu}$ and
$X\leftrightarrow\chi{\chi}$ we know the chemical potentials of these species
are related to be $2\mu_{\nu}^{\rm eq}=2\mu_{\chi}^{\rm eq}=\mu_{X}^{\rm
eq}$555Note that in the global case considered in the main text the relevant
processes are $X\leftrightarrow\nu{\nu}$ and $X\leftrightarrow\nu{\chi}$.
Although they are different they also lead to the condition $2\mu_{\nu}^{\rm
eq}=2\mu_{\chi}^{\rm eq}=\mu_{X}^{\rm eq}$. Similarly, while the relevant
processes in the gauge case are $X\leftrightarrow\chi{\chi}$,
$X\leftrightarrow\nu{\chi}$, and $X\leftrightarrow\chi{\chi}$, we end up
having the same result..

In addition to energy conservation, we know that the processes we consider
respect a conservation law for the number density of the various species:

$\displaystyle n_{\nu}(T_{\nu},0)=n_{\nu}(T_{\rm eq},\mu_{\nu}^{\rm
eq})+N_{\chi}(T_{\rm eq},\mu_{\chi}^{\rm eq})+2n_{X}(T_{\rm eq},\mu_{X}^{\rm
eq})\,,$ (B.2)

where here the 2 comes from the fact that a pair of neutrinos or $\chi$
particles can be interchanged by an $X$ state. Solving these two equations
simultaneously provides $T_{\rm eq}$ and $\mu_{\nu}^{\rm eq}$ from a given
$T_{\nu}$.

After thermalisation happens the Universe expands and eventually the $X$ boson
decays away at $T\ll m_{X}$. Since we are in thermal equilibrium, in this
system entropy density is conserved and number density is also conserved as
before. This explicitly means:

$\displaystyle\\!\\!\\!\\!a^{3}(s_{\nu}(T,\mu)+s_{\chi}(T,\mu)+s_{X}(T,2\mu))=\text{constant}\,,$
(B.3)
$\displaystyle\\!\\!\\!\\!a^{3}(n_{\nu}(T,\mu)+N_{\chi}(T,\mu)+2n_{X}(T,2\mu))=\text{constant}\,.$
(B.4)

Thus, solving eqs. (B.3) and (B.4) simultaneously using as initial condition
$T_{\rm eq}$ and $\mu_{\rm eq}$ allows one to find $T_{f}$ and $\mu_{f}$.
Namely, the final temperature and chemical potential of both neutrinos and
massless $\chi$ species as relevant for CMB observations once the $X$ particle
has disappeared from the plasma.

It is important to highlight that in order to obtain the decrease in the
neutrino number density after this process one does not actually need to solve
these equations explicitly. The number density conservation in eq. (B.2) and
in eq. (B.4) imply

$\displaystyle
n_{\nu}(T_{\nu},0)=n_{\nu}(T_{f},\mu_{f})+N_{\chi}(T_{f},\mu_{f})\,.$ (B.5)

Since the $\chi$ particles are fermions we can easily see that666While in the
main text we consider ultralight fermions the same consideration can be done
for bosons. In such a case a compact formula cannot be obtained. A numerical
fit to the result is: $n_{\nu}(T_{f},\mu_{f})=n_{\nu}^{\rm
SM}/\left({1+g_{\chi}N_{\chi}/6+0.07}\right)$ which is accurate to better than
0.5% for any $N_{\chi}$.:

$\displaystyle n_{\nu}(T_{f},\mu_{f})=n_{\nu}^{\rm
SM}\frac{1}{1+g_{\chi}N_{\chi}/6}\,,$ (B.6)

as highlighted by eq. (2.5).

Obtaining $N_{\rm eff}$, however, does require solving explicitly for these
four equations simultaneously. In fig. 5 we show its value as a function of
$N_{\chi}$. We can appreciate that unless $N_{\chi}$ is small,
$N_{\chi}\lesssim 3$ the contribution to $N_{\rm eff}$ at the time of
recombination should be undetectable even with an ultrasensitive Stage-IV CMB
experiment.

Finally, we note that actually the processes we consider lead to a significant
production of entropy because a large number of new relativistic species are
produced. Importantly this does not change the number density of photons and
therefore have not relevant impact on what the number density of baryons or
dark matter particles should be. For completeness, however we report an
approximate fitting function of $g_{\star S}$ which reads:

$\displaystyle g_{\star S}\simeq 3.5+\frac{6}{10}\log(N_{\chi}^{7/10})\,,$
(B.7)

and that we see depends only logarithmically upon $N_{\chi}$.

Figure 5: $N_{\rm eff}$ as relevant for CMB observations as a function of the
number of massless species that eventually thermalise with neutrinos
$N_{\chi}$. In blue when the mediator is a vector boson with $g_{X}=3$ and in
red when we have a scalar mediator with $g_{X}=1$. Dashed lines correspond to
the case of $\chi$ being bosons and solid lines to the case of $\chi$ being
bosons. The black dotted line shows the SM prediction of $N_{\rm eff}=3.044$.

##  Appendix C Decay Rates

Here we outline relevant formulas used to compute the various of the
constraints in the main text. The thermally averaged rate of the process
$a+b\leftrightarrow X$ with $X$ being either a scalar, pseudoscalar or a
vector boson and $a,\,b$ being two light fermions is given by Escudero Abenza
(2020) :

$\langle\Gamma\rangle=\Gamma_{X}\,\dfrac{g_{X}}{4}\left(\dfrac{m_{X}}{T_{a}}\right)^{2}\,K_{1}\left(\dfrac{m_{X}}{T_{a}}\right),$
(C.1)

where $\Gamma_{X}$ is the decay rate of $X$ in vacuum, $g_{X}$ is the number
of degrees of freedom of $X$ (1 for the scalar/pseudoscalar and 3 for the
$Z^{\prime}$), $m_{X}$ is the mass of $X$, $T_{a}$ is the temperature of the
fermions (in case of neutrinos, after decoupling we will approximate
$T_{\gamma}/T_{\nu}\simeq 1.4$) and $K_{1}$ is a modified Bessel function of
the first kind. The decay rate in vacuum for the scalar and pseudoscalar is

$\Gamma_{X}=\dfrac{\lambda^{2}}{16\pi}m_{X},$ (C.2)

and for the vector boson

$\Gamma_{X}=\dfrac{\lambda^{2}}{8\pi}m_{X},$ (C.3)

with $\lambda$ the corresponding coupling, see eqs. (3.15) and (3.18). For the
cross sections of 2 to 2 processes with interaction coupling $g$ and mediator
$X$ we have used

$\sigma\sim\dfrac{g^{4}}{16\pi^{2}}\dfrac{1}{m_{X}^{4}}T^{2},$ (C.4)

when they are mediated by a massive ($m_{X}\gtrsim T$) mediator and

$\sigma\sim\dfrac{g^{4}}{16\pi^{2}}\dfrac{1}{T^{2}},$ (C.5)

when the mediator is light ($m_{X}\ll T$).

## References

  * Dolgov (2002) A. D. Dolgov, Phys. Rept. 370, 333 (2002), arXiv:hep-ph/0202122 .
  * Lesgourgues and Pastor (2006) J. Lesgourgues and S. Pastor, Phys. Rept. 429, 307 (2006), arXiv:astro-ph/0603494 .
  * Aghanim  _et al._ (2020) N. Aghanim _et al._ (Planck), Astron. Astrophys. 641, A6 (2020), arXiv:1807.06209 [astro-ph.CO] .
  * Di Valentino  _et al._ (2021) E. Di Valentino, S. Gariazzo, and O. Mena, Phys. Rev. D 104, 083504 (2021), arXiv:2106.15267 [astro-ph.CO] .
  * Aghamousa  _et al._ (2016) A. Aghamousa _et al._ (DESI), (2016), arXiv:1611.00036 [astro-ph.IM] .
  * Amendola  _et al._ (2018) L. Amendola _et al._ , Living Rev. Rel. 21, 2 (2018), arXiv:1606.00180 [astro-ph.CO] .
  * Brinckmann  _et al._ (2019) T. Brinckmann, D. C. Hooper, M. Archidiacono, J. Lesgourgues, and T. Sprenger, JCAP 01, 059 (2019), arXiv:1808.05955 [astro-ph.CO] .
  * Chudaykin and Ivanov (2019) A. Chudaykin and M. M. Ivanov, JCAP 11, 034 (2019), arXiv:1907.06666 [astro-ph.CO] .
  * Aker  _et al._ (2019) M. Aker _et al._ (KATRIN), Phys. Rev. Lett. 123, 221802 (2019), arXiv:1909.06048 [hep-ex] .
  * Aker  _et al._ (2022) M. Aker _et al._ (KATRIN), Nature Phys. 18, 160 (2022), arXiv:2105.08533 [hep-ex] .
  * Agostini  _et al._ (2022) M. Agostini, G. Benato, J. A. Detwiler, J. Menéndez, and F. Vissani, (2022), arXiv:2202.01787 [hep-ex] .
  * Gando  _et al._ (2016) A. Gando _et al._ (KamLAND-Zen), Phys. Rev. Lett. 117, 082503 (2016), [Addendum: Phys. Rev. Lett.117,no.10,109903(2016)], arXiv:1605.02889 [hep-ex] .
  * Agostini  _et al._ (2020) M. Agostini _et al._ (GERDA), Phys. Rev. Lett. 125, 252502 (2020), arXiv:2009.06079 [nucl-ex] .
  * Adams  _et al._ (2021) D. Q. Adams _et al._ (CUORE), (2021), arXiv:2104.06906 [nucl-ex] .
  * Anton  _et al._ (2019) G. Anton _et al._ (EXO-200), Phys. Rev. Lett. 123, 161802 (2019), arXiv:1906.02723 [hep-ex] .
  * Giuliani  _et al._ (2019) A. Giuliani, J. J. Gomez Cadenas, S. Pascoli, E. Previtali, R. Saakyan, K. Schäffner, and S. Schönert (APPEC Committee), (2019), arXiv:1910.04688 [hep-ex] .
  * Ettengruber  _et al._ (2022) M. Ettengruber, M. Agostini, A. Caldwell, P. Eller, and O. Schulz, Phys. Rev. D 106, 073004 (2022), arXiv:2208.09954 [hep-ph] .
  * Chacko  _et al._ (2020) Z. Chacko, A. Dev, P. Du, V. Poulin, and Y. Tsai, JHEP 04, 020 (2020), arXiv:1909.05275 [hep-ph] .
  * Escudero  _et al._ (2020) M. Escudero, J. Lopez-Pavon, N. Rius, and S. Sandner, JHEP 12, 119 (2020), arXiv:2007.04994 [hep-ph] .
  * Escudero and Fairbairn (2019) M. Escudero and M. Fairbairn, Phys. Rev. D 100, 103531 (2019), arXiv:1907.05425 [hep-ph] .
  * Barenboim  _et al._ (2021) G. Barenboim, J. Z. Chen, S. Hannestad, I. M. Oldengott, T. Tram, and Y. Y. Y. Wong, JCAP 03, 087 (2021), arXiv:2011.01502 [astro-ph.CO] .
  * Chacko  _et al._ (2021) Z. Chacko, A. Dev, P. Du, V. Poulin, and Y. Tsai, Phys. Rev. D 103, 043519 (2021), arXiv:2002.08401 [astro-ph.CO] .
  * Abellán  _et al._ (2021) G. F. Abellán, Z. Chacko, A. Dev, P. Du, V. Poulin, and Y. Tsai, (2021), arXiv:2112.13862 [hep-ph] .
  * Chen  _et al._ (2022) J. Z. Chen, I. M. Oldengott, G. Pierobon, and Y. Y. Y. Wong, Eur. Phys. J. C 82, 640 (2022), arXiv:2203.09075 [hep-ph] .
  * Dvali and Funcke (2016) G. Dvali and L. Funcke, Phys. Rev. D 93, 113002 (2016), arXiv:1602.03191 [hep-ph] .
  * Dvali  _et al._ (2021) G. Dvali, L. Funcke, and T. Vachaspati, (2021), arXiv:2112.02107 [hep-ph] .
  * Lorenz  _et al._ (2019) C. S. Lorenz, L. Funcke, E. Calabrese, and S. Hannestad, Phys. Rev. D 99, 023501 (2019), arXiv:1811.01991 [astro-ph.CO] .
  * Lorenz  _et al._ (2021) C. S. Lorenz, L. Funcke, M. Löffler, and E. Calabrese, Phys. Rev. D 104, 123518 (2021), arXiv:2102.13618 [astro-ph.CO] .
  * Esteban and Salvado (2021) I. Esteban and J. Salvado, JCAP 05, 036 (2021), arXiv:2101.05804 [hep-ph] .
  * Farzan and Hannestad (2016) Y. Farzan and S. Hannestad, JCAP 02, 058 (2016), arXiv:1510.02201 [hep-ph] .
  * Renk  _et al._ (2021) J. J. Renk _et al._ (GAMBIT Cosmology Workgroup), JCAP 02, 022 (2021), arXiv:2009.03286 [astro-ph.CO] .
  * Cuoco  _et al._ (2005) A. Cuoco, J. Lesgourgues, G. Mangano, and S. Pastor, Phys. Rev. D 71, 123501 (2005), arXiv:astro-ph/0502465 .
  * Oldengott  _et al._ (2019) I. M. Oldengott, G. Barenboim, S. Kahlen, J. Salvado, and D. J. Schwarz, JCAP 04, 049 (2019), arXiv:1901.04352 [astro-ph.CO] .
  * Alvey  _et al._ (2022a) J. Alvey, M. Escudero, and N. Sabti, JCAP 02, 037 (2022a), arXiv:2111.12726 [astro-ph.CO] .
  * Alvey  _et al._ (2022b) J. Alvey, M. Escudero, N. Sabti, and T. Schwetz, Phys. Rev. D 105, 063501 (2022b), arXiv:2111.14870 [hep-ph] .
  * Chun  _et al._ (1995) E. J. Chun, A. S. Joshipura, and A. Y. Smirnov, Phys. Lett. B 357, 608 (1995), arXiv:hep-ph/9505275 .
  * Barry  _et al._ (2011) J. Barry, W. Rodejohann, and H. Zhang, JHEP 07, 091 (2011), arXiv:1105.3911 [hep-ph] .
  * Zhang (2012) H. Zhang, Phys. Lett. B 714, 262 (2012), arXiv:1110.6838 [hep-ph] .
  * Heeck and Zhang (2013) J. Heeck and H. Zhang, JHEP 05, 164 (2013), arXiv:1211.0538 [hep-ph] .
  * Ballett  _et al._ (2019) P. Ballett, M. Hostert, and S. Pascoli, Phys. Rev. D 99, 091701 (2019), arXiv:1903.07590 [hep-ph] .
  * Escudero Abenza (2020) M. Escudero Abenza, JCAP 05, 048 (2020), arXiv:2001.04466 [hep-ph] .
  * Akita and Yamaguchi (2020) K. Akita and M. Yamaguchi, JCAP 08, 012 (2020), arXiv:2005.07047 [hep-ph] .
  * Froustey  _et al._ (2020) J. Froustey, C. Pitrou, and M. C. Volpe, JCAP 12, 015 (2020), arXiv:2008.01074 [hep-ph] .
  * Bennett  _et al._ (2021) J. J. Bennett, G. Buldgen, P. F. De Salas, M. Drewes, S. Gariazzo, S. Pastor, and Y. Y. Y. Wong, JCAP 04, 073 (2021), arXiv:2012.02726 [hep-ph] .
  * Allahverdi  _et al._ (2020) R. Allahverdi _et al._ ,  (2020), 10.21105/astro.2006.16182, arXiv:2006.16182 [astro-ph.CO] .
  * Pospelov and Pradler (2010) M. Pospelov and J. Pradler, Ann. Rev. Nucl. Part. Sci. 60, 539 (2010), arXiv:1011.1054 [hep-ph] .
  * Taule  _et al._ (2022) P. Taule, M. Escudero, and M. Garny, Phys. Rev. D 106, 063539 (2022), arXiv:2207.04062 [astro-ph.CO] .
  * Georgi and Glashow (2000) H. Georgi and S. L. Glashow, Phys. Rev. D 61, 097301 (2000), arXiv:hep-ph/9808293 .
  * Ellis and Lola (1999) J. R. Ellis and S. Lola, Phys. Lett. B 458, 310 (1999), arXiv:hep-ph/9904279 .
  * Casas  _et al._ (1999) J. A. Casas, J. R. Espinosa, A. Ibarra, and I. Navarro, Nucl. Phys. B 556, 3 (1999), arXiv:hep-ph/9904395 .
  * de Medeiros Varzielas  _et al._ (2009) I. de Medeiros Varzielas, G. G. Ross, and M. Serna, Phys. Rev. D 80, 073002 (2009), arXiv:0811.2226 [hep-ph] .
  * Boudjemaa and King (2009) S. Boudjemaa and S. F. King, Phys. Rev. D 79, 033001 (2009), arXiv:0808.2782 [hep-ph] .
  * Grimus and Lavoura (2000) W. Grimus and L. Lavoura, JHEP 11, 042 (2000), arXiv:hep-ph/0008179 .
  * Gherghetta  _et al._ (2019) T. Gherghetta, J. Kersten, K. Olive, and M. Pospelov, Phys. Rev. D 100, 095001 (2019), arXiv:1909.00696 [hep-ph] .
  * Bauer and Foldenauer (2022) M. Bauer and P. Foldenauer, Phys. Rev. Lett. 129, 171801 (2022), arXiv:2207.00023 [hep-ph] .
  * Janka (2012) H.-T. Janka, Ann. Rev. Nucl. Part. Sci. 62, 407 (2012), arXiv:1206.2503 [astro-ph.SR] .
  * Raffelt (1996) G. G. Raffelt,  _Stars as Laboratories for Fundamental Physics: The Astrophysics of Neutrinos, Axions, and Other Weakly Interacting Particles_ (1996).
  * Fiorillo  _et al._ (2022) D. F. G. Fiorillo, G. G. Raffelt, and E. Vitagliano, (2022), arXiv:2209.11773 [hep-ph] .
  * Barbieri and Dolgov (1991) R. Barbieri and A. Dolgov, Nucl. Phys. B 349, 743 (1991).
  * Abazajian (2006) K. Abazajian, Phys. Rev. D 73, 063506 (2006), arXiv:astro-ph/0511630 .
  * Pisanti  _et al._ (2021) O. Pisanti, G. Mangano, G. Miele, and P. Mazzella, JCAP 04, 020 (2021), arXiv:2011.11537 [astro-ph.CO] .
  * Böser  _et al._ (2020) S. Böser, C. Buck, C. Giunti, J. Lesgourgues, L. Ludhova, S. Mertens, A. Schukraft, and M. Wurm, Prog. Part. Nucl. Phys. 111, 103736 (2020), arXiv:1906.01739 [hep-ex] .
  * Dentler  _et al._ (2018) M. Dentler, A. Hernández-Cabezudo, J. Kopp, P. A. N. Machado, M. Maltoni, I. Martinez-Soler, and T. Schwetz, JHEP 08, 010 (2018), arXiv:1803.10661 [hep-ph] .
  * Forastieri  _et al._ (2019) F. Forastieri, M. Lattanzi, and P. Natoli, Phys. Rev. D 100, 103526 (2019), arXiv:1904.07810 [astro-ph.CO] .
  * Abada  _et al._ (2022) A. Abada, P. Escribano, X. Marcano, and G. Piazza, (2022), arXiv:2208.13882 [hep-ph] .
  * Buchmuller  _et al._ (2005) W. Buchmuller, P. Di Bari, and M. Plumacher, Annals Phys. 315, 305 (2005), arXiv:hep-ph/0401240 .
  * Davidson  _et al._ (2008) S. Davidson, E. Nardi, and Y. Nir, Phys. Rept. 466, 105 (2008), arXiv:0802.2962 [hep-ph] .
  * Davidson and Ibarra (2002) S. Davidson and A. Ibarra, Phys. Lett. B 535, 25 (2002), arXiv:hep-ph/0202239 .
  * Arkani-Hamed  _et al._ (2016) N. Arkani-Hamed, T. Cohen, R. T. D’Agnolo, A. Hook, H. D. Kim, and D. Pinner, Phys. Rev. Lett. 117, 251801 (2016), arXiv:1607.06821 [hep-ph] .
  * Buchmuller  _et al._ (2007) W. Buchmuller, K. Hamaguchi, O. Lebedev, S. Ramos-Sanchez, and M. Ratz, Phys. Rev. Lett. 99, 021601 (2007), arXiv:hep-ph/0703078 .
  * Ellis and Lebedev (2007) J. R. Ellis and O. Lebedev, Phys. Lett. B 653, 411 (2007), arXiv:0707.3419 [hep-ph] .
  * Hinshaw  _et al._ (2013) G. Hinshaw _et al._ (WMAP), Astrophys. J. Suppl. 208, 19 (2013), arXiv:1212.5226 [astro-ph.CO] .
  * Lattanzi and Gerbino (2018) M. Lattanzi and M. Gerbino, Front. in Phys. 5, 70 (2018), arXiv:1712.07109 [astro-ph.CO] .
  * Beacom  _et al._ (2004) J. F. Beacom, N. F. Bell, and S. Dodelson, Phys. Rev. Lett. 93, 121302 (2004), arXiv:astro-ph/0404585 .
  * Hannestad (2005) S. Hannestad, JCAP 02, 011 (2005), arXiv:astro-ph/0411475 .