††institutetext: Department of Physics, Royal Holloway, University of London,
Egham, Surrey, TW20 0EX, United Kingdom

# Dark Matter and Neutrino Masses in a Portalino-like Model

Dugald Hepburn  and Stephen M. West<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We explore a Portalino-like model of dark matter and neutrino masses in which
right-handed neutrino fields connect gauge neutral operators from the Standard
Model and Hidden Sector. Neutrino masses are generated via a seesaw-like
mechanism that can explain the light active neutrino masses. The model
includes a “Portalino” state that connects the two sectors via the neutrino
portal. Dark Matter in this model consists of a hidden sector Dirac fermion
that dominantly freezes-out via resonant annihilations into other hidden
sector states, which ultimately results in a population of Portalinos. Due to
small mixing in the extended neutrino sector these Portalinos tend to be
cosmologically long lived, decaying into Standard Model particles leading to
constraints on the model from Big Bang Nucleosynthesis and measurements of the
Cosmic Microwave Background radiation. Combining these limits with direct
constraints on the size of the Portalino-neutrino mixing and the assumptions
of the model the viable mass ranges for the Portalino states are found to be
$$0.02\text{\,}\mathrm{eV}$\lesssim m_{n}\lesssim$6.4\text{\,}\mathrm{eV}$$ or
$$489\text{\,}\mathrm{MeV}$\lesssim m_{n}\lesssim$ TeV. Indirect dark matter
signals in the form of highly boosted, mono-energetic Portalinos produced in
Dark Matter annihilations provide a target for neutrino telescopes.

##  1 Introduction

The question of how the dark sector interacts with the visible sector, if it
does at all, underpins the uncertainty surrounding the nature and origin of
Dark Matter (DM). Many proposals have been made for how a connection can be
established through so called “portals”, including the Higgs portal see e.g.
Schabinger:2005ei ; Patt:2006fw ; MarchRussell:2008yu , through the Kinetic
mixing portal Holdom:1985ag ; Dienes:1996zr ; DelAguila:1993px ; Babu:1996vt ,
neutrino portal Lindner:2010rr ; Falkowski:2009yz ; GonzalezMacias:2015rxl ;
Batell:2017rol ; Batell:2017cmf ; Blennow:2019fhy , axion portal Nomura:2008ru
, or perhaps there is no portal at all in which case the dark sector evolves
independently but may still have observable effects March-Russell:2020nun .

In this paper we focus on the neutrino portal, and in particular examine a
model inspired by the Portalino scenario in which a singlet fermion field
connects gauge neutral fermion operators from the Standard Model (SM) and
hidden sector Schmaltz:2017oov .

In a simple realisation of the Portalino framework introduced in
Schmaltz:2017oov  the SM is supplemented by two additional gauge singlet
fermions and a complex scalar singlet. One of the fermion states plays the
role of the right-handed neutrino, $\nu_{R}$, and couples to the gauge
invariant combination of the SM Higgs and Lepton doublets generating a Dirac
like neutrino mass term after electroweak symmetry breaking. This right-handed
neutrino state also couples to a second gauge invariant operator composed of
the second singlet fermion, which we call $\psi$, and the complex scalar
field, call it $\Phi$. If there is a dark $U(1)$ under which $\Phi$ and $\psi$
both transform then we can construct Yukawa interactions that lead to Dirac
masses after the spontaneous symmetry breaking of the dark $U(1)$ such that

$\displaystyle\mathcal{L}\subset\lambda_{1}\nu_{R}^{\dagger}H_{0}\nu_{L}+\lambda_{2}\nu_{R}^{\dagger}\Phi\;\psi=m_{d}\nu_{R}^{\dagger}(\sin\theta\;\nu_{L}+\cos\theta\;\psi),$
(1)

where $H_{0}$ is the neutral component of the SM Higgs doublet, and the linear
combination of $\psi$ and $\nu_{L}$ forms a massive Dirac state with
$\nu_{R}$. As $\psi$ has vector interactions with the dark gauge sector the
light neutrino (zero) mass eigenstate,
$\nu=\cos\theta\;\nu_{L}-\sin\theta\;\psi$ inherits these interactions albeit
suppressed by a factor of $\sin\theta$. As pointed out in Schmaltz:2017oov ,
this scenario is a specific version of a $Z^{\prime}$ model in which the only
interactions between the new hidden sector $U(1)$ and the SM is via the
neutrinos.111Although, given the introduction of a new scalar field, the Higgs
portal also connects the two sectors.

Introducing DM into the hidden sector is straightforward. For example, in
Schmaltz:2017oov  a Yukawa interaction involving the scalar state $\Phi$ and a
new Dirac fermion, call it $X$, was included. The dark sector dominantly
interacts with the neutrino sector potentially leading to the $X$ DM states
freezing-out via annihilation to neutrinos. This removes, or greatly
suppresses, the usual modes for probing DM in direct and indirect detection
detection experiments allowing for models that consider a wider range of
potentially viable DM masses. On the other hand this makes the model harder to
probe.

The simple model outlined above however requires modification in order to
include neutrino masses. There are a number of choices we can make to do this.
One possibility is to add a Majorana mass term for the $\psi$ field leading to
a model along the lines of the inverse see-saw model, see e.g. Dias:2012xp .
In Liu:2019ixm , it was suggested that it may be possible to produce non-zero
neutrino masses in a Zee-type model including two Higgs doublets via a
$(l_{i}h)(hl_{j})$ term generated at loop level.

An alternative is to change particle content by introducing further
generations of the singlet fermion fields. In this paper, we extend the model
to include two more generations of right-handed neutrino and introduce
associated large Majorana mass terms for these states. This set-up generates
masses for two of the three generations of light neutrinos, with the heavy
Majorana masses suppressing the mass scale of these two mass eigenstates
through a seesaw-like mechanism. Without the heavy Majorana masses, the light
neutrinos will be Dirac states with Dirac neutrino masses of
$\mathcal{O}\left(\lambda_{\nu}v_{h}\right)$. Although with sufficiently small
Yukawa couplings this is in principle a viable model, we choose instead to
adopt the Majorana case.

The introduction of the large Majorana mass scale leads to small mixing angles
in this combined neutrino-hidden state sector, which in turn generates
suppressed couplings for the more massive hidden sector states. This leads to
relatively long lifetimes for these states, giving rise to interesting
cosmological implications and constraints on the model.

In Section 2 we describe the model in full, including the detailed properties
of the putative DM candidate. In Section 3 we outline the model’s predictions
for neutrino masses and mixings, and how the experimentally observed values
can be accommodated. In Section 4.1 we specify the viable parameter space
capable of generating the correct DM abundance. We explore the phenomenology
of - and constraints on - the new hidden sector states, which can have
lifetimes up to and exceeding the age of the universe, in Section 5.

##  2 The Model

The model consists of the SM supplemented by a number of SM singlet fields.
These include three generations of right-handed neutrino, $\nu_{R_{\alpha}}$
($\alpha=$ 1, 2, 3), a complex scalar, $\Phi$, and three Weyl fermions,
$\psi$, $X_{L}$ and $X_{R}$. The $X_{L}$, $X_{R}$ fields will combine to form
a Dirac fermion state and will be our DM candidate. We further introduce a new
abelian gauge symmetry, $U(1)_{\rm d}$, under which $\Phi$, $\psi$ and $X_{R}$
transform each with charge $1/2$. The right-handed neutrinos, $X_{L}$ and all
other SM states are uncharged under the new symmetry. Additionally both
$X_{L}$ and $X_{R}$ are charged under a separate $Z_{3}$ symmetry uncharged.
The role of this $Z_{3}$ is two-fold, firstly this forbids an explicit
Majorana mass term for $X_{L}$ and secondly it ensures the stability of the
$X$ DM state. A summary of these charges is displayed in Table 1.

Table 1: Charge assignments of the field content in the hidden sector under $U(1)_{\rm d}$ and $Z_{3}$. All fields in the table are Standard Model singlets Field | $\nu_{R_{\alpha}}$ | $\Phi$ | $\psi$ | $X_{L}$ | $X_{R}$
---|---|---|---|---|---
$U(1)_{\rm d}$ | 0 | 1/2 | 1/2 | 0 | 1/2
$Z_{3}$ | + | + | + | $e^{i\pi/3}$ | $e^{i\pi/3}$

Given this particle content and charge assignment, the Lagrangian for the
model reads

$\displaystyle\mathcal{L}=\left(-\sqrt{2}\lambda^{\nu}_{\alpha\beta}L^{\dagger}_{\alpha}H\nu_{R\beta}-\sqrt{2}\lambda^{\psi}_{\alpha}\psi^{\dagger}\Phi\nu_{R\alpha}+\frac{i}{2}M_{R\alpha\beta}\nu_{R\alpha}^{T}\sigma_{2}\nu_{R\beta}-\sqrt{2}\lambda_{X}X_{R}^{\dagger}\Phi
X_{L}+{\rm h.c.}\right)$
$\displaystyle+\mu_{H}^{2}|H|^{2}-\lambda_{H}|H|^{4}+\mu_{\Phi}^{2}|\Phi|^{2}-\lambda_{\Phi}|\Phi|^{4}-\lambda_{H,\Phi}|H|^{2}|\Phi|^{2}\ldots,$

where the ellipsis represents the SM Lagrangian terms and all kinetic terms
for the new states including all relevant gauge interactions with the
$U(1)_{\rm d}$ gauge boson, $\omega_{\mu}$, and we specify that
$\mu_{\Phi}^{2},\mu_{H}^{2}>0$. In principle, we may expect a kinetic mixing
term that mixes the field strengths associated with the new $U(1)_{\rm d}$ and
SM hypercharge $U(1)_{Y}$. However, we assume, for simplicity, that this term
is sufficiently small that it does not impact the phenomenology of the
model222Following Gherghetta:2019coi , the leading contribution to the loop
induced kinetic mixing term is $\sim
10^{-7}\left(\lambda^{\nu}\lambda^{\psi}\right)^{2}$, which arises from a
3-loop diagram. Even with the couplings $\lambda^{\nu}$ and $\lambda^{\psi}$
of order 1 the size of the induced kinetic mixing term is not relevant for the
model phenomenology and evades bounds on the size of the kinetic mixing
parameter, Gherghetta:2019coi ..

The form of the potential leads to the spontaneous breaking of $SU(2)_{\rm
L}\times U(1)_{\rm y}\rightarrow U(1)_{\rm em}$ and $U(1)_{\rm d}\rightarrow$
nothing. We parameterise both $\Phi$ and $H$ in terms of excitations,
$\phi^{\prime}$ and $h^{\prime}$ respectively, around the corresponding vacuum
expectation values, expressed in the unitary gauge as
$H=\frac{1}{\sqrt{2}}\begin{pmatrix}0\\\ v_{h}+h^{\prime}\end{pmatrix}$ and
$\Phi=\frac{1}{\sqrt{2}}\left(v_{\phi}+\phi^{\prime}\right)$, where the
expectation values are given by

$\displaystyle
v_{h}^{2}=\frac{2\mu_{H}^{2}\lambda_{\Phi}-\mu_{\Phi}^{2}\lambda_{H,\Phi}}{4\lambda_{H}\lambda_{\Phi}-\lambda_{H,\Phi}^{2}},\qquad
v_{\phi}^{2}=\frac{2\mu_{\Phi}^{2}\lambda_{H}-\mu_{H}^{2}\lambda_{H,\Phi}}{4\lambda_{\Phi}\lambda_{H}-\lambda_{H,\Phi}^{2}}.$
(3)

The Lagrangian after spontaneous symmetry breaking reads

$\displaystyle\mathcal{L}=\left(-M_{d\alpha\beta}\nu_{l\alpha}^{\dagger}\nu_{R\beta}-M_{\psi_{\alpha}}\psi^{\dagger}\nu_{R\alpha}+\frac{i}{2}M_{R\alpha\beta}\nu_{R\alpha}^{T}\sigma_{2}\nu_{R\beta}-m_{X}X_{R}^{\dagger}X_{L}+{\rm
h.c.}\right)$
$\displaystyle+\left(-\lambda^{\nu}_{\alpha\beta}\nu_{l\alpha}^{\dagger}\nu_{R\beta}h^{\prime}-\lambda^{\psi}_{\alpha}\psi^{\dagger}\phi^{\prime}\nu_{R\alpha}-\lambda_{X}X_{R}^{\dagger}\phi^{\prime}X_{L}+{\rm
h.c.}\right)-V(h^{\prime},\phi^{\prime})+\ldots,$ (4)

where
$M_{d\alpha\beta}=\lambda^{\nu}_{\alpha\beta}v_{h},\;\;M_{\psi\alpha}=\lambda^{\psi}_{\alpha}v_{\phi},\;\;m_{X}=\lambda_{X}v_{\phi}$,
and where we have assumed $\lambda_{X}$ is real. The scalar potential now
reads

$\displaystyle V(h^{\prime},\phi^{\prime})=\lambda_{H}v_{h}^{2}h^{\prime
2}+\lambda_{\Phi}v_{\phi}^{2}\phi^{\prime
2}+\lambda_{H,\Phi}v_{h}v_{\phi}h^{\prime}\phi^{\prime}$
$\displaystyle\qquad\qquad\qquad+\lambda_{H}v_{h}h^{\prime
3}+\frac{\lambda_{H}}{4}h^{\prime 4}+\lambda_{\Phi}v_{\phi}\phi^{\prime
3}+\frac{\lambda_{\Phi}}{4}\phi^{\prime 4}$
$\displaystyle\qquad\qquad\qquad+\frac{\lambda_{H,\Phi}v_{h}}{2}h^{\prime}\phi^{\prime
2}+\frac{\lambda_{H,\Phi}v_{\phi}}{2}h^{\prime
2}\phi^{\prime}+\frac{\lambda_{H,\Phi}}{4}h^{\prime 2}\phi^{\prime 2}.$ (5)

The ellipsis in Equation 2 again include the rest of the SM Lagrangian with
the addition of all the beyond the SM kinetic terms and interactions of the
states charged under $U(1)_{\rm d}$ with the associated gauge boson, $\omega$,
whose mass is given by $m_{\omega}=(v_{\phi}\tilde{g})/2$ after symmetry
breaking.

The scalar sector is diagonalised via the rotation defined by

$\displaystyle\begin{pmatrix}h\\\
\phi\end{pmatrix}=\begin{pmatrix}\cos\theta&\sin\theta\\\
-\sin\theta&\cos\theta\end{pmatrix}\begin{pmatrix}h^{\prime}\\\
\phi^{\prime}\end{pmatrix},\qquad\mbox{where}\qquad\tan
2\theta=-\frac{\lambda_{H,\Phi}v_{h}v_{\phi}}{\lambda_{\Phi}v_{\phi}^{2}-\lambda_{H}v_{h}^{2}}.$
(6)

The measured value of the couplings of the SM gauge bosons to the Higgs are
very close to that predicted by the SM and consequently the mixing angle
$\theta$ must be small - in the region of interest
($v_{\phi}\gtrsim$\mathrm{TeV}$$), the limit is approximated by
Lebedev:2021xey ; Falkowski:2015iwa :

$\displaystyle|\sin\theta|\lesssim\frac{0.3}{\sqrt{1+\log\left(\frac{m_{\phi}}{$\mathrm{TeV}$}\right)}}.$
(7)

This can be achieved by insisting $v_{\phi}\gg v_{h}$ and by assuming that the
coupling $\lambda_{H,\phi}$ is moderately suppressed compared with the other
dimensionless couplings in the scalar potential. Suppressing
$\lambda_{H,\phi}$ also has the effect of shutting off the Higgs Portal as a
channel for DM annihilation, see Section 4.1 for details. The $v_{\phi}\gg
v_{h}$ hierarchy is also necessary for achieving light neutrinos with
phenomenologically viable masses.

In this limit the scalar mass eigenstates read

$\displaystyle m_{h}^{2}$
$\displaystyle=2\lambda_{H}v_{h}^{2}\left(1-\frac{\lambda_{H,\Phi}^{2}}{4\lambda_{\Phi}\lambda_{H}}+\mathcal{O}\left(\left(\frac{\lambda_{H,\Phi}v_{h}}{\lambda_{\Phi}v_{\phi}}\right)^{2}\right)\right),$
(8) $\displaystyle m_{\phi}^{2}$
$\displaystyle=2\lambda_{\Phi}v_{\phi}^{2}\left(1+\left(\frac{\lambda_{H,\Phi}v_{h}}{\lambda_{\Phi}v_{\phi}}\right)^{2}+\mathcal{O}\left(\left(\frac{\lambda_{H,\Phi}v_{h}}{\lambda_{\Phi}v_{\phi}}\right)^{2}\frac{\lambda_{H}v_{h}^{2}}{\lambda_{\Phi}v_{\phi}^{2}}\right)\right).$
(9)

