# Muon $g-2$ in a two-Higgs-doublet model with a type-II seesaw mechanism

Chuan-Hung Chen<EMAIL_ADDRESS>Department of Physics, National
Cheng-Kung University, Tainan 70101, Taiwan  Physics Division, National Center
for Theoretical Sciences, Taipei 10617, Taiwan  Cheng-Wei Chiang
<EMAIL_ADDRESS>Department of Physics and Center for Theoretical
Physics, National Taiwan University, Taipei 10617, Taiwan  Physics Division,
National Center for Theoretical Sciences, Taipei 10617, Taiwan  Takaaki Nomura
<EMAIL_ADDRESS>College of Physics, Sichuan University, Chengdu 610065,
China

###### Abstract

We study the two-Higgs-doublet model with type-II seesaw mechanism. In view of
constraints from the Higgs data, we consider the aligned two-Higgs-doublet
scheme and its effects on muon anomalous magnetic dipole moment, $a_{\mu}$,
including both one-loop and two-loop Barr-Zee type diagrams. Thanks to a
sizable trilinear scalar coupling, the Barr-Zee type diagrams mediated by the
Higgs triplet fields have a dominant effect on $a_{\mu}$. In particular,
unlike the usual two-Higgs-doublet models that require exotic Higgs bosons
light in mass, the masses of the corresponding particles in the model are of
${\cal O}(100)$ GeV. The doubly-charged Higgs boson presents a different decay
pattern from the usual Higgs triplet model and thus calls for a new collider
search strategy, such as multi-$\tau$ searches at the LHC.

##  I Introduction

A long-standing anomaly in particle physics is the muon anomalous magnetic
dipole moment (dubbed the muon $g-2$ anomaly) denoted by
$a_{\mu}\equiv(g-2)_{\mu}/2$, where the data and Standard Model (SM) show an
over $3\sigma$ disagreement. The E821 experiment at Brookhaven National Lab
(BNL) has presented a precision measurement of:

$a^{\rm exp}_{\mu}=116592089(63)\times 10^{-11}\,,$ (1)

with an uncertainty of $0.54$ ppm Bennett:2006fi . The current theoretical
estimate of $a_{\mu}$ within the SM has also reached a comparable precision of
$0.369$ ppm, and is shown to be Aoyama:2020ynm :

$\displaystyle a^{\rm SM}_{\mu}=116591810(43)\times 10^{-11}\,.$ (2)

The deviation between the experiment and the SM is $\Delta a_{\mu}=a^{\rm
exp}_{\mu}-a^{\rm SM}_{\mu}=279(76)\times 10^{-11}$ with an achievement of
$3.7\sigma$. The new muon $g-2$ measurement performed in the E989 Run 1
experiment at Fermilab, designed to have a precision of $0.14$ ppm, reports
its first measurement as Abi:2021gix :

$a^{\rm FNAL}_{\mu}=116592040(54)\times 10^{-11}\,.$ (3)

Combining all available measurements on the quantity, we now have a
$4.2\sigma$ deviation between experiment and SM expectation111The latest
lattice QCD calculation for the leading hadronic vacuum polarization from the
BMW collaboration is obtained as $a^{\rm LO-HVP}_{\mu}=707.5(5.5)\times
10^{-10}$, which leads to a larger $a_{\mu}$, can be found in Borsanyi:2020mff
., accentuating the muon $g-2$ anomaly.

On the other hand, since the discovery of Higgs boson at the LHC in 2012
summer, measurements of the Higgs signal strengths, commonly used as a measure
of deviations from the SM, have been improving over the years. They are found
to be quite consistent with the SM expectations and, hence, models with
extensions in the scalar sector are severely constrained. One possibility for
a new physics (NP) model to achieve such a good agreement with the SM in the
Higgs couplings while having exotic Higgs bosons of mass at ${\cal O}(100)$
GeV scale is when the model shows the so-called alignment limit Gunion:2002zf
; Carena:2013ooa ; Bernon:2015qea .

In this work, we study the contributions of a model with an extended scalar
sector to the muon $g-2$ when the relevant theoretical and experimental
constraints are taken into account. One purpose is to revisit the two-Higgs-
doublet models (2HDMs), where the earlier studies can be found in Refs.
Leveille:1977rc ; Haber:1978jt ; Krawczyk:1996sm ; Dedes:2001nx ;
Cheung:2001hz ; Larios:2001ma ; Chen:2001kn ; Arhrib:2001xx ; Krawczyk:2002df
; Gunion:2008dg ; Broggio:2014mna ; Ilisie:2015tra ; Han:2015yys ; Li:2018aov
; Rose:2020nxm ; Sabatta:2019nfg ; Jana:2020pxx ; Li:2020dbg ; Chen:2021rnl ;
Han:2018znu ; Botella:2020xzf . It is known that to explain the muon $g-2$ in
this framework, the new scalar or pseudoscalar boson are required to be as
light as ${\cal O}(10)$ GeV. Although such a parameter space is still allowed
by the current data, it is of interest to probe the scenarios where the new
scalar masses can be more relaxed and of $\sim{\cal O}(100)$ GeV by further
extending the scalar sector. More importantly, such a new extension should
also address some other unsolved issues, such as the origin of neutrino mass,
that the simple 2HDMs cannot accommodate.

To achieve the above-mentioned goals, we consider the 2HDM with type-II seesaw
mechanism Chen:2014xva ; Ouazghour:2018mld . In addition to the SM Higgs
doublet, the scalar sector contains another complex doublet and a complex
triplet. Moreover, we will consider the so-called aligned two-Higgs-doublet
scheme (A2HDS), where the Yukawa couplings of the two Higgs doublets to the SM
fermions are proportional to each other and one of the neutral physical Higgs
boson is SM-like. The A2HDS has the interesting feature that it reduces to
various 2HDM types by taking proper limits on the alignment parameters. With a
small vacuum expectation value (VEV) induced by electroweak symmetry breaking
from the two Higgs doublets, the Higgs triplet in the model provides Majorana
mass to neutrinos through the so-called type-II seesaw mechanism Magg:1980ut ;
Lazarides:1980nt ; Schechter:1980gr ; Cheng:1980qt ; Mohapatra:1980yp ;
Bilenky:1980cx .

It is found that rather than a simple combination of the 2HDM and the type-II
seesaw model (also called the Higgs triplet model or HTM), the model presents
several interesting features:

  1. 1.

The coupling between the heavier neutral Higgs boson in the 2HDM and the
doubly-charged Higgs boson in the HTM can significantly enhance the muon $g-2$
through two-loop Barr-Zee type diagrams Barr:1990vd ; Chang:2000ii , even when
the heavier neutral Higgs mass is $\sim{\cal O}(100)$ GeV.

  2. 2.

The Higgs triplet VEV is now determined by three lepton number-violating
parameters instead of just one in the simple HTM. As a result of the extra
freedom, these parameters are not necessarily of the same order as the Higgs
triplet VEV Chen:2014xva .

  3. 3.

With a sizable Higgs triplet VEV, the doubly-charged Higgs boson shows a
richer decay pattern. As a result, the doubly-charged Higgs boson can evade
the recent ATLAS lower bound of 350 GeV in pair production Aad:2021lzu . In
addition to the like-sign diboson channel, the doubly-charged Higgs boson can
also be probed via channels involving the light charged Higgs boson.

The paper is organized as follows. In Sec. II, we derive the Yukawa couplings
in the A2HDS and show the relations between the scheme and the various types
of 2HDMs with $\mathbb{Z}_{2}$ symmetry. The mass-square relations of the
triplet Higgs bosons are discussed, and the CP-even neutral Higgs couplings
with the charged Higgses are given. In Sec. III, we discuss the results of
one-loop and the dominant two-loop Barr-Zee type diagrams. Using the bounded
parameters, we present the detailed numerical analysis and discussion in Sec.
IV. Sec. V summarizes our findings in this work. The full scalar mass matrices
and their approximations in the limit of neglecting $v_{\Delta}$ are given in
appendix A.

##  II Model and Interactions

We consider a model where the scalar sector is extended with a doublet with
$Y=1/2$ and a complex triplet with $Y=1$. In the following, we discuss the
general Yukawa interactions and scalar potential in this model.

###  II.1 Scalar potential and the trilinear scalar couplings

Since the scalar sector is an extension of 2HDM or of type-II seesaw, in the
following, we briefly discuss the essential parts for our analysis. First, as
we will assume negligibly small mixing between the doublet fields and the
triplet field, it is useful to go to the Higgs basis in the usual 2HDM,
defined by:

$\begin{pmatrix}H_{1}\\\ H_{2}\\\
\end{pmatrix}=\begin{pmatrix}c_{\beta}&s_{\beta}&\\\
-s_{\beta}&c_{\beta}&\end{pmatrix}\begin{pmatrix}\Phi_{1}\\\ \Phi_{2}\\\
\end{pmatrix}\,,$ (4)

where $v_{i}$ is the VEV of $\Phi_{i}$ ($i=1,2$),
$c_{\beta}~{}(s_{\beta})=\cos\beta~{}(\sin\beta)$, $\tan\beta=v_{2}/v_{1}$ and
$v=\sqrt{v^{2}_{1}+v^{2}_{2}}\simeq 246$ GeV. Written in terms of field
components, the Higgs doublets $H_{1,2}$ and triplet $\Delta$ are:

$\displaystyle H_{1}$ $\displaystyle=\begin{pmatrix}G^{+}\\\
(v+H^{0}_{1}+iG^{0})/\sqrt{2}\\\
\end{pmatrix}\,,~{}~{}H_{2}=\begin{pmatrix}H^{+}\\\
(H^{0}_{2}+iA^{0})/\sqrt{2}\\\ \end{pmatrix}\,,$ $\displaystyle\Delta$
$\displaystyle=\begin{pmatrix}\delta^{+}/\sqrt{2}&\delta^{++}\\\
(v_{\Delta}+\delta^{0}+i\eta^{0})/\sqrt{2}&~{}-\delta^{+}/\sqrt{2}\\\
\end{pmatrix}\,.$ (5)

In the conventional CP-conserving 2HDM, $G^{\pm(0)}$ are the Goldstone bosons,
and $H^{\pm}$ and $A^{0}$ are the charged Higgs boson and the CP-odd
pseudoscalar physical states. In addition, the CP-even scalars $H^{0}_{1}$ and
$H^{0}_{2}$ mix to give their mass eigenstates via:

$\begin{pmatrix}H\\\ h\\\
\end{pmatrix}=\begin{pmatrix}c_{\beta-\alpha}&-s_{\beta-\alpha}\\\
s_{\beta-\alpha}&c_{\beta-\alpha}\end{pmatrix}\begin{pmatrix}H^{0}_{1}\\\
H^{0}_{2}\\\ \end{pmatrix}\,,$ (6)

where $h$ is the 125-GeV SM-like Higgs boson,
$c_{\beta-\alpha}=\cos(\beta-\alpha)$, $s_{\beta-\alpha}=\sin(\beta-\alpha)$,
and $\alpha$ is the mixing angle of $\Phi^{0}_{1}$ and $\Phi^{0}_{2}$.
Although $\delta^{0}$, $\eta^{0}$ and $\delta^{\pm}$ generally mix with
$(H,h)$, $(G^{0},A^{0})$, and $(G^{\pm},H^{\pm})$, respectively, such mixings
are small and phenomenologically negligible when $v_{\Delta}\ll 1$ GeV, as is
the case considered in this work. Hence, it is a good approximation to take
$h$, $H$, $A^{0}$, and $H^{\pm}$ as the physical states.

The scalar potential of two-Higgs-doublet fields and the Higgs triplet field
under the $SU(2)_{L}\times U(1)_{Y}$ gauge symmetry is given by:

$V=V(\Phi_{1},\Phi_{2})+V(\Delta)+V(\Phi_{1},\Phi_{2},\Delta)\,,$ (7)

where each term is more explicitly given by Chen:2014xva ; Ouazghour:2018mld

$\displaystyle V(\Phi_{1},\Phi_{2})=$ $\displaystyle
m^{2}_{1}\Phi^{\dagger}_{1}\Phi_{1}+m^{2}_{2}\Phi^{\dagger}_{2}\Phi_{2}-m^{2}_{12}(\Phi^{\dagger}_{1}\Phi_{2}+\mbox{H.c.})+\frac{1}{2}\lambda_{1}(\Phi^{\dagger}_{1}\Phi_{1})^{2}$
$\displaystyle+\frac{1}{2}\lambda_{2}(\Phi^{\dagger}_{2}\Phi_{2})^{2}+\lambda_{3}\Phi^{\dagger}_{1}\Phi_{1}\Phi^{\dagger}_{2}\Phi_{2}+\lambda_{4}\Phi^{\dagger}_{1}\Phi_{2}\Phi^{\dagger}_{2}\Phi_{1}+\left[\frac{1}{2}\lambda_{5}(\Phi^{\dagger}_{1}\Phi_{2})^{2}\right.$
$\displaystyle\left.+\lambda_{6}(\Phi^{\dagger}_{1}\Phi_{1})(\Phi^{\dagger}_{2}\Phi_{1})+\lambda_{7}(\Phi^{\dagger}_{2}\Phi_{2})(\Phi^{\dagger}_{2}\Phi_{1})+\mbox{H.c.}\right]\,,$
(8) $\displaystyle V(\Delta)=$ $\displaystyle
m^{2}_{\Delta}Tr\Delta^{\dagger}\Delta+\lambda_{\Delta
1}(Tr\Delta^{\dagger}\Delta)^{2}+\lambda_{\Delta
2}Tr(\Delta^{\dagger}\Delta)^{2}\,,$ (9) $\displaystyle
V(\Phi_{1},\Phi_{2},\Delta)=$
$\displaystyle\left(\mu_{1}\Phi^{T}_{1}i\tau_{2}\Delta^{\dagger}\Phi_{1}+\mu_{2}\Phi^{T}_{2}i\tau_{2}\Delta^{\dagger}\Phi_{2}+\mu_{3}\Phi^{T}_{1}i\tau_{2}\Delta^{\dagger}\Phi_{2}+\mbox{H.c.}\right)$
$\displaystyle+\left(\lambda_{8}\Phi^{\dagger}_{1}\Phi_{1}+\lambda_{9}\Phi^{\dagger}_{2}\Phi_{2}+(\lambda_{12}\Phi^{\dagger}_{1}\Phi_{2}+\mbox{H.c.})\right)Tr\Delta^{\dagger}\Delta$
$\displaystyle+\lambda^{\prime}_{8}\Phi^{\dagger}_{1}\Delta\Delta^{\dagger}\Phi_{1}+\lambda^{\prime}_{9}\Phi^{\dagger}_{2}\Delta\Delta^{\dagger}\Phi_{2}+(\lambda^{\prime}_{12}\Phi^{\dagger}_{1}\Delta\Delta^{\dagger}\Phi_{2}+\mbox{H.c.})\,.$
(10)

In terms of the Higgs basis, the minimal conditions for the VEVs of $H_{1,2}$
and $\Delta$ can be obtained as:

