diagram

Lepton universality in a model with three generations of sterile Majorana
neutrinos

M.N. Dubinina, D.M. Kazarkinb111corresponding author, e-mail:
<EMAIL_ADDRESS>
aSkobeltsyn Institute of Nuclear Physics, Lomonosov Moscow State University
119991, Moscow, Russia

bPhysics Department, Lomonosov Moscow State University
119991, Moscow, Russia

Abstract

> The extension of the Standard Model lepton sector by three right-handed
> Majorana neutrinos with masses up to 102 MeV scale is considered. While the
> lightest heavy neutral lepton is the dark matter particle with mass of the
> order of 10 keV, the remaining two leptons ensure standard (active) neutrino
> mass generation by means of the see-saw type I mechanism. Besides lepton
> flavor violation these heavy neutral leptons can violate lepton universality
> in the charged meson decays. Deviations from universality in the meson
> decays may be experimentally observable. Their magnitude within two natural
> mixing scenarios taking into account the three-particle decays in
> calculation of the meson total width is analysed. The strong dependence of
> the violation of lepton universality on the mixing scenario is discussed.

##  1 Introduction. The main features of the model.

Extension of the Standard Model (SM) leptonic sector by heavy neutral leptons
of right helicity (HNL, also referred to as sterile Majorana neutrinos) known
for a long time [1, 2] has been analyzed multilaterally recently due to the
attractive general features of such an extension, within the framework of
which the symmetry between right and left neutrinos is restored, new large
energy scales are not necessarily introduced, neutrino oscillations and their
masses generation by means of the see-saw mechanism are successfully explained
and a number of important cosmological applications of the model, such as
baryon asymmetry of the Universe generation, description of the inflationary
stage of the early Universe and its accelerated expansion at the present time
are successfully interpreted [3, 4]. New useful features are realized in the
construction of the so-called minimal neutrino standard model $\nu$MSM [5, 6],
where the HNL masses do not exceed the electroweak scale. Cosmological
observations impose significant limitations on the model parameter space,
which lead to at least three HNL and establish a strict upper limit on the
mass of the lightest active neutrino $min(\nu_{e,\mu,\tau})$. The lightest
heavy lepton $N_{1}$ with mass of the order of 10 keV, the lifetime is more
than $\tau_{Universe}=10^{17}$ sec and mixing parameter is of the order of
10-13 – 10-7, plays a role of the dark matter (DM) particle in such an
extension [7, 8]. The direct method of $N_{1}$ DM detection is due to the
possibility of observation of the one-loop decay process of
$N_{1}\to\gamma\nu$ in galactic media [9]. The two remaining heavy leptons
ensure the mechanism of mass generation of active neutrinos, their masses can
vary in a wide range of values up to multiGeV scale, however, enough baryon
asymmetry can be generated even if they are of the order of MeV and their mass
splitting is rather small [10].

Mixing of light enough HNL states with standard (or active) neutrino states
could lead to observable HNL production in charged meson decays such as
$\pi^{+}\to e^{+}N_{2,3}$ and $K^{+}\to l^{+},N_{2,3}$, $l=e,\mu$ which could
violate lepton universality principle demonstrating departures from the SM
ratio $R_{\pi,K}=\Gamma(\pi,K\to e\nu)/\Gamma(\pi,K\to\mu\nu)$ [11, 12, 14,
15, 16] which is a quantity stable with respect to radiative corrections and
hadronization uncertainties.

In this paper we estimate the possible departure of the lepton universality
parameter

$\Delta r_{M}=\frac{R_{M}}{R^{SM}_{M}}-1$ (1)

from zero value due to HNL contributions in the $\nu MSM$-like model where an
explicit form of mixing for the three lepton generations is used. The
Lagrangian of extension has the form

$\mathcal{L}=\mathcal{L}_{SM}+i\overline{\nu}_{R}\partial_{\mu}\gamma^{\mu}\nu_{R}-\left(F\;\overline{l}_{L}\nu_{R}\tilde{H}+\frac{M_{M}}{2}\overline{\nu^{c}}_{R}\nu_{R}+h.c\right),$
(2)

where $l_{L}=(\nu_{L},\,e_{L})^{T}$ is the left lepton doublet, $\nu_{R}$ are
HNL flavor states, $(\nu_{R})^{c}\equiv C\overline{\nu}_{R}^{T}$
($C=i\gamma_{2}\gamma_{0}$),
$\overline{\nu}_{R}\equiv\nu_{R}^{\dagger}\gamma^{0}$, $H$ is the Higgs
doublet ($\tilde{H}=\epsilon_{ij}H^{\dagger}$), $F$ is the Yukawa matrix and
$M_{M}$ is a Majorana mass matrix. After spontaneous symmetry breaking
$M_{D}=F\langle H\rangle=Fv$ ($v=174$ GeV) is the matrix of Yukawa term. The
full 6$\times$6 mass matrix

$\displaystyle\frac{1}{2}(\overline{\nu}_{L}\vspace{2mm}\overline{\nu^{c}}_{R})\mathcal{M}\left(\begin{array}[]{c}\nu_{L}^{c}\\\
\nu_{R}\end{array}\right)+h.c.\,=\,\frac{1}{2}(\overline{\nu}_{L}\vspace{2mm}\overline{\nu^{c}}_{R})\left(\begin{array}[]{cc}0&M_{D}\\\
M_{D}^{T}&M_{M}\end{array}\right)\left(\begin{array}[]{c}\nu_{L}^{c}\\\
\nu_{R}\end{array}\right)+h.c.,$ (9)

where the flavor states $(\nu_{L})_{\alpha},(\nu_{R})_{I}$ and the mass states
$\nu_{k},N_{I}$ ($\alpha=e,\mu,\tau$, $k,I=1,2,3$) are connected by the
unitary transformation

$\displaystyle\left(\begin{array}[]{c}\nu_{L}\\\
\nu_{R}^{c}\end{array}\right)=\mathcal{U}P_{L}\left(\begin{array}[]{c}\upnu\\\
N\end{array}\right),\hskip 22.76219pt\mathcal{U}=\mathcal{W}\cdot{\rm
diag}(U_{\nu},U_{N}^{*})$ (14)

where $P_{L}$ is the left projector and $U_{\nu},U_{N}$ are unitary $3\times
3$ matrices. Block-diagonal form of the mass matrix (9) looks as

$\displaystyle\mathcal{U}^{\dagger}\mathcal{M}\mathcal{U}^{*}=\left(\begin{array}[]{cc}U_{\nu}^{\dagger}&0\\\
0&U_{N}^{T}\end{array}\right)\mathcal{W}^{\dagger}\mathcal{M}\mathcal{W}^{*}\left(\begin{array}[]{cc}U_{\nu}^{*}&0\\\
0&U_{N}\end{array}\right)=\left(\begin{array}[]{cc}U_{\nu}^{\dagger}m_{\nu}U_{\nu}^{*}&0\\\
0&U_{N}^{T}M_{N}U_{N}\end{array}\right)=\left(\begin{array}[]{cc}\hat{m}&0\\\
0&\hat{M}\end{array}\right),$ (23)

where $\hat{m}={\rm diag}(m_{1},m_{2},m_{3})$, $\hat{M}={\rm
diag}(M_{1},M_{2},M_{3})$,
$\mathcal{W}^{\dagger}\mathcal{M}\mathcal{W}=diag(m_{\nu},\;M_{N})$,
$\mathcal{M}$ is defined by Eq.(9). In the following diagonalization procedure
[17] the unitary $\mathcal{W}$-matrix is represented as an exponent of an
antihermitian matrix

$\mathcal{W}=\exp\left(\begin{array}[]{cc}0&\theta\\\
-\theta^{\dagger}&0\end{array}\right)$ (24)

and decomposed

$\displaystyle\mathcal{W}=\left(\begin{array}[]{cc}1-\frac{1}{2}\theta\theta^{\dagger}+O(\theta^{4})&\theta+O(\theta^{3})\\\
-\theta^{\dagger}+O(\theta^{3})&1-\frac{1}{2}\theta^{\dagger}\theta+O(\theta^{4})\end{array}\right).$
(27)

The flavor states are related to the mass states in the following form

$\displaystyle\nu_{L}$ $\displaystyle\simeq$
$\displaystyle\left(1-\frac{1}{2}\theta\theta^{\dagger}\right)U_{\nu}P_{L}\upnu+\theta
U_{N}^{*}P_{L}N,$ (28) $\displaystyle\nu^{c}_{R}$ $\displaystyle\simeq$
$\displaystyle-\theta^{\dagger}U_{\nu}P_{L}\upnu+\left(1-\frac{1}{2}\theta^{\dagger}\theta\right)U_{N}^{*}P_{L}N.$
(29)

For the left neutrino the main contribution in the flavor basis is given by
the first term in (28) which corresponds to the the well-known
phenomenological relation $\nu_{L\alpha}=\sum_{\alpha}(U_{\rm PMNS})_{\alpha
j}P_{L}\upnu_{j}$ [18]. Deviation from unitarity for the PMNS matrix in the
approximations $\mathcal{W}\sim\mathcal{O}(\theta^{2})$ (and also
$\mathcal{W}\sim\mathcal{O}(\theta^{3})$) is given by $U_{\rm
PMNS}\simeq\left(1-\frac{1}{2}\theta\theta^{\dagger}\right)U_{\nu}$ and
defined by the $\eta=-\frac{1}{2}\theta\theta^{\dagger}$-matrix. The
lagrangian terms for HNL currents interaction with $W^{\pm},Z$ bosons have the
form