Moving to the fermionic content of the model, the first two mass mixing terms
of Equation 2 encode the Portalino-like mixing, as detailed in Equation 1. The
picture is necessarily complicated by the Majorana mass term for the $\nu_{R}$
fields and this is what leads to non-zero light neutrino masses. In this work
we do not propose a full flavour model, instead we assume that there are no
significant hierarchies within the components of
$\lambda^{\nu}_{\alpha\beta}$, $\lambda^{\psi}_{\alpha}$ or $M_{R\alpha}$.
Under this assumption, we define

$\displaystyle\lambda^{\nu}_{\alpha\beta}\equiv\lambda^{\nu}F^{\nu}_{\alpha\beta}\equiv\frac{m_{d}}{v_{u}}F^{\nu}_{\alpha\beta},\;\;\lambda^{\psi}_{\alpha}\equiv\lambda^{\psi}F^{\psi}_{\alpha}\equiv\frac{m_{\psi}}{v_{\phi}}F^{\psi}_{\alpha},\;\;M_{R\alpha\beta}\equiv
m_{R}F_{R\alpha\beta}$ (10)

where the parameters without flavour indices, which we define to be real, will
be used to indicate the typical size of the entries of each term leaving the
precise flavour dependence to the objects labelled $F$.

In order to obtain the correct mass spectrum, we require that the $\nu_{R}$
Majorana mass is much larger than its mass mixing with either $\psi$ or the
active neutrinos $\nu_{L}$, and that the mixing with $\psi$ is much larger
than the mixing with $\nu_{L}$, that is
$\lambda^{\nu}v_{h}\ll\lambda^{\psi}v_{\phi}\ll m_{R}$ or equivalently,
$m_{d}\ll m_{\psi}\ll m_{R}$.

Given this hierarchy of scales, the mass matrix mixing the states
$\psi,\nu_{l}$ and $\nu_{R}$ can be approximately diagonalised via the
following transformations

$\displaystyle\nu_{i}$ $\displaystyle\sim
U_{i\alpha}^{(\nu\nu_{l})}\begin{pmatrix}\nu_{l\alpha}\\\
i\sigma_{2}\nu_{l\alpha}^{*}\end{pmatrix}+\frac{m_{d}}{m_{\psi}}\;U_{i}^{(\nu\psi)}\begin{pmatrix}\psi\\\
i\sigma_{2}\psi^{*}\end{pmatrix}+\frac{m_{d}}{m_{R}}\;U_{i\alpha}^{(\nu\nu_{R})}\begin{pmatrix}-i\sigma_{2}\nu_{R_{\alpha}}^{*}\\\
\nu_{R_{\alpha}}\end{pmatrix},$ $\displaystyle n$
$\displaystyle\sim\frac{m_{d}}{m_{\psi}}U_{\alpha}^{(n\nu_{l})}\begin{pmatrix}\nu_{l_{\alpha}}\\\
i\sigma_{2}\nu_{l_{\alpha}}^{*}\end{pmatrix}+U^{(n\psi)}\begin{pmatrix}\psi\\\
i\sigma_{2}\psi^{*}\end{pmatrix}+\frac{m_{\psi}}{m_{R}}U_{\alpha}^{(n\nu_{R})}\;\begin{pmatrix}-i\sigma_{2}\nu_{R_{\alpha}}^{*}\\\
\nu_{R_{\alpha}}\end{pmatrix},$ $\displaystyle N_{i}$
$\displaystyle\sim\frac{m_{d}}{m_{R}}U_{i\alpha}^{(N\nu_{l})}\begin{pmatrix}\nu_{l_{\alpha}}\\\
i\sigma_{2}\nu_{l_{\alpha}}^{*}\end{pmatrix}+\frac{m_{\psi}}{m_{R}}\;U_{i}^{(N\psi)}\begin{pmatrix}\psi\\\
i\sigma_{2}\psi^{*}\end{pmatrix}+U_{i\alpha}^{(N\nu_{R})}\begin{pmatrix}-i\sigma_{2}\nu_{R_{\alpha}}^{*}\\\
\nu_{R_{\alpha}}\end{pmatrix},$ (11)

where $i,\alpha=1,2,3$ and in the above the definitions in Equation 10 have
been used to factor out the leading order behaviour while the various factors
of $U$ contain all the detailed flavour mixing.

The full diagonalisation of the $(\nu_{l\alpha},\psi,\nu_{R\alpha})$ system is
presented in Appendix A and includes the full expression for the unitary
matrix that Equation 2 derives from, including the explicit form of the $U$
factors.

To leading order the three light neutrinos, $\nu_{i}$, have masses

$\displaystyle m_{\nu_{1}}=0,\quad m_{\nu_{2,3}}\sim\frac{m_{d}^{2}}{m_{R}},$
(12)

and the three heavy neutrinos, $N_{i}$, have masses

$\displaystyle m_{N_{i}}\sim m_{R}.$ (13)

The field $n$, which we choose to call the Portalino333It is not entirely
clear which of our states is the analogue of the Portalino from the earlier
example outlined in Equation 1, where $\nu_{R}$ was the Portalino. It should
perhaps, morally speaking, be the fields $N_{i}$ that should take on the
Portlino title given that their largest component comes from the $\nu_{R}$
fields. We prefer however to adopt the naming conventions from neutrino mass
models where the $N_{i}$s are the heavy neutrinos, the $\nu_{i}$ the light
neutrinos, leaving the $n$ field which we will refer to as the Portalino., has
a mass suppressed relative to the mass scale $m_{\psi}$ given by

$\displaystyle m_{n}\sim\frac{m_{\psi}^{2}}{m_{R}}.$ (14)

With $m_{\psi}=\lambda_{\psi}v_{\phi}$, it is clear that under the hierarchy
assumption of $m_{R}\gg\lambda_{\psi}v_{\phi}$, it must be true that $m_{n}\ll
v_{\phi}$. The full Lagrangian in the mass eigenbasis is given in Appendix B.

##  3 Reconstructing the PMNS matrix

In Appendix A the full masses and mixings of the
$(\nu_{l\alpha},\psi,\nu_{R\alpha})$ system are calculated and presented as
approximate analytic expressions following the hierarchy in masses scales
$m_{d}\ll m_{\psi}\ll m_{R}$. As stated above, we are assuming that there are
no significant hierarchies between the flavours of the individual masses.

Due to the additional states mixing with the left-handed neutrinos, the PMNS
matrix will no longer be unitary. Assuming no mixing in the charged lepton
sector the PMNS matrix is determined by the mixing in the extended neutrino
sector only. The allowed 3$\sigma$ ranges on the entries PMNS matrix (once the
assumption of unitary is dropped) are Parke:2015goa :

$\displaystyle|V|=\begin{pmatrix}0.76\rightarrow 0.85&0.50\rightarrow
0.60&0.13\rightarrow 0.16\\\ 0.21\rightarrow 0.54&0.42\rightarrow
0.70&0.61\rightarrow 0.79\\\ 0.18\rightarrow 0.58&0.38\rightarrow
0.72&0.40\rightarrow 0.78\end{pmatrix}.$ (15)

In addition to the mixing, the masses of the light neutrinos must fall within
the following ranges (the lightest neutrino is massless in this model) -
assuming normal ordering Esteban:2020cvm :

$\displaystyle
m_{2}\in\left[$8.2\text{\,}\mathrm{meV}$,$9.0\text{\,}\mathrm{meV}$\right],m_{3}\in\left[$49.0\text{\,}\mathrm{meV}$,$50.9\text{\,}\mathrm{meV}$\right].$
(16)

A flavour model for the structure of $\lambda^{\nu},\lambda^{\psi}$ and
$M_{R}$ is beyond the scope of this work, and without such a model the task of
finding values for the components of these matrices that satisfy the mixing
and mass constraints is an under-constrained problem.

##  4 Dark Matter Phenomenology

###  4.1 Dark Matter Abundance

Moving now to the DM phenomenology in this model. In our numerical analysis
below we have used micrOMEGAs, Belanger:2018ccd , to compute the freeze-out
abundance for a range of parameter values. We can eliminate a number of
parameters in favour of the measured values of the Higgs mass, $m_{h}$, and
the masses of the SM gauge bosons. The DM phenomenology is not sensitive to
the relative sizes of the neutrino masses and mixings. In order to ensure that
we consider parameter values that can reproduce light neutrino masses we
set444This mass sets the scale for the light neutrinos, the precise masses and
mixings are determined by other flavour parameters that do not play a leading
role in the determination of the DM abundance. In order to numerically
calculate the DM abundance we do need to input some structure by hand and we
assume a simple parameterised form of the components of the full $7\times 7$
neutrino mixing matrix, these are listed in Appendix C.
$(\frac{m_{d}}{m_{\psi}})^{2}m_{n}=$20\text{\,}\mathrm{meV}$$ in order to fix
$\lambda_{\nu}$. The remaining relevant masses and couplings are determined,
at least to leading order, by seven parameters:
$v_{\phi},\lambda_{\psi},\lambda_{\Phi},\lambda_{H,\Phi},m_{n}$, $\tilde{g}$,
and $m_{X}$.

Figure 1: Left: DM abundance against $m_{X}$, for
$v_{\phi}=$2\text{\,}\mathrm{TeV}$$,
$m_{n}=$100\text{\,}\mathrm{keV}$,\lambda_{\psi}=1,\tilde{g}=1,\lambda_{\Phi}=1,\lambda_{H,\Phi}=0.1$.
The horizontal dotted line indicates the observed DM abundance. The first
trough corresponds to the resonant $X\bar{X}\rightarrow nn$ process with
$\omega$ in the s-channel with $m_{X}\approx m_{\omega}/2$
($m_{\omega}=\frac{\tilde{g}}{2}v_{\phi}=$1\text{\,}\mathrm{TeV}$$) and the
second, shallower, trough corresponds to the resonant
$X\bar{X}\rightarrow\omega\omega$ process with $\phi$ in the s-channel with
$m_{X}\approx m_{\phi}/2$
($m_{\phi}=\sqrt{2\lambda_{\phi}}v_{\phi}\approx$2.8\text{\,}\mathrm{TeV}$$).
Middle (Right): Relative (absolute) contributions of each channel to
$\left(\Omega h^{2}\right)^{-1}$. The line labelled ‘Higgs Resonance’ includes
several channels which are only significant near the $h$ resonance at
$m_{X}\approx m_{h}/2$. These are dominated by $X\overline{X}\rightarrow
b\overline{b}$; the next largest contributions come from
$X\overline{X}\rightarrow GG$, $X\overline{X}\rightarrow\tau^{+}\tau^{-}$ and
$X\overline{X}\rightarrow c\overline{c}$.

In the left-hand panel of Figure 1 we demonstrate how the DM abundance behaves
as a function of the DM mass, $m_{X}$, with other parameters fixed at
$v_{\phi}=$2\text{\,}\mathrm{TeV}$,m_{n}=$100\text{\,}\mathrm{keV}$,\lambda_{\psi}=1,\tilde{g}=1,\lambda_{\Phi}=1,\lambda_{H,\Phi}=0.1$
(unless stated otherwise these are the parameter values used in all plots in
this section). The dynamics of the DM freeze-out is rather insensitive to the
exact value of $m_{n}$ provided $m_{n}<m_{X}$ by at least a factor of 10.

The horizontal dotted line in Figure 1 indicates the observed DM abundance. In
the middle and right-hand panel of Figure 1 we show the absolute and relative
contributions of different DM annihilation channels to $\left(\Omega
h^{2}\right)^{-1}$ respectively.

For $m_{X}<m_{\omega}$ the dominant annihilation process is
$X\bar{X}\rightarrow nn$, which proceeds via s-channel exchange of the hidden
sector gauge boson, $\omega$. For $m_{X}>m_{\omega}$, DM annihilation into
pairs of $\omega$ gauge bosons is kinematically possible and becomes the
dominant channel for $m_{X}\sim m_{\phi}/2$ and above.

The annihilation cross section for $X\bar{X}\rightarrow nn$ expanded in powers
of relative velocity, $v$, reads

$\displaystyle\sigma(X\bar{X}\rightarrow nn)v$ $\displaystyle\approx$
$\displaystyle\frac{\tilde{g}^{4}m_{X}^{2}\absolutevalue{U^{n\psi}}^{4}}{128\pi\left[\left(4m_{X}^{2}-m_{\omega}^{2}\right)^{2}+\Gamma^{2}_{\omega}m^{2}_{\omega}\right]},$
(17)

where $\Gamma_{\omega}$ is the total decay width of $\omega$ and we have
assumed555This assumption is only made in the analytic expressions for the
cross section, the mass of the Portalinos is included in the numerical
calculations. $m_{X}\gg m_{n}$. For $X\bar{X}\rightarrow\omega\omega$, the
cross section has the form

$\displaystyle\sigma(X\bar{X}\rightarrow\omega\omega)v$ $\displaystyle\approx$
$\displaystyle\frac{\tilde{g}^{4}\left(m_{X}^{2}-m_{\omega}^{2}\right)^{3/2}}{256\pi
m_{\omega}^{2}m_{X}\left(2m_{X}^{2}-m_{\omega}^{2}\right)}+\frac{v^{2}\tilde{g}^{2}\;F(m_{X},m_{\omega},\tilde{g},\lambda_{X},\theta)}{\left[\left(4m_{X}^{2}-m_{\phi}^{2}\right)^{2}+m_{\phi}^{2}\Gamma_{\phi}^{2}\right]}$
(18)

where $\Gamma_{\phi}$ is the total decay width of $\phi$ and
$F(m_{X},m_{\omega},\tilde{g},\lambda_{X},\theta)$ is a complicated expression
detailed in Appendix D. The $\mathcal{O}(v^{2})$ term in Equation 18 includes
the s-channel diagram with $\phi$ in the intermediate state, that, although
$p$-wave, will dominate the $\sigma(X\bar{X}\rightarrow\omega\omega)$ around
$m_{X}\sim m_{\phi}/2$. The s-wave term in Equation 18 comes from a diagram
with $X$ in the t-channel.

Looking more closely at the left panel in Figure 1, the structure of the plot
is dominated by two resonances, one in each of the channels described above.
The first with an on-shell $\omega$ in the s-channel appearing at $m_{X}\sim
m_{\omega}/2$, and the second with an on-shell $\phi$ appearing at $m_{X}\sim
m_{\phi}/2$.

The middle and right panels of Figure 1 demonstrate over what mass range the
two processes dominantly contribute to the determination of the DM abundance.
In the middle panel we plot the relative contributions of all channels with
more than a $1\%$ contribution and it is clear that the hidden
sector/Portalino only channels dominate.

There are some contributions from SM model channels, all of which are enabled
by the Higgs Portal via the mixing between the SM Higgs and hidden sector
$\phi$. For example, contributions from the $W^{+}W^{-}$, $ZZ$, $hh$ final
state channels are present due to resonant s-channel exchange of the $\phi$
field. The size of these SM channel contributions is ultimately controlled by
the parameter $\lambda_{H,\Phi}$. Suppressing this parameter or even setting
it to zero shuts off the Higgs Portal and removes the contributions from the
SM channels in Figure 1. On one hand this may be desirable as it means the DM
abundance is determined entirely by hidden sector/Portalino physics. The
usually close link in freeze-out models between the annihilation process
determining the abundance and the predicted signal rate in direct and indirect
DM detection experiments is then decoupled. There are however still potential
constraints on this model from the phenomenology of the Portalinos described
in Section 5 and potential signals from indirect detection described in
Section 5.4.

Conversely, if the Higgs Portal is activated by increasing the size of
$\lambda_{H,\Phi}$ the role of the SM states in both generating the DM
abundance (mainly around the $\phi$ resonance) and in constraining the model
become more important and can lead to interesting signals, for example in
indirect detection signals where the DM states are annihilating to SM final
states. These processes are however p-wave and therefore velocity suppressed
and do not trouble current limits.

For sufficiently large masses the process $X\bar{X}\rightarrow\phi\omega$ can
play a role, with a modest dip in the abundance towards
$m_{X}\sim$2\text{\,}\mathrm{TeV}$$. At smaller masses, the Higgs s-channel
resonance can also contribute but only in a very narrow range, as can be seen
in the middle panel of Figure 1. This latter channel is again only present due
to the Higgs Portal and will be reduced if $\lambda_{H,\Phi}$ is further
suppressed below the value of 0.1 used here.

Figure 2:  DM abundance against $m_{X}$ for different values of $\tilde{g}$
(left panel) and $\lambda_{\Phi}$ (right panel). Apart from the values of
$\tilde{g}$ and $\lambda_{\Phi}$ indicated in the plots, both panels display
results for parameter values as stated for the left panel of Figure 1. The
effect of varying the values of both these parameters is seen in the position
and shape of the troughs in the DM abundance (see text for details).

In Figure 2, we demonstrate the dependence of the DM abundance on $\tilde{g}$
(left panel) and $\lambda_{\Phi}$ (right panel). In particular, the way in
which these parameters determine the position and shape of the troughs in the
abundance. With reference to the left panel, $\tilde{g}$ modifies the
abundance in three ways. For $m_{X}<$1.25\text{\,}\mathrm{TeV}$$, the process
$X\bar{X}\rightarrow nn$ dominates the determination of the abundance. With
the cross section for this process going as $\sim\tilde{g}^{4}$, reducing the
value of $\tilde{g}$ increases the abundance, which can be seen in the left
panel of Figure 2.