$\displaystyle\left(\frac{c^{4}_{\beta}\lambda_{1}+s^{4}_{\beta}\lambda_{2}}{2}+c^{2}_{\beta}s^{2}_{\beta}\lambda_{345}+2c^{3}_{\beta}s_{\beta}\lambda_{6}+2c_{\beta}s^{3}_{\beta}\lambda_{7}\right)v^{2}=-c^{2}_{\beta}m^{2}_{1}-s^{2}_{\beta}m^{2}_{2}+2c_{\beta}s_{\beta}m^{2}_{12}$
$\displaystyle+\sqrt{2}v_{\Delta}\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\right)-\frac{v^{2}_{\Delta}}{2}\left(c^{2}_{\beta}\bar{\lambda}_{8}+s^{2}_{\beta}\bar{\lambda}_{9}\right)\,,$
(11a) $\displaystyle
c_{\beta}s_{\beta}(m^{2}_{2}-m^{2}_{1})-m^{2}_{12}c_{2\beta}+\frac{v^{2}}{2}\left(-c^{3}_{\beta}s_{\beta}\lambda_{1}+c_{\beta}s^{3}_{\beta}\lambda_{2}\right)+\frac{\lambda_{345}v^{2}}{2}c_{\beta}s_{\beta}c_{2\beta}+\frac{\lambda_{6}v^{2}}{2}(-3c^{2}_{\beta}+c^{4}_{\beta})$
$\displaystyle+\frac{\lambda_{7}v^{2}}{2}(3c^{2}_{\beta}-s^{4}_{\beta})=-\sqrt{2}c_{\beta}s_{\beta}v_{\Delta}(\mu_{1}-\mu_{2})+\frac{c_{2\beta}v_{\Delta}\mu_{3}}{\sqrt{2}}+\frac{\bar{\lambda}_{8}+\bar{\lambda}_{9}}{2}c_{\beta}s_{\beta}v^{2}_{\Delta}\,,$
(11b)
$\displaystyle\left(m^{2}_{\Delta}+\frac{\bar{\lambda}_{8}}{2}c^{2}_{\beta}v^{2}+\frac{\bar{\lambda}_{9}}{2}s^{2}_{\beta}v^{2}+(\lambda_{\Delta
1}+\lambda_{\Delta
2})v^{3}_{\Delta}\right)v_{\Delta}=\frac{v^{2}}{2}\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}+s_{2\beta}\bar{\lambda}_{12}\right)\,,$
(11c)

where $\bar{\lambda}_{8}=\lambda_{8}+\lambda^{\prime}_{8}$,
$\bar{\lambda}_{9}=\lambda_{9}+\lambda^{\prime}_{9}$,
$\bar{\lambda}_{12}=\lambda_{12}+\lambda^{\prime}_{12}$, and the VEV of
$\Delta$ is denoted by $v_{\Delta}$. These relations are useful to simplify
the expressions of scalar masses and trilinear scalar couplings. We note that
since no discrete symmetry is imposed in the 2HDM, both $H_{1,2}$ Higgs
doublets are indistinguishable. It is simpler to directly use the Higgs basis
in the scalar potential. To compare our results with those given in Ref.
Ouazghour:2018mld , here we employ the generic Higgs flavor basis, which is
used in Ref. Ouazghour:2018mld . Nevertheless, we show the more compact
expressions in appendix B.

If we drop the small effect from $(\lambda_{\Delta 1}+\lambda_{\Delta
2})v^{3}_{\Delta}$, the Higgs triplet VEV can be obtained as:

$v_{\Delta}\approx\frac{v^{2}}{\sqrt{2}}\frac{c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}}{\tilde{M}^{2}_{\Delta}}\,,$
(12)

with
$\tilde{M}^{2}_{\Delta}=m^{2}_{\Delta}+(\bar{\lambda}_{8}c^{2}_{\beta}v^{2}+\bar{\lambda}_{9}s^{2}_{\beta}v^{2}+\bar{\lambda}_{12}s_{2\beta})v^{2}/2$.
Since $v_{\Delta}$ is bounded by the electroweak precision measurement, and
with the exception of the neutrino mass, its effect is irrelevant to the
current study. Precision measurement of the electroweak $\rho$ parameter gives
a constraint that $v_{\Delta}\lesssim 8$ GeV PDG . To illustrate the
importance of trilinear scalar couplings between the two Higgs doublets and
the Higgs triplet on the muon $g-2$, we take $v_{\Delta}\sim
O(10^{-3}-10^{-4})$ GeV. The considered parameter region can be easily
achieved. For instance, using
$c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\approx
10^{-3}$ GeV and $\tilde{M}_{\Delta}=500$ GeV, we obtain $v_{\Delta}\approx
1.7\times 10^{-4}$ GeV. As a result, the mass mixings of scalars,
pseudoscalars, and charged scalars between $\Phi_{i}$ and $\Delta$ are
phenomenologically negligible, justifying our earlier assumption. Due to the
doublet-triplet coupling terms in Eq. (10), the $G^{0,\pm}$ and
$H^{\pm}(A^{0})$ in the Higgs basis of 2HDM are not the Goldstone modes and
the physical states, and the $(G^{0},A^{0},\eta^{0})$ and
$(G^{\pm},H^{\pm},\delta^{\pm})$ states will respectively mix. The only
nonmixing state is the doubly-charged Higgs, where from Eq. (9), its mass can
be expressed as:

$m^{2}_{\delta^{\pm\pm}}=m^{2}_{\Delta}+\lambda_{\Delta
1}v^{2}_{\Delta}+\frac{v^{2}}{2}\left(c^{2}_{\beta}\lambda_{8}+s^{2}_{\beta}\lambda_{9}+\lambda_{12}s_{2\beta}\right)\,.$
(13)

It can be seen that the new doublet-triplet couplings terms shift the
$\delta^{\pm\pm}$ mass.

The detailed discussions for the scalar, pseudoscalar, and charged scalar mass
matrix are given in appendix A. We summarized the characteristics as follows:
from Eqs. (49) and (54), it can be seen that $m^{2}_{G^{0}}$ and
$m^{2}_{G^{\pm}}$ are proportional to $v^{2}_{\Delta}$. With $v_{\Delta}\sim
O(10^{-4})$ GeV, their values can be dropped. If $\mu_{3}=\tan
2\beta(\mu_{1}-\mu_{2})$ is required, we can find that the mixing matrix
elements of $G^{0}A^{0}$ and $G^{-}H^{+}$ become $O(v^{2}_{\Delta})$, and the
$G^{0}(A^{0})\eta^{0}$, $G^{-}(H^{-})\delta^{+}$, and $h(H)\delta^{0}$ matrix
elements are of $O(v_{\Delta})$. In comparison with the mass-square elements
of other massive particles, their mixing effects are small. Although the small
mixing effects can have important influence on some processes, e.g.
$\delta^{\pm\pm}\to W^{\pm}H^{\pm}$ can be induced, their influence on the
muon $g-2$ can be indeed neglected. Hence, when we numerically estimate the
muon $g-2$, we take $h(H)$, $H^{\pm}(\delta^{\pm})$ and $A^{0}$ as the
physical states; however, for other processes, one can take the mixing effects
into account if necessary.

The new doublet-triplet couplings can cause the triplet scalar mass
splittings, and the mass differences can be found as:

$\displaystyle m^{2}_{\delta^{\pm}}-m^{2}_{\delta^{\pm\pm}}$
$\displaystyle=\frac{v^{2}_{\Delta}}{2}\lambda_{\Delta
2}+\frac{v^{2}}{4}\left(c^{2}_{\beta}\lambda^{\prime}_{8}+s^{2}_{\beta}\lambda^{\prime}_{9}+\lambda^{\prime}_{12}s_{2\beta}\right)\,,$
$\displaystyle m^{2}_{\delta^{0}}-m^{2}_{\delta^{\pm\pm}}$
$\displaystyle=(2\lambda_{\Delta 1}+3\lambda_{\Delta
2})v^{2}_{\Delta}+\frac{v^{2}}{2}\left(c^{2}_{\beta}\lambda^{\prime}_{8}+s^{2}_{\beta}\lambda^{\prime}_{9}+\lambda^{\prime}_{12}s_{2\beta}\right)\,.$
(14)

It can be seen that the mass split can be or be less than $O(100)$ GeV .

The trilinear interactions among a neutral Higgs boson and two charged Higgs
bosons are given by:

${\cal
L}_{H^{0}_{i}SS}=-v\left[\lambda_{H^{0}_{i}\delta^{--}\delta^{++}}H^{0}_{i}\delta^{--}\delta^{++}+\lambda_{H^{0}_{i}\delta^{-}\delta^{+}}H^{0}_{i}\delta^{-}\delta^{+}+\lambda_{H^{0}_{i}H^{-}H^{+}}H^{0}_{i}H^{-}H^{+}\right]\,,$
(15)

where the couplings are written as:

$\displaystyle\begin{split}\lambda_{H^{0}_{1}\delta^{--}\delta^{++}}=&\lambda_{8}c^{2}_{\beta}+\lambda_{9}s^{2}_{\beta}+\lambda_{12}s_{2\beta}\,,\\\
\lambda_{H^{0}_{2}\delta^{--}\delta^{++}}=&(-\lambda_{8}+\lambda_{9})c_{\beta}s_{\beta}+\lambda_{12}c_{2\beta}\,,\\\
\lambda_{H^{0}_{1}\delta^{-}\delta^{+}}=&\left(\lambda_{8}+\frac{\lambda^{\prime}_{8}}{2}\right)c^{2}_{\beta}+\left(\lambda_{9}+\frac{\lambda^{\prime}_{9}}{2}\right)s^{2}_{\beta}+\left(\lambda_{12}+\frac{\lambda^{\prime}_{12}}{2}\right)s_{2\beta}\,,\\\
\lambda_{H^{0}_{2}\delta^{-}\delta^{+}}=&-\left(\lambda_{8}-\lambda_{9}+\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2}\right)+\left(\lambda_{12}+\frac{\lambda^{\prime}_{12}}{2}\right)c_{2\beta}\,,\\\
\lambda_{H^{0}_{1}H^{-}H^{+}}=&[\lambda_{1}+\lambda_{2}-2(\lambda_{4}+\lambda_{5})]c^{2}_{\beta}s^{2}_{\beta}+\lambda_{3}(c^{4}_{\beta}+s^{4}_{\beta})v-2(\lambda_{6}-\lambda_{7})c_{\beta}s_{\beta}c_{2\beta}~{},\\\
\lambda_{H^{0}_{2}H^{-}H^{+}}=&-\lambda_{1}c_{\beta}s^{3}_{\beta}+\lambda_{2}c^{3}_{\beta}s_{\beta}-\lambda_{345}c_{\beta}s_{\beta}c_{2\beta}-\lambda_{6}(3c^{2}_{\beta}s^{2}_{\beta}-s^{4}_{\beta})+\lambda_{7}\left(c^{4}_{\beta}-3c^{2}_{\beta}s^{2}_{\beta}\right)\,.\end{split}$
(16)

In the alignment limit of $c_{\beta-\alpha}=0$, the $h$ and $H$ trilinear
terms can be easily obtained by the replacement of $H^{0}_{1}\to h$ and
$H^{0}_{2}\to-H$. Therefore, the corresponding trilinear couplings have the
relations:

$\displaystyle\lambda_{h\delta^{--(-)}\delta^{++(+)}}$
$\displaystyle=\lambda_{H^{0}_{1}\delta^{--(-)}\delta^{++(+)}}\,,~{}~{}\lambda_{hH^{-}H^{+}}=\lambda_{H^{0}_{1}H^{-}H^{+}}\,,$
$\displaystyle\lambda_{H\delta^{--(-)}\delta^{++(+)}}$
$\displaystyle=-\lambda_{H^{0}_{2}\delta^{--(-)}\delta^{++(+)}}\,,~{}~{}\lambda_{HH^{-}H^{+}}=-\lambda_{H^{0}_{2}H^{-}H^{+}}\,.$
(17)

We note that the pseudoscalar $A^{0}$ does not couple to the charged scalars
in the CP-conserving case.

###  II.2 Yukawa interactions

The most general Yukawa couplings in the model are given by:

$\displaystyle-{\cal L}_{Y}$
$\displaystyle=\bar{Q}_{L}Y^{d}_{1}D_{R}\Phi_{1}+\bar{Q}_{L}Y^{d}_{2}D_{R}\Phi_{2}+\bar{Q}_{L}Y^{u}_{1}U_{R}\tilde{\Phi}_{1}+\bar{Q}_{L}Y^{u}_{2}U_{R}\tilde{\Phi}_{2}$
$\displaystyle+\bar{L}Y^{\ell}_{1}\ell_{R}\Phi_{1}+\bar{L}Y^{\ell}_{2}\ell_{R}\Phi_{2}+\frac{1}{2}L^{T}C{\bf
y}^{\nu}\,i\tau_{2}\,\Delta L+{\rm H.c.}\,,$ (18)

where the flavor indices are suppressed, ${\bf y^{\nu}}$ is a symmetric
matrix, $Q_{L}$ $(L)$ denotes the quark (lepton) doublets, $q_{R}$
$(\ell_{R})$ denotes the quark (lepton) singlets, $Y^{f}_{1,2}$ with
$f=u,d,\ell$ are respectively the Yukawa matrices for the up-type quarks,
down-type quarks, and charged leptons, $C$ is the charge conjugation operator,
and $\tilde{\Phi}_{i}\equiv i\tau_{2}\Phi^{*}_{i}$ with $\tau_{2}$ being the
Pauli matrix.

Since $\Phi_{1}$ and $\Phi_{2}$ simultaneously couple to each type of
fermions, flavor-changing neutral currents (FCNCs) naturally arise at tree
level. The FCNC effects are usually suppressed by introducing, for example, a
$Z_{2}$ discrete symmetry Glashow:1976nt . In this case, the 2HDM can be
categorized into Type-I Haber:1978jt ; Hall:1981bc , Type-II Donoghue:1978cj ;
Hall:1981bc , Type-X, and Type-Y Barger:1989fj ; Savage:1991qh ;
Grossman:1994jb ; Akeroyd:1994ga . See Ref. Branco:2011iw  for a detailed
review. In addition to the above-mentioned schemes in 2HDM, the tree-level
FCNCs can also be avoided by imposing a certain relation between $Y_{1}^{f}$
and $Y_{2}^{f}$, where $f=u$, $d$, and $\ell$. The A2HDS assumes the relation
$Y^{f}_{2}=\xi_{f}Y^{f}_{1}$, where $\xi_{f}$ is a proportionality constant
Pich:2009sp . Alternatively, one may also impose the condition
$Y^{f}_{2}=N_{I}Y^{f}_{1}N^{\dagger}_{I}$ Ahn:2010zza ; Chen:2011wp ;
Chen:2016xju , where the possible $N_{I}$ matrices can be found in Ref.
Ahn:2010zza . In this work, we are considering the A2HDS.

With the assumed VEVs of $\Phi_{i}$, $Y^{f}_{1}$ and $Y^{f}_{2}$ in Eq. (18)
can be linearly combined to form two matrices:

$\displaystyle X_{f}$ $\displaystyle=c_{\beta}Y^{f}_{1}+s_{\beta}Y^{f}_{2}\,,$
$\displaystyle Z_{f}$
$\displaystyle=-s_{\beta}Y^{f}_{1}+c_{\beta}Y^{f}_{2}\,,$ (19)

so that $X_{f}(Z_{f})$ is associated with the doublet $H_{1(2)}$, and the
fermion mass matrix can be obtained as $M_{f}=X_{f}v/\sqrt{2}$. Moreover,
$M_{f}$ can be diagonalized by the unitary matrices $U^{f}_{L,R}$ in the way
$M^{\rm dia}_{f}=U^{f}_{L}M_{f}U^{f{\dagger}}_{R}$. If $Y^{f}_{1}$ and
$Y^{f}_{2}$ are two linearly independent matrices and cannot be diagonalized
simultaneously, then tree-level FCNCs can arise due to the $Z_{f}$ couplings
because its off-diagonal elements cannot be removed when $X_{f}$ is
diagonalized. When the A2HDS relation $Y^{f}_{2}=\xi_{f}Y^{f}_{1}$ is taken,
the Yukawa matrix can be related to the mass matrix as:

$\displaystyle Y^{f}_{1}$
$\displaystyle=\frac{\sqrt{2}}{c_{\beta}v}\frac{1}{1+\xi_{f}t_{\beta}}M_{f}\,.$
(20)

As a result, both $Y^{f}_{1}$ and $Y^{f}_{2}$ now can be diagonalized
simultaneously, and the $H^{0}_{1,2}$ FCNCs are suppressed at the tree level.

For simplicity, we only concentrate on the CP-conserving case and assume
$\xi_{f}$ to be real, though they can generally be complex. Using Eqs. (4) and
$(\ref{eq:2HDM_mass_Basis})$ and the notations used in Pich:2009sp , the mass
terms and Yukawa interactions with $h$, $H$, $A^{0}$, and $H^{\pm}$ are found
to be:

$\displaystyle-{\cal L}_{Y}=$
$\displaystyle\sum_{f=u,d,\ell}\left[\bar{f}_{L}M^{\rm
dia}_{f}f_{R}+\frac{c_{\beta-\alpha}-\zeta_{f}s_{\beta-\alpha}}{v}\bar{f}_{L}M^{\rm
dia}_{f}f_{R}H+\frac{s_{\beta-\alpha}+\zeta_{f}c_{\beta-\alpha}}{v}\bar{f}_{L}M^{\rm
dia}_{f}f_{R}h\right]$
$\displaystyle+\sum_{f=u,d,\ell}\frac{\zeta_{f}s_{f}}{v}\bar{f}_{L}M^{\rm
dia}_{f}f_{R}(iA^{0})+\frac{\sqrt{2}}{v}\bar{u}\left(V_{\rm
CKM}\zeta_{d}M^{\rm dia}_{d}P_{R}-M^{\rm dia}_{u}\zeta^{*}_{u}V_{\rm
CKM}P_{L}\right)dH^{+}$ $\displaystyle+\frac{\sqrt{2}}{v}\bar{\nu}\left(V_{\rm
PMNS}\zeta_{\ell}M^{\rm dia}_{\ell}P_{R}\right)\ell H^{+}+\mbox{H.c.}\,,$ (21)

where $s_{d,\ell}=+1$, $s_{u}=-1$; $P_{R(L)}=(1\pm\gamma_{5})/2$ are the
chirality projection operators, $V_{\rm CKM}=V^{u}_{L}V^{d{\dagger}}_{L}$ is
the Cabibbo-Kobayashi-Maskawa (CKM) matrix Cabibbo:1963yz ; Kobayashi:1973fv ,
$V_{\rm PMNS}=V^{\nu}_{L}V^{\ell{\dagger}}_{L}$ is the Pontecorvo-Maki-
Nakagawa-Sakata (PMNS) matrix Pontecorvo:1957cp ; Maki:1962mu , and

$\zeta_{f}=\frac{\xi_{f}-t_{\beta}}{1+\xi_{f}t_{\beta}}\,,$ (22)

with $t_{\beta}=\tan\beta$. In general, $\zeta_{f}$ can be complex numbers as
$\xi_{f}$ Pich:2009sp , and their magnitudes can be large without requiring a
large $\tan\beta$. In this study, we only focus on the CP-conserving case. The
$A^{0}$ and $H^{\pm}$ Yukawa couplings do not depend on $c_{\beta-\alpha}$
$(s_{\beta-\alpha})$. For comparison, we show in Table 1 the vanishing and
non-vanishing $Y^{f}_{1,2}$ for various 2HDM types and the associated
$\zeta_{f}$. In particular, $\zeta_{f}$ in Type-I, -II, -X, and -Y can be
obtained from the A2HDS by taking an appropriate limit of $\xi_{f}$ and thus
$\zeta_{f}$:

$\displaystyle Y^{f}_{2}$
$\displaystyle=0:~{}\xi_{f}=0\,~{}~{}\rightarrow~{}\zeta_{f}=-t_{\beta}\,,$
$\displaystyle Y^{f}_{1}$
$\displaystyle=0:~{}\xi_{f}=\infty~{}\rightarrow~{}\zeta_{f}=t^{-1}_{\beta}\,.$
(23)

Table 1:  Vanishing (mark by 0) and non-vanishing (marked by $\times$) Yukawa matrices of various 2HDM types and the corresponding $\zeta_{f}$.  |  $Y^{d}_{1}$ |  $Y^{u}_{1}$ |  $Y^{\ell}_{1}$ |  $Y^{d}_{2}$ |  $Y^{u}_{2}$ |  $Y^{\ell}_{2}$ |  $\zeta_{u}$ |  $\zeta_{d}$ |  $\zeta_{\ell}$
---|---|---|---|---|---|---|---|---|---
Type I | 0 | 0 | 0 | $\times$ | $\times$ | $\times$ | $t^{-1}_{\beta}$ | $t^{-1}_{\beta}$ | $t^{-1}_{\beta}$
Type II | $\times$ | $0$ | $\times$ | $0$ | $\times$ | $0$ | $t^{-1}_{\beta}$ | $-t_{\beta}$ | $-t_{\beta}$
Type X | $0$ | $0$ | $\times$ | $\times$ | $\times$ | $0$ | $t^{-1}_{\beta}$ | $t^{-1}_{\beta}$ | $-t_{\beta}$
Type Y | $\times$ | $0$ | $0$ | $0$ | $\times$ | $\times$ | $t^{-1}_{\beta}$ | $-t_{\beta}$ | $t^{-1}_{\beta}$
A2HDS | $\times$ | $\times$ | $\times$ | $\times$ | $\times$ | $\times$ | $\displaystyle\frac{\xi_{u}-t_{\beta}}{1+\xi_{u}t_{\beta}}$ | $\displaystyle\frac{\xi_{d}-t_{\beta}}{1+\xi_{d}t_{\beta}}$ | $\displaystyle\frac{\xi_{\ell}-t_{\beta}}{1+\xi_{\ell}t_{\beta}}$

Since the SM-like Higgs couplings generally depend on $c_{\beta-\alpha}$
$(s_{\beta-\alpha})$ and $\zeta_{f}$, the current Higgs production and decay
measurements put stringent constraints on the value of $c_{\beta-\alpha}$.
Here we simply take the alignment limit with $\beta-\alpha\approx\pi/2$
Gunion:2002zf , i.e., $c_{\beta-\alpha}\to 0$ ($s_{\beta-\alpha}\to 1$). As a
result, the $H$ and $A^{0}$ couplings to the SM fermions have the same
magnitude and are dictated by $\zeta_{f}$. In this work, we demonstrate how a
large $\zeta_{\ell}$ can affect the muon $g-2$ when $m_{H}>m_{h}$ and the
$H\delta^{++}\delta^{--}$ and $H\delta^{+}\delta^{-}$ couplings are present.

Using the component fields of the Higgs triplet shown in Eq. (5), the neutrino
mass and lepton Yukawa interactions with the triplet fields are given by:

$\displaystyle{\cal
L}_{Y}\supset\frac{1}{2}\overline{\nu^{C}_{L}}M_{\nu}\nu_{L}+\frac{1}{2}\overline{\nu^{C}_{L}}M_{\nu}\nu_{L}\frac{\delta^{0}+i\eta^{0}}{v_{\Delta}}-\overline{\nu^{C}_{L}}\frac{{\bf
y^{\nu}}}{\sqrt{2}}\ell_{L}\delta^{+}-\frac{1}{2}\overline{\ell^{C}_{L}}{\bf
y^{\nu}}\ell_{L}\delta^{++}+\mbox{H.c.}\,,$ (24)

where $f^{C}=C\gamma^{0}f^{*}$ and $M_{\nu}={\bf y^{\nu}}v_{\Delta}/\sqrt{2}$
is the neutrino mass matrix. In order to fit the neutrino data, the values of
$(M_{\nu})_{ij}$ has to be of $O(10^{-3}-10^{-2})$ eV Chen:2014xva ;
Chen:2017gvf ; Chen:2020ptg . In the type-II seesaw model with the assumed
triplet VEV $v_{\Delta}\sim{\cal O}(10^{-3}-10^{-4})$ GeV, the Yukawa
couplings ${\bf y}^{\nu}_{ij}$ are very small, $\lesssim{\cal O}(10^{-7})$.
Therefore, $\delta^{\pm}$ and $\delta^{\pm\pm}$ of ${\cal O}(10^{2})$ GeV mass
have negligible effects on most lepton processes.

As we will numerically show below, the Yukawa couplings as well as the
trilinear scalar couplings $H\delta^{--(-)}\delta^{++(+)}$ and $HH^{-}H^{+}$,
arising from the scalar potential, play important roles in producing a sizable
correction to the muon $g-2$.

##  III One- and two-loop muon $g-2$

The electromagnetic interaction of a lepton can be written as:

$\overline{\ell}(p^{\prime})\Gamma^{\mu}\ell(p)=\overline{\ell}(p^{\prime})\left[\gamma^{\mu}F_{1}(k^{2})+\frac{i\sigma^{\mu\nu}k_{\nu}}{2m_{\ell}}F_{2}(k^{2})\right]\ell(p)\,.$
(25)

The lepton anomalous magnetic dipole moment is then defined by

$a_{\ell}=\frac{g_{\ell}-2}{2}=F_{2}(0)\,.$ (26)

Since the magnetic moment is associated with dipole operator, the lepton $g-2$
originates from radiative quantum corrections. In the model, the one-loop
corrections from new physics are induced by the mediation of $H$, $A^{0}$, and
$H^{\pm}$, where the associated Feynman diagrams are shown in Fig. 1(a) and
(b). Moreover, it is known that the two-loop Barr-Zee type diagrams can have
important contributions to the magnetic dipole moment due to a large coupling
enhancement Barr:1990vd ; Chang:2000ii . The potentially large two-loop
diagrams mediated by heavy fermions, including top, bottom, and $\tau$, are
shown in Fig. 1(c). The essential mechanism contributing to the muon $g-2$ in
the model is the two-loop with Barr-Zee type diagram mediated by the charged
scalars, including $\delta^{++}$, $\delta^{+}$, and $H^{+}$, as shown in Fig.
1(d). In addition to the lepton Yukawa coupling, such diagrams further enjoy
the enhancement of the electric charges associated with the charged scalars.

Figure 1: One-loop and two-loop Barr-Zee type Feynman diagrams for the muon
$g-2$, where $f$ in plot (c) includes top (bottom) quark and $\tau$ lepton.

The one-loop corrections to the anomalous magnetic dipole moment in the 2HDM
have been studied long time ago Leveille:1977rc ; Haber:1978jt ;
Krawczyk:1996sm ; Dedes:2001nx . Using the Yukawa couplings shown in Eq. (21),
the muon $g-2$ from Fig. 1(a) and (b) can be expressed as:

$\displaystyle\Delta a^{1,H/A^{0}}_{\mu}$
$\displaystyle=\frac{m^{2}_{\mu}}{8\pi^{2}v^{2}}\left[(y^{H}_{\ell})^{2}\int^{1}_{0}\frac{r^{\mu}_{H}x^{2}(2-x)}{1-x(1-r^{\mu}_{H}x)}-(y^{A^{0}}_{\ell})^{2}\int^{1}_{0}\frac{r^{\mu}_{A^{0}}x^{3}}{1-x(1-r^{\mu}_{A^{0}}x)}\right]\,,$
$\displaystyle\Delta a^{1,H^{\pm}}_{\mu}$
$\displaystyle=\frac{m^{2}_{\mu}\zeta^{2}_{\ell}}{8\pi^{2}v^{2}}\int^{1}_{0}\frac{r^{\mu}_{H^{\pm}}x^{2}(1-x)}{1-x(1-r^{\mu}_{H^{\pm}}x)}\,,$
(27)

where $r^{f}_{B}=m^{2}_{f}/m^{2}_{B}$ with $B=H,~{}A^{0},~{}H^{\pm}$ and the
Yukawa couplings $y^{H,A^{0}}_{f}$ are defined as:

$y^{H}_{f}=s_{\beta-\alpha}\zeta_{f}-c_{\beta-\alpha}\,,~{}~{}y^{A^{0}}_{f}=-s_{f}\zeta_{f}\,.$
(28)

From the expressions, it can be seen that the induced muon $g-2$ is
proportional to $m^{4}_{\mu}$. One of the four factors of $m_{\mu}$ comes from
the definition in Eq. (25), another comes from the mass insertion for
chirality flip, and the rest two enter through the two Yukawa interaction
vertices, each of which is proportional to the muon mass. Thus, if the
intermediate scalar mass is of ${\cal O}(100)$ GeV, the resulting muon $g-2$
is far below $10^{-9}$. To get $\Delta a_{\mu}$ up to $\sim 10^{-9}$, the
mediating particle has to be as light as tens of GeV. This was the observation
previously found for the 2HDM in the literature.

Following the results shown in Ref. Ilisie:2015tra , the two-loop Barr-Zee
type diagrams with fermion and charged scalars can be written as:

$\displaystyle\Delta a^{2,f}_{\mu}$
$\displaystyle=\sum_{f=t,b,\tau}\frac{\alpha_{\rm
em}m^{2}_{\mu}N^{f}_{C}Q^{2}_{f}}{4\pi^{3}v^{2}}\left[y^{H}_{f}y^{H}_{\ell}J^{H}_{f}\left(r^{f}_{H}\right)+y^{A^{0}}_{f}y^{A^{0}}_{\ell}J^{A^{0}}_{f}\left(r^{f}_{A^{0}}\right)\right]\,,$
$\displaystyle\Delta a^{2,S}_{\mu}$
$\displaystyle=\sum_{S=\delta^{++},\delta^{+},H^{+}}\frac{\alpha_{\rm
em}Q^{2}_{S}r^{\mu}_{H}}{8\pi^{3}}\zeta_{\ell}\lambda_{HS^{*}S}\,J_{S}(r^{S}_{H})\,,$
(29)

where $N^{f}_{C}$ is the number of color for the fermion $f$, $Q_{P}$
($P=f,~{}S$) is the electric charge of the particle, and the loop functions
are given by:

$\displaystyle J^{H}_{f}(z)$
$\displaystyle=\frac{z}{2}\int^{1}_{0}\frac{2x(1-x)-1}{z-x(1-x)}\ln\left(\frac{z}{x(1-x)}\right)\,,$
$\displaystyle J^{A^{0}}_{f}(z)$
$\displaystyle=\frac{z}{2}\int^{1}_{0}\frac{1}{z-x(1-x)}\ln\left(\frac{z}{x(1-x)}\right)\,,$
$\displaystyle J_{S}(z)$
$\displaystyle=\frac{1}{2}\int^{1}_{0}\frac{x(1-x)}{z-x(1-x)}\ln\left(\frac{z}{x(1-x)}\right)\,.$
(30)

The two-loop results are proportional to $m^{2}_{\mu}$ because there is only
one muon Yukawa coupling involved. It can be seen that when $y^{H/A}_{u,d}$
are strictly bounded by the experimental data, their contributions become
subleading, and $\Delta a^{2,S}_{\mu}$ is the dominant effect.

##  IV Numerical Analysis

In this section, we present how we choose the parameters in our model, how
they affect $\Delta a_{\mu}$ at one-loop and partial two-loop levels, and how
the doubly-charged Higgs phenomenology at the LHC is modified.

###  IV.1 Parameter choice

Among the parameters in the Yukawa and scalar sectors, most relevant ones for
the muon $g-2$ are combinations of the Yukawa matrix elements, the quartic
scalar couplings, and $t_{\beta}$ that appear in various couplings. More
explicitly, the relevant parameters are:
$c_{\beta-\alpha}~{}(s_{\beta-\alpha})$, $\zeta_{u,d,\ell}$,
$m_{H,A^{0},H^{\pm}}$, $m_{\delta^{\pm\pm},\delta^{\pm}}$,
$\lambda_{H\delta^{--}\delta^{--}}$, $\lambda_{H\delta^{-}\delta^{-}}$, and
$\lambda_{HH^{-}H^{+}}$. We will show how they contribute the muon $g-2$.

Before a numerical analysis, we first need to find the allowed parameter space
for the model. All potential constraints from experimental measurements,
including various flavor physics processes, Higgs data, and electroweak
precision observables, and theoretical bounds, such as perturbative unitarity
and positivity of the scalar potential, have to be taken into account.
Recently, such a global fit, considering the theoretical constraints has been
done in the A2HDS Eberhardt:2020dat . In this work, we will follow the global
fit results in Ref. Eberhardt:2020dat  when the parameter values are taken for
the numerical estimations.

Two scenarios, the light scenario and the heavy scenario, are discussed by in
Ref. Eberhardt:2020dat , where the former refers to the case with
$m_{H}>m_{h}$ and the latter has $m_{H}<m_{h}$. Since we are interested in the
heavy scalar boson contribution to $\Delta a_{\mu}$, we will concentrate on
the light scenario.