$\displaystyle\mathcal{L}_{NC}^{\nu}$ $\displaystyle=$
$\displaystyle-\frac{g}{2c_{W}}\gamma^{\mu}\overline{\upnu}_{L}U_{\rm
PMNS}^{\dagger}U_{\rm PMNS}\upnu_{L}Z_{\mu},$ (30)
$\displaystyle\mathcal{L}_{CC}^{\nu}$ $\displaystyle=$
$\displaystyle-\frac{g}{\sqrt{2}}\overline{l}_{L}\gamma^{\mu}U_{\rm
PMNS}\upnu_{L}W_{\mu}^{-}+h.c.,$ $\displaystyle\mathcal{L}_{NC}^{N}$
$\displaystyle=$
$\displaystyle-\frac{g}{2c_{W}}\overline{N}_{L}\gamma^{\mu}U_{N}^{T}\theta^{\dagger}\theta
U_{N}^{*}N_{L}Z_{\mu}$ $\displaystyle-$
$\displaystyle\left[\frac{g}{2c_{W}}\overline{\upnu}_{L}\gamma^{\mu}U_{\nu}^{\dagger}\theta
U_{N}^{*}N_{L}Z_{\mu}+h.c.\right],$ $\displaystyle\mathcal{L}_{CC}^{N}$
$\displaystyle=$
$\displaystyle-\frac{g}{\sqrt{2}}\overline{l}_{L}\gamma^{\mu}\theta
U_{N}^{*}N_{L}W_{\mu}^{-}+h.c.$

In the following consideration we are keeping only the first and the second
order terms in $\theta$, as it is customary to do in the available literature.
Some details of the diagonalization in the third order in $\theta$ are
discussed in [19]. The HNL mixing parameter is defined in the approximation
$\mathcal{W}\sim\mathcal{O}(\theta^{2})$ as $\Theta\equiv\theta U_{N}^{*}$.
The standard set of active neutrino masses is defined in the framework of the
$\mathcal{O}(\theta^{2})$ scenario as a solution of see-saw type I equation

$m_{\nu}\simeq-M_{D}\theta^{T}\simeq-M_{D}M^{-1}_{M}M_{D}^{T}$ (31)

with ambiguous definition of $M_{D}$ by means of the $U_{N}$ mass matrix in
the HNL sector [17, 20]

$M_{D}=\pm iU_{\rm PMNS}\sqrt{\hat{m}}\Omega\sqrt{\hat{M}}U_{N}^{\dagger},$
(32)

where $\Omega$ is an arbitrary orthogonal matrix, $\Omega\Omega^{T}=$0\.
Discussion of ambiguity of the choice of the type of matrix $\Omega$ can be
found in [20], the most significant of them are related to the processes of
lepton flavor violation [21]. In the following we use the
$\mathcal{W}\sim\mathcal{O}(\theta^{2})$ analysis of [22], where the ”minimal
parametric mixing” $\Omega=$I is discussed. The mixing matrix has the form

$\displaystyle\Theta_{min}=\theta
U^{*}_{N}=\left(\begin{array}[]{ccc}iU_{e1}\sqrt{\frac{m_{1}}{M_{1}}}&iU_{e2}\sqrt{\frac{m_{2}}{M_{2}}}&iU_{e3}\sqrt{\frac{m_{3}}{M_{3}}}\\\
iU_{\mu 1}\sqrt{\frac{m_{1}}{M_{1}}}&iU_{\mu
2}\sqrt{\frac{m_{2}}{M_{2}}}&iU_{\mu 3}\sqrt{\frac{m_{3}}{M_{3}}}\\\ iU_{\tau
1}\sqrt{\frac{m_{1}}{M_{1}}}&iU_{\tau 2}\sqrt{\frac{m_{2}}{M_{2}}}&iU_{\tau
3}\sqrt{\frac{m_{3}}{M_{3}}}\end{array}\right),$ (36)

where $U_{\alpha I}$ are the elements of $U_{\rm PMNS}$. In further
consideration, this case of $\Omega$=I with ”minimal parametric mixing” (36)
is designated as mixing scenario 1. Cosmological consequences of such a
scenario coming from the lifetime of the DM particle $N_{1}$ and the
restriction on the energy fraction carried by the DM particles have been
analysed in [19, 22] where it is shown that the case of moderate $N_{1}$
masses up to 10 keV with normal neutrino mass hierarchy is preferred. Inverse
hierarchy is almost excluded, NHL signals in collider experiments are rather
suppressed due to the mass ratio $(m_{i}/M_{i})^{1/2}$ in the matrix elements
of $\Theta_{min}$. Other forms of mixing have been considered in the existing
literature. In the limiting case of the third HNL decoupling $M_{3}\to\infty$
(two sterile neutrino generations)the 3$\times$3 $\Omega$-matrix proposed in
[23] reads in our notation

$\Omega_{\rm NH}=\left(\begin{array}[]{ccc}1&0&0\\\
0&\cos(\omega)&-\sin(\omega)\\\ 0&\xi\sin(\omega)&\xi\cos(\omega)\\\
\end{array}\right)\hskip 28.45274pt\Omega_{\rm
IH}=\left(\begin{array}[]{ccc}0&\cos(\omega)&-\sin(\omega)\\\
0&\xi\sin(\omega)&\xi\cos(\omega)\\\ 1&0&0\\\ \end{array}\right)$ (37)

for normal hierarchy (NH) and inverted hierarchy (IH), respectively. The case
of only two HNL generations in comparison with the case of three generations
has been analysed within the cosmological framework in [5]. Significant
differences from the case of $\Omega=$I appear with complex-valued $\omega$
parameter which lead to the factors $X_{\omega}=exp(Im(\omega))$ in the mixing
matrix $\Theta$. Three new parameters are introduced in this case, $\xi=\pm
1$, $Re(\omega)$ and $Im(\omega)$. Assuming the framework of thermal
leptogenesis, when the lepton asymmetry is formed due to decays of the right
sterile neutrinos [24], the CP violation parameter reaches sufficiently large
values. Detailed and accurate phenomenological analysis of active and sterile
neutrino mixing in [25] showed that a phenomenologically consistent hierarchy
of mixings $\Theta_{eI}$, $\Theta_{\mu I}$ and $\Theta_{\tau I}$ with
suppressed $\Theta_{eI}$ relative to other matrix elements can be achieved in
a wide interval of $X_{\omega}$ independently on the values of HNL masses.
Translating the experimental upper bounds on $\Theta_{\alpha I}$ from the
shortest possible lifetimes of $N_{2,3}$ from $\pi^{\pm}$ and $K^{\pm}$ meson
decays into the upper bound on $X_{\omega}$ one obtains at the HNL mass scale
102 MeV $Im(\omega)=$4.5 for the lifetime of the order of 1 sec and
$Im(\omega)\sim$7 for the lifetime of the order of 0.01 sec in the case of NH.
Values of the omega parameter greater than seven can lead to a large mixing
parameters of the $W^{\pm},Z$-neutrino interaction not consistent with the
data. In the case of IH the bounds on $Im(\omega)$ are stronger, see details
in [25]. Complex-valued parametrization of $\Omega$ was also used for the
study of HNL properties at the TeV scale [26, 27], see also [28]. In further
consideration, the case of Eq.37 is designated as mixing scenario 2. The
analysis of lepton universality within scenario 2 assuming the $\nu MSM$
framework in the approximation of the two-particle decays was performed in
[16].

For further analysis of the lepton universality within these two scenarios,
calculations are reproduced for two-particle semileptonic HNL decays (Section
2) and calculations are made for three-particle leptonic HNL decays (Section
3). They are used for the lifetime restrictions on the mixing parameter space.
Bounds on the characteristic lepton universality parameter, Eq.1, in the
decays of $K^{\pm}$ and $\pi^{\pm}$ are considered in the framework of both
mixing scenarios in Section 3. In both scenarios defined above the
contribution of $N_{1}$, the DM particle, to the lepton universality parameter
is negligibly small and the degree of violation depends on $N_{2,3}$.

##  2 Semileptonic HNL decays

HNL decays with mesons in the final state can be divided in four groups [29]:

  * •

pseudoscalar neutral meson
$M^{0}_{ps}=\pi^{0},\eta,\eta^{\prime},K^{0},D^{0},B^{0}$ in the final state
with the decay width

$\displaystyle\Gamma(N_{I}\to\nu_{\alpha}M^{0}_{ps})=\frac{G_{F}^{2}M_{I}^{3}}{32\pi}f_{M^{0}_{ps}}|V_{q_{1}q_{2}}|^{2}|\Theta_{\alpha
I}|^{2}\left(1-x_{M^{0}_{ps}}^{2}\right)^{2}.$ (38)

  * •

pseudoscalar charged meson
$M^{\pm}_{ps}=\pi^{\pm},K^{\pm},D^{\pm},D_{s}^{\pm},B^{\pm},B_{s}^{\pm},B_{c}^{\pm}$

$\Gamma(N_{I}\to
l^{-}_{\alpha}M^{+}_{ps})=\frac{G_{F}^{2}M_{I}^{3}}{16\pi}f_{M^{+}_{ps}}|V_{q_{1}q_{2}}|^{2}|\Theta_{\alpha
I}|^{2}\left(\left(1-x_{\alpha}^{2}\right)^{2}-x_{M^{+}_{v}}^{2}\left(1+x_{\alpha}^{2}\right)\right)\times\\\
\times\sqrt{\left(1-(x_{M^{+}_{v}}-x_{\alpha})^{2}\right)\left(1-(x_{M^{+}_{v}}+x_{\alpha})^{2}\right)}.$
(39)

  * •

vector neutral meson $M^{0}_{v}=\rho^{0},\varphi,\omega$

$\Gamma(N_{I}\to\nu_{\alpha}M^{0}_{v})=\frac{G_{F}^{2}M_{I}^{3}}{16\pi}\frac{g^{2}_{M^{0}_{v}}}{m^{2}_{M^{0}_{v}}}|\Theta_{\alpha
I}|^{2}(1+2x^{2}_{M^{0}_{v}})\left(1-x_{M^{0}_{v}}^{2}\right)^{2},$ (40)

  * •

vector charged meson $M^{\pm}_{v}=\rho^{\pm}$ in the final state with the
decay width

$\Gamma(N_{I}\to
l^{-}_{\alpha}M^{+}_{v})=\frac{G_{F}^{2}M_{I}^{3}}{8\pi}\frac{g^{2}_{M^{0}_{v}}}{m^{2}_{M^{0}_{v}}}|\Theta_{\alpha
I}|^{2}\left(\left(1-x_{\alpha}^{2}\right)^{2}+x_{M^{+}_{v}}^{2}\left(1-2x_{M^{+}_{v}}^{2}+x_{\alpha}^{2}\right)\right)\times\\\
\times\sqrt{\left(1-(x_{M^{+}_{v}}-x_{\alpha})^{2}\right)\left(1-(x_{M^{+}_{v}}+x_{\alpha})^{2}\right)}.$
(41)

