# Neutrino masses and mixing from milli-charged dark matter

Sudip Jana<EMAIL_ADDRESS>Max-Planck-Institut für Kernphysik,
Saupfercheckweg 1, 69117 Heidelberg, Germany  Michael Klasen
<EMAIL_ADDRESS>Institut für Theoretische Physik, Universität
Münster, Wilhelm-Klemm-Straße 9, 48149 Münster, Germany  Vishnu P.K.
<EMAIL_ADDRESS>Institut für Theoretische Physik, Universität
Münster, Wilhelm-Klemm-Straße 9, 48149 Münster, Germany  Luca Paolo Wiggering
<EMAIL_ADDRESS>Institut für Theoretische Physik, Universität
Münster, Wilhelm-Klemm-Straße 9, 48149 Münster, Germany

###### Abstract

We propose a simple extension to the Standard Model, wherein neutrinos
naturally attain small Majorana masses through a one-loop radiative mechanism
featuring particles within the loops characterized by milli-charges. Unlike
the conventional scotogenic model, our approach avoids imposing a discrete
symmetry or expanding the gauge sector. The minuscule electric charges ensure
the stability of the lightest particle within the loop as a viable dark matter
candidate. Our investigation systematically scrutinizes the far-reaching
phenomenological implications arising from these minuscule charges.

††preprint: MS-TP-24-12

##  I Introduction

The electric charges of all known elementary particles are integer multiples
of $Q_{e}/3$, where $Q_{e}$ represents the electron charge. The apparent
quantization of electric charge remains a persistent puzzle in particle
physics Marinelli:1983nd ; Dylla:1973zz ; Gahler:1982xi ; Kyuldjiev:1984kz ;
Davidson:2000hf . Developing a theoretical understanding of this quantization
is crucial for comprehending the interactions of elementary particles.
However, it is theoretically possible for particles to have electric charges
of $\epsilon Q_{e}$, where $\epsilon$ is any real number less than 1. These
are referred to as milli-charged particles. Milli-charged particles can be
incorporated into the Standard Model through several mechanisms. One well-
studied approach introduces a new dark gauge symmetry $U(1)^{\prime}$, where
mixing between the photon and the dark photon allows particles charged under
$U(1)^{\prime}$ to have a small coupling to the photon, resulting in a milli-
charge Holdom:1985ag . Another method allows Standard Model neutrinos to have
tiny charges consistent with non-standard electric charges, given by
$Q_{st}+\epsilon(L_{i}-L_{j})$, where $i,j=e,\mu,\tau\ (i\neq j)$, $Q_{st}$ is
the standard electric charge, $L_{i}$ is a family-lepton number, and
$\epsilon$ is an arbitrary parameter. These non-standard charges are
compatible with the Standard Model Lagrangian and preserve the gauge anomaly
cancellation of $SU(3)_{c}\otimes SU(2)_{L}\otimes U(1)_{Y}$ Foot:1990uf ;
Babu:1992sw .

In contrast, without introducing any additional symmetries and preserving the
charges of Standard Model particles, we extend the SM by adding vector-like
fermions or scalars with a small fractional charge $\epsilon$. The anomaly
cancellation of the SM is maintained because the fermions are vector-like.
These vector-like fermions and scalars will manifest as milli-charged
particles within the theory. Our focus is on determining whether such milli-
charged particles can address some of the significant unresolved issues in
particle physics. Specifically, we are interested in exploring whether these
particles could serve as candidates for dark matter (DM), given the compelling
evidence supporting the existence of DM as a strong indication for physics
beyond the Standard Model Bertone:2004pz ; Klasen:2015uma . Additionally,
considering the well-established phenomenon of neutrino flavor oscillations
confirmed by numerous experiments Super-Kamiokande:1998kpq ; SNO:2001kpb ;
SNO:2002tuh , we investigate whether milli-charged DM could play a crucial
role in neutrino mass generation.

Our aim is to address the two major challenges of neutrino mass generation and
DM stability within a unified framework employing milli-charged particles. The
scotogenic model provides an elegant solution by introducing a symmetry that
enables both the suppression of neutrino masses through loop processes and the
stabilization of DM Tao:1996vb ; Ma:2006km . Traditionally, this approach
necessitates additional discrete symmetries Tao:1996vb ; Krauss:2002px ;
Ma:2006km ; Aoki:2008av ; Ma:2008cu ; Restrepo:2013aga ; Gustafsson:2012vj ;
Fiaschi:2018rky ; Esch:2018ccs ; Jana:2020joi , which are vulnerable to
gravitational anomalies, or requires the extension of the SM gauge group
Ma:2013yga ; Ho:2016aye ; deBoer:2023phz ; Jana:2019mez ; Jana:2019mgj . We
propose a minimal extension to the Standard Model in which neutrinos naturally
acquire small Majorana masses via a one-loop radiative seesaw mechanism
involving milli-charged particles. Unlike the conventional scotogenic model,
our method does not require discrete symmetries or the expansion of the gauge
sector. The tiny electric charges ensure the stability of the lightest
particle in the loop, making it a viable candidate for DM.

In our model, we propose that the lightest milli-charged particle, a vector-
like fermion $F$, can serve as a viable DM candidate. Its minuscule electric
charge ensures its stability and allows for DM-photon coupling at tree level.
However, a purely photon-portal milli-charged DM model is inconsistent with
current observations McDermott:2010pa . The Planck satellite’s CMB data rules
it out because milli-charged DM particles would couple with the baryon-photon
plasma even at low temperatures, affecting CMB acoustic peaks and damping DM
density fluctuations. Considering the freeze-in mechanism, where $F$ is
initially out of equilibrium with the SM thermal bath but produced later
through feeble interactions, generally works for sub-GeV masses
Dvorkin:2019zdi . However, if milli-charged DM also contributes to neutrino
mass generation, significant Yukawa interactions are needed, which immediately
establish thermal equilibrium between dark matter and the SM bath given a
negligible dark matter abundance after reheating, making the freeze-in
mechanism incompatible with our model. Therefore, while milli-charged DM can
generate neutrino masses and mixing successfully, it must annihilate into SM
particles via gauge or leptonic portals through the freeze-out scenario to
accurately reproduce the relic density. We analyze this scenario for milli-
charged DM from sub-GeV to TeV mass regimes, discussing various
phenomenological consequences such as lepton flavor-violating observables and
electroweak precision constraints. Additionally, the presence of milli-charged
DM modifies the amplitude of $0\nu\beta\beta$ decay, with implications for
current and future experiments.

The structure of the paper is as follows: First, in Sec. II, we present the
model. Subsequently, in Sec. III, we analyze the milli-charged DM scenarios
for both freeze-in and freeze-out mechanisms, discussing the current
constraints and various phenomenological implications. Following this, we
introduce in Sec. IV the DM scenario specific to our model, analyze its
phenomenology demonstrating its capability to generate neutrino masses and
mixing successfully. The resulting implications on the neutrinoless double
beta decay are presented in Sec. V. We conclude in Sec. VI.

##  II The Model

We extend the SM with a single generation of a milli-charged fermion $F$ and
two electroweak (EW) scalar doublets $\phi_{1,2}$. No additional symmetries
are imposed. The quantum numbers of these multiplets under the SM gauge
symmetry $SU(3)_{c}\otimes SU(2)_{L}\otimes U(1)_{Y}$ are given below

$\displaystyle F_{R,L}\sim(1,1,\epsilon),$ (1)
$\displaystyle\phi_{1}=\begin{pmatrix}\phi_{1}^{1+\epsilon}\\\
\phi_{1}^{\epsilon}\end{pmatrix}\sim(1,2,\dfrac{1}{2}+\epsilon),$ (2)
$\displaystyle\phi_{2}=\begin{pmatrix}\phi_{2}^{1-\epsilon}\\\
\phi_{2}^{-\epsilon}\end{pmatrix}\sim(1,2,\dfrac{1}{2}-\epsilon)\,.$ (3)

Among the five BSM states, three of them
$\\{F,\phi_{1}^{\epsilon},\phi_{2}^{\epsilon}\\}$ carry electric charge
$\epsilon$. The lightest of these milli-charged particles can be a potential
candidate for DM, where its stability is guaranteed by its minuscule electric
charge.

The Yukawa interactions of $\\{F,\phi_{1},\phi_{2}\\}$ and bare mass of $F$
are given by

$\displaystyle-\mathcal{L}_{\rm Yuk}\supset$ $\displaystyle
Y_{1}\overline{F}_{R}\widetilde{\phi}_{1}^{\dagger}\ell_{L}+Y_{2}\overline{\ell^{c}_{L}}(i\sigma_{2})\phi_{2}F_{L}+m_{F}\overline{F}_{L}F_{R}+h.c.$
(4)

Here, $\ell_{L}\equiv(\nu_{L},e_{L})^{T}$ is the SM lepton doublet, while the
corresponding charge conjugate field is denoted as $\ell^{c}_{L}$,
$\widetilde{\phi}_{1}\equiv i\sigma_{2}\phi^{*}_{1}$, and $\sigma_{2}$ denotes
the second Pauli matrix.

