Lepton flavor violating $Z$ and Higgs decays in the scotogenic model

Raghavendra Srikanth Hundi

Department of Physics, Indian Institute of Technology Hyderabad,

Kandi - 502 284, India.

E-mail address<EMAIL_ADDRESS>
###### Abstract

In this work, we have studied lepton flavor violating (LFV) decays of $Z$
gauge boson and Higgs boson ($H$) in the scotogenic model. We have computed
branching ratios for the decays $Z\to\ell_{\alpha}\ell_{\beta}$ and
$H\to\ell_{\alpha}\ell_{\beta}$ in this model. Here, $\ell_{\alpha}$ and
$\ell_{\beta}$ are different charged lepton fields. After fitting to the
neutrino oscillation observables in the scotogenic model, we have found that
the branching ratios for the LFV decays of $Z$ and $H$ can be as large as
$\sim 10^{-8}$ and $\sim 10^{-3}$ respectively. However, after satisfying the
constraints due to non-observation of $\ell_{\alpha}\to\ell_{\beta}\gamma$
decays, the above mentioned branching ratio results are found to be suppressed
by a factor of $\sim 10^{-7}$.

##  1 Introduction

Physics beyond the standard model [1] can be probed by searching for lepton
flavor violating (LFV) [2] processes in experiments. So far no LFV signal has
been observed in experiments, and as result, upper bounds exist on various LFV
processes [3]. In the standard model these experimental limits are satisfied,
since LFV processes are highly suppressed due to Glashow-Iliopoulos-Maiani
cancellation mechanism. On the other hand, in a beyond standard model, the
branching ratios for these processes can be appreciably large and the model
can be constrained by experiments.

Scotogenic model [4] is an extension of standard model, which explains the
neutrino mass and dark matter problems, which are briefly described below.
Neutrino masses are found to be tiny [5], and hence, in order to explain the
smallness of neutrino masses a different mechanism should be proposed for them
[6]. Regarding the dark matter problem, it is known that the universe consists
of nearly 25$\%$ of energy in the form of non-baryonic matter [7], which
cannot be explained by the standard model. In the scotogenic model, the origin
of neutrino masses are explained by a radiative mechanism by proposing an
extra scalar doublet ($\eta$), three right-handed Majorana neutrinos ($N_{k}$)
and an additional $Z_{2}$ symmetry. Under $Z_{2}$ symmetry, which is unbroken,
$\eta$ and $N_{k}$ are odd and all the standard model fields are even. As a
result of this, the lightest among the neutral $Z_{2}$-odd particles can be a
candidate for the dark matter.

Various phenomenological consequences of scotogenic model have been studied in
relation to LFV, dark matter, matter-antimatter asymmetry and colliders [8, 9,
10]. In the studies on LFV in the scotogenic model, the following processes
have been analyzed: $\ell_{\alpha}\to\ell_{\beta}\gamma$, $\ell_{\alpha}\to
3\ell_{\beta}$ and conversion of $\mu$ to $e$ [8, 9]. In a related direction,
see Ref. [11], for a study on LFV in the supersymmetric scotogenic model [12].
In contrast to above mentioned studies, in this work, we analyze the LFV
decays of $Z$ and Higgs boson in the scotogenic model [4]. The decays
$Z\to\ell_{\alpha}\ell_{\beta}$ and $H\to\ell_{\alpha}\ell_{\beta}$ are driven
at 1-loop level by $\eta^{\pm}$ and $N_{k}$, where $\eta^{\pm}$ is the charged
component of $\eta$. We compute branching ratios for these decays, which we
find to be dependent on the Yukawa couplings and masses of $\eta^{\pm}$ and
$N_{k}$. By varying the parameters of the model, we study on the reach of the
values of the above mentioned branching ratios.

The current experimental bounds on the branching ratios of
$Z\to\ell_{\alpha}\ell_{\beta}$ and $H\to\ell_{\alpha}\ell_{\beta}$ are as
follows.

$\displaystyle{\rm Br}(Z\to e\mu)$ $\displaystyle<$ $\displaystyle 7.5\times
10^{-7}~{}\cite[cite]{[\@@bibref{}{zemu}{}{}]},$ $\displaystyle{\rm Br}(Z\to
e\tau)$ $\displaystyle<$ $\displaystyle 9.8\times
10^{-6}~{}\cite[cite]{[\@@bibref{}{zetau}{}{}]},$ $\displaystyle{\rm
Br}(Z\to\mu\tau)$ $\displaystyle<$ $\displaystyle 1.2\times
10^{-5}~{}\cite[cite]{[\@@bibref{}{zmutau}{}{}]}.$ (1) $\displaystyle{\rm
Br}(H\to e\mu)$ $\displaystyle<$ $\displaystyle 6.1\times
10^{-5}~{}\cite[cite]{[\@@bibref{}{hemu}{}{}]},$ $\displaystyle{\rm Br}(H\to
e\tau)$ $\displaystyle<$ $\displaystyle 4.7\times
10^{-3}~{}\cite[cite]{[\@@bibref{}{hetau}{}{}]},$ $\displaystyle{\rm
Br}(H\to\mu\tau)$ $\displaystyle<$ $\displaystyle 2.5\times
10^{-3}~{}\cite[cite]{[\@@bibref{}{hmutau}{}{}]}.$ (2)

In future, LFV decays of $Z$ and $H$ will be probed. For instance, in the
upcoming $e^{+}e^{-}$ collider such as the FCC-ee, the following sensitivites
can be probed for the LFV decays of $Z$ [19].

$\displaystyle{\rm Br}(Z\to e\mu)$ $\displaystyle\sim$ $\displaystyle
10^{-10}-10^{-8},$ $\displaystyle{\rm Br}(Z\to e\tau)$ $\displaystyle\sim$
$\displaystyle 10^{-9},$ $\displaystyle{\rm Br}(Z\to\mu\tau)$
$\displaystyle\sim$ $\displaystyle 10^{-9}.$ (3)

Similarly, the bounds on LFV decays of Higgs boson, given in Eq. (2), may be
reduced in future by the LHC. Since in future experiments, there is an
interest to probe LFV decays of $Z$ and $H$, it is worth to compute the
branching ratios of these decays in the scotogenic model. It is also
interesting to analyze the status of the above mentioned decays in this model,
in relation to the present and future bounds on them.

As already stated, the LFV decays of $Z$ and $H$ are mediated at 1-loop by
$\eta^{\pm},N_{k}$ in the scotogenic model. The same mediating particles, in
this model, can also drive $\ell_{\alpha}\to\ell_{\beta}\gamma$ at 1-loop
level. As a result of this, there exist a correlation between branching ratios
of $Z,H\to\ell_{\alpha}\ell_{\beta}$ and that of
$\ell_{\alpha}\to\ell_{\beta}\gamma$. Since stringent bounds exists on the
non-observation of $\ell_{\alpha}\to\ell_{\beta}\gamma$ [3], we have studied
the implications of those bounds on the branching ratios of
$Z,H\to\ell_{\alpha}\ell_{\beta}$. For related studies on LFV decays of $Z$
and $H$, see Refs.[20, 21].

The neutrino masses in the scotogenic model are generated at 1-loop level
through the mediation of neutral components of $\eta$ and $N_{k}$. As a result
of this, neutrino masses in this model depend on neutrino Yukawa couplings and
masses of neutral components of $\eta$ and $N_{k}$. As already stated before,
the branching ratios for $Z,H\to\ell_{\alpha}\ell_{\beta}$ also depend on the
neutrino Yukawa couplings and masses of $\eta^{\pm}$ and $N_{k}$. One can
notice that there exist a correlation between branching ratios of
$Z,H\to\ell_{\alpha}\ell_{\beta}$ and neutrino masses and mixing angles. We
have explored this correlation and we have found that the branching ratios of
$Z\to\ell_{\alpha}\ell_{\beta}$ can reach as high as $10^{-8}$ by satisfying
the perturbativity limits on the parameters of the scotogenic model. On the
other hand, the branching ratios for $H\to\ell_{\alpha}\ell_{\beta}$ can reach
as high as $10^{-3}$. However, the above mentioned results are obtained
without imposing the constraints due to non-observation of
$\ell_{\alpha}\to\ell_{\beta}\gamma$. After imposing the constraints due to
$\ell_{\alpha}\to\ell_{\beta}\gamma$, we have found that the above mentioned
results on the branching ratios are suppressed by a factor of $10^{-7}$. As a
result of this, the decay $H\to\mu\tau$ is found to have the highest branching
ratio of $\sim 10^{-10}$, in our analysis on the LFV decays of $Z$ and $H$.

In this work, although we study LFV decays of both $Z$ and $H$, only the LFV
decays of $H$ have been studied in Ref. [22]. Our method of computing the
branching ratios for LFV decays of $H$ is different from that of Ref. [22].
Moreover, only an estimation on the branching ratio of $H\to\mu\tau$ has been
made in Ref. [22], in the context of scotogenic model. Whereas, we have
studied branching ratios for all LFV Higgs decays in more details here. We
compare our results with that of Ref. [22] towards the end of this paper. See
Ref. [23] for some discussion on LFV decays of $Z$ and $H$ in the context of
generalized scotogenic model.