where $G_{F}$ is the Fermi constant, $V_{qq}$ is the Cabibbo-Kobayashi-Maskawa
(CKM) matrix element, $f_{M}$ and $g_{M}$ are the corresponding meson decay
constants [30], $x_{M}=m_{M}/M_{HNL}$, $x_{\alpha}=m_{l}/M_{HNL}$. These two-
particle widths, see Fig.1 and Fig.2, introduce essential contributions to the
total $N_{2}$ width starting from the $\pi^{0}$ threshold.

Figure 1: Branching ratios for semileptonic two-body decays of $N_{2}$ with
charged pseudoscalar or vector mesons and charged lepton in the final state.
Upper plot – scenario 1 for both the normal and the inverted hierarchy. Middle
plot – scenario 2 with $\Omega=\Omega_{NH}(\xi,\omega)$ where $\xi=+1$,
$Re(\omega)=0$ and $Im(\omega)=7$. This parameter set comes from the Big Bang
nucleosynthesis (BBN) limits, see Fig.5. Bottom plot – the same as the middle
plot but for the case of IH. Solid lines correspond to $e^{-}$ in final state,
long-dashed lines – to $\mu^{-}$, dotted lines – to $\tau^{-}$. Blue lines –
decay with $\pi^{+}$, green – $K^{+}$, orange – $\rho^{+}$, purple – $D^{+}$,
red – $B^{+}$. Figure 2: Branching ratios for semileptonic two-body decays of
$N_{2}$ with neutral pseudoscalar or vector mesons and active neutrino in the
final state. Curves for NH and IH practically coincide. Upper plot, middle
plot and bottom plot are for the same scenarios as in Fig.1. Black lines –
decay with $\pi^{0}$, red – $K^{0}$, green – $\eta$, blue – $\eta^{\prime}$,
orange – $D^{0}$, purple – $B^{0}$, dashed black – $\rho^{0}$, dashed gray –
$\varphi$.

##  3 Leptonic HNL decays

In this section we calculate squared amplitudes and HNL decay widths for the
four-fermion effective interaction terms when $M_{1,2,3}\ll M_{W},M_{Z}$.
Three different decay amplitudes with respect to the mixing factors can be
distinguished, $N_{I}\to\upnu_{i},\upnu_{j},\upnu_{j};$
$N_{I}\to\upnu_{i},l^{+}_{\beta},l^{-}_{\beta}$ and
$N_{I}\to\upnu_{i},l^{+}_{\alpha},l^{-}_{\beta}$. Three-particle decay widths
are calculated symbolically ($p=k_{1}+k_{2}+k_{3}$, $p^{2}=M_{I}^{2}$,
$k_{i}^{2}=m_{i}^{2}$, i=1,2,3) keeping all masses of leptons nonzero by
integrating in invariant variables over the Dalitz plot. The details of these
calculations are given in Appendix B. In the model under consideration, all
six neutrinos are Majorana neutrinos, which requires careful determination of
the signs of the interference terms, see [31], also [38]. Three cases for
amplitudes take place:

Case 1. Three active neutrino mass states appear in the final state:
$N_{I}(p)\to\upnu_{i}(k_{1}),\upnu_{j}(k_{2}),\upnu_{j}(k_{3})$, $i,j=1,2,3$.
Insofar as the approximation $\mathcal{O}(\theta^{2})$ is used, the PMNS
matrix can be operated as a unitary one $(U^{\dagger}U)=\delta_{ij}$ and there
are no $Z\upnu_{j}\upnu_{k\neq j}$ vertices.

(120,120) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
plain,label=$\upnu_{i}(k_{1})$, label.side=right, label.dist=0.1cmo4,v2
plaino1,v1,o3 plain,tension=3,label=$\upnu_{j}(k_{3})$, label.side=right,
label.dist=0.1cmo1,v1 plain,label=$\upnu_{j}(k_{2})$, label.side=right,
label.dist=0.1cmv1,o3 photon, tension=2,label=$Z$, label.side=right,
label.dist=0.2cmv2,v1

Diagram for the case 1. HNL decay into three active Majorana neutrinos.

Neglecting the active neutrino masses for the sake of clarity we get the
squared amplidute of the following structure

$|M_{3\upnu}|^{2}=64G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}~{}\left[(pk_{3})(k_{1}k_{2})+(pk_{2})(k_{1}k_{3})\right]$

Details of the calculations can be found in Appendix A. The mixing factor
$|(U^{\dagger}\Theta)_{iI}|^{2}$ is reduced taking the sum over the active
neutrino mass states

$\sum\limits_{i}|(U^{\dagger}\Theta)_{iI}|^{2}=\sum\limits_{i}\sum\limits_{\alpha,\beta}U^{*}_{\alpha
i}\Theta_{\alpha I}U_{\beta i}\Theta^{*}_{\beta
I}=\sum\limits_{\alpha,\beta}(UU^{\dagger})_{\alpha\beta}\Theta_{\alpha
I}\Theta^{*}_{\beta
I}=\sum\limits_{\alpha,\beta}\delta_{\alpha\beta}\Theta_{\alpha
I}\Theta^{*}_{\beta I}=\sum\limits_{\alpha}|\Theta_{\alpha I}|^{2}$

and the decay width

$\Gamma(N_{I}\to\sum\limits_{i}\upnu_{i},\upnu_{j},\upnu_{j})=\frac{G_{F}^{2}M_{I}^{5}}{192\pi^{3}}\sum\limits_{\alpha=e,\mu,\tau}|\Theta_{\alpha
I}|^{2}.$ (42)

Case 2. Two leptons of different flavors and a neutrino mass state appear in
the final state:
$N_{I}(p)\to\upnu_{j}(k),l_{\alpha\neq\beta}^{+}(p_{1}),l_{\beta}^{-}(p_{2})$.
In this charged current amplitude with a neutral Majorana particle decay the
interference term appears between diagrams with intermediate $W^{+}$ and
$W^{-}$,

(100,100) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=right, label.dist=0.1cmo4,v2
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=right,
label.dist=0.1cmv1,o3 photon, tension=2,label=$W^{+}$, label.side=right,
label.dist=0.2cmv2,v1 (100,100) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=left, label.dist=0.1cmv2,o4
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=left,
label.dist=0.1cmo3,v1 photon, tension=2,label=$W^{-}$, label.side=right,
label.dist=0.2cmv2,v1

Diagrams for the case 2. HNL decay into Dirac charged lepton and antilepton of
different flavors associated with one active Majorana neutrino

$|M_{W}|^{2}=128G_{F}^{2}\left[|\Theta_{\alpha I}|^{2}|U_{\beta
i}|^{2}(pp_{2})(p_{1}k)+|\Theta_{\beta I}|^{2}|U_{\alpha
i}|^{2}(pp_{1})(p_{2}k)\right].$

After summing by the active neutrino mass states and using $UU^{\dagger}=$I
due to approximate unitarity of PMNS matrix, the decay width takes the form

$\Gamma(N_{I}\to\sum\limits_{i=1,2,3}\upnu_{i}~{}l^{+}_{\alpha\neq\beta}l^{-}_{\beta})=\frac{G_{F}^{2}M_{I}^{5}}{192\pi^{3}}\left(|\Theta_{\alpha
I}|^{2}+|\Theta_{\beta I}|^{2}\right)\mathcal{G}(r_{\alpha},r_{\beta}),$ (43)

where

$\displaystyle\mathcal{G}(x,y)$ $\displaystyle=$
$\displaystyle(1-7x-7x^{2}+x^{3}+12xy-7y-7y^{2}+y^{3}-7x^{2}y-7xy^{2})R+$
$\displaystyle+12(y^{2}+x^{2}y^{2}-2x^{2})\ln\left(\frac{1+x-y+R}{2}\right)+12x^{2}(1-y^{2})\ln(\frac{1}{x})+$
$\displaystyle+12y^{2}(1-x^{2})\ln\left(\frac{1-x-y+R}{1-x+y-R}\right),~{}~{}~{}~{}R=\sqrt{1-2x+x^{2}-2y+y^{2}-2xy}$

In the limiting case $m_{e}\ll m_{\mu}$, $m_{\mu}\ll m_{\tau}$ Eq.43 is
reduced to

$\Gamma(N_{I}\to\sum\limits_{i=1,2,3}\upnu_{i}~{}l^{+}_{\alpha\neq\beta}l^{-}_{\beta})\Big{|}_{m_{\beta}\to
0}=\Gamma(N_{I}\to\sum\limits_{i=1,2,3}\upnu_{i}~{}l^{+}_{\alpha\neq\beta}l^{-}_{\beta})\Big{|}_{m_{\alpha}\to
0}=\\\ =\frac{G_{F}^{2}M_{I}^{5}}{192\pi^{3}}\left(|\Theta_{\alpha
I}|^{2}+|\Theta_{\beta
I}|^{2}\right)\left(1-8r+8r^{3}-r^{4}-12r^{2}\ln(r)\right)$ (44)

where $r=\frac{m^{2}}{M_{I}^{2}}$, $m=max\\{m_{\alpha},m_{\beta}\\}$. This
result is the same as the one obtained in [29].

To account for charge conjugated final states, the above width must be
multiplied by a factor of two.

Case 3. Two leptons of the same flavor and a neutrino mass state appear in the
final state:
$N_{I}(p)\to\upnu_{j}(k),l_{\alpha}^{+}(p_{1}),l_{\alpha}^{-}(p_{2})$. The
amplitude contains both charged and neutral currents and includes in this case
three interfering diagrams