A second variation arises due to the fact that the value of $\tilde{g}$
determines the mass of the vector boson $\omega$ for fixed $v_{\phi}$ and
hence shifts the position of the resonance in $m_{X}$ and in turn shifts where
the tough appears in the abundance. Decreasing the value of $\tilde{g}$
therefore shifts the trough in the abundance to lower DM masses.

Finally, the width of the trough/resonance depends on $\tilde{g}$ via the
decay rate of $\omega$. A smaller value of $\tilde{g}$ decreases
$\Gamma_{\omega}$ producing a more narrow trough/resonance. The second trough
remains largely unchanged.

In the right panel of Figure 2, the dependence of the DM abundance on
$\lambda_{\Phi}$ is demonstrated. The first trough is unchanged as this is
dominantly determined by the $nn$ final state channel, but the decrease in
$\lambda_{\Phi}$ shifts the second trough to lower $m_{X}$ due to the decrease
in $m_{\phi}$. The width of the $\phi$ resonance/trough is narrower for
smaller $\lambda_{X}$.

We note that the modest dip in the abundance at large $m_{X}$ is no longer
visible in the left panel of Figure 2 when $\tilde{g}$ is decreased. The
reason is that the $X\bar{X}\rightarrow\omega\omega$ annihilation process will
dominate in this mass range due to having a dominant contribution that goes
like $\sim\tilde{g}^{2}\lambda_{X}$ compared with the leading contribution for
the $\omega\phi$ channel, which goes like $\sim\tilde{g}^{4}$. In the left
panel of Figure 2 however, the dip is clearly visible and appears at a lower
value of $m_{X}$ for $\lambda_{\Phi}=0.5$ owing to the reduced value of
$m_{\phi}$.

Finally, we summarise the dependence on the remaining free parameters. For
fixed $v_{\phi}$ the DM abundance doesn’t depend on $\lambda_{\nu}$,
$\lambda^{\psi}$ or $m_{n}$, as can be seen from Equations 17 and 18. There is
a degeneracy in these parameters whereby a change in one can be compensated by
another with no effect on the DM abundance. In Section 5, however, we show
that there are constraints on the Portalino that constrain these parameters of
the model.

For fixed $\tilde{g},\lambda_{X},\lambda_{\psi}$, increasing $v_{\phi}$
increases the masses of the hidden sector states. For DM masses around the
$\omega$ resonance the correct abundance can still be achieved up to
$v_{\phi}\sim 100\;$TeV, whereas the correct abundance for DM masses around
the $\phi$ resonance can be achieved up to $v_{\phi}\sim 7\;$TeV. However, in
both these extreme cases this is only possible if we are precisely on
resonance. Given the assumption that $m_{n}\ll v_{\phi}$ from Section 2 the
viability of the DM model limits the maximum mass of $m_{n}$ to be $\sim$TeV,
the precise limit depending on the mass of the DM state and the degree of
tuning to the resonance one can tolerate.

In summary, we have shown that it is relatively straightforward to reconstruct
the correct DM abundance in this model with the $X$ states freezing-out
dominantly via the annihilation channels: $XX\rightarrow nn$ and
$XX\rightarrow\omega\omega$. There is an important question however about the
fate of the Portalino, $n$, states. It is expected that there is a significant
number density of these states left after the DM states have frozen-out and
all other dark sector states have decayed. The Portalino states are unstable
with potentiality long lifetimes and may disrupt, for example, Big Bang
Nucleosynthesis (BBN) or the Cosmic Microwave Background Radiation (CMBR) as
they decay to SM particles. Constraints coming from the Portalino
phenomenology are discussed in the Section 5. They will also play an important
role in indirect detection as discussed in the Section 5.4.

###  4.2 Direct Detection

Direct DM detection signals can be generated if the Higgs portal is active,
that is the parameter $\lambda_{H,\Phi}$ is non-zero. The dominant
contribution to the direct detection signal comes from Higgs exchange with
scattering cross section per nucleon approximately given by Kurylov:2003ra

$\displaystyle\sigma\sim\frac{m_{\rm
r}^{2}}{2\pi}\left(\frac{\lambda_{X}\sin{2\theta}}{v_{h}m_{h}^{2}}\right)^{2}\;f_{p}^{2},$
(19)

where $m_{\rm r}$ is the reduced mass of the DM-proton given by $m_{\rm
r}=m_{X}m_{p}/(m_{X}+m_{p})$ and we have assumed that the interactions with
protons and neutrons are the same with

$\displaystyle
f_{p}=m_{p}\left[\sum_{u,d,s}f_{T_{q}}+\frac{6}{27}f_{T_{G}}\right]\sim
0.30\;m_{p},$ (20)

where, following Ellis:2018dmb , we have used
$(f_{T_{u}},f_{T_{d}},f_{T_{s}},f_{T_{G}})=(0.018,0.027,0.037,0.917)$.

Assuming a small mixing angle $\theta$ and applying
$\lambda_{\Phi}v_{\phi}^{2}\gg\lambda_{H}v_{h}^{2}$ to
$\tan\left(2\theta\right)$ in Equation 6, we find

$\displaystyle\sigma\sim 5\times 10^{-46}\;{\rm
cm}^{-2}\left(\frac{\lambda_{H,\Phi}}{0.1}\right)^{2}\left(\frac{2\;{\rm
TeV}}{v_{\phi}}\right)^{4}\left(\frac{m_{X}}{2\;{\rm
TeV}}\right)^{2}\left(\frac{1}{\lambda_{\Phi}}\right)^{2}.$ (21)

This value is just below the constraint from LUX-ZEPLIN (LZ), LUX-
ZEPLIN:2022qhg , at 2 TeV. For smaller values of $m_{X}$, the direct detection
limit decreases linearly with decreasing mass (until around 30 GeV where it
flattens off) in contrast, the predicted cross section from Equation 21 with
fixed values of $\lambda_{H,\Phi}$, $v_{\phi}$ and $\lambda_{\Phi}$ decreases
with $m_{X}^{2}$. As a result, masses below 2 TeV are allowed for
$\lambda_{H,\Phi}$, $v_{\phi}$ and $\lambda_{\Phi}$ fixed at the values
indicated in Equation 21.

To get a more general understanding of the direct detection limits Equation 21
can be compared to a linear approximation of the LZ bound (which holds for
$m_{X}\gtrsim$40\text{\,}\mathrm{GeV}$)$ and reads

$\displaystyle\sigma_{\text{max}}\sim 5.5\times 10^{-46}\;{\rm
cm}^{-2}\left(\frac{m_{X}}{2\;{\rm TeV}}\right).$ (22)

Using this, we can write:

$\displaystyle
m_{X}\lesssim$2.2\text{\,}\mathrm{TeV}$\left(\frac{0.1}{\lambda_{H,\Phi}}\right)^{2}\left(\frac{v_{\phi}}{2\;{\rm
TeV}}\right)^{4}\left(\frac{\lambda_{\Phi}}{1}\right)^{2}.$ (23)

Focusing now on parameter values where the observed DM abundance is correctly
reproduced in the model, it is clear from Figure 1 that we need to be near one
of the troughs corresponding to the $\omega$ or $\phi$ resonances. These occur
at $m_{X}=m_{\omega}/2$ and $m_{X}=m_{\phi}/2$ respectively, or equivalently
at $m_{X}/v_{\phi}=\tilde{g}/2$ and $m_{X}/v_{\phi}=\sqrt{\lambda_{\phi}/2}$.
Comparing these to Equation 23, the troughs will be allowed by direct
detection limits if

$\displaystyle\lambda_{H,\Phi}$ $\displaystyle\lesssim$ $\displaystyle
0.15\;\lambda_{\Phi}\;\tilde{g}^{-1/2}\left(\frac{v_{\phi}}{2\;{\rm
TeV}}\right)^{\frac{3}{2}}\qquad(\omega\;\;{\rm resonance}),$ (24)
$\displaystyle\lambda_{H,\Phi}$ $\displaystyle\lesssim$ $\displaystyle
0.12\;\lambda_{\Phi}^{3/4}\;\left(\frac{v_{\phi}}{2\;{\rm
TeV}}\right)^{\frac{3}{2}}\qquad\;\;\;\;\,(\phi\;\;{\rm resonance}).$ (25)

In summary, direct detection can play a role in limiting the allowed parameter
space, but it is always possible to suppress the predicted signal by reducing
the size of $\lambda_{H,\Phi}$. Reducing this parameter has no significant
impact on whether the correct abundance can be achieved. This parameter,
however, cannot be arbitrarily small as it provides the interaction that keeps
the dark sector in thermal equilibrium.

##  5 Portalino Phenomenology

The Portalino mass and dark sector masses are all related to $v_{\phi}$, as a
result the scale of the Portalino mass can be linked to the dark sector
masses. In particular, the mass of the DM particle $X$ can be written in terms
of the Portalino mass as

$\displaystyle m_{X}$
$\displaystyle=\lambda_{X}v_{\phi}=\left(\frac{\lambda^{\nu}\lambda_{X}}{\lambda^{\psi}}\right)\left(\frac{m_{n}}{m_{\nu}}\right)^{\frac{1}{2}}v_{h}$
$\displaystyle\approx$500\text{\,}\mathrm{GeV}$\left(\frac{\lambda_{X}}{0.25}\right)\left(\frac{m_{n}}{$\mathrm{GeV}$}\right)^{\frac{1}{2}}\left(\frac{1}{\lambda^{\psi}}\right)\left(\frac{\lambda^{\nu}}{3.6\times
10^{-5}}\right).$ (26)

Comparing with Figure 1, the observed DM abundance can be produced even in
scenarios with a relatively heavy Portalino. Note though that this does
require a relatively small value of $\lambda_{\nu}$. This is due to the fact
that increasing the Portalino mass $m_{n}$ in Equation 14 while holding
$v_{\phi}$ ($m_{\psi}$) constant is only possible via decreasing the Majorana
mass $m_{R}$. In turn, a decrease in $m_{R}$ necessitates a smaller value of
$\lambda_{\nu}$ in order to obtain the correct neutrino mass scale in Equation
12.

Beyond the phenomenological role Portalinos play in dark matter freeze-out,
indirect detection and neutrino mixing, their presence in the Early Universe
may also lead to significant constraints due to their potentially very long
lifetime. The decay modes for the Portalino are to SM neutrinos; a neutrino
plus neutral meson; a neutrino plus charged lepton pair; or a charged lepton
plus charged meson, the first two via a SM $Z$ boson and the latter two via
the SM $W^{\pm}$. The mixing between the Portalino and neutrino can vary from
the first order approximations described by Equation 2 depending on the
underlying flavour parameters. To account for this in a simplified manner we
include an overall scaling parameter $\eta$ in the Portalino-neutrino mixings
(see Appendix C) and neglect further details of any potential flavour
structure. Within this parametrisation, the lifetime of the Portalino is given
approximately by

$\displaystyle\tau_{n}\sim$3\text{\,}\mathrm{s}$\left(\frac{1}{\eta}\right)^{2}\left(\frac{$\mathrm{GeV}$}{m_{n}}\right)^{4}C(m_{n}),$
(27)

where $C(m_{n})$ accounts for the different decay modes possible for a given
mass, $m_{n}$. The value of $C(m_{n})$ is plotted in Figure 3 and ranges from
$\sim 1$ for $m_{n}<2m_{e}$ to $\sim 0.01$ for a Portalino mass just below the
$W^{\pm}$ mass. The steps down in the plot correspond to mass thresholds of
Standard Models particles into which the Portalino can decay.

Figure 3:  The value of $C(m_{n})$ from Equation 27 plotted as a function of
the Portalino mass $m_{n}$.

The important thing to note is that these decays are cosmologically important,
as any decays into the neutrino or photon sector which occur after neutrino
decoupling affect the neutrino-photon temperature ratio (and hence
$N_{\text{eff}}$). Furthermore there are constraints on long lived decaying
particles from the primordial abundance of light elements set during BBN,
Kawasaki:2017bqm .

This relationship between the Portalino mass and its coupling to the SM, in
particular to the SM leptons, determines the Portalino decoupling temperature.
Portalinos are primarily held in thermal equilibrium by processes such as
$en\leftrightarrow e\nu$, which depends on a coupling of order
$\mathcal{O}\left(m_{d}/m_{\psi}\right)$ or in terms of the physical mass
states $\;\mathcal{O}\left(\sqrt{m_{\nu}/m_{n}}\right)$. The more massive the
Portalino state is the weaker its coupling with the SM becomes and the earlier
it will decouple, for example in the limit where the Portalino decouples
relativistically, the rate of $en\leftrightarrow e\nu$ is given approximately
by