The paper is organized as follows. In the next section, we briefly describe
the scotogenic model. In Sec. 3, we present analytic expressions on the
branching ratios of $Z\to\ell_{\alpha}\ell_{\beta}$ and
$H\to\ell_{\alpha}\ell_{\beta}$ in the scotogenic model. In Sec. 4, we analyze
these branching ratios and present our numerical results on them. We conclude
in the last section.

##  2 Scotogenic model

Scotogenic model [4] is an extension of the standard model, where the
additional fields are one $SU(2)$ scalar doublet
$\eta=(\eta^{+},\eta^{0})^{T}$ and three singlet right-handed neutrinos
$N_{k}$. This model has an additional discrete $Z_{2}$ symmetry, under which
$\eta,N_{k}$ are odd and all the standard model fields are even. To construct
the invariant Lagrangian of this model, we can choose a basis where the Yukawa
couplings for charged leptons and the masses of right-handed neutrinos are
diagonal. In such a basis, the Lagrangian of this model in the lepton sector
is [4]

$-{\cal L}_{Y}=f_{\alpha}\bar{L}_{L\alpha}\phi\ell_{R\alpha}+h_{\alpha
k}\bar{L}_{L\alpha}\eta^{c}N_{k}+\frac{M_{k}}{2}\overline{N^{c}_{k}}N_{k}+h.c.$
(4)

Here, $\alpha=e,\mu,\tau$ and $k=1,2,3$.
$L_{L\alpha}=(\nu_{L\alpha},\ell_{L\alpha})^{T}$ is a left-handed lepton
doublet, $\ell_{R\alpha}$ is a right-handed singlet charged lepton,
$\phi=(\phi^{+},\phi^{0})^{T}$ is the scalar Higgs doublet and
$\eta^{c}=i\sigma_{2}\eta^{*}$, where $\sigma_{2}$ is a Pauli matrix. $\phi$
and $\eta$ are the only two scalar fields of this model. The scalar potential
between these two fields is given below [4].

$\displaystyle V$ $\displaystyle=$ $\displaystyle
m_{1}^{2}\phi^{\dagger}\phi+m_{2}^{2}\eta^{\dagger}\eta+\frac{1}{2}\lambda_{1}(\phi^{\dagger}\phi)^{2}+\frac{1}{2}\lambda_{2}(\eta^{\dagger}\eta)^{2}+\lambda_{3}(\phi^{\dagger}\phi)(\eta^{\dagger}\eta)+\lambda_{4}(\phi^{\dagger}\eta)(\eta^{\dagger}\phi)$
(5) $\displaystyle+\frac{1}{2}\lambda_{5}[(\phi^{\dagger}\eta)^{2}+h.c.].$

Here, $\lambda_{5}$ is chosen to be real, without loss of generality. Since
$Z_{2}$ is an exact symmetry of this model, we should have $m_{1}^{2}<0$ and
$m_{2}^{2}>0$ so that only $\phi$ acquires vacuum expectation value (VEV),
whereas $\eta$ does not acquire VEV. Since only $\phi$ acquires VEV, the
physical fields in the neutral components of $\phi$ and $\eta$ can be written
as

$\phi^{0}=\frac{H}{\sqrt{2}}+v,\quad\eta^{0}=\frac{1}{\sqrt{2}}(\eta_{R}+i\eta_{I})$
(6)

Here, $H$ is the Higgs boson and $v\approx$ 174 GeV. Now, after the
electroweak symmetry breaking, the physical components of $\phi$ and $\eta$
acquire masses, whose expressions in the form of mass-squares are given below
[4].

$\displaystyle m^{2}(H)\equiv m_{H}^{2}$ $\displaystyle=$ $\displaystyle
2\lambda_{1}v^{2},$ $\displaystyle m^{2}(\eta^{\pm})\equiv m_{\eta^{\pm}}^{2}$
$\displaystyle=$ $\displaystyle m_{2}^{2}+\lambda_{3}v^{2},$ $\displaystyle
m^{2}(\eta_{R})\equiv m_{R}^{2}$ $\displaystyle=$ $\displaystyle
m_{2}^{2}+(\lambda_{3}+\lambda_{4}+\lambda_{5})v^{2}=m_{0}^{2}+\lambda_{5}v^{2},$
$\displaystyle m^{2}(\eta_{I})\equiv m_{I}^{2}$ $\displaystyle=$
$\displaystyle
m_{2}^{2}+(\lambda_{3}+\lambda_{4}-\lambda_{5})v^{2}=m_{0}^{2}-\lambda_{5}v^{2}$
(7)

Here, $m_{0}^{2}=m_{2}^{2}+(\lambda_{3}+\lambda_{4})v^{2}$.

After the electroweak symmetry breaking, the first term of Eq. (4) give masses
to charged leptons, whose expressions can be written as

$m_{\ell_{\alpha}}=f_{\alpha}v$ (8)

On the other hand, since $\eta$ does not acquire VEV, the second term of Eq.
(4) do not generate Dirac masses for neutrinos. As a result of this, neutrinos
are massless at tree level. However, at 1-loop level, neutrinos acquire masses
through the mediation of neutral components of $\eta$ and $N_{k}$ [4]. By
taking $\Lambda={\rm diag}(\Lambda_{1},\Lambda_{2},\Lambda_{3})$, the mass
expressions for neutrinos at 1-loop level can be written as follows [4].

$\displaystyle(M_{\nu})_{\alpha\beta}$ $\displaystyle=$
$\displaystyle(h\Lambda h^{T})_{\alpha\beta}=\sum_{k=1}^{3}h_{\alpha
k}h_{\beta k}\Lambda_{k},$ $\displaystyle\Lambda_{k}$ $\displaystyle=$
$\displaystyle\frac{M_{k}}{16\pi^{2}}\left[\frac{m_{R}^{2}}{m_{R}^{2}-M_{k}^{2}}\ln\frac{m_{R}^{2}}{M_{k}^{2}}-\frac{m_{I}^{2}}{m_{I}^{2}-M_{k}^{2}}\ln\frac{m_{I}^{2}}{M_{k}^{2}}\right]$
(9)

Using the Casas-Ibarra parametrization [24], the matrix containing Yukawa
couplings $h_{\alpha k}$ can be parametrized as

$h=U_{PMNS}^{*}\sqrt{m_{\nu}}R\sqrt{\Lambda}^{-1}$ (10)

Here, $U_{PMNS}$ is the Pontecorvo-Maki-Nakagawa-Sakata matrix, which can be
parametrized [3] in terms of the three neutrino mixing angles, one $CP$
violating Dirac phase and two Majorana phases. $m_{\nu}$ is a diagonal matrix
containing the neutrino mass eigenvalues, which can be written as
$m_{\nu}={\rm diag}(m_{1},m_{2},m_{3})$. $R$ is a complex orthogonal matrix
which satisfies $RR^{T}=I=R^{T}R$. Using the parametrization of Eq. (10), one
can notice that

$M_{\nu}=U_{PMNS}^{*}m_{\nu}U_{PMNS}^{\dagger}$ (11)

From the above equation, we can see that the unitary matrix which diagonalize
$M_{\nu}$ is $U_{PMNS}$. Hence, the mixing pattern in the neutrino sector of
the scotogenic model can be explained by parametrizing the Yukawa couplings as
given by Eq. (10).

As described in Sec. 1, the aim of this work is to analyze LFV decays of $Z$
and $H$. One can notice that the LFV processes in the scotogenic model are
driven by the off-diagonal Yukawa couplings of the second term of Eq. (4). In
the next section, we explicitly show that the branching ratios of the LFV
decays for $Z$ and $H$ are proportional to off-diagonal elements of $h_{\alpha
k}$. As a result of this, the above mentioned branching ratios are
unsuppressed if $h_{\alpha k}\sim 1$. On the other hand, $h_{\alpha k}$ also
determine neutrino masses from Eq. (9). As already pointed in Sec. 1, masses
of neutrinos are very small. Hence, in order to explain the smallness of
neutrino masses along with $h_{\alpha k}\sim 1$, one can make $\Lambda_{k}$
very small. The above statement is possible if one takes $m_{R}^{2}$ and
$m_{I}^{2}$ to be nearly degenerate, which is described below. In this work,
we take the masses of the components of $\eta$ and $M_{k}$ to be around few
hundred GeV. Now, after using $\lambda_{5}\ll 1$ in the expressions for
$m_{R}^{2}$ and $m_{I}^{2}$, up to first order in $\lambda_{5}$, we get

$\Lambda_{k}=\frac{M_{k}}{8\pi^{2}}\frac{\lambda_{5}v^{2}}{m_{0}^{2}-M_{k}^{2}}\left[1-\frac{M_{k}^{2}}{m_{0}^{2}-M_{k}^{2}}\ln\frac{m_{0}^{2}}{M_{k}^{2}}\right]$
(12)

Using the above equation, one can notice that the smallness of neutrino masses
in the scotogenic model can be explained by suppressing the $\lambda_{5}$
coupling. For this choice of $\lambda_{5}$, the Yukawa couplings $h_{\alpha
k}$ are ${\cal O}(1)$, which can lead to unsuppressed decay rates for LFV
processes in the scotogenic model.

##  3 Analytic expressions for the branching ratios of
$Z\to\ell_{\alpha}\ell_{\beta}$ and $H\to\ell_{\alpha}\ell_{\beta}$

In the scotogenic model, the LFV decays of $Z$ and $H$ are dominantly driven
by $\eta^{\pm}$ and $N_{k}$, which are shown in Fig. 1.