(100,100) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=right, label.dist=0.1cmo4,v2
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=right,
label.dist=0.1cmv1,o3 photon, tension=2,label=$W^{+}$, label.side=right,
label.dist=0.2cmv2,v1 (100,100) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=left, label.dist=0.1cmv2,o4
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=left,
label.dist=0.1cmo3,v1 photon, tension=2,label=$W^{-}$, label.side=right,
label.dist=0.2cmv2,v1 (100,100) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
plain,label=$\upnu_{k}(k)$, label.side=right, label.dist=0.1cmo4,v2
fermiono1,v1,o3 fermion,tension=3,label=$l_{\alpha}^{+}(p_{1})$,
label.side=right, label.dist=0.2cmo1,v1 fermion,label=$l_{\alpha}^{-}(p_{2})$,
label.side=right, label.dist=0.2cmv1,o3 photon, tension=2,label=$Z$,
label.side=right, label.dist=0.2cmv2,v1

Diagrams for the case 3. HNL decay into Dirac charged lepton and antilepton
with the same flavor and one active Majorana neutrino.

The mixing factor appearing in the squared amplitude is reduced by summing
over the states of active neutrinos222The squared amplitude not summed up by
neutrino mass states $\upnu_{1,2,3}$ and details of calculation can be found
in Appendix A.

$\sum\limits_{i}\Theta_{\alpha I}U^{*}_{\alpha
i}(U^{\dagger}\Theta)_{iI}^{*}=\sum\limits_{i,\beta}\Theta_{\alpha
I}U^{*}_{\alpha i}U_{\beta i}\Theta_{\beta
I}=\sum\limits_{\beta}(UU^{\dagger})_{\alpha\beta}\Theta_{\alpha
I}\Theta^{*}_{\beta I}=\\\
=\sum\limits_{\beta}\delta_{\alpha\beta}\Theta_{\alpha I}\Theta^{*}_{\beta
I}=|\Theta_{\alpha I}|^{2}.$ (45)

Integration of the amplitude

