Radiative Dirac Neutrino Mass with Dark Matter and it's implication to $0\nu 4\beta$ in the $U(1)_{B-L}$ extension of the Standard Model

TL;DR

This paper proposes a $U(1)_{B-L}$ extended Standard Model with residual $Z_4$ symmetry that explains Dirac neutrino masses, stabilizes dark matter, and predicts an observable neutrinoless quadruple beta decay ($0 u 4eta$).

Contribution

It introduces a novel anomaly-free $U(1)_{B-L}$ extension with residual $Z_4$ symmetry that links Dirac neutrino mass generation, dark matter stability, and predicts $0 u 4eta$ decay.

Findings

Absence of neutrinoless double beta decay ($0 u 2eta$).

Predicted enhancement of neutrinoless quadruple beta decay ($0 u 4eta$).

Compatibility with dark matter detection and collider constraints.

Abstract

The Standard Model gauge symmetry is extended by $U(1)_{B-L}$ which when spontaneously broken leads to residual $\mathbb{Z}_4$ symmetry. $U(1)_{B-L}$ gauge symmetry made anomaly free by introducing exotic SM singlets with corresponding $U(1)_{B-L}$ charges of $13$, $-14$, and $15$. $\mathbb{Z}_4$ symmetry ensures the Dirac nature of neutrinos, simultaneously stabilizing dark matter. Dirac neutrino mass is generated through scotogenic scenario. Dark matter, direct detection, cosmological constraints, and collider constraints analysis is performed. $\mathbb{Z}_4$ symmetry predicts the exact absence of neutrinoless double beta decay ($0\nu 2\beta$) and gives a prediction for an enhanced neutrinoless quadruple beta decay ($0\nu 4\beta$) via which this model can be tested. Model allows for Majorana dark matter as well as for long-lived dark matter candidates.

Tables

Table 1. Table 1 : Model particle content. All fields are given as Left-handed chiral fields. Bold-faced numbers represent non-Abelian group irreducible representations.

Field	SU(3)c	SU(2)L	U(1)Y	U(1)B-L	Flavor
Q	3	2	$\frac{1}{6}$	$\frac{1}{3}$	3
uc	$\overline{\mathrm{𝟑}}$	1	$-\frac{2}{3}$	$-\frac{1}{3}$	3
dc	$\overline{\mathrm{𝟑}}$	1	$\frac{1}{3}$	$-\frac{1}{3}$	3
L	1	2	$-\frac{1}{2}$	$-1$	3
ec	1	1	$1$	1	3
${\nu }^c$	1	1	$0$	5	3
N	1	1	$0$	$-6$	3
Nc	1	1	$0$	2	3
${\Psi }_{\text{1}}$	1	1	$0$	$13$	3
${\Psi }_{\text{2}}$	1	1	$0$	$-14$	$2\times 3$
${\Psi }_{\text{3}}$	1	1	$0$	$15$	3
H	1	2	$\frac{1}{2}$	$0$	$1$
$\eta$	1	2	$\frac{1}{2}$	$-1$	$1$
$\chi$	1	1	$0$	$1$	$1$
S	1	1	$0$	$2$	$1$
S4	1	1	$0$	$4$	$1$

Table 2. Table 2 : Global U ​ ( 1 ) ( B − L ) 𝑈 subscript 1 𝐵 𝐿 U(1)_{(B-L)} and ℤ 4 subscript ℤ 4 \mathbb{Z}_{4} transformations of fields. a = 1 , 2 , 3 𝑎 1 2 3 a=1,2,3 is the flavor index and i = 1 , 2 𝑖 1 2 i=1,2 .

${\mathbb{Z}}_4$	Fields
$𝟙$	$H$, $S_4$
$w^{\ast }$	${\nu }^a$, $e^a$, ${\Psi }_{1L}^{a\ast }$, ${\Psi }_{3R}^a$, ${\eta }^0$, ${\eta }^{+}$, ${\chi }^{\ast }$
$w^2$	$N_{L,R}^a$, ${\Psi }_{i2L}^a$, $S$
Global $U{\left(1\right)}_{\left(B-L\right)}$	Fields
$\frac{1}{3}$	$u^a$, $d^a$

Table 3. Table 4: Some relevant phenomenological bounds.

Limit	Condition	Reference
	$\forall M_{Z^{\prime }}$	(55)
$M_{Z^{\prime }}>4.2$TeV	$L=36.1$fb-1;$\sqrt{s}=13$TeV	(56)
$M_{Z^{\prime }}>4.4$TeV	$\sigma \times Br\left(pp\to Z^{\prime }\to e\mu \right)\wedge Br\left(Z^{\prime }\to e\mu \right)=0.10$, Model independent	CMS CR-2018-371 (57)
$g_{BL}$<0.236(0.127)	$Q_{BL}^{Max}=15\left(28\right)$	Perturbativity bound
		(58)
	$m_{A^{\prime }}\sim O\left(10-100\text{MeV}/c^2\right)$	(59)
		(60; 61)
	$m_{A^{\prime }}\approx O\left(\text{GeV}\right)$	(62; 63)

Equations

$$\frac{f_{ij}}{2\Lambda }\left({\nu }_i{\varphi }^0-l_i{\varphi }^{+}\right)\left({\nu }_j{\varphi }^0-l_j{\varphi }^{+}\right)+\text{h.c.},$$

$${\mathcal{L}}_{\text{Yuk}}$$

$$-{\bar{N}}_{aR}{Y}_{ND}^{ab}{N}_{bL}{S}_4-N_{aR}{Y}_{NM}^{ab}{N}_{bR}{S}_4$$

$$-{\Psi }_{1aL}{Y}_{12}^{ab}{\Psi }_{2bL}\chi -{\Psi }_{2aL}{Y}_{23}^{ab}{\Psi }_{3bL}{\chi }^{\ast }-{\Psi }_{1aL}{Y}_{13}^{ab}{\Psi }_{3bL}\frac{{\left(S_4^{\ast }\right)}^7}{{\Lambda }^6}-{\Psi }_{2aL}{Y}_{22}^{ab}{\Psi }_{2bL}\frac{S_4^7}{{\Lambda }^6}+\text{H.c.}$$

$${\sum Q_{B-L}|}_{\text{Grav.}}$$

$${\mathcal{L}}_{Gauge}$$

$${\epsilon }_f$$

$${\epsilon }_s$$