The scalar potential of the model can be written as

$\displaystyle V(H,\phi_{1},$
$\displaystyle\phi_{2})=-\mu_{H}^{2}H^{\dagger}H+\mu_{\phi_{1}}^{2}\phi_{1}^{\dagger}\phi_{1}+\mu_{\phi_{2}}^{2}\phi_{2}^{\dagger}\phi_{2}$
$\displaystyle+\frac{\lambda_{H}}{2}(H^{\dagger}H)^{2}+\frac{\lambda_{\phi_{1}}}{2}(\phi_{1}^{\dagger}\phi_{1})^{2}+\frac{\lambda_{\phi_{2}}}{2}(\phi_{2}^{\dagger}\phi_{2})^{2}$
$\displaystyle+\lambda_{H\phi_{1}}(H^{\dagger}H)(\phi_{1}^{\dagger}\phi_{1})+\lambda_{H\phi_{1}}^{\prime}(H^{\dagger}\phi_{1})(\phi_{1}^{\dagger}H)$
$\displaystyle+\lambda_{H\phi_{2}}(H^{\dagger}H)(\phi_{2}^{\dagger}\phi_{2})+\lambda_{H\phi_{2}}^{\prime}(H^{\dagger}\phi_{2})(\phi_{2}^{\dagger}H)$
$\displaystyle+\left\\{\lambda_{H\phi_{1}\phi_{2}}(H^{\dagger}\phi_{1})(H^{\dagger}\phi_{2})+h.c.\right\\}$
$\displaystyle+\lambda_{\phi_{1}\phi_{2}}(\phi_{1}^{\dagger}\phi_{1})(\phi_{2}^{\dagger}\phi_{2})+\lambda_{\phi_{1}\phi_{2}}^{\prime}(\phi_{1}^{\dagger}\phi_{2})(\phi_{2}^{\dagger}\phi_{1}).$
(5)

Here, $H\equiv(0,\frac{v+h}{\sqrt{2}})^{T}$ represents the SM Higgs doublet
after EW symmetry braking with $v\simeq 246\,\mathrm{GeV}$ denoting the EW
vacuum expectation value and $h$ is the SM Higgs. The quartic coupling
$\lambda_{H\phi_{1}\phi_{2}}$ leads to mixing between the milli-charged
scalars $\phi_{1}^{\epsilon}$ and $\phi_{2}^{\epsilon}$, resulting in the
following mass squared matrix:

$\displaystyle\mathcal{M}^{2}_{\phi^{{\epsilon}}}=\begin{pmatrix}\mu^{2}_{\phi_{1}}+\frac{(\lambda_{H\phi_{1}}+\lambda^{\prime}_{H\phi_{1}})}{2}v^{2}&\frac{\lambda_{H\phi_{1}\phi_{2}}}{2}v^{2}\\\
\frac{\lambda_{H\phi_{1}\phi_{2}}}{2}v^{2}&\mu^{2}_{\phi_{2}}+\frac{(\lambda_{H\phi_{2}}+\lambda^{\prime}_{H\phi_{2}})}{2}v^{2}\end{pmatrix}.$
(6)

The corresponding mass eigenstates $S_{1,2}^{{\epsilon}}$ with masses
$m_{S_{1,2}}$ are related to $\phi_{1,2}^{{\epsilon}}$ through the following
mixing matrix

$\begin{pmatrix}S_{1}^{{\epsilon}}\\\
S_{2}^{{\epsilon}}\end{pmatrix}=\begin{pmatrix}\cos{\theta}&-\sin{\theta}\\\
\sin{\theta}&\cos{\theta}\end{pmatrix}\begin{pmatrix}\phi_{1}^{{\epsilon}}\\\
\phi_{2}^{{\epsilon}}\end{pmatrix}\,.$ (7)

Here, $\theta$ is the mixing angle, which can be related to the quartic
coupling $\lambda_{H\phi_{1}\phi_{2}}$ via

$\sin{2\theta}=\dfrac{\lambda_{H\phi_{1}\phi_{2}}v^{2}}{m^{2}_{S_{2}}-m^{2}_{S_{1}}}.$
(8)

The masses of the scalars $\phi_{1}^{1+{\epsilon}}$ and
$\phi_{2}^{1-{\epsilon}}$ are denoted as $m_{\phi_{1}}$ and $m_{\phi_{2}}$,
respectively, and are given by

$\displaystyle
m^{2}_{\phi_{1}}=\mu^{2}_{\phi_{1}}+\frac{\lambda_{H\phi_{1}}}{2}v^{2},\quad
m^{2}_{\phi_{2}}=\mu^{2}_{\phi_{2}}+\frac{\lambda_{H\phi_{2}}}{2}v^{2}.$ (9)

Figure 1: Feynman diagram for neutrino masses and mixings generated by milli-
charged particles.

The simultaneous presence of the Lagrangian terms in Eq. (4) and the quartic
coupling term
$\lambda_{H\phi_{1}\phi_{2}}(H^{\dagger}\phi_{1})(H^{\dagger}\phi_{2})+h.c$
leads to lepton number violation by two units, which results in radiative
neutrino mass generation at the one-loop level via the diagram depicted in
Fig. 1. The corresponding neutrino mass matrix is given by

$\displaystyle\mathcal{M}^{\nu}$
$\displaystyle=m_{0}\left(Y_{1}.Y_{2}^{T}+Y_{2}.Y_{1}^{T}\right),$ (10)

where

$\displaystyle m_{0}$ $\displaystyle=\frac{m_{F}\sin
2\theta}{32\pi^{2}}\left\\{\frac{m_{S_{2}}^{2}\ln\left(\frac{m_{S_{2}}^{2}}{m^{2}_{F}}\right)}{m_{S_{2}}^{2}-m_{F}^{2}}-\frac{m_{S_{1}}^{2}\ln\left(\frac{m_{S_{1}}^{2}}{m^{2}_{F}}\right)}{m_{S_{1}}^{2}-m_{F}^{2}}\right\\}.$
(11)

The above mass matrix can be diagonalized by a unitary transformation of the
form $U^{T}\mathcal{M}^{\nu}U$, where $U$ denotes the PMNS matrix. Note that
$\mathcal{M}^{\nu}$ is of rank 2, which implies that the lightest neutrino is
massless in our model. Moreover, the particular form of $\mathcal{M}^{\nu}$
allows one to parameterize the Yukawa couplings $Y_{1,2}$ in terms of the
elements of the PMNS matrix, two non-zero neutrino masses $m_{i}$, and masses
and mixings of the BSM states that enter in the loop diagram shown in Fig. 1.
The corresponding parametrization for the normal hierarchy (NH) of neutrino
masses is given by Cordero-Carrion:2018xre ; Cordero-Carrion:2019qtu ;
Hagedorn:2018spx

$\displaystyle
Y_{1}=\sqrt{\frac{k}{2m_{0}}}\left(\begin{array}[]{c}\sqrt{m_{3}}U_{13}^{*}+i\sqrt{m_{2}}U_{12}^{*}\\\
\sqrt{m_{3}}U_{23}^{*}+i\sqrt{m_{2}}U_{22}^{*}\\\
\sqrt{m_{3}}U_{33}^{*}+i\sqrt{m_{2}}U_{32}^{*}\\\ \end{array}\right),$ (15)
$\displaystyle
Y_{2}=\sqrt{\frac{1}{2m_{0}k}}\left(\begin{array}[]{c}\sqrt{m_{3}}U_{13}^{*}-i\sqrt{m_{2}}U_{12}^{*}\\\
\sqrt{m_{3}}U_{23}^{*}-i\sqrt{m_{2}}U_{22}^{*}\\\
\sqrt{m_{3}}U_{33}^{*}-i\sqrt{m_{2}}U_{32}^{*}\\\ \end{array}\right),$ (19)

where $k$ signifies the relative hierarchy between the Yukawa couplings
$Y_{1}$ and $Y_{2}$. For the inverted hierarchy (IH) of neutrino masses, the
parametrization takes the following form

$\displaystyle
Y_{1}=\sqrt{\frac{k}{2m_{0}}}\left(\begin{array}[]{c}\sqrt{m_{1}}U_{11}^{*}+i\sqrt{m_{2}}U_{12}^{*}\\\
\sqrt{m_{1}}U_{21}^{*}+i\sqrt{m_{2}}U_{22}^{*}\\\
\sqrt{m_{1}}U_{31}^{*}+i\sqrt{m_{2}}U_{32}^{*}\\\ \end{array}\right),$ (23)
$\displaystyle
Y_{2}=\sqrt{\frac{1}{2m_{0}k}}\left(\begin{array}[]{c}\sqrt{m_{1}}U_{11}^{*}-i\sqrt{m_{2}}U_{12}^{*}\\\
\sqrt{m_{1}}U_{21}^{*}-i\sqrt{m_{2}}U_{22}^{*}\\\
\sqrt{m_{1}}U_{31}^{*}-i\sqrt{m_{2}}U_{32}^{*}\\\ \end{array}\right).$ (27)

For our analysis, we use the above parametrizations for $Y_{1}$ and $Y_{2}$.