The values of parameters used in our numerical analysis are described below.
Using the experimental data at $2\sigma$ errors, the global fit gives
$|c_{\beta-\alpha}|<0.04$. Thus, we will take the alignment limit of
$c_{\beta-\alpha}=0$. Under this limit, the $HW^{-}W^{+}$ and $HZZ$ couplings
vanish identically. It is found that $\zeta_{u,d}$ have to be of the same sign
and their values are restricted to small-value regions when $|\zeta_{\ell}|$
approaches the boundary of maximum, i.e., $|\zeta_{\ell}|=100$. Moreover, the
sign of $\zeta_{\ell}$ cannot be determined by the global fit, and it always
appears in the product along with other parameters in $\Delta a_{\mu}$, e.g.,
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$. We can thus fix the value of
$\zeta_{\ell}$ and let the associated parameter vary. In numerical
calculations, we take:

$\displaystyle\zeta_{u}=0.1,~{}\zeta_{d}=10,~{}\zeta_{\ell}=-100\,.$ (31)

Since the maximally allowed value of $|\zeta_{\ell}|$ in the negative region
is larger than that in positive region, we assume $\zeta_{\ell}$ to be
negative. The signs of $\zeta_{u,d}$ are taken to fit the positive $\Delta
a_{\mu}$.

From Eq. (16), it is seen that the trilinear couplings have involved relations
with the parameters in the scalar potential. If we assume that the 2HDM with a
$Z_{2}$ symmetry contributes little to the $HH^{-}H^{+}$ coupling, the
$\lambda_{6,7}$ parameters in A2HDS become the dominant source, i.e.,

$\lambda_{HH^{-}H^{+}}\simeq\lambda_{6}\left(3c^{2}_{\beta}s^{2}_{\beta}-s^{4}_{\beta}\right)-\lambda_{7}\left(c^{4}_{\beta}-3c^{2}_{\beta}s^{2}_{\beta}\right)\,.$
(32)

Since the constrained $\lambda_{6,7}$ values allow $\lambda_{6}\simeq-3.5$ and
$\lambda_{7}\simeq-2.5$, we will take $\lambda_{HH^{-}H^{+}}\approx 1.5$ to
estimate the muon $g-2$.

According to the results shown in Ref. Ouazghour:2018mld , the allowed values
for $|\lambda_{8,9}|$ and $|\lambda^{\prime}_{8,9}|$ can be of ${\cal
O}(10)$222The parameters $\lambda_{8,9}$ and $\lambda^{\prime}_{8,9}$ used in
this paper correspond respectively to the parameters $\lambda_{6,7,8,9}$ in
Ref. Ouazghour:2018mld . , where the Higgs data and the theoretical
constraints have been imposed. With $c_{\beta}=s_{\beta}=1/\sqrt{2}$, it is
expected that $\lambda_{H\delta^{--(-)}\delta^{++(+)}}\lesssim 5$ can be
conservative upper bounds. The upper bound is consistent with that used in
Ref. Ilisie:2014hea  for the $HH^{-}H^{+}$ coupling due to the perturbative
requirement. For simplicity, we take
$\lambda_{H\delta^{-}\delta^{+}}=\lambda_{H\delta^{--}\delta^{++}}$. Since
$\zeta_{\ell}$ and $\lambda_{H\delta^{--}\delta^{++}}$ show up together in
$\Delta a_{\mu}$, we can consider
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$ as a single variable because
$\zeta_{\ell}=-100$ is fixed.

Global fits cannot determine the masses of the involved new scalar particles.
Nevertheless, their mass differences are strongly correlated and constrained.
In our numerical analysis, we take $m_{\delta^{\pm}}=m_{\delta^{\pm\pm}}+100$
GeV according to Eq. (14). When discussing the CP-even or CP-odd scalar
effects, we take $m_{H(A^{0})}$ as a free parameter. When combining the
effects of $H$ and $A^{0}$ together, we take $m_{H}=m_{H^{\pm}}$ and
$m_{A^{0}}=m_{H}\pm 50$ GeV.

The parameters $|\zeta_{\ell}|$ and $m_{H^{\pm}}$ can also be bounded by the
Michel parameters Michel:1949qe  in the $\tau$ decays Abe:2015oca ;
Crivellin:2015hha ; Chun:2015hsa ; Enomoto:2015wbn ; Tobe:2016qhz . The Michel
parameters in the leptonic $\tau$ decay are defined as:

$\displaystyle\frac{d\Gamma_{\tau\to\ell\bar{\nu}_{\ell}\nu_{\tau}}}{dx}$
$\displaystyle\propto
f_{0}(x)+\rho_{\ell}f_{1}(x)+\eta_{\ell}\frac{m_{\ell}}{m_{\tau}}f_{2}(x)-P_{\tau}\xi_{\ell}\left(g_{1}(x)+\delta_{\ell}g_{2}(x)\right)\,,$
(33)

where $x=E_{\ell}/E_{\ell_{\rm max}}$ and $E_{\ell_{\rm
max}}=m_{\tau}(1+m_{\ell}^{2}/m_{\tau}^{2})/2$, $P_{\tau}$ is the
$\tau$-lepton polarization, and the explicit expressions of functions $f_{i}$
can be found in PDG . In the SM, the Michel parameters are predicted as
$\rho_{\ell}=3/4$, $\eta_{\ell}=0$, $\xi_{\ell}=1$, and $\delta_{\ell}=3/4$,
whereas the current experimental values are: $\rho^{\rm exp}_{\ell}=0.749\pm
0.008$, $\eta^{\rm exp}_{\ell}=0.015\pm 0.021$, $\xi^{\rm exp}_{\ell}=0.981\pm
0.031$, and $(\delta_{\ell}\xi_{\ell})^{\rm exp}=0.79\pm 0.04$. We will
concentrate on the $\eta_{\ell}$ and $\xi_{\ell}$ parameters because they are
sensitive to the scalar couplings.

The general transition matrix element for
$\tau\to\ell\bar{\nu}_{\ell}\nu_{\tau}$ can be written as Fetscher:1986uj :

${\cal M}=4\frac{G_{F}}{\sqrt{2}}\sum_{\begin{subarray}{c}\kappa=S,V,T\\\
\epsilon,\lambda=R,T\end{subarray}}g^{\kappa}_{\epsilon\lambda}\langle\bar{\ell}_{\epsilon}|\Gamma^{\kappa}|(\nu_{\ell})_{n}\rangle\langle(\bar{\nu}_{\tau})_{m}|\Gamma_{\kappa}|\tau_{\lambda}\rangle\,,$
(34)

where $\kappa=S,V,T$ denotes the type of interaction, $\epsilon(\lambda)=R,L$
is the lepton chirality, and the chirality of $m(n)$ can be determined when
$\kappa$ and $\epsilon(\lambda)$ are fixed. In the SM, due to the $V-A$
interaction, we only have $g^{V}_{LL}=-1$. Since the involved couplings in the
A2HDM are scalar and vector types and the $H^{\pm}$-Yukawa coupling is
proportional to the lepton mass, we only need to consider the muon mode and
the effective couplings $g^{S}_{RR}$ and $g^{V}_{LL}$. Thus, the Michel
parameters of $\eta_{\mu}$ and $\xi_{\mu}$ are expressed as:

$\displaystyle\eta_{\mu}$ $\displaystyle=\frac{1}{2}{\rm
Re}(g^{S}_{RR}\,g^{V*}_{LL})\,,$ $\displaystyle\xi_{\mu}$
$\displaystyle=|g^{V}_{LL}|^{2}-\frac{1}{4}|g^{S}_{RR}|^{2}\,,$ (35)

where $g^{S}_{RR}$ and $g^{V}_{LL}$ with one-loop corrections
Crivellin:2015hha  in the model are given by:

$\displaystyle g^{S}_{RR}$
$\displaystyle=\frac{m_{\mu}}{m_{\tau}}\left(\frac{m_{\tau}\zeta_{\ell}}{m_{H^{\pm}}}\right)^{2}\,,$
$\displaystyle g^{V}_{LL}$
$\displaystyle=-1-\frac{\zeta_{\ell}^{2}}{32\pi^{2}}\frac{m^{2}_{\tau}}{v^{2}}\left(F\left(\frac{m^{2}_{A^{0}}}{m^{2}_{H^{\pm}}}\right)+F\left(\frac{m^{2}_{H}}{m^{2}_{H^{\pm}}}\right)\right)\,,$
(36)

with the loop function $F$ defined by

$F(a)=\frac{1}{2}+\frac{1+a}{4(1-a)}\ln a\,.$ (37)

To illustrate the constraints from the measured Michel parameters, we show the
contours of $\eta_{\mu}$ (left panel) and $\xi_{\mu}$ (right panel) in the
$|\zeta_{\ell}|$-$m_{H^{\pm}}$ plane in Fig. 2, where $m_{A^{0}}=50$ GeV and
$m_{H}=100$ GeV are used. The dashed line in the left plot corresponds to the
$2\sigma$ lower bound of $\eta^{\rm exp}_{\mu}$, pointing to the lower right
region as more favorable parameter space. For example, $m_{H^{\pm}}\lesssim
180$ GeV is excluded when $|\zeta_{\ell}|=100$. The right plot, on the other
hand, does not show much constraining powers as the entire parameter space
gives values consistent with $\xi^{\rm exp}_{\mu}$ at the $2\sigma$ level.
Therefore, to avoid the constraint from $\eta^{\rm exp}_{\mu}$, we take
$m_{H^{\pm}}\gtrsim 180$ GeV in the following numerical analysis.

Figure 2:  Contours for the Michel parameters $\eta_{\mu}$ (left panel) and
$\xi_{\mu}$ (right panel) as a function of $|\zeta_{\ell}|$ and $m_{H^{\pm}}$,
where $m_{A^{0}}=50$ GeV and $m_{H}=100$ GeV are used.

###  IV.2  Muon $g-2$

In the following, we divide our discussion of muon $g-2$ into three
contributing parts.

####  IV.2.1 One-loop contribution

According to Eq. (III), $\Delta a^{1,H/A^{0}}_{\mu}$ and $\Delta
a^{1,H^{\pm}}_{\mu}$ with $c_{\beta-\alpha}=0$ depend only on the scalar boson
masses and the parameter $\zeta_{\ell}$. Thus, we show $\Delta a^{1,X}_{\mu}$
as a function of $m_{X}$ for $X=H,A^{0},H^{\pm}$ in Fig. 3(a), (b), and (c),
respectively, where the curves correspond to $|\zeta_{\ell}|=(30,50,100)$.
Clearly, the pseudoscalar $A^{0}$ and the charged-Higgs $H^{\pm}$
contributions are always negative at the one-loop level. The $H^{\pm}$
contribution is small and can be neglected. Although $H$ makes a positive
contribution to $\Delta a_{\mu}$, $\Delta a^{1,H}_{\mu}>10^{-9}$ is possible
only when $m_{H}$ is lighter than about 50 GeV. According to the global fit
results in the A2HDS with $m_{H}<m_{h}$ scenario Eberhardt:2020dat , such
light CP-even scalar is still allowed. However, the associated $m_{A^{0}}$
parameter can be of $O(100)$ GeV. To demonstrate the correlation between
$m_{H}$ and $m_{A^{0}}$, the contours of the combined $\Delta a_{\mu}$ as a
function of $m_{H}$ and $m_{A^{0}}$ are shown in Fig. 3(d), where
$\zeta_{\ell}=-100$ and $m_{H^{\pm}}=180$ GeV are used. When $A^{0}$ with the
allowed mass of $O(10^{2})$ GeV is included, comparing to the case without
$A^{0}$ contribution, the $m_{H}$ value required for $\Delta a_{\mu}>10^{-9}$
has to be shifted downward. Hence, in the region of $m_{H}>100$ GeV, the one-
loop contribution to the muon $g-2$ within the 2HDM is far below $10^{-9}$.

Figure 3:  One-loop muon $g-2$ within the 2HDM, induced by the mediation of
(a) $H$, (b) $A^{0}$, and (c) $H^{\pm}$, as a function of respectively
$m_{H,A^{0},H^{\pm}}$, where the curves correspond to
$|\zeta_{\ell}|=(30,50,100)$. (d) Contours of $\Delta a_{\mu}$ with the
combination of $H^{0}$, $A^{0}$, and $H^{\pm}$, where $\zeta_{\ell}=-100$ and
$m_{H^{\pm}}=180$ GeV are used.

####  IV.2.2 Barr-Zee contribution within the 2HDM

Before discussing the new two-loop effects on $\Delta a_{\mu}$ in our model,
we first focus on the analysis within the 2HDM. As discussed before, the main
two-loop Barr-Zee diagram contributions are usually from the fermion and
charged Higgs loops shown in Figs. 1(c) and 1(d), respectively. In order to
understand the influence of these effects on $\Delta a_{\mu}$, we separately
show the fermion and $H^{\pm}$ contributions in Fig. 4, where plot (a) [(b)]
is mediated by $H~{}[A^{0}]$; the dashed, dotted, dot-dashed, and dot-dot-
dashed curves are the top, bottom, $\tau$-lepton, and $H^{\pm}$ contributions,
respectively, and the solid curves are the combined fermion- and
$H^{\pm}$-loop Barr-Zee contributions and the one-loop results. In the
numerical estimates, we take $\zeta_{u}=0.1$, $\zeta_{b}=10$, and
$\zeta_{\ell}=-100$, consistent with the numerical results given in Ref.
Eberhardt:2020dat . To estimate the $H^{\pm}$ contribution, we use the
conservative value of $\lambda_{HH^{-}H^{+}}=2$ and set $m_{H^{\pm}}=m_{H}$.
In the $H$-mediated part given in Fig. 4(a), it can be seen that due to the
strict bounds on $\zeta_{u,d}$ the top- and bottom-quark contributions are
smaller than $2\times 10^{-10}$, although an enhancement factor of
$\zeta_{\ell}$ is already applied. The $H$-mediated $\tau$-loop contribution
is always negative and sizable in magnitude. Intriguingly, using
$\lambda_{HH^{-}H^{+}}=2$, it is found that the $H^{\pm}$-loop in the mass
region of $m_{H,H^{\pm}}>m_{h}$ gives the dominant effect, and it can overcome
the negative $\tau$-loop contribution, so that $\Delta a_{\mu}$ is positive.

Figure 4:  Two-loop muon $g-2$ from Figs. 1(c) and 1(d) induced by the
mediation of (a) $H$ and (b) $A^{0}$ as a function of respectively
$m_{H,A^{0}}$, where the dashed, dotted, dot-dashed, and dot-dot-dashed curves
are the top, bottom, $\tau$-lepton, and $H^{\pm}$ contributions, respectively.
The solid curves contain both one-loop and two-loop effects.

In the $A^{0}$-mediated part given in Fig. 4(b), similar to the case mediated
by $H$, the top-quark and bottom-quark effects are small. However, the bottom-
quark and $\tau$-lepton contributions interchange sign, and the latter becomes
the dominant effect. The sign difference arises from the loop functions
$J^{H}_{\tau}(z)$ and $J^{A^{0}}_{\tau}(z)$ shown in Eq. (III). For the region
of $m_{A^{0}}>m_{h}$, it can be seen from Fig. 4(b) that its two-loop effect
on $\Delta a_{\mu}$ is smaller than $10^{-9}$. Nevertheless, when $A^{0}$ is
lighter than the SM-like Higgs, its contribution to $\Delta a_{\mu}$ increases
significantly. When the negative one-loop contribution is included, one
observes that the $A^{0}$ contribution to $\Delta a_{\mu}$ can reach $2\times
10^{-9}$ when $m_{A^{0}}\approx 55$ GeV. Following the global analysis in the
A2HDS Eberhardt:2020dat , when $m_{H}>m_{h}$, such a light pseudoscalar boson
is not excluded by the current experimental data. Hence, the conclusion is
consistent with that obtained in the 2HDM type-X Broggio:2014mna ; Han:2015yys
.

####  IV.2.3 Barr-Zee contribution from the triplet extension