$\displaystyle\sum\limits_{i}|M_{\rm WZ}|^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}\left[\left(\mathcal{C}_{1}^{2}+\mathcal{C}_{2}^{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}+\left(1-2\mathcal{C}_{1}\right)|\Theta_{\alpha
I}|^{2}\right]~{}\left[(pp_{2})(p_{1}k)+(pp_{1})(p_{2}k)\right]+$
$\displaystyle+128G_{F}^{2}\left[\left(4\mathcal{C}_{1}\mathcal{C}_{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}-4\mathcal{C}_{2}^{2}|\Theta_{\alpha
I}|^{2}\right]~{}\frac{m_{\alpha}^{2}}{2}(pk)$

gives the decay width

$\displaystyle\Gamma(N_{I}\to\sum\limits_{i=1,2,3}\upnu_{i}~{}l^{+}_{\alpha}l^{-}_{\alpha})$
$\displaystyle=$
$\displaystyle\frac{G_{F}^{2}M_{I}^{5}}{96\pi^{3}}\left[\left(\mathcal{C}_{1}^{2}+\mathcal{C}_{2}^{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}+\left(1-2\mathcal{C}_{1}\right)|\Theta_{\alpha
I}|^{2}\right]\mathcal{F}_{1}(r)+$ (46)
$\displaystyle+\left[\left(2\mathcal{C}_{1}\mathcal{C}_{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}-2\mathcal{C}_{2}^{2}|\Theta_{\alpha
I}|^{2}\right]~{}\mathcal{F}_{2}(r)$

where $\mathcal{C}_{1}=s_{W}^{2}-\frac{1}{2}$, $\mathcal{C}_{2}=s_{W}^{2}$,
$r=\frac{m_{\alpha}^{2}}{M_{I}^{2}}$ ($s_{W}=\sin\theta_{W}$ where
$\theta_{W}$ is the Weinberg angle),

$\displaystyle\mathcal{F}_{1}(r)$ $\displaystyle=$
$\displaystyle\left(1-14r-2r^{2}-12r^{3}\right)\sqrt{1-4r}~{}+~{}12r^{2}\left(1-r^{2}\right)\ln\left(\frac{1-3r+(1-r)\sqrt{1-4r}}{r(1-\sqrt{1-4r})}\right),$
$\displaystyle\mathcal{F}_{2}(r)$ $\displaystyle=$
$\displaystyle(2r+10r^{2}-12r^{3})\sqrt{1-4r}~{}-~{}(6r^{2}-12r^{3}+12r^{4})\ln\left(\frac{1-3r+(1-r)\sqrt{1-4r}}{r(1-\sqrt{1-4r})}\right).$

Case 4. Decay channels with HNL in the final state $N_{2,3}\to N_{1}+...$ are
suppressed by the factor of $\sim\mathcal{O}(\Theta^{2})$ in comparison with
the channels described above. They are insignificant for the following
analysis.

Figure 3: Branching ratios for pure leptonic three-body decays of $N_{2}$.
Upper plot, middle plot and bottom plot describe the same scenarios as in
Fig.1 and Fig.2, respectively. Gray lines – $\upnu\upnu\upnu$ decay mode,
green lines – $\upnu e^{-}e^{+}$, purple lines – $\upnu\mu^{-}\mu^{+}$, dashed
purple lines – $e^{-}\mu^{+}$, dashed blue – $\upnu e^{-}\tau^{+}$, dashed
orange – $\upnu\mu^{-}\tau^{+}$, red lines – $\upnu\tau^{-}\tau^{+}$. mode.
There is a very weak sensitivity of partial widths in relation to the choice
of the mixing scenario. Figure 4:  The lifetime of $N_{2}$ as a function of
its mass in the case of normal (solid green line) and inverted (dashed green
line) hierarchies, scenario 1. Horizontal lines correspond to BBN constraints
of $\tau<1$ sec, solid black line and $\tau<0.1$ sec, dashed black line. These
bounds are used in [16]. The bound $\tau<0.02$ sec, dotted gray line, is taken
from [13]. Figure 5: The lifetime of $N_{2}$ as a function of its mass in
scenario 2 for the cases of normas and inverted hierarchies, case (1)
$\Omega=\Omega_{NH}(\xi=\pm 1,\omega)$; case (2) $\Omega=\Omega_{IH}(\xi=\pm
1,\omega)$ with $Re(\omega)=0$ (green lines) and $Re(\omega)=\pi/2$ (brown
lines). Solid lines correspond to the border of the strong BBN limit
$\tau<0.02$ sec, long-dashed – $\tau=0.1$ sec and short-dashed – $\tau=1$ sec.
The domains from the lower left corner to the corresponding border line are
excluded.

##  4 Restrictions on the lepton universality parameter in the decays of
$\pi^{\pm}$ and $K^{\pm}$ mesons

Limitations on the lepton universality parameter are imposed by the
restrictions on the $N_{2}$ and $N_{3}$ lifetime from the Big Bang scenario
and experimental restrictions from the resonance observation in meson decays
by PS191, E949 and NA62 collaborations [32]. Limitations on the part of
collaborations require careful interpretation when brought to a common
denominator, since the so-called ”model-independent approach” (see e.g. [28])
is used in data processing, when the mixing parameter and the HNL mass are
independent variables. This assumption seems to be quite strong even within
the framework of consideration with one generation of leptons, see for example
[33]. Due to these circumstances, we do not specify exactly the lower boundary
of the experimentally allowed region in Fig.5, limiting ourselves only to the
remark that it lies at the level $Im(\omega)>$7 at $M_{HNL}<$100 MeV in the
case of NH, slowly going down at higher masses. In the case of IH this bound
is stronger. Cosmological considerations imposing restrictions on the lifetime
of $N_{2}$ and $N_{3}$ on the level of $\tau_{N_{2,3}}<$0.1 – 1 sec [34] were
recently improved in [13] to the minimal level of 0.02 sec. These bounds are
model dependent and obtained in the framework of a rather specific scenarios
of the Big Bang nucleosynthesis.

In the scenario 1 framework the minimmal mixing matrix (36) does not contain
redundant parameters, so the lifetime dependence on the HNL mass is
unambiguous, see Fig.4. In the scenario 2 first the allowed lifetime domain on
the $|Im(\omega)|-M$ plane is found, see Fig.5, which is then translated to
the allowed domain for the lepton universality parameter $\Delta r_{M}$, see
Fig.6 and Fig.7.

Although the calculations were performed taking into account all the
abovementioned decay modes, it is instructive to write out simple formulas for
the case of leptonic two-particle decays. In the effective four-fermion
approximation the width of the scalar meson $M=\pi^{\pm},K^{\pm}$

$\Gamma(M^{+}\to
l_{\alpha}^{+},N_{I})=\frac{G_{F}^{2}f_{M}^{2}}{4\pi}|\Theta_{\alpha
I}|^{2}m_{M}^{4}\lambda^{1/2}(1,r_{I},r_{\alpha})\left[r_{I}+r_{\alpha}-(r_{I}-r_{\alpha})^{2}\right]$
(47)

where $G_{F}$ is Fermi constant and $f_{M}$ is meson formfactor,
$\alpha=e,\mu,\tau$ (only channels allowed by energy-momentum conservation are
admitted), $I=1,2,3$ number of HNL generation, $j=1,2,3$ active neutrino mass
states; $r_{\alpha}=m^{2}_{\alpha}/m^{2}_{M}$, $r_{I}=M_{I}^{2}/m_{M}^{2}$,
$\lambda(a,b,c)=a^{2}+b^{2}+c^{2}-2ab-2bc-2ac$.
In the case of only one active neutrino in the final state the mass
corrections are neglected and

$\Gamma(M^{+}\to
l_{\alpha}^{+},\upnu_{j})=\frac{G_{F}^{2}f_{M}^{2}}{4\pi}|U_{\alpha
i}|^{2}m_{M}^{4}r_{\alpha}(1-r_{\alpha})^{2}.$ (48)

Convenient variable to analyse BSM deviations from lepton universality is the
ratio [11]

$R_{M}=\frac{\sum_{i}\Gamma(M^{+}\to e^{+}\upnu_{i})+\sum_{j}\Gamma(M^{+}\to
e^{+}N_{j})}{\sum_{i}\Gamma(M^{+}\to\mu^{+}\upnu_{i})+\sum_{j}\Gamma(M^{+}\to\mu^{+}N_{j})}$
(49)

or its derivative demonstrating the deviation of the ratio from zero

$\Delta
r_{M}\equiv\frac{R_{M}}{R^{SM}_{M}}\,-1=\frac{\sum_{i}|U_{ei}|^{2}+\sum_{I}|\Theta_{eI}|^{2}G^{M}_{eI}}{\sum_{i}|U_{\mu
i}|^{2}+\sum_{I}|\Theta_{\mu I}|^{2}G^{M}_{\mu I}}\,-1.$ (50)

where masses of active neutrinos are neglected and only decays which are
allowed kinematically $M_{I}<m_{M}-m_{\alpha}$ are accounted for. If
$M_{I}>m_{M}-m_{\alpha}$ then $G_{\alpha I}^{M}=0$. The SM parameter and the
BSM parameter are

$\displaystyle R_{M}^{(0)}$ $\displaystyle=$
$\displaystyle\frac{r_{\mu}(1-r_{\mu})^{2}}{r_{e}(1-r_{e})^{2}},$
$\displaystyle G^{M}_{\alpha I}$ $\displaystyle=$
$\displaystyle\frac{\lambda^{1/2}(1,r_{I},r_{\alpha})\left[r_{I}+r_{\alpha}-(r_{\alpha}-r_{I})^{2}\right]}{r_{\alpha}(1-r_{\alpha})^{2}}.$

Figure 6: Lepton universality parameters $\Delta r_{\pi}$, red line, and
$\Delta r_{K}$, blue line, scenario 1, as a function of the HNL mass in the
case of quasidegenerate $M_{2}\simeq M_{3}\equiv M$, scenario 1. First
generation HNL $N_{1}$ with the mass $M_{1}\simeq 5$ keV is the dark matter
particle. Active neutrino masses for the normal hierarchy (NO) are $m_{1}=0$,
$m_{2}=\sqrt{\Delta m^{2}_{sol}}$, $m_{3}=\sqrt{\Delta m^{2}_{atm}+\Delta
m^{2}_{sol}}$, for the inverted hierarchy (IO) $m_{3}=0$, $m_{2}=\sqrt{\Delta
m^{2}_{atm}+\Delta m^{2}_{sol}}$, $m_{1}=\sqrt{\Delta m^{2}_{atm}}$.
Horisontal dotted lines demonstrate possible experimental accuracy of $\Delta
r_{\pi}$ (blue) and $\Delta r_{K}$ (red) measurement, see [16].  Figure 7:
Lepton universality parameters $\Delta r_{\pi}$, red line, and $\Delta r_{K}$,
blue line, scenario 2, as a function of the HNL mass in the case of
quasidegenerate $M_{2}\simeq M_{3}\equiv M$. The choice of $M_{1}$ and the
active neutrino masses for both normal and inverted hierarchies are the same
as in Fig.6. The mixing matrix corresponds to Eq.(37). The domains lower than
contours are excluded by the lifetime restrictions $\tau<1$ sec (weak BBN
constraint used in [16]) and $\tau<0.02$ sec (improved BBN constraint from
[13]).

Using the unitarity condition for the full $6\times 6$ mixing matrix

$\sum\limits_{i=1}^{3}|U_{\alpha i}|^{2}+\sum\limits_{I=1}^{3}|\Theta_{\alpha
I}|^{2}=1$ (51)

one can rewrite Eq.50 in the form

$\Delta
r_{M}=\frac{1+\sum\limits_{I=1}^{3}|\Theta_{eI}|^{2}(G^{M}_{eI}-1)}{1+\sum\limits_{I=1}^{3}|\Theta_{\mu
I}|^{2}(G^{M}_{\mu I}-1)}-1$ (52)

If the HNL mass is below the $m_{M}-m_{\mu}$ threshold, $\Delta r_{M}$ is
positive and proportional to the sum of
$|\Theta_{lI}|^{2}(M_{HNL}/m_{l})^{2}\,(1-M^{2}_{HNL}/M^{2}_{M})^{2}$,
$I=2,3$, $l=e$ or $\mu$, the mass ratio squared factor provides a strong
enhancement factor with respect to the mixing factor constrained from above by
the data. If the NHL mass is above the $m_{M}-m_{e}$ threshold, both leptonic
modes are closed kinematically and nonzero $\Delta r_{M}$ appears due to
nonunitarity of the PMNS mixing matrix, see Eq.51. In this case $\Delta r_{M}$
is not enhanced by the mass ratio squared and defined by the sum of
$|\Theta_{eI}|^{2}-|\Theta_{\mu I}|^{2}$, $I=2,3$. Experimental upper limit on
the mixing matrix elements was found only for some intervals of the HNL mass
less than 450 MeV (PS191 and NA62 collaborations). If $M_{HNL}>$450 MeV, there
is no experimental limitation and the BBN lifetime restrictions play a major
role.

##  5 Summary

Limitations of the possible form of the mixing matrix in the see-saw type I
models, which are of undoubted interest, can be obtained by combining
constraints on the matrix elements imposed from above by the absence of
signals on colliders and beam dump experiments, while constraints on the
mixing from below are provided by baryogenesis scenarios within the framework
of the Big Bang concept, as well as restrictions on flavor oscillations.
However, the numerical values for the boundaries strongly depend on the mixing
scenario within a particular model which, from a technical point of view, is
ambiguously implemented by a certain choice of the $\Omega$-matrix in the
$M_{D}-U_{PMNS}-U_{N}$ connection, Eq.12. The lepton universality violation
parameter in $K^{\pm}$ and $\pi^{\pm}$ decays is sensitive to the neutral
heavy leptons due to their additional contributions to meson decay widths,
dependent on the mixing matrix.

Two mixing scenarios analysed above are rather different, in the ”minimal
parametric mixing”, scenario 1, the mixing matrix depends on the mass ratio
$(m_{i}/M_{i})^{1/2}$ of active neutrino and HNL, which suppresses any HNL
production or decay channel and imposes rather strict restrictions on $N_{1}$
from the cosmological lifetime and DM energy fraction requirements, not
affecting $N_{2}$ and $N_{3}$ masses. In the scenario 2 with three additional
parameters for the $N_{2}$ and $N_{3}$ mixing, intensively discussed in the
literature, a sort of tuning is needed to generate small masses of active
neutrinos by means of see-saw type I and the baryon asymmetry of the Universe
by means of $N_{2}-N_{3}$ flavor oscillation mechanism. HNL production is
enhanced by the mixing parameter $exp(Im(\omega))$ and stronger restrictions
on the $N_{2}$ and $N_{3}$ masses ($N_{2}\sim N_{3}$) from below are imposed
not affecting to large extent $N_{1}$ which has a negligible mixing. Strong
hierarchy of $\Theta_{eI}$ and $\Theta_{\mu I}$, $\Theta_{\tau I}$ elements of
the mixing matrix naturally appears in scenario 1 and can be easily configured
for the NH case in scenario 2. The cosmological upper bound on the HNL
lifetime is critical for determining the bound on the masses of $N_{2}$ and
$N_{3}$, so their lifetime was carefully evaluated taking into account the
three-particle leptonic decay channels which are dominant below the thresholds
of the two-particle channels.

The lepton universality violation parameter $\Delta r_{M}$,
$M=\pi^{\pm},K^{\pm}$ does not exceed ${\cal O}(10^{-5})$ at $M_{2}$ and
$M_{3}$ masses of the order of 102 MeV in the scenario 1 both for normal and
inverted hierarchy, demonstrating very low sensitivity to the cosmological
lifetime restrictions which are essential only in the $N_{2}$ and $N_{3}$
masses range above 1 GeV. The situation is more involved in the three-
parametric scenario 2 where the lifetime bound is given by an exclusion
contours on the mixing – mass, $|Im(\omega)|-M$ plane which are strongly
parameter dependent, especially for the IH case when $Im(\omega)$ parameter
value in comparison with the NH case parameter value can be several orders of
magnitude different. Similar NH/IH dependency behavior is observed only if
$|Im(\omega)|>$5\. Strong mass dependence of the universality parameter
$\Delta r_{M}$ is observed, allowing to impose restrictions on the HNL masses
of the second and third generation. Using the theoretical SM value for
$R_{\pi}$ from [35] and the experimental value from Particle Data Group (PDG)
[36], the lepton universality parameter for $\pi^{\pm}$ is found to be $\Delta
r_{\pi}=(-4\pm 3)\times 10^{-4}$, so $M_{2,3}>$133 MeV bound is consistent
with $\Delta r_{\pi}<10^{-4}$ at 2$\sigma$ level, see [16], demonstrating weak
dependence on the neutrino hierarchy and the $N_{2,3}$ lifetime restrictions.
Using the theoretical SM value for $R_{K}$ [37, 35] and the experimental PDG
value for $K^{\pm}$ parameter, $\Delta r_{K}$ is found to be $\Delta
r_{K}=(4\pm 4)\times 10^{-3}$, so $M_{2,3}$ bound varies from 170 MeV to 320
MeV significantly depending on the hierarchy and the HNL lifetime
restrictions. It follows that the currently achieved level of experimental
accuracy is quite moderate, as a result of which future experiments of high
precision have a great potential for the discovery of BSM physics lepton
mixing scenarios.

Acknowlegement  This work was supported by the Russian Science Foundation
Grant No. 22-12-00152.

## Appendix A. Amplitudes for the three-particle leptonic decays

Models with Majorana fermions have a number of features of the use of diagram
technique, see [31, 38]. Diagram technique for Majorana fermions proposed in
[38] where the charge conjugation matrix was explicitly included in the
Feynman rules, has been modified in [31] where standard propagators are used,
vertices do not include the charge conjugation matrix and a specific fermion
flow defined for fermion lines is introduced. This Appendix provides the
specifics of calculations in the cases under consideration.

Trace/combination of traces |  Result of calculation
---|---
$Tr[\hat{A}\gamma_{\alpha}P_{L(R)}\hat{B}\gamma^{\beta}P_{L(R)}\hat{C}\gamma^{\alpha}P_{L(R)}\hat{D}\gamma_{\beta}P_{L(R)}~{}]$ |  $-16(AC)(BD)$
$Tr[\hat{A}\gamma_{\alpha}P_{L(R)}\hat{B}\gamma_{\beta}P_{L(R)}]\cdot Tr[\hat{C}\gamma^{\alpha}P_{L(R)}\hat{D}\gamma^{\beta}P_{L(R)}]$ |  $16(AC)(BD)$
$Tr[\hat{A}\gamma_{\alpha}P_{L}\hat{B}\gamma_{\beta}P_{L}]Tr[\hat{C}\gamma^{\alpha}P_{R}\hat{D}\gamma^{\beta}P_{R}]$ |  $16(AD)(BC)$
$Tr[\hat{A}\gamma_{\alpha}P_{R(L)}\hat{B}P_{L(R)}\hat{C}\gamma^{\alpha}P_{L(R)}\hat{D}P_{R(L)}]$ |  $8(AD)(BC)$
$Tr[\hat{A}\gamma_{\alpha}P_{L(R)}\hat{B}\gamma_{\beta}P_{L(R)}]~{}Tr[\gamma^{\alpha}P_{L(R)}\gamma^{\beta}P_{R(L)}]$ |  $-8(AB)$
$Tr[\hat{A}\gamma_{\alpha}P_{L(R)}\hat{B}\gamma_{\beta}P_{L(R)}]~{}Tr[\gamma^{\alpha}P_{R(L)}\gamma^{\beta}P_{L(R)}]$ |  $-8(AB)$
$Tr[\hat{A}\gamma_{\alpha}P_{L(R)}\gamma_{\beta}P_{R(L)}\gamma^{\alpha}P_{L(R)}\hat{D}\gamma^{\beta}P_{L(R)}]$ |  $8(AD)$

Table 1: Formulas for all nonzero traces (combination of traces) of gamma
matrices arising during the calculation in cases 1,2,3

.

(70,70) i1 o1,o2 photon,label.side=left,label=$Z$i1,w1 plain,label.side=right,label=$\upnu_{j}$w1,o1 plain,label.side=left,label=$N_{I}$w1,o2 lab=$\circlearrowright$w1 | $~{}~{}~{}-i\frac{g}{2c_{W}}\gamma^{\mu}\left[(U^{\dagger}\Theta)^{*}_{jI}P_{L}-(U^{\dagger}\Theta)_{jI}P_{R}\right]$
---|---
(70,70) i1 o1,o2 photon,label.side=right,label=$Z$i1,w1 plain,label.side=right,label=$\upnu_{i}$w1,o1 plain,label.side=left,label=$\upnu_{i}$w1,o2 lab=$\circlearrowright$w1 | $-i\frac{g}{2c_{W}}\gamma^{\mu}(P_{L}-P_{R})$
(70,70) i1 o1,o2 photon,label.side=left,label=$W^{+}$i1,w1 plain,label.side=right,label=$N_{I}$w1,o1 fermion,label.side=left,label=$l_{\alpha}^{-}$w1,o2 lab=$\circlearrowright$w1 | $-i\frac{g}{\sqrt{2}}\Theta_{\alpha I}\gamma^{\mu}P_{L}$
(70,70) i1 o1,o2 photon,label.side=left,label=$W^{-}$i1,w1 plain,label.side=right,label=$N_{I}$w1,o1 fermion,label.side=right,label=$l_{\alpha}^{+}$o2,w1 lab=$\circlearrowright$w1 | $+i\frac{g}{\sqrt{2}}\Theta_{\alpha I}^{*}\gamma^{\mu}P_{R}$
(70,70) i1 o1,o2 photon,label.side=left,label=$Z$i1,w1 fermion,label.side=left,label=$l_{\alpha}^{+}$o1,w1 fermion,label.side=left,label=$l_{\alpha}^{-}$w1,o2 lab=$\circlearrowright$w1 | $~{}~{}-i\frac{g}{2c_{W}}\gamma^{\mu}\left((2s_{W}^{2}-1)P_{L}+2s_{W}^{2}P_{R}\right)$
Table 2: Feynman rules for Majorana fermions used in calculations of pure
leptonic decay amplitudes. The direction of the fermion flow is chosen from
the bottom up, as indicated by an arrow near the vertex,
$P_{R,L}=(1\pm\gamma_{5})/2$, $s_{W}=\sin\theta_{W}$.

Case 1: Here we choose the direction of the fermion flow (see [31]) as shown
on diagram A1, then

(120,120) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
plain,label=$\upnu_{i}(k_{1})$, label.side=right, label.dist=0.1cmo4,v2
plaino1,v1,o3 plain,tension=3,label=$~{}\upnu_{j}(k_{3})$, label.side=right,
label.dist=0.1cmo1,v1 plain,label=$~{}\upnu_{j}(k_{2})$, label.side=right,
label.dist=0.1cmv1,o3 photon, tension=2,label=$Z$, label.side=right,
label.dist=0.2cmv2,v1 lab=$\leftarrow$v2 lab=$\circlearrowleft$v1

Diagram A1. Diagram for the case 1 with with an explicitly specified choice of
fermion flow direction indicated by arrows.

$\displaystyle
M_{3\upnu}=-i\sqrt{2}G_{F}(\overline{v}_{N}\gamma_{\mu}\left[(U^{\dagger}\Theta)^{*}_{iI}P_{L}-(U^{\dagger}\Theta)_{iI}P_{R}\right]v_{1})(\overline{u}_{3}\gamma^{\mu}(P_{L}-P_{R})v_{2})$
(53)

Denoting the terms of the full decay amplitude as

$\displaystyle M_{LL}$ $\displaystyle=$
$\displaystyle-i\sqrt{2}G_{F}(U^{\dagger}\Theta)_{iI}^{*}(\overline{v}_{N}\gamma_{\mu}P_{L}v_{1})(\overline{u}_{3}\gamma^{\mu}P_{L}v_{2})$
$\displaystyle M_{LR}$ $\displaystyle=$
$\displaystyle-i\sqrt{2}G_{F}(U^{\dagger}\Theta)_{iI}^{*}(\overline{v}_{N}\gamma_{\mu}P_{L}v_{1})(\overline{u}_{3}\gamma^{\mu}P_{L}v_{2})$
$\displaystyle M_{RL}$ $\displaystyle=$
$\displaystyle-i\sqrt{2}G_{F}(U^{\dagger}\Theta)_{iI}^{*}(\overline{v}_{N}\gamma_{\mu}P_{L}v_{1})(\overline{u}_{3}\gamma^{\mu}P_{L}v_{2})$
$\displaystyle M_{RR}$ $\displaystyle=$
$\displaystyle-i\sqrt{2}G_{F}(U^{\dagger}\Theta)_{iI}^{*}(\overline{v}_{N}\gamma_{\mu}P_{L}v_{1})(\overline{u}_{3}\gamma^{\mu}P_{L}v_{2}),$

one can observe that there is no interference in the approximation of zero
masses of active neutrinos in the final state

$\displaystyle
Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{L(R)}\hat{k}_{1}\gamma_{\nu}P_{R(L)}]=0,$
$\displaystyle
Tr[\hat{k}_{3}\gamma^{\mu}P_{L(R)}\hat{k}_{2}\gamma^{\nu}P_{R(L)}]=0.$

For the amplitude terms we get

$\displaystyle|M_{LL}|^{2}$ $\displaystyle=$ $\displaystyle
2G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{L}\hat{k}_{1}\gamma_{\nu}P_{L}]~{}Tr[\hat{k}_{3}\gamma^{\mu}P_{L}\hat{k}_{2}\gamma^{\nu}P_{L}]=$
$\displaystyle=$ $\displaystyle
32G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}(pk_{3})(k_{1}k_{2}),$
$\displaystyle|M_{RL}|^{2}$ $\displaystyle=$ $\displaystyle
2G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{R}\hat{k}_{1}\gamma_{\nu}P_{R}]~{}Tr[\hat{k}_{3}\gamma^{\mu}P_{L}\hat{k}_{2}\gamma^{\nu}P_{L}]=$
$\displaystyle=$ $\displaystyle
32G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}(pk_{2})(k_{1}k_{3}),$
$\displaystyle|M_{LR}|^{2}$ $\displaystyle=$ $\displaystyle
2G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{L}\hat{k}_{1}\gamma_{\nu}P_{L}]~{}Tr[\hat{k}_{3}\gamma^{\mu}P_{R}\hat{k}_{2}\gamma^{\nu}P_{R}]=$
$\displaystyle=$ $\displaystyle
32G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}(pk_{2})(k_{1}k_{3}),$
$\displaystyle|M_{RR}|^{2}$ $\displaystyle=$ $\displaystyle
2G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{R}\hat{k}_{1}\gamma_{\nu}P_{R}]~{}Tr[\hat{k}_{3}\gamma^{\mu}P_{R}\hat{k}_{2}\gamma^{\nu}P_{R}]=$
$\displaystyle=$ $\displaystyle
32G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}(pk_{3})(k_{1}k_{2}),$