$$F:$$

$$S:$$

$$\text{Total}:$$

$${\epsilon }_0$$

$${\epsilon }^{\prime }$$

$$m_{\nu }^{ab}$$

$$F(x)$$

$${\mathcal{L}}_N$$

$$m_{N_{1,2}}$$

$$\left(\begin{array}{c}{N}_1\\ N_2\end{array}\right)$$

$$\text{tan}(2{\theta }_N)$$

$${\bar{N}}_R{Y}_{R2}{\Psi }_2\frac{S_4^3}{{\Lambda }^2}+N_L{Y}_{L2}{\Psi }_2\frac{S_4^5}{{\Lambda }^4}+\text{H.c.},$$

$${\bar{\Psi }}_3{Y}_{3R}{\nu }_R\frac{S_4^5}{{\Lambda }^4}+{\Psi }_1{Y}_{1L}LH\frac{S_4^{\star 3}}{{\Lambda }^3}+\text{H.c.}$$

$$V\left(H,\eta ,\chi ,S,S_4\right)=$$

$$+m_{\chi }^2{\chi }^{\ast }\chi +\frac{1}{2}{\lambda }_{\chi }{\left({\chi }^{\ast }\chi \right)}^2+m_S^2{S}^{\ast }S+\frac{1}{2}{\lambda }_S{\left(S^{\ast }S\right)}^2$$

$$-m_{S_4}^2{S}_4^{\ast }{S}_4+\frac{1}{2}{\lambda }_{S_4}{\left(S_4^{\ast }{S}_4\right)}^2+{\lambda }_{H\eta }\left(H^{†}H\right)\left({\eta }^{†}\eta \right)$$

$$+{\lambda }_{H\eta }^{\prime }\left(H^{†}\eta \right)\left({\eta }^{†}H\right)+{\lambda }_{H\chi }\left(H^{†}H\right)\left({\chi }^{\ast }\chi \right)$$

$$+{\lambda }_{HS}\left(H^{†}H\right)\left(S^{\ast }S\right)+{\lambda }_{H{S}_4}\left(H^{†}H\right)\left(S_4^{\ast }{S}_4\right)$$

$$+{\lambda }_{\eta \chi }\left({\eta }^{†}\eta \right)\left({\chi }^{\ast }\chi \right)+{\lambda }_{\eta S}\left({\eta }^{†}\eta \right)\left(S^{\ast }S\right)$$

$$+{\lambda }_{\eta S_4}\left({\eta }^{†}\eta \right)\left(S_4^{\ast }{S}_4\right)+{\lambda }_{\chi S}\left({\chi }^{\ast }\chi \right)\left(S^{\ast }S\right)+{\lambda }_{\chi S_4}\left({\chi }^{\ast }\chi \right)\left(S_4^{\ast }{S}_4\right)$$

$$+{\lambda }_{S{S}_4}\left(S^{\ast }S\right)\left(S_4^{\ast }{S}_4\right)+[{\mu }_H\left(H^{†}\eta \right)\chi +{\mu }_S{\chi }^2{S}^{\ast }+{\mu }_{S_4}{S}^2{S}_4^{\ast }$$

$$+{\lambda }_1\left(H^{†}\eta \right){\chi }^{\ast }S+{\lambda }_2{\chi }^{2}S{S}_4^{\ast }+\text{H.c.}]$$

$${\lambda }_H{v}^2+{\lambda }_{H{S}_4}{v}_4^2-2m_H^2=0,$$

$${\lambda }_{S_4}{v}_4^2+{\lambda }_{H{S}_4}{v}^2-2m_{S_4}^2=0.$$

$$m^2$$

$$m_{1,2}^2$$

$$\text{tan}\left(2{\theta }_{H^0{S}_4}\right)$$

$$\text{with }\left(\begin{array}{c}{s}_h\\ s_H\end{array}\right)$$

$$m^2$$

$$m_{1,2}^2$$

$$\pm \frac{1}{2}\sqrt{{\left[m_{\eta }^2-m_{\chi }^2+\left({\lambda }_{H\eta }+{\lambda }_{H\eta 2}-{\lambda }_{H\chi }\right)\frac{v^2}{2}+\left({\lambda }_{\eta S_4}-{\lambda }_{\chi S_4}\right)\frac{v_4^2}{2}\right]}^2+4{\mu }_H^2\frac{v^2}{2}},$$

$$\text{tan}\left(2{\theta }_{\xi R/I}\right)$$

Full text

scotogenic, Dirac, neutrino mass, neutrinoless quadruple beta decay, dark matter, B$-$L

††preprint: arXiv:1903.xxxx††preprint: Prepared for submission to Phys. Rev. D

1. I Introduction
2. II Model
3. III Neutrino masses
4. IV Fermion sector
5. V Scalar sector
6. VI Dark Matter
7. VII Neutrinoless Quadruple Beta Decay
8. VIII Results and discussion
9. IX Conclusion
10. A Note on removing phases from Yukawa terms in $N_{L,R}$ sector
11. B $C_{0}(s_{1},s_{12},s_{2};m_{0},m_{1},m_{2})$
12. C N annihilation diagrams

I Introduction

The Standard Model (SM) of strong and electroweak interactions has proven to be very successful so far with the last remaining piece experimentally discovered on July, 4’th 2012 (1; 2).Nevertheless there are experimental observations that require new physics beyond Standard Model (BSM). One of these problems is the experimental observation of neutrino oscillations (3; 4; 5; 6; 7; 8; 9; 10; 11; 12) back in 1990’s. Theoretical explanation of neutrino masses requires addition of new particles BSM. The most minimalistic and simplest realizations of this are the seesaw mechanism of type I (13; 14; 15; 16) which adds a fermion singlet to SM. Next would be seesaw of type II (17; 18; 19; 20) which extends the SM by a scalar triplet. Last of this kind of realizations is seesaw of type III in which SM is extended by a fermionic electroweak triplet. All these tree level realizations of naturally small neutrino masses require either a small couplings or heavy new physics in order to explain the smallness of neutrino masses and they are lead to unique dimension-five effective operator

[TABLE]