$\Gamma\sim\left\\{\begin{array}[]{ll}\eta^{2}\left(\frac{m_{\nu}}{\pi^{3}m_{n}}\right)G_{F}^{2}T^{5},&\mbox{for}\;\;T<m_{Z},\\\\[8.53581pt]
\eta^{2}\left(\frac{m_{\nu}}{\pi^{3}m_{n}}\right)G_{F}^{2}m_{Z}^{2}T^{3}&\mbox{for}\;\;T>m_{Z},\end{array}\right.$
(28)

where we have neglected the masses of the electron and neutrino and have
quoted the result for decoupling temperatures above and below the mass of the
$Z$ (and $W$) SM gauge boson.

By comparing these rates to the Hubble parameter we can approximately
determine the decoupling temperature of the Portalinos as

$\displaystyle T_{n,{\rm
decouple}}\sim\left\\{\begin{array}[]{ll}$23\text{\,}\mathrm{GeV}$\left(\frac{1}{\eta}\right)^{\frac{2}{3}}\left(\frac{m_{n}}{$1\text{\,}\mathrm{GeV}$}\right)^{\frac{1}{3}},&\mbox{for}\;\;T_{n,{\rm
decouple}}<m_{Z},\\\\[8.53581pt]
$165\text{\,}\mathrm{GeV}$\left(\frac{0.1}{\eta}\right)^{2}\left(\frac{m_{n}}{$1\text{\,}\mathrm{GeV}$}\right)&\mbox{for}\;\;T_{n,{\rm
decouple}}>m_{Z},\end{array}\right.$ (31)

where we have set the total number of effectively massless degrees of freedom
$g_{*}\sim 100$ and $m_{\nu}=0.2\;$meV.

For these example parameter values the Portalino decouples while relativistic
with a significant energy density. Any Portalino decays producing neutrinos
that occur after neutrino decoupling and before/during BBN or recombination
would affect $N_{\text{eff}}$ during these times. Any shift in
$N_{\text{eff}}$ is tightly constrained Fields:2019pfx ; Aghanim:2018eyx . In
addition, sufficiently long lived Portalinos decaying into SM particles during
BBN directly impact the abundance of light elements, Kawasaki:2017bqm ,
further constraining the Portalino parameter space. The constraints from BBN
leave us with two options: the Portalinos must decay before BBN and neutrino
decoupling, or after recombination. The latter possibility can be further
split into two scenarios: one in which the Portalinos decay after
recombination, and another in which the Portalinos don’t decay within the
lifetime of the universe. These scenarios will be discussed in the following
sections.

###  5.1 The Heavy Portalino Case

We first focus on the heavy Portalino case in which Portalinos decay most
rapidly. For a given mass, the value of $\eta$ will determine the decoupling
temperature and decay time of the Portalino and the condition that the
Portalinos must decay before BBN corresponds to a minimum value of $\eta$. The
parameter $\eta$ also feeds into the mixing between the Portalino and the
active neutrinos (see Appendix A for details). Limits on the maximum size of
this mixing comes from electroweak precision tests, collider searches for the
direct production of Portalinos, beam dump experiments, and measurements of
meson decays such that combined with the BBN constraints limits $\eta$ to a
range of allowed values. The extent of this range narrows with decreasing
Portalino mass and shrinks to zero at $m_{n}=$489\text{\,}\mathrm{MeV}$$.

To show this, we first evaluate the combined constraints on $|V_{e4}|^{2}$ for
Portalino masses between $0.1\text{\,}\mathrm{GeV}$ and
$100\text{\,}\mathrm{GeV}$, a summary of which is shown in Figure 4. The
Portalino is not expected to mix more strongly with any particular neutrino so
only constraints for $|V_{e4}|^{2}$ are shown, as these are strongest.

Figure 4: Current and future bounds on heavy Portalino mixing with the
electron neutrino, combined with the constraint that the Portalino must decay
before neutrino decoupling. See text for details of constraints such as
Collider, EWPT etc. The line labelled $|V_{e4}(m_{n})|^{2}$ indicates an
approximate expected size of the Portalino-neutrino mixing with $\eta=1$
leaving $|V_{e4}(m_{n})|^{2}\sim m_{\nu}/m_{n}$. The four colour graded lines
(for
$m_{n}=$489\text{\,}\mathrm{MeV}$,$500\text{\,}\mathrm{MeV}$,$1\text{\,}\mathrm{GeV}$$
and $2\text{\,}\mathrm{GeV}$), projected into the $m_{n}$-$|V_{e4}|^{2}$ plane
are lines on which the model lives as a function of $\eta$.

Electroweak precision tests: The Portalino mixing with active neutrinos can
affect several electroweak observables such as the invisible $Z$ decay width.
These effects are primarily dependent on the size of the Portalino-neutrino
mixing $V_{n\nu}$, but there is some mass dependence for lower values of
$m_{n}$. Global fits of sterile neutrino mixing have been performed on
electroweak precision data (taken from Deppisch:2015qwa , which draws from
delAguila:2008pw ; Akhmedov:2013hec ; Basso:2013jka ; deBlas:2013gla ;
Antusch:2015mia ), and these bounds can be applied to the Portalino. This
constraint is displayed in Figure 4 as a green coloured region labelled
‘EWPT’.

Collider searches: Portalinos can be produced directly (e.g. via
$e^{+}e^{-}\rightarrow n\nu$ or $pp\rightarrow W^{*}\rightarrow nl^{\pm}$) or
via $Z$-boson decays. They could then decay into visible products, possibly
with detectable displaced vertices. Searches for such decays have been carried
out using data from LEP L3:1992xaz ; DELPHI:1996qcc , ATLAS ATLAS:2019kpx ,
and CMS CMS:2018jxx ; CMS:2018iaf ; CMS:2022fut . Limits have also been
projected for future experiments such as MATHUSLA Curtin:2018mvb , FCC-ee
Antusch:2017pkq  and ILC Antusch:2017pkq . These constraints are displayed in
Figure 4 as a red region labelled ‘Collider’, along with projected limits for
future experiments labelled ‘MATHUSLA’, ‘FCC-ee’ and ‘ILC’.

Beam dump experiments: Portalinos with a mass of around
$1\text{\,}\mathrm{GeV}$ can have significant lifetimes and so may decay at
some distance from the production site. Visible decay products can be searched
for in beam dump experiments with the detector positioned at a distance from
the production site. Many such experiments have been carried out NOMAD:2001eyx
; Ruchayskiy:2011aa ; Atre:2009rg ; Castro:2013jsn ; Yuan:2013yba ;
Wang:2014lda ; LHCb:2014osd . These constraints are displayed in Figure 4 as a
blue region labelled ‘Beam dump’, along with projected limits for future
experiments labelled ‘DUNE’ LBNE:2013dhi  and ‘SHiP’ SHiP:2018xqw . Note that
CHARM and PS191 bounds have been adjusted to account for the Majorana nature
of the Portalino, where the bounds are twice as strong in this case
SHiP:2015gkj .

Meson decays: The Portalino may take part in charged meson decays such as
$X^{\pm}\rightarrow l^{\pm}n$, with a branching ratio proportional to
$|V_{n\nu}|^{2}$. This would manifest as an additional peak in the charged
meson decay spectrum. Constraints from decays such as $\pi^{+}\rightarrow
e^{+}\nu$ are compiled in Bryman:2019ssi ; Deppisch:2015qwa , these
constraints are displayed in Figure 4 labelled as ‘$\pi\rightarrow e\nu$’.
Additionally the Belle experiment, which searched for the decay $B\rightarrow
XlN$ or $B\rightarrow lN$ followed by $N\rightarrow l\pi$ (where $N$ is a
sterile neutrino and X is a meson), would also place constraints on the
neutrino-Portalino mixing Belle:2013ytx . This constraint is the dark pink
region labelled ‘Belle’.

Lepton number violation in meson decays: The Majorana mass term violates
lepton number. Hence in the Portalino model lepton number violating (LNV)
processes such as $K^{+}\rightarrow l^{+}l^{+}\pi^{-}$ may take place. Many
searches for LNV processes have been carried out (e.g. Atre:2009rg ). However,
the bounds from lepton number violation are weaker than other limits and are
not shown in Figure 4.

BBN and $N_{\text{eff}}$: Starting with the impact on $N_{\text{eff}}$ from
Portalino decays into neutrinos. Neutrinos decouple from the rest of the
thermal bath before electron-positron annihilation, and hence the entropy from
electrons and positrons is transferred into the photons alone, raising the
photon temperature relative to the neutrino temperature. In the standard case
this leads to the ratio
$T_{\nu}\approx\left(\frac{4}{11}\right)^{\frac{1}{3}}T_{\gamma}$. However,
Portalinos can decay into neutrinos and charged leptons, so if the Portalinos
decay after neutrino decoupling they will alter the neutrino-photon
temperature ratio. A convenient way to parameterise this is as a constraint on
the effective number of neutrino species $N_{\text{eff}}$: $\Delta
N_{\text{eff}}=N_{\text{eff}}-N^{{}^{\prime}}_{\text{eff}}$, which is
constrained to be less than $0.16$ at BBN Fields:2019pfx , and less than
$0.33$ at recombination Aghanim:2018eyx  (where
$N^{{}^{\prime}}_{\text{eff}}=3.046$ is the SM prediction deSalas:2016ztq ).

The form of $N_{\text{eff}}$ can be defined via the relationship between the
total radiation energy density and the energy density in photons:

$\displaystyle\rho_{r}$
$\displaystyle=\rho^{{}^{\prime}}_{\gamma}+\rho^{{}^{\prime}}_{\nu}+\Delta_{\rho}$
(32)
$\displaystyle=\rho_{\gamma}\left(1+\frac{7}{8}N_{\text{eff}}\left(\frac{4}{11}\right)^{\frac{4}{3}}\right).$
(33)

where $\Delta\rho$ is the energy density due to Portalino decay products, and
the ′ superscript refers to quantities ignoring any Portalino contributions.

The size of $\Delta\rho$ depends on whether the Portalinos decouple
relativistically or remain in thermal equilibrium long enough such that they
freeze-out non-relativistically with a Boltzmann suppressed abundance. The
latter scenario where the Portalinos freeze-out non-relativistically does not
lead to a modification of $N_{\text{eff}}$ during BBN (or later) as the
increased size of the coupling required to keep the Portalinos in thermal
equilibrium long enough leads to a short Portalino lifetime meaning all states
will have decayed well before BBN.

For Portalinos decoupling while relativistic (and assuming they decay at a
temperature $T_{n\text{,decay}}<m_{n}$), the energy density due to the decay
products at (photon thermal bath) temperatures $T$ is given by

$\displaystyle\Delta\rho(T)$
$\displaystyle=m_{n}n_{n}\left(T_{n\text{,decay}}\right)\left(\frac{T}{T_{n\text{,decay}}}\right)^{4}$
$\displaystyle=\frac{3m_{n}\zeta(3)}{2\pi^{2}}T_{n\text{,decouple}}^{3}\left(\frac{a_{n\text{,decouple}}}{a_{n\text{,decay}}}\right)^{3}\left(\frac{T}{T_{n\text{,decay}}}\right)^{4},$
(34)

where $n_{n}$ is the Portalino number density and the Riemann zeta function
$\zeta(3)\approx 1.2$. Applying conservation of entropy we find

$\displaystyle\frac{\Delta\rho(T)}{T^{4}}=\frac{3\zeta(3)}{2\pi^{2}}\frac{g_{*}(T_{n\text{,decay}})}{g_{*}(T_{n\text{,decouple}})}\frac{m_{n}}{T_{n\text{,decay}}}.$
(35)

Assuming the Portalinos instantaneously decay at (thermal bath) temperature
$T_{n\text{,decay}}$, they deposit energy densities $\beta\Delta\rho$ and
$(1-\beta)\Delta\rho$ into the neutrino and photon sectors respectively. The
resulting change in $N_{\text{eff}}$ reads

$\displaystyle\Delta
N_{\text{eff}}=\frac{15\left(\beta\frac{8}{7}\left(\frac{11}{4}\right)^{\frac{4}{3}}-(1-\beta)N^{{}^{\prime}}_{\text{eff}}\right)}{\pi^{2}+15(1-\beta)\frac{\Delta\rho}{(T_{n\text{,decay}})^{4}}}\frac{\Delta\rho}{(T_{n\text{,decay}})^{4}}.$
(36)

For simplicity in our analysis, we make the conservative choice of $\beta=1$
when producing constraints as the Portalino abundance is so large that $\Delta
N_{\text{eff}}\gg 0.16$ whenever Portalinos decay after neutrino decoupling,
regardless of the value of $\beta$ (as long as $\beta\gtrsim 0.4$, below which
$\Delta N_{\text{eff}}$ becomes negative). The constraint is therefore that
any decay occurring after neutrino decoupling is ruled out. This constraint is
the grey region labelled “$N_{\text{eff}}$" in Figure 4.

Further to the constraints on changes to $N_{\text{eff}}$, Portalinos decaying
into SM states during BBN may directly impact the abundance of light elements.
The yield of Portalinos that decouple while relativistic will be equal to the
equilibrium yield, $Y_{n}^{\rm EQ}\sim 0.4/g_{*}(T_{\rm n,decouple})$. Even
for very high decoupling temperatures, $g_{*}(T_{\rm n,decouple})$ will at
most be $\sim 10^{2}$ meaning the yield of decaying Portalinos will be large.
For values of $m_{n}Y_{n}^{\rm EQ}\gtrsim 10^{-8}$, the lifetime of Portalinos
decaying and producing quark-antiquark pairs, for example, is restricted to
less than $\sim$ 0.03 seconds, Kawasaki:2017bqm . As with the $N_{\text{eff}}$
constraint, Portalinos decoupling when non-relativistic with a Boltzmann
suppressed abundance have short lifetimes that means they will have decayed
well before BBN. The constraint on relativistically decoupling Portalinos is
given by the pink region labelled “BBN” in Figure 4, where this constraint
continues “behind” the grey region towards the bottom left hand side of the
plot.

In Figure 4 the model lives on the vertical multi coloured lines. Each line
corresponds to a different Portalino mass with values 2 GeV, 1 GeV, 500 MeV
and 489 MeV plotted. The colour gradient on these lines represents the
changing values of $\eta$ moving from large values at the top of the figure
down to small values at the bottom. As the mass of the Portalino is lowered
the vertical model line moves towards the left and for $m_{n}=489\;$MeV the
full line is completed excluded meaning that we require Portalino masses with
$m_{n}>489\;$GeV in this scenario.

Figure 5: Portalino decoupling temperature against Portalino lifetime, for a
range of Portalino masses. Marked on the plot are: the decoupling temperature
as a function of lifetime, $T_{n,\text{decouple}}(\tau)$, parameterised by
$\eta$ (the enhancement/suppression of the Portalino-neutrino mixings); the
region in pink and to the right of this is ruled out by the impact of the
Portalino decays on the abundance of light elements during BBN and the grey
region indicates where the Portalino decay temperature $T_{n,\text{decay}}$ is
equal to or below the neutrino decoupling temperature
$T_{\nu,\text{decouple}}$, within this region the Portalino decays impact on
the effective relativistic degrees of freedom $N_{\text{eff}}$; and the region
where the Portalino-neutrino mixing is larger than allowed by constraints is
indicated by the light blue shading. The allowed region lies on the
multicoloured line $T_{n,\text{decouple}}(\tau)$, between the blue shading and
the pink region.

Figure 5 illustrates the same model constraints but now in terms of the
Portalino decoupling temperature and lifetime. The model parameter space lies
along the diagonal multicoloured line labelled
$T_{n\text{,decouple}}\left(\tau\right)$. This line was calculated numerically
using the exact tree-level cross section expression for the process
$en\leftrightarrow e\nu$ using FeynCalc, Shtabovenko:2020gxv , in order
evaluate the rate $\Gamma(T,\eta)=n_{n}\langle\sigma v\rangle(T,\eta)$, where
$n_{n}$ is the number density of Portlainos.

The decoupling temperature for each value of $\eta$ was found by equating this
rate to the Hubble parameter. The lifetime of the Portalino is also evaluated
numerically as a function of $\eta$ to produce the multicoloured line, where
the gradient of colours represents the size of $\eta$. Moving to the left on
this model line the size of $\eta$ increases and as a result so does the
corresponding Portalino-neutrino mixing. The blue shading around the line
indicates values of $\eta$ ruled out by the constraints on
$\absolutevalue{V_{e4}}^{2}$.

On the right hand side of all plots in Figure 5 the constraints on the
Portalino lifetime stemming from $\Delta N_{\text{eff}}$ (grey region,
labelled “$\Delta N_{\rm eff}>0.16$”) and the abundance of light elements
during BBN (pink region, labelled “BBN”) are shown. The BBN constraints
continue to longer lifetimes behind the grey region.

From the fourth plot in Figure 5 the limits from $\absolutevalue{V_{e4}}^{2}$
meet those from BBN ruling out all values of $\eta$ for Portalino masses equal
to or less than 489 MeV indicating that for the heavy Portalino case we have
viable parameter space for $m_{n}>489$MeV for $\mathcal{O}(1)$ values of the
flavour parameter $\eta$.

###  5.2 Intermediate Portalino

Decreasing the mass of the Portalino (and/or decreasing $\eta$) allows for
their decays to occur after recombination. This means that they don’t affect
the neutrino-photon temperature ratio at recombination, and hence they evade
constraints on $N_{\text{eff}}$ at this point, potentially opening up an
additional region of parameter space. However, we will show that this set-up
tends to lead to an early extra period of (Portalino) matter domination, and a
universe which, at the present day temperature, has an energy density and
expansion rate that is too high.

Firstly, $\eta$ increases the decay rate of the Portalino. The condition that
the Portalinos must decay after recombination can be recast as a condition on
$\eta$:

$\displaystyle\textbf{Condition 1:
}\tau>t_{\text{recombination}}\implies\eta\lesssim
160\left(\frac{$100\text{\,}\mathrm{keV}$}{m_{n}}\right)^{2}.$

Next, we can consider the time dependence of the expansion of the universe.
Similarly to the heavy case, in this scenario the Portalino tends to decouple
while relativistic and with a significant number density. If
$m_{n}\gtrsim$100\text{\,}\mathrm{eV}$$ the energy density in Portalinos comes
to dominate the universe until they decay. This allows us to place a lower
limit on $\eta$ given that the smaller the value of $\eta$ the longer lived
the Portalinos are. An increase in the Portalino lifetime increases the length
of the period of early matter domination, and leads to a lower temperature (or
equivalently, a larger scale factor $a(t)$) at the point that Portalinos
decay. Under the assumption that $\tau>t_{\text{recombination}}$, this doesn’t
leave enough time to reach matter domination between Portalino decays and the
point when the temperature of the universe reaches $T_{0}$ (i.e. the present
day temperature). This can be seen from the following approximation for the
scale factor at matter-radiation equality (where radiation includes Portalino
decay products):

$\displaystyle a_{\rm MRE}\approx 5.8\times
10^{3}\left(\frac{61.75}{g_{*}\left(T_{n\text{,decouple}}\right)}\right)^{\frac{2}{9}}\eta^{-\frac{4}{3}}\left(\frac{$100\text{\,}\mathrm{keV}$}{m_{n}}\right)^{\frac{7}{3}},$
(37)

where $a_{\rm MRE}>1$ would mean that the present day is radiation-dominated,
with a higher expansion rate $H_{0}$ than observed. The condition that $a_{\rm
MRE}<1$ (i.e. that matter-radiation equality is reached before the present
day) can be translated into a constraint on $\eta$:

$\displaystyle\textbf{Condition 2: }a_{\rm MRE}<1\implies\eta\gtrsim
670\left(\frac{61.75}{g_{*}\left(T_{n\text{,decouple}}\right)}\right)^{\frac{1}{6}}\left(\frac{$100\text{\,}\mathrm{keV}$}{m_{n}}\right)^{\frac{7}{4}}.$

The final condition that must be taken into account is that $|V_{e4}|^{2}<1$
and so $\eta$ cannot be too large:

$\displaystyle\textbf{Condition 3: }|V_{e4}|^{2}<1\implies\eta\lesssim
2.2\times
10^{3}\left(\frac{m_{n}}{$100\text{\,}\mathrm{keV}$}\right)^{\frac{1}{2}}.$

Conditions 1, 2 and 3 cannot be simultaneously satisfied, for any value of
$m_{n}$. This is equivalent to the statement that the existence of a long-
lived Portalino ($\tau>t_{\text{recombination}}$) which comes to dominate the
universe inevitably leads to a current-day universe which is dominated by
Portalino decay products (or Portalinos themselves), and is growing more
quickly than we observe. Hence this scenario is ruled out.

This only leaves the possibility that the initial Portalino density is so low
(via low mass and/or density) that they never come to dominate the universe,
which brings us on to the next section - the light Portalino.

###  5.3 Light Portalino

The final possibility is a very light
($m_{n}\lesssim$10\text{\,}\mathrm{eV}$$) Portalino. Similarly to the above
cases, DM freezes out at a temperature $T\sim\text{few hundred
}$\mathrm{GeV}$$ following this, the Portalinos decouple (possibly long)
before the QCD phase transition. Again, there is still a significant
population of Portalinos after decoupling. However, unlike in either of the
above cases the light Portalinos never come to dominate the energy density of
the universe, and tend not to decay within the lifetime of the universe. They
will however contribute $\Delta N_{\text{eff}}$ at BBN and will behave like
light sterile neutrinos and will be constrained by measurements of the CMB by
Planck Aghanim:2018eyx .

There are several other bounds for this scenario coming from the Portalino-
neutrino mixing, e.g. those that arise from beta decay experiments (see for
example deGouvea:2015euy ). However these bounds are far weaker than the
constraints from $\Delta N_{\text{eff}}$ and Planck.

As the light Portalino will be relativistic until well after BBN, its
contribution to the energy density at BBN will follow

$\displaystyle\rho_{n}=\frac{7\pi^{2}}{120}T_{n}^{4}=\frac{7\pi^{2}}{120}\left(\frac{g_{*s}(T_{\nu\text{,decouple}})}{g_{*s}(T)}\right)T^{4},$
(38)

where $g_{*s}(T_{\nu\text{,decouple}})$ is the number of relativistic degrees
of freedom at neutrino decoupling and $T$ is the temperature of the Universe.
As above,

$\Delta
N_{\text{eff}}=\frac{120\rho_{n}}{7\pi^{2}T_{\nu}^{4}}=\left(\frac{T_{n}}{T_{\nu}}\right)^{4}=\left(\frac{g_{*s}(T_{\nu\text{,decouple}})}{g_{*s}(T_{n\text{,decouple}})}\right)^{\frac{4}{3}}.$

Imposing $\Delta N_{\text{eff}}<0.16$ Fields:2019pfx  (and inserting
$g_{*s}(T_{\nu\text{,decouple}})=43/4$), implies that
$g_{*s}(T_{n\text{,decouple}})>42.5$, or equivalently
$T_{n\text{,decouple}}\gtrsim$150\text{\,}\mathrm{MeV}$$ for the light
Portalino. Comparing this to Equation 31 this implies

$\displaystyle\eta\lesssim
0.042\left(m_{n}/$10\text{\,}\mathrm{eV}$\right)^{\frac{1}{2}}.$ (39)

Figure 6: Constraints on the light Portalino scenario with the blue region
ruled out by constraints on $\Delta N_{\rm eff}$ at BBN and the grey region
ruled out by Planck’s analysis of the CMB anisotropies Aghanim:2018eyx . The
vertical dashed line indicates the scale of neutrino masses and is included to
highlight the limitations of the model assumption that the Portalino is more
massive than the neutrinos.

The second constraint on this scenario comes from Planck’s determination of
cosmological parameters from measurements of the CMB anisotropies, which
combines data from temperature and polarisation maps with lensing and Baryon
Acoustic Oscillation (BAO) measurements. In particular, the light
($m_{n}<$10\text{\,}\mathrm{eV}$$) Portalino is constrained by the Planck
TT,TE,EE+lowE+lensing+BAO analysis limit on an effective sterile neutrino
mass, $m_{\rm eff}$, where
$m_{\text{eff}}=\Omega_{\text{sterile}}h^{2}\left($94.1\text{\,}\mathrm{eV}$\right)$,
is constrained to be less than $0.65\text{\,}\mathrm{eV}$ Aghanim:2018eyx .

Applying this to our case the light Portalino abundance reads

$\displaystyle\Omega_{n}=\frac{4\zeta(3)GT_{0}^{3}}{\pi
H_{0}^{2}}\frac{g_{*S}(T_{0})}{g_{*S}(T_{n\text{,decouple}})}m_{n},$ (40)

where $G$ is the gravitational constant. Hence this limit is almost entirely a
constraint on $m_{n}$ alone, with a small adjustment depending on
$T_{n\text{,decouple}}$. For example, if the Portalinos decouple extremely
early, before top quarks, $m_{\text{eff}}<$0.65\text{\,}\mathrm{eV}$$
translates to $m_{n}\lesssim$6.4\text{\,}\mathrm{eV}$$.

In Figure 6 the combined limits on the light Portalino scenario is mapped out
as a function of $\eta$ and $m_{n}$ with the coloured regions ruled out. The
blue region is ruled out due to the Portalino contributing too much to $\Delta
N_{\rm eff}$ at BBN with the shape of the region determined by Equation 39.
The grey region represents the parameter region ruled out by the Planck
constraint on $m_{\rm eff}$, where the stepped shape comes from the
temperature dependence of $g_{*S}$.

###  5.4 Portalinos and the Indirect Detection of Dark Matter

As outlined above there are two DM mass regions of interest corresponding to
resonant annihilation processes in which the model can generate the observed
DM relic abundance. For DM masses around $m_{\omega}/2$ the resonant
annihilation into $nn$ pairs dominates and for masses around $m_{\phi}/2$ the
resonant annihilation into $\omega\omega$ pairs dominates. The $\omega$ states
decay quickly into $nn$ pairs and as a result for DM states with masses around
the $\phi$ resonance, DM annihilation results in the production of four $n$
states.

Resonant DM annihilation into $nn$ pairs is an s-wave process and is therefore
DM velocity independent. This means that the DM annihilation rate in the
Galactic Centre or Dwarf Spheroidal Galaxies will be the same as that during
freeze-out. In the case where the $\phi$ resonance dominates, the resonant
part of the annihilation cross section to $\omega\omega$ is p-wave and with
the velocity in astrophysical environments, such as the Galactic Centre, at
$v\sim 10^{-3}$ (or lower in Dwarf Spheroidal Galaxies), the indirect
detection signals coming from such a process will be velocity suppressed and
play no role in constraining this mass region. There is an s-wave contribution
to the $\omega\omega$ channel but this is not resonant, with an annihilation
rate for our parameter values well below the ‘thermal relic cross section’
($\sim 2\times 10^{-26}{\rm cm}^{3}{\rm s}^{-1}$) that is required to generate
the observed relic abundance and so does not lead to constraints from indirect
detection searches.

Focusing for the rest of this Section on the s-wave DM annihilation into $nn$
pairs, the first important point is that the Portalinos produced are mono-
energetic, each with an energy around the mass of the DM particle. The
indirect detection phenomenology is then determined by whether or not the
Portalinos decay before they reach Earth. In the case of heavy
($m_{n}\gtrsim$489\text{\,}\mathrm{MeV}$$) Portalinos produced in DM
annihilations, the allowed lifetime can be up to $0.03\;$seconds (see Section
5.1), with the Portalino travelling on average a distance,

$d\sim 10^{10}\;{\rm m}\left(\frac{m_{X}}{500\;{\rm
GeV}}\right)\left(\frac{489\;{\rm MeV}}{m_{n}}\right),$

where we have assumed the Portalinos are ultra relativistic, that is,
$m_{n}\ll m_{X}$. For Portalino masses closer to the DM mass the distance
travelled will be significantly shorter. Hence, for the heavy Portalino case,
the Portalinos will decay before reaching Earth, even if the DM annihilation
has taken place in the Sun rather than the Galactic Centre or Dwarf Spheroidal
Galaxy.

The Portalinos decay via virtual SM gauge bosons, $n\rightarrow\nu_{i}Z\;{\rm
or}\;l^{-}W^{+}$, with the $Z$ and $W$ decaying either hadronically or
leptonically leading to indirect DM signals, in particular in gamma rays.
Fermi-LAT searched for gamma rays produced in DM annihilations in Milky Way
Dwarf Spheroidal Galaxies, Fermi-LAT:2015att  and taking the most constrained
case of DM annihilating into $\bar{b}b$ final states only they were able to
rule out the thermal relic cross section for DM masses below $100\;$GeV, see
e.g. Figures 1 and 2 of Fermi-LAT:2015att . Above this mass the limit weakens
and the thermal relic cross section is not constrained. Limits are also
presented for DM that dominantly annihilates into $\tau^{+}\tau^{-}$ pairs,
the resulting lower bound on the DM mass is only marginally lower, still of
order 100 GeV.

For $m_{n}>2m_{b}$ DM in the Portalino model does not dominantly annihilate
into $\bar{b}b$ pairs, but we can use the Fermi-LAT limit to place a
conservative bound on the model parameter space. For example, the lower bound
on the DM mass of 100 GeV leads to $\lambda_{X}\gtrsim 0.1\;({\rm
TeV}/v_{\phi})$ and given that in order to produce the correct relic abundance
we must be near the $\omega$ resonance, such that $m_{X}\sim m_{\omega}/2$, we
have that $\tilde{g}\gtrsim 0.4\;({\rm TeV}/v_{\phi})$. For Portalino masses
between the bottom mass and the $\tau$ mass we will find similar constraints
but for masses below the mass of the $\tau$, the constraints weaken
considerably with the lower bound on the DM mass dropping to around 10 GeV,
see for example MAGIC:2016xys  where constraints are derived for the case
where the $\mu^{+}\mu^{-}$ annihilation channel dominates. Due to the
Portalinos being produced on shell it is their mass that determines the
spectrum of SM states produced rather than the mass of the DM, in contrast to
the usual scenario of the DM annihilating directly to SM states. This offers a
way to evade the indirect limits.

The projected limits from the Cherenkov Telescope Array (CTA) will reach the
thermal relic cross section (and below) for DM masses between a few 100 GeV
and just over a TeV assuming DM annihilates dominantly into $W^{+}W^{-}$ or
$\bar{b}b$ final states, CTA:2020hii . These projected limits will be weakened
for DM annihilating into decaying Portalinos that produce a spectrum of
different SM final states or the particular the case of Portalinos with masses
below the $b$-quark mass. A dedicated study of the gamma ray flux produced in
this mode is required to understand both the current limits and the future
constraints in detail. This includes the impact of the Portalinos travelling
significant distances from the DM annihilation site before decaying into SM
states.

Beyond the signal from gamma rays, neutrinos produced in Portalino decays can
also be searched for via neutrino telescopes. The DM annihilation cross
section limits from IceCube IceCube:2015rnn  and ANTARES ANTARES:2022aoa  are
at least an order of magnitude above the thermal relic cross section for all
DM masses we consider and so do not constrain the model. For example,
IceCube’s most stringent limit for DM with masses between 30 GeV and 10 TeV
annihilating into a pair of SM neutrinos is $\sim 4\times 10^{-24}{\rm
cm}^{3}{\rm s}^{-1}$, IceCube:2015rnn .

For the lighter Portalino case ($m_{n}\lesssim$6.4\text{\,}\mathrm{eV}$$), the
Portalinos do not decay at all and can travel all the way to the Earth.
Neutrino search experiments are potentially capable of detecting Portalinos
e.g. IceCube Abbasi:2021miz  and ANTARES ANTARES:2022aoa , but with reduced
sensitivity as the Portalinos will have a coupling to the SM suppressed by
$\sqrt{m_{\nu}/m_{n}}$ compared with neutrino interactions and for light
Portalinos of around an eV this ratio is $\sim 10^{-2}$. As with the heavy
Portalino case, these Portalino states will be mono-energetic, the energy of
each Portalino being equal to the DM mass, which in our model means they are
produced with TeV energies. These Portalinos will produce a line in the
neutrino spectrum, which provides an intriguing target for neutrino telescopes
with the location of the line providing a way to measure the DM mass. As
discussed above the current sensitivity of neutrino telescopes, such as
IceCube, is at least an order of magnitude too weak to detect neutrinos from
DM annihilations with thermal relic cross sections and with the potential
suppression of Portalinos compared with SM neutrinos this signal will be
challenging even with the development of the next generation of neutrino
telescope. However, this does provide an exciting target to focus on given the
potential for measuring both the DM mass and its annihilation cross section.

##  6 Conclusions and Discussion

It seems fairly natural that if a singlet right-handed neutrino does exist (as
is the case in many models of neutrino mass) that it may have interactions
with the hidden sector and may play a role in DM dynamics, creating a
connection or portal between the SM and dark sectors.

In this work we have expanded on the Portalino model outlined in
Schmaltz:2017oov  to include neutrino masses and an expanded dark sector
including a DM state. A dark $U(1)$ gauge symmetry was introduced and
spontaneously broken generating the dark sector masses including for the Dirac
Fermion DM state whose stability is ensured via an unbroken $Z_{3}$ symmetry.
The observed DM abundance in this model can be generated by freeze-out via DM
self-annihilations to either the Portalino states or the now massive dark
sector gauge bosons. Following the freeze-out of DM, a population of Portalino
states is produced. The Portalino lifetime is cosmologically relevant and as a
result provides the main constraints on the properties of these states and the
viability of the model. Portalinos can also potentially play an important role
in the indirect detection of DM in this model. DM can annihilate into high
energy mono-energetic Portalinos which in the case where they travel to the
Earth can then be search for in neutrino telescopes or, if they decay, their
SM decay products will lead to signals in gamma rays. Although current
indirect limits do not rule these out, it is possible future experiments will
née sensitive to these signals.

We have considered three qualitatively different scenarios categorised in
terms of the mass of the Portalino: an intermediate case, which is not
cosmologically viable; a heavy and light case, the former with allowed
parameter space for a Portalino mass $m_{n}\gtrsim$489\text{\,}\mathrm{MeV}$$
and the latter viable for $m_{n}\lesssim$6.4\text{\,}\mathrm{eV}$$ provided
Portalinos decouple before top quarks. In the heavy Portalino case there is an
upper limit of $m_{n}\lesssim$TeV due to the restriction that $m_{n}\ll
v_{\phi}$, where the requirement on the successful generation of the observed
DM abundance limits the maximum size of $v_{\phi}$ to the multi TeV mass scale
at most.

Throughout this work we have only considered including one Portalino (with
three heavy singlet neutrinos), but we could consider models with multiple
Portalinos (and/or a different number of heavy neutrinos). We can put concrete
restrictions on which configurations are viable by imposing that they must
give rise to at most one massless active neutrino. In the absence of specific
flavour symmetries, the number of massless states is given by $n_{0}={\rm
min}(0,3-n_{R}+n_{p})$ where $n_{R}$ is the number of heavy neutrinos and
$n_{p}$ is the number of Portalino states. Hence for a model with $n_{p}$
Portalinos the number of heavy right-handed neutrinos required is at least
$2+n_{p}$.

The inclusion of a full neutrino flavour model was beyond the scope of this
work. The details of such a flavour model will feed into the Portalino
phenomenology in a more complicated way compared with the parameterisation
used here in terms of $\eta$ and may lead to ways to widen the allowed
parameter range found in this work.

Another variation on what has been presented here is to remove the Majorana
mass term $M_{R}$. In order to generate the light active neutrino masses the
neutrino Yukawa couplings $\lambda_{\nu}$ would need to be small ($\sim
10^{-13}$). The structure of this model is significantly different: for
example the Portalino mass isn’t suppressed relative to the other dark states
and the Portalino-neutrino mixing is significantly suppressed. Interestingly,
if the $h-\phi$ interactions were turned off or also significantly suppressed,
the Portalino-neutrino mixing could be the strongest interaction between the
dark and visible sectors, and could potentially lead to the freeze-in
production McDonald:2001vt ; Hall:2009bx  of the Portalino or other dark
sector states.

The Portalino can in principle provide explanations for some anomalies.
Firstly, a decaying sterile state (which could be the heavy Portalino) has
been proposed as a solution to short baseline anomalies Aguilar-
Arevalo:2018gpe ; Aguilar:2001ty ; Fischer:2019fbw . An $\mathrm{eV}$ sterile
neutrino has also been mooted as a solution to these anomalies, however the
most straightforward case of an $\mathrm{eV}$ Portalino with sufficiently
strong mixing ($V_{n\nu}\sim 0.1$) with SM neutrinos would be ruled out by
cosmological constraints Aguilar-Arevalo:2018gpe ; Hamann:2011ge .

##  Appendix A Neutrino masses and mixing

In this Appendix, a detailed presentation of the diagonalisation of the
$(\nu_{L\alpha},\psi,\nu_{R\alpha}$) sector is outlined. The objective is to
evaluate the mass eigenstates and mixings of this seven by seven system with
mass matrix given by

$M=\left(\begin{array}[]{cc}\mbox{\Large 0}&\begin{matrix}M_{d}^{T}\\\
\mathbf{M}_{\psi}^{T}\end{matrix}\\\
\begin{matrix}M_{d}&\mathbf{M}_{\psi}\end{matrix}&M_{R}\end{array}\right),$

where $M_{d}$ and $M_{R}$ are $3\times 3$ matrices and $\mathbf{M}_{\psi}$ is
a three vector.

First it is noted that the mass matrix, $M$, has a zero eigenvalue,
$m_{\nu_{1}}$, with eigenvector:

$\mathbf{e}_{1}=N\begin{pmatrix}-M_{d}^{-1}\mathbf{M}_{\psi}\\\ 1\\\ 0\\\ 0\\\
0\end{pmatrix}$

where $N=1/\sqrt{1+\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}}$. We define an
orthonormal basis which includes this zero eigenvector:

$\mathbf{e}_{1},\;\mathbf{e}_{2}=\begin{pmatrix}\mathbf{x}_{1}\\\ 0\\\ 0\\\
0\\\ 0\end{pmatrix},\mathbf{e}_{3}=\begin{pmatrix}\mathbf{x}_{2}\\\ 0\\\ 0\\\
0\\\
0\end{pmatrix},\mathbf{e}_{4}=\frac{N}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\begin{pmatrix}M_{d}^{-1}\mathbf{M}_{\psi}\\\
\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}\\\ 0\\\ 0\\\
0\end{pmatrix},(\mathbf{e}_{i})_{j}=\delta_{ij},$