and finally for the full squared amplitude

$|M_{3\upnu}|^{2}=64G_{F}^{2}|(U^{\dagger}\Theta)_{iI}|^{2}~{}\left[(pk_{3})(k_{1}k_{2})+(pk_{2})(k_{1}k_{3})\right].$

Note that for Dirac active neutrino terms $LR$ and $RR$ do not arise and the
factor in front of right-hand side will be 32. It is also necessary to take
into account the charge conjugate mode by multiplying by two. There is no
charge conjugated mode for Majorana neutrinos.

Case 2. Acting by the same rules in combination with Fiertz transformations,
see Table 2, the amplitudes $M_{1}$ and $M_{2}$ for the diagram with
intermediate $W^{+}$ and the diagram with intermediate $W^{-}$ can be written
as

(120,120) i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=right, label.dist=0.1cmo4,v2
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=right,
label.dist=0.1cmv1,o3 photon, tension=2,label=$W^{+}$, label.side=right,
label.dist=0.2cmv2,v1 lab=$\circlearrowleft$v1 lab=$\leftarrow$v2 (120,120)
i1,i2,id1,id2,i3,i4,i5 o1,o2,od1,od2,o3,o4,o5
plain,tension=1.5,label=$N_{I}(p)$, label.side=right, label.dist=0.1cmv2,i4
fermion,label=$l_{\beta}^{-}(p_{2})$, label.side=left, label.dist=0.1cmv2,o4
plaino1,v1,o3 plain,tension=3,label=$\upnu_{k}(k)$, label.side=right,
label.dist=0.1cmo1,v1 fermion,label=$l_{\alpha}^{+}(p_{1})$, label.side=left,
label.dist=0.1cmo3,v1 photon, tension=2,label=$W^{-}$, label.side=right,
label.dist=0.2cmv2,v1 lab=$\circlearrowleft$v1 lab=$\leftarrow$v2