known as Weinberg operator (21). In order to avoid the requirement of heavy new physics or small couplings, for instance neutrino masses can be generated radiatively at one-loop order. Examples of this realizations include (22) the Zee model from 1980, the canonical scotogenic model (23) (scotos from Greek meaning darkness) from 2006, and radiative inverse seesaw model (24). Since neutrinos are neutral and colorless they can be of Dirac or Majorana type. Currently there is no experimental evidence toward any direction. But if neutrinos are Dirac in nature there must be a symmetry (conserved quantity) responsible for the absence of Majorana mass of neutrinos. This issue was systematically studied in (25). The symmetry for the Dirac nature of neutrinos can be the lepton number already present in SM as an accidental symmetry. Experiments such as COURE, GERDA $0\nu 2\beta$, NEMO 3, The MAJORANA Neutrinoless Double-beta Decay Experiment, etc. looking for neutrinoless double and quadruple beta decay can solve this problem in the near future. On the other hand, if $(0\nu 2\beta)$ experiment sees no positive results this could hint in the direction of Dirac neutrinos. But there exist even more exotic scenarios like was explained in (26). In the case of absence of positive results from neutrinoless double beta decay and confirmation of neutrinoless quadruple beta decay, one needs to find a theoretical explanation for this kind of experimental observation.

In our work we present a UV complete model where lepton number is gauged in $U(1)_{B-L}$ symmetry. Spontaneously breaking $U(1)_{B-L}$ to residual $\mathbb{Z}_{4}$ discrete symmetry allows for Dirac neutrinos which obtain their masses radiatively via scotogenic scenario. Model naturally predicts neutrinoless quadruple beta decay whereas neutrinoless double beta decay is exactly absent. Furthermore, model allows for stable dark matter, fermionic or bosonic, and leptogenesis. Similar works done on $U(1)_{B-L}$ extension of SM are (27; 28; 29; 30; 31; 32; 33; 34; 35; 36).

The paper is organized as follows: in Sec. II the model is introduced and the cancellation of chiral anomalies is explained; Sec. III demonstrates how radiative Dirac neutrino masses are generated; Secs. IV and V give the fermion and scalar mass spectrum, respectively; Sec. VI discusses dark matter candidates; in Sec. VII we go over the neutrinoless quadruple beta decay prediction; Sec. VIII presents the results and discusses relevant constraints for our model; and Sec. IX concludes.

II Model

SM gauge symmetry is extended to $SU(3)_{c}\times SU(2)_{L}\times U(1)_{Y}\times U(1)_{B-L}$ and the field content of the model is shown in Tab. 1. All fields are given as left-chiral fields and "c" in superscript denotes the charge conjugation. The Yukawa part of the Lagrangian of the model is given below

[TABLE]

The model is constructed as follows: $\nu_{R}$ is introduced as the Dirac partner for Left-handed neutrinos, $N_{L,R}$ fermions, $\eta$, and $\chi$ scalars are introduced to complete the loop for radiative neutrino mass generation, i.e. scotogenic scenario, $\Psi_{i}$ are introduced for anomaly cancellation, and lastly $S_{4}$ is needed for spontaneous symmetry breaking (SSB) of $U(1)_{B-L}$ to residual $\mathbb{Z}_{4}$ discrete symmetry in the leptonic sector. Here, the residual $\mathbb{Z}_{4}$ symmetry is given by $e^{\imath(B-L)2\pi/4}=w^{(B-L)}$, where $w=e^{\imath 2\pi/4}$ with $w^{4}=1$.

$H$ serves the role of Standard Model (SM) Higgs field that couples and gives masses to SM quarks and charged leptons. $S_{4}$ is introduced to spontaneously break $U(1)_{B-L}$ gauge symmetry. Since $SU(2)_{L}$ (same for $SU(3)_{c}$) group is orthogonal to $U(1)_{B-L}$, $SU(2)_{L}$ irreducible representation’s (irrep’s) components must transform identically under $U(1)_{B-L}$. Therefore, for quarks, i.e. $u_{L,R}$, $d_{L,R}$, $U(1)_{(B-L)}$ gauge symmetry is broken to global $U(1)_{(B-L)}$ symmetry. Whereas for $\eta^{+},\eta^{0},e_{L,R}$, and $\nu_{L,R}$, $U(1)_{(B-L)}$ gauge symmetry is broken to $\mathbb{Z}_{4}$ symmetry under which they transform as $w^{*}$. All field transformations under residual $\mathbb{Z}_{4}$ and $U(1)_{B-L}$ global symmetries are summarized in Tab. 2. Fermions that transform as $w^{2}$ are of Majorana type, i.e. $N_{L,R}$ and $\Psi_{2L}$ in this case and fermions that transform as $w$ and $w^{*}$ under $\mathbb{Z}_{4}$ residual symmetry arrange themselves into Dirac pairs, $\nu_{L,R},\Psi_{1,3L}$ are among them. When electroweak symmetry is broken by $H$ the $\eta^{0}$ and $\chi^{\star}$ scalars mix through $H$, which is needed for the neutrino mass generation. Furthermore, when $U(1)_{B-L}$ symmetry is broken by $S_{4}$ vacuum expectation value (VEV), the mass eigenstates of $(\eta^{0},\chi^{\star})$(call them $\xi_{1,2}$) obtain an effective operator $\xi_{i}^{4}$+H.c., which is invariant under $Z_{4}$ and generates the neutrinoless quadruple beta decay. Neutrinoless double beta decay is forbidden by $\mathbb{Z}_{4}$ symmetry, therefore neutrinoless quadruple beta decay will be dominant. More on this in Secs. VII and V.

Chiral anomalies

Model is chiral anomaly free and cancellation of anomalies per family is shown in Tab. 3. $U(1)_{Y}$ gravitational anomaly is cancelled like in SM and $U(1)_{B-L}$ gravitational anomaly is cancelled as follows:

[TABLE]

Abelian kinetic mixing

Model Lagrangian must be augmented with renormalizable Abelian kinetic mixing(KM) $\epsilon_{0}$ counter-term, since the one-loop corrections has singular contribution, same as in (38)

[TABLE]

where $B$ and $B^{\prime}$ are strength tensors for hypercharge($U(1)_{Y}$) and $B-L$($U(1)_{B-L}$) gauge groups. respectively. $\varepsilon_{0}$ represents the bare Abelian kinetic mixing counter-term which must be included to renormalize the divergent one-loop corrections (39; 40)

[TABLE]

Fermions that contribute Abelian KM are $Q,u^{c},d^{c},L,e^{c}$ and scalar that contributes is $\eta$. Their contributions to the divergence are given as