Figure 1: Feynman diagrams representing the decays
$Z,H\to\ell_{\alpha}\ell_{\beta}$. In these diagrams, wavy line corresponds to
either $Z$ gauge boson or Higgs boson.

The amplitudes from the individual diagrams of Fig. 1 can have divergences.
But the sum of the amplitudes from the diagrams of Fig. 1 is finite. For
computing the amplitudes from the diagrams of Fig. 1, we have followed the
work of Ref. [25]. In the individual diagrams of Fig. 1, we assign the
momentum $p$ to the incoming $Z$ or $H$. We assign momentum $p_{1}$ and
$p_{2}$ to the outgoing charged leptons $\ell_{\alpha}$ and $\ell_{\beta}$,
respectively. In the next two subsections, we present analytic results for the
branching ratios of $Z,H\to\ell_{\alpha}\ell_{\beta}$.

###  3.1 Branching ratios of $Z\to\ell_{\alpha}\ell_{\beta}$

In all the diagrams of Fig. 1, we can see propagators due to $\eta^{\pm}$ and
$N_{k}$. Hence, it is convenient to define the following quantities

$D_{k}=q^{2}-M_{k}^{2},\quad D_{1\eta}=(q+p_{1})^{2}-m_{\eta^{\pm}}^{2},\quad
D_{2\eta}=(q-p_{2})^{2}-m_{\eta^{\pm}}^{2}$ (13)

Here, $q$ is a 4-momentum. While computing the amplitudes from the diagrams of
Fig. 1, one come across the following integrals [26], through which we define
the quantities $b_{1,2}^{k}$, $c_{1,2}^{k}$, $d_{1,2}^{k}$, $f^{k}$ and
$u^{k}$.

$\displaystyle\int\frac{d^{d}q}{(2\pi)^{d}}\frac{q^{\mu}}{D_{k}D_{1\eta}}=-b_{1}^{k}p_{1}^{\mu},\quad\int\frac{d^{d}q}{(2\pi)^{d}}\frac{q^{\mu}}{D_{k}D_{2\eta}}=b_{2}^{k}p_{2}^{\mu},$
$\displaystyle\int\frac{d^{d}q}{(2\pi)^{d}}\frac{q^{\mu}}{D_{k}D_{1\eta}D_{2\eta}}=-c_{1}^{k}p_{1}^{\mu}+c_{2}^{k}p_{2}^{\mu},$
$\displaystyle\int\frac{d^{d}q}{(2\pi)^{d}}\frac{q^{\mu}q^{\nu}}{D_{k}D_{1\eta}D_{2\eta}}=d_{1}^{k}p_{1}^{\mu}p_{1}^{\nu}+d_{2}^{k}p_{2}^{\mu}p_{2}^{\nu}-f^{k}(p_{1}^{\mu}p_{2}^{\nu}+p_{2}^{\mu}p_{1}^{\nu})+u^{k}g^{\mu\nu}$
(14)

The above integrals are expressed in $d$-dimensions and at the end of the
calculations we take $d\to 4$. From these integrals, we can notice that
$b_{1,2}^{k}$ and $u^{k}$ are divergent quantities. On the other hand,
$c_{1,2}^{k}$, $d_{1,2}^{k}$ and $f^{k}$ are finite. Using the integrals of
Eq. (14), one can obtain the following relations

$\displaystyle
b_{1}^{k}-b_{2}^{k}=(d_{1}^{k}-d_{2}^{k})m_{Z}^{2}+(\kappa_{1}^{k}+\kappa_{2}^{k})(m_{\ell_{\alpha}}^{2}-m_{\ell_{\beta}}^{2}),$
(15) $\displaystyle
m_{\ell_{\alpha}}^{2}b_{1}^{k}-m_{\ell_{\beta}}^{2}b_{2}^{k}=(m_{\ell_{\alpha}}^{2}d_{1}^{k}-m_{\ell_{\beta}}^{2}d_{2}^{k})m_{Z}^{2}+(m_{\ell_{\alpha}}^{2}-m_{\ell_{\beta}}^{2})[2u^{k}-f^{k}m_{Z}^{2}+\kappa_{1}^{k}m_{\ell_{\alpha}}^{2}+\kappa_{2}^{k}m_{\ell_{\beta}}^{2}],$
(16)
$\displaystyle\kappa_{1}^{k}=d_{1}^{k}+f^{k}-c_{1}^{k},\quad\kappa_{2}^{k}=d_{2}^{k}+f^{k}-c_{2}^{k}$
(17)

Here, $m_{Z}$ is the mass of $Z$ gauge boson.

All the diagrams in Fig. 1 give divergent amplitudes for the case of
$Z\to\ell_{\alpha}\ell_{\beta}$. However, one can notice that the sum of the
amplitudes from these diagrams is finite, after using Eqs. (15) and (16). For
the decay $Z\to\ell^{+}_{\alpha}\ell^{-}_{\beta}$, we have found the total
amplitude from the diagrams of Fig. 1 as

$\displaystyle-i{\cal M}_{Z}$ $\displaystyle=$
$\displaystyle\bar{u}(p_{2})[A_{1}^{L}\gamma^{\mu}P_{L}+A_{1}^{R}\gamma^{\mu}P_{R}+A_{2}^{L}i\sigma^{\mu\nu}p_{\nu}P_{L}+A_{2}^{R}i\sigma^{\mu\nu}p_{\nu}P_{R}]v(p_{1})\epsilon_{\mu}(p),$
$\displaystyle
P_{L(R)}=\frac{1\mp\gamma_{5}}{2},\quad\sigma^{\mu\nu}=\frac{i}{2}[\gamma^{\mu},\gamma^{\nu}],$
$\displaystyle A_{1}^{L}$ $\displaystyle=$
$\displaystyle\sum_{k=1}^{3}\frac{g}{c_{W}}(s_{W}^{2}-\frac{1}{2})h_{\alpha
k}^{*}h_{\beta k}(d_{Z}^{k}-f_{Z}^{k})m_{Z}^{2},\quad
A_{1}^{R}=\sum_{k=1}^{3}\frac{g}{c_{W}}h_{\alpha k}^{*}h_{\beta
k}\kappa_{Z}^{k}m_{\ell_{\alpha}}m_{\ell_{\beta}},$ $\displaystyle A_{2}^{L}$
$\displaystyle=$
$\displaystyle\sum_{k=1}^{3}\frac{g}{c_{W}}(s_{W}^{2}-\frac{1}{2})h_{\alpha
k}^{*}h_{\beta k}\kappa_{Z}^{k}m_{\ell_{\beta}},\quad
A_{2}^{R}=\sum_{k=1}^{3}\frac{g}{c_{W}}(s_{W}^{2}-\frac{1}{2})h_{\alpha
k}^{*}h_{\beta k}\kappa_{Z}^{k}m_{\ell_{\alpha}},$ $\displaystyle d_{Z}^{k}$
$\displaystyle=$ $\displaystyle
d_{1}^{k}=d_{2}^{k}=\frac{-i}{16\pi^{2}}\int_{0}^{1}dx\int_{0}^{1-x}dy\frac{y^{2}}{-y(1-x-y)m_{Z}^{2}+xM_{k}^{2}+(1-x)m_{\eta^{\pm}}^{2}},$
$\displaystyle f_{Z}^{k}$ $\displaystyle=$
$\displaystyle\frac{-i}{16\pi^{2}}\int_{0}^{1}dx\int_{0}^{1-x}dy\frac{y(1-x-y)}{-y(1-x-y)m_{Z}^{2}+xM_{k}^{2}+(1-x)m_{\eta^{\pm}}^{2}},$
$\displaystyle c_{Z}^{k}$ $\displaystyle=$ $\displaystyle
c_{1}^{k}=c_{2}^{k}=\frac{-i}{16\pi^{2}}\int_{0}^{1}dx\int_{0}^{1-x}dy\frac{y}{-y(1-x-y)m_{Z}^{2}+xM_{k}^{2}+(1-x)m_{\eta^{\pm}}^{2}},$
$\displaystyle\kappa_{Z}^{k}$ $\displaystyle=$
$\displaystyle\kappa_{1}^{k}=\kappa_{2}^{k}=d_{Z}^{k}+f_{Z}^{k}-c_{Z}^{k}.$
(18)

Here, $s_{W}(c_{W})=\sin\theta_{W}(\cos\theta_{W})$, where $\theta_{W}$ is the
weak-mixing angle. $g$ is the coupling strength of $SU(2)$ gauge group of the
standard model. From the above amplitude, notice that, except $A_{1}^{L}$,
rest of the form factors of it are proportional to charged lepton masses.
Since $\frac{m_{\ell_{\alpha}}^{2}}{m_{Z}^{2}}\ll 1$, the form factors
$A_{1}^{R}$ and $A_{2}^{L,R}$ give subleading contributions to the branching
ratio of $Z\to\ell^{+}_{\alpha}\ell^{-}_{\beta}$. As a result of this, the
leading contribution to the branching ratio of $Z\to\ell_{\alpha}\ell_{\beta}$
is found to be