It has been shown that when the CP-even and CP-odd scalar masses are heavier
than the observed Higgs mass, the $H$-mediated and $A^{0}$-mediated effects in
the 2HDM become ineffective to accommodate the measured $\Delta a_{\mu}$. In
the following analysis, we discuss the new contributions from the doubly- and
singly-charged Higgs bosons derived from the Higgs triplet extension. In the
analysis, we focus on the $m_{H}>m_{h}$ scenario, following the parameter
constraints given in Ref. Eberhardt:2020dat .

According to Eq. (III), in addition to the mass factor $\Delta a^{2,Y}_{\mu}$
($Y=\delta^{++},\delta^{+}$) further depends on the product
$\zeta_{\ell}\lambda_{HY^{*}Y}$. Since the $\delta^{+}$-loop effect is similar
to the $H^{+}$-loop effect, its contribution is expected to be of ${\cal
O}(10^{-10})$. However, the doubled electric charge of $\delta^{++}$ results
in a factor of 4 enhancement. We show contours of $\Delta
a^{2,\delta^{++}+\delta^{+}}_{\mu}$ in the
$m_{H}$-$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$ plane in Fig. 5(a),
where $\lambda_{H\delta^{-}\delta^{+}}=\lambda_{H\delta^{--}\delta^{++}}$,
$m_{\delta^{\pm}}=m_{\delta^{\pm\pm}}+100$ GeV, and $m_{\delta^{\pm\pm}}=350$
GeV. Fixing $\zeta_{\ell}=-100$, we treat $\lambda_{HY^{*}Y}$ as a variable
and set $\lambda_{HY^{*}Y}\lesssim 5$ to satisfy the perturbativity bound
Ilisie:2014hea . We observe that $\Delta
a^{2,\delta^{++}+\delta^{+}}_{\mu}\approx 1.3\times 10^{-9}$ can be achieved
in the model when $m_{H}\approx 200$ GeV and
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}\approx-320$. Even when using
the maximal value of $|\zeta_{\ell}\lambda_{H\delta^{--}\lambda^{++}}|=500$,
$\Delta a^{2,\delta^{++}+\delta^{+}}_{\mu}$ can still reach $10^{-9}$ at
$m_{H}\approx 500$ GeV. In Fig. 5(b), we show contours of $\Delta
a^{2,\delta^{++}+\delta^{+}}_{\mu}$ in the
$m_{\delta^{\pm\pm}}$-$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$ plane,
where $m_{H}=200$ GeV is taken. Clearly, with
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}\approx-500$, we can have
$\Delta a^{2,\delta^{++}+\delta^{+}}_{\mu}\approx 10^{-9}$ at
$m_{\delta^{\pm\pm}}\approx 500$ GeV. These results demonstrate that the
measured muon $g-2$ can be readily achieved even when the exotic Higgs bosons
in our model have mass of a few $\times 100$ GeV.

Figure 5:  Contours of $\Delta a_{\mu}$ from the 2-loop Barr-Zee type diagrams
of $\delta^{\pm\pm}$ and $\delta^{\pm}$ in the plane of (a) $m_{H}$ and
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$ with $m_{\delta^{\pm\pm}}=350$
GeV and (b) $m_{\delta^{\pm\pm}}$ and
$\zeta_{\ell}\lambda_{H\delta^{--}\delta^{++}}$ with $m_{H}=200$ GeV.

Besides the enhancement from the two units of electric charge, it is
interesting to note the other enhancement factor associated with the doubly-
charged Higgs boson by comparing the result with that induced from the
$\tau$-loop. Because the $H$ and $A^{0}$ couplings to muon are the same in the
alignment limit, if we further set $m_{H}=m_{A^{0}}=m_{X}$ the only different
factors come from the couplings, $A^{0}\tau\tau$ and $HH^{-}H^{+}$, and from
the loop integrals, $J^{A^{0}}_{f}(z)$ and $J_{S}(z)$ defined in Eq. (III).
For simplicity, we use $z_{\chi}J_{\chi}$ to represent the effect for the
$\tau$-loop and the $\delta^{++}$-loop, where $z_{\chi}$ and $J_{\chi}$ are
the associated coupling factor and integral function, respectively. We show
$z_{\chi}J_{\chi}$ as a function $m_{X}$ in Fig. 6, where $\zeta_{\ell}=-100$,
$m_{\delta^{++}}=350$ GeV, and $\lambda_{H\delta^{--}\delta^{++}}=3$ are
applied. It can be seen that once $m_{X}>160$ GeV, the $\delta^{++}$-loop
contribution is larger than the $\tau$-loop.

Figure 6:  Product of coupling factor and loop-integral function for the
$\delta^{++}$-loop and $\tau$-loop.

So far, we have just analyzed the effects of the Higgs triplet field on the
muon $g-2$. Although the 2HDM contribution to $\Delta a_{\mu}$ becomes less
significant as $m_{A^{0}}>m_{h}$, its effects are part of the model and we
should combine them with the Higgs triplet effects altogether. According to
the global fit analysis presented in Ref. Eberhardt:2020dat , the masses of
$H$, $A^{0}$, and $H^{\pm}$ are strongly correlated and, in the $m_{H}>m_{h}$
scenario, $A^{0}$ and $H^{\pm}$ can be heavier or lighter than $h$. For the
purpose of illustration, we show the cases with $m_{H^{\pm}}=m_{H}$ for
$m_{A^{0}}=m_{H}+50$ GeV and for $m_{A^{0}}=m_{H}-50$ GeV in Fig. 7(a) and
Fig. 7(b), respectively. The dashed, dotted, and dot-dashed curves give
respectively the fermion-loop, $H^{\pm}$, and $\delta^{++,+}$ contributions.
The dot-dot-dashed curves are the two-loop 2HDM results, where the fermion-
and $H^{\pm}$-loop contributions are summed up. The solid curves combine all
the above-mentioned contributions, including the one-loop effects. To show the
maximal contribution from the Higgs triplet field, we take
$\lambda_{H\delta^{--}\delta^{++}}=5$ in the plots, and the other parameter
are taken to be the same as those used in the earlier plots. In order to see
the effects of $\zeta_{\ell}$ and $m_{\delta^{\pm\pm}}$ on $\Delta a_{\mu}$,
we show the contours of $\Delta a_{\mu}$ (in units of $10^{-9}$) in the plane
of $\zeta_{\ell}$ and $m_{\delta^{\pm\pm}}$ in Fig. 8, where the parameter
values are taken to be the same as those used in Fig. 7 with the exception of
$m_{H,H^{\pm}}=180$ GeV and $m_{A^{0}}=120$ GeV.

Figure 7:  Combined results of one-loop and two-loop contributions in the
model for (a) $m_{H^{\pm}}=m_{H}$ and $m_{A^{0}}=m_{H}+50$ GeV and for (b)
$m_{H^{\pm}}=m_{H}$ and $m_{A^{0}}=m_{H}-50$ GeV, where
$m_{\delta^{\pm\pm}}=350$ GeV is used and the other parameters are the same as
in plot (a).  Figure 8:  Contours for the combined one- and two-loop
contributions to $\Delta a_{\mu}$ (in units of $10^{-9})$ in the
$\zeta_{\ell}$-$m_{\delta^{\pm\pm}}$ plane.

Based on the analyses, we summarize the results as follows:

  1. 1.

The contribution of Barr-Zee type diagram with the $H^{\pm}$-loop can only be
of ${\cal O}(10^{-10})$.

  2. 2.

When the $H$\- and $A^{0}$-mediated two-loop diagrams are combined, $A^{0}$
with a mass of $70$ GeV can lead to $\Delta a_{\mu}\sim 2\times 10^{-9}$,
similar to the situation of 2HDM type-X Broggio:2014mna . However, when
$m_{A^{0}}>m_{h}$, the 2HDM contribution is below $10^{-9}$.

  3. 3.

The $H$-mediated $\delta^{++(+)}$-loop contribution is independent of
$m_{A^{0}}$, and it can play an important role on $\Delta a_{\mu}$ in a wide
range of $m_{H}$ and $m_{\delta^{\pm\pm}}$, particularly when the exotic Higgs
masses are a few hundred GeV.

Before closing this subsection, we make a brief remark about the implication
of the model on the anomalous magnetic dipole moment of the electron (electron
$g-2$). Applying the accurate measurements of the fine structure constant from
133Cs and 87Rb to the theoretical calculations Aoyama:2014sxa ; ATN2019 , the
differences in the electron $g-2$ between the experiments and the SM
expectation are found to be:

$\displaystyle\Delta a_{e}(^{133}\text{Cs})$ $\displaystyle=-(8.8\pm
3.6)\times 10^{-13}~{}~{}\text{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{PYZEM}{\@@citephrase{(}}{\@@citephrase{)}}}}\,,$
$\displaystyle\Delta a_{e}(^{87}\text{Rb})$ $\displaystyle=(4.8\pm 3.0)\times
10^{-13}~{}~{}\text{\cite[cite]{\@@bibref{Authors
Phrase1YearPhrase2}{MYCG}{\@@citephrase{(}}{\@@citephrase{)}}}}\,,$ (38)

i.e., having $-2.4\sigma$ and $1.6\sigma$ deviations, respectively. Their
weighted average is $\Delta a_{e}=-0.8\pm 2.3$. In spite of the inconclusively
experimental results, the model predicts a concrete correlation in the
corrections for muon $g-2$ and electron $g-2$. Since the dominant contribution
comes from the Barr-Zee type diagrams in the model and the lepton Yukawa
couplings to $H(A^{0})$ are proportional to $m_{\ell}$ according to Eq. (III),
we have the ratio $\Delta a_{e}/\Delta a_{\mu}=(m_{e}/m_{\mu})^{2}\simeq
2.36\times 10^{-5}$. Hence, $\Delta a_{e}$ and $\Delta a_{\mu}$ should have
the same sign in the model.

###  IV.3 Doubly-charged Higgs decays

In this subsection we discuss the decay branching ratios of $\delta^{\pm\pm}$
in our scenario and its related collider signature. With our choice of
$v_{\Delta}>\mathcal{O}(10^{-4})$ GeV, $\delta^{\pm\pm}$ dominantly decays
into charged scalar and/or $W$ bosons, evading the stringent constraint of LHC
searches for same-sign dileptons. The interactions relevant to the dominant
decays of $\delta^{\pm\pm}$ are obtained from the scalar potential and kinetic
term as follows

$\displaystyle\mathcal{L}\supset$
$\displaystyle-(\mu_{1}s_{\beta}^{2}+\mu_{2}c_{\beta}^{2}-\mu_{3}s_{\beta}c_{\beta})\delta^{++}H^{-}H^{-}-\frac{v}{2}\left(s_{2\beta}\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2}-c_{2\beta}\lambda^{\prime}_{12}\right)\delta^{++}H^{-}\delta^{-}$
$\displaystyle-\frac{g^{2}}{\sqrt{2}}v_{\Delta}\delta^{++}W^{-\mu}W^{-}_{\mu}-igW^{-\mu}(\partial_{\mu}\delta^{-}\delta^{++}-\delta^{-}\partial_{\mu}\delta^{++})+\mbox{h.c.}\
,$ (39)

where $g$ is the $SU(2)_{L}$ gauge coupling. It can be seen that the
$\delta^{++}H^{-}H^{-}$ and $\delta^{++}W^{-}W^{-}$ couplings are dictated by
the factor of $v_{\Delta}$. That is, for the doubly-charged-Higgs decays, the
mixing of ${\cal O}(v_{\Delta})$ among $G^{\pm}$, $H^{\pm}$, and
$\delta^{\pm}$ has to be taken into account though its effects are small in
the muon $g-2$. If we take $\mu_{3}=t_{2\beta}(\mu_{1}-\mu_{2})$ and
$\lambda_{8}c^{2}_{\beta}+\lambda_{9}s^{2}_{\beta}+\lambda^{\prime}_{12}s_{2\beta}\sim-2(2m_{\Delta}/v)^{2}$
in Eq. (54), it is found that $m^{2}_{G^{-}\delta^{+}}\sim{\cal
O}(v^{2}_{\Delta})$ and $m^{2}_{H^{-}\delta^{+}}\sim{\cal O}(v_{\Delta})$.
With such parameters and dropping ${\cal O}(v^{2}_{\Delta})$ terms, the
relevant charged scalar mass-square matrix can be simplified as a $2\times 2$
matrix:

$\begin{pmatrix}H^{-}\delta^{-}\end{pmatrix}\begin{pmatrix}m^{2}_{H^{\pm}}&m^{2}_{H^{-}\delta^{+}}\\\
m^{2}_{H^{-}\delta^{+}}&m^{2}_{\delta^{\pm}}\\\
\end{pmatrix}\begin{pmatrix}H^{+}\\\ \delta^{+}\\\ \end{pmatrix}\,,$ (40)

where the matrix can be diagonalized by an $SO(2)$ transformation, similar to
the expression given in Eq. (55) but using the $\theta_{\pm}$ mixing angle
instead of $\alpha$. Due to the fact that $v_{\Delta}\ll v$, we can ignore the
influence of $m^{2}_{H^{-}\delta^{+}}$ on $m^{2}_{H^{\pm}}$ and
$m^{2}_{\delta^{\pm}}$, and the physical states and the small mixing angle can
be simply expressed as:

$\displaystyle\begin{split}H^{\pm}_{1}&\simeq
H^{\pm}+\theta_{\pm}\,\delta^{\pm}~{},\\\
H^{\pm}_{2}&\simeq-\theta_{\pm}\,H^{\pm}+\delta^{\pm}~{},\\\
\theta_{\pm}&\simeq\frac{vv_{\Delta}}{2\sqrt{2}(m^{2}_{\delta^{\pm}}-m^{2}_{H^{\pm}})}\left(s_{2\beta}\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2}+s^{2}_{\beta}\lambda^{\prime}_{12}\right)~{}.\end{split}$
(41)

Taking $\lambda^{\prime}_{8}-\lambda^{\prime}_{9}=1$,
$\lambda^{\prime}_{12}=1$, $m_{\delta^{\pm}}=350$ GeV, $m_{H^{\pm}}=180$ GeV,
$v_{\Delta}=5\times 10^{-4}$ and $c_{\beta}=s_{\beta}=1/\sqrt{2}$, the mixing
angle value is estimated to be $\theta_{\pm}\simeq 4.8\times 10^{-7}$. Since
$H^{\pm}_{1(2)}$ only carries a tiny component of $\delta^{\pm}(H^{\pm})$, in
the following analysis we keep using $H^{\pm}(\delta^{\pm})$ instead of
$H^{\pm}_{1}(H^{\pm}_{2})$.

Because of the introduction of $\theta_{\pm}$, in addition to the
$W^{\pm}W^{\pm}$ and $H^{\pm}H^{\pm}$ modes, the doubly-charged Higgs can also
decay into $H^{\pm}W^{\pm}$. In the case an on-shell decay to $H^{\pm}H^{\pm}$
is impossible, we consider the decay channel $H^{\pm}H^{\pm*}$, i.e., one of
$H^{\pm}$ is off-shell, with the assumption that the $H^{\pm}\to\tau^{\pm}\nu$
is dominant. The partial decay widths of these dominant modes are explicitly
expressed as:

$\displaystyle\Gamma(\delta^{\pm\pm}\to W^{\pm}W^{\pm})\simeq$
$\displaystyle\frac{g^{4}v_{\Delta}^{2}m_{\delta^{\pm\pm}}^{3}}{16\pi
m_{W}^{4}}\left(3r_{W}^{2}-r_{W}+\frac{1}{4}\right)\sqrt{1-4r_{W}}~{},$ (42)
$\displaystyle\Gamma(\delta^{\pm\pm}\to H^{\pm}W^{\pm})\simeq$
$\displaystyle\frac{g^{2}m^{3}_{\delta^{\pm\pm}}}{16\pi
m_{W}^{2}}\theta_{\pm}^{2}\lambda(1,r_{W},r_{\pm})^{\frac{2}{3}}~{},$ (43)
$\displaystyle\Gamma(\delta^{\pm\pm}\to H^{\pm}H^{\pm})\simeq$
$\displaystyle\frac{v^{2}\sqrt{1-4r_{\pm}}}{4\pi
m_{\delta^{\pm\pm}}}|\lambda_{\delta^{++}H^{-}H^{-}}|^{2}~{},$ (44)
$\displaystyle\Gamma(\delta^{\pm\pm}\to H^{\pm}H^{\pm*})\simeq$
$\displaystyle\frac{1}{2^{6}\pi^{3}m_{\delta^{\pm\pm}}}\left|\sqrt{2}m_{\tau}\zeta_{\ell}\lambda_{\delta^{++}H^{-}H^{-}}\right|^{2}$
$\displaystyle\times\int^{x_{\rm max}}_{x_{\rm
min}}dx\frac{x}{(x-r_{\pm})^{2}}\sqrt{\lambda\left(1,r_{\pm},x\right)}~{},$
(45)

where $r_{W}=m^{2}_{W}/m^{2}_{\delta^{\pm\pm}}$,
$r_{\pm}=m_{H^{\pm}}^{2}/m_{\delta^{\pm\pm}}^{2}$;
$x=q^{2}/m^{2}_{\delta^{\pm\pm}}$ and $q^{2}$ is the invariant mass of
$\tau\nu_{\tau}$ in the $H^{+}\to\tau^{+}\nu_{\tau}$ decay; $x_{\rm
max}=(1-\sqrt{r_{\pm}})^{2}$, $x_{\rm
min}=m^{2}_{\tau}/m^{2}_{\delta^{\pm\pm}}$;
$\lambda(1,x,y)=1+x^{2}+y^{2}-2xy-2x-2y$, and

$\lambda_{\delta^{++}H^{-}H^{-}}=\frac{1}{v}\left(\mu_{1}s_{\beta}^{2}+\mu_{2}c_{\beta}^{2}-\mu_{3}s_{\beta}c_{\beta}\right)+\frac{1}{2}\left(s_{2\beta}\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2}-c_{2\beta}\lambda^{\prime}_{12}\right)\theta_{\pm}\,.$
(46)

For the purpose of comparison, we include the
$\delta^{++}\to\ell^{+}_{i}\ell^{+}_{j}$ channels. According to the
$\delta^{++}$ Yukawa couplings to $\ell^{-}_{i}\ell^{-}_{j}$ shown in Eq.
(24), the decay rate is given by:

$\Gamma(\delta^{++}\to\ell^{+}_{i}\ell^{+}_{j})=\frac{m_{\delta^{\pm\pm}}}{8\pi(1+\delta_{ij})}\left|\frac{(M_{\nu})_{ij}}{v_{\Delta}}\right|^{2}\,,$
(47)

where $\delta_{ij}$ is the Kronecker delta. We use the $(M_{\nu})_{ij}$ values
determined by a global fit to the neutrino data and given in deSalas:2017kay .

For illustration purposes, we show in Fig. 9 their branching ratios as
functions of $\mu_{1}$ with the parameter choice of $t_{\beta}=1$,
$|\zeta_{\ell}|=100$, $\epsilon_{\Delta}=v_{\Delta}=5\times 10^{-4}$ GeV,
$m_{\delta^{\pm\pm}}=300~{}(400)$ GeV and $m_{H^{\pm}}=180$ GeV for the left
(right) plot where Eq. (A3) is applied to fix $\mu_{2}$ and $\mu_{3}$. For the
neutrino mass matrix element in Eq. (47), we apply the dominant ones
$(M_{\nu})_{\mu\mu,\tau\tau,\mu\tau}\sim 2\times 10^{-2}$ eV according to the
global fit in deSalas:2017kay  for the normal ordering of neutrino masses.
Note that one of $H^{+}$ is off-shell in the left plot while both $H^{+}$ are
on-shell in the right plot. In the former case, we find that the
$H^{+}H^{+*}(\to\tau^{+}\nu)$ mode can be dominant when $\mu_{1}\gg
v_{\Delta}$, even though it is a three-body decay. In the latter case, we find
that the $H^{+}H^{+}$ mode is dominant for $\mu_{1}>\mathcal{O}(v_{\Delta})$
and that the ratio of branching ratio is ${\rm Br}(H^{+}H^{+}):{\rm
Br}(W^{+}W^{+})\simeq 0.07:1$ when $\mu_{1}\lesssim v_{\Delta}$ except for the
region $\mu_{1}\sim v_{\Delta}/2$. The suppression of ${\rm Br}(H^{+}H^{+})$
at around $\mu_{1}=v_{\Delta}/2$ is due to a cancellation in the coupling. The
cancellation occurs when Eq. (A3) is applied. With $t_{\beta}=1$, we have
$\mu_{1}=\mu_{2}$ and $\mu_{3}=2(\epsilon_{\Delta}-\mu_{1})$,
$\mu_{1}s^{2}_{\beta}+\mu_{2}c^{2}_{\beta}-\mu_{3}c_{\beta}s_{\beta}\simeq
2\mu_{1}-\epsilon_{\Delta}$. The mixing angle $\theta_{\pm}$ could induce a
sizable effect if $\lambda^{\prime}_{8}-\lambda^{\prime}_{9}=2$ is used. Note
also that $\lambda^{\prime}_{12}$ does not contribute to the coupling
$\lambda_{\delta^{++H^{-}H^{-}}}$ when $t_{\beta}=1$ is taken. In Fig. 10, we
also show the branching ratios of $\delta^{++}$ as functions of
$m_{\delta^{++}}$ with $\mu_{1}=0.01$ GeV, where the other parameter values
are the same as those shown in Fig. 9. We find that the $H^{+}H^{+(*)}$ mode
becomes dominant when $m_{\delta^{\pm\pm}}>2m_{H^{+}}$ and that the branching
ratios of $W^{+}W^{+}$ and $H^{+}W^{+}$ slightly increase with
$m_{\delta^{\pm\pm}}$.

Figure 9:  Branching ratios of $\delta^{++}$ as functions of $\mu_{1}$, with
$t_{\beta}=1$, $\epsilon_{\Delta}=v_{\Delta}=5\times 10^{-4}$ GeV,
$m_{\delta^{\pm\pm}}=300~{}(400)$ GeV and $m_{H^{\pm}}=180$ GeV for left
(right) plots.  Figure 10:  Branching ratios of $\delta^{++}$ as functions of
$m_{\delta^{\pm\pm}}$, with $t_{\beta}=1$,
$\epsilon_{\Delta}=v_{\Delta}=5\times 10^{-4}$ GeV, $\mu_{1}=0.01$ GeV and
$m_{H^{\pm}}=180$ GeV.

Finally we discuss signals from doubly-charged Higgs boson production at the
LHC. The $\delta^{++}\delta^{--}$ pair can be produced via electroweak
interactions in proton-proton collision process. In our scenario, the produced
$\delta^{\pm\pm}$ dominantly decay into $W^{\pm}W^{\pm}$ and/or
$H^{\pm}H^{\pm(*)}$ mode, depending on the value of $\mu_{1}$ as discussed
above. Here we focus on the $H^{\pm}H^{\pm(*)}$ mode as it is a special
channel in our model. The scenario when the $W^{\pm}W^{\pm}$ mode is more
dominant is equivalent to the Higgs triplet model with
$v_{\Delta}>\mathcal{O}(10^{-4})$ GeV. The singly-charged Higgs $H^{\pm}$
dominantly decays into $\tau\nu$ mode since the Yukawa interaction between
$H^{\pm}$ and leptons is enhanced by the large $\zeta_{\ell}$ factor that is
required for a sizable muon $g-2$ contribution. Therefore, the signature from
$\delta^{++}\delta^{--}$ production is $4\tau+\not{E}_{T}$ in our scenario.
Note that we can relax the bound on $m_{\delta^{\pm\pm}}$ of $\sim 350$ GeV
Aad:2021lzu  because the analysis assumes that the $W^{\pm}W^{\pm}$ mode is
dominant and considers only the muons/electrons in the final state. It would
be difficult to reconstruct the doubly-charged Higgs mass because of the
missing transverse energy carried away by neutrinos from $H^{\pm}$ decays. Our
signature could be tested in multi-tau searches in future LHC experiments.

##  V Summary

In this work, we have studied an extension of the Standard Model only in the
scalar sector, with the addition of one Higgs doublet and one complex Higgs
triplet, rending a two-Higgs-doublet model (2HDM) with the type-II seesaw
mechanism. For the 2HDM part, we consider the aligned two-Higgs-doublet scheme
(A2HDS) to avoid undesired flavor-changing neutral currents induced by the two
Higgs doublet fields and to satisfy the current Higgs data constraints. The
Higgs triplet field obtains a small vacuum expectation value (VEV) induced by
the electroweak symmetry breaking and gives Majorana mass to neutrinos through
Yukawa couplings.

We have examined how the model can accommodate the measured muon $g-2$
deviation. Simple 2HDMs usually require CP-even and -odd Higgs bosons ($H$ and
$A^{0}$) to be sufficiently light (about a few $\times{\cal O}(10)$ GeV) and
rely on the contributions of Barr-Zee type diagrams to account for the muon
$g-2$ anomaly, $\Delta a_{\mu}$. In our model, the Barr-Zee type diagrams get
additional contributions from a large $H\delta^{++(+)}\delta^{--(-)}$
coupling, an enhanced coupling between charged leptons and the charged Higgs
boson ($H^{\pm}$) in 2HDM, and the electric charges of the charged Higgs
bosons ($\delta^{\pm}$ and $\delta^{\pm\pm}$) from the Higgs triplet field,
independent of the mass of CP-odd Higgs boson. In fact, the mass of the exotic
Higgs bosons is allowed to have a wider range, up to a few hundred GeV.

Owing to the new interactions with the other charged Higgs and $W$ bosons, the
doubly-charged Higgs boson presents a different decay pattern than the usual
Higgs triplet model. With the assumed Higgs triplet VEV, $v_{\Delta}\sim
5\times 10^{-4}$ GeV, the doubly-charged Higgs boson may dominantly decay into
like-sign charged Higgs bosons in the 2HDM rather than like-sign $W$ bosons,
when the magnitudes of the trilinear couplings $\mu_{1,2,3}$ are greater than
$10^{-3}~{}(10^{-2})$ GeV for both (one of the) charged Higgs bosons being on-
shell. Therefore, pair productions of the doubly-charged Higgs bosons will
lead to the signature of 4 $\tau$-leptons and missing energy at the LHC.

## Acknowledgments

This work is supported in part by the Ministry of Science and Technology
(MOST) of Taiwan under Grant Nos. MOST-108-2112-M-006-003-MY2 (CHC) and
MOST-108-2112-M-002-005-MY3 (CWC).

##  Appendix A Mass matrices for neutral and charged Higgs bosons

In this appendix, we show the full mass matrices at tree level for scalar,
pseudoscalar, and charged Higgs bosons. Since $A^{0}$ and $H^{\pm}$ are the
physical states in the 2HDM, it is useful to show the scalar mass matrices in
terms of the $(G^{0},A^{0},\eta^{0})$ and $(G^{+},H^{+},\delta^{+})$ bases
when the Higgs triplet field $\Delta$ is introduced. Since the CP-even
$H^{0}_{1}$ and $H^{0}_{2}$ scalars mix, we show the mass matrix in the basis
of the $(h,H)$ states, defined in Eq. (6).

From the scalar potentials given in Eqs. (8)$-$(10) and the Higgs basis in Eq.
(4), the mass matrix for the CP-odd components $G^{0}$, $A^{0}$ and $\eta^{0}$
is given by

$\frac{1}{2}\begin{pmatrix}G^{0}~{}A^{0}~{}\eta^{0}\end{pmatrix}\begin{pmatrix}m^{2}_{G^{0}G^{0}}&m^{2}_{G^{0}A^{0}}&m_{G^{0}\eta^{0}}^{2}\\\
m^{2}_{G^{0}A^{0}}&m_{A^{0}A^{0}}^{2}&m_{A^{0}\eta^{0}}^{2}\\\
m_{G^{0}\eta^{0}}^{2}&m_{A^{0}\eta^{0}}^{2}&m_{\eta^{0}\eta^{0}}^{2}\end{pmatrix}\begin{pmatrix}G^{0}\\\
A^{0}\\\ \eta^{0}\end{pmatrix},$ (48)

where the mass matrix elements

$\displaystyle m^{2}_{G^{0}G^{0}}=$
$\displaystyle\frac{3v_{\Delta}}{\sqrt{2}}\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\right)\,,$
$\displaystyle m^{2}_{G^{0}A^{0}}=$
$\displaystyle-\sqrt{2}v_{\Delta}\left[s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3}\right]\,,$
$\displaystyle m^{2}_{G^{0}\eta^{0}}=$
$\displaystyle-\sqrt{2}v\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\right)\,,$
$\displaystyle m_{A^{0}A^{0}}^{2}\equiv$ $\displaystyle
m^{2}_{A^{0}}=\frac{m^{2}_{12}}{s_{\beta}c_{\beta}}-\lambda_{5}v^{2}-\frac{v^{2}}{2}\left(\lambda_{6}t^{-1}_{\beta}+\lambda_{7}t_{\beta}\right)+2\sqrt{2}s^{2}_{\beta}v_{\Delta}\mu_{1}$
$\displaystyle+\frac{c^{2}_{\beta}v_{\Delta}\mu_{3}}{\sqrt{2}t_{\beta}}(1-t^{2}_{\beta})^{2}-\frac{v^{2}_{\Delta}\bar{\lambda}_{12}}{2}(s_{2\beta}+c^{3}_{\beta}+s^{3}_{\beta})\,,$
$\displaystyle m_{A^{0}\eta^{0}}^{2}=$
$\displaystyle\frac{v}{\sqrt{2}}\left[s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3}\right]\,,$
$\displaystyle m_{\eta^{0}\eta^{0}}^{2}\equiv$ $\displaystyle
m^{2}_{\eta^{0}}=m^{2}_{\Delta}+\left(\lambda_{\Delta 1}+\lambda_{\Delta
2}\right)v^{2}_{\Delta}+\frac{\bar{\lambda}_{8}c^{2}_{b}+\bar{\lambda}_{9}s^{2}_{\beta}}{2}v^{2}+\frac{\bar{\lambda}_{12}s_{2\beta}}{2}v^{2}\,,$
(49)

with $\bar{\lambda}_{8}=\lambda_{8}+\lambda^{\prime}_{8}$ and
$\bar{\lambda}_{9}=\lambda_{9}+\lambda^{\prime}_{9}$. From Eq. (12), it is
known that
$c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\sim{\cal
O}(v_{\Delta})$. Neglecting term of ${\cal O}(v_{\Delta})$, one can see that
$m^{2}_{G^{0}G^{0}}$, $m^{2}_{G^{0}A^{0}}$, and $m^{2}_{G^{0}\eta^{0}}$ are
negligibly small. Thus, to a good approximation, $G^{0}$ represents the
neutral Goldstone boson. Moreover, if we further demand that the factor
$s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3}$ in $m^{2}_{A^{0}\eta^{0}}$
vanish, $A^{0}$ and $\eta^{0}$ decouple from each other. To understand the
correlations among $\mu_{1}$, $\mu_{2}$, and $\mu_{3}$ under the conditions:

$\displaystyle\begin{cases}c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}&=\epsilon_{\Delta}\,,\\\
s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3}&=0\,,\end{cases}$ (50)

where $\epsilon_{\Delta}$ is an parameter of ${\cal O}(v_{\Delta})$, we solve
and obtain:

$\mu_{2}=\mu_{1}t^{-2}_{\beta}+\epsilon_{\Delta}\left(1-t^{-2}_{\beta}\right)\,,~{}~{}\mu_{3}=-2t^{-1}_{\beta}(\mu_{1}-\epsilon_{\Delta})\,.$
(51)

For $\tan\beta=1$, we have $\mu_{2}\simeq\mu_{1}$ and
$\mu_{3}\simeq-2\mu_{1}$. For large $\tan\beta$, they can be approximated as:

$\mu_{2}\simeq\epsilon_{\Delta}+\mu_{1}t^{-2}_{\beta}\,,~{}~{}\mu_{3}\simeq-2\mu_{1}t^{-1}_{\beta}\,.$
(52)

It is seen that the $\mu_{2,3}$ scale is determined by $\mu_{1}$ and
$t_{\beta}$.

Similarly, the mass matrix for $G^{+}$, $H^{+}$ and $\delta^{+}$ is given by:

$(G^{-}H^{-}\delta^{-})\left(\begin{array}[]{ccc}m^{2}_{G^{-}G^{+}}&m^{2}_{G^{-}H^{+}}&m^{2}_{G^{-}\delta^{+}}\\\
m^{2}_{G^{-}H^{+}}&m^{2}_{H^{-}H^{+}}&m^{2}_{H^{-}\delta^{+}}\\\
m^{2}_{G^{-}\delta^{+}}&m^{2}_{H^{-}\delta^{+}}&m^{2}_{\delta^{-}\delta^{+}}\\\
\end{array}\right)\left(\begin{array}[]{c}G^{+}\\\ H^{+}\\\ \delta^{+}\\\
\end{array}\right)\,,$ (53)

where the mass matrix elements

$\displaystyle m^{2}_{G^{-}G^{+}}$
$\displaystyle=\sqrt{2}v_{\Delta}\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\right)-\frac{v^{2}_{\Delta}}{2}(c^{2}_{\beta}\lambda^{\prime}_{8}+s^{2}_{\beta}\lambda^{\prime}_{9}+\lambda^{\prime}_{12}s_{2\beta})\,,$
$\displaystyle m^{2}_{G^{-}H^{+}}$
$\displaystyle=-\frac{v_{\Delta}}{\sqrt{2}}(s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3})+\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2}c_{\beta}s_{\beta}v^{2}_{\Delta}-\frac{\lambda^{\prime}_{12}}{2}c_{2\beta}v^{2}_{\Delta}\,,$
$\displaystyle m^{2}_{G^{-}\delta^{+}}$
$\displaystyle=-v\left(c^{2}_{\beta}\mu_{1}+s^{2}_{\beta}\mu_{2}+c_{\beta}s_{\beta}\mu_{3}\right)+\frac{v_{\Delta}v}{2\sqrt{2}}\left(c^{2}_{\beta}\lambda^{\prime}_{8}+s^{2}_{\beta}\lambda^{\prime}_{9}+s_{2\beta}\lambda^{\prime}_{12}\right)\,,$
$\displaystyle m^{2}_{H^{-}H^{+}}$ $\displaystyle\equiv
m^{2}_{H^{\pm}}=\frac{m^{2}_{12}}{c_{\beta}s_{\beta}}-\frac{\lambda_{4}+\lambda_{5}}{2}v^{2}-\frac{\lambda_{6}v^{2}}{2}t^{-1}_{\beta}-\frac{\lambda_{7}v^{2}}{2}t_{\beta}\,,$
$\displaystyle+\sqrt{2}v_{\Delta}\left(s^{2}_{\beta}\mu_{1}-c^{2}_{\beta}\mu_{2}\right)+\frac{\mu_{3}v_{\Delta}}{\sqrt{2}c_{\beta}s_{\beta}}\left(c^{4}_{\beta}+s^{4}_{\beta}\right)-\frac{v^{2}_{\Delta}}{2}\left(\lambda^{\prime}_{8}s^{2}_{\beta}+\lambda^{\prime}_{9}c^{2}_{\beta}\right)$
$\displaystyle-\frac{v^{2}_{\Delta}}{2}\left(s_{2\beta}\lambda_{12}+(c^{3}_{\beta}+s^{3}_{\beta})\bar{\lambda}_{12}\right)\,,$
$\displaystyle m^{2}_{H^{-}\delta^{+}}$
$\displaystyle=\frac{v}{2}\left(s_{2\beta}(\mu_{1}-\mu_{2})-c_{2\beta}\mu_{3}\right)-\frac{\lambda^{\prime}_{8}-\lambda^{\prime}_{9}}{2\sqrt{2}}c_{\beta}s_{\beta}vv_{\Delta}-\frac{\lambda^{\prime}_{12}}{2\sqrt{2}}s^{2}_{\beta}vv_{\Delta}\,,$
$\displaystyle m^{2}_{\delta^{-}\delta^{+}}$ $\displaystyle\equiv
m^{2}_{\delta^{\pm}}=m^{2}_{\Delta}+\frac{2\lambda_{9}+\lambda^{\prime}_{9}}{4}s^{2}_{\beta}v^{2}+\frac{2\lambda_{8}+\lambda^{\prime}_{8}}{4}c^{2}_{\beta}v^{2}+\frac{2\lambda_{12}+\lambda^{\prime}_{12}}{4}s_{2\beta}v^{2}\,.$
(54)

Analogous to the case of CP-odd scalar mass matrix, $m^{2}_{G^{-}G^{+}}$,
$m^{2}_{G^{-}H^{+}}$, and $m^{2}_{G^{-}\delta^{+}}$ vanish if we drop the
${\cal O}(v_{\Delta})$ terms. Therefore, $G^{\pm}$ can be approximated as the
charged Goldstone bosons. If we further demand
$\mu_{3}=t_{2\beta}(\mu_{1}-\mu_{2})$, $m^{2}_{H^{-}\delta^{+}}$ also
vanishes, and $H^{\pm}$ and $\delta^{\pm}$ decouple each other and are
approximately the physical states.

Next, we discuss the CP-even scalars. In 2HDMs, the physical $(H,h)$ states
and the $(\Phi^{0}_{1},\Phi^{0}_{2})$ states are related by:

$\begin{pmatrix}H\\\ h\end{pmatrix}=\begin{pmatrix}c_{\alpha}&s_{\alpha}\\\
-s_{\alpha}&c_{\alpha}\end{pmatrix}\begin{pmatrix}\Phi^{0}_{1}\\\
\Phi^{0}_{2}\end{pmatrix}\,,$ (55)

where $\alpha$ is the mixing angle for CP-even scalars. To derive the mass
matrix in the $(h,H,\delta^{0})$ states, it is more convenient to start from
the basis of $(\Phi^{0}_{1},\Phi^{0}_{2})$. Once the mass eigenvalues and
eigenstates of $h$ and $H$ in the 2HDM are obtained, we then include
$\delta^{0}$ to form a three-component basis $(h,H,\delta^{0})$. As we will
explicitly see below, these will be approximately the physical states as long
as ${\cal O}(v_{\Delta})$ terms are ignored.

The $2\times 2$ mass matrix for $(\Phi^{0}_{1},\Phi^{0}_{2})$ is:

$\frac{1}{2}\begin{pmatrix}\Phi^{0}_{1}~{}\Phi^{0}_{2}\end{pmatrix}\begin{pmatrix}m^{2}_{\phi_{1}\phi_{1}}&m^{2}_{\phi_{1}\phi_{2}}\\\
m^{2}_{\phi_{1}\phi_{2}}&m^{2}_{\phi_{2}\phi_{2}}\end{pmatrix}\begin{pmatrix}\Phi^{0}_{1}\\\
\Phi^{0}_{2}\end{pmatrix}\,,$ (56)

where the matrix elements

$\displaystyle m^{2}_{\phi_{1}\phi_{1}}$
$\displaystyle=\left(m^{2}_{12}+\frac{v_{\Delta}\mu_{3}}{\sqrt{2}}\right)t_{\beta}+\lambda_{1}v^{2}c^{2}_{\beta}+\frac{3\lambda_{6}v^{2}}{2}c_{\beta}s_{\beta}-\frac{\lambda_{7}v^{2}}{2}s^{2}_{\beta}t_{\beta}\,,$
$\displaystyle m^{2}_{\phi_{1}\phi_{2}}$
$\displaystyle=-\left(m^{2}_{12}+\frac{v_{\Delta}\mu_{3}}{\sqrt{2}}\right)+\lambda_{345}v^{2}c_{\beta}s_{\beta}+\frac{3v^{2}}{2}\left(c^{2}_{\beta}\lambda_{6}+s^{2}_{\beta}\lambda_{7}\right)\,,$
$\displaystyle m^{2}_{\phi_{2}\phi_{2}}$
$\displaystyle=\left(m^{2}_{12}+\frac{v_{\Delta}\mu_{3}}{\sqrt{2}}\right)t^{-1}_{\beta}+\lambda_{2}v^{2}s^{2}_{\beta}-\frac{\lambda_{6}v^{2}}{2}c^{2}_{\beta}t^{-1}_{\beta}+\frac{3\lambda_{7}v^{2}}{2}c_{\beta}s_{\beta}\,.$
(57)

Using the parametrization in Eq. (55), we obtain the eigenvalues of Eq. (55)
as:

$\displaystyle
m^{2}_{h,H}=\frac{m^{2}_{\phi_{1}\phi_{1}}+m^{2}_{\phi_{2}\phi_{2}}}{2}\pm\frac{1}{2}\left[\left(m^{2}_{\phi_{1}\phi_{1}}-\phi^{2}_{\phi_{2}\phi_{2}}\right)^{2}+4(m^{2}_{\phi_{1}\phi_{2}})^{2}\right]^{1/2}\,,$
(58)

and the mixing angle is determined by:

$\tan{2\alpha}=\frac{2m^{2}_{\phi_{1}\phi_{2}}}{m^{2}_{\phi_{1}\phi_{1}}-m^{2}_{\phi_{2}\phi_{2}}}\,.$
(59)

Using the rotational matrix,

$\begin{pmatrix}-s_{\alpha}&c_{\alpha}&0\\\ c_{\alpha}&s_{\alpha}&0\\\
0&0&1\end{pmatrix}\,,$ (60)

we can transform the basis from $(\Phi^{0}_{1},\Phi^{0}_{2},\delta^{0})$ to
$(h,H,\delta^{0})$, and the transformed mass matrix is:

$\frac{1}{2}\left(h~{}H~{}\delta^{0}\right)\begin{pmatrix}m^{2}_{h}&0&m_{h\delta^{0}}^{2}\\\
0&m_{H}^{2}&m_{H\delta^{0}}^{2}\\\
m_{h\delta^{0}}^{2}&m_{H\delta^{0}}^{2}&m_{\delta^{0}\delta^{0}}^{2}\end{pmatrix}\begin{pmatrix}h\\\
H\\\ \delta^{0}\end{pmatrix}\,,$ (61)

where the additional elements

$\displaystyle m_{h\delta^{0}}^{2}=$
$\displaystyle-\frac{v}{\sqrt{2}}\left[(\mu_{1}+\mu_{2})s_{\beta-\alpha}-(\mu_{1}-\mu_{2})s_{\beta+\alpha}+\mu_{3}c_{\beta+\alpha}\right]$
$\displaystyle-\frac{v_{\Delta}v}{2}\left[(\bar{\lambda}_{8}+\bar{\lambda}_{9})s_{\beta-\alpha}+(\bar{\lambda}_{8}-\bar{\lambda}_{9})s_{\beta+\alpha}-2\bar{\lambda}_{12}c_{\beta+\alpha}\right]\,,$
$\displaystyle m_{H\delta^{0}}^{2}=$
$\displaystyle-\frac{v}{\sqrt{2}}\left[(\mu_{1}+\mu_{2})c_{\beta-\alpha}+(\mu_{1}-\mu_{2})c_{\beta+\alpha}+\mu_{3}s_{\beta+\alpha}\right]$
$\displaystyle+\frac{v_{\Delta}v}{2}\left[(\bar{\lambda}_{8}+\bar{\lambda}_{9})c_{\beta-\alpha}+(\bar{\lambda}_{8}-\bar{\lambda}_{9})c_{\beta+\alpha}+2\bar{\lambda}_{12}s_{\beta+\alpha}\right]\,,$
$\displaystyle m_{\delta^{0}\delta^{0}}^{2}=$ $\displaystyle
m^{2}_{\Delta}+3(\lambda_{\Delta 1}+\lambda_{\Delta
2})v^{2}_{\Delta}+\frac{v^{2}}{2}\left(c^{2}_{\beta}\bar{\lambda}_{8}+s^{2}_{\beta}\bar{\lambda}_{9}+s_{2\beta}\bar{\lambda}_{12}\right)\,.$
(62)

In the model, $m^{2}_{h\delta^{0}}$ and $m^{2}_{H\delta^{0}}$ are generally
not small. Nevertheless, if we apply the conditions in Eq. (51),
$m^{2}_{h\delta^{0}}$ and $m^{2}_{H\delta^{0}}$ can be rewritten as:

$\displaystyle m^{2}_{h\delta^{0}}=$
$\displaystyle-\sqrt{2}v\epsilon_{\Delta}s_{\beta-\alpha}-\frac{v_{\Delta}v}{2}\left[(\bar{\lambda}_{8}+\bar{\lambda}_{9})s_{\beta-\alpha}+(\bar{\lambda}_{8}-\bar{\lambda}_{9})s_{\beta+\alpha}-\bar{\lambda}_{12}c_{\beta+\alpha}\right]\,,$
$\displaystyle m_{H\delta^{0}}^{2}=$
$\displaystyle-\sqrt{2}v\epsilon_{\Delta}s_{\beta+\alpha}+\frac{v_{\Delta}v}{2}\left[(\bar{\lambda}_{8}+\bar{\lambda}_{9})c_{\beta-\alpha}+(\bar{\lambda}_{8}-\bar{\lambda}_{9})c_{\beta+\alpha}+\bar{\lambda}_{12}s_{\beta+\alpha}\right]\,,$
(63)

both being $\sim O(v_{\Delta})$. Comparing to the dominant diagonal elements
$m^{2}_{h,H,\delta^{0}}$, the mixing effects among $h$, $H$, and $\delta^{0}$
are thus small. Therefore, to the leading order in $v_{\Delta}$, we will
neglect such mixing and $h$ and $H$ decouple from $\delta^{0}$.

##  Appendix B Scalar potential and trilinear couplings in the Higgs basis

The scalar potential in terms of the Higgs basis $H_{1,2}$ can be written as:

$V=V(H_{1},H_{2})+V(\Delta)+V(H_{1},H_{2},\Delta)\,,$ (64)

where each term is more explicitly given by

$\displaystyle V(H_{1},H_{2})=$ $\displaystyle
Y_{1}H^{\dagger}_{1}H_{1}+Y_{2}H^{\dagger}_{2}H_{2}+Y_{3}(H^{\dagger}_{1}H_{2}+\mbox{H.c.})+\frac{1}{2}\Lambda_{1}(H^{\dagger}_{1}H_{1})^{2}$
$\displaystyle+\frac{1}{2}\Lambda_{2}(H^{\dagger}_{2}H_{2})^{2}+\Lambda_{3}H^{\dagger}_{1}H_{1}H^{\dagger}_{2}H_{2}+\Lambda_{4}H^{\dagger}_{1}H_{2}H^{\dagger}_{2}H_{1}+\left[\frac{1}{2}\Lambda_{5}(H^{\dagger}_{1}H_{2})^{2}\right.$
$\displaystyle\left.+\Lambda_{6}(H^{\dagger}_{1}H_{1})(H^{\dagger}_{2}H_{1})+\Lambda_{7}(H^{\dagger}_{2}H_{2})(H^{\dagger}_{2}H_{1})+\mbox{H.c.}\right]\,,$
(65) $\displaystyle V(\Delta)=$ $\displaystyle
M^{2}_{\Delta}Tr\Delta^{\dagger}\Delta+\Lambda_{\Delta
1}(Tr\Delta^{\dagger}\Delta)^{2}+\Lambda_{\Delta
2}Tr(\Delta^{\dagger}\Delta)^{2}\,,$ (66) $\displaystyle
V(H_{1},H_{2},\Delta)=$
$\displaystyle\left(W_{1}H^{T}_{1}i\tau_{2}\Delta^{\dagger}H_{1}+W_{2}H^{T}_{2}i\tau_{2}\Delta^{\dagger}H_{2}+W_{3}H^{T}_{1}i\tau_{2}\Delta^{\dagger}H_{2}+\mbox{H.c.}\right)$
$\displaystyle+\left[\Lambda_{8}H^{\dagger}_{1}H_{1}+\Lambda_{9}H^{\dagger}_{2}H_{2}+(\Lambda_{12}H^{\dagger}_{1}H_{2}+\mbox{H.c.})\right]Tr\Delta^{\dagger}\Delta$
$\displaystyle+\Lambda^{\prime}_{8}H^{\dagger}_{1}\Delta\Delta^{\dagger}H_{1}+\Lambda^{\prime}_{9}H^{\dagger}_{2}\Delta\Delta^{\dagger}H_{2}+(\Lambda^{\prime}_{12}H^{\dagger}_{1}\Delta\Delta^{\dagger}H_{2}+\mbox{H.c.})\,.$
(67)