[TABLE]

In order to regularize the divergence $\varepsilon_{0}$ must be given by

[TABLE]

where the tree level finite piece is denoted by $\varepsilon_{0}^{\text{finite}}$ and the one-loop corrected finite contribution to Abelian KM is given by

[TABLE]

III Neutrino masses

Neutrino tree level mass is forbidden by $U(1)_{B-L}$ symmetry. This is the $S$ symmetry from Ref. (25) and the neutrino mass is generated via first scenario of one-loop radiative case from Ref. (25). Neutrino masses are obtained via a diagram shown in Fig. 1. Neutrinos transform as $w^{*}$ under residual $\mathbb{Z}_{4}$ symmetry, therefore $\mathbb{Z}_{4}$ guarantees the Dirac nature of neutrinos in our model. Interesting feature of this model is that Dirac neutrino masses are generated through the Majorana dark sector $N_{1,2}$ fermions which transform as $w^{2}$ under $\mathbb{Z}_{4}$ symmetry and are allowed to have Majorana masses. Other interesting property is that the $\mathbb{Z}_{4}$ residual symmetry which originated from gauged $U(1)_{(B-L)}$ symmetry is responsible both for Diracness of the neutrinos as well as for the stability of dark matter in our model. It is actually the $\mathbb{Z}_{4}$ plus the Lorentz symmetry that stabilizes the dark matter.

Neutrino radiative mass is given as

[TABLE]

where $F(x)$ is defined as

[TABLE]

and mixing angles of $\left(\eta,\chi\right)$ and $\left(N_{L},N_{R}^{c}\right)$ are given in Eqs. 29 and 17, respectively.

IV Fermion sector

The SM fermions generate their masses in a usual way. Since neutrinos transform as $w^{*}$ under residual $\mathbb{Z}_{4}$ symmetry, their masses are of Dirac type and were given in Sec. III. $N_{L,R}$ transform as $w^{2}$ under $\mathbb{Z}_{4}$ therefore they obtain Majorana masses through the seesaw-I texture matrix form

[TABLE]

In general Yukawas here can be complex but the Majorana phases of $N_{L,R}$ can be used to remove this phases, so they are not physical. If the $(1,1)$ component of the $\left(\bar{N}_{L}^{c},\bar{N}_{R}\right)$ mass matrix was non-zero this would not be the case. See App. A for more details on this. Eigenvalues and eigenvectors are given by

[TABLE]