where $\mathbf{x}_{1,2}$ are chosen such that
$\mathbf{x}_{1,2}^{T}M_{d}^{-1}\mathbf{M}_{\psi}=\mathbf{x}_{1}^{T}\mathbf{x}_{2}=0$
and $\|\mathbf{x}_{1,2}\|^{2}=1$.

After rotating away the zero eigenstate, we obtain a 7$\times$7 matrix with a
non-zero 6$\times$6 sub-matrix with a seesaw-type structure:

$P^{-1}MP=\left(\begin{array}[]{cc}\mbox{\Large
0}_{4,4}&\begin{matrix}0_{1,3}\\\ M_{D}\end{matrix}\\\
\begin{matrix}0_{3,1}&M_{D}^{T}\end{matrix}&M_{R}\end{array}\right),$

where for the sake of clarity we have indicated the dimensions of the zero
matrices (e.g. $0_{n,m}$ is a $n\times m$ zero matix) and where

$M_{D}=\begin{pmatrix}\mathbf{x}_{1}^{T}M_{d}^{T}\\\
\mathbf{x}_{2}^{T}M_{d}^{T}\\\
\frac{\mathbf{M}_{\psi}^{T}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\end{pmatrix}$

and

$P=\left(\begin{array}[]{cc}\begin{matrix}-NM_{d}^{-1}\mathbf{M}_{\psi}&\mathbf{x}_{1}&\mathbf{x}_{2}&\frac{NM_{d}^{-1}\mathbf{M}_{\psi}}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\\\\[5.69054pt]
N&0&0&N\|M_{d}^{-1}\mathbf{M}_{\psi}\|\end{matrix}&\mbox{\Large
0}_{4,3}\\\\[5.69054pt] \;\;\;\;\mbox{\Large
0}_{3,4}&I_{3}\end{array}\right).$

.

The resulting matrix can be approximately block diagonalised, assuming the
hierarchy of masses $m_{d},m_{\psi}\ll m_{R}$:

$\displaystyle Q^{-1}P^{-1}MPQ=\begin{pmatrix}0&0_{1,3}&0_{1,3}\\\\[2.84526pt]
0_{3,1}&\;\;\;-M_{D}M_{R}^{-1}M_{D}^{T}+\mathcal{O}(M_{D}^{3}M_{R}^{-2})&0_{3,3}\\\\[2.84526pt]
0_{3,1}&0_{3,3}&\;\;M_{R}+\mathcal{O}(M_{D})\end{pmatrix}$ (41)

where

$Q=\begin{pmatrix}1&0_{1,3}&0_{1,3}\\\\[5.69054pt]
0_{3,1}&\;\;I_{3}-\frac{1}{2}M_{D}M_{R}^{-1}M_{R}^{-1}M_{D}^{T}+\mathcal{O}(m_{D}^{3}M_{R}^{-3})&M_{D}M_{R}^{-1}+\mathcal{O}(M_{D}^{3}M_{R}^{-3})\\\\[5.69054pt]
0_{3,1}&-M_{R}^{-1}M_{D}^{T}+\mathcal{O}(M_{D}^{3}M_{R}^{-3})&\;\;I_{3}-\frac{1}{2}M_{R}^{-1}M_{D}^{T}M_{D}M_{R}^{-1}+\mathcal{O}(M_{D}^{3}M_{R}^{-3})\end{pmatrix}.$

We note that the eigenvalues of the $M_{R}$ mass matrix will correspond to the
masses of the three heavy neutrino states, labelled $N_{i}$ in Section 2.

The remaining three mass eigenvalues contained within the central $3\times 3$
mass matrix block in Equation 41 are identified as the remaining two light
neutrino masses, along with the Portalino mass, $m_{n}$ in Section 2. The
explicit form of this mass matrix is given by

$-M_{D}M_{R}^{-1}M_{D}^{T}=-\begin{pmatrix}\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{1}&\;\;\;\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}&\;\;\frac{\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}\mathbf{M}_{\psi}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\\\\[5.69054pt]
\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{1}&\;\;\;\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}&\;\;\frac{\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}\mathbf{M}_{\psi}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\\\\[5.69054pt]
\frac{\mathbf{M}_{\psi}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{1}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}&\frac{\mathbf{M}_{\psi}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}&\;\;\frac{\mathbf{M}_{\psi}^{T}M_{R}^{-1}\mathbf{M}_{\psi}}{N^{2}\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}}\end{pmatrix}$