Using the representations of $H_{1,2}$ and $\Delta$ in Eq. (5), the minimum
conditions for the the VEVs of $H_{1,2}$ and $\Delta$ can be obtained as:

$\displaystyle Y^{2}_{1}+\frac{\Lambda^{2}_{1}}{2}v^{2}$
$\displaystyle=\sqrt{2}W_{1}v_{D}-\frac{\bar{\Lambda}_{8}}{2}v^{2}_{D}\,,$
$\displaystyle Y^{2}_{3}+\frac{\Lambda^{2}_{6}}{2}v^{2}$
$\displaystyle=\frac{W_{3}}{\sqrt{2}}v_{D}-\frac{\bar{\Lambda}_{12}}{2}v^{2}_{D}\,,$
$\displaystyle M^{2}_{\Delta}+\frac{\bar{\Lambda}_{8}}{2}v^{2}$
$\displaystyle=\frac{v^{2}}{\sqrt{2}}\frac{W_{1}}{v_{\Delta}}-(\Lambda_{\Delta
1}+\Lambda_{\Delta 2})v^{2}_{\Delta}\,,$ (68a)

with $\bar{\Lambda}_{8}=\Lambda_{8}+\Lambda^{\prime}_{8}$ and
$\bar{\Lambda}_{12}=\Lambda_{12}+\Lambda^{\prime}_{12}$.

From Eqs. (65)–(67), the doubly-charged Higgs mass is obtained as:

$m^{2}_{\delta^{\pm\pm}}=m^{2}_{\Delta}+\frac{\Lambda_{8}v^{2}}{2}+\Lambda_{\Delta
1}v^{2}_{\Delta}\,.$ (69)

The mass splittings in the Higgs triplet are:

$\displaystyle m^{2}_{\delta^{\pm}}-m^{2}_{\delta^{\pm\pm}}$
$\displaystyle=\frac{\Lambda_{\Delta
2}v^{2}_{\Delta}}{2}+\frac{\Lambda^{\prime}_{8}v^{2}}{4}\,,$ $\displaystyle
m^{2}_{\delta^{0}}-m^{2}_{\delta^{\pm\pm}}$ $\displaystyle=(2\Lambda_{\Delta
1}+3\Lambda_{\Delta 2})v^{2}_{\Delta}+\frac{\Lambda^{\prime}_{8}v^{2}}{2}\,.$
(70)

The trilinear couplings of the neutral scalars to the charged Higgs scalars
with $c_{\beta-\alpha}=0$ are given by:

${\cal
L}_{H_{i}SS}=-v\left[\lambda_{H_{i}\delta^{--}\delta^{++}}H_{i}\delta^{--}\delta^{++}+\lambda_{H_{i}\delta^{-}\delta^{+}}H_{i}\delta^{-}\delta^{+}+\lambda_{H_{i}H^{-}H^{+}}H_{i}H^{-}H^{+}\right]\,,$
(71)

where the couplings:

$\displaystyle\lambda_{h\delta^{--}\delta^{++}}$
$\displaystyle=\Lambda_{8}\,,~{}\Lambda_{H\delta^{--}\delta^{++}}=-\Lambda_{12}\,,$
$\displaystyle\lambda_{h\delta^{-}\delta^{+}}$
$\displaystyle=\Lambda_{8}+\frac{\Lambda^{\prime}_{8}}{2}\,,~{}\lambda_{H\delta^{-}\delta^{+}}=-\Lambda_{12}-\frac{\Lambda^{\prime}_{12}}{2}\,,$
$\displaystyle\lambda_{hH^{-}H^{+}}$
$\displaystyle=\Lambda_{3}\,,~{}\lambda_{HH^{-}H^{+}}=-\Lambda_{7}\,.$ (72)

## References

  * (1) G. W. Bennett et al. [Muon g-2], Phys. Rev. D 73, 072003 (2006) [arXiv:hep-ex/0602035 [hep-ex]].
  * (2) T. Aoyama, N. Asmussen, M. Benayoun, J. Bijnens, T. Blum, M. Bruno, I. Caprini, C. M. Carloni Calame, M. Cè and G. Colangelo, et al. Phys. Rept. 887, 1-166 (2020) [arXiv:2006.04822 [hep-ph]].
  * (3) B. Abi et al. [Muon g-2], Phys. Rev. Lett. 126, 141801 (2021) [arXiv:2104.03281 [hep-ex]].
  * (4) S. Borsanyi et al., Nature (2021), 2002.12347. [arXiv:2002.12347 [hep-lat]].
  * (5) J. P. Leveille, Nucl. Phys. B 137, 63-76 (1978).
  * (6) H. E. Haber, G. L. Kane and T. Sterling, Nucl. Phys. B 161, 493-532 (1979)
  * (7) M. Krawczyk and J. Zochowski, Phys. Rev. D 55, 6968-6974 (1997) [arXiv:hep-ph/9608321 [hep-ph]].
  * (8) A. Dedes and H. E. Haber, JHEP 05, 006 (2001) [arXiv:hep-ph/0102297 [hep-ph]].
  * (9) K. m. Cheung, C. H. Chou and O. C. W. Kong, Phys. Rev. D 64, 111301 (2001) [arXiv:hep-ph/0103183 [hep-ph]].
  * (10) F. Larios, G. Tavares-Velasco and C. P. Yuan, Phys. Rev. D 64, 055004 (2001) [arXiv:hep-ph/0103292 [hep-ph]].
  * (11) C. H. Chen and C. Q. Geng, Phys. Lett. B 511, 77-84 (2001) [arXiv:hep-ph/0104151 [hep-ph]].
  * (12) A. Arhrib and S. Baek, Phys. Rev. D 65, 075002 (2002) [arXiv:hep-ph/0104225 [hep-ph]].
  * (13) M. Krawczyk, Acta Phys. Polon. B 33, 2621-2634 (2002) [arXiv:hep-ph/0208076 [hep-ph]].
  * (14) J. F. Gunion, JHEP 08, 032 (2009) [arXiv:0808.2509 [hep-ph]].
  * (15) A. Broggio, E. J. Chun, M. Passera, K. M. Patel and S. K. Vempati, JHEP 11, 058 (2014) [arXiv:1409.3199 [hep-ph]].
  * (16) V. Ilisie, JHEP 04, 077 (2015) [arXiv:1502.04199 [hep-ph]].
  * (17) T. Han, S. K. Kang and J. Sayre, JHEP 02, 097 (2016) [arXiv:1511.05162 [hep-ph]].
  * (18) S. P. Li, X. Q. Li and Y. D. Yang, Phys. Rev. D 99, no.3, 035010 (2019) [arXiv:1808.02424 [hep-ph]].
  * (19) X. F. Han, T. Li, L. Wang and Y. Zhang, Phys. Rev. D 99 (2019) no.9, 095034 [arXiv:1812.02449 [hep-ph]].
  * (20) F. J. Botella, F. Cornet-Gomez and M. Nebot, Phys. Rev. D 102 (2020) no.3, 035023 [arXiv:2006.01934 [hep-ph]].
  * (21) L. Delle Rose, S. Khalil and S. Moretti, [arXiv:2012.06911 [hep-ph]].
  * (22) D. Sabatta, A. S. Cornell, A. Goyal, M. Kumar, B. Mellado and X. Ruan, Chin. Phys. C 44, no.6, 063103 (2020) [arXiv:1909.03969 [hep-ph]].
  * (23) S. Jana, V. P. K. and S. Saad, Phys. Rev. D 101, no.11, 115037 (2020) [arXiv:2003.03386 [hep-ph]].
  * (24) S. P. Li, X. Q. Li, Y. Y. Li, Y. D. Yang and X. Zhang, JHEP 01, 034 (2021) [arXiv:2010.02799 [hep-ph]].
  * (25) N. Chen, B. Wang and C. Y. Yao, [arXiv:2102.05619 [hep-ph]].
  * (26) C. H. Chen and T. Nomura, Phys. Rev. D 90, no.7, 075008 (2014) [arXiv:1406.6814 [hep-ph]].
  * (27) B. A. Ouazghour, A. Arhrib, R. Benbrik, M. Chabab and L. Rahili, Phys. Rev. D 100, no.3, 035031 (2019) [arXiv:1812.07719 [hep-ph]].
  * (28) A. Pich and P. Tuzon, Phys. Rev. D 80, 091702 (2009) [arXiv:0908.1554 [hep-ph]].
  * (29) M. Jung, A. Pich and P. Tuzon, JHEP 11, 003 (2010) [arXiv:1006.0470 [hep-ph]].
  * (30) A. Celis, V. Ilisie and A. Pich, JHEP 07, 053 (2013) [arXiv:1302.4022 [hep-ph]].
  * (31) A. Celis, V. Ilisie and A. Pich, JHEP 12, 095 (2013) [arXiv:1310.7941 [hep-ph]].
  * (32) V. Ilisie and A. Pich, JHEP 09, 089 (2014) [arXiv:1405.6639 [hep-ph]].
  * (33) M. Magg and C. Wetterich, Phys. Lett. B 94, 61 (1980).
  * (34) G. Lazarides, Q. Shafi and C. Wetterich, Nucl. Phys. B 181, 287 (1981).
  * (35) J. Schechter and J. W. F. Valle, Phys. Rev. D 22, 2227 (1980).
  * (36) T. P. Cheng and L. -F. Li, Phys. Rev. D 22, 2860 (1980).
  * (37) R. N. Mohapatra and G. Senjanovic, Phys. Rev. D 23, 165 (1981).
  * (38) S. M. Bilenky, J. Hosek and S. T. Petcov, Phys. Lett. B 94, 495 (1980).
  * (39) S. M. Barr and A. Zee, Phys. Rev. Lett. 65, 21-24 (1990) [erratum: Phys. Rev. Lett. 65, 2920 (1990)].
  * (40) D. Chang, W. F. Chang, C. H. Chou and W. Y. Keung, Phys. Rev. D 63, 091301 (2001) [arXiv:hep-ph/0009292 [hep-ph]].
  * (41) G. Aad et al. [ATLAS], arXiv:2101.11961 [hep-ex].
  * (42) P. A. Zyla et al. (Particle Data Group), Prog. Theor. Exp. Phys. 2020, 083C01 (2020).
  * (43) S. L. Glashow and S. Weinberg, Phys. Rev. D 15, 1958 (1977).
  * (44) L. J. Hall and M. B. Wise, Nucl. Phys. B 187, 397-408 (1981).
  * (45) J. F. Donoghue and L. F. Li, Phys. Rev. D 19, 945 (1979).
  * (46) V. D. Barger, J. L. Hewett and R. J. N. Phillips, Phys. Rev. D 41, 3421-3441 (1990).
  * (47) M. J. Savage, Phys. Lett. B 266, 135-141 (1991).
  * (48) Y. Grossman, Nucl. Phys. B 426, 355-384 (1994) [arXiv:hep-ph/9401311 [hep-ph]].
  * (49) A. G. Akeroyd and W. J. Stirling, Nucl. Phys. B 447, 3-17 (1995).
  * (50) G. C. Branco, P. M. Ferreira, L. Lavoura, M. N. Rebelo, M. Sher and J. P. Silva, Phys. Rept. 516, 1-102 (2012) [arXiv:1106.0034 [hep-ph]].
  * (51) Y. H. Ahn and C. H. Chen, Phys. Lett. B 690, 57-61 (2010) [arXiv:1002.4216 [hep-ph]].
  * (52) C. H. Chen, C. W. Chiang, T. Nomura and Y. Fu-Sheng, Phys. Rev. D 85, 075018 (2012) [arXiv:1105.2870 [hep-ph]].
  * (53) C. H. Chen and T. Nomura, Phys. Lett. B 767, 443-449 (2017) [arXiv:1609.01874 [hep-ph]].
  * (54) N. Cabibbo, Phys. Rev. Lett. 10, 531-533 (1963).
  * (55) M. Kobayashi and T. Maskawa, Prog. Theor. Phys. 49, 652-657 (1973).
  * (56) B. Pontecorvo, Sov. Phys. JETP 6, 429 (1957) [Zh. Eksp. Teor. Fiz. 33, 549 (1957)].
  * (57) Z. Maki, M. Nakagawa and S. Sakata, Prog. Theor. Phys. 28, 870 (1962).
  * (58) R. Benbrik, C. H. Chen and T. Nomura, Phys. Rev. D 93, no.9, 095004 (2016) [arXiv:1511.08544 [hep-ph]].
  * (59) J. F. Gunion and H. E. Haber, Phys. Rev. D 67, 075019 (2003) [arXiv:hep-ph/0207010 [hep-ph]].
  * (60) M. Carena, I. Low, N. R. Shah and C. E. M. Wagner, JHEP 04, 015 (2014) [arXiv:1310.2248 [hep-ph]].
  * (61) J. Bernon, J. F. Gunion, H. E. Haber, Y. Jiang and S. Kraml, Phys. Rev. D 92, no.7, 075004 (2015) doi:10.1103/PhysRevD.92.075004 [arXiv:1507.00933 [hep-ph]].
  * (62) C. H. Chen and T. Nomura, Nucl. Phys. B 940, 292-311 (2019) [arXiv:1705.10620 [hep-ph]].
  * (63) C. H. Chen and T. Nomura, [arXiv:2001.07515 [hep-ph]].
  * (64) O. Eberhardt, A. P. Martínez and A. Pich, [arXiv:2012.09200 [hep-ph]].
  * (65) L. Michel, Proc. Phys. Soc. A 63, 514-531 (1950).
  * (66) T. Abe, R. Sato and K. Yagyu, JHEP 07, 064 (2015) [arXiv:1504.07059 [hep-ph]].
  * (67) A. Crivellin, J. Heeck and P. Stoffer, Phys. Rev. Lett. 116, no.8, 081801 (2016) [arXiv:1507.07567 [hep-ph]].
  * (68) E. J. Chun, Z. Kang, M. Takeuchi and Y. L. S. Tsai, JHEP 11 (2015), 099 [arXiv:1507.08067 [hep-ph]].
  * (69) T. Enomoto and R. Watanabe, JHEP 05 (2016), 002 [arXiv:1511.05066 [hep-ph]].
  * (70) K. Tobe, JHEP 10, 114 (2016) [arXiv:1607.04447 [hep-ph]].
  * (71) W. Fetscher, H. J. Gerber and K. F. Johnson, Phys. Lett. B 173, 102-106 (1986).
  * (72) T. Aoyama, M. Hayakawa, T. Kinoshita and M. Nio, Phys. Rev. D 91, no. 3, 033006 (2015) Erratum: [Phys. Rev. D 96, no. 1, 019901 (2017)] [arXiv:1412.8284 [hep-ph]];
  * (73) T. Aoyama, T. Kinoshita, and M. Nio, Atoms 7, 28 (2019).
  * (74) R. H. Parker, C. Yu, W. Zhong, B. Estey, and H. Müeller, Science 360, 191 (2018).
  * (75) L. Morel, Z. Yao, P. Clad$\rm\acute{e}$, and S. Guellati-Kh$\rm\acute{e}$lifa, Science 588, 61 (2020).
  * (76) P. F. de Salas, D. V. Forero, C. A. Ternes, M. Tortola and J. W. F. Valle, Phys. Lett. B 782, 633-640 (2018) [arXiv:1708.01186 [hep-ph]].