In order for $\Psi_{1,2,3}$ to get their masses, for instance, a SM singlet scalar with $U(1)_{B-L}\sim 28$ (say $S_{28}$) can be introduced. When $S_{28}$ obtains a non-zero VEV $\left\langle(S_{28}\right\rangle\sim v_{28}$ which is allowed by $\mathbb{Z}_{4}$ symmetry, $\Psi_{i}$ masses would be generated through Lagrangian terms $\Psi_{1L}^{a}Y_{13}^{ab}\Psi_{3L}^{b}S_{28}^{*}$ and $\Psi_{2L}^{a}Y_{22}^{ab}\Psi_{2L}^{b}S_{28}$. But in our case, we generate $\Psi_{1,2,3}$ masses through effective dimension-ten Lagrangian terms given in Eq. II. $\Lambda$ scale is associated with mass scale involved in generation of $m_{\Psi_{i}}$. Since $\Psi_{1}^{*},\Psi_{3}\sim w$ and $\Psi_{2i}\sim w^{2}$ under $\mathbb{Z}_{4}$, $\Psi_{1}$ and $\Psi_{3}$ form a Dirac fermion and $\Psi_{2i}$ form two Majorana fermions. Their masses are given by $m_{\Psi_{1,3}}^{ab}=Y_{13}^{ab}\frac{\left\langle S_{8}\right\rangle^{7}}{\Lambda^{6}}$ and $m_{\Psi_{22}}^{\alpha\beta}=Y_{22}^{\alpha\beta}\frac{\left\langle S_{8}\right\rangle^{7}}{\Lambda^{6}}$ with $a,b=1-3$ and $\alpha,\beta=1-6$, respectively. For the number estimates: if $Y_{13,22}\sim O(1)$ and $m_{13,22}\sim O(1\text{TeV})$ we get $\log(\Lambda)\approx(7\log v_{4}-3)/6$ relation between $v_{4}$ and $\Lambda$ scales. If $v_{4}\approx 10^{5}(10^{11})$ GeV then $\Lambda\approx 2\cdot 10^{5}(10^{12})$ GeV, respectively.

Remark regarding the mixing of $N_{L,R}$ and $\nu_{L,R}$ with $\Psi_{2i}$ and $\Psi_{1,3}$, respectively. Since there is no symmetry distinguishing $N_{R,L}$ from $\Psi_{2}$ and similarly for $\nu_{L,R}$ and $\Psi_{1,3}$, these will mix via dimension-6, 7, and 8 effective operators given by

[TABLE]

On the other hand if $S_{28}$ is included, as was explained above, then $\Psi_{1,2,3}$ do generate their masses. But since scalars with charge $\sim 20$ and $\sim 12$ under $U(1)_{B-L}$ are not included, there will be no mixing between $\Psi$ sector and $(N,\nu)$ sector. Which means there is an inherent $\mathbb{Z}_{2}$ symmetry induced under which $\Psi$ sector is odd and all other particles are even(trivial). If this is the case, then lightest of $\Psi$ eigenstates and $\mathbb{Z}_{4}$ LSP can be DM candidates which would give a multi-component DM scenario.

V Scalar sector

Most general scalar potential is given as

[TABLE]

Potential minimization conditions are

[TABLE]

Due to $\mathbb{Z}_{4}$ residual symmetry mass eigenstates can be divided into three groups: trivially transforming under $\mathbb{Z}_{4}$ (singlet representation), transforming as $w$ or $w^{*}$ under $\mathbb{Z}_{4}$ (complex irreducible representation), and transforming as $w^{2}$ under $\mathbb{Z}_{4}$ (real irreducible representation). Scalars transforming trivially under $\mathbb{Z}_{4}$ are those that obtain none-zero VEVs $S_{4}^{0}$, $H^{0}$ and $H^{0}$’s charged multiplet partner $H^{\pm}$. Their mass matrices are given as

[TABLE]

for the $\left(Re[H^{0}],Re[S_{4}]\right)$ basis. Here $\left\langle H^{0}\right\rangle=\frac{v}{\sqrt{2}}$ and $\left\langle S_{4}\right\rangle=\frac{v_{4}}{\sqrt{2}}$. Corresponding mass eigenvalues are

[TABLE]

and mixing angle is given by

[TABLE]

$Im[H^{0}]$ and $Im[S_{4}]$ correspond to would-be Nambu-Goldstone bosons of $Z$ Standard Model gauge boson corresponding to weak neutral current and $Z^{\prime}$ corresponding to spontaneously broken $U(1)_{B-L}$, hence they get eaten-up and have zero mass matrix. Similarly for $H^{\pm}$, it is a would-be Nambu-Goldstone boson and corresponds to SM $W^{\pm}$. Here $H^{\pm}$ does not mix with $\eta^{\pm}$ due to $\mathbb{Z}_{4}$ residual discrete symmetry, the former transforms trivially and the later transforms as $w^{*}$ under $\mathbb{Z}_{4}$ residual discrete symmetry. The mass of $\eta^{\pm}$ is given as $m^{2}_{\eta^{\pm}}=m_{\eta}^{2}+\lambda_{H\eta}v^{2}/2+\lambda_{\eta S_{4}}v_{4}^{2}/2$. Scalars transforming as $w$ or $w^{*}$ under $\mathbb{Z}_{4}$ also mix, their corresponding mass matrix is given by

[TABLE]

with plus sign corresponding to the $\left(Re[\eta],Re[\chi]\right)$ basis and minus sign corresponding to the $\left(Im[\eta],Im[\chi]\right)$ basis. Corresponding mass eigenvalues are the same for both scalar and pseudo-scalar parts, since they transform as complex $w$ representation under $\mathbb{Z}_{4}$, and are given as

[TABLE]

and mixing angles are given by

[TABLE]

Lastly, scalars transforming as $w^{2}$ under $\mathbb{Z}_{4}$ consist only of $S$. Corresponding scalar and pseuso-scalar mass eigenvalues are

[TABLE]

Here the mass splitting is due to $\mu_{S_{4}}$ term in Eq. V which is allowed by $\mathbb{Z}_{4}$ residual symmetry, since $S\sim w^{2}$.

VI Dark Matter

Dark matter is stabilized by same $\mathbb{Z}_{4}$ symmetry that ensures the Dirac nature of the neutrinos. Since our model has beyond Standard Model (BSM) fields that transform as $w$ and $w^{2}$ under $\mathbb{Z}_{4}$, Dirac as well as Majorana type DM candidates are possible. For Dirac type DM, either $\left(\Psi_{1}\Psi_{3}\right)$ Dirac fermion or lighter of $\xi_{1,3}$ mass eigenstates is possible. Whereas if DM of Majorana type, lighter mass eigenstate of $N_{L,R}$ or lighter mass eigenstate of $\Psi_{2i}$ is a viable candidate.

Now in order to calculate the relic abundance of the particle dark matter which was in thermal equilibrium, we would need to calculate the Boltzmann equation

[TABLE]

where $Y_{DM}=n_{DM}/s$, $n_{DM}$ is the number density of the dark matter and $s$ is the entropy density. $H$ is the Hubble expansion, $z=M_{DM}/T$ where $T$ is the background temperature and $\langle\sigma v\rangle$ is the thermally averaged cross-section of the dark matter annihilation process given as

[TABLE]

We can write the partial wave expansion $\langle\sigma v\rangle=a+bv^{2}$. Now, the solution of the above Boltzmann equation in terms of this expansion can be given as

[TABLE]

where $M_{pl}=2.4\times 10^{18}$ GeV and $g_{*}$ is the number of relativistic degrees of freedom at the time of freeze-out. The freeze-out tempertaure can be calculated by the following expression

[TABLE]

which in turn derived from the equality condition of rate of expansion of the Universe $H\approx g^{1/2}_{*}T^{2}/M_{pl}$.

Now, since in our case we have additional particles with mass differences close the dark matter, then they can be thermally accessible during the freeze-out. This will eventually give rise to many additional channels through which the dark matter can co-annihilate and give Standard Model (SM) particles in the final states. The effective cross-section in this case would be as follows

[TABLE]

where $x_{f}=M_{DM}/T$, $\Delta_{i}=M_{i}/M_{DM}-1$ and

[TABLE]

And the thermally averaged cross-section is given as

[TABLE]

One remarkable thing here is that the symmetry that stabilizes DM is the same symmetry that makes neutrinos of Dirac type. The consequence of this is that neutrinos transform non-trivially under DM symmetry, $\mathbb{Z}_{4}$ in this case. Therefore, any field that transforms as $w^{2}$ under $\mathbb{Z}_{4}$ and is in tensor irrep of Poincare symmetry will always decay to pair of neutrinos. On the other hand, fields that transform as $w^{2}$ and are in spinor irrep of Poincare symmetry will not be able to decay to only neutrinos, therefore the lightest can be DM candidate.

$\Psi_{i}$ will not be considered for DM candidate since, as can be seen from eq. II, they do not participate directly in neutrino mass and $0\nu 4\beta$  generations and will not lead to interesting phenomenology. Main candidates to consider are $\xi$, $N$, $S$. $\xi$ has a mixing with the neutral component of the $\eta$ doublet, therefore it will have a direct detection channel mediated by $Z$ SM gauge boson and is severely constrained (41; 42; 43). $N$ for which $m_{N}\approx-v_{4}\frac{Y_{D}^{2}}{Y_{M}}$ is the best DM candidate, since this $N$ is naturally LSP as required by the smallness of neutrino mass and enhancement of $0\nu 4\beta$. The only neutral $\mathbb{Z}_{4}$ non-trivial particle that is lighter than $N$ is neutrino, but $N$ decay to neutrinos is forbidden by $\mathbb{Z}_{4}$ and Poincare symmetry. Decay to the other $\mathbb{Z}_{4}$ non-trivial particles is forbidden by $U(1)_{em}\times\mathbb{Z}_{4}$. The annihilation channels for $N$ as a DM candiate are shown in fig.17. And since the dominant channel will be near resonance i.e ($m_{Z^{\prime}}=2m_{N}$), we have imposed the resonance condition while doing the analysis. The allowed parameter region to satisfy the relic is shown in fig. 2.

From the plot in Fig. 2 we infer that in order for $N$ to be a plausible dark matter candidate the mass of the lightest $N$ has to be between $2.2\sim 7.8$ TeV and the coupling $g_{B-L}$ to be between $0.1\sim 0.02$, Fig. 3. For the analysis we have implemented the model into SARAH 4 (44) and then we took the output to SPheno 3.1 (45) to calculate the mass spectrum. Finally for the dark matter analysis we used MicrOmega 4.3 (46), using the mass spectrum from SPheno 3.1.

Now we focus on $S$ being DM candidate. For $S$ to be a viable DM candidate we assume the following particle mass hierarchy: $m_{\Psi_{i}},m_{N_{j}},m_{\xi_{k}},m_{s4},m_{\eta^{+}}>m_{s}>m_{h}>m_{w},m_{z}>m_{e,u,d}>m_{\nu}$. Since $S$ is a neutral scalar boson that transforms as $w^{2}$ under residual $\mathbb{Z}_{4}$ symmetry, $\mathbb{Z}_{4}$, $U(1)_{em}$, and Poincare symmetries allow $S$ to decay only to $\nu$’s. Assuming all BSM $\mathbb{Z}_{4}$ non-singlets are heavier than $S$, the decay of $S$ to neutrino pair is radiative and shown in Fig. 4.

The amplitude of the diagram in Fig. 4 is given by

[TABLE]

where $\mu_{ij}^{x}$ is given in eqs. 45 and 46, $x,y$ are spinors in 2 component notation, $s(c)_{N}$ is the mixing angle of $N$ states given in eq. 17, $s(c)_{\xi}$ is the mixing of $(\eta^{0},\chi^{\ast})$ states given in eq. 29, and $C_{0}$ is given in eq. B. Then the decay width is given by

[TABLE]

where

[TABLE]

We assume the $S_{R,I}$ mass scale is above EW scale ($v=246$GeV) but below $U(1)_{B-L}$ spontaneous breaking ($v<m_{s}<v_{4}$), so at the moment of freeze-out of $S_{R,I}$ EW symmetry is conserved whereas $U(1)_{B-L}$ symmetry is broken to $\mathbb{Z}_{4}$. Then annihilation of $S$ to SM particles will proceed through the Feynman diagrams shown in Fig. 5.

The inelastic scatterring of ($S_{R},(n,p)\rightarrow S_{I},(n,p)$) assuming $S_{R}$ as DM via t$-$channel $Z^{\prime}$ mediator can be avoided using same trick as was used in (47), namely by making $\left|\Delta m_{s}^{2}\right|=2\sqrt{2}|\mu_{S_{4}}|v_{4}>O(10^{2-4}keV^{2})$. Furthermore, $Z$ SM gauge boson mediated t$-$channel DD is absent since $S_{R,I}$ couples only to $Z^{\prime}$ and $Z-Z^{\prime}$ mixing appears only at the one-loop order.

$S_{R,I}$’s lifetime, $\tau_{s}$, dependence on $v_{4}$, $m_{S}$, and $\mu^{x}$ (coupling between $\xi_{i}\xi_{j}$ and $s_{R,I}$; eqs. 45 and 46) is shown in Fig. 6.

The age of the Universe is $4.35\times 10^{17}$s, but the bound on decaying dark matter (48) is much greater, $\tau>10^{25}$, from the cosmic microwave background (CMB) constraint.

Important remark regarding Fig. 6 and $S_{R,I}$ being a viable DM candidate is that to make $S_{R,I}$ long-lived, $\tau_{s}>10^{25}$s, $\mu^{x}$ must be tiny ($\approx 10^{-5}$eV). As will be explained in sec. VII, in order to have enhanced $0\nu 4\beta$  $\mu^{x}\sim\mu_{S}\approx 10^{6}$TeV is required. So, for $S_{R,I}$ to be a viable DM candidate means strongly suppressed $0\nu 4\beta$. There are two ways to make $\mu^{x}$ tiny: either $\left|\sqrt{2}\mu_{s}+\lambda_{2}v_{4}\right|<2\times 10^{-15}$TeV (strong fine-tuning), which will allow for observable $0\nu 4\beta$  via the other $S$ component ($S_{I}$) or $\mu_{S},\lambda_{2}v_{4}<10^{-15}$TeV in which can $0\nu 4\beta$  will be strongly suppressed.

We assume that the mass splitting $\Delta m^{2}_{S}$ between $S_{R}$ and $S_{I}$ is small, therefore both $S_{R}$ and $S_{I}$ freeze-out simultaneously (with $S_{I}$ decaying to $S_{R}$ for $m_{S_{R}}<m_{S_{I}}$). Diagrams shown in Fig. 5 contribute to $\sigma(S_{R,I}S_{R,I}\rightarrow XX)$ $S$ annihilation cross-section in order to get the correct relic abundance for $S_{R,I}$, $\Omega_{S}h^{2}=0.120$ (49). The contact diagram annihilation to Higgs pair is dominant since the $Z^{\prime}$ s$-$channel diagram is suppressed due to large $m_{Z^{\prime}}>4.2$TeV(Sec. VIII). Even at the resonance, $m_{Z^{\prime}\approx 2m_{S}}$, the $Z^{\prime}$ s$-$channel diagram is sub-dominant due to $g_{B-L}<0.127$(Sec. VIII). Therefore, $S_{R}$ relic abundance and effective annihilation cross-section for $S_{R,I}$ as a function of DM mass($m_{S_{R}}$) and coupling $\lambda_{HS}$ with other parameters fixed is plotted in Fig. 7.

The correlation between $m_{S_{R}}$, $\lambda_{HS}$, and $g_{B-L}$ is shown in Fig. 8 for the range $0.119<\Omega_{S}h^{2}<0.121$ and $\frac{m_{S_{I}}}{m_{S_{R}}}-1=10^{-2}$, $M_{Z^{\prime}}=4.2$TeV fixed.

As can be seen $S$ is a viable long-lived DM candidate that also allows for correct neutrino masses to be satisfied but will simultaneously lead to highly suppressed $0\nu 4\beta$  signal. In this case $0\nu 2\beta$ decay is forbidden by Dirac nature of neutrino masses, whereas $0\nu 4\beta$  signal is highly suppressed. As was shown above, the situation with $N$ is quite different!

VII Neutrinoless Quadruple Beta Decay

In our construction of the model, by design, due to $\mathbb{Z}_{4}$ residual symmetry neutrinoless double beta decay ($0\nu 2\beta$) is exactly absent. Therefore the dominant multipole will be neutrinoless quadruple beta decay ($0\nu 4\beta$). Contribution to neutrinoless quadruple beta decay is shown in Fig. 9. There will be also a diagram mediated by $\nu_{R}$ right-chiral neutrinos with $N_{R}$ replaced by $N_{L}$ in Fig. 9. But due to suppression with neutrino mass at every leg and seesaw suppressed $N_{L}=\left(\text{cos}\theta N_{1}-\text{sin}\theta N_{2}\right)$ Majorana mass (Eq. 15), contribution mediated by $\nu_{R}$ can be safely ignored.

Reference (50) is the first paper to study experimental side of $0\nu 4\beta$ with $B-L$ breaking to $\mathbb{Z}_{2n}$ where $n=2$ naturally leading to neutrinoless quadruple beta decay.

Neutrinoless quadruple beta decay has been searched for and experimentally studied by NEMO$-$3 collaboration in Refs. (51; 52). Another study was performed using ${}^{150}Nd$ (53) nuclei at Kimballton Underground Research Facility setting upper limit for half life-time for $0\nu 4\beta$.

The diagram in Fig. 10 effectively gives the $\mathbb{Z}_{4}$ invariant vertex

[TABLE]

The relation between interaction eigenstates $(\eta^{0},\chi)$ and mass eigenstates $\xi_{i}$ is given in Eq. 29 and $\Delta m_{s}^{2}=m_{s_{R}}^{2}-m_{s_{I}}^{2}=2\sqrt{2}\mu_{S_{4}}v_{4}$ is due to $\mu_{S_{4}}$ term in the scalar potential Eq. V. The coefficients $\mu_{ij}^{R,I}$ in the basis $(\xi_{1},\xi_{2})$ are given as

[TABLE]

where $s$ and $c$ stand for sin$\theta_{\xi}$ and cos$\theta_{\xi}$, respectively, and $\theta_{\xi}$ was defined in Eq. 29. $0\nu 4\beta$ can be calculated as two one-loop diagrams. Neutrinoless quadruple beta decay is given by

[TABLE]

where $Q^{abcd}$ represents quadruple strength, $\Lambda$ is the new physics scale relevant for the neutrinoless quadruple beta decay. $Q^{abcd}/\Lambda^{2}$ explicitly is given by

[TABLE]

where the sum over repeated indices is assumed and the Majorana $N$ mass represents the $\Lambda$ scale in Eq. 47. $a,b,c,d,\alpha,\beta,\gamma,\delta$ are flavor indices and take values $1-3$. $v_{S}$ and $M_{S}$ are given as

[TABLE]

$\theta_{\xi}$ is mixing angle between $\eta$ and $\chi$ scalars and was given in Eq. 29. $F(x,y,z)$ is the loop function and is given by

[TABLE]

In Eq. VII, $x_{ij}$ and $y_{ij}$ are given by

[TABLE]

where $m_{i}$ is the mass eigenstate of $\xi_{i}$ given in Eq. V and $m_{N_{i}}$ is the Majorana mass eigenstate of $N_{i}$ given in Eq. 15. $s_{N}$ and $c_{N}$ stand for the sine and cosine of the mixing angle of the $N_{L,R}$ fermions and are given in Eq. 17. Lastly, $c_{4}$ is the combinatorics factor and is given as

[TABLE]

Important remark regarding eq. VII is the presence of the $\lambda_{2}v_{4}$ cross term in the last line. If $\lambda_{2}$ was absent (forbidden) in the model then $0\nu 4\beta$ would be proportional to the splitting of $S$ scalar and pseudo-scalar masses, which is controlled by $\mu_{S_{4}}$ term in eq. V. Neutrino mass suppression factors like $Y_{L}$, $\theta_{\xi}$, $\Delta m_{\xi}$, $v_{4}$ also suppress $0\nu 4\beta$  but $\mu_{S}$ freedom can be used to control the enhancement of $0\nu 4\beta$. In the case if $\mu_{S}\gg v_{4}\sim O(10^{2-3}\text{TeV})$ the $\mu_{S}^{2}$ term will dominate and $0\nu 4\beta$  will scale as $\mu_{S}^{2}$.

Below numerical calculation of $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ is performed using pySecDec (54) software tool. Diagrams that have dominant contribution to $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ are the ones with $\nu_{L}$ legs and are shown in Fig. 11. There are also diagrams with $\nu_{L}$ replaced by $\nu_{R}$ but they are suppressed by a factor of $\frac{m_{\nu}}{p_{\nu}}$ for each $\nu_{L}\rightarrow\nu_{R}$ leg replacement.

Diagrams in Fig. 11 produce loop integrals

[TABLE]

with $\sum_{i}p_{i}=0$ and $Y_{L}^{\prime}$ given below.

[TABLE]

where

[TABLE]

After using pySecDec python code to calculate these integrals numerically, we compare numerical results with analytically obtained results in eq. VII and plot both in Figs. 12 and 13.

$\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ dependence on $m_{\xi}$, $v_{4}$, and $\mu_{S}$ for the analytical result from eq. VII is plotted in fig. 12 with the other parameters fixed. As can be seen from eq. VII, for $\lambda_{2}v_{4}\ll\mu_{S}$ $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}\propto\mu_{S}^{2}$ and the $\mu_{S}$ can be used to enhance the $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ for possible detection in the upcoming $0\nu 4\beta$  experiments.