.

We use the remaining freedom to choose $\mathbf{x}_{1}$ (or equivalently
$\mathbf{x}_{2}$) to aid in further diagonalising. For example, choose
$\mathbf{x}_{1}\propto\left(M_{d}^{T}M_{R}^{-1}\mathbf{M}_{\psi}\right)\times\left(M_{d}^{-1}\mathbf{M}_{\psi}\right)$
(note that if this is zero then $M_{d}^{T}M_{R}^{-1}\mathbf{M}_{\psi}\propto
M_{d}^{-1}\mathbf{M}_{\psi}$ and hence we can choose $\mathbf{x}_{1}$ and
$\mathbf{x}_{2}$ such that the 1, 3 part is already block diagonalised - so
assume this isn’t the case), then:

$\displaystyle-M_{D}M_{R}^{-1}M_{D}^{T}=-\begin{pmatrix}c&d&0\\\ d&e&b\\\
0&b&a\end{pmatrix}$ (42)

where

$\displaystyle a$ $\displaystyle=$
$\displaystyle\frac{\mathbf{M}_{\psi}^{T}M_{R}^{-1}\mathbf{M}_{\psi}}{N^{2}\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}}=\mathcal{O}\left(\frac{m_{\psi}^{2}}{m_{R}}\right)$
$\displaystyle b$ $\displaystyle=$
$\displaystyle\frac{\mathbf{M}_{\psi}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}=\mathcal{O}\left(\frac{m_{d}m_{\psi}}{m_{R}}\right)$
$\displaystyle c$ $\displaystyle=$
$\displaystyle\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{1}=\mathcal{O}\left(\frac{m_{d}^{2}}{m_{R}}\right)$
$\displaystyle d$ $\displaystyle=$
$\displaystyle\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}=\mathcal{O}\left(\frac{m_{d}^{2}}{m_{R}}\right)$
$\displaystyle e$ $\displaystyle=$
$\displaystyle\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}M_{d}\mathbf{x}_{2}=\mathcal{O}\left(\frac{m_{d}^{2}}{m_{R}}\right),$

where we have used the definitions in Equation 10 to write the leading order
behaviour of these expressions assuming the hierarchy of masses $m_{d}\ll
m_{\psi}\ll m_{R}$. Utilising this hierarchy further we can apply a rotation,
$R_{24}$, to the mass matrix in Equation 42 such that

$-R_{24}^{-1}M_{D}M_{R}^{-1}M_{D}^{T}R_{24}=-\begin{pmatrix}c&d&0\\\
d&\;\;e-\frac{|b|^{2}}{a}&0\\\
0&0&a+\frac{|b|^{2}}{a}\end{pmatrix}+\mathcal{O}\left(\frac{m_{d}^{3}}{m_{R}m_{\psi}}\right),$

where

$R_{24}=\begin{pmatrix}1&0&0\\\ 0&\cos\theta_{24}&\sin\theta_{24}\\\
0&\sin\theta_{24}&\cos\theta_{24}\end{pmatrix},\;\;\theta_{24}=-\frac{|b|}{a}+\mathcal{O}\left(\frac{m^{3}_{d}}{m^{3}_{\psi}}\right).$

This leaves a final 2$\times$2 matrix to diagonalise. All elements are of the
same order, and $\mathbf{x}_{2}$ is already fixed by the orthogonality
constraints. A final rotation leaves the system diagonal:

$-R_{23}^{-1}R_{24}^{-1}M_{D}M_{R}^{-1}M_{D}^{T}R_{24}R_{23}=-\begin{pmatrix}m_{\nu_{2}}&0&0\\\
0&m_{\nu_{3}}&0\\\
0&0&m_{n}\end{pmatrix}+\mathcal{O}\left(\frac{m_{d}^{3}}{m_{R}m_{\psi}}\right)$

where

$R_{23}=\begin{pmatrix}\cos\theta_{23}&\sin\theta_{23}&0\\\
-\sin\theta_{23}&\cos\theta_{23}&0\\\ 0&0&1\end{pmatrix},$
$\displaystyle\cos\theta_{23}$ $\displaystyle=$
$\displaystyle\frac{\text{sign}(d)}{\sqrt{2}}\sqrt{1-\frac{\left(c+\frac{|b|^{2}}{a}-e\right)}{\sqrt{\left(c+\frac{|b|^{2}}{a}-e\right)^{2}+4|d|^{2}}}},$
$\displaystyle\sin\theta_{23}$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\sqrt{1+\frac{\left(c+\frac{|b|^{2}}{a}-e\right)}{\sqrt{\left(c+\frac{|b|^{2}}{a}-e\right)^{2}+4|d|^{2}}}}.$

The three masses, $m_{\nu_{1,2}}$ and $m_{n}$ read

$\displaystyle m_{\nu_{1}}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left(c+e-\frac{|b|^{2}}{a}-\sqrt{\left(c+\frac{|b|^{2}}{a}-e\right)^{2}+4|d|^{2}}\right)=\mathcal{O}\left(\frac{m_{d}^{2}}{m_{R}}\right)$
$\displaystyle m_{\nu_{2}}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\left(c+e-\frac{|b|^{2}}{a}+\sqrt{\left(c+\frac{|b|^{2}}{a}-e\right)^{2}+4|d|^{2}}\right)=\mathcal{O}\left(\frac{m_{d}^{2}}{m_{R}}\right)$
$\displaystyle m_{n}$ $\displaystyle=$ $\displaystyle
a+\frac{|b|^{2}}{a}=\mathcal{O}\left(\frac{m_{\psi}^{2}}{m_{R}}\right).$

Summarising the above, the 7$\times 7$ unitary matrix that diagonalises the
mass matrix, $M$, is given to leading order by

$V=\begin{pmatrix}i&0&0\\\ 0&iI_{3}&0\\\
0&0&I_{3}\end{pmatrix}PQ\begin{pmatrix}1&0&0\\\ 0&R_{24}&0\\\
0&0&I_{3}\end{pmatrix}\begin{pmatrix}1&0&0\\\ 0&R_{23}&0\\\
0&0&I_{3}\end{pmatrix}\begin{pmatrix}-i&0&0\\\ 0&-iI_{3}&0\\\
0&0&I_{3}\end{pmatrix}$

$=\left(\begin{array}[]{ccc}\begin{matrix}-NM_{d}^{-1}\mathbf{M}_{\psi}\\\\[8.53581pt]
N\end{matrix}&\begin{matrix}c_{23}\mathbf{x}_{1}-s_{23}c_{24}\mathbf{x}_{2}&s_{23}\mathbf{x}_{1}+c_{23}c_{24}\mathbf{x}_{2}&\frac{c_{24}NM_{d}^{-1}\mathbf{M}_{\psi}}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|}-s_{24}\mathbf{x}_{2}\\\\[5.69054pt]
-s_{23}s_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|&c_{23}s_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|&c_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|\end{matrix}&\begin{matrix}A\\\\[8.53581pt]
-i\mathbf{M}_{\psi}^{T}M_{R}^{-1}\end{matrix}\\\\[17.07164pt]
0_{3,1}&iM_{R}^{-1}M_{D}^{T}R_{24}R_{23}&I_{3}-\frac{1}{2}M_{R}^{-1}M_{D}^{T}M_{D}M_{R}^{-1}\end{array}\right),$

where

$A=-i\frac{\left(M_{d}^{-1}\mathbf{M}_{\psi}\right)\left(M_{R}^{-1}\mathbf{M}_{\psi}\right)^{T}}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}}-i\left(\mathbf{x}_{1}\right)\left(M_{R}^{-1}M_{d}\mathbf{x}_{1}\right)^{T}-i\left(\mathbf{x}_{2}\right)\left(M_{R}^{-1}M_{d}\mathbf{x}_{2}\right)^{T}$

and $c_{23}=\cos\theta_{23}$ etc. The order of the terms in the mixing matrix
$V$ are

$V\sim\begin{pmatrix}\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)\\\
\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)\\\
\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)\\\
\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{\psi}}\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)\\\
0&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{\psi}^{2}}{m^{2}_{R}}\right)&\mathcal{O}\left(\frac{m^{2}_{\psi}}{m^{2}_{R}}\right)\\\
0&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)&\mathcal{O}\left(\frac{m^{2}_{\psi}}{m^{2}_{R}}\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{\psi}}{m^{2}_{R}}\right)\\\
0&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{R}}\right)&\mathcal{O}\left(\frac{m_{\psi}}{m_{R}}\right)&\mathcal{O}\left(\frac{m^{2}_{\psi}}{m^{2}_{R}}\right)&\mathcal{O}\left(\frac{m^{2}_{\psi}}{m^{2}_{R}}\right)&\mathcal{O}\left(1\right)\end{pmatrix}.$

$V\sim\left(\begin{array}[]{ccc:c:ccc}&&&&&&\\\\[3.0pt]
&\left(U^{(\nu\nu_{l})}\right)^{T}&&\frac{m_{d}}{m_{\psi}}\mathbf{U}^{(n\nu_{l})}&&\frac{m_{d}}{m_{R}}\left(U^{(N\nu_{l})}\right)^{T}&\\\\[3.0pt]
&&&&&&\\\\[3.0pt]
\hdashline&\frac{m_{d}}{m_{\psi}}\left(\mathbf{U}^{(\nu\psi)}\right)^{T}&&U^{(n\psi)}&&\frac{m_{\psi}}{m_{R}}\left(\mathbf{U}^{(N\psi)}\right)^{T}&\\\\[3.0pt]
\hdashline&&&&&&\\\\[3.0pt]
&\frac{m_{d}}{m_{R}}\left(U^{(\nu\nu_{R})}\right)^{T}&&\frac{m_{\psi}}{m_{R}}\mathbf{U}^{(n\nu_{R})}&&\left(U^{(N\nu_{R})}\right)^{T}&\\\\[3.0pt]
&&&&&&\\\\[3.0pt] \end{array}\right).$

In Equation 2, the flavour mixing was parameterised by factors of $U$. These
nine factors have the following forms

$\displaystyle U^{(\nu\nu_{l})}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}-N\left(M_{d}^{-1}\mathbf{M}_{\psi}\right)^{T}\\\
c_{23}\mathbf{x}_{1}^{T}-s_{23}c_{24}\mathbf{x}_{2}^{T}\\\
s_{23}\mathbf{x}_{1}^{T}+c_{23}c_{24}\mathbf{x}_{2}^{T}\end{array}\right),\;\;\;\frac{m_{d}}{m_{\psi}}\mathbf{U}^{(\nu\psi)}=\left(\begin{array}[]{c}N\\\
-s_{23}s_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|\\\
c_{23}s_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|\end{array}\right),$ (49)
$\displaystyle\frac{m_{d}}{m_{R}}U^{(\nu\nu_{R})}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}0_{1,3}\\\
c_{23}\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}+s_{23}c_{24}\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}+s_{23}s_{24}\frac{M_{R}^{-1}\mathbf{M}_{\psi}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\\\
-s_{23}\mathbf{x}_{1}^{T}M_{d}^{T}M_{R}^{-1}+c_{23}c_{24}\mathbf{x}_{2}^{T}M_{d}^{T}M_{R}^{-1}+c_{23}s_{24}\frac{M_{R}^{-1}\mathbf{M}_{\psi}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|}\end{array}\right),$
(53) $\displaystyle\frac{m_{d}}{m_{\psi}}\mathbf{U}^{(n\nu_{l})}$
$\displaystyle=$
$\displaystyle\frac{c_{24}NM_{d}^{-1}\mathbf{M}_{\psi}}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|}-s_{24}\mathbf{x}_{2},\;\;\;\;\;U^{(n\psi)}=c_{24}N\|M_{d}^{-1}\mathbf{M}_{\psi}\|,$
$\displaystyle\frac{m_{\psi}}{m_{R}}\mathbf{U}^{(n\nu_{R})}$ $\displaystyle=$
$\displaystyle
s_{24}\frac{M_{R}^{-1}\mathbf{M}_{\psi}}{N\|M_{d}^{-1}\mathbf{M}_{\psi}\|},$
$\displaystyle\frac{m_{d}}{m_{R}}U^{(N\nu_{l})}$ $\displaystyle=$
$\displaystyle-i\frac{\left(M_{R}^{-1}\mathbf{M}_{\psi}\right)\left(M_{d}^{-1}\mathbf{M}_{\psi}\right)^{T}}{\|M_{d}^{-1}\mathbf{M}_{\psi}\|^{2}}-i\left(M_{R}^{-1}M_{d}\mathbf{x}_{1}\right)\left(\mathbf{x}_{1}\right)^{T}-i\left(M_{R}^{-1}M_{d}\mathbf{x}_{2}\right)\left(\mathbf{x}_{2}\right)^{T},$
$\displaystyle\frac{m_{\psi}}{m_{R}}\mathbf{U}^{(N\psi)}$ $\displaystyle=$
$\displaystyle-
iM_{R}^{-1}\mathbf{M}_{\psi},\;\;\;\;U^{(N\nu_{R})}=I_{3}-\frac{1}{2}M_{R}^{-1}M_{D}^{T}M_{D}M_{R}^{-1}$