##  III Constraints on the millicharged dark matter

Before investigating the phenomenological implications of our model, we
discuss the constraints on the standard milli-charged DM scenario, where the
DM couples to the SM particles only via photons. These constraints are
summarized in Fig. 2.

Figure 2:  Current status of milli-charged DM in $m_{F}$-$\epsilon$ plane. The
bounds are imposed by the following constraints: CMB+Lyman-$\alpha$
constraints Xu:2018efh ; Planck:2015fie ; SDSS:2004kjl  (light-orange region),
DM halo constraints from the bullet cluster McDermott:2010pa ; Ibe:2009mk
(light-purple region) and elliptical galaxies McDermott:2010pa ;
Markevitch:2003at  (orange region), MilliQAn demonstrator Citron:2021woe ;
Ball:2020dnx  (light-green region), SLACmQ Prinz:1998ua ; Davidson:2000hf
(purple region), LEP + beam-dump experiments Davidson:2000hf  (green region),
LSND Magill:2018tbb ; LSND:1996jxj  (brown region), MiniBooNE Magill:2018tbb ;
MiniBooNE:2018esg  (maroon region), ArgoNeuT ArgoNeuT:2019ckq  (blue region),
and DM direct detection constraints from DarkSide-50 DarkSide:2018bpj ;
Belanger:2020gnr  (blue shaded region), Xenon1T XENON:2018voc ;
Belanger:2020gnr  (red shaded region), and Lux-Zeplin LZ:2022lsv  (green
shaded region). The dashed and dot-dashed curves represent the DM relic
density $\Omega h^{2}=0.12$ for the scenarios of freeze-out and freeze-in
mechanisms, respectively. See text for more details.

The interactions of milli-charged DM with the SM thermal plasma during the
recombination epoch are constrained by the precise measurements of CMB
observables. Limits on such interactions were derived in Ref. Xu:2018efh
using Planck satellite data Planck:2015fie  and the Sloan Digital Sky Survey
SDSS:2004kjl , which are indicated in Fig. 2 through the light-orange shaded
region.

There are strong constraints on the electric charge of DM from direct
detection experiments since the charged DM can scatter off the nuclei at tree
level via photon exchange. The bounds from DarkSide-50 DarkSide:2018bpj  and
Xenon1T XENON:2018voc  derived in Ref. Belanger:2020gnr , are depicted in Fig.
2 by the blue-shaded region and the red-shaded region, respectively. Using
these, we recast the limits from Lux-Zeplin LZ:2022lsv , which are represented
in Fig. 2 by the green-shaded region.

Figure 3: Different portals available for DM annihilation/production in our
model.

At collider or fixed target experiments, milli-charged DM can be produced via
the Drell-Yan process or via the decay of various mesons. The null results
from the following experiments impose stringent constraints on the millicharge
and are depicted in Fig. 2: MilliQAn demonstrator Citron:2021woe ;
Ball:2020dnx  (light-green colored region), SLACmQ Prinz:1998ua ;
Davidson:2000hf  (purple colored region), LEP + beam-dump experiments
Davidson:2000hf  (green colored region), LSND Magill:2018tbb ; LSND:1996jxj
(brown colored region), MiniBooNE Magill:2018tbb ; MiniBooNE:2018esg  (maroon
colored region), and ArgoNeuT ArgoNeuT:2019ckq  (blue colored region).

The self-interaction cross section of DM is constrained through observations
of dark matter halos like galaxy cluster mergers and the ellipticity of
galaxies. The resulting limits from the bullet cluster Ibe:2009mk  and the
elliptical galaxy NGC 720 Markevitch:2003at  on the electromagnetic charge of
dark matter are derived in Ref. McDermott:2010pa  and shown in Fig. 2 in
light-purple and orange, respectively.

##  IV Dark matter Phenomenology

In our model, the lightest milli-charged particle can be a viable candidate
for DM, which we choose to be the vector-like fermion $F$. The minuscule
electric charge of $F$ ensures its stability, which also results in tree-level
DM-photon coupling. In addition to its gauge interactions, $F$ also has
leptonic interactions via the Yukawa couplings given in Eq. (4). Hence, the
milli-charged DM candidate in our model is not pure gauge portal DM like in
Refs. Dvorkin:2019zdi ; McDermott:2010pa . As discussed earlier, these Yukawa
couplings are crucial for explaining neutrino masses and mixings. We will show
later in this section that they also play a vital role in DM phenomenology.

We investigate both the freeze-in and the freeze-out mechanisms for explaining
the observed amount of DM. The dominant channels that contribute to the relic
abundance are shown in Fig. 3.

### Freeze-in scenario

First, we look into the freeze-in scenario. In this case, the DM candidate $F$
is assumed to a have a negligible abundance after reheating but is produced
from the SM bath at later times.

In the case of pure gauge coupling-induced processes, $F$ is mainly produced
via the annihilation of the SM particles (see Fig. 3). The amplitude of these
processes depends on the mass of $F$ and its electric charge $\epsilon$. Using
micrOMEGAs 5.3.41 Belanger:2018ccd , we calculated the relic abundance for
this scenario. The corresponding parameter space consistent with the Planck
relic density constraint Planck:2018vyg  $\Omega h^{2}=0.12$ is shown as a
dot-dashed line in Fig. 2. As one can see, $\epsilon\sim 10^{-11}$ is required
to generate the correct relic abundance. Such values of $\epsilon$ are also
consistent with the out-of-equilibrium condition that is required for the
successful working of the freeze-in scenario.

However, once the Yukawa couplings are turned on as per Eq. (19) or Eq. (27),
a random parameter scan revealed that the DM immediately thermalizes with the
SM plasma, rendering the neutrino mass mechanism and the freeze-in mechanism
incompatible with each other in our model.

To illustrate this, first, we calculate the limit on the Yukawa couplings
imposed by the out-of-equilibrium condition, which read as

$\displaystyle\Gamma(\phi\xrightarrow{}F\ell)\lesssim H|_{T\simeq
m_{\phi}}\Longrightarrow|Y_{1,2}|\lesssim
10^{-8}\sqrt{\dfrac{m_{\phi_{1,2}}}{100\,\mathrm{GeV}}}.$ (28)

Here $H\simeq\frac{T^{2}}{M_{\mathrm{Pl}}}$ stands for the Hubble expansion
rate, $M_{\mathrm{Pl}}=2.4\times 10^{18}\,\,\mathrm{GeV}$ is the reduced
Planck mass and $\Gamma$ denotes the production rate of $F$ via the decay of
BSM scalars. With this limit on the Yukawa couplings, we estimate the neutrino
masses for a scenario where $m_{F}\ll m_{S_{1}}\sim m_{\phi_{1}},m_{S_{2}}\sim
m_{\phi_{2}}$ giving

$\displaystyle M_{\nu}$
$\displaystyle\simeq\frac{Y_{1}Y_{2}\sin{2\theta}}{16\pi^{2}}m_{F}\ln\left(\frac{m_{S_{2}}^{2}}{m_{S_{!}}^{2}}\right),\,$
$\displaystyle\simeq\begin{cases}10^{-9}\,\mathrm{eV}\left(\frac{m_{F}}{1\,\mathrm{GeV}}\right)\left(\frac{\sin
2\theta}{1}\right)&m_{S_{1,2}}\sim\mathcal{O}(100)\,\mathrm{GeV},\\\
10^{-11}\,\mathrm{eV}\left(\frac{m_{F}}{1\,\mathrm{GeV}}\right)\left(\frac{\lambda_{H\phi_{1}\phi_{2}}}{2}\right)&m_{S_{1,2}}\sim\mathcal{O}(100)\,\mathrm{TeV},\end{cases}$
(29)

As one can see, the resulting neutrino masses are about eight to ten orders of
magnitude smaller than the observed values.

### Freeze-out scenario

Next, we turn our attention to the freeze-out mechanism. In this case, the
relic abundance of DM is mainly governed by its (co)annihilation cross-section
into SM final states ($\sigma v$). As aforementioned, in our model, $F$
annihilates into SM particles either via the gauge or the leptonic portal.

The contributions from the gauge portal processes depend on $M_{F}$ and
$\epsilon$. In Fig. 2, we show $\epsilon$ as a function of $M_{F}$ for $\Omega
h^{2}=0.12$ as a dashed line. As one can see, this parameter space is excluded
by the CMB constraints Xu:2018efh . Moreover, to be consistent with the
various constraints on the electric charge of DM, $\epsilon$ is required to be
very small, $\epsilon\ll 1$. Hence, the contribution from gauge portal
processes is negligible in our scenario such that we set $\epsilon=10^{-12}$
in the remainder of the paper.

On the other hand, the contributions that arise from the leptonic portal can
be significant. The $s$-wave components of the corresponding cross sections
are given by Herms:2023cyy