Current half-life lower limit on $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ is given as $\tau_{1/2}^{0\nu 4\beta}>3.2\times 10^{21}$years (51). The relation between half-life and $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ is given as

[TABLE]

where $G_{0\nu 4\beta}$ is the four particle phase space factor and $A_{0\nu 4\beta}$ is the matrix element for $0\nu 4\beta$  process. (51) and (53) use 150Nd$\rightarrow^{150}$Gd which has $Q=2.079-2.084$MeV. $q$ can be estimated as $p_{\nu}=|q|\approx 1$fm${}^{-1}\approx 100$MeV. Then $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ can be estimated from

[TABLE]

where the last factor was inserted for dimensional matching. Using this estimate and half-life lower limit of $\tau_{1/2}^{0\nu 4\beta}>3.2\times 10^{21}$ we get

[TABLE]

Fig. 13 shows the comparison of numerical results from pySecDec with approximate analytical expression from eq. VII. As can be seen from the plot, the coupled loop(Fig. 11(b)) is relevant at $v_{4}<10^{6}$TeV scales, where it of the order of the decoupled loop(Fig. 11(a)) and can interfere destructively ($10^{3}\text{TeV}<v_{4}<10^{6}\text{TeV}$). For $v_{4}>10^{6}$TeV coupled loop becomes irrelevant. Important thing to notice is that $\mu_{S}$ plays crucial role at enhancing $0\nu 4\beta$  at $v_{4}<10^{6}$TeV scales (for $\lambda_{2}=1$). So, this model predicts a possibility for an enhanced $\frac{Q_{0\nu 4\beta}}{\Lambda^{2}}$ which can be probed in the future $0\nu 4\beta$  experiments.