.

The PMNS matrix, once again assuming there is no contribution from the charged
lepton sector, is then

$\begin{pmatrix}\nu_{e}\\\ \nu_{\mu}\\\ \nu_{\tau}\\\ \psi\end{pmatrix}=V_{\rm
PMNS}^{4\times 4}\begin{pmatrix}\nu_{1}\\\ \nu_{2}\\\ \nu_{3}\\\
n\end{pmatrix}+\begin{pmatrix}\mathcal{O}\left(\frac{m_{n}}{M_{R}}\right)\\\
\mathcal{O}\left(\frac{m_{d}}{M_{R}}\right)\\\
\mathcal{O}\left(\frac{m_{d}}{M_{R}}\right)\\\
\mathcal{O}\left(\frac{m_{d}}{M_{R}}\right)\end{pmatrix},$

where

$V_{\rm PMNS}^{4\times 4}=\\\
\begin{pmatrix}-NM_{d}^{-1}\mathbf{m}_{n}&\cos\theta_{23}\mathbf{x}_{1}-\sin\theta_{23}\cos\theta_{24}\mathbf{x}_{2}&\sin\theta_{23}\mathbf{x}_{1}+\cos\theta_{23}\cos\theta_{24}\mathbf{x}_{2}&\frac{\cos\theta_{24}NM_{d}^{-1}\mathbf{m}_{n}}{\|M_{d}^{-1}\mathbf{m}_{n}\|}-\sin\theta_{24}\mathbf{x}_{2}\\\
N&-\sin\theta_{23}\sin\theta_{24}N\|M_{d}^{-1}\mathbf{m}_{n}\|&\cos\theta_{23}\sin\theta_{24}N\|M_{d}^{-1}\mathbf{m}_{n}\|&\cos\theta_{24}N\|M_{d}^{-1}\mathbf{m}_{n}\|\end{pmatrix}$
(54)
$=\begin{pmatrix}\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)\\\
\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)\\\
\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(1\right)&\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)\\\
\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)&\mathcal{O}\left(\frac{m_{d}}{m_{n}}\right)&\mathcal{O}\left(1\right)\end{pmatrix}.$

##  Appendix B Full Lagrangian in mass eigenbasis

In this appendix we detail the dominant contributions to interactions in the
mass eigenbasis Lagrangian. For some components more than one term is included
if the dominant contribution depends on relative sizes of couplings. The full
Lagrangian in the mass eignestate basis reads

$\mathcal{L}=\mathcal{L}_{\rm matter-scalar}+\mathcal{L}_{\rm gauge-
matter}+\mathcal{L}_{h-\phi},$

where

$\displaystyle\mathcal{L}_{\rm matter-scalar}\supset$
$\displaystyle-\frac{m_{d}}{m_{R}}\cos\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta
j}^{(\nu\nu_{R})}}{2}\right)\overline{\nu}_{i}\nu_{j}h$
$\displaystyle-\left[\frac{m_{d}}{m_{R}}\cos\theta\left(\frac{U_{\alpha}^{(n\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta}^{(n\nu_{R})}}{2}\right)+\frac{m_{\psi}}{m_{R}}\sin\theta\left(\frac{U^{(n\psi)}\lambda^{\psi}_{\alpha}U_{\alpha}^{(n\nu_{R})}}{2}\right)\right]\overline{n}nh$
$\displaystyle-\bigg{[}\frac{m_{d}}{m_{R}}\cos\theta\left(\frac{U_{\alpha
i}^{(N\nu_{l})*}\lambda^{\nu*}_{\alpha\beta}U_{\beta
j}^{(N\nu_{R})*}}{2}\right)+\frac{m_{\psi}}{m_{R}}\sin\theta\left(\frac{U_{i}^{(N\psi)*}\lambda^{\psi*}_{\alpha}U_{\alpha
j}^{(N\nu_{R})*}}{2}\right)\bigg{]}\overline{N}_{i}N_{j}h$
$\displaystyle-\frac{m_{\psi}}{m_{R}}\cos\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta}^{(n\nu_{R})}}{2}\right)\overline{\nu_{i}}nh-\cos\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta
j}^{(N\nu_{R})}}{2}\right)\overline{\nu}_{i}N_{j}h$
$\displaystyle-\left[\frac{m_{d}}{m_{\psi}}\cos\theta\left(\frac{U_{\alpha}^{(n\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta
i}^{(N\nu_{R})}}{2}\right)+\sin\theta\left(\frac{U^{(n\psi)}\lambda^{\psi}_{\alpha}U_{\alpha
i}^{(N\nu_{R})}}{2}\right)\right]\overline{n}N_{i}h$
$\displaystyle-\bigg{[}\frac{m_{d}^{2}}{m_{R}m_{\psi}}\cos\theta\left(\frac{U_{i}^{(\nu\psi)}\lambda^{\psi}_{\alpha}U_{\alpha
j}^{(\nu\nu_{R})}}{2}\right)-\frac{m_{d}}{m_{R}}\sin\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta
j}^{(\nu\nu_{R})}}{2}\right)\bigg{]}\overline{\nu}_{i}\nu_{j}\phi$
$\displaystyle-\frac{m_{\psi}}{m_{R}}\cos\theta\left(\frac{U^{(n\psi)}\lambda^{\psi}_{\alpha}U_{\alpha}^{(n\nu_{R})}}{2}\right)\overline{n}n\phi-\frac{m_{\psi}}{m_{R}}\cos\theta\left(\frac{U_{i}^{(N\psi)*}\lambda^{\psi*}_{\alpha}U_{\alpha
j}^{(N\nu_{R})*}}{2}\right)\overline{N}_{i}N_{j}\phi$
$\displaystyle-\bigg{[}\frac{m_{d}}{m_{R}}\cos\theta\left(\frac{U_{i}^{(\nu\psi)}\lambda^{\psi}_{\alpha}U_{\alpha}^{(n\nu_{R})}}{2}\right)-\frac{m_{d}}{m_{R}}\cos\theta\left(\frac{U^{(n\psi)}\lambda^{\psi}_{\alpha}U_{\alpha
i}^{(\nu^{\prime}nu_{R})}}{2}\right)$
$\displaystyle\qquad\qquad\qquad-\frac{m_{\psi}}{m_{R}}\sin\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta}^{(n\nu_{R})}}{2}\right)\bigg{]}\overline{\nu_{i}}n\phi$
$\displaystyle-\left[\frac{m_{d}}{m_{\psi}}\cos\theta\left(\frac{U_{i}^{(\nu\psi)}\lambda^{\psi}_{\alpha}U_{\alpha
j}^{(N\nu_{R})}}{2}\right)-\sin\theta\left(\frac{U_{\alpha
i}^{(\nu\nu_{l})}\lambda^{\nu}_{\alpha\beta}U_{\beta
j}^{(N\nu_{R})}}{2}\right)\right]\overline{\nu}_{i}N_{j}\phi$
$\displaystyle-\cos\theta\left(\frac{U^{(n\psi)}\lambda^{\psi}_{\alpha}U_{\alpha
i}^{(N\nu_{R})}}{2}\right)\overline{n}N_{i}\phi+{\rm
h.c.}-\lambda_{X}\cos\theta\overline{X}X\phi+\lambda_{X}\sin\theta\overline{X}Xh,$
(55) $\displaystyle\mathcal{L}_{\rm gauge-matter}\supset$
$\displaystyle\frac{g}{\sqrt{2}}\bigg{[}\left(U_{\alpha
i}^{(\nu\nu_{l})}\overline{\nu}_{i}+\frac{m_{d}}{m_{\psi}}U_{\alpha}^{(n\nu_{l})}\overline{n}+\frac{m_{d}}{m_{R}}U_{\alpha
j}^{(N\nu_{l})}\overline{N}_{j}\right)\gamma^{\mu}\frac{1}{2}(1-\gamma^{5})e_{\alpha}W_{\mu}^{+}$
$\displaystyle\qquad\qquad\qquad+{\rm h.c.}\bigg{]}$
$\displaystyle-\frac{\sqrt{g^{2}+g^{\prime 2}}}{2}\bigg{[}\left(U_{\alpha
i}^{(\nu\nu_{l})}\overline{\nu}_{i}+\frac{m_{d}}{m_{\psi}}U_{\alpha}^{(n\nu_{l})}\overline{n}+\frac{m_{d}}{m_{R}}U_{\alpha
j}^{(N\nu_{l})}\overline{N}_{j}\right)\gamma^{\mu}\frac{1}{2}\gamma^{5}$
$\displaystyle\qquad\qquad\qquad\left(U_{\alpha
i}^{(\nu\nu_{l})*}\nu_{i}+\frac{m_{d}}{m_{\psi}}U_{\alpha}^{(n\nu_{l})*}n+\frac{m_{d}}{m_{R}}U_{\alpha
j}^{(N\nu_{l})*}N_{j}\right)Z_{\mu}\bigg{]}$
$\displaystyle+\frac{\tilde{g}}{2}\bigg{[}\left(\frac{m_{d}}{m_{\psi}}U_{i}^{(\nu\psi)}\overline{\nu}_{i}+U^{(n\psi)}\overline{n}+\frac{m_{\psi}}{m_{R}}U_{j}^{(N\psi)}\overline{N}_{j}\right)\gamma^{\mu}\frac{1}{2}\gamma^{5}$
$\displaystyle\qquad\qquad\qquad\left(\frac{m_{d}}{m_{\psi}}U_{i}^{(\nu\psi)*}\nu_{i}+U^{(n\psi)*}n+\frac{m_{\psi}}{m_{R}}U_{j}^{(N\psi)*}N_{j}\right)\omega_{\mu}\bigg{]}$
$\displaystyle-\frac{\tilde{g}}{4}\overline{X}\gamma^{\mu}\left(1+\gamma^{5}\right)X\omega_{\mu}$
$\displaystyle+\left(\frac{g^{2}}{2}W^{+}_{\mu}W^{-\mu}+\left(\frac{g^{2}+g^{\prime
2}}{4}\right)Z_{\mu}Z^{\mu}\right)$ $\displaystyle\qquad\times\left(\cos\theta
v_{h}h+\frac{\cos^{2}\theta}{2}h^{2}-\cos\theta\sin\theta h\phi-\sin\theta
v_{h}\phi+\frac{\sin^{2}\theta}{2}\phi^{2}\right)$
$\displaystyle+\frac{\tilde{g}^{2}}{4}\omega_{\mu}\omega^{\mu}\left(\cos\theta
v_{\phi}\phi+\frac{\cos^{2}\theta}{2}\phi^{2}+\cos\theta\sin\theta
h\phi+\sin\theta v_{\phi}h+\frac{\sin^{2}\theta}{2}h^{2}\right),$
$\displaystyle{\rm and}$ $\displaystyle\mathcal{L}_{h-\phi}\supset$
$\displaystyle\enspace\left(\lambda_{H}v_{h}\cos\theta-\frac{\lambda_{H,\Phi}v_{\phi}}{2}\sin\theta\right)\cos^{2}\theta
h^{3}+\frac{\lambda_{H}}{4}\cos^{4}\theta
h^{4}+\lambda_{\Phi}v_{\phi}\cos^{3}\theta\phi^{3}$
$\displaystyle+\frac{\lambda_{\Phi}}{4}\cos^{4}\theta\phi^{4}+\frac{\lambda_{H,\Phi}v_{\phi}}{2}\cos^{3}\theta
h^{2}\phi$
$\displaystyle+\left(\frac{\lambda_{H,\Phi}v_{h}}{2}\cos\theta-3\lambda_{\Phi}v_{\phi}\sin\theta+\lambda_{H,\Phi}v_{\phi}\sin\theta\right)\cos^{2}\theta
h\phi^{2}$
$\displaystyle+\left(\lambda_{H}-\frac{\lambda_{H,\Phi}}{2}\right)\cos^{3}\theta\sin\theta
h^{3}\phi+\frac{\lambda_{H,\Phi}}{4}\cos^{4}\theta h^{2}\phi^{2}$
$\displaystyle-\left(\lambda_{\Phi}-\frac{\lambda_{H,\Phi}}{2}\right)\cos^{3}\theta\sin\theta
h\phi^{3}.$

##  Appendix C Parameterisation of the mixing matrix for numerical evaluation.

Without a flavour model we have no guidance for what form the components (that
is, the factors of $U$ in Equation 2) of the $7\times 7$ unitary matrix, $V$,
will take. In order to evaluate the DM phenomenology we take the following
assignments

$\displaystyle U^{(\nu\nu_{l})}$
$\displaystyle=\sqrt{\left(1-\left(\frac{m_{d}}{m_{\psi}}\right)^{2}\right)\left(1-\left(\frac{m_{\psi}}{m_{R}}\right)^{2}\right)}N_{\rm
PMNS}^{\dagger},\;\;\;\mathbf{U}^{(\nu\psi)}=\frac{1}{\sqrt{3}}\left(\begin{array}[]{c}1\\\
1\\\ 1\end{array}\right)\eta,$ (59) $\displaystyle U^{(\nu\nu_{R})}$
$\displaystyle=\frac{1}{\sqrt{3}}\left(\begin{array}[]{ccc}1&1&1\\\ 1&1&1\\\
1&1&1\end{array}\right),\;\;\;\mathbf{U}^{(n\nu_{l})}=\frac{1}{\sqrt{3}}\left(\begin{array}[]{c}1\\\
1\\\ 1\end{array}\right)\eta,$ (66) $\displaystyle U^{(n\psi)}$
$\displaystyle=\sqrt{\left(1-\left(\frac{m_{d}}{m_{\psi}}\right)^{2}\right)\left(1-\left(\frac{m_{\psi}}{m_{R}}\right)^{2}\right)},\;\;\;\mathbf{U}^{(n\nu_{R})}=\frac{1}{\sqrt{3}}\left(\begin{array}[]{c}1\\\
1\\\ 1\end{array}\right),$ (70) $\displaystyle U^{(N\nu_{l})}$
$\displaystyle=\frac{1}{\sqrt{3}}\left(\begin{array}[]{ccc}1&1&1\\\ 1&1&1\\\
1&1&1\end{array}\right),\;\;\;\mathbf{U}^{(N\psi)}=\frac{1}{\sqrt{3}}\left(\begin{array}[]{c}1\\\
1\\\
1\end{array}\right),\;\;\;U^{(N\nu_{R})}=\sqrt{\left(1-\left(\frac{m_{d}}{m_{\psi}}\right)^{2}\right)\left(1-\left(\frac{m_{\psi}}{m_{R}}\right)^{2}\right)}I_{3},$
(77)

where $N_{\rm PMNS}$ is the experimentally observed PMNS matrix Gonzalez-
Garcia:2021dve .

##  Appendix D Dark matter annihilation cross section to hidden vectors

The annihilation rate for DM states into hidden sector gauge bosons reads

$\displaystyle\sigma(X\bar{X}\rightarrow\omega\omega)v=\frac{\tilde{g}^{4}\left(m_{X}^{2}-m_{\omega}^{2}\right)^{3/2}}{256\pi
m_{\omega}^{2}m_{X}\left(2m_{X}^{2}-m_{\omega}^{2}\right)}+\frac{v^{2}\tilde{g}^{2}\;F(m_{X},m_{\omega},\tilde{g},\lambda_{X},\theta)}{\left[\left(4m_{X}^{2}-m_{\phi}^{2}\right)^{2}+m_{\phi}^{2}\Gamma_{\phi}^{2}\right]}$
$\displaystyle\equiv\frac{\tilde{g}^{4}\left(m_{X}^{2}-m_{\omega}^{2}\right)^{3/2}}{256\pi
m_{\omega}^{2}m_{X}\left(2m_{X}^{2}-m_{\omega}^{2}\right)}+\frac{v^{2}\tilde{g}^{2}\sqrt{m_{X}^{2}-m_{\omega}^{2}}}{6144\pi
m_{X}\left[\left(4m_{X}^{2}-m_{\phi}^{2}\right)^{2}+m_{\phi}^{2}\Gamma_{\phi}^{2}\right]\left(m_{\omega}^{3}-2m_{\omega}m_{X}^{2}\right)^{4}}$
$\displaystyle\times\left[2K_{1}\lambda_{X}m_{\omega}\left(m_{\omega}^{2}-2m_{X}^{2}\right)^{2}+18K_{2}\lambda_{X}^{2}m_{\omega}^{2}\left(m_{\omega}^{2}-2m_{X}^{2}\right)^{4}+\tilde{g}^{2}K_{4}\left(\left(m_{\phi}^{2}-4m_{X}^{2}\right)^{2}+m_{\phi}^{2}\Gamma_{\phi}^{2}\right)\right.$
$\displaystyle\hskip
28.45274pt\left.+8K_{3}\;\tilde{g}\lambda_{X}m_{\omega}m_{X}\left(m_{\phi}^{2}-4m_{X}^{2}\right)\left(m_{\omega}^{2}-2m_{X}^{2}\right)^{2}\right],$