$\displaystyle(\sigma
v)_{S_{1}}^{\bar{F}F\to\bar{\nu}\nu}=\frac{m_{F}^{2}}{32\pi}\left(\frac{(|Y_{1}|^{2}c_{\theta}^{2}+|Y_{2}|^{2}s_{\theta}^{2})^{2}}{(m_{F}^{2}+m_{S_{1}}^{2})^{2}}\right),$
(30) $\displaystyle(\sigma
v)_{S_{2}}^{\bar{F}F\to\bar{\nu}\nu}=\frac{m_{F}^{2}}{32\pi}\left(\frac{(|Y_{1}|^{2}s_{\theta}^{2}+|Y_{2}|^{2}c_{\theta}^{2})^{2}}{(m_{F}^{2}+m_{S_{2}}^{2})^{2}}\right),$
(31) $\displaystyle(\sigma
v)_{\phi_{1}}^{\bar{F}F\to\ell^{+}\ell^{-}}=\frac{m_{F}^{2}}{32\pi}\left(\frac{|Y_{1}|^{4}}{(m_{F}^{2}+m_{\phi_{1}}^{2})^{2}}\right),$
(32) $\displaystyle(\sigma
v)_{\phi_{2}}^{\bar{F}F\to\ell^{+}\ell^{-}}=\frac{m_{F}^{2}}{32\pi}\left(\frac{|Y_{2}|^{4}}{(m_{F}^{2}+m_{\phi_{2}}^{2})^{2}}\right),$
(33)

where the subscripts denote the scalar mediators appearing in the underlying
$t$-channel processes. It is worth mentioning that the annihilation processes
mediated by the milli-charged scalars
$\\{S_{1}^{{\epsilon}},S_{2}^{{\epsilon}}\\}$ result in neutrinos as
annihilation products, whereas the processes mediated by
$\\{\phi_{1}^{1+{\epsilon}},\phi_{2}^{1-{\epsilon}}\\}$ result in charged
leptons.

In the following, we explore separately the low-mass and the high-mass regimes
of DM, which we define through the conditions $M_{F}<1\,\mathrm{GeV}$ and
$10\,\mathrm{GeV}<M_{F}$, respectively.

_Light dark matter_.– In this case, in order to generate an adequately large
contribution to $(\sigma v)$, at least one of the mediator states of
annihilation processes needs to be light Okawa:2020jea ; Herms:2023cyy . Since
the masses of the scalars
$\\{\phi_{1}^{1+{\epsilon}},\phi_{2}^{1-{\epsilon}}\\}$ are constrained to be
above $110$ GeV to satisfy collider bounds Babu:2019mfe ; Iguro:2022tmr ;
Jana:2020pxx , one of the milli-charged scalars
$\\{S_{1}^{{\epsilon}},S_{2}^{{\epsilon}}\\}$ has to be light to fulfill this
condition. Moreover, the sum $m_{S_{1}}+m_{S_{2}}$ needs to be larger than
$m_{Z}$ to avoid $Z$-boson decays into
$S_{1}^{\pm\epsilon}S_{2}^{\mp\epsilon}$ which are severely constrained by
measurements of its decay width Electroweak:2003ram ; Tanabashi:2018oca .
Without loss of generality, $S_{1}^{{\epsilon}}$ is therefore chosen to be the
light mediator state while the mass of $S_{2}^{{\epsilon}}$ is selected such
that $m_{S_{1}}+m_{S_{2}}>m_{Z}$. The requirement of
$m_{S_{1}}<1\,\mathrm{GeV}$ and
$m_{S_{2}},m_{\phi_{1}},m_{\phi_{2}}\sim\mathcal{O}(100)\,\,\mathrm{GeV}$
imposes additional constraints on the model parameter space, which we address
next.

Figure 4: Parameter space of the light-DM scenario consistent with the Planck
relic density constraint in terms of $(m_{F},|Y_{1}|)$ for both NH (left
panel) and IH (right panel) of neutrino masses. The red-colored data points
are excluded by cLFV constraints given in Eq. (IV). The data points
represented by green (blue) color can (not) be probed in future cLFV
experiments.

First of all, the mixing angle $\theta$ is constrained from the $Z$-decay
width measurements Electroweak:2003ram ; Tanabashi:2018oca , as it induces the
invisible decay of $Z$ into $S_{1}^{{\epsilon}}S_{1}^{{-\epsilon}}$ whose
amplitude is directly proportional to $\cos 2\theta$. To satisfy this
constraint, we choose the mixing angle to be $\theta\simeq\pi/4$ Herms:2023cyy
.

Since the mass splitting between $S_{1}^{{\epsilon}}$ and the other scalar
states
$\\{S_{2}^{{\epsilon}},\phi_{1}^{1+{\epsilon}},\phi_{2}^{1-{\epsilon}}\\}$ is
required to be large, the new contributions to the EW-precision observables,
particularly to the $T$-parameter Peskin:1990zt ; Peskin:1991sw  could be
significant in this scenario. In the limit $\epsilon\ll 1$, the corresponding
contributions to the $T$-parameter is given by Herms:2023cyy :

$\displaystyle T$ $\displaystyle\simeq\dfrac{1}{8\pi^{2}\alpha_{\rm
em}(M_{Z})v^{2}}\bigg{\\{}c_{\theta}^{2}\mathcal{F}(m_{\phi_{1}^{+}}^{2},m_{S_{1}}^{2})$
$\displaystyle+s_{\theta}^{2}\mathcal{F}(m_{\phi_{2}^{+}}^{2},m_{S_{1}}^{2})+s_{\theta}^{2}\mathcal{F}(m_{\phi_{1}^{+}}^{2},m_{S_{2}}^{2})$
$\displaystyle+c_{\theta}^{2}\mathcal{F}(m_{\phi_{2}^{+}}^{2},m_{S_{2}}^{2})-s^{2}_{2\theta}\mathcal{F}(m_{S_{1}}^{2},m_{S_{2}}^{2})\bigg{\\}}\,,$
(34)

where $\mathcal{F}$ is defined as

$\mathcal{F}(m_{1}^{2},m_{2}^{2})\ \equiv\
\frac{1}{2}(m_{1}^{2}+m_{2}^{2})-\frac{m_{1}^{2}m_{2}^{2}}{m_{1}^{2}-m_{2}^{2}}\ln\left(\frac{m_{1}^{2}}{m_{2}^{2}}\right)\,.$
(35)

For a mass hierarchy of the form $m_{S_{1}}\ll m_{S_{2}}\simeq
m_{\phi_{1}}\simeq m_{\phi_{2}}$, one obtains $T\propto\cos^{2}{2\theta}$,
which is naturally suppressed for the chosen value of $\theta\simeq\pi/4$.

Now, we investigate the implications of the quartic couplings. With the choice
of the mass hierarchy $m_{S_{1}}\ll m_{S_{2}}\simeq m_{\phi_{1}}\simeq
m_{\phi_{2}}$ and $\theta\simeq\pi/4$ , one obtains the following form for

$\displaystyle\lambda^{\prime}_{H\phi_{1,2}}\simeq-\lambda_{H\phi_{1}\phi_{2}}\simeq\frac{m^{2}_{S_{1}}-m^{2}_{S_{2}}}{v^{2}}.$
(36)

On the other hand, the combination of the quartic couplings,
$(\lambda_{H\phi_{1}}+\lambda_{H\phi_{2}}+\lambda^{\prime}_{H\phi_{1}}+\lambda^{\prime}_{H\phi_{2}}-2\lambda_{H\phi_{1}\phi_{2}})$
is required to be small, as it induces the invisible decay of the SM Higgs
into $S_{1}^{{\epsilon}}S_{1}^{{-\epsilon}}$, which is constrained by the
current bounds on the branching ratio of this decay process

$\displaystyle\mathrm{Br}_{\rm
inv}\leq\begin{cases}0.145&\mathrm{ATLAS}\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{ATLAS:2022yvh}{\@@citephrase{(}}{\@@citephrase{)}}}},\\\
0.18&\mathrm{CMS}\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{CMS:2022qva}{\@@citephrase{(}}{\@@citephrase{)}}}}.\end{cases}$
(37)

This leads to the following relation

$\displaystyle\lambda_{H\phi_{1}}+\lambda_{H\phi_{2}}$ $\displaystyle\simeq
4\lambda_{H\phi_{1}\phi_{2}}\simeq\frac{4(m^{2}_{S_{2}}-m^{2}_{S_{1}})}{v^{2}}.$
(38)

However, this condition has non-trivial consequences on the
$h\rightarrow\gamma\gamma$ process, as the new physics contributions to this
decay mode depend on $\lambda_{H\phi_{1}}$ and $\lambda_{H\phi_{2}}$ in our
model. To see this, the relevant interactions pertaining to this decay mode
are given below

$\displaystyle
V\supset\lambda_{H\phi_{1}}v(h|\phi_{1}^{1+{\epsilon}}|^{2})+\lambda_{H\phi_{2}}v(h|\phi_{2}^{1-{\epsilon}}|^{2}).$
(39)