Diagrams A2. Diagrams for case 2 with an explicitly specified choice of
fermion flow direction (thin arrows).

$\displaystyle M_{1}$ $\displaystyle=$ $\displaystyle
i2\sqrt{2}G_{F}\Theta^{*}_{\alpha I}U_{\beta
i}(\overline{v}_{N}\gamma_{\mu}P_{L}v_{1})(\overline{u}_{\nu}\gamma^{\mu}P_{R}v_{2}),$
$\displaystyle M_{2}$ $\displaystyle=$ $\displaystyle
i2\sqrt{2}G_{F}\Theta_{\beta I}U^{*}_{\alpha
i}(\overline{v}_{N}\gamma_{\mu}P_{R}v_{2})(\overline{u}_{\nu}\gamma^{\mu}P_{L}v_{1}),$

the squared terms are

$\displaystyle|M_{1}|^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}|\Theta^{*}_{\alpha I}|^{2}|U_{\beta i}|^{2}(pp_{2})(p_{1}k),$
$\displaystyle|M_{2}|^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}|\Theta_{\beta i}|^{2}|U_{\beta i}|^{2}(pp_{1})(p_{2}k),$
$\displaystyle M_{1}M_{2}^{\dagger}$ $\displaystyle\sim$ $\displaystyle
Tr[(\hat{p}-M_{I})\gamma_{\mu}P_{L}(\hat{p}_{1}-m_{1})\gamma^{\nu}P_{L}\hat{k}\gamma^{\nu}P_{R}(\hat{p}_{2}-m_{2})\gamma_{\nu}P_{R}]=0,$

and the squared amplitude $|M_{W}|^{2}=|M_{1}+M_{2}|^{2}$ for the case 2 has
the form

$|M_{W}|^{2}=128G_{F}^{2}\left[|\Theta_{\alpha I}|^{2}|U_{\beta
i}|^{2}(pp_{2})(p_{1}k)+|\Theta_{\beta I}|^{2}|U_{\alpha
i}|^{2}(pp_{1})(p_{2}k)\right].$

The decay width for all lepton masses nonzero is given by Eq.21.

Case 3. For the three diagrams in this case, the notation $M_{3_{XY}}$,
$XY={LL,LR,RL,RR}$ is used for the neutral current diagram and the notation
$M_{1},M_{2}$ is used for $W^{+}$ and $W^{-}$ exchange diagrams, respectively.
Then the full amplitude $|M_{\rm WZ}|$ contains 21 terms: 6 terms for the
squared diagrams and $\frac{36-6}{2}=15$ terms for the interferences.
Moreover, we need to change the fermion flow depending on each specific
interference term to build a trace of gamma matrices. For instance, we need to
change the fermion flow for $Z\upnu l$-vertex in the diagram with $W^{+}$ in
order to calculate $M_{1}M_{3_{XY}}^{\dagger}$. Full interference term looks
as

$M_{\rm
interf}^{2}\equiv\sum\limits_{XY}2Re\left[M_{1}{({M_{3}}^{\prime}_{XY})}^{\dagger}\right]+\sum\limits_{XY}2Re\left[M_{2}^{\prime}{({M_{3}}^{\prime\prime}_{XY})}^{\dagger}\right].$

Symbols ”$\prime$” and ”$\prime\prime$” mark those diagrams in which it was
necessary to change the fermion flow for calculation,

$\displaystyle|M_{\rm WZ}|^{2}=|M_{Z}|^{2}+|M_{W}|_{\alpha=\beta}^{2}+M_{\rm
interf}^{2},$

where

$\displaystyle|M_{Z}|^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}\left(\mathcal{C}_{1}^{2}+\mathcal{C}_{2}^{2}\right)~{}|(U\Theta)_{iI}|^{2}~{}\left[(pp_{2})(p_{1}k)+(pp_{1})(p_{2}k)\right]+$
$\displaystyle+128G_{F}^{2}(4\mathcal{C}_{1}\mathcal{C}_{2})~{}|(U^{\dagger}\Theta)_{iI}|^{2}~{}\frac{m_{\alpha}^{2}}{2}(pk),$
$\displaystyle|M_{W}|_{\alpha=\beta}^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}~{}|\Theta_{\alpha I}|^{2}|U_{\alpha
i}|^{2}\cdot\left[(pp_{2})(p_{1}k)+(pp_{1})(p_{2}k)\right],$ $\displaystyle-
M_{\rm interf}^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}\mathcal{C}_{1}\cdot 2\left[(U^{\dagger}\Theta)_{iI}\Theta_{\alpha
I}^{*}U_{\alpha
i}~{}(pp_{2})(p_{1}k)+(U^{\dagger}\Theta)^{*}_{iI}\Theta_{\alpha I}U_{\alpha
i}^{*}~{}(pp_{1})(p_{2}k)\right]+$
$\displaystyle+128G_{F}^{2}\mathcal{C}_{2}\cdot
2\left[(U^{\dagger}\Theta)_{iI}\Theta_{\alpha I}^{*}U_{\alpha
i}+(U^{\dagger}\Theta)^{*}_{iI}\Theta_{\alpha I}U_{\alpha
i}^{*}\right]\cdot\frac{m_{\alpha}^{2}}{2}(pk),$

and we use notation ($s_{W}=\sin\theta_{W}$ where $\theta_{W}$ is the Weinberg
angle)

$\displaystyle\mathcal{C}_{1}=\left(s_{W}^{2}-\dfrac{1}{2}\right),~{}~{}~{}\mathcal{C}_{2}=s_{W}^{2}.$

The square of the full amplitude is simplified after summing by the neutrino
mass states