where

$\displaystyle K_{1}$ $\displaystyle=3K_{2}\lambda_{X}(4\cos{2\theta}+\cos
4\theta)\,m_{\omega}\left(m_{\omega}^{2}-2m_{X}^{2}\right)^{2}+4K_{3}\cos{2\theta}\,\tilde{g}\,m_{X}\left(m_{\phi}^{2}-4m_{X}^{2}\right)$
$\displaystyle K_{2}$
$\displaystyle=4m_{X}^{4}-4m_{\omega}^{2}m_{X}^{2}+3m_{\omega}^{4},$
$\displaystyle K_{3}$
$\displaystyle=10m_{\omega}^{4}m_{X}^{2}-24m_{\omega}^{2}m_{X}^{4}+m_{\omega}^{6}+16m_{X}^{6},$
$\displaystyle K_{4}$
$\displaystyle=3m_{\omega}^{8}m_{X}^{2}-84m_{\omega}^{6}m_{X}^{4}+152m_{\omega}^{4}m_{X}^{6}-80m_{\omega}^{2}m_{X}^{8}+7m_{\omega}^{10}+32m_{X}^{10}.$

## References

  * (1) R.M. Schabinger and J.D. Wells, _A Minimal spontaneously broken hidden sector and its impact on Higgs boson physics at the large hadron collider_ , _Phys. Rev. D_ 72 (2005) 093007 [hep-ph/0509209].
  * (2) B. Patt and F. Wilczek, _Higgs-field portal into hidden sectors_ , hep-ph/0605188.
  * (3) J. March-Russell, S.M. West, D. Cumberbatch and D. Hooper, _Heavy Dark Matter Through the Higgs Portal_ , _JHEP_ 07 (2008) 058 [0801.3440].
  * (4) B. Holdom, _Two U(1)’s and Epsilon Charge Shifts_ , _Phys. Lett. B_ 166 (1986) 196.
  * (5) K.R. Dienes, C.F. Kolda and J. March-Russell, _Kinetic mixing and the supersymmetric gauge hierarchy_ , _Nucl. Phys. B_ 492 (1997) 104 [hep-ph/9610479].
  * (6) F. Del Aguila, _The Physics of z-prime bosons_ , _Acta Phys. Polon. B_ 25 (1994) 1317 [hep-ph/9404323].
  * (7) K.S. Babu, C.F. Kolda and J. March-Russell, _Leptophobic U(1) $s$ and the R($b$) - R($c$) crisis_, _Phys. Rev. D_ 54 (1996) 4635 [hep-ph/9603212].
  * (8) M. Lindner, A. Merle and V. Niro, _Enhancing Dark Matter Annihilation into Neutrinos_ , _Phys. Rev. D_ 82 (2010) 123529 [1005.3116].
  * (9) A. Falkowski, J. Juknevich and J. Shelton, _Dark Matter Through the Neutrino Portal_ , 0908.1790.
  * (10) V. Gonzalez Macias and J. Wudka, _Effective theories for Dark Matter interactions and the neutrino portal paradigm_ , _JHEP_ 07 (2015) 161 [1506.03825].
  * (11) B. Batell, T. Han and B. Shams Es Haghi, _Indirect Detection of Neutrino Portal Dark Matter_ , _Phys. Rev. D_ 97 (2018) 095020 [1704.08708].
  * (12) B. Batell, T. Han, D. McKeen and B. Shams Es Haghi, _Thermal Dark Matter Through the Dirac Neutrino Portal_ , _Phys. Rev. D_ 97 (2018) 075016 [1709.07001].
  * (13) M. Blennow, E. Fernandez-Martinez, A. Olivares-Del Campo, S. Pascoli, S. Rosauro-Alcaraz and A. Titov, _Neutrino Portals to Dark Matter_ , _Eur. Phys. J. C_ 79 (2019) 555 [1903.00006].
  * (14) Y. Nomura and J. Thaler, _Dark Matter through the Axion Portal_ , _Phys. Rev. D_ 79 (2009) 075008 [0810.5397].
  * (15) J. March-Russell, H. Tillim and S.M. West, _Reproductive freeze-in of self-interacting dark matter_ , _Phys. Rev. D_ 102 (2020) 083018 [2007.14688].
  * (16) M. Schmaltz and N. Weiner, _A Portalino to the Dark Sector_ , _JHEP_ 02 (2019) 105 [1709.09164].
  * (17) A.G. Dias, C.A. de S. Pires, P.S. Rodrigues da Silva and A. Sampieri, _A Simple Realization of the Inverse Seesaw Mechanism_ , _Phys. Rev. D_ 86 (2012) 035007 [1206.2590].
  * (18) D. Liu and N. Weiner, _A Portalino to the Twin Sector_ , 1905.00861.
  * (19) T. Gherghetta, J. Kersten, K. Olive and M. Pospelov, _Evaluating the price of tiny kinetic mixing_ , _Phys. Rev. D_ 100 (2019) 095001 [1909.00696].
  * (20) O. Lebedev, _The Higgs portal to cosmology_ , _Prog. Part. Nucl. Phys._ 120 (2021) 103881 [2104.03342].
  * (21) A. Falkowski, C. Gross and O. Lebedev, _A second Higgs from the Higgs portal_ , _JHEP_ 05 (2015) 057 [1502.01361].
  * (22) S. Parke and M. Ross-Lonergan, _Unitarity and the three flavor neutrino mixing matrix_ , _Phys. Rev. D_ 93 (2016) 113009 [1508.05095].
  * (23) I. Esteban, M. Gonzalez-Garcia, M. Maltoni, T. Schwetz and A. Zhou, _The fate of hints: updated global analysis of three-flavor neutrino oscillations_ , _JHEP_ 09 (2020) 178 [2007.14792].
  * (24) G. Bélanger, F. Boudjema, A. Goudelis, A. Pukhov and B. Zaldivar, _micrOMEGAs5.0 : Freeze-in_ , _Comput. Phys. Commun._ 231 (2018) 173 [1801.03509].
  * (25) A. Kurylov and M. Kamionkowski, _Generalized analysis of weakly interacting massive particle searches_ , _Phys. Rev. D_ 69 (2004) 063503 [hep-ph/0307185].
  * (26) J. Ellis, N. Nagata and K.A. Olive, _Uncertainties in WIMP Dark Matter Scattering Revisited_ , _Eur. Phys. J. C_ 78 (2018) 569 [1805.09795].
  * (27) LUX-ZEPLIN collaboration, _First Dark Matter Search Results from the LUX-ZEPLIN (LZ) Experiment_ , 2207.03764.
  * (28) M. Kawasaki, K. Kohri, T. Moroi and Y. Takaesu, _Revisiting Big-Bang Nucleosynthesis Constraints on Long-Lived Decaying Particles_ , _Phys. Rev. D_ 97 (2018) 023502 [1709.01211].
  * (29) B.D. Fields, K.A. Olive, T.-H. Yeh and C. Young, _Big-Bang Nucleosynthesis After Planck_ , _JCAP_ 03 (2020) 010 [1912.01132].
  * (30) Planck collaboration, _Planck 2018 results. VI. Cosmological parameters_ , _Astron. Astrophys._ 641 (2020) A6 [1807.06209].
  * (31) F.F. Deppisch, P.S. Bhupal Dev and A. Pilaftsis, _Neutrinos and Collider Physics_ , _New J. Phys._ 17 (2015) 075019 [1502.06541].
  * (32) F. del Aguila, J. de Blas and M. Perez-Victoria, _Effects of new leptons in Electroweak Precision Data_ , _Phys. Rev. D_ 78 (2008) 013010 [0803.4008].
  * (33) E. Akhmedov, A. Kartavtsev, M. Lindner, L. Michaels and J. Smirnov, _Improving Electro-Weak Fits with TeV-scale Sterile Neutrinos_ , _JHEP_ 05 (2013) 081 [1302.1872].
  * (34) L. Basso, O. Fischer and J.J. van der Bij, _Precision tests of unitarity in leptonic mixing_ , _EPL_ 105 (2014) 11001 [1310.2057].
  * (35) J. de Blas, _Electroweak limits on physics beyond the Standard Model_ , _EPJ Web Conf._ 60 (2013) 19008 [1307.6173].
  * (36) S. Antusch and O. Fischer, _Testing sterile neutrino extensions of the Standard Model at future lepton colliders_ , _JHEP_ 05 (2015) 053 [1502.05915].
  * (37) L3 collaboration, _Search for isosinglet neutral heavy leptons in Z0 decays_ , _Phys. Lett. B_ 295 (1992) 371.
  * (38) DELPHI collaboration, _Search for neutral heavy leptons produced in Z decays_ , _Z. Phys. C_ 74 (1997) 57.
  * (39) ATLAS collaboration, _Search for heavy neutral leptons in decays of $W$ bosons produced in 13 TeV $pp$ collisions using prompt and displaced signatures with the ATLAS detector_, _JHEP_ 10 (2019) 265 [1905.09787].
  * (40) CMS collaboration, _Search for heavy Majorana neutrinos in same-sign dilepton channels in proton-proton collisions at $\sqrt{s}=13$ TeV_, _JHEP_ 01 (2019) 122 [1806.10905].
  * (41) CMS collaboration, _Search for heavy neutral leptons in events with three charged leptons in proton-proton collisions at $\sqrt{s}=$ 13 TeV_, _Phys. Rev. Lett._ 120 (2018) 221801 [1802.02965].
  * (42) CMS collaboration, _Search for long-lived heavy neutral leptons with displaced vertices in proton-proton collisions at $\sqrt{\mathrm{s}}$ =13 TeV_, _JHEP_ 07 (2022) 081 [2201.05578].
  * (43) D. Curtin et al., _Long-Lived Particles at the Energy Frontier: The MATHUSLA Physics Case_ , _Rept. Prog. Phys._ 82 (2019) 116201 [1806.07396].
  * (44) S. Antusch, E. Cazzato, M. Drewes, O. Fischer, B. Garbrecht, D. Gueter et al., _Probing Leptogenesis at Future Colliders_ , _JHEP_ 09 (2018) 124 [1710.03744].
  * (45) NOMAD collaboration, _Search for heavy neutrinos mixing with tau neutrinos_ , _Phys. Lett. B_ 506 (2001) 27 [hep-ex/0101041].
  * (46) O. Ruchayskiy and A. Ivashko, _Experimental bounds on sterile neutrino mixing angles_ , _JHEP_ 06 (2012) 100 [1112.3319].
  * (47) A. Atre, T. Han, S. Pascoli and B. Zhang, _The Search for Heavy Majorana Neutrinos_ , _JHEP_ 05 (2009) 030 [0901.3589].
  * (48) G.L. Castro and N. Quintero, _Bounding resonant Majorana neutrinos from four-body B and D decays_ , _Phys. Rev. D_ 87 (2013) 077901 [1302.1504].
  * (49) H. Yuan, T. Wang, G.-L. Wang, W.-L. Ju and J.-M. Zhang, _Lepton-number violating four-body decays of heavy mesons_ , _JHEP_ 08 (2013) 066 [1304.3810].
  * (50) Y. Wang, S.-S. Bao, Z.-H. Li, N. Zhu and Z.-G. Si, _Study Majorana Neutrino Contribution to B-meson Semi-leptonic Rare Decays_ , _Phys. Lett. B_ 736 (2014) 428 [1407.2468].
  * (51) LHCb collaboration, _Search for Majorana neutrinos in $B^{-}\to\pi^{+}\mu^{-}\mu^{-}$ decays_, _Phys. Rev. Lett._ 112 (2014) 131802 [1401.5361].
  * (52) LBNE collaboration, _The Long-Baseline Neutrino Experiment: Exploring Fundamental Symmetries of the Universe_ , in _Snowmass 2013: Workshop on Energy Frontier_ , 7, 2013 [1307.7335].
  * (53) SHiP collaboration, _Sensitivity of the SHiP experiment to Heavy Neutral Leptons_ , _JHEP_ 04 (2019) 077 [1811.00930].
  * (54) SHiP collaboration, _A facility to Search for Hidden Particles (SHiP) at the CERN SPS_ , 1504.04956.
  * (55) D.A. Bryman and R. Shrock, _Improved Constraints on Sterile Neutrinos in the MeV to GeV Mass Range_ , _Phys. Rev. D_ 100 (2019) 053006 [1904.06787].
  * (56) Belle collaboration, _Search for heavy neutrinos at Belle_ , _Phys. Rev. D_ 87 (2013) 071102 [1301.1105].
  * (57) P.F. de Salas and S. Pastor, _Relic neutrino decoupling with flavour oscillations revisited_ , _JCAP_ 07 (2016) 051 [1606.06986].
  * (58) V. Shtabovenko, R. Mertig and F. Orellana, _FeynCalc 9.3: New features and improvements_ , _Comput. Phys. Commun._ 256 (2020) 107478 [2001.04407].
  * (59) A. de Gouvêa and A. Kobach, _Global Constraints on a Heavy Neutrino_ , _Phys. Rev. D_ 93 (2016) 033005 [1511.00683].
  * (60) Fermi-LAT collaboration, _Searching for Dark Matter Annihilation from Milky Way Dwarf Spheroidal Galaxies with Six Years of Fermi Large Area Telescope Data_ , _Phys. Rev. Lett._ 115 (2015) 231301 [1503.02641].
  * (61) MAGIC, Fermi-LAT collaboration, _Limits to Dark Matter Annihilation Cross-Section from a Combined Analysis of MAGIC and Fermi-LAT Observations of Dwarf Satellite Galaxies_ , _JCAP_ 02 (2016) 039 [1601.06590].
  * (62) CTA collaboration, _Sensitivity of the Cherenkov Telescope Array for probing cosmology and fundamental physics with gamma-ray propagation_ , _JCAP_ 02 (2021) 048 [2010.01349].
  * (63) IceCube collaboration, _Search for Dark Matter Annihilation in the Galactic Center with IceCube-79_ , _Eur. Phys. J. C_ 75 (2015) 492 [1505.07259].
  * (64) ANTARES collaboration, _Search for secluded dark matter towards the Galactic Centre with the ANTARES neutrino telescope_ , _JCAP_ 06 (2022) 028 [2203.06029].
  * (65) R. Abbasi et al., _Indirect search for dark matter in the Galactic Centre with IceCube_ , _PoS_ ICRC2021 (2021) 524.
  * (66) J. McDonald, _Thermally generated gauge singlet scalars as selfinteracting dark matter_ , _Phys. Rev. Lett._ 88 (2002) 091304 [hep-ph/0106249].
  * (67) L.J. Hall, K. Jedamzik, J. March-Russell and S.M. West, _Freeze-In Production of FIMP Dark Matter_ , _JHEP_ 03 (2010) 080 [0911.1120].
  * (68) MiniBooNE collaboration, _Significant Excess of ElectronLike Events in the MiniBooNE Short-Baseline Neutrino Experiment_ , _Phys. Rev. Lett._ 121 (2018) 221801 [1805.12028].
  * (69) LSND collaboration, _Evidence for neutrino oscillations from the observation of $\bar{\nu}_{e}$ appearance in a $\bar{\nu}_{\mu}$ beam_, _Phys. Rev. D_ 64 (2001) 112007 [hep-ex/0104049].
  * (70) O. Fischer, A. Hernández-Cabezudo and T. Schwetz, _Explaining the MiniBooNE excess by a decaying sterile neutrino with mass in the 250 MeV range_ , _Phys. Rev. D_ 101 (2020) 075045 [1909.09561].
  * (71) J. Hamann, S. Hannestad, G.G. Raffelt and Y.Y.Y. Wong, _Sterile neutrinos with eV masses in cosmology: How disfavoured exactly?_ , _JCAP_ 09 (2011) 034 [1108.4136].
  * (72) M.C. Gonzalez-Garcia, M. Maltoni and T. Schwetz, _NuFIT: Three-Flavour Global Analyses of Neutrino Oscillation Experiments_ , _Universe_ 7 (2021) 459 [2111.03086].