Through these interactions, the scalars
$\\{\phi_{1}^{1+{\epsilon}},\phi_{2}^{1-{\epsilon}}\\}$ contribute to the
$h\rightarrow\gamma\gamma$ process at the one-loop level. Such new physics
contributions alter the Higgs signal strength to $\gamma\gamma$, which is
defined as
$R_{\gamma\gamma}=\frac{\mathrm{Br}(h\rightarrow\gamma\gamma)}{\mathrm{Br}(h\rightarrow\gamma\gamma)_{\rm
SM}}$. In the limit $\epsilon\ll 1$ and $m_{\phi_{1}}\simeq
m_{\phi_{2}}\simeq\mathcal{O}(100)\,\,\mathrm{GeV}$, we obtain
$R_{\gamma\gamma}\simeq 0.8$. This is consistent with experimental
measurements at the 3 sigma level ATLAS:2022tnm ; CMS:2022wpo .

In addition to these requirements on the quartic couplings, they also need to
satisfy various conditions to ensure the stability of the scalar potential.
The corresponding conditions are given by Hagedorn:2018spx

$\displaystyle\lambda_{H\phi_{1}}\geq
0,\quad\lambda_{H\phi_{1}}+\lambda^{\prime}_{H\phi_{1}}\geq 0,$ (40)
$\displaystyle\lambda_{H\phi_{2}}\geq
0,\quad\lambda_{H\phi_{2}}+\lambda^{\prime}_{H\phi_{2}}\geq 0,$ (41)
$\displaystyle\sqrt{(\lambda_{H\phi_{1}}+\lambda^{\prime}_{H\phi_{1}})(\lambda_{H\phi_{2}}+\lambda^{\prime}_{H\phi_{2}})}-\lambda_{H\phi_{1}\phi_{2}}\geq
0.$ (42)

We checked whether the above conditions are consistent with Eqs. (36) and (38)
and found that they are compatible.

Figure 5: Left (right) panel: parameter space in $\Delta m_{S_{1}}\equiv
m_{S_{1}}-m_{F}$ vs. $|Y_{1}|$ $(m_{F}$ vs. $\Delta m_{S_{1}})$ plane
consistent with the Planck relic density constraint. Here, we set
$S_{1}^{\epsilon}$ to be the coannihilation partner of DM. The horizontal
gray-shaded band denotes the region excluded by the perturbativity bound on
the Yukawa coupling $|Y_{1}|<2$. The other color codes represent the same as
in Fig. 4.

The scalars $\\{\phi_{1}^{1+{\epsilon}},\phi_{2}^{1-{\epsilon}}\\}$ contribute
to charged lepton flavor violation (cLFV) processes via the Yukawa couplings
$Y_{1,2}$. In the limit $\epsilon\ll 1$, the corresponding contributions to
$\ell_{\alpha}\to\ell_{\beta}\gamma$ can be written as Herms:2023cyy ;
Hagedorn:2018spx

$\displaystyle\frac{\mathrm{BR}\left(\ell_{\alpha}\to\ell_{\beta}\gamma\right)}{\mathrm{BR}\left(\ell_{\alpha}\to\ell_{\beta}\nu_{\alpha}\overline{\nu}_{\beta}\right)}\simeq\frac{3\alpha_{\mathrm{EM}}}{64\pi
G^{2}_{F}}\times$ $\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\bigg{\\{}\left|\frac{Y_{1,\beta}^{*}Y_{1,\alpha}}{m^{2}_{\phi_{1}}}f\left(\frac{m_{F}^{2}}{m^{2}_{\phi_{1}}}\right)+\frac{Y_{2,\beta}^{*}Y_{2,\alpha}}{m^{2}_{\phi_{2}}}f\left(\frac{m_{F}^{2}}{m^{2}_{\phi_{2}}}\right)\right|^{2}\bigg{\\}},$
(43)

with the auxiliary function

$\displaystyle f(x)=\frac{1-6x+3x^{2}+2x^{3}-6x^{2}\ln x}{6(1-x)^{4}}\,.$ (44)

The current bounds (and future sensitivities) on the relevant branching ratios
at $90\text{\,}\mathrm{\char 37\relax}$ CL are given by

$\displaystyle\mathrm{BR}(\mu\to e\gamma)<3.1\times 10^{-13}\,\,(6.0\times
10^{-14})\,\,\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{MEGII:2023ltw}{\@@citephrase{(}}{\@@citephrase{)}}}},$
$\displaystyle\mathrm{BR}(\tau\to e\gamma)<3.3\times
10^{-8}\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{BaBar:2009hkt}{\@@citephrase{(}}{\@@citephrase{)}}}}\,\,(3\times
10^{-9}\,\mbox{\cite[cite]{\@@bibref{Authors Phrase1YearPhrase2}{Belle-
II:2018jsg}{\@@citephrase{(}}{\@@citephrase{)}}}}),$
$\displaystyle\mathrm{BR}(\tau\to\mu\gamma)<4.2\times
10^{-8}\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{Belle:2021ysv}{\@@citephrase{(}}{\@@citephrase{)}}}}\,\,(10^{-9}\,\mbox{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{Belle-II:2018jsg}{\@@citephrase{(}}{\@@citephrase{)}}}}),$
(45)

and impose stringent limits on the Yukawa couplings.

Taking into account all these constraints, we explore the parameter space that
can explain the correct relic abundance for DM. For this, we perform a random
scan over the ranges
$m_{F}<m_{S_{1}}\in(10\,\,\mathrm{MeV}-10\,\,\mathrm{GeV})$ and
$m_{S_{2}}=m_{\phi_{1}}=m_{\phi_{2}}\in(110\,\,\mathrm{GeV}-500\,\,\mathrm{GeV})$.
The Yukawa couplings are fixed through the parametrization given in Eq. (19)
for NH and Eq. (27) for IH. The entries of the PMNS matrix are determined by
varying over the $3\sigma$ range of neutrino observables taken from NuFIT 5.3
Esteban:2020cvm  under the assumption of flat priors. We also fix the
hierarchy of these couplings to be $|Y_{1}|>|Y_{2}|$. Since one requires
sizeable Yukawa couplings to generate a sufficiently large contribution to
$\sigma v$ Herms:2023cyy , we vary $k$ in the range $(10^{2},10^{8})$. From
this data sample, we select the data points that satisfy the DM relic density
constraint Planck:2018vyg . The results of this analysis are shown in the DM
mass ($m_{F}$) versus Yukawa coupling ($|Y_{1}|$) plane for both NH (left
panel) and IH (right panel) in Fig. 4. In this figure, the red-colored data
points are excluded by the cLFV constraints, whereas the green (blue)- colored
data points are allowed and can (not) be probed in future cLFV experiments.
Notice that, for larger DM masses, sizeable values of Yukawa couplings are
required to be consistent with the DM relic density constraint. However, such
large values of Yukawa couplings are excluded by the cLFV constraints for
$m_{F}\gtrsim 0.8\,\,\mathrm{GeV}$. Moreover, it is also worth mentioning that
these constraints are stronger for the IH of neutrino masses as already DM
masses $m_{F}\gtrsim 0.5\,\,\mathrm{GeV}$ are ruled out.

Figure 6: Left (right) panel: parameter space in $\Delta m_{\phi_{1}}\equiv
m_{\phi_{1}}-m_{F}$ vs. $|Y_{1}|$ $(m_{F}$ vs. $\Delta m_{\phi_{1}})$ plane
consistent with the Planck relic density constraint. Here, we set
$\phi_{1}^{1+\epsilon}$ to be the coannihilation partner of DM. The color
codes represent the same as in Fig. 4.

_Heavy dark matter_.– As discussed above, for larger DM masses, its
annihilation into SM leptons via the $t$-channel processes is excluded through
the cLFV constraints. However, the coannihilations with the new scalars are
less severely constrained by these constraints, as they only depend on the
square of the Yukawa couplings, while coannihilations of the new scalars into,
e.g., gauge bosons are independent of the Yukawas. Hence, for the heavy DM
scenario, the coannihilation processes can become important for explaining the
correct relic abundance. To generate a sufficiently large contribution to
$\sigma v$ through these processes, the mass splitting between the DM and the
coannihilation partner is required to be small Griest:1990kh . In our
framework, such small mass splittings are allowed. For the numerical study, we
consider two scenarios. In the first scenario, the milli-charged scalar
$S_{1}^{\epsilon}$ is considered to be the next-to-lightest new particle,
whereas in the second, $\phi_{1}^{1+\epsilon}$ is considered to be the next-
to-lightest particle. For the analysis, we use the following toolchain.
minimal-lagrangians May:2020sod  is used to create the model input for SARAH
4.15.1 Staub:2013tta  which in turn generates the SPheno 4.0.5 Porod:2011nf
module as well as the CalcHEP Belyaev:2012qa  model used to calculate the
particle spectrum and the cLFV observables (Eq. (IV)), the $S$, $T$ and $U$
parameters and the files for micrOMEGAs 5.3.41 Belanger:2018ccd  used to
calculate the relic density.

In Fig. 5 (left panel), we show the parameter space for the first scenario
consistent with the Planck relic density constraint Planck:2018vyg  in the
$\Delta m_{S_{1}}\equiv m_{S_{1}}-m_{F}$ versus $|Y_{1}|$ plane using the same
color codes as in Fig. 4 for the NH of neutrino masses. The associated scan is
over the ranges