VIII Results and discussion

Concerning collider constraints on the model: as can be seen from Tab. 4 $g_{B-L}$ has a upper bound of $0.25$ from collider searches and upper bound of $0.127$ from perturbativity constraints (due to large $Q_{B-L}$ charges in the model). This together with a lower bound on $M_{Z^{\prime}}$ give a lower bound on $U(1)_{B-L}$ breaking scale $v_{4}$. Constraints on $v_{4}$ are shown in Fig. 14.

Neutrino mass scale dependency on $v_{4}$, $m_{\xi}$, $\frac{Y_{D}}{Y_{M}}$, and $\frac{m_{\xi_{1}}}{m_{\xi_{2}}}$ is shown in Fig. 15.

The correlation between $v_{4}$ scale, $m_{\xi}$, and $\frac{Y_{D}}{Y_{M}}$ ratio for a fixed $\frac{m_{\xi_{1}}}{m_{\xi_{2}}}-1=10^{-2}$ mass splitting that produce neutrino mass of the order $O(1-0.1\text{eV})$ is shown in Fig. 16.

Neutrino mass can be made small by the following ways: loop suppression, small $Y_{L,R,D}$ Yukawas, large $v_{4}$ $U(1)_{B-L}$ breaking scale, and small mass splitting between $\xi_{i}$ mass eigenstates. In order to suppress neutrino mass but keep $0\nu 4\beta$  large the following parameter choices were made: $Y_{R}$, $Y_{D}\ll 1$ to suppress $m_{\nu}$; $Y_{L}\sim O(1)$, $\theta_{\xi}=\pi/4$(max. $\xi$ mixing) since $0\nu 4\beta$  has quartic dependence on them; $Y_{M}\sim O(1)$ since $0\nu 4\beta$  depends quadratically on it; $\mu_{S}\gg v_{4}$ since it is used to enhance the $0\nu 4\beta$.