$\displaystyle{\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ $\displaystyle=$
$\displaystyle\frac{\Gamma(Z\to\ell^{+}_{\alpha}\ell^{-}_{\beta})+\Gamma(Z\to\ell^{-}_{\alpha}\ell^{+}_{\beta})}{\Gamma_{Z}}$
(19) $\displaystyle=$
$\displaystyle\left(\frac{g}{c_{W}}\right)^{2}\left(s_{W}^{2}-\frac{1}{2}\right)^{2}\frac{m_{Z}^{5}}{12\Gamma_{Z}}\left|\sum_{k=1}^{3}h_{\alpha
k}^{*}h_{\beta k}(d_{Z}^{k}-f_{Z}^{k})\right|^{2}$

Here, $\Gamma_{Z}$ is the total decay width of $Z$ gauge boson. In our
numerical analysis, which is presented in the next section, we have taken
$\Gamma_{Z}=$ 2.4952 GeV [3].

###  3.2 Branching ratios of $H\to\ell_{\alpha}\ell_{\beta}$

While computing the amplitude for $H\to\ell_{\alpha}\ell_{\beta}$, we can
define the integrals of Eq. (14). Moreover, the relations in Eqs. (15) and
(16) are also valid in this case after replacing $m_{Z}^{2}$ with $m_{H}^{2}$
in these equations. Now, for the case of $H\to\ell_{\alpha}\ell_{\beta}$, the
top two diagrams of Fig. 1 give divergent amplitudes, whereas, the bottom
diagram of this figure give finite contribution. Hence, only the analog of Eq.
(15) is sufficient to see the cancellation of divergences between the top two
diagrams of Fig. 1. Now, after summing the amplitudes from the diagrams of
Fig. 1, for the decay $H\to\ell^{+}_{\alpha}\ell^{-}_{\beta}$, we have found
the total amplitude as

$\displaystyle i{\cal M}_{H}$ $\displaystyle=$
$\displaystyle\bar{u}(p_{2})[A_{H}^{L}P_{L}+A_{H}^{R}P_{R}]v(p_{1})$
$\displaystyle A_{H}^{L}$ $\displaystyle=$
$\displaystyle\sqrt{2}\sum_{k=1}^{3}h_{\alpha k}^{*}h_{\beta
k}(\lambda_{3}c_{H}^{k}+\frac{m_{\ell_{\alpha}}^{2}}{v^{2}}\kappa_{H}^{k})vm_{\ell_{\beta}},$
$\displaystyle A_{H}^{R}$ $\displaystyle=$
$\displaystyle\sqrt{2}\sum_{k=1}^{3}h_{\alpha k}^{*}h_{\beta
k}(\lambda_{3}c_{H}^{k}+\frac{m_{\ell_{\beta}}^{2}}{v^{2}}\kappa_{H}^{k})vm_{\ell_{\alpha}}$
(20)

The expressions for $c_{H}^{k}$ and $\kappa_{H}^{k}$ are respectively same as
that for $c_{Z}^{k}$ and $\kappa_{Z}^{k}$, after replacing $m_{Z}^{2}$ with
$m_{H}^{2}$ in these expressions. The first term in $A_{H}^{L,R}$ is arising
due to the bottom diagram of Fig. 1. On the other hand, the top two diagrams
of Fig. 1 contribute to the second term in $A_{H}^{L,R}$. One can see that for
$\lambda_{3}\sim 1$, the second term in $A_{H}^{L,R}$ gives negligibly small
contribution. In our numerical analysis, we consider $\lambda_{3}\sim 1$.
Hence, for a case like this, the branching ratio for
$H\to\ell_{\alpha}\ell_{\beta}$ is found to be

$\displaystyle{\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ $\displaystyle=$
$\displaystyle\frac{\Gamma(H\to\ell^{+}_{\alpha}\ell^{-}_{\beta})+\Gamma(H\to\ell^{-}_{\alpha}\ell^{+}_{\beta})}{\Gamma_{H}}$
(21) $\displaystyle=$
$\displaystyle\frac{m_{H}}{4\pi\Gamma_{H}}(\lambda_{3}v)^{2}(m_{\ell_{\alpha}}^{2}+m_{\ell_{\beta}}^{2})\left|\sum_{k=1}^{3}h_{\alpha
k}^{*}h_{\beta k}c_{H}^{k}\right|^{2}$

Here, $\Gamma_{H}$ is the total Higgs decay width.

In our numerical analysis, which is presented in the next section, we have
taken $m_{H}=$ 125.1 GeV [3] and $\Gamma_{H}=4.08\times 10^{-3}$ GeV [27].
This value of $\Gamma_{H}$ is same as that for the Higgs boson of standard
model. We have taken this value for $\Gamma_{H}$ in order to simplify our
numerical analysis. The above mentioned value of $\Gamma_{H}$ has an
implication that the Higgs boson should not decay into $Z_{2}$-odd particles
of the scotogenic model. We comment further about this later.

##  4 Numerical analysis

From the analytic expressions given in the previous section, we can see that
the branching ratios of $Z,H\to\ell_{\alpha}\ell_{\beta}$ are proportional to
the Yukawa couplings $h_{\alpha k}$. The same Yukawa couplings also drive
neutrino masses which are described in Sec. 2. It is worth to explore the
correlation between neutrino oscillation observables and the branching ratios
of $Z,H\to\ell_{\alpha}\ell_{\beta}$. Here, our objective is to fit the
neutrino oscillation observables in the scotogenic model in such a way that
the branching ratios for $Z,H\to\ell_{\alpha}\ell_{\beta}$ can become maximum
in this model. It is explained in Sec. 2 that the above objective can be
achieved by taking $h_{\alpha k}\sim 1$ and $\Lambda_{k}$ very small. Below we
describe the procedure in order to achieve this objective.

The neutrino oscillation observables can be explained in the scotogenic model
by parametrizing the Yukawa couplings as given in Eq. (10). In this equation,
$R$ is an orthogonal matrix, whose elements can have a magnitude of ${\cal
O}(1)$. To simplify our numerical analysis we take $R$ to be a unit matrix. In
such a case we get

$h=U^{*}_{PMNS}\cdot{\rm
diag}\left(\sqrt{\frac{m_{1}}{\Lambda_{1}}},\sqrt{\frac{m_{2}}{\Lambda_{2}}},\sqrt{\frac{m_{3}}{\Lambda_{3}}}\right)$
(22)

In our analysis we have parametrized $U_{PMNS}$ as [3]

$U_{PMNS}=\left(\begin{array}[]{ccc}c_{12}c_{13}&s_{12}c_{13}&s_{13}e^{-i\delta_{CP}}\\\
-s_{12}c_{23}-c_{12}s_{23}s_{13}e^{i\delta_{CP}}&c_{12}c_{23}-s_{12}s_{23}s_{13}e^{i\delta_{CP}}&s_{23}c_{13}\\\
s_{12}s_{23}-c_{12}c_{23}s_{13}e^{i\delta_{CP}}&-c_{12}s_{23}-s_{12}c_{23}s_{13}e^{i\delta_{CP}}&c_{23}c_{13}\end{array}\right)$
(23)

Here, $c_{ij}=\cos\theta_{ij}$, $s_{ij}=\sin\theta_{ij}$ and $\delta_{CP}$ is
the $CP$ violating Dirac phase. We have taken Majorana phases to be zero in
$U_{PMNS}$. Shortly below we describe the numerical values for neutrino masses
and mixing angles. Using these values, we can see that the elements of
$U_{PMNS}$ can have a magnitude of ${\cal O}(1)$. Hence, we need to make
$\frac{m_{k}}{\Lambda_{k}}\sim 1$ for $k=1,2,3$ in order to get $h_{\alpha
k}\sim 1$. Since neutrino mass eigenvalues $m_{k}$ are very small,
$\Lambda_{k}$ should be proportionately small in order to achieve $h_{\alpha
k}\sim 1$. It is described in Sec. 2 that $\Lambda_{k}$ can be made very small
by suppressing the $\lambda_{5}$ parameter.

From the global fits to neutrino oscillation data the following mass-square
differences among the neutrino fields are found [5].

$m_{s}^{2}=m_{2}^{2}-m_{1}^{2}=7.5\times 10^{-5}~{}{\rm eV}^{2},\quad
m_{a}^{2}=\left\\{\begin{array}[]{c}m_{3}^{2}-m_{1}^{2}=2.55\times
10^{-3}~{}{\rm eV}^{2}~{}~{}{\rm(NO)}\\\ m_{1}^{2}-m_{3}^{2}=2.45\times
10^{-3}~{}{\rm eV}^{2}~{}~{}{\rm(IO)}\end{array}\right..$ (24)

Here, NO(IO) represents normal(inverted) ordering. In the above equation we
have given the best fit values. In order to fit the above mass-square
differences, we take the neutrino mass eigenvalues as

$\displaystyle{\rm NO}:\quad m_{1}=0.1m_{s},\quad
m_{2}=\sqrt{m_{s}^{2}+m_{1}^{2}},\quad m_{3}=\sqrt{m_{a}^{2}+m_{1}^{2}}.$
$\displaystyle{\rm IO}:\quad m_{3}=0.1m_{s},\quad
m_{1}=\sqrt{m_{a}^{2}+m_{3}^{2}},\quad m_{2}=\sqrt{m_{s}^{2}+m_{1}^{2}}.$ (25)

The above neutrino mass eigenvalues satisfy the cosmological upper bound on
the sum of neutrino masses, which is 0.12 eV [28]. Apart from neutrino masses,
neutrino mixing angles are also found from the global fits to neutrino
oscillation data [5]. The best fit and 3$\sigma$ ranges for these variables
are given in Table 1.

parameter | best fit | 3$\sigma$ range
---|---|---
$\sin^{2}\theta_{12}/10^{-1}$ | 3.18 | 2.71 - 3.69
$\sin^{2}\theta_{13}/10^{-2}$ (NO) | 2.200 | 2.000 - 2.405
$\sin^{2}\theta_{13}/10^{-2}$ (IO) | 2.225 | 2.018 - 2.424
$\sin^{2}\theta_{23}/10^{-1}$ (NO) | 5.74 | 4.34 - 6.10
$\sin^{2}\theta_{23}/10^{-1}$ (IO) | 5.78 | 4.33 - 6.08
$\delta_{CP}/{\rm o}$ (NO) | 194 | 128 - 359
$\delta_{CP}/{\rm o}$ (IO) | 284 | 200 - 353
Table 1: Best fit and 3$\sigma$ ranges for the neutrino mixing angles and $CP$
violating Dirac phase, which are obtained from the global fits to neutrino
oscillation data [5].

In the next two subsections, we present numerical results on the branching
ratios of $Z,H\to\ell_{\alpha}\ell_{\beta}$. From the analytic expressions
given in the previous section, we can see that the above mentioned branching
ratios can become maximum for large values of Yukawa couplings and
$\lambda_{3}$ parameter. In order to satisfy the perturbativity limits on
these variables, we apply the following constraints on the Yukawa couplings
and the $\lambda$ parameters of the scotogenic model.

$|h_{\alpha k}|\leq\sqrt{4\pi},\quad|\lambda_{i}|\leq 4\pi$ (26)

###  4.1 $Z\to\ell_{\alpha}\ell_{\beta}$

As explained previously, to satisfy perturbativity limit, $|h_{\alpha k}|$ can
be as large as $\sqrt{4\pi}$. Since the magnitude of the elements of
$U_{PMNS}$ are less than about one, from Eq. (22) we can see that
$\frac{m_{k}}{\Lambda_{k}}$ can be as large as 4$\pi$ in order to satisfy the
above mentioned perturbativity limit. $\Lambda_{k}$ depends on $M_{k}$,
$m_{0}$ and $\lambda_{5}$. We have plotted $\frac{m_{k}}{\Lambda_{k}}$ versus
$\lambda_{5}$ in Fig. 2.

Figure 2: Plots between $\frac{m_{k}}{\Lambda_{k}}$ and $\lambda_{5}$. Red,
blue and green lines are for $\frac{m_{1}}{\Lambda_{1}}$,
$\frac{m_{2}}{\Lambda_{2}}$ and $\frac{m_{3}}{\Lambda_{3}}$ respectively.
Horizontal line indicates the value 4$\pi$. Left- and right-hand side plots
are for NO and IO respectively. In both the plots, we have taken $m_{0}$ = 150
GeV, $M_{1}$ = 100 GeV, $M_{2}=M_{1}+50$ GeV and $M_{3}=M_{2}+50$ GeV.

In these plots, we have chosen masses for right-handed neutrinos to be between
100 to 200 GeV. The reason for such a choice is that, for these low masses of
right-handed neutrinos ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ can become
maximum. Results related to ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ will be
presented shortly later. In the plots of Fig. 2, for the case of NO, all the
lines are distinctly spaced because of the fact that the neutrino masses are
hierarchical in this case. On the other hand, the neutrino mass eigenvalues
$m_{1}$ and $m_{2}$ are nearly degenerate for the case of IO. As a result of
this, red and blue lines in the right-hand side plot of Fig. 2 are close to
each other. From this figure, we can see that $\frac{m_{k}}{\Lambda_{k}}$
increases when $\lambda_{5}$ is decreasing. This follows from the fact that in
the limit $\lambda_{5}\to 0$, $m_{R}^{2}$ and $m_{I}^{2}$ are degenerate, and
hence, $\Lambda_{k}$ becomes vanishingly small. From Fig. 2, in the case of
NO, for $\lambda_{5}=3\times 10^{-3}$ we get $\frac{m_{3}}{\Lambda_{3}}\approx
4\pi$. Hence, for $\lambda_{5}<3\times 10^{-3}$ and for the values of
$m_{0},M_{k}$ taken in Fig. 2, the perturbativity limit for Yukawa couplings,
which is given in Eq. (26), can be violated. Similarly, from the right-hand
side plot of Fig. 2, we can see that the above mentioned perturbativity limit
can be violated for $\lambda_{5}<3.7\times 10^{-3}$, in the case of IO.

From Fig. 2, we have obtained the minimum value of $\lambda_{5}$ through which
the perturbativity limit on the Yukawa couplings can be satisfied. Using this
minimum value of $\lambda_{5}$ we have plotted branching ratios for
$Z\to\ell_{\alpha}\ell_{\beta}$ in Fig. 3 for the case of NO.

Figure 3: Plots between ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ and
$\delta_{CP}$ for the case of NO, without applying the constraints due to non-
observation of $\ell_{\alpha}\to\ell_{\beta}\gamma$. Numerical values for
neutrino masses are taken from Eq. (25). Neutrino mixing angles are taken to
be the best fit values, which are given in Table 1. In both the plots, we have
taken $\lambda_{5}=3\times 10^{-3}$, $m_{0}$ = 150 GeV, $m_{\eta^{\pm}}$ = 170
GeV, $M_{1}$ = 100 GeV, $M_{2}=M_{1}+50$ GeV and $M_{3}=M_{2}+50$ GeV. In the
left-hand side plot, red and blue lines are for $e\mu$ and $e\tau$ modes
respectively.

In the plots of this figure, we have taken $m_{\eta^{\pm}}$ to be as low as
170 GeV. One can understand that by increasing this value, branching ratios
for $Z\to\ell_{\alpha}\ell_{\beta}$ decreases. The plots in Fig. 3 are made
after fitting to the neutrino oscillation observables in the scotogenic model.
We can see that the branching ratios for $Z\to\ell_{\alpha}\ell_{\beta}$, in
this model, can be as large as $10^{-8}-10^{-9}$. These branching ratio values
are lower than the current experimental limits on them, which are given in Eq.
(1). On the other hand, these values can be probed in the future FCC-ee
collider, which can be seen in Eq. (3). However, as will be described below,
the above mentioned branching ratio values will be suppressed, if constraints
due to non-observation of $\ell_{\alpha}\to\ell_{\beta}\gamma$ are applied. We
have also made the analog plots of Fig. 3, for the case of IO, by taking
$\lambda_{5}=3.7\times 10^{-3}$. We have found that, in the case of IO, the
branching ratios for $Z\to\ell_{\alpha}\ell_{\beta}$ are slightly higher than
that of plots in Fig. 3. But, otherwise, the shape of the curves for ${\rm
Br}(Z\to\ell_{\alpha}\ell_{\beta})$, in the case of IO, are same as that of
Fig. 3.

Regarding the shape of the curves in Fig. 3, we can notice that the shapes of
${\rm Br}(Z\to e\mu)$ and ${\rm Br}(Z\to\mu\tau)$, with respect to
$\delta_{CP}$, are similar. On the other hand, the shapes of ${\rm Br}(Z\to
e\mu)$ and ${\rm Br}(Z\to e\tau)$, with respect to $\delta_{CP}$, are opposite
to each other. We have found that the shapes of the curves for ${\rm Br}(Z\to
e\mu)$ and ${\rm Br}(Z\to e\tau)$, with respect to $\delta_{CP}$, do not
change by changing the values for neutrino mixing angles. On the other hand,
the shape of the curve for ${\rm Br}(Z\to\mu\tau)$, with respect to
$\delta_{CP}$, changes with $s_{23}^{2}$. For $s_{23}^{2}>0.5$, which is the
case considered in Fig. 3, the shape of the curves for ${\rm Br}(Z\to e\mu)$
and ${\rm Br}(Z\to\mu\tau)$ are found to be similar. In contrast to this, for
$s_{23}^{2}<0.5$, the shape of the curve for ${\rm Br}(Z\to\mu\tau)$ is found
to be similar to that of ${\rm Br}(Z\to e\tau)$. Whereas, for
$s_{23}^{2}=0.5$, the shape of the curve for ${\rm Br}(Z\to\mu\tau)$ has no
resemblance with either to that of ${\rm Br}(Z\to e\mu)$ and ${\rm Br}(Z\to
e\tau)$. The shapes of the above mentioned branching ratios with respect to
$\delta_{CP}$ depend on the Yukawa couplings, which in our case is given in
Eq. (22). After using the form of these Yukawa couplings in the branching
ratio expressions of Eq. (19), one can understand the above described shapes
with respect to $\delta_{CP}$.

Plots in Fig. 3 are made for a minimum value of $\lambda_{5}$ for which the
Yukawa couplings can be close to a value of $\sqrt{4\pi}$. However, the Yukawa
couplings $h_{\alpha k}$ can also drive the decays
$\ell_{\alpha}\to\ell_{\beta}\gamma$, whose branching ratios in the scotogenic
model are as follows [8].

$\displaystyle{\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ $\displaystyle=$
$\displaystyle\frac{3\alpha_{EW}}{64\pi
G_{F}^{2}m_{\eta^{\pm}}^{4}}\left|\sum_{k=1}^{3}h^{*}_{\alpha k}h_{\beta
k}F_{2}\left(\frac{M_{k}^{2}}{m_{\eta^{\pm}}^{2}}\right)\right|^{2},$
$\displaystyle F_{2}(x)$ $\displaystyle=$
$\displaystyle\frac{1-6x+3x^{2}+2x^{3}-6x^{2}\ln x}{6(1-x)^{4}}$ (27)

Here, $\alpha_{EW}$ and $G_{F}$ are fine-structure and Fermi constants,
respectively. The decays $\ell_{\alpha}\to\ell_{\beta}\gamma$ are not observed
in experiments. Hence, the branching ratios for these decays are constrained
as follows.

$\displaystyle{\rm Br}(\mu\to e\gamma)$ $\displaystyle<$ $\displaystyle
4.2\times 10^{-13}~{}\cite[cite]{[\@@bibref{}{meg}{}{}]},$ $\displaystyle{\rm
Br}(\tau\to e\gamma)$ $\displaystyle<$ $\displaystyle 3.3\times
10^{-8}~{}\cite[cite]{[\@@bibref{}{babar}{}{}]},$ $\displaystyle{\rm
Br}(\tau\to\mu\gamma)$ $\displaystyle<$ $\displaystyle 4.4\times
10^{-8}~{}\cite[cite]{[\@@bibref{}{babar}{}{}]}$ (28)

After comparing Eqs. (19) and (27), we can see that the same set of model
parameters which determine ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ also
determine ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. For the set of model
parameters taken in Fig. 3, we have found that the branching ratios for
$\ell_{\alpha}\to\ell_{\beta}\gamma$ exceed the experimental bounds of Eq.
(28). The reason for this is as follows. In the plots of Fig. 3, the Yukawa
couplings are close to $\sqrt{4\pi}$ and the masses of mediating particles are
between 100 to 200 GeV. For such large Yukawa couplings and low masses, the
branching ratios for $\ell_{\alpha}\to\ell_{\beta}\gamma$ are quite large that
they do not respect the bounds of Eq. (28). Hence, the plots in Fig. 3 give us
the maximum values that the branching ratios of
$Z\to\ell_{\alpha}\ell_{\beta}$ can reach in the scotogenic model, without
applying constraints due to non-observation of
$\ell_{\alpha}\to\ell_{\beta}\gamma$.

Now, it is our interest to know the branching ratios of
$Z\to\ell_{\alpha}\ell_{\beta}$ after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. One can notice that ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ depends on Yukawa couplings, masses
of right-handed neutrinos and $\eta^{\pm}$. Hence, to satisfy the bounds on
${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$, one has to suppress Yukawa
couplings and increase the masses for right-handed neutrinos and $\eta^{\pm}$.
The mass of $\eta^{\pm}$ can be written as
$m_{\eta^{\pm}}=\sqrt{m_{0}^{2}-\lambda_{4}v^{2}}$. To satisfy the
perturbativity limit on $\lambda_{4}$, we choose $\lambda_{4}=-4\pi$. With
this choice, the mass of $\eta^{\pm}$ can take maximum value, for a fixed
value of $m_{0}$. Now, the Yukawa couplings depend on $m_{0}$, $\lambda_{5}$
and masses of right-handed neutrinos, apart from neutrino oscillation
observables. Hence, for the above mentioned choice, ${\rm
Br}(Z\to\ell_{\alpha}\ell_{\beta})$ and ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ depend on $m_{0}$, $\lambda_{5}$ and
masses of right-handed neutrinos, apart from neutrino oscillation observables.
In Fig. 4, we have plotted branching ratios of $Z\to\ell_{\alpha}\ell_{\beta}$
after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$.

Figure 4: Plots between ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ and
$\lambda_{5}$ for the case of NO, after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. In these plots, solid lines are
allowed and dotted lines are excluded by the constraints due to ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. Numerical values for neutrino masses
are taken from Eq. (25). Neutrino mixing angles and $\delta_{CP}$ are taken to
be the best fit values, which are given in Table 1. We have taken $m_{0}$ =
150 GeV, $m_{\eta^{\pm}}=\sqrt{m_{0}^{2}+4\pi v^{2}}$, $M_{1}$ = 1000 GeV,
$M_{2}=M_{1}+100$ GeV and $M_{3}=M_{2}+100$ GeV.

In Fig. 4, we have varied $\lambda_{5}$ up to 0.61. The reason for this is
explained below. For the parametric values of Fig. 4, we can see that the
lightest $Z_{2}$-odd particle in the scotogenic model is $\eta_{I}$. The mass
of $\eta_{I}$ decreases with $\lambda_{5}$. At $\lambda_{5}=0.61$, we get
$m_{I}\approx$ 63.5 GeV. Since the Higgs boson mass is 125.1 GeV, for
$\lambda_{5}>0.61$ there is a possibility that the Higgs can decay into a pair
of $\eta_{I}$. It is described in the previous section that the total decay
width for the Higgs boson in our analysis is taken to be the same as that in
the standard model. Hence, to avoid the above mentioned decay, we have varied
$\lambda_{5}$ up to 0.61 in Fig. 4.

From Fig. 4, we can see that the branching ratios for
$Z\to\ell_{\alpha}\ell_{\beta}$ vary in the range of $10^{-17}-10^{-15}$.
These values are suppressed by about $10^{-7}$ as compared that in Fig. 3. The
reason for this suppression is due to the fact that the $\lambda_{5}$ and the
masses of right-handed neutrinos and $\eta^{\pm}$ are large as compared to
those in Fig. 3. As already stated before, the masses of right-handed
neutrinos and $\eta^{\pm}$ should be taken large, otherwise, the constraints
on ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ cannot be satisfied. The
mass of lightest right-handed neutrino in Fig. 4 is taken to be 1 TeV. We have
found that, for the case of $M_{2}=M_{1}+100$ GeV and $M_{3}=M_{2}+100$ GeV,
$M_{1}$ should be at least around 500 GeV in order to satisfy the constraints
from ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. However, in such a case,
the allowed range for $\lambda_{5}$ becomes narrower than that in Fig. 4 and
the allowed ranges for ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ are found to
be nearly same as that in Fig. 4. Although the right-handed neutrino masses
are taken to be non-degenerate in Fig. 4, the plots in this figure do not vary
much with degenerate right-handed neutrinos of 1 TeV masses. It is stated
above that another reason for the suppression of ${\rm
Br}(Z\to\ell_{\alpha}\ell_{\beta})$ in Fig. 4 is due to the fact that
$\lambda_{5}$ is large. This suppression is happening because Yukawa couplings
reduce with increasing $\lambda_{5}$. This fact can be understood with the
plots of Fig. 2 and also with Eq. (22).

In the plots of Fig. 4, we have fixed $m_{0}$ to 150 GeV. By increasing this
value to 500 GeV, we have found that ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$
reduces as compared to that in Fig. 4. This is happening because
$m_{\eta^{\pm}}$ increases. Another difference we have noticed is that, for
$m_{0}=$ 500 GeV and right-handed neutrino masses to be same as in Fig. 4, the
allowed range for $\lambda_{5}$ is found to be $\sim 1.5-8.0$. This is
happening because, by increasing $m_{0}$, one has to increase $\lambda_{5}$ in
order to suppress the Yukawa couplings and thereby satisfy the constraints on
${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$.

We have plotted ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ for the case of IO,
which are presented in Fig. 5.

Figure 5: Plots between ${\rm Br}(Z\to\ell_{\alpha}\ell_{\beta})$ and
$\lambda_{5}$ for the case of IO, after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. In these plots, solid lines are
allowed and dotted lines are excluded by the constraints due to ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. Numerical values for neutrino masses
are taken from Eq. (25). Neutrino mixing angles and $\delta_{CP}$ are taken to
be the best fit values, which are given in Table 1. We have taken $m_{0}$ =
150 GeV, $m_{\eta^{\pm}}=\sqrt{m_{0}^{2}+4\pi v^{2}}$, $M_{1}$ = 2000 GeV,
$M_{2}=M_{1}+100$ GeV and $M_{3}=M_{2}+100$ GeV.

In this case, for $m_{0}=$ 150 GeV we have found that $M_{1}$ should be at
least around 1.7 TeV in order to satisfy the constraints on ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. As a result of this, in the plots of
Fig. 5, we have taken $M_{1}=$ 2 TeV. Comparing the plots of Figs. 4 and 5, we
can conclude the following points. In both the cases of NO and IO, ${\rm
Br}(Z\to\mu\tau)$ is larger than that for the other LFV decays of $Z$ gauge
boson. In the case of NO, ${\rm Br}(Z\to e\mu)$ is one order less than ${\rm
Br}(Z\to e\tau)$. On the other hand, in the case of IO, ${\rm Br}(Z\to e\mu)$
is slightly larger than ${\rm Br}(Z\to e\tau)$.

###  4.2 $H\to\ell_{\alpha}\ell_{\beta}$

In this subsection, we present numerical results on the branching ratios of
$H\to\ell_{\alpha}\ell_{\beta}$. After comparing Eqs. (19) and (21), we can
see that a common set of parameters determine both ${\rm
Br}(Z\to\ell_{\alpha}\ell_{\beta})$ and ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$. Apart from this common set of parameters,
$\lambda_{3}$ is an additional parameter which determine ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$. In our analysis, we have taken
$\lambda_{3}=4\pi$ in order to satisfy the perturbativity limit and also to
maximize ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$. Apart from the above
mentioned parameters, ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ also depends
on the charged lepton masses. We have taken these masses to be the best fit
values, which are given in Ref. [3].

First we present the results on ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$
after fitting to the neutrino oscillation observables, but without satisfying
the constraints from ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. These
results are given in Fig. 6 for the case of NO.

Figure 6: Plots between ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ and
$\delta_{CP}$ for the case of NO, without applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. See the caption of Fig. 3, for
parametric values and neutrino oscillation observables, which are used in
these plots.

One can compare the branching ratios in this figure with the current limits on
them, which are given in Eq. (2). We can see that the values for ${\rm
Br}(H\to e\mu)$ and ${\rm Br}(H\to e\tau)$ from this figure are marginally
lower than the current experimental limits on them. Whereas, the values for
${\rm Br}(H\to\mu\tau)$ are just below the current experimental limit on this.
However, in the plots of Fig. 6, we have taken $\lambda_{5}=3\times 10^{-3}$
and the masses of right-handed neutrinos and $\eta^{\pm}$ are chosen to be
between 100 to 200 GeV. For this choice of parameters, as already explained in
the previous subsection, the Yukawa couplings can be large, and hence, ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$ can become maximum. Plots in Fig. 6 are
made for the case of NO. We have plotted ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$ for the case of IO by taking
$\lambda_{5}=3.7\times 10^{-3}$ and for the mass parameters which are
described above. In this case, we have found a slight enhancement in the
values of ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ as compared to that of
Fig. 6. But otherwise, in the case of IO, the shape of the curves for ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$ are found to be the same as that in Fig.
6.

In the plots of Fig. 6, constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ are not applied. After applying the
constraints from ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$, branching
ratios for $H\to\ell_{\alpha}\ell_{\beta}$ are given in Fig. 7 for the case of
NO.

Figure 7: Plots between ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ and
$\lambda_{5}$ for the case of NO, after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. In these plots, solid lines are
allowed and dotted lines are excluded by the constraints due to ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. See the caption of Fig. 4, for
parametric values and neutrino oscillation observables, which are used in
these plots.

One can see that the branching ratios in this figure are suppressed by a
factor of about $10^{-7}$ as compared to that in Fig. 6. The reason for this
suppression, which can be understood from the reasoning’s given around Fig. 4,
is due to the fact that $\lambda_{5}$ and masses of right-handed neutrinos and
$\eta^{\pm}$ are large as compared that in Fig. 6. The mass of lightest right-
handed neutrino is 1 TeV in Fig. 7. As already pointed around Fig. 4, the
value of $M_{1}$ should be at least around 500 GeV in order to satisfy the
constraints from ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ for the case
of Fig. 7. Even with $M_{1}=$ 500 GeV, we have found the allowed ranges for
${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ are nearly same as that of Fig. 7.
Although the right-handed neutrino masses are non-degenerate in Fig. 7, with
degenerate right-handed neutrinos with masses of 1 TeV we have found that the
allowed ranges for ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ are similar to
that in Fig. 7. In this figure, among the three LFV decays of $H$, the
branching ratios of $H$ into $\tau$ mode are large, since these branching
ratios are proportional to $m_{\tau}^{2}$.

We have plotted ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$, after applying the
constraints from ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$, for the case
of IO. These plots are given in Fig. 8.

Figure 8: Plots between ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$ and
$\lambda_{5}$ for the case of IO, after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. In these plots, solid lines are
allowed and dotted lines are excluded by the constraints due to ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$. See the caption of Fig. 5, for
parametric values and neutrino oscillation observables, which are used in
these plots.

The masses for right-handed neutrinos are different in this figure as compared
to that in Fig. 7. Nevertheless, the allowed range of values for ${\rm
Br}(H\to\ell_{\alpha}\ell_{\beta})$ are found to be nearly same in Figs. 7 and
8.

Among the LFV decays of $Z$ and $H$, after applying the constraints from ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$, $H\to\mu\tau$ is found to have the
largest branching ratio, which is around $10^{-10}$. This indicates that
probing LFV decays of Higgs boson in experiments is one possible way to test
the scotogenic model. However, in our analysis of LFV decays of $H$, we have
taken $\lambda_{3}=4\pi$, which is the maximum possible value for this
parameter. In this model, the $\lambda_{3}$ coupling can also drive the decay
$H\to\gamma\gamma$. In the LHC experiment, it is found that there is no
enhancement in the signal strength of this decay as compared to the standard
model prediction [3]. As a result of this, one can expect some constraints on
$\lambda_{3}$ parameter. Apart from this, the model parameters of the
scotogenic model can get additional constraints due to precision electroweak
observables and relic abundance of dark matter. One may expect that the above
mentioned constraints can lower the allowed ranges for the branching ratios of
LFV decays of $Z$ and $H$ in this model.

As stated in Sec. 1, in the context of scotogenic model, branchig ratio for
$H\to\mu\tau$ has been estimated as ${\rm Br}(H\to\mu\tau)\mathrel{\raise
2.58334pt\hbox{$<$}\kern-7.7778pt\lower
2.79857pt\hbox{$\sim$}}10^{-7}\lambda_{3}^{2}$ [22], after applying the
constraint from ${\rm Br}(\tau\to\mu\gamma)$. In our analysis, we have applied
constraints due to non-observation of all LFV decays of the form
$\ell_{\alpha}\to\ell_{\beta}\gamma$ and we have found that ${\rm
Br}(H\to\mu\tau)$ can be as large as $\sim 10^{-10}$, even with
$\lambda_{3}=4\pi$. Hence, our result on ${\rm Br}(H\to\mu\tau)$ is more
stringent than the above mentioned estimation of Ref. [22].

##  5 Conclusions

In this work, we have studied LFV decays of $Z$ gauge boson and Higgs boson in
the scotogenic model. After deriving analytic expressions for the branching
ratios of the above mentioned decays, numerically we have studied how large
they can be in this model. The above mentioned numerical study has been done
by satisfying the following quantities: fit to neutrino oscillation
observables, constraints on ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ and
perturbativity limits on the parameters of the model. If we satisfy only the
fit to neutrino oscillation observables and the perturbativity limits on the
model parameters, we have found the following maximum values for the branching
ratios of LFV decays of $Z$ and $H$: ${\rm Br}(Z\to e\mu,e\tau)\sim 10^{-9}$,
${\rm Br}(Z\to\mu\tau)\sim 10^{-8}$, ${\rm Br}(H\to e\mu)\sim 10^{-7}$, ${\rm
Br}(H\to e\tau)\sim 10^{-4}$, ${\rm Br}(H\to\mu\tau)\sim 10^{-3}$. However, in
addition to satisfying the above mentioned quantities, after satisfying
constraints on ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$, the above
mentioned results on the branching ratios get an additional suppression of
about $10^{-7}$. If the scotogenic model is true, results obtained in this
work can give indication about future results on LFV decays of $Z$ and $H$ in
the upcoming experiments.

Note added: While this manuscript was under preparation, Ref. [31] had
appeared where LFV decays of Higgs boson were studied in the scotogenic model.
The method of computing the branching ratios for these decays and numerical
study done on them in Ref. [31] are found to be different from what we have
done in this work. After comparing ${\rm Br}(H\to\ell_{\alpha}\ell_{\beta})$
versus ${\rm Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ in Ref. [31], it is
shown that the allowed values for ${\rm
Br}(\ell_{\alpha}\to\ell_{\beta}\gamma)$ are suppressed to around $10^{-34}$.
Moreover, the branching ratio for $H\to\mu\tau$ is also shown to be suppressed
to around $10^{-37}$. The above mentioned results are different from what we
have presented here.

## References

  * [1] C. Quigg, [arXiv:hep-ph/0404228 [hep-ph]]; J. Ellis, Nucl. Phys. A 827, 187C-198C (2009) doi:10.1016/j.nuclphysa.2009.05.034 [arXiv:0902.0357 [hep-ph]].
  * [2] T. Mori, eConf C060409, 034 (2006) [arXiv:hep-ex/0605116 [hep-ex]];
  * [3] P. A. Zyla et al. [Particle Data Group], PTEP 2020, no.8, 083C01 (2020) doi:10.1093/ptep/ptaa104.
  * [4] E. Ma, Phys. Rev. D 73, 077301 (2006) doi:10.1103/PhysRevD.73.077301 [arXiv:hep-ph/0601225 [hep-ph]].
  * [5] P. F. de Salas, D. V. Forero, S. Gariazzo, P. Martínez-Miravé, O. Mena, C. A. Ternes, M. Tórtola and J. W. F. Valle, JHEP 02, 071 (2021) doi:10.1007/JHEP02(2021)071 [arXiv:2006.11237 [hep-ph]].
  * [6] M. C. Gonzalez-Garcia and M. Maltoni, Phys. Rept. 460, 1-129 (2008) doi:10.1016/j.physrep.2007.12.004 [arXiv:0704.1800 [hep-ph]].
  * [7] G. Bertone, D. Hooper and J. Silk, Phys. Rept. 405, 279-390 (2005) doi:10.1016/j.physrep.2004.08.031 [arXiv:hep-ph/0404175 [hep-ph]].
  * [8] J. Kubo, E. Ma and D. Suematsu, Phys. Lett. B 642, 18-23 (2006) doi:10.1016/j.physletb.2006.08.085 [arXiv:hep-ph/0604114 [hep-ph]].
  * [9] T. Toma and A. Vicente, JHEP 01, 160 (2014) doi:10.1007/JHEP01(2014)160 [arXiv:1312.2840 [hep-ph]].
  * [10] D. Aristizabal Sierra, J. Kubo, D. Restrepo, D. Suematsu and O. Zapata, Phys. Rev. D 79, 013011 (2009) doi:10.1103/PhysRevD.79.013011 [arXiv:0808.3340 [hep-ph]]; D. Suematsu, T. Toma and T. Yoshida, Phys. Rev. D 79, 093004 (2009) doi:10.1103/PhysRevD.79.093004 [arXiv:0903.0287 [hep-ph]]; A. Adulpravitchai, M. Lindner and A. Merle, Phys. Rev. D 80, 055031 (2009) doi:10.1103/PhysRevD.80.055031 [arXiv:0907.2147 [hep-ph]]; A. Vicente and C. E. Yaguna, JHEP 02, 144 (2015) doi:10.1007/JHEP02(2015)144 [arXiv:1412.2545 [hep-ph]]; A. Ahriche, A. Jueid and S. Nasri, Phys. Rev. D 97, no.9, 095012 (2018) doi:10.1103/PhysRevD.97.095012 [arXiv:1710.03824 [hep-ph]]; T. Hugle, M. Platscher and K. Schmitz, Phys. Rev. D 98, no.2, 023020 (2018) doi:10.1103/PhysRevD.98.023020 [arXiv:1804.09660 [hep-ph]]; S. Baumholzer, V. Brdar and P. Schwaller, JHEP 08, 067 (2018) doi:10.1007/JHEP08(2018)067 [arXiv:1806.06864 [hep-ph]]; D. Borah, P. S. B. Dev and A. Kumar, Phys. Rev. D 99, no.5, 055012 (2019) doi:10.1103/PhysRevD.99.055012 [arXiv:1810.03645 [hep-ph]]; A. Ahriche, A. Arhrib, A. Jueid, S. Nasri and A. de La Puente, Phys. Rev. D 101, no.3, 035038 (2020) doi:10.1103/PhysRevD.101.035038 [arXiv:1811.00490 [hep-ph]]; S. Baumholzer, V. Brdar, P. Schwaller and A. Segner, JHEP 09, 136 (2020) doi:10.1007/JHEP09(2020)136 [arXiv:1912.08215 [hep-ph]].
  * [11] R. S. Hundi, Phys. Rev. D 93, 015008 (2016) doi:10.1103/PhysRevD.93.015008 [arXiv:1510.02253 [hep-ph]].
  * [12] E. Ma, Annales Fond. Broglie 31, 285 (2006) [arXiv:hep-ph/0607142 [hep-ph]].
  * [13] G. Aad et al. [ATLAS], Phys. Rev. D 90, no.7, 072010 (2014) doi:10.1103/PhysRevD.90.072010 [arXiv:1408.5774 [hep-ex]].
  * [14] R. Akers et al. [OPAL], Z. Phys. C 67, 555-564 (1995) doi:10.1007/BF01553981.
  * [15] P. Abreu et al. [DELPHI], Z. Phys. C 73, 243-251 (1997) doi:10.1007/s002880050313.
  * [16] G. Aad et al. [ATLAS], Phys. Lett. B 801, 135148 (2020) doi:10.1016/j.physletb.2019.135148 [arXiv:1909.10235 [hep-ex]].
  * [17] G. Aad et al. [ATLAS], Phys. Lett. B 800, 135069 (2020) doi:10.1016/j.physletb.2019.135069 [arXiv:1907.06131 [hep-ex]].
  * [18] A. M. Sirunyan et al. [CMS], JHEP 06, 001 (2018) doi:10.1007/JHEP06(2018)001 [arXiv:1712.07173 [hep-ex]].
  * [19] M. Dam, SciPost Phys. Proc. 1, 041 (2019) doi:10.21468/SciPostPhysProc.1.041 [arXiv:1811.09408 [hep-ex]].
  * [20] J. I. Illana and T. Riemann, Phys. Rev. D 63, 053004 (2001) doi:10.1103/PhysRevD.63.053004 [arXiv:hep-ph/0010193 [hep-ph]]; E. O. Iltan and I. Turan, Phys. Rev. D 65, 013001 (2002) doi:10.1103/PhysRevD.65.013001 [arXiv:hep-ph/0106068 [hep-ph]]; J. I. Illana and M. Masip, Phys. Rev. D 67, 035004 (2003) doi:10.1103/PhysRevD.67.035004 [arXiv:hep-ph/0207328 [hep-ph]]; J. Cao, Z. Xiong and J. M. Yang, Eur. Phys. J. C 32, 245-252 (2004) doi:10.1140/epjc/s2003-01391-1 [arXiv:hep-ph/0307126 [hep-ph]]; E. O. Iltan, Eur. Phys. J. C 56, 113-118 (2008) doi:10.1140/epjc/s10052-008-0644-0 [arXiv:0802.1277 [hep-ph]]; M. J. Herrero, X. Marcano, R. Morales and A. Szynkman, Eur. Phys. J. C 78, no.10, 815 (2018) doi:10.1140/epjc/s10052-018-6281-3 [arXiv:1807.01698 [hep-ph]]; V. Cirigliano, K. Fuyuto, C. Lee, E. Mereghetti and B. Yan, JHEP 03, 256 (2021) doi:10.1007/JHEP03(2021)256 [arXiv:2102.06176 [hep-ph]].
  * [21] E. Arganda, A. M. Curiel, M. J. Herrero and D. Temes, Phys. Rev. D 71, 035011 (2005) doi:10.1103/PhysRevD.71.035011 [arXiv:hep-ph/0407302 [hep-ph]]; E. Arganda, M. J. Herrero, X. Marcano and C. Weiland, Phys. Rev. D 91, no.1, 015001 (2015) doi:10.1103/PhysRevD.91.015001 [arXiv:1405.4300 [hep-ph]]; E. Arganda, M. J. Herrero, X. Marcano, R. Morales and A. Szynkman, Phys. Rev. D 95, no.9, 095029 (2017) doi:10.1103/PhysRevD.95.095029 [arXiv:1612.09290 [hep-ph]]; N. H. Thao, L. T. Hue, H. T. Hung and N. T. Xuan, Nucl. Phys. B 921, 159-180 (2017) doi:10.1016/j.nuclphysb.2017.05.014 [arXiv:1703.00896 [hep-ph]]; Q. Qin, Q. Li, C. D. Lü, F. S. Yu and S. H. Zhou, Eur. Phys. J. C 78, no.10, 835 (2018) doi:10.1140/epjc/s10052-018-6298-7 [arXiv:1711.07243 [hep-ph]]; A. Vicente, Front. in Phys. 7, 174 (2019) doi:10.3389/fphy.2019.00174 [arXiv:1908.07759 [hep-ph]]; Z. N. Zhang, H. B. Zhang, J. L. Yang, S. M. Zhao and T. F. Feng, Phys. Rev. D 103, no.11, 115015 (2021) doi:10.1103/PhysRevD.103.115015 [arXiv:2105.09799 [hep-ph]].
  * [22] J. Herrero-Garcia, N. Rius and A. Santamaria, JHEP 11, 084 (2016) doi:10.1007/JHEP11(2016)084 [arXiv:1605.06091 [hep-ph]].
  * [23] C. Hagedorn, J. Herrero-García, E. Molinaro and M. A. Schmidt, JHEP 11, 103 (2018) doi:10.1007/JHEP11(2018)103 [arXiv:1804.04117 [hep-ph]].
  * [24] J. A. Casas and A. Ibarra, Nucl. Phys. B 618, 171-204 (2001) doi:10.1016/S0550-3213(01)00475-8 [arXiv:hep-ph/0103065 [hep-ph]].
  * [25] W. Grimus and L. Lavoura, Phys. Rev. D 66, 014016 (2002) doi:10.1103/PhysRevD.66.014016 [arXiv:hep-ph/0204070 [hep-ph]].
  * [26] G. ’t Hooft and M. J. G. Veltman, Nucl. Phys. B 44, 189-213 (1972) doi:10.1016/0550-3213(72)90279-9; G. Passarino and M. J. G. Veltman, Nucl. Phys. B 160, 151-207 (1979) doi:10.1016/0550-3213(79)90234-7.
  * [27] S. Heinemeyer et al. [LHC Higgs Cross Section Working Group], doi:10.5170/CERN-2013-004 [arXiv:1307.1347 [hep-ph]].
  * [28] N. Aghanim et al. [Planck], Astron. Astrophys. 641, A6 (2020) [erratum: Astron. Astrophys. 652, C4 (2021)] doi:10.1051/0004-6361/201833910 [arXiv:1807.06209 [astro-ph.CO]].
  * [29] A. M. Baldini et al. [MEG], Eur. Phys. J. C 76, no.8, 434 (2016) doi:10.1140/epjc/s10052-016-4271-x [arXiv:1605.05081 [hep-ex]].
  * [30] B. Aubert et al. [BaBar], Phys. Rev. Lett. 104, 021802 (2010) doi:10.1103/PhysRevLett.104.021802 [arXiv:0908.2381 [hep-ex]].
  * [31] M. Zeleny-Mora, J. L. Díaz-Cruz and O. Félix-Beltrán, [arXiv:2112.08412 [hep-ph]].