$\displaystyle
m_{F}\in[$10\text{\,}\mathrm{GeV}$,${10}^{5}\text{\,}\mathrm{GeV}$]\,,$
$\displaystyle
m_{S_{1}}-m_{F}<$100\text{\,}\mathrm{GeV}$\,\mathrm{with}\,m_{S_{2}}\in[m_{S_{1}},${10}^{5}\text{\,}\mathrm{GeV}$]\,,$
$\displaystyle
m_{\phi_{1}},m_{\phi_{2}}\in[\max(m_{S_{1}},$110\text{\,}\mathrm{GeV}$),${10}^{5}\text{\,}\mathrm{GeV}$]\,,$
$\displaystyle\lambda_{H\phi_{i}}\in[\max(-\lambda^{\prime}_{H\phi_{i}},10^{-10}),\sqrt{4\pi}]\,,$
$\displaystyle\lambda_{H\phi_{1}\phi_{2}}\in[10^{-8},\sqrt{4\pi}]\,\mathrm{and}\,k\in[1,10^{10}]$
(46)

under the assumption of logarithmic priors and using the three-flavor neutrino
oscillation parameters taken from the global fit of the NuFIT 5.2
collaboration Esteban:2020cvm . The vacuum stability conditions in Eqs. (40)
to (42) are enforced as well. Notice that, for the larger mass splittings,
larger values of the Yukawa couplings are required to generate the correct
relic abundance, which are severely constrained by the cLFV constraints. On
the other hand, for lower mass splittings, the Yukawa couplings are allowed to
be much smaller since the coannihilation processes become the dominant
annihilation mode in this regime. In the figure shown in the right panel, we
show the parameter space for the same scenario in the DM mass versus $\Delta
m_{S_{1}}$ plane. Notice that the region where $m_{F}\gtrsim
40\,\,\mathrm{GeV}$ is consistent with all the constraints including the cLFV
bounds. However, DM masses below $40\,\,\mathrm{GeV}$ are not favored from
these constraints, since both coannihilation and annihilation of
$S_{1}^{\epsilon}$ are either suppressed or kinematically forbidden in this
region.

In Fig. 6, the parameter space for the second scenario involving dominantly
$\phi_{1}^{1+\epsilon}$ coannihilations is shown. The corresponding scan
ranges are identical to those in Eq. (46) except with the roles of
$S_{1}^{\epsilon}$ and $\phi_{1}^{1+\epsilon}$ being interchanged, i.e.,

$\displaystyle m_{\phi_{1}}-\max(110,m_{F})<$110\text{\,}\mathrm{GeV}$$
$\displaystyle m_{S_{1}}\in[m_{\phi_{1}},${10}^{5}\text{\,}\mathrm{GeV}$]$
(47)

As one can see from this figure, DM masses below 100 GeV are excluded by the
cLFV constraints since the coannihilation processes are suppressed in this
region because of $m_{\phi_{1}}\gtrsim 110\,\,\mathrm{GeV}$.

##  V Implications on Neutrinoless double beta decay

Our framework predicts neutrinoless double beta decay at an observable level.
To demonstrate this, first, we calculate the contribution to the effective
neutrino mass $m_{\beta\beta}$ in our model. Since the lightest neutrino is
massless in our setup, $m_{\beta\beta}$ takes the following form Herms:2023cyy
:

$\displaystyle m_{\beta\beta}^{NH}$
$\displaystyle=\left|m_{2}s_{12}^{2}c_{13}^{2}+m_{3}s_{13}^{2}e^{i\alpha}\right|,$
$\displaystyle m_{\beta\beta}^{IH}$
$\displaystyle=\left|m_{1}c_{12}^{2}c_{13}^{2}+m_{2}s_{12}^{2}c_{13}^{2}e^{i\alpha}\right|,$
(48)

where $s_{ij}\equiv\sin{\theta_{ij}}$ ($c_{ij}\equiv\cos{\theta_{ij}}$), and
$\theta_{ij}$ and $\alpha$ denote the neutrino mixing angles and Majorana
phase, respectively. As a consequence of neutrino masses being generated
radiatively in our framework, $m_{\beta\beta}$ given in Eq. (48) will be
modified by the loop effect Rodejohann:2019quz . In the limit where the scale
of neutrino self-energy (denoted as $\Lambda$) is greater than the momentum
transfer of the $0\nu\beta\beta$ process (denoted as $p$), the modification of
$m_{\beta\beta}$ is given by Rodejohann:2019quz ; Herms:2023cyy

$m_{\beta\beta}\longrightarrow
m_{\beta\beta}\left(1+\frac{p^{2}}{\Lambda^{2}}\right),$ (49)

where $\Lambda$ is given by Herms:2023cyy

$\displaystyle\Lambda^{2}=2\frac{\mathcal{G}[m^{2}_{S_{2}}]-\mathcal{G}[m^{2}_{S_{1}}]}{\mathcal{H}[m^{2}_{S_{2}}]-\mathcal{H}[m^{2}_{S_{1}}]}.$
(50)

Here, the loop functions are defined as

$\displaystyle\mathcal{G}[m^{2}_{S}]\equiv\frac{m_{S_{i}}^{2}\ln\left(\frac{m_{S}^{2}}{m^{2}_{F}}\right)}{m_{S}^{2}-m_{F}^{2}},$
(51)
$\displaystyle\mathcal{H}[m^{2}_{S}]\equiv\frac{m^{4}_{F}-m^{4}_{S}+2m^{2}_{S}m^{2}_{F}\ln{\dfrac{m^{2}_{S}}{m^{2}_{F}}}}{(m^{2}_{S}-m^{2}_{F})^{3}}.$
(52)

Notice from Eq. (49) that the loop effect becomes significant when
$\Lambda\simeq\mathcal{O}(p)$. This implies that the loop effect is more
significant for the light DM scenario compared to the heavy DM case since the
typical scale of momentum transfer of $0\nu\beta\beta$ process is
$\mathcal{O}(100)\,\,\mathrm{MeV}$ Simkovic:2018hiq ; Shimizu:2017qcy ;
Graf:2023dzf .

In Fig. 7, we show the parameter space in DM mass versus the effective
neutrino mass plane consistent with the Planck relic density constraint and
the cLFV constraints. For concreteness, we fix $\alpha=0$ and
$p=200\,\,\mathrm{MeV}$. Notice that from this figure, the contribution to
$m_{\beta\beta}$ is enhanced for the lower DM masses, whereas for the higher
masses, the effect is negligible.

Figure 7: Parameter space in DM mass ($m_{F}$) vs. effective neutrino mass
($m_{\beta\beta}$) for both NH and IH of neutrino masses. For concreteness, we
fix $p=200\,\mathrm{MeV}$ and $\alpha=0$. The blue and green color codes
denote the same as in Fig. 4.

##  VI Conclusions

In this work, we have presented a radiative neutrino mass mechanism in which
the particles within the loops carry small electric charges. Contrary to the
conventional scotogenic setup, this scheme does not require the ad-hoc
imposition of any new symmetry for stabilizing the DM candidate since its
stability is ensured through its minuscule electric charge. The minimal UV
completion of this scenario requires one generation of a milli-charged fermion
$F$ and two EW-scalar doublets $\phi_{1}$ and $\phi_{2}$ in addition to the SM
particles.

We have investigated the DM phenomenology associated with the lightest milli-
charged particle $F$, for both freeze-in and freeze-out scenarios finding that
the freeze-in mechanism is not compatible with the neutrino mass mechanism in
our model while the freeze-out mechanism works successfully. We explored the
thermal relic scenario from sub-GeV to TeV DM masses. In the sub-GeV mass
regime, the DM candidate mainly self-annihilates into neutrinos via the
leptonic portal whereas, for the heavy DM scenario, coannihilation processes
set the relic abundance.

We have also discussed various other phenomenological implications of our
model, such as lepton flavor-violating observables, electroweak precision
observables, and neutrinoless double beta decay. We found that a vast region
of parameter space can be probed in future cLFV experiments. Additionally, our
model also predicts neutrinoless double beta decay at an observable level.

###### Acknowledgements.

Acknowledgments.– We would like to thank A. Pukhov for correspondence on the
implementation of milli-charged DM into micrOMEGAs. This work has been
supported by the DFG through the Research Training Network 2149 “Strong and
Weak Interactions - from Hadrons to Dark Matter”.