Detailed study of phenomenology of the $U(1)_{B-L}$ model is done in (36). Now, in our model the right handed neutrinos can have a strongly hierarchical neutrino Yukawa structure. Which can create leptonic asymmetry through the decays right handed neutrino as shown in (64; 65).

IX Conclusion

The $U(1)_{B-L}$ extension of the SM was presented which is then spontaneously broken to residual $\mathbb{Z}_{4}$ symmetry. The $\mathbb{Z}_{4}$ symmetry is both responsible for the Dirac nature of neutrinos as well as for the stability of DM, a unique feature for this type of construction. Neutrino masses are generated radiatively through scotogenic scenario. Since the neutrinos are of Dirac type the neutrinoless double beta decay is exactly absent, but the $\mathbb{Z}_{4}$ symmetry allows for non-zero neutrinoless quadruple beta decay, which is despite being an experimentally tiny effect is the dominant of the neutrinoless multipole beta decays. If future experiments on $0\nu 2n\beta$ see no positive results in $0\nu 2\beta$ but do observe non-zero $0\nu 4\beta$, this will be a strong indication toward neutrinos of Dirac type while still violating lepton number by 4 units and will hint toward this type of model. $\mathbb{Z}_{4}$ allows for several WIMP like DM candidates in our model: best DM candidate is Majorana $N$ which allows for small neutrino masses of $O(0.1\text{eV})$ scale, enhanced $0\nu 4\beta$ decay, $U(1)_{B-L}$ breaking scale as low as $O(10\text{TeV})$, and DM masses of $O(1\text{TeV})$; other possible DM candidate is $S$ real scalar field which has a radiative decay to neutrinos and is suitable long-lived DM candidate, making $S$ long-lived also suppresses $0\nu 4\beta$ decay, so it predicts no observable $0\nu 4\beta$ in current or future $0\nu 2n\beta$ experiments without fine-tuning. In many models like this one, $0\nu 4\beta$ might be predicted to be non-zero but even in that case it is expected to be well below the sensitivity of current and future experiments looking for $0\nu 2n\beta$ decays. Model presented here allows for arbitrary enhanced $0\nu 4\beta$ decay which can be made as large as $10^{16-19}$. The prize we pay for this is the introduction of $S$ field which gives us a freedom of the enhancement of $0\nu 4\beta$ without effecting neutrino mass generation and DM related processes (for the $N$ DM case). We have also shown that the model can satisfy all required collider constraints. More detailed collider phenomenology will be presented elsewhere. Here we focused on demonstrating a way for Dirac scotogenic neutrinos with $\Delta L=4$ and dominant $0\nu 4\beta$ decay multipole where Baryon and Lepton number symmetries and violations are obtained from $U(1)_{B-L}$ gauge symmetry.

Acknowledgements.

This work was supported by the National Research Foundation of Korea Grants No. 2009-0083526, No. 2017K1A3A7A09016430, and No. 2017R1A2B4006338.

Appendix A Note on removing phases from Yukawa terms in $N_{L,R}$ sector

[TABLE]

with

[TABLE]

where $c=$cos$\theta$, $s=$sin$\theta$, and in the last equality we have set $\Delta\alpha=2\Delta$. As can be seen $\Delta\alpha$ and $\alpha_{N}$ Majorana phases can be used to freely adjust $\phi$ and $\Delta$ phases in the unitary transformation. Next,

[TABLE]

where

[TABLE]

As can be seen, Majorana phases of $N_{L,R}$ fermion fields can be used to remove phases from the mass matrix.

Appendix B $C_{0}(s_{1},s_{12},s_{2};m_{0},m_{1},m_{2})$

[TABLE]

where Kallen $\lambda$ is defined as

[TABLE]

Appendix C $N$ annihilation diagrams