$\displaystyle\sum\limits_{i}|M_{\rm WZ}|^{2}$ $\displaystyle=$ $\displaystyle
128G_{F}^{2}\left[\left(\mathcal{C}_{1}^{2}+\mathcal{C}_{2}^{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}+\left(1-2\mathcal{C}_{1}\right)|\Theta_{\alpha
I}|^{2}\right]~{}\left[(pp_{2})(p_{1}k)+(pp_{1})(p_{2}k)\right]+$
$\displaystyle+128G_{F}^{2}\left[\left(4\mathcal{C}_{1}\mathcal{C}_{2}\right)\sum\limits_{\beta}|\Theta_{\beta
I}|^{2}-4\mathcal{C}_{2}|\Theta_{\alpha
I}|^{2}\right]~{}\frac{m_{\alpha}^{2}}{2}(pk),$

the decay width for lepton nonzero masses is given by Eq.46. Note that for
Dirac neutrinos, other signs of interference terms above may appear.

## Appendix B. Kinematics of $1\to 3$ decay with mass terms

This Appendix contains the details of integration by the Dalitz plot in the
general case of non-zero masses using invariant variables. For the process
$N_{I}(p)\to a(p_{1}),b(p_{2}),c(p_{3})$ invariant kinematic variables are
defined as $S_{1}=(p-p_{3})^{2}=(p_{1}+p_{2})^{2}$,
$S_{2}=(p-p_{1})^{2}=(p_{2}+p_{3})^{2}$,
$S_{3}=(p-p_{2})^{2}=(p_{3}+p_{1})^{2}$ which are then redefined to
dimensionless variables for convenience, $s_{i}\equiv\frac{S_{i}}{M_{I}^{2}}$,
$r_{i}\equiv\frac{m_{i}^{2}}{M_{I}^{2}}$ where $m_{i}^{2}=p_{i}^{2}$,
$i=1,2,3$ and the width

$\displaystyle\Gamma(N_{I}\to a,b,c)$ $\displaystyle=$
$\displaystyle\frac{1}{256\pi^{3}M_{I}^{3}}\int
dS_{1}dS_{2}|M(S_{1},S_{2})|^{2}~{}\Uptheta\left[-G(S_{1},S_{2},M_{I}^{2},m_{2}^{2},m_{1}^{2},m_{3}^{2})\right]=$
$\displaystyle=$ $\displaystyle\frac{M_{I}^{5}}{256\pi^{3}}\int
ds_{1}ds_{2}|M(s_{1},s_{2})|^{2}~{}\Uptheta\left[-G(s_{1},s_{2},1,r_{2},r_{1},r_{3})\right]$

where $\Uptheta[x]$ is the step function and $G$ is the Gram determinant

$\displaystyle G(x,y,z,a,b,c)$ $\displaystyle=$ $\displaystyle
x^{2}y+xy^{2}-xy(z+a+b+c)-$ $\displaystyle-bc(x+y+z+a)-za(x+y+b+c)+$
$\displaystyle+xzc+xab+yzc+yac+z^{2}a+za^{2}+b^{2}c+bc^{2}.$

The equations for the boundaries of the physical region in terms of invariants
are obtained requiring $G(s_{1},s_{2},1,r_{2},r_{1},r_{3})=0$, so for $s_{1}$
they are

$\displaystyle
s_{1}^{\pm}=1+r_{3}-\frac{1}{2s_{2}}\left((1-r_{1}+s_{2})(s_{2}-r_{2}+r_{3})\mp\lambda^{1/2}(s_{2},1,r_{1})\lambda^{1/2}(s_{2},r_{2},r_{3})\right)$

and for $s_{2}$ they are given by

$(m_{2}+m_{3})^{2}\leq
S_{2}\leq(M_{I}-m_{1})^{2}~{}~{}~{}\text{or}~{}~{}~{}(\sqrt{r_{2}}+\sqrt{r_{3}})^{2}\leq
s_{2}\leq(1-\sqrt{r_{1}})^{2}$

so the width

$\Gamma(N_{I}\to
a,b,c)=\frac{M_{I}^{5}}{256\pi^{3}}\int\limits_{(\sqrt{r_{2}}+\sqrt{r_{3}})^{2}}^{(1-\sqrt{r_{1}})^{2}}ds_{2}\int\limits_{s_{1}^{-}}^{s_{1}^{+}}ds_{1}|M(s_{1},s_{2})|^{2}$

## References

  * [1] R.N. Mohapatra and G. Senjanovic, Neutrino mass and spontaneous parity violation, Phys. Rev. Lett. 44 (1980) 912
  * [2] J. Schechter and J.M.F. Valle, Neutrino masses in SU(2) $\times$ U(1) theories, Phys. Rev. D22 (1980) 2227
  * [3] G. Bellini, L. Ludhova, G. Ranucci and F. Villante, Neutrino oscillations, Adv. High Energy Phys. 2014 (2014) 191960
  * [4] R. Adhikari et al., White paper on keV sterile neutrino dark matter, JCAP 01 (2017) 025
  * [5] T. Asaka, S. Blanchet and M. Shaposhnikov, The $\nu$MSM, dark matter and neutrino masses, Phys. Lett. B631 (2005) 151
  * [6] M. Shaposhnikov, A possible symmetry of the $\nu$MSM, Nucl. Phys. B763 (2007) 49
  * [7] A. Boyarsky, O. Ruchayskiy and Shaposhnikov M., The role of sterile neutrinos in cosmology and astrophysics, Ann. Rev. Nucl. Part. Sci. 59 (2009) 191
  * [8] A. Merle, KeV neutrino model building, Int. J. Mod. Phys. D22 (2013) 1330020
  * [9] A. Boyarsky, O. Ruchayskiy, D. Iakubovskyi and J. Franse, Unidentiﬁed Line in X-Ray Spectra of the Andromeda Galaxy and Perseus Galaxy Cluster, Phys. Rev. Lett. 113 (2014) 251301 (arXiv:1402.4119 [astro-ph.CO]).
  * [10] T. Asaka and M. Shaposhnikov, The nuMSM, dark matter and baryon asymmetry of the universe, Phys. Lett. B620 (2005) 17 (hep-ph/0505013)
  * [11] R.E. Shrock, General Theory of Weak Leptonic and Semileptonic Decays. 1. Leptonic Pseudoscalar Meson Decays, with Associated Tests For, and Bounds on, Neutrino Masses and Lepton Mixing, Phys. Rev. D24 (1981) 1232
  * [12] R.E. Shrock, General Theory of Weak Processes Involving Neutrinos. 2. Pure Leptonic Decays, Phys. Rev. D24 (1981) 1275
  * [13] A. Boyarsky, M. Ovchinnikov, O. Ruchayskiy and V. Syvolap, Improved big bang nucleosynthesis constraints on heavy neutral leptons, Phys.Rev. D104 (2021) 023517 (arXiv:2008.00749 [hep-ph])
  * [14] A. Abada, D. Das, A.M. Teixeira, A. Vicente and C. Weiland, Tree-level lepton universality violation in the presence of sterile neutrinos: impact for $R_{K}$ and $R_{\pi}$, JHEP 1302 (2013) 048 (arXiv:1211.3052 [hep-ph])
  * [15] A. Abada, D. Das, A.M. Teixeira, A. Vicente and C. Weiland, Sterile neutrinos in leptonic and semileptonic decays, JHEP 1402 (2014) 091 (arXiv:1311.2830 [hep-ph])
  * [16] T. Asaka, S. Eijima and K. Takeda, Lepton universality in the $\nu$MSM, Phys. Lett. B742 (2015) 303 [arXiv:1410.0432 [hep-ph]].
  * [17] A. Ibarra, E. Molinaro and S. Petcov, TeV TeV scale see-saw mechanisms of neutrino mass generation, the Majorana nature of the heavy singlet neutrinos and $(\beta\beta)_{0\nu}$ decay, JHEP 09 (2010) 108 (arXiv:1007.2378 [hep-ph])
  * [18] B. Pontecorvo, Inverse beta processes and nonconservation of lepton charge, Sov. Phys. JETP, 7 (1958) 172
Z. Maki, M. Nakagawa and S. Sakata, Remarks on the unified model of elementary
particles, Prog. Theor. Phys. 28 (1962) 870

  * [19] M.N. Dubinin, E.Yu. Fedotova and D.M. Kazarkin, Cosmological constraints for three generations see-saw mechanism of keV sterile neutrinos, Phys. Part. Nucl. Lett. 20 (2023), to appear
  * [20] J. Casas and A. Ibarra, Oscillating neutrinos and $\mu\to e,\gamma$, Nucl.Phys. B618 (2001) 171
  * [21] J. Ellis, M. E. Gomez, G. K. Leontaris, S. Lola and D. V. Nanopoulos, Charged lepton flavor violation in the light of the Super-Kamiokande data, Eur. Phys. J. C14 (2000) 319 (arXiv:hep-ph/9911459) R. Barbieri, L. Hall and A. Strumia, Violations of lepton flavor and CP in supersymmetric unified theories, Nucl. Phys. B445 (1995) 219 (arXiv:hep-ph/9501334) S. Bilenky, S. Petcov and B. Pontecorvo, Lepton mixing, $\mu\to e\gamma$ decay and neutrino oscillations, Phys. Lett. B67 (1977) 309
  * [22] M.N. Dubinin and D.M. Kazarkin, Improved cosmological bounds for a fine-tuned see-saw mechanism of keV sterile neutrinos, arXiv:2206.05186 [hep-ph]
  * [23] A. Ibarra and G.G. Ross, Neutrino phenomenology: the case of two right-handed neutrinos, Phys. Lett. B591 (2004) 285 (arXiv:hep-ph/0312138)
  * [24] W. Buchmuller, P. Di Bari and M. Plumacher, The neutrino mass window for baryogenesis, Nucl. Phys. B665 (2003) 445 (arXiv:hep-ph/0302092)
  * [25] T. Asaka, S. Eijima and H. Ishida, Mixing of Active and Sterile Neutrinos, JHEP 04 (2011) 011 (arXiv:1101.1382 [hep-ph])
  * [26] J. Kersten and A. Smirnov, Right-handed neutrinos at LHC and the mechanism of neutrino mass generation, Phys. Rev. D76 (2007) 073005 (arXiv:0705.3221 [hep-ph])
  * [27] A. Ibarra, E. Molinaro and S. Petcov, Low energy signatures of the TeV scale see-saw mechanism, Phys. Rev. D84 (2011) 013005 (arXiv:1103.6217 [hep-ph])
  * [28] S. Alekhin et al, A facility to Search for Hidden Particles at the CERN SPS: the SHiP physics case , Rept. Prog. Phys. 79 (2016) 124201 (arXiv:1504.04855 [hep-ph])
  * [29] K. Bondarenko, A. Boyarsky, D. Gorbunov and O. Ruchayskiy, Phenomenology of GeV-scale heavy neutral leptons, JHEP 11 (2018) 032 (arXiv:1805.08567 [hep-ph]) D. Gorbunov and M. Shaposhnikov, How to find neutral leptons of the $\nu$MSM, JHEP 10 (2007) 015 (arXiv:0705.1729 [hep-ph]).
  * [30] J. L. Rosner, S. Stone and R. S. Van de Water, Leptonic decays of charged pseudoscalar mesons - 2015, arXiv:1509.02220 [hep-ph]
  * [31] A. Denner, H. Eck, O. Hahn, and J. Kublbeck, Feynman rules for fermion number violating interactions, Nucl. Phys. B387 (1992) 467
  * [32] Eduardo Cortina Gil et al (NA62 Collaboration), Search for heavy neutral lepton production in $K^{+}$ decays, Phys.Lett. B778 (2018) 137 (arXiv:1712.00297 [hep-ex]) A.V. Artamonov et al (E949 Collaboration), Search for heavy neutrinos in $K^{+}\to\mu^{+}\nu$ decay, Phys.Rev. D91 (2015) 052001, erratum: Phys.Rev. D91 (2015) 059903, (arXiv: 1411.3963 [hep-ex]) A.T. Shaikhiev and Yu.G. Kudenko, New results on sterile neutrinos searches , Phys.Atom.Nucl. 78 (2015) 1591 G. Bernardi et al, Further limits on heavy neutrino couplings, Phys. Lett. B203 (1988) 332 G. Bernardi et al, Search for neutrino decay, Phys. Lett. B166 (1986) 479
  * [33] M. Gronau, C. Leung and J. Rosner, Extending Limits on Neutral Heavy Leptons, Phys. Rev. D29 (1984) 2539
  * [34] A.D. Dolgov, S.H. Hansen, G. Raffelt and D. Semikoz, Cosmological and astrophysical bounds on a heavy sterile neutrino and the KARMEN anomaly, Nucl. Phys. B580 (2000) 331 (arXiv:hep-ph/0002223) A.D. Dolgov, S.H. Hansen, G. Raffelt and D. Semikoz, Heavy sterile neutrinos: Bounds from big bang nucleosynthesis and SN1987A, Nucl. Phys. B590 (2000) 562 (arXiv:hep-ph/0008138)
  * [35] V. Cirigliano and I. Rosell, Two-loop effective theory analysis of $\pi(K)\to e\nu_{e}(\gamma)$ branching ratios, Phys. Rev. Lett. 99 (2007) 231801 (arXiv:0707.3439 [hep-ph])
  * [36] R.L. Workman et al., Particle Data Group, Prog. Theor. Exp. Phys. 2022 (2022) 083C01
  * [37] M. Finkmeier, Radiative corrections to $\pi$(l2) and K(l2) decays, Phys. Lett. B387 (1996) 391 (hep-ph/9505434)
  * [38] H. Haber and G. Kane, The search for supersymmetry: probing physics beyond the Standard Model, Phys. Repts. 117 (1985) 75