## References

  * (1) M. Marinelli and G. Morpurgo, “The Electric Neutrality of Matter: A Summary,” Phys. Lett. B 137 (1984) 439–442.
  * (2) H. F. Dylla and J. G. King, “Neutrality of Molecules by a New Method,” Phys. Rev. A 7 (1973) 1224–1229.
  * (3) R. Gahler, J. Kalus, and W. Mampe, “EXPERIMENTAL LIMIT FOR THE CHARGE OF THE FREE NEUTRON,” Phys. Rev. D 25 (1982) 2887–2894.
  * (4) A. V. Kyuldjiev, “Searching for Effects of Neutrino Magnetic Moments at Reactors and Accelerators,” Nucl. Phys. B 243 (1984) 387–397.
  * (5) S. Davidson, S. Hannestad, and G. Raffelt, “Updated bounds on millicharged particles,” JHEP 05 (2000) 003, arXiv:hep-ph/0001179.
  * (6) B. Holdom, “Two U(1)’s and Epsilon Charge Shifts,” Phys. Lett. B 166 (1986) 196–198.
  * (7) R. Foot, G. C. Joshi, H. Lew, and R. R. Volkas, “Charge quantization in the standard model and some of its extensions,” Mod. Phys. Lett. A 5 (1990) 2721–2732.
  * (8) K. S. Babu and R. R. Volkas, “Bounds on minicharged neutrinos in the minimal Standard Model,” Phys. Rev. D 46 (1992) R2764–R2768, arXiv:hep-ph/9208260.
  * (9) G. Bertone, D. Hooper, and J. Silk, “Particle dark matter: Evidence, candidates and constraints,” Phys. Rept. 405 (2005) 279–390, arXiv:hep-ph/0404175.
  * (10) M. Klasen, M. Pohl, and G. Sigl, “Indirect and direct search for dark matter,” Prog. Part. Nucl. Phys. 85 (2015) 1–32, arXiv:1507.03800 [hep-ph].
  * (11) Super-Kamiokande Collaboration, Y. Fukuda et al., “Evidence for oscillation of atmospheric neutrinos,” Phys. Rev. Lett. 81 (1998) 1562–1567, arXiv:hep-ex/9807003.
  * (12) SNO Collaboration, Q. R. Ahmad et al., “Measurement of the rate of $\nu_{e}+d\to p+p+e^{-}$ interactions produced by 8B solar neutrinos at the Sudbury Neutrino Observatory,” Phys. Rev. Lett. 87 (2001) 071301, arXiv:nucl-ex/0106015.
  * (13) SNO Collaboration, Q. R. Ahmad et al., “Direct evidence for neutrino flavor transformation from neutral current interactions in the Sudbury Neutrino Observatory,” Phys. Rev. Lett. 89 (2002) 011301, arXiv:nucl-ex/0204008.
  * (14) Z.-j. Tao, “Radiative seesaw mechanism at weak scale,” Phys. Rev. D 54 (1996) 5693–5697, arXiv:hep-ph/9603309.
  * (15) E. Ma, “Verifiable radiative seesaw mechanism of neutrino mass and dark matter,” Phys. Rev. D 73 (2006) 077301, arXiv:hep-ph/0601225.
  * (16) L. M. Krauss, S. Nasri, and M. Trodden, “A Model for neutrino masses and dark matter,” Phys. Rev. D 67 (2003) 085002, arXiv:hep-ph/0210389.
  * (17) M. Aoki, S. Kanemura, and O. Seto, “Neutrino mass, Dark Matter and Baryon Asymmetry via TeV-Scale Physics without Fine-Tuning,” Phys. Rev. Lett. 102 (2009) 051805, arXiv:0807.0361 [hep-ph].
  * (18) E. Ma and D. Suematsu, “Fermion Triplet Dark Matter and Radiative Neutrino Mass,” Mod. Phys. Lett. A 24 (2009) 583–589, arXiv:0809.0942 [hep-ph].
  * (19) D. Restrepo, O. Zapata, and C. E. Yaguna, “Models with radiative neutrino masses and viable dark matter candidates,” JHEP 11 (2013) 011, arXiv:1308.3655 [hep-ph].
  * (20) M. Gustafsson, J. M. No, and M. A. Rivera, “Predictive Model for Radiatively Induced Neutrino Masses and Mixings with Dark Matter,” Phys. Rev. Lett. 110 no. 21, (2013) 211802, arXiv:1212.4806 [hep-ph]. [Erratum: Phys.Rev.Lett. 112, 259902 (2014)].
  * (21) J. Fiaschi, M. Klasen, and S. May, “Singlet-doublet fermion and triplet scalar dark matter with radiative neutrino masses,” JHEP 05 (2019) 015, arXiv:1812.11133 [hep-ph].
  * (22) S. Esch, M. Klasen, and C. E. Yaguna, “A singlet doublet dark matter model with radiative neutrino masses,” JHEP 10 (2018) 055, arXiv:1804.03384 [hep-ph].
  * (23) S. Jana, P. K. Vishnu, W. Rodejohann, and S. Saad, “Dark matter assisted lepton anomalous magnetic moments and neutrino masses,” Phys. Rev. D 102 no. 7, (2020) 075003, arXiv:2008.02377 [hep-ph].
  * (24) E. Ma, I. Picek, and B. Radovčić, “New Scotogenic Model of Neutrino Mass with $U(1)_{D}$ Gauge Interaction,” Phys. Lett. B 726 (2013) 744–746, arXiv:1308.5313 [hep-ph].
  * (25) S.-Y. Ho, T. Toma, and K. Tsumura, “Systematic $U(1)_{B–L}$ extensions of loop-induced neutrino mass models with dark matter,” Phys. Rev. D 94 no. 3, (2016) 033007, arXiv:1604.07894 [hep-ph].
  * (26) T. de Boer, M. Klasen, and S. Zeinstra, “Anomaly-free dark matter models with one-loop neutrino masses and a gauged U(1) symmetry,” JHEP 01 (2024) 013, arXiv:2309.06920 [hep-ph].
  * (27) S. Jana, P. K. Vishnu, and S. Saad, “Minimal dirac neutrino mass models from $\hbox{U}(1)_{\mathrm{R}}$ gauge symmetry and left–right asymmetry at colliders,” Eur. Phys. J. C 79 no. 11, (2019) 916, arXiv:1904.07407 [hep-ph].
  * (28) S. Jana, P. K. Vishnu, and S. Saad, “Minimal realizations of Dirac neutrino mass from generic one-loop and two-loop topologies at $d=5$,” JCAP 04 (2020) 018, arXiv:1910.09537 [hep-ph].
  * (29) S. D. McDermott, H.-B. Yu, and K. M. Zurek, “Turning off the Lights: How Dark is Dark Matter?,” Phys. Rev. D 83 (2011) 063509, arXiv:1011.2907 [hep-ph].
  * (30) C. Dvorkin, T. Lin, and K. Schutz, “Making dark matter out of light: freeze-in from plasma effects,” Phys. Rev. D 99 no. 11, (2019) 115009, arXiv:1902.08623 [hep-ph]. [Erratum: Phys.Rev.D 105, 119901 (2022)].
  * (31) I. Cordero-Carrión, M. Hirsch, and A. Vicente, “Master Majorana neutrino mass parametrization,” Phys. Rev. D 99 no. 7, (2019) 075019, arXiv:1812.03896 [hep-ph].
  * (32) I. Cordero-Carrión, M. Hirsch, and A. Vicente, “General parametrization of Majorana neutrino mass models,” Phys. Rev. D 101 no. 7, (2020) 075032, arXiv:1912.08858 [hep-ph].
  * (33) C. Hagedorn, J. Herrero-García, E. Molinaro, and M. A. Schmidt, “Phenomenology of the Generalised Scotogenic Model with Fermionic Dark Matter,” JHEP 11 (2018) 103, arXiv:1804.04117 [hep-ph].
  * (34) W. L. Xu, C. Dvorkin, and A. Chael, “Probing sub-GeV Dark Matter-Baryon Scattering with Cosmological Observables,” Phys. Rev. D 97 no. 10, (2018) 103530, arXiv:1802.06788 [astro-ph.CO].
  * (35) Planck Collaboration, P. A. R. Ade et al., “Planck 2015 results. XIII. Cosmological parameters,” Astron. Astrophys. 594 (2016) A13, arXiv:1502.01589 [astro-ph.CO].
  * (36) SDSS Collaboration, P. McDonald et al., “The Lyman-alpha forest power spectrum from the Sloan Digital Sky Survey,” Astrophys. J. Suppl. 163 (2006) 80–109, arXiv:astro-ph/0405013.
  * (37) M. Ibe and H.-b. Yu, “Distinguishing Dark Matter Annihilation Enhancement Scenarios via Halo Shapes,” Phys. Lett. B 692 (2010) 70–73, arXiv:0912.5425 [hep-ph].
  * (38) M. Markevitch, A. H. Gonzalez, D. Clowe, A. Vikhlinin, L. David, W. Forman, C. Jones, S. Murray, and W. Tucker, “Direct constraints on the dark matter self-interaction cross-section from the merging galaxy cluster 1E0657-56,” Astrophys. J. 606 (2004) 819–824, arXiv:astro-ph/0309303.
  * (39) M. D. Citron, “Search for millicharged particles at the LHC with the milliQan prototype,” PoS ICHEP2020 (2021) 232.
  * (40) A. Ball et al., “Search for millicharged particles in proton-proton collisions at $\sqrt{s}=13$ TeV,” Phys. Rev. D 102 no. 3, (2020) 032002, arXiv:2005.06518 [hep-ex].
  * (41) A. A. Prinz et al., “Search for millicharged particles at SLAC,” Phys. Rev. Lett. 81 (1998) 1175–1178, arXiv:hep-ex/9804008.
  * (42) G. Magill, R. Plestid, M. Pospelov, and Y.-D. Tsai, “Millicharged particles in neutrino experiments,” Phys. Rev. Lett. 122 no. 7, (2019) 071801, arXiv:1806.03310 [hep-ph].
  * (43) LSND Collaboration, C. Athanassopoulos et al., “The Liquid scintillator neutrino detector and LAMPF neutrino source,” Nucl. Instrum. Meth. A 388 (1997) 149–172, arXiv:nucl-ex/9605002.
  * (44) MiniBooNE Collaboration, A. A. Aguilar-Arevalo et al., “Significant Excess of ElectronLike Events in the MiniBooNE Short-Baseline Neutrino Experiment,” Phys. Rev. Lett. 121 no. 22, (2018) 221801, arXiv:1805.12028 [hep-ex].
  * (45) ArgoNeuT Collaboration, R. Acciarri et al., “Improved Limits on Millicharged Particles Using the ArgoNeuT Experiment at Fermilab,” Phys. Rev. Lett. 124 no. 13, (2020) 131801, arXiv:1911.07996 [hep-ex].
  * (46) DarkSide Collaboration, P. Agnes et al., “Low-Mass Dark Matter Search with the DarkSide-50 Experiment,” Phys. Rev. Lett. 121 no. 8, (2018) 081307, arXiv:1802.06994 [astro-ph.HE].
  * (47) G. Belanger, A. Mjallal, and A. Pukhov, “Recasting direct detection limits within micrOMEGAs and implication for non-standard Dark Matter scenarios,” Eur. Phys. J. C 81 no. 3, (2021) 239, arXiv:2003.08621 [hep-ph].
  * (48) XENON Collaboration, E. Aprile et al., “Dark Matter Search Results from a One Ton-Year Exposure of XENON1T,” Phys. Rev. Lett. 121 no. 11, (2018) 111302, arXiv:1805.12562 [astro-ph.CO].
  * (49) LZ Collaboration, J. Aalbers et al., “First Dark Matter Search Results from the LUX-ZEPLIN (LZ) Experiment,” Phys. Rev. Lett. 131 no. 4, (2023) 041002, arXiv:2207.03764 [hep-ex].
  * (50) G. Bélanger, F. Boudjema, A. Goudelis, A. Pukhov, and B. Zaldivar, “micrOMEGAs5.0 : Freeze-in,” Comput. Phys. Commun. 231 (2018) 173–186, arXiv:1801.03509 [hep-ph].
  * (51) Planck Collaboration, N. Aghanim et al., “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641 (2020) A6, arXiv:1807.06209 [astro-ph.CO]. [Erratum: Astron.Astrophys. 652, C4 (2021)].
  * (52) J. Herms, S. Jana, V. P. K., and S. Saad, “Light neutrinophilic dark matter from a scotogenic model,” Phys. Lett. B 845 (2023) 138167, arXiv:2307.15760 [hep-ph].
  * (53) S. Okawa and Y. Omura, “Light mass window of lepton portal dark matter,” JHEP 02 (2021) 231, arXiv:2011.04788 [hep-ph].
  * (54) K. S. Babu, P. S. B. Dev, S. Jana, and A. Thapa, “Non-Standard Interactions in Radiative Neutrino Mass Models,” JHEP 03 (2020) 006, arXiv:1907.09498 [hep-ph].
  * (55) S. Iguro, S. Okawa, and Y. Omura, “Light lepton portal dark matter meets the LHC,” JHEP 03 (2023) 010, arXiv:2208.05487 [hep-ph].
  * (56) S. Jana, V. P. K., and S. Saad, “Resolving electron and muon $g-2$ within the 2HDM,” Phys. Rev. D 101 no. 11, (2020) 115037, arXiv:2003.03386 [hep-ph].
  * (57) LEP, ALEPH, DELPHI, L3, OPAL, LEP Electroweak Working Group, SLD Electroweak Group, SLD Heavy Flavor Group Collaboration, t. S. Electroweak, “A Combination of preliminary electroweak measurements and constraints on the standard model,” arXiv:hep-ex/0312023.
  * (58) Particle Data Group Collaboration, M. Tanabashi et al., “Review of Particle Physics,” Phys. Rev. D 98 no. 3, (2018) 030001.
  * (59) M. E. Peskin and T. Takeuchi, “A New constraint on a strongly interacting Higgs sector,” Phys. Rev. Lett. 65 (1990) 964–967.
  * (60) M. E. Peskin and T. Takeuchi, “Estimation of oblique electroweak corrections,” Phys. Rev. D 46 (1992) 381–409.
  * (61) ATLAS Collaboration, G. Aad et al., “Search for invisible Higgs-boson decays in events with vector-boson fusion signatures using 139 fb-1 of proton-proton data recorded by the ATLAS experiment,” JHEP 08 (2022) 104, arXiv:2202.07953 [hep-ex].
  * (62) CMS Collaboration, A. Tumasyan et al., “Search for invisible decays of the Higgs boson produced via vector boson fusion in proton-proton collisions at s=13 TeV,” Phys. Rev. D 105 no. 9, (2022) 092007, arXiv:2201.11585 [hep-ex].
  * (63) ATLAS Collaboration, “Measurement of the properties of Higgs boson production at $\sqrt{s}=13$ TeV in the $H\to\gamma\gamma$ channel using $139$ fb-1 of $pp$ collision data with the ATLAS experiment,” arXiv:2207.00348 [hep-ex].
  * (64) CMS Collaboration, “Measurement of the Higgs boson inclusive and differential fiducial production cross sections in the diphoton decay channel with pp collisions at $\sqrt{s}$ = 13 TeV,” arXiv:2208.12279 [hep-ex].
  * (65) MEG II Collaboration, K. Afanaciev et al., “A search for $\mu^{+}\rightarrow\textrm{e}^{+}\gamma$ with the first dataset of the MEG II experiment,” Eur. Phys. J. C 84 no. 3, (2024) 216, arXiv:2310.12614 [hep-ex].
  * (66) BaBar Collaboration, B. Aubert et al., “Searches for Lepton Flavor Violation in the Decays tau+- —$>$ e+- gamma and tau+- —$>$ mu+- gamma,” Phys. Rev. Lett. 104 (2010) 021802, arXiv:0908.2381 [hep-ex].
  * (67) Belle-II Collaboration, W. Altmannshofer et al., “The Belle II Physics Book,” PTEP 2019 no. 12, (2019) 123C01, arXiv:1808.10567 [hep-ex]. [Erratum: PTEP 2020, 029201 (2020)].
  * (68) Belle Collaboration, A. Abdesselam et al., “Search for lepton-flavor-violating tau-lepton decays to $\ell\gamma$ at Belle,” JHEP 10 (2021) 19, arXiv:2103.12994 [hep-ex].
  * (69) I. Esteban, M. C. Gonzalez-Garcia, M. Maltoni, T. Schwetz, and A. Zhou, “The fate of hints: updated global analysis of three-flavor neutrino oscillations,” JHEP 09 (2020) 178, arXiv:2007.14792 [hep-ph].
  * (70) K. Griest and D. Seckel, “Three exceptions in the calculation of relic abundances,” Phys. Rev. D 43 (1991) 3191–3203.
  * (71) S. May, “minimal-lagrangians: Generating and studying dark matter model Lagrangians with just the particle content,” Comput. Phys. Commun. 261 (2021) 107773, arXiv:2003.08037 [hep-ph].
  * (72) F. Staub, “SARAH 4 : A tool for (not only SUSY) model builders,” Comput. Phys. Commun. 185 (2014) 1773–1790, arXiv:1309.7223 [hep-ph].
  * (73) W. Porod and F. Staub, “SPheno 3.1: Extensions including flavour, CP-phases and models beyond the MSSM,” Comput. Phys. Commun. 183 (2012) 2458–2469, arXiv:1104.1573 [hep-ph].
  * (74) A. Belyaev, N. D. Christensen, and A. Pukhov, “CalcHEP 3.4 for collider physics within and beyond the Standard Model,” Comput. Phys. Commun. 184 (2013) 1729–1769, arXiv:1207.6082 [hep-ph].
  * (75) W. Rodejohann and X.-J. Xu, “Loop-enhanced rate of neutrinoless double beta decay,” JHEP 11 (2019) 029, arXiv:1907.12478 [hep-ph].
  * (76) F. Šimkovic, A. Smetana, and P. Vogel, “$0\nu\beta\beta$ nuclear matrix elements, neutrino potentials and $\mathrm{SU}(4)$ symmetry,” Phys. Rev. C 98 no. 6, (2018) 064325, arXiv:1808.05016 [nucl-th].
  * (77) N. Shimizu, J. Menéndez, and K. Yako, “Double Gamow-Teller Transitions and its Relation to Neutrinoless $\beta$$\beta$ Decay,” Phys. Rev. Lett. 120 no. 14, (2018) 142502, arXiv:1709.01088 [nucl-th].
  * (78) L. Gráf, S. Jana, O. Scholer, and N. Volmer, “Neutrinoless Double Beta Decay without Vacuum Majorana Neutrino Mass,” arXiv:2312.15016 [hep-ph].