# Dual electroweak phase transition in the two-Higgs-doublet model with the
$S_{3}$ discrete symmetry

Vo Quoc Phonga,b<EMAIL_ADDRESS>aDepartment of Theoretical Physics,
University of Science, Ho Chi Minh City 70000, Vietnam
bVietnam National University, Ho Chi Minh City 70000, Vietnam  Nguyen Minh
Anha,b<EMAIL_ADDRESS>aDepartment of Theoretical Physics, University of
Science, Ho Chi Minh City 70000, Vietnam
bVietnam National University, Ho Chi Minh City 70000, Vietnam  Hoang Ngoc
Longc,d<EMAIL_ADDRESS>(corresponding author) c Subatomic Physics
Research Group, Science and Technology Advanced Institute, Van Lang
University, Ho Chi Minh City 70000, Vietnam
d Faculty of Applied Technology, School of Technology, Van Lang University, Ho
Chi Minh City 70000, Vietnam

###### Abstract

In this work, dual electroweak phase transition (EWPT) consisting of two
phases, is carefully studied in the two-Higgs-doublet model with the $S_{3}$
discrete symmetry. The role of $S_{3}$ here is to further separate the stages
of the electroweak phase transition, compared to that of the original two-
Higgs-doublet model (2HDM). The strength of the electroweak phase transition
$(S)$ in the model under consideration is large enough for the first-order
EWPT, specifically $1<S<2.8$. The ratio between the two vacuum expectation
values (VEVs), $\tan\beta=v_{2}/v_{1}$, is proven to have no effects on the
strength of the phase transition. This ratio only affects the mass domain that
causes the first-order phase transition. Furthermore, in this paper we will
show clearly that when studying the EWPT in models of more than one scalar
field that generates masses, one needs to analyze the problem of phase
transition under multiple stages. In other words, the effect of the first
stage of symmetry breaking to the second one, is to simplify by suggestion
that vacuum expectation value of the Higgs boson responsible for the initial
stage is proportional to that of the field for the next stage.

###### pacs:

11.15.Ex, 12.60.Fr, 98.80.Cq

Keywords: Spontaneous breaking of gauge symmetries, Extensions of electroweak
Higgs sector, Particle-theory models (Early Universe)

###### Contents

  1. I INTRODUCTION
  2. II Review on the Higgs potential and comments on EWPT in the 2HDM
    1. II.1 The Higgs potential in the 2HDM
    2. II.2 Higgs potential at the tree level
    3. II.3 The masses of gauge bosons
    4. II.4 Remarks on EWPT structure in the 2HDM
    5. II.5 Comments on EWPT in the 2HDM
  3. III Review on the 2HDM-$S_{3}$
    1. III.1 Particle content
    2. III.2 Higgs potential
    3. III.3 The soft breaking of $S_{3}$ group
  4. IV Electroweak phase transition in the 2HDM-$S_{3}$
    1. IV.1 A vital role of $S_{3}$
    2. IV.2 Structure of EWPT
    3. IV.3 The effective potential
    4. IV.4 Probing the independence of EWPT strength on $\tan\beta$
    5. IV.5 The true critical temperatures
    6. IV.6 Searching the first-order EWPT and the role of $\tan\beta$
  5. V CONCLUSION AND OUTLOOKS
  6. A Effective action of multiscalar field models

##  I INTRODUCTION

The Standard Model (SM), an outstanding achievement of physics in particular,
and a memorable milestone for the scientific community in general is a
systematic theory of elementary particles and their interactions. The model
predicted the results of many experiments, the existence of the Higgs
particle; together with the Higgs mechanism, it shows us the nature of
subatomic particles. However, the model still has some shortcomings such as
not being able to unify gravity, describe dark matter or small neutrino mass,
etc.

One of the significant phenomena in cosmology that cannot be explained by the
SM is baryon asymmetry, also known as matter-antimatter asymmetry. This
problem is explains why there is an imbalance between matter and antimatter in
the universe. For a strong first-order electroweak phase transition, the third
of Sakharov’s three conditions sakharov , plays an important role in
explaining this asymmetry. It indeed important because this condition not only
explains the thermal imbalance but also provides a link between the violation
of B and CP and the other two conditions of Sakharov.

The thermal imbalance is expressed through a first-order electroweak phase
transition (EWPT) which should be considered first. The SM does not have
enough triggers for a first-order phase transition mkn ; SME ; SMEb ; SMEc ;
SMEd ; michela . Therefore, in the beyond SM, this problem must be considered
(see, for example, Refs. plv ; 2b ; 2c ; BSM ; BSMb ; majorana ; majoranab ;
thdm ; thdmb ; ESMCO ; elptdm ; elptdma ; elptdmb ; elptdmc ; elptdmd ;
phonglongvan ; phonglongvanb ; phonglongvan2 ; SMS ; dssm ; munusm ; lr ;
singlet ; singletb ; singletc ; singletd ; mssm1 ; mssm1b ; mssm1c ; twostep ;
twostepb ; twostepc ; 1101.4665 ; 1101.4665b ; 1101.4665c ; jjgb ; jjgc ; jjgd
; Ahriche1 ; Ahriche2 ; Ahriche2b ; Ahriche3 ; Ahriche3b ; Ahriche3c ;
Ahriche3d ; Ahriche3e ; Fuyuto ; Fuyutob ; Fuyutoc ; span ; chr ; cde ;
kusenko ). The different scenarios that can be enumerated in these references
are as follows: doing high-temperature effective potential, analysis of the
trigger roles of new particles, the decoupling conditions, bubble nucleations,
and sphalerons.

The triggers for the EWPT can be new particles (beyond SM) or parameter
corrections in SM plv ; 2b ; 2c ; BSM ; BSMb ; majorana ; majoranab ; thdm ;
thdmb ; ESMCO ; elptdm ; elptdma ; elptdmb ; elptdmc ; elptdmd ; phonglongvan
; phonglongvanb ; phonglongvan2 ; SMS ; dssm ; munusm ; lr ; singlet ;
singletb ; singletc ; singletd ; mssm1 ; mssm1b ; mssm1c ; twostep ; twostepb
; twostepc ; 1101.4665 ; 1101.4665b ; 1101.4665c ; jjgb ; jjgc ; jjgd ;
Ahriche2 ; Ahriche2b ; Ahriche3 ; Ahriche3b ; Ahriche3c ; Ahriche3d ;
Ahriche3e ; Fuyutob ; Fuyutoc ; span ; chr ; cde ; kusenko . In the SM, it is
contrary to the experiment that the mass of Higgs boson must be less than
$125$ GeV, for a strong first-order EWPT mkn ; SME ; SMEb ; SMEc ; SMEd ;
michela . In a model, if the new particles are the cause of the violent EWPT,
then that model can have more than one Higgs field plv ; 2b ; 2c ; majorana ;
majoranab ; thdm ; thdmb ; ESMCO ; elptdm ; elptdma ; elptdmb ; elptdmc ;
elptdmd ; phonglongvan ; phonglongvanb ; phonglongvan2 ; epjc ; zb ; singlet ;
singletb ; singletc ; singletd ; mssm1 ; mssm1b ; mssm1c ; twostep ; twostepb
; twostepc ; chiang3 . Another interesting point but consistent with the
physical nature, the strength of EWPT is gauge-independent zb ; 1101.4665 ;
1101.4665b ; 1101.4665c ; Arefe . The self-energy term or daisy loops cause a
problem for effective potentials at high temperatures. However, it is not the
main trigger for EWPT and it reduces the strength of EWPT r23 .

By making sure that the C and CP violations exist, the third condition of
baryogenesis is given as
$\Gamma_{sph}\sim\mathcal{A}(T)\textrm{exp}\left(\frac{-E_{sph}}{T}\right)\ll
H_{rad}$ spha-huble ; decoupling ; decouplingb ; decouplingc  in the context
of topological transitions, where $\Gamma_{sph}$ and $E_{sph}$ are the
sphaleron rate and energy, respectively, $H_{rad}$ is the Hubble expansion
rate in the radiation-dominated period, and $\mathcal{A}(T)\approx T^{4}$.
This is referred to as the sphaleron decoupling situation. This condition is
frequently written as $S=v_{c}/T_{c}>1$ in the SM using the approximation
$E_{sph}(T)\approx[v(T)/v]E_{sph}(T=0)$ 47 ; Ahriche1 ; Ahriche2 ; Ahriche2b .
However, this approximation should be used with caution in models that go
beyond the SM.

Through the above brief summaries, ones have another aspect of EWPT survey.
Currently, there are two scenarios as follows: the first one is that the EWPT
process has only one stage; the second is that this process has two or three
stages. Some of the theoretical models that study one stage of EWPT are the
SM, Zee-Babu zb , SMEFT pkll , SM with corrections to the Yukawa interaction
for quarks yukawa1 ; yukawa2 . Models that study more than one stage of EWPT
are the ones that consist of more than one Higgs field. Some of them as
follows: 3-3-1 models plv ; phonglongvan2 ; ptl , 2-2-1 model pa .

However, there is one quite special model, the 2HDM has two vacuum expectation
values (VEVs) having values/ranges in the electroweak scale. In this paper,
the ways of studying the EWPT with only one stage of this model will be
reconsidered. Then the strengths and weaknesses of it will be analyzed. At the
same time, there are many interesting versions of the model that go beyond
itself, in which there is a model of two Higgs doublets with symmetry of
$S_{3}$. This symmetry can account for quark mixing bpal1 ; bpal2 ; bpal3 .
The $S_{3}$ symmetry has been proposed as the basic flavor symmetry in various
frameworks. This discrete symmetry in the lepton sector is to produce a
$\mu-\tau$ symmetry 28t ; 29t ; 30t ; 31t  or a tribimaximal neutrino mixing
matrix 32t ; 33t . In the quark sector, this symmetry can produce Fritzsch and
Fritzsch mass textures 34t . In addition, the nearest neighbor interactions
(NNI) mass texture is hidden in a $S_{3}$ flavor symmetry 35t . Therefore, the
model with $S_{3}$ is interested, since $S_{3}$ would simplify the Higgs
sector, which is very important in studying EWPT. Furthermore, since there are
only two VEVs, our work of studying the multi-stages of the 2HDM-$S_{3}$ would
be relatively easier compared to that of the SM.

The two-Higgs-doublet model with $S_{3}$ symmetry (2HDM-$S_{3}$) 1601 , one of
the extended versions the 2HDM, that has the potential to ”possess” a strongly
first-order electroweak phase transition, because of the following factors:
the model has the heavy Higgs boson, as well as the charged Higgs boson; at
the same time, it obeys the smallest non-Abelian discrete symmetry group.

More specifically, we will consider whether the first-order phase transition
in the 2HDM-$S_{3}$ is strong or not. And when it is strong, the range of
values of the phase transition strength and the mass of the new particles and
related parameters will be investigated.

The paper has the following structure. Except for the Introduction (Sec. I)
and the Conclusion and Outlooks (Sec. V), Appendix A and Sec. II give a quick
review of the effective potentials, as well as some comments and remarks on
the electroweak phase transition in the 2HDM. In Sec. III and Sec. IV, the
electroweak phase transition in the 2HDM-$S_{3}$ is studied. More
specifically, the effective potential of the dual electroweak phase transition
will be studied, the strength of first-order phase transition, the mass domain
for the first-order phase transition, all of which will be given with a
parameter $a$ that will be introduced to replace the popular parameter
$\tan\beta=v_{2}/v_{1}$.

##  II Review on the Higgs potential and comments on EWPT in the 2HDM

###  II.1 The Higgs potential in the 2HDM

The fermion and scalar spectrum with their assignments under the
$SU(2)_{L}\times U(1)_{Y}$ gauge group are given by davidson ; Herrero-
Garcia:2017xdu :

$\displaystyle L_{a}$ $\displaystyle=\begin{pmatrix}\nu_{aL}\\\
e_{aL}\end{pmatrix}\sim(2,-1)\,,\;\quad Q_{aL}=\begin{pmatrix}u_{aL}\\\
d_{aL}\end{pmatrix}\sim\left(2,\dfrac{1}{3}\right)\,,$ $\displaystyle e_{aR}$
$\displaystyle\sim(1,-2),\quad u_{aR}\sim\left(1,\dfrac{4}{3}\right)\,,\quad
d_{aR}\sim\left(1,-\dfrac{2}{3}\right),\;a=1,2,3\,,$
$\displaystyle\Phi_{i}=\begin{pmatrix}\phi^{+}_{i}\\\
\phi^{0}_{i}\end{pmatrix}\sim(2,1)\,,\;\langle\Phi_{i}\rangle=\frac{1}{\sqrt{2}}\begin{pmatrix}0\\\
v_{i}\end{pmatrix}\sim(2,1),i=1,2\,.$ (1)

For details of the quark sector of different types of the 2HDM, the reader is
referred to Refs. davidson ; Herrero-Garcia:2017xdu . In the 2HDM, the charge
operator is defined as

$Q=T_{3}+\frac{Y}{2}\,.$ (2)

The most general form of the effective potential in the 2HDM contains $14$
parameters, and there may exist the CP conservation (charge and mirror
symmetry), CP violation, and charge violation. When expressing the form of the
potential, we must be careful in defining the quantities and distinguishing
the parameters, since when applying the group rotations, the physics can be
changed. However, in the studies of the phenomenology of the 2HDM, assumptions
are often made to simplify the calculations. For example, CP is often assumed
to be conserved in the Higgs fields (only then can we distinguish the scalar
field from the scalar pseudo-field). Now, the CP discrete symmetry will
eliminate all the fourth power terms that contain the odd number of one of two
Higgs fields from the potential (for instance,
$\Phi_{1}^{\dagger}\Phi_{1}\Phi_{1}^{\dagger}\Phi_{2}$). We can also assume
that all of the parameters corresponding to the fourth power terms are of real
number, including the term added to break the symmetries.

With the above assumptions, putting in the symmetry $Z_{2}$
($\Phi_{1}\rightarrow\Phi_{1}$, $\Phi_{2}\rightarrow-\Phi_{2}$) in order not
to have FCNC at the tree level. After this, the most general form for the
scalar potential for the two doublets $\Phi_{1}$ and $\Phi_{2}$ with
supercharge $+1$ would have the following form tdhm :

$\displaystyle\begin{split}V=&\,\,\mu_{11}^{2}\Phi_{1}^{2}+\mu_{22}^{2}\Phi_{2}^{2}-\mu_{12}^{2}\left(\Phi_{1}.\Phi_{2}+\Phi_{2}.\Phi_{1}\right)+\dfrac{\lambda_{1}}{2}\Phi_{1}^{4}+\dfrac{\lambda_{2}}{2}\Phi_{2}^{4}\\\
&+\lambda_{3}\Phi_{1}^{2}\Phi_{2}^{2}+\lambda_{4}(\Phi_{1}.\Phi_{2})(\Phi_{2}.\Phi_{1})+\dfrac{\lambda_{5}}{2}\left(\Phi_{1}^{4}+\Phi_{2}^{4}\right),\end{split}$
(3)

where we have denoted $\Phi_{i}^{2}=\Phi_{i}^{\dagger}\Phi_{i},\,=1,2$ and
$\Phi_{i}.\Phi_{j}=\Phi_{i}^{\dagger}\Phi_{j},\,i\neq j$. All the parameters
$\mu_{11}$, $\mu_{22}$, $\lambda_{i}$ $(i=1,\dots,5)$ are all real and the
term containing $\mu_{12}$ “softly” breaks the symmetry $Z_{2}$. Following
this, there are two complex scalar doublets $SU(2)$ containing eight fields,

$\displaystyle\Phi_{a}=\begin{pmatrix}\Delta_{a}^{+}\\\
\frac{1}{\sqrt{2}}(v_{a}+\rho_{a}+i\eta_{a})\end{pmatrix},\qquad a=1,2.$ (4)

Three of them are Goldstone bosons eaten by the massive gauge bosons $W^{\pm}$
and $Z^{0}$ to generate their masses. The other five are physical scalar Higgs
fields, including one scalar carrying charge, two neutral scalars, and one
pseudo scalar.

Averaging over the whole space, VEVs read

$\displaystyle\langle\Phi_{1}\rangle_{0}=\frac{1}{\sqrt{2}}\begin{pmatrix}0\\\
v_{1}\end{pmatrix},\langle\Phi_{2}\rangle_{0}=\frac{1}{\sqrt{2}}\begin{pmatrix}0\\\
v_{2}\end{pmatrix},$ (5)

One of the most important parameters of the model is as follows:

$\displaystyle\tan\beta=\frac{s_{\beta}}{c_{\beta}}=\dfrac{v_{2}}{v_{1}}\,,$
(6)

where we have used the notations
$s_{\beta}\equiv\sin\beta,\,c_{\beta}\equiv\cos\beta$. Here $\beta$ is the
rotational angle when normalizing the squared mass matrix of the charged and
pseudoscalars.

Two general neutral even CP states $\rho_{1}$ and $\rho_{2}$ are not physical
mass states. The mass matrix corresponding to them can be diagonalized by a
rotation of a mixed angle of $\rho_{1}$ and $\rho_{2}$.

Having two minima as above inserted into the Lagrangian which contains the
mass terms corresponding to the scalar fields. First, the term that
corresponds to the charged scalar field is

$\displaystyle\mathcal{L}_{\phi^{\pm},mass}=\left[\mu_{12}^{2}-(\lambda_{4}+\lambda_{5})v_{1}v_{2}\right]\begin{pmatrix}\Delta_{1}^{-}&\Delta_{2}^{-}\end{pmatrix}\begin{pmatrix}\tan\beta&-1\\\
-1&\coth\beta\end{pmatrix}\begin{pmatrix}\Delta_{1}^{+}\\\
\Delta_{2}^{+}\end{pmatrix}.$ (7)

After diagonalizing the above mass matrix, the squared mass of the charged
Higgs particle,

$\displaystyle
m_{H^{\pm}}^{2}=\left(\dfrac{\mu_{12}^{2}}{v_{1}v_{2}}-\lambda_{4}-\lambda_{5}\right)(v_{1}^{2}+v_{2}^{2})=\left(\dfrac{\mu_{12}^{2}}{v_{1}v_{2}}-\lambda_{4}-\lambda_{5}\right)v^{2}.$
(8)

Next, the mass term that corresponds to the pseudoscalar field,

$\displaystyle\mathcal{L}_{\eta,mass}=\left(\dfrac{\mu_{12}^{2}}{v_{1}v_{2}}-\lambda_{5}\right)\begin{pmatrix}\eta_{1}&\eta_{2}\end{pmatrix}\begin{pmatrix}v_{2}^{2}&-v_{1}v_{2}\\\
-v_{1}v_{2}&v_{1}^{2}\end{pmatrix}\begin{pmatrix}\eta_{1}\\\
\eta_{2}\end{pmatrix}.$ (9)

The physical squared mass of the pseudoscalar after normalization,

$\displaystyle
m_{A}^{2}=\left(\dfrac{\mu_{12}^{2}}{v_{1}v_{2}}-2\lambda_{5}\right)(v_{1}^{2}+v_{2}^{2}).$
(10)

Finally, the mass term that corresponds to the remaining two scalar fields is

$\displaystyle\mathcal{L}_{\rho,mass}=-\begin{pmatrix}\rho_{1}&\rho_{2}\end{pmatrix}\begin{pmatrix}\mu_{12}^{2}\frac{v_{2}}{v_{1}}+\lambda_{1}v_{1}^{2}&-\mu_{12}^{2}+\lambda_{345}v_{1}v_{2}\\\
-\mu_{12}^{2}+\lambda_{345}v_{1}v_{2}&\mu_{12}^{2}\frac{v_{1}}{v_{2}}+\lambda_{2}v_{2}^{2}\end{pmatrix}\begin{pmatrix}\rho_{1}\\\
\rho_{2}\end{pmatrix}.$ (11)

With $\lambda_{345}=\lambda_{3}+\lambda_{4}+\lambda_{5}$, after normalizing
the matrix in Eq.(11), the squared masses of the light $(h)$ and the heavy
$(H)$ Higgs particles are, respectively:

$\displaystyle\begin{split}&m_{h}^{2}=\dfrac{1}{2}\left[(\lambda_{1}v_{1}^{2}+\lambda_{2}v_{2}^{2})+\mu_{12}^{2}\dfrac{v^{2}}{v_{1}v_{2}}\right]\\\
&\qquad-\sqrt{\left[\dfrac{\lambda_{1}v_{1}^{2}-\lambda_{2}v_{2}^{2}}{2}-\dfrac{\mu_{12}^{2}}{2v_{1}v_{2}}(v_{1}^{2}-v_{2}^{2})\right]^{2}+(\lambda_{345}v_{1}v_{2}-\mu_{12}^{2})^{2}},\\\
&m_{H}^{2}=\dfrac{1}{2}\left[(\lambda_{1}v_{1}^{2}+\lambda_{2}v_{2}^{2})+\mu_{12}^{2}\dfrac{v^{2}}{v_{1}v_{2}}\right]\\\
&\qquad+\sqrt{\left[\dfrac{\lambda_{1}v_{1}^{2}-\lambda_{2}v_{2}^{2}}{2}-\dfrac{\mu_{12}^{2}}{2v_{1}v_{2}}(v_{1}^{2}-v_{2}^{2})\right]^{2}+(\lambda_{345}v_{1}v_{2}-\mu_{12}^{2})^{2}},\end{split}$
(12)

From the Eqs. (8), (10), (12), the squared masses of the Higgs particles
always contain in them troublesome mixing terms of VEVs. A quantity is
introduced as follows:

$M^{2}_{Higgs}(v_{1},v_{2})=m^{2}_{H}H^{2}+m^{2}_{h}h^{2}+m^{2}_{A}A^{2}+m^{2}_{H^{\pm}}(H^{\pm})^{+}H^{\pm}.$
(13)

###  II.2 Higgs potential at the tree level

From the Higgs potential given in Eq. (3), $V_{0}$ in a form that is dependent
on the VEVs as follows:

$V_{0}(v_{1},v_{2})=\frac{\mu_{11}^{2}}{2}v_{1}^{2}+\frac{\mu_{22}^{2}}{2}v^{2}_{2}-\frac{\mu_{12}^{2}}{2}v_{1}v_{2}+\left(\dfrac{\lambda_{1}}{8}+\dfrac{\lambda_{5}}{8}\right)v_{1}^{4}+\left(\dfrac{\lambda_{2}}{8}+\dfrac{\lambda_{5}}{8}\right)v_{2}^{4}+\left(\frac{\lambda_{3}}{4}+\frac{\lambda_{4}}{4}\right)v_{2}^{2}v_{1}^{2}.$
(14)

$V_{0}(v_{1},v_{2})$ has a quartic form like in the SM. On the other hand, by
developing the Higgs potential Eq. (3), two minimum equations which permit us
to transform the mixing between $v_{1}$ and $v_{2}$,

$\displaystyle\frac{\partial V}{\partial v_{1}}$ $\displaystyle=$
$\displaystyle
0,\mu_{11}^{2}v_{1}-\mu_{12}^{2}v_{2}+\left(\lambda_{1}+\lambda_{5}\right)\frac{v^{3}_{1}}{2}+\left(\lambda_{3}+\lambda_{4}\right)\frac{v_{2}^{2}v_{1}}{2}=0.$
$\displaystyle\frac{\partial V}{\partial v_{2}}$ $\displaystyle=$
$\displaystyle
0,\mu_{22}^{2}v_{2}-\mu_{12}^{2}v_{1}+\left(\lambda_{2}+\lambda_{5}\right)\frac{v^{3}_{2}}{2}+\left(\lambda_{3}+\lambda_{4}\right)\frac{v_{2}v_{1}^{2}}{2}=0.$
(15)

From Eq. (15), a relationship between VEVs, such as

$\mu_{22}^{2}v^{2}_{2}+\left(\lambda_{2}+\lambda_{5}\right)\frac{v^{4}_{2}}{2}=\mu_{12}^{2}v_{1}v_{2}-\left(\lambda_{3}+\lambda_{4}\right)\frac{v_{2}^{2}v_{1}^{2}}{2}.$
(16)

Substituting Eqs. (16) into Eq. (14) yields

$V_{0}(v_{1},v_{2})=\frac{\mu_{11}^{2}}{2}v_{1}^{2}+\left(\dfrac{\lambda_{1}}{8}+\dfrac{\lambda_{5}}{8}\right)v_{1}^{4}-\left(\dfrac{\lambda_{2}}{8}+\dfrac{\lambda_{5}}{8}\right)v_{2}^{4}=V_{0}(v_{1})+V_{0}(v_{2}),$
(17)

where
$V_{0}(v_{1})=\frac{\mu_{11}^{2}}{2}v_{1}^{2}+\left(\dfrac{\lambda_{1}}{8}+\dfrac{\lambda_{5}}{8}\right)v_{1}^{4}$
and
$V_{0}(v_{2})=\left(-\dfrac{\lambda_{2}}{8}-\dfrac{\lambda_{5}}{8}\right)v_{2}^{4}$
are in the quartic form. In addition, there are alternative ways to arrive Eq.
(17) which has other forms but $V_{0}(v_{1})$ and $V_{0}(v_{2})$ are still in
the quartic form.

If the potential at the tree level had the quartic form of each vacuum
expectation value, that is there are no mixing terms in it, the job of
calculating the effective potential for each VEV would be much easier. This
will be made clear in the following sections.

###  II.3 The masses of gauge bosons

In order to find the gauge boson masses, we starting from the kinetic term of
the Higgs fields. In the 2HDM-$S_{3}$, there are two components in the kinetic
term for the two Higgs doublets,

$\mathcal{L}^{GB}_{mass}=\left(\mathcal{D}_{\mu}\langle\Phi_{1}\rangle\right)^{\dagger}\left(\mathcal{D}^{\mu}\langle\Phi_{1}\rangle\right)+\left(\mathcal{D}_{\mu}\langle\Phi_{2}\rangle\right)^{\dagger}\left(\mathcal{D}^{\mu}\langle\Phi_{2}\rangle\right)\equiv
A+B,$ (18)

in which, the covariant derivatives act on $\Phi_{1}$ and $\Phi_{2}$. as
follows

$\mathcal{D}^{\mu}=\partial^{\mu}-igW^{\mu}_{i}T_{i}-\frac{i}{2}g^{\prime}B^{\mu}Y\,.$
(19)

Note that the gauge fields ($W^{\mu}_{i},B^{\mu}$) inside the covariant
derivatives of $A$ and $B$ are the same. So after diagonalizing, the gauge
fields in $A$ and $B$ are the same, and gauge bosons $\gamma$, $Z$, $W^{\pm}$
are inferred.

From the term $A$, one obtains the mass components of the physical gauge
bosons only depends on $v_{1}$,

$\begin{split}M^{A}_{bosons}=m^{2}_{W^{\pm}}(v_{1})W^{+}_{\mu}W^{-\mu}+m^{2}_{Z}(v_{1})Z_{\mu}Z^{\mu}.\end{split}$
(20)

From the term $B$, one obtains the mass components of the physical gauge
bosons that only depend on $v_{2}$,

$\begin{split}M^{B}_{bosons}=m^{2}_{W^{\pm}}(v_{2})W^{+}_{\mu}W^{-\mu}+m^{2}_{Z}(v_{2})Z_{\mu}Z^{\mu}.\end{split}$
(21)

Therefore the bosons masses can be split into two parts,

$m^{2}_{gauge-boson}(v_{1},v_{2})=m^{2}_{boson}(v_{1})+m^{2}_{boson}(v_{2}).$
(22)

Similar arguments for the problem can be found in Ref. phonglongvan2 .

###  II.4 Remarks on EWPT structure in the 2HDM

The Higgs and gauge boson sectors from the full Higgs Lagrangian,

$\mathcal{L}=\mathcal{L}^{GB}_{mass}-V(\Phi_{1},\Phi_{2}),$ (23)

where $V(\Phi_{1},\Phi_{2})$ is given by Eq. (3).

Expanding the Higgs fields $\Phi_{1}$ and $\Phi_{2}$ around their VEVs which
are $v_{1}$, $v_{2}$, yields

$\displaystyle\begin{split}\mathcal{L}=&\frac{1}{2}\partial^{\mu}v_{1}\partial_{\mu}v_{1}+\frac{1}{2}\partial^{\mu}v_{2}\partial_{\mu}v_{2}-V_{0}(v_{1},v_{2})+M^{A}_{gauge-
bosons}+M^{B}_{gauge-bosons}\\\ &+M^{2}_{Higgs}(v_{1},v_{2})+[m_{top-
quark}(v_{1})+m_{top-quark}(v_{2})]t\bar{t}.\end{split}$ (24)

Therefore, from the Lagrangian in Eq. (24), two motion equations according to
$v_{1}$ and $v_{2}$ are calculated,

$\displaystyle\partial^{\mu}v_{1}\partial_{\mu}v_{1}-\frac{\partial
V_{0}(v_{1})}{\partial v_{1}}+\sum\frac{\partial
m^{2}_{bosons}(v_{1})}{\partial v_{1}}W^{\mu}W_{\mu}+\frac{\partial
m_{top}(v_{1})}{\partial v_{1}}t\bar{t}+\frac{\partial
M^{2}_{Higgs}(v_{1},v_{2})}{\partial v_{1}}$ $\displaystyle=$ $\displaystyle
0,$ (25) $\displaystyle\partial^{\mu}v_{2}\partial_{\mu}v_{2}-\frac{\partial
V_{0}(v_{2})}{\partial v_{2}}+\sum\frac{\partial
m^{2}_{bosons}(v_{2})}{\partial v_{2}}W^{\mu}W_{\mu}+\frac{\partial
m_{top}(v_{2})}{\partial v_{2}}t\bar{t}+\frac{\partial
M^{2}_{Higgs}(v_{1},v_{2})}{\partial v_{2}}$ $\displaystyle=$ $\displaystyle
0,$ (26)

where $W$ runs over all gauge fields.

If $M^{2}_{Higgs}(v_{1},v_{2})$ is like Eq. (13), in which there are no mixing
terms of VEVs, this term can be separated into two terms such that each of the
new terms only depends on one VEV. However, the fifth Higgs particles’ masses
in the 2HDM all have mixing terms.

Next, there are an important observation, that when the universe was cooling
down to the value of $v_{2}$ after the big bang, the field $\Phi_{2}$ broke
the electroweak symmetry, and after that when the universe continued to cool
down to the value of $v_{1}$, the field $\Phi_{1}$ continued to break the
electroweak symmetry once again. The process of this electroweak symmetry
breaking must be sequential. Hence, we cannot combine the two stages to study.

Therefore the rules for generating the masses of the particles through two
symmetry breakings as follows:

Remark 1: At stage 1, $\Phi_{2}$ breaks the symmetry or $v_{2}\neq 0$, but now
$\Phi_{1}$ has not yet broken the symmetry so $v_{1}$ is still equal to 0.
Hence, all the Higgs particles’ masses only contain $v_{2}$.

Remark 2: When the breaking symmetry occurs at $\Phi_{1}$, the interactions
between $\Phi_{2}$ and $\Phi_{1}$ are turned on and $v_{1}$ would not be 0. In
this stage, the further generated masses can not only depend on $v_{2}$, but
also on $v_{1}$ and the mixing of $v_{2}$ and $v_{1}$.

Remark 3: With the above two remarks, through the mixing terms in the masses
of the Higgs particles, the effects of the first stage has on the second
stage. But they also make it difficult for investigating the phase transition
at the later stage.

Remark 4: In order to view the two phase transition stages with the separated
effective potentials, we can apply the following approximation to the mixing
terms: $v_{1}.v_{2}\sim\kappa.v_{1}^{2}$. Since at this time, $v_{2}$ can
still change as the temperature decreases, it can consider the change of
$v_{2}$ is now equal to $\kappa v_{1}$. This remark is actually a mathematical
treatment like the approximation $v^{2}=v_{1}^{2}+v^{2}_{2}$. But when it is
combined with the above third remarks, they make physical sense in the
analysis of EWPT.

Therefore, from four remarks, all the Higgs particles’ masses can be split
into two different components,

$M^{2}_{Higgs}(v_{1},v_{2})=m^{2}_{Higgs}(v_{1})+m^{2}_{Higgs}(v_{2}).$ (27)

Also the squared masses of the gauge and Higgs particles all can be split into
two separate components at the tree level. From Eqs. (25) and (26), averaging
over space and using Bose-Einstein and Fermi-Dirac distributions respectively
for bosons and fermions to average over space, the one-loop effective
potential can be obtained at high temperatures. Also according to the analysis
of Appendix A, the analysis of the Lagrangian of 2HDM into two separate
components (as shown in Secs. II.1, II.2, II.3), the total effective potential
in the 2HDM model can be rewritten as

$V^{2HDM}_{eff}=V_{eff}^{2HDM}(v_{1})+V_{eff}^{2HDM}(v_{2}).$ (28)

For further clarity, we restate the calculating process of the effective
potential from the contributions of one-loop diagrams. The process of
calculating the one-loop effective potential is the process of calculating the
contribution of 1-loop diagrams with $n$ external lines that are Higgs scalar
fields (fields that act as mass generators). In the 2HDM, there are two Higgs
fields ($h,H$) that act as such, corresponding to two VEVs ($v_{1},v_{2}$).
One-loop diagrams are shown in Figs. 1, 2, 3, 4.

$\varphi_{1}$$\varphi_{1}$

\+ $\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$ \+
$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$+…

Figure 1: The 1-loop contributions of the scalar fields

$\varphi_{1}$$\varphi_{1}$

\+ $\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$ \+
$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$+…

Figure 2: The 1-loop contributions of the fermion fields

$\varphi_{1}$$\varphi_{1}$

\+ $\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$ \+
$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$
+…

Figure 3: The 1-loop contributions of the gauge fields

Similarly, also having diagrams where the external lines are $\varphi_{2}$.
The mixing $\varphi_{1}-\varphi_{2}$ diagrams as Fig. 4.

$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{2}$$\varphi_{2}$$\varphi_{2}$$\varphi_{2}$

\+
$\varphi_{1}$$\varphi_{2}$$\varphi_{2}$$\varphi_{1}$$\varphi_{2}$$\varphi_{1}$\+
$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{1}$$\varphi_{2}$$\varphi_{2}$$\varphi_{2}$$\varphi_{2}$
+…

Figure 4: The diagrams of one-loop contribution when there is a mix of
$v_{1}-v_{2}$

In the above diagrams, $\varphi_{1}$ and $\varphi_{2}$ correspond to the terms
of Higgs fields that only contain $v_{1}$ or $v_{2}$. Since the 2HDM consists
of two VEVs, having the mixing diagrams as in Fig. 4.

For the $\phi^{4}$ theory, calculating the contributions of the diagrams as in
Figs. 1,2 and 3 is really easy, and this has been shown in Ref. quiros .
However, the difficulty is to calculate the contributions of the diagrams in
Fig. 4. As in Appendix A, we explicitly study the first diagram in Fig. 4 that
translates to the following expression:

$\displaystyle\Gamma_{n}(v_{1},v_{2})=$ $\displaystyle
i\frac{1}{2n}\int\frac{d^{4}p}{(2\pi)^{4}}\left(\frac{i}{p^{2}-m^{2}(v_{1},v_{2})+i\varepsilon}\right)^{n}(-i2\lambda_{i})^{n_{1}/2}v_{1}^{n_{1}}(-i2\lambda_{2})^{n_{2}/2}v_{2}^{n_{2}},$
(29)

with $n_{1}+n_{2}=2n$. Summing over all $\Gamma_{n}$ with $n$ runs from $0$ to
infinity, and hence $n_{1}$ and $n_{2}$ will also run from $0$ to infinity.
Hence calculating the integrals $\Gamma_{n}$ and the infinite sums are very
tricky. Since for each value of $n$, there are a sum that runs with $n_{1}$ or
$n_{2}$.

To quickly calculate $\Gamma_{n}$, from remarks 2 and 4,
$\Gamma_{n}(v_{1},v_{2})\equiv\Gamma_{n}(v_{1})$,

$\displaystyle\Gamma_{n}(v_{1},v_{2})=$ $\displaystyle
i\frac{1}{2n}\int\frac{d^{4}p}{(2\pi)^{4}}\left(\frac{i}{p^{2}-m^{2}(v_{1})+i\varepsilon}\right)^{n}v_{1}^{n_{1}}v_{1}^{n_{2}}\kappa^{n_{2}}(-i2\lambda_{1})^{n_{1}/2}(-i2\lambda_{2})^{n_{2}/2}$
(30)
$\displaystyle=i\frac{1}{2n}\int\frac{d^{4}p}{(2\pi)^{4}}\left(-i2\lambda^{\prime}_{n}\frac{i}{p^{2}-m^{2}(v_{1})+i\varepsilon}\right)^{n}v_{1}^{2n}.$
(31)

Here, $\lambda^{\prime}_{n}$ must have some very small values in order for the
series to converge. Therefore $\lambda^{\prime}_{n}\sim\lambda^{\prime}$.
Hence, $\Gamma_{n}(v_{1},v_{2})\equiv\Gamma_{n}(v_{1})$,

$\displaystyle\Gamma_{n}(v_{1})=i\frac{1}{2n}\int\frac{d^{4}p}{(2\pi)^{4}}\left(2\lambda^{\prime}\frac{v^{2}_{1}}{p^{2}-m^{2}(v_{1})+i\varepsilon}\right)^{n}.$
(32)

Therefore the contributions from these mixing diagrams can be combined into
the contribution from the diagrams that the external lines are just all
$\varphi_{1}$. In other words, by calculating the contributions from the
diagrams, and applying Remarks 2, and 4, the mixing diagrams can be processed
to turn the effective potential into two clearly separated components, one
depends only on $v_{1}$, while the other depends only on $v_{2}$. In other
words, the effective potential will be expressed as Eq. (28).

In the previous studies, to overcome the difficulties of dealing with the
mixing terms of VEVs and to investigate the process of electroweak symmetry
breaking similar to the SM, the authors have converted $v_{1}$ and $v_{2}$ to
$v$ through $\tan\beta$ tdhm ; Fuyutob . This is a very clever technique, but
in terms of physics, it needs to be interpreted with care. Since $v_{1}\neq
v_{2}$, the symmetry-breaking stages must be sequential. The fact that we
write the same effective potential for $v$ (shortly denoted as $V_{eff}(v)$)
to calculate the strength of phase transition is not wrong, but we should only
write it when $v_{1}\ll v_{2}$, to clearly show the nature of the physics in
the 2HDM. The limitations and utilities of this technique will be analyzed in
the following sections where investigating the 2HDM-$S_{3}$.

###  II.5 Comments on EWPT in the 2HDM

First, we have summarized of experimental as well as theoretical calculations
leading to parameter regions in 2HDM models:

  * •

Since the 125 GeV Higgs boson observed at the LHC, the model becomes
consistent with the LHC Higgs data when the model provides such a Higgs
particle bpal1 . $1<\tan\beta$ can be precisely determined from the
requirement of the light mass of the up- and down-quarks bpal1 . From here
also admit the scenario $\mu_{1},\mu_{2},l_{1},l_{2}\neq 0$ and having the
soft-breaking potential. Also according to Ref. bpal2 , the mass of the
nonstandard particles are less than 1 TeV, $0.3<\tan\beta<17$.

  * •

Some references used the data from decay channels in LHC and investigated the
value range of $\tan\beta$. In Ref. beta1 , the authors removed the
circumstance $\tan\beta<1$ for four types of 2HDM. In Ref. beta2 , the authors
saw that in the decay channel $pp\rightarrow W\rightarrow hH^{\pm}$, the
scattering amplitude of that channel is lasting over the change of $\tan\beta$
in the range of $1<\tan\beta<10$.

  * •

The FCNC structure exists in the model. It is therefore compatible with
current experimental data on quarks bpal1 ; 36t . The mass of exotic particle
below 190 or 300 GeV has been already excluded by the data from the LHC Run-II
and the HL-LHC, the most of the parameter region would be explored mayumi .
Since having a significant amount of the $b\bar{b}$ branching ratios for the
additional Higgs bosons mayumi .

  * •

The model under consideration contains the pseudoscalar field $A$ being
attached subject for recent experimental study CMS2 .

  * •

Next, in general, there are four types of 2HDM models of $Z_{2}$ symmetry. Two
of the four types were investigated in Ref. dori  and also with parameter
domains consistent with the above conclusions.

The above conclusions lead to an instruction for surveying 2HDM-$S_{3}$ also
with a parameter region, $1\sim\tan\beta<17$ and the masses of additional
bosons must be larger than 200 GeV.

More importantly, the following comments about the EWPT in 2HDM models after
the observation of the 125 GeV Higgs at LHC (2012 are in order):

  * •

According to Ref. 1305.6610 , with LHC data and decay channel
$h^{0}\longrightarrow\gamma\gamma$, for a strong first order EWPT in 2HDM,
$m_{A}>400$ GeV, a mass hierarchy $m_{H^{\pm}}<m_{H}<m_{A}$ and
$1<\tan\beta<10$. This does not define the upper bounds of the masses of the
particles, but the lower bounds are about 400 GeV.

  * •

In Ref. 1504.05949 , when analyzing the inert 2HDM model, for $500\text{
GeV}<m_{H}<afew$ TeV, for a first order EWPT, boson $H$ could be a candidate
for dark matter.

  * •

The authors in Ref. PRL  performed a nonperturbation study of EWPT in 2HDM. To
have a first order EWPT, the condition is $m_{A}>m_{H}+m_{Z}$. This suggests
that the mass of additional Higgs bosons must be larger mass than one of SM-
like Higgs boson.

  * •

In Ref. 2111.13079J. High Energy Phys , for a first order EWPT and combined
with LHC data, the masses of additional Higgs bosons are typically $300-400$
GeV. The triple Higgs boson coupling is predicted to be $35-55$ % larger than
the standard model value.

  * •

In Ref. dori , analyzing decay channels $A,H,H^{\pm}\longrightarrow tt,tb$,
combining with HL-LHC signal and gravitational wave observations at LISA. The
2HDM model for a first order EWPT. It also shows that these decay channels can
be key channels to authenticate the first order EWPT in 2HDM.

  * •

In particular, in Ref. dori , Fig. 1 shows that the ratio between two phase
transition strengths of the two model types is almost independent of
$\tan\beta$. Different scenarios between the masses of the additional Higgs
particles were analyzed in the EWPT problem, such as $m_{H^{\pm}}=m_{A}$ or
$m_{H}=m_{H^{\pm}}$.

The summaries of the EWPT results in the 2HDM model are important indicators
for the parameter space in the calculation of EWPT in the 2HDM-$S_{3}$. It
also shows the effect of the $S_{3}$ symmetry, which will be discussed in the
following sections.

##  III Review on the 2HDM-$S_{3}$

###  III.1 Particle content

To solve with FCNC for the 2HDM, ones can realize by implement of $S_{3}$
symmetry S32dM . The particle contents and their charge assignment are given
in Table 1.

Particle | $Q_{a}$ | $L_{a}$ | $L_{\tau}$ | $u_{aR}$ | $d_{iR}$ | $e_{aR}$ | $\tau_{R}$ | $\Phi_{1}$ | $\Phi_{2}$
---|---|---|---|---|---|---|---|---|---
$SU(2)_{L}$ | $\bm{2}$ | $\bm{2}$ | $\bm{2}$ | $\bm{1}$ | $1$ | $\bm{1}$ | $\bm{1}$ | $\bm{2}$ | $\bm{2}$
$U(1)_{Y}$ | $\frac{1}{6}$ | $-\frac{1}{2}$ | $-\frac{1}{2}$ | $\frac{2}{3}$ | $-\frac{1}{3}$ | $-1$ | $-1$ | $\frac{1}{2}$ | $\frac{1}{2}$
$\mathbb{S}_{3}$ | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{1}$ | $\mathbf{1^{\prime}}$ |  $\mathbf{1}$ or $\mathbf{1^{\prime}}$ | $\mathbf{2}$ |  $\mathbf{1}$ or $\mathbf{1^{\prime}}$ | $\mathbf{1}$ | $\mathbf{1^{\prime}}$
Table 1: The particle contents and their charge assignment of the
$SU(2)_{L}\times U(1)_{Y}\times S_{3}$ symmetry.

There are two kinds of representations for $S_{3}$: real and complex; and it
is easier to work with complex representation ma .

###  III.2 Higgs potential

The generic scalar potential of 2HDM-$S_{3}$ 1601  can be written as

$\displaystyle\begin{split}V_{H}=\quad&m_{11}^{2}\Phi_{1}^{2}+m_{22}^{2}\Phi_{2}^{2}-m_{12}^{2}\Phi_{1}.\Phi_{2}-(m_{12}^{2})^{*}\Phi_{2}.\Phi_{1}\\\
&+\dfrac{\lambda_{1}}{2}\Phi_{1}^{4}+\dfrac{\lambda_{2}}{2}\Phi_{2}^{4}+\lambda_{3}\Phi_{1}^{2}\Phi_{2}^{2}+\lambda_{4}(\Phi_{1}.\Phi_{2})(\Phi_{2}.\Phi_{1})\\\
&+\left[\dfrac{\lambda_{5}}{2}(\Phi_{1}.\Phi_{2})^{2}+(\lambda_{6}\Phi_{1}^{2}+\lambda_{7}\Phi_{2}^{2})(\Phi_{1}.\Phi_{2})+\textrm{H.
c.}\right],\end{split}$ (33)

where any couplings other than $m^{2}_{12}$, $\lambda_{5}$, $\lambda_{6}$, and
$\lambda_{7}$ are real. Using $v_{1}/\sqrt{2}$ and $v_{2}/\sqrt{2}$ to denote
the VEVs, but ignore the factor $\sqrt{2}$ (however, it was still included in
the calculations).

If choosing $\Phi\sim(\mathbf{1},\mathbf{1^{\prime}})$ to be the
representations of $S_{3}$, then all the odd terms only containing $\Phi_{2}$,
such as $m_{12}^{2}\Phi_{1}.\Phi_{2}$, $(m_{12}^{2})^{*}\Phi_{2}.\Phi_{1}$ and
$(\lambda_{6}\Phi_{1}^{2}+\lambda_{7}\Phi_{2}^{2})(\Phi_{1}.\Phi_{2})$ must be
terminated for the Lagrangian to be invariant under the $S_{3}$ group
transformation, since the representation $\mathbf{1^{\prime}}$ changes the
sign of the fields with odd permutation. In this case, getting a $Z_{2}$
symmetric potential, in which $m^{2}_{12}=\lambda_{6}=\lambda_{7}=0$.

In general, we can assume that $\Phi\sim s$, where $s$ denotes either of the
two alternatives $(\mathbf{1},\mathbf{1})$ or
$(\mathbf{1},\mathbf{1^{\prime}})$. Because the labels 1 and 2 were selected
at random, the case $\Phi\sim(\mathbf{1^{\prime}},\mathbf{1})$ is also
included. For its simplicity, a complicated representation can be selected to
work with. The most common Higgs potential of a $S_{3}$ doublet is of the
following form where the two scalars $\Phi=(\phi_{1},\phi_{2})^{T}$ transform
as a doublet in a complex representation 1601 :

$\displaystyle\begin{split}V_{C}=\quad&\mu_{1}^{2}(\phi_{2}^{2}+\phi_{1}^{2})+\frac{1}{2}l_{1}(\phi_{2}^{2}+\phi_{1}^{2})^{2}+\frac{1}{2}l_{2}(\phi_{2}^{2}-\phi_{1}^{2})^{2}\\\
&+l_{3}(\phi_{1}.\phi_{2})(\phi_{2}.\phi_{1}).\end{split}$ (34)

Note that both complex and real representations provide the similar result
given in (34) which coincides with generic one $V_{H}$ in (33) for the
conditions 1601

$m^{2}_{11}=m^{2}_{22}\equiv\mu^{2}_{1}\,,\,m^{2}_{12}=0\,,\lambda_{1}=\lambda_{2}\,,\,\lambda_{5}=0\,.$
(35)

###  III.3 The soft breaking of $S_{3}$ group

To break $S_{3}$ softly, Ma and Melic melic  include a soft potential by hand
to the full one, while still preserving the $\phi_{1}\leftrightarrow\phi_{2}$
symmetry,

$\displaystyle
V_{\text{soft}}=-\mu_{2}^{2}\left(\phi_{1}^{\dagger}\phi_{2}+\phi_{2}^{\dagger}\phi_{1}\right).$
(36)

Soft breaking terms here mean that they violate the original symmetry in the
Higgs potential or the Lagrangian. They are ’soft’ because the couplings
associated with those terms are small. Without the above mentioned term, we
will face the trouble of massless pseudoscalar $A$. Realistically, there
should be some breaking terms which will take care of this problem.

In a spontaneously broken case, you break the symmetry of the ground state and
it naturally breaks the symmetry in the Lagrangian. On the other hand, in this
case, those terms must be added by hand to break the symmetry. The potential
then becomes $V_{C}+V_{\text{soft}}$.

The term containing $V_{soft}$ makes the mass of the pseudoscalar Higgs
particle $A$ always nonzero in all cases. If there is no $V_{soft}$ term, the
Higgs particle A will have a mass of 0, when $v_{1}=v_{2}$. Besides, the term
$V_{soft}$ must exist, since it represents the direct interaction between
$\phi_{1}$ and $\phi_{2}$ Therefore, we must study the Higgs potential that
contains the term $V_{soft}$.

##  IV Electroweak phase transition in the 2HDM-$S_{3}$

###  IV.1 A vital role of $S_{3}$

First, the function of the $S_{3}$ group in the 2HDM-$S_{3}$ can be shown by
comparing the Higgs potential of the two models before and after adding
$S_{3}$ symmetry. Let us consider the Higgs potential of 2HDM which has the
following form

$\displaystyle V_{\text{Higgs}}^{\text{2HDM}}$
$\displaystyle=m_{11}^{2}\Phi_{1}^{\dagger}\Phi_{1}+m_{22}^{2}\Phi_{2}^{\dagger}\Phi_{2}-m_{12}^{2}\left(\Phi_{1}^{\dagger}\Phi_{2}+\Phi_{2}^{\dagger}\Phi_{1}\right)+\dfrac{\lambda_{1}}{2}\left(\Phi_{1}^{\dagger}\Phi_{1}\right)^{2}$
(37)
$\displaystyle+\dfrac{\lambda_{2}}{2}\left(\Phi_{2}^{\dagger}\Phi_{2}\right)^{2}+\lambda_{3}\Phi_{1}^{\dagger}\Phi_{1}\Phi_{1}^{\dagger}\Phi_{2}+\lambda_{4}\Phi_{1}^{\dagger}\Phi_{2}\Phi_{2}^{\dagger}\Phi_{1}+\dfrac{\lambda_{5}}{2}\left[\left(\Phi_{1}^{\dagger}\Phi_{1}\right)^{2}+\left(\Phi_{2}^{\dagger}\Phi_{2}\right)^{2}\right].$

The above Higgs potential contains 8 parameters, and the squared masses of
Higgs particles are given by

$\displaystyle\begin{split}&m_{H^{\pm},2HDM}^{2}=\left(\dfrac{m_{12}^{2}}{v_{1}v_{2}}-\lambda_{4}-\lambda_{5}\right)v^{2},\\\
&m_{A,2HDM}^{2}=\left(\dfrac{m_{12}^{2}}{v_{1}v_{2}}-2\lambda_{5}\right)(v_{1}^{2}+v_{2}^{2}),\\\
&m_{h,2HDM}^{2}=\dfrac{1}{2}\left[(\lambda_{1}v_{1}^{2}+\lambda_{2}v_{2}^{2})+m_{12}^{2}\dfrac{v^{2}}{v_{1}v_{2}}\right]\\\
&\qquad-\sqrt{\left[\dfrac{\lambda_{1}v_{1}^{2}-\lambda_{2}v_{2}^{2}}{2}-\dfrac{m_{12}^{2}}{2v_{1}v_{2}}(v_{1}^{2}-v_{2}^{2})\right]^{2}+(\lambda_{345}v_{1}v_{2}-m_{12}^{2})^{2}},\\\
&m_{H,2HDM}^{2}=\dfrac{1}{2}\left[(\lambda_{1}v_{1}^{2}+\lambda_{2}v_{2}^{2})+m_{12}^{2}\dfrac{v^{2}}{v_{1}v_{2}}\right]\\\
&\qquad+\sqrt{\left[\dfrac{\lambda_{1}v_{1}^{2}-\lambda_{2}v_{2}^{2}}{2}-\dfrac{m_{12}^{2}}{2v_{1}v_{2}}(v_{1}^{2}-v_{2}^{2})\right]^{2}+(\lambda_{345}v_{1}v_{2}-m_{12}^{2})^{2}}.\end{split}$
(38)

Looking at the formulas of masses in (38), we see that they contain the very
annoying mixing terms of $v_{1}$ and $v_{2}$. Hence, at nonzero temperatures,
these mixing terms make the job of calculating the contributions from the
particles to the effective potential very difficult. They turn the symmetry
breaking process into the process of breaking the ambiguous mixings of
$\Phi_{1}$ and $\Phi_{2}$.

Meanwhile, in the complex representation, the Higgs potential of 2HDM-$S_{3}$
has a simpler form:

$\displaystyle V_{\text{Higgs}}^{\text{2HDM}\otimes S_{3}}\equiv$
$\displaystyle\,V(\phi_{1},\phi_{2})$
$\displaystyle=\mu_{1}^{2}(\phi_{2}^{2}+\phi_{1}^{2})+\frac{1}{2}l_{1}(\phi_{2}^{2}+\phi_{1}^{2})^{2}+\frac{1}{2}l_{2}(\phi_{2}^{2}-\phi_{1}^{2})^{2}$
$\displaystyle\qquad+l_{3}(\phi_{1}.\phi_{2})(\phi_{2}.\phi_{1})-\mu_{2}^{2}\left(\phi_{1}.\phi_{2}+\phi_{2}.\phi_{1}\right).$
(39)

Seeing the structure of those Higgs potentials containing $S_{3}$ symmetry,
ones conclude that the Higgs potential becomes simpler, with fewer parameters
(from 8 reduced to 5). Thanks to that, the forms of Higgs mass are also
simpler. By replacing $v_{1}$ and $v_{2}$ by $c_{\beta}v$ and $s_{\beta}v$, in
expressions of $H$ and $h$ masses, ones get a form of $v$-dependent.

$\displaystyle\begin{split}&m_{H^{\pm},S_{3}}^{2}=-l_{2}v^{2},\\\
&m_{A,S_{3}}^{2}=-\dfrac{1}{2}(2l_{2}-l_{3})v^{2},\\\
&m_{H,S_{3}}^{2}=\dfrac{1}{4}v^{2}\bigg{[}2l_{1}+l_{3}+\sqrt{16l_{1}l_{2}(c_{\beta}^{2}-s_{\beta}^{2})^{2}-8l_{2}l_{3}+l_{3}^{2}-4l_{1}l_{3}(c_{\beta}^{4}-6c_{\beta}^{2}s_{\beta}^{2}+s_{\beta}^{4})}\bigg{]}=\frac{1}{4}f_{H}.v^{2},\\\
&m_{h,S_{3}}^{2}=\dfrac{1}{4}v^{2}\bigg{[}2l_{1}+l_{3}-\sqrt{16l_{1}l_{2}(c_{\beta}^{2}-s_{\beta}^{2})^{2}-8l_{2}l_{3}+l_{3}^{2}-4l_{1}l_{3}(c_{\beta}^{4}-6c_{\beta}^{2}s_{\beta}^{2}+s_{\beta}^{4})}\bigg{]}=\frac{1}{4}f_{h}.v^{2}\,.\end{split}$
(40)

The $S_{3}$ symmetry has removed the mixing parts of the two VEVs inside the
mass of two charged Higgs bosons $H^{\pm}$ and one neutral Higgs boson $A$.
Otherwise, $h$ and $H$ still have the mixing parts hiding in $c_{\beta}$ and
$s_{\beta}$. However, this mixing would be simpler in the 2HDM, since there
are no such ratios as $v_{1}/v_{2}$ or $v^{2}/(v_{1}.v_{2})$, their mass
formulas then have fewer parameters and are simpler than the ones of 2HDM
without $S_{3}$.

Therefore, according to the remarks for the 2HDM, $S_{3}$ had made the process
of the electroweak phase transition occur in each VEV. This will be discussed
in the next section.

###  IV.2 Structure of EWPT

The procedure to describe the structure of electroweak phase transition in
this model is similar to that in the SM, whereas Higgs and gauge bosons are
the main contributors in the breaking symmetry process. For that reason,
determination mass can also affect the phase transition.

The 2HDM with $S_{3}$ symmetry has Higgs Lagrangian with kinetic and potential
elements as

$\displaystyle\mathcal{L}_{\text{Higgs}}=(D_{\mu}\phi_{1})^{\dagger}(D^{\mu}\phi_{1})+(D_{\mu}\phi_{2})^{\dagger}(D^{\mu}\phi_{2})-V(\phi_{1},\phi_{2})\,.$
(41)

Averaging all over the space, then replacing fields with VEVs, the Higgs
Lagrangian with variables $v_{1}$ and $v_{2}$ (with $v_{1}\neq v_{2}$) has the
following form:

$\displaystyle\begin{split}\mathcal{L}_{\text{Higgs}}&=\dfrac{1}{2}\partial^{\mu}v_{1}\partial_{\mu}v_{1}+\dfrac{1}{2}\partial^{\mu}v_{2}\partial_{\mu}v_{2}-V_{0}(v_{1},v_{2})\\\
&\qquad\sum_{i=\text{vector
boson}}m_{i}^{2}(v_{1},v_{2})W^{\mu}W_{\mu}+\sum_{j=\text{Higgs
boson}}m_{j}^{2}(v_{1},v_{2})H^{2},\end{split}$ (42)

whereas $W$ and $H$ are the vector boson and scalar fields, respectively.

Table 2 contains the squared mass of the particles contributing to the EWPT,
in the form of depending VEVs; $n$ is the degree of freedom of the fields. The
masses of known particles are in generic form and at 0K, shown in Table 3.

Particles | $m^{2}(v_{1},v_{2})$ | $m^{2}(v_{1})$ | $m^{2}(v_{2})$ | $n$
---|---|---|---|---
$m^{2}_{W^{\pm}}$ | $\frac{g^{2}v^{2}}{4}$ | $\frac{g^{2}v_{1}^{2}}{4}$ | $\frac{g^{2}v_{2}^{2}}{4}$ | $6$
$m^{2}_{Z}$ | $(g^{2}+g^{\prime 2})\frac{v^{2}}{4}$ | $(g^{2}+g^{\prime 2})\frac{v_{1}^{2}}{4}$ | $(g^{2}+g^{\prime 2})\frac{v_{2}^{2}}{4}$ | $3$
$m^{2}_{h}$ | $\frac{1}{4}f_{h}v^{2}$ | $\frac{1}{4}f_{h}v_{1}^{2}$ | $\frac{1}{4}f_{h}v_{2}^{2}$ | $1$
$m^{2}_{H}$ | $\frac{1}{4}f_{H}v^{2}$ | $\frac{1}{4}f_{H}v_{1}^{2}$ | $\frac{1}{4}f_{H}v_{2}^{2}$ | $1$
$m^{2}_{A}$ | $-\dfrac{1}{2}(2l_{2}-l_{3})v^{2}$ | $-\dfrac{1}{2}(2l_{2}-l_{3})v_{1}^{2}$ | $-\dfrac{1}{2}(2l_{2}-l_{3})v_{2}^{2}$ | $1$
$m^{2}_{H^{\pm}}$ | $-l_{2}v^{2}$ | $-l_{2}v_{1}^{2}$ | $-l_{2}v_{2}^{2}$ | $2$
$m_{t}^{2}$ | $f_{t}^{2}v^{2}$ | $f_{t}^{2}v_{1}^{2}$ | $f_{t}^{2}v_{2}^{2}$ | $-12$
Table 2: Squared mass of the gauge bosons and scalar bosons in 2HDM-$S_{3}$; whereas mass of the $W^{\pm}$, $Z$ and $t$ is the same as the one in SM; $v^{2}=v_{1}^{2}+v_{2}^{2}$. Particles | $m_{W^{\pm}}(v_{0})$ | $m_{Z}(v_{0})$ | $m_{h}(v_{0})$ | $m_{t}(v_{0})$
---|---|---|---|---
$m(v_{0})~{}[\textrm{GeV}]$ | $80.442$ | $91.18$ | $125$ | $173.1$
Table 3: Mass of particles (GeV) at $0$K in 2HDM-$S_{3}$

Table 2 shows us that all the particles in the model depend on two VEVs. But
$v_{1}$ and $v_{2}$ depend on each other, $v_{2}/v_{1}=\tan\beta$. The forms
of mass could be changed into one-VEV-depended
($v=\sqrt{v_{1}^{2}+v_{2}^{2}}$) by replacing $v_{1}=vc_{\beta}$ and
$v_{2}=vs_{\beta}$. Therefore, in this model, assuming the remark 4, the
electroweak phase transition can be considered as a dual transition, with two
VEVs accomplished to the condition $v=\sqrt{v_{1}^{2}+v_{2}^{2}}$ and at 0K,
$v_{0}=\sqrt{v_{01}^{2}+v_{02}^{2}}=246$ GeV.

Here in remark 4, it is also a note that the coefficient $\kappa$ is now equal
to $\tan\beta$. Carefully observing Table 2, although the masses of $h$ and
$H$ can be split, there still exists the coefficient $\kappa$ in $f_{h}$ and
$f_{H}$. These constants are only meaningful in that the masses of $h$ and $H$
can be split into terms that each one of them depends only on one VEV. In
other words, the contributions of the mixing of $v_{1}$ and $v_{2}$ are all
brought back to only one VEV is $v_{1}$ or $v_{2}$, and the difference between
$v_{1}$ and $v_{2}$ is put into the constants $f_{h}$ and $f_{H}$. Hence, the
investigation of phase transition is somewhat relatively easier. But in the
end, the replacement is really not that important, since in the next sections,
this coefficient in fact will not have any effects on the strength of phase
transition.

###  IV.3 The effective potential

This dual-phase transition has the participation of new particles as two
charged Higgs $H^{\pm}$, one neutral CP-odd $A$, and one neutral CP-even $H$.
More importantly, there are also the presence of SM particles as one neutral
CP-even Higgs boson $h$, two charged gauge bosons $W^{\pm}$, one neutral boson
$Z$ and top quark $t$.

The effective potential for each stage can be calculated in two ways. The
effective potentials only contain the contributions from the particles outside
of the SM and the gauge bosons, SM-like Higgs boson and top quark. The other
particles have small values of mass so they can be just ignored. The process
of calculating the effective potential is in detail given in Ref. pkll .

The effective potential of one phase transition without daisy loops has the
form:

$\displaystyle\begin{split}V_{eff}(\mathcal{V},T)=&V_{0}(\mathcal{V})+\dfrac{1}{64\pi^{2}}\bigg{[}6m_{W^{\pm}}^{4}(\mathcal{V})\ln\dfrac{m_{W^{\pm}}^{2}(\mathcal{V})}{Q^{2}}+3m_{Z}^{4}(\mathcal{V})\ln\dfrac{m_{Z}^{2}(\mathcal{V})}{Q^{2}}+2m_{H^{\pm}}^{4}(\mathcal{V})\ln\dfrac{m_{H^{\pm}}^{2}(\mathcal{V})}{Q^{2}}\\\
&+m_{h}^{4}(\mathcal{V})\ln\dfrac{m_{h}^{2}(\mathcal{V})}{Q^{2}}+m_{H}^{4}(\mathcal{V})\ln\dfrac{m_{H}^{2}(\mathcal{V})}{Q^{2}}+m_{A}^{4}(\mathcal{V})\ln\dfrac{m_{A}^{2}(\mathcal{V})}{Q^{2}}-12m_{t}^{4}(\mathcal{V})\ln\dfrac{m_{t}^{2}(\mathcal{V})}{Q^{2}}\bigg{]}\\\
&+\dfrac{T^{4}}{4\pi^{2}}\bigg{[}6F_{-}\left(\dfrac{m_{W^{\pm}}(\mathcal{V})}{T}\right)+3F_{-}\left(\dfrac{m_{Z}(\mathcal{V})}{T}\right)+2F_{-}\left(\dfrac{m_{H^{\pm}}(\mathcal{V})}{T}\right)+F_{-}\left(\dfrac{m_{h}(\mathcal{V})}{T}\right)\\\
&+F_{-}\left(\dfrac{m_{H}(\mathcal{V})}{T}\right)+F_{-}\left(\dfrac{m_{A}(\mathcal{V})}{T}\right)+12F_{+}\left(\dfrac{m_{t}(\mathcal{V})}{T}\right)\bigg{]},\end{split}$
(43)

whereas

$\displaystyle
F_{\pm}\left(\dfrac{m_{\phi}}{T}\right)=\int_{0}^{\frac{m_{\phi}}{T}}\alpha
J_{\pm}^{(1)}(\alpha,0)d\alpha,$ (44) $\displaystyle
J_{\pm}^{(1)}(\alpha,0)=2\int_{\alpha}^{\infty}\dfrac{(x^{2}-\alpha^{2})^{\nu/2}}{e^{x}\pm
1}dx.$ (45)

Then,

$\displaystyle\begin{cases}&J_{-}^{(1)}(\alpha,0)=\dfrac{\pi^{2}}{3}-\pi\alpha-\dfrac{\alpha^{2}}{2}\left(\ln\dfrac{\alpha}{4\pi}+C-\dfrac{1}{2}\right)+\mathcal{O}(\alpha^{2}),\\\
&J_{+}^{(1)}(\alpha,0)=\dfrac{\pi^{2}}{6}-\dfrac{\alpha^{2}}{2}\left(\ln\dfrac{\alpha}{\pi}+C-\dfrac{1}{2}\right)+\mathcal{O}(\alpha^{2}).\\\
\end{cases}$ (46)

This process has the contribution from five Higgs particles in total into the
effective potential. However, there are only two scalar Higgs particles $h,H$,
which are associated with the two nonzero vacuum expectation values
$v_{1},v_{2}$. Therefore, the minimum conditions then are

$\displaystyle V_{eff}(\mathcal{V}_{0},0)=0,\qquad\dfrac{\partial
V_{eff}(\mathcal{V},0)}{\partial\mathcal{V}}\bigg{|}_{\mathcal{V}=\mathcal{V}_{0}}=0,$
(47)
$\displaystyle\dfrac{\partial^{2}V_{eff}(\mathcal{V},0)}{\partial\mathcal{V}^{2}}\bigg{|}_{\mathcal{V}=\mathcal{V}_{0}}=\left[m_{h}^{2}(\mathcal{V})+m_{H}^{2}(\mathcal{V})\right]_{\mathcal{V}=\mathcal{V}_{0}}.$
(48)

With the minimum conditions, expanding the functions of $J_{\pm}$, the
effective potential can be rewritten as

$\displaystyle
V_{eff}(\mathcal{V})=\dfrac{\lambda_{T}}{4}\mathcal{V}^{4}-\theta
T\mathcal{V}^{3}+\gamma(T^{2}-T_{0}^{2})\mathcal{V}^{2},$ (49)

where,

$\displaystyle\begin{split}&\lambda_{T}=\dfrac{m_{h}^{2}(\mathcal{V}_{0})+m_{H}^{2}(\mathcal{V}_{0})}{2\mathcal{V}_{0}^{2}}\bigg{\\{}1+\dfrac{1}{8\pi^{2}\mathcal{V}_{0}^{2}\left[m_{h}^{2}(\mathcal{V}_{0})+m_{H}^{2}(\mathcal{V}_{0})\right]}\times\\\
&\quad\bigg{[}6m_{W^{\pm}}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{W^{\pm}}^{2}(\mathcal{V}_{0})}+3m_{Z}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{Z}^{2}(\mathcal{V}_{0})}+2m_{H^{\pm}}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{H^{\pm}}^{2}(\mathcal{V}_{0})}\\\
&\quad+m_{h}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{h}^{2}(\mathcal{V}_{0})}+m_{H}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{H}^{2}(\mathcal{V}_{0})}+m_{A}^{4}(\mathcal{V}_{0})\ln\dfrac{bT^{2}}{m_{A}^{2}(\mathcal{V}_{0})}-12m_{t}^{4}(\mathcal{V}_{0})\ln\dfrac{b_{F}T^{2}}{m_{t}^{2}(\mathcal{V}_{0})}\bigg{]}\bigg{\\}},\\\
&\theta=\dfrac{1}{12\pi\mathcal{V}_{0}^{3}}\bigg{[}6m_{W^{\pm}}^{3}(\mathcal{V}_{0})+3m_{Z}^{3}(\mathcal{V}_{0})+2m_{H^{\pm}}^{3}(\mathcal{V}_{0})+m_{h}^{3}(\mathcal{V}_{0})+m_{H}^{3}(\mathcal{V}_{0})+m_{A}^{3}(\mathcal{V}_{0})\bigg{]},\\\
&\gamma=\dfrac{1}{24\mathcal{V}_{0}^{2}}\bigg{[}6m_{W^{\pm}}^{2}(\mathcal{V}_{0})+3m_{Z}^{2}(\mathcal{V}_{0})+2m_{H^{\pm}}^{2}(\mathcal{V}_{0})+m_{h}^{2}(\mathcal{V}_{0})+m_{H}^{2}(\mathcal{V}_{0})+m_{A}^{2}(\mathcal{V}_{0})+6m_{t}^{2}(\mathcal{V}_{0})\bigg{]},\\\
&T_{0}^{2}=\dfrac{1}{4\gamma}\bigg{\\{}m_{h}^{2}(\mathcal{V}_{0})+m_{H}^{2}(\mathcal{V}_{0})-\dfrac{1}{8\pi^{2}\mathcal{V}_{0}^{2}}\bigg{[}6m_{W^{\pm}}^{4}(\mathcal{V}_{0})+3m_{Z}^{4}(\mathcal{V}_{0})+2m_{H^{\pm}}^{4}(\mathcal{V}_{0})\\\
&\qquad\qquad+m_{h}^{4}(\mathcal{V}_{0})+m_{H}^{4}(\mathcal{V}_{0})+m_{A}^{4}(\mathcal{V}_{0})-12m_{t}^{4}(\mathcal{V}_{0})\bigg{]}\bigg{\\}}.\end{split}$
(50)

The critical temperature $T_{c}$ is given by

$T_{c}=\frac{T_{0}}{\sqrt{1-\theta^{2}/[\gamma\lambda_{T_{c}}}]},$ (51)

and the critical VEV can be derived as

$\mathcal{V}_{c}=\frac{2\theta T_{c}}{\lambda_{T_{c}}}.$ (52)

Therefore, the strength of EWPT is

$S=\frac{\mathcal{V}_{c}}{T_{c}}=\frac{2\theta}{\lambda_{T_{c}}}.$ (53)

Next, taking into account daisy loops, the effective potential will have the
form:

$\displaystyle V^{daisy}_{eff}=V_{eff}(\mathcal{V})-V^{d}(\mathcal{V}),$ (54)

in which the second component on the right-hand side of Eq. (54) is the
contribution of daisy loops carrington ; curtin ; katz  (especially the
appendix A in Ref. katz ). Here, degrees of freedom are given by:
$g_{Z}=3,g_{W}=6,g_{h}=g_{A}=g_{H}=1,g_{H^{\pm}}=2$ and

$\displaystyle
V^{d}(\mathcal{V})=\frac{T}{12\pi}\sum_{i=h,W,Z,A,H,H^{\pm}}g_{i}\left\\{\left[\frac{m_{i}^{2}(\mathcal{V}_{0})\mathcal{V}^{2}}{\mathcal{V}^{2}_{0}}+\Pi_{i}(T)\right]^{3/2}-\frac{m_{i}^{3}(\mathcal{V}_{0})\mathcal{V}^{3}}{\mathcal{V}^{3}_{0}}\right\\},$
(55)
$\displaystyle\Pi_{W}(T)=\frac{22}{3}\frac{m_{W}^{2}(\mathcal{V}_{0})}{\mathcal{V}^{2}_{0}}T^{2},$
$\displaystyle\Pi_{Z}(T)=\frac{22}{3}\frac{(m_{Z}^{2}(\mathcal{V}_{0})-m_{W}^{2}(\mathcal{V}_{0}))}{\mathcal{V}^{2}_{0}}T^{2},$
$\displaystyle\Pi_{h}(T)=\frac{2m_{W}^{2}(\mathcal{V}_{0})+m_{Z}^{2}(\mathcal{V}_{0})+m_{h}^{2}(\mathcal{V}_{0})+2m_{t}^{2}(\mathcal{V}_{0})}{4\mathcal{V}^{2}_{0}}T^{2}+(\Lambda_{hH}+\Lambda_{hA}+\Lambda_{hH^{\pm}}).T^{2}.$
(56)

As for exotic Higgs particles then

$\displaystyle\Pi_{h-exotic}=(\Lambda_{hH}+\Lambda_{hA}+\Lambda_{hH^{\pm}}).T^{2}$
(57)
$\displaystyle\Pi_{A}(T)\sim\frac{m^{2}_{A}(\mathcal{V}_{0})}{\mathcal{V}_{0}^{2}}T^{2},$
$\displaystyle\Pi_{H}(T)\sim\frac{m^{2}_{H}(\mathcal{V}_{0})}{\mathcal{V}_{0}^{2}}T^{2},$
$\displaystyle\Pi_{H^{\pm}}(T)\sim\frac{m^{2}_{H^{\pm}}(\mathcal{V}_{0})}{\mathcal{V}_{0}^{2}}T^{2}.$
(58)

$\Lambda_{hH}+\Lambda_{hA}+\Lambda_{hH^{\pm}}$ are coefficients representing
the contribution from exotic Higgs daisy loops to SM-like Higgs boson. $A$,
$H$, $H^{\pm}$ are called exotic Higgs for short.

The daisy loops of exotic Higgs boson can be omitted, since these masses of
particles are large and $m(\mathcal{V})/T\sim 1$ carrington . This can be
explained in Sec. IV.6.

Note that $\mathcal{V}$ can be $v_{1},v_{2}$ or $v$. Let $v^{2}=a.v^{2}_{2}$,
the relation between $\tan\beta$ and $a$ is

$\tan\beta=\sqrt{1/(a-1)}.$ (59)

From this, the masses of the particles in terms of $a$ are given in Table 4.

Particles | $m^{2}(v_{1},v_{2})$ | $m^{2}(v_{2})$ | $m^{2}(v_{1})$
---|---|---|---
$m^{2}_{W^{\pm}}$ | $\frac{g^{2}v^{2}}{4}$ | $m^{2}_{W^{\pm}}/a$ | $m^{2}_{W^{\pm}}(v_{2}).(a-1)$
$m^{2}_{Z}$ | $(g^{2}+g^{\prime 2})\frac{v^{2}}{4}$ | $m^{2}_{Z}/a$ | $m^{2}_{Z}(v_{2}).(a-1)$
$m^{2}_{h}$ | $\frac{1}{4}f_{h}v^{2}$ | $m^{2}_{h}/a$ | $m^{2}_{h}(v_{2}).(a-1)$
$m^{2}_{H}$ | $\frac{1}{4}f_{H}v^{2}$ | $m^{2}_{H}/a$ | $m^{2}_{H}(v_{2}).(a-1)$
$m^{2}_{A}$ | $-\dfrac{1}{2}(2l_{2}-l_{3})v^{2}$ | $m^{2}_{A}/a$ | $m^{2}_{A}(v_{2}).(a-1)$
$m^{2}_{H^{\pm}}$ | $-l_{2}v^{2}$ | $m^{2}_{H^{\pm}}/a$ | $m^{2}_{H^{\pm}}(v_{2}).(a-1)$
$m_{t}^{2}$ | $f_{t}^{2}v^{2}$ | $m^{2}_{t}/a$ | $m^{2}_{t}(v_{2}).(a-1)$
Table 4: Squared mass of the gauge bosons and scalar bosons in 2HDM-$S_{3}$.

With $v_{1}$ different from $v_{2}$, this model has two stages of phase
transition. We assume that $v_{1}<v_{2}$, which means $1<a<2$ or
$1<\tan\beta$. As the above sections pointed out, in 2HDM-$S_{3}$, the
particles’ masses can be changed such that there are no mixing terms between
$v_{1}$ and $v_{2}$. So the correct effective potential for this model is

$V_{eff}^{S_{3}}=V_{eff}^{S_{3}}(v_{1})+V_{eff}^{S_{3}}(v_{2}).$ (60)

With the above formula for the potential, the phase transition’s strength does
not depend on $a$.

###  IV.4 Probing the independence of EWPT strength on $\tan\beta$

Note that the EWPT in this model occurs in two stages and the mass components
of the involved particles are given in Table 4. Let us assume that the phase
transition’s strength of the first stage has already been calculated,
$S_{1}>1$ for $\mathcal{V}=v_{2}$.

It would like to prove that the second phase transition’s strength ($S_{2}$),
corresponds to $\mathcal{V}^{2}=v_{1}^{2}=v_{2}^{2}(a-1)$ will actually not
change, that it is still equal to $S_{1}$. Or in other words, it does not
depend on $a$.

To do this, the functions
$\lambda_{T_{c}}(v_{1}),\theta(v_{1}),\gamma(v_{1}),T_{c}(v_{1})$ must be
indicated the independent of $a$ (equal to themselves when calculated with
$v_{2}$). First, let us consider the function $\theta$ correspond to
$\mathcal{V}=v_{1}$,

$\theta(v_{1})=\dfrac{1}{12\pi
v_{1}^{3}}\bigg{[}6m_{W^{\pm}}^{3}(v_{1})+3m_{Z}^{3}(v_{1})+2m_{H^{\pm}}^{3}(v_{1})+m_{h}^{3}(v_{1})+m_{H}^{3}(v_{1})+m_{A}^{3}(v_{1})\bigg{]}.$
(61)

The masses in the bracket and $v_{1}$ all have the same power of 3. Hence they
are all proportional to $(a-1)^{3/2}$. By extracting this factor out of the
masses and canceling it with the exact same factor from $v_{1}$ in the
denominator, $\theta(v_{1})$ is independent of $a$, or
$\theta(v_{1})=\theta(v_{2})$. Similarly $\gamma(v_{1})$ does not depend on
$a$, but $T^{2}_{0}(v_{1})$ depends on $(a-1)$, specifically
$T^{2}_{0}(v_{1})=(a-1)T^{2}_{0}(v_{2})$. Because of this dependence on
$(a-1)$ of $T_{0}^{2}(v_{1})$, $\lambda_{T_{c}}(v_{1})$ will not depend on
$a$. The proof is as follows. Consider the function $\lambda_{T_{c}}$
corresponds to $v_{1}$,

$\displaystyle\begin{split}\lambda_{T_{c}}(v_{1})&=\dfrac{m_{h}^{2}(v_{1})+m_{H}^{2}(v_{1})}{2v_{1}^{2}}\bigg{\\{}1+\dfrac{1}{8\pi^{2}v_{1}^{2}\left[m_{h}^{2}(v_{1})+m_{H}^{2}(v_{1})\right]}\times\\\
&\quad\bigg{[}6m_{W^{\pm}}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{W^{\pm}}^{2}(v_{1})}+3m_{Z}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{Z}^{2}(v_{1})}+2m_{H^{\pm}}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{H^{\pm}}^{2}(v_{1})}\\\
&\quad+m_{h}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{h}^{2}(v_{1})}+m_{H}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{H}^{2}(v_{1})}+m_{A}^{4}(v_{1})\ln\dfrac{b{T_{c}}^{2}(v_{1})}{m_{A}^{2}(v_{1})}\\\
&-12m_{t}^{4}(v_{1})\ln\dfrac{b_{F}{T_{c}}^{2}(v_{1})}{m_{t}^{2}(v_{1})}\bigg{]}\bigg{\\}}.\end{split}$
(62)

With $T_{c}$ is given by Eq. (51) in the no daisy loop case. By the same
reasoning from above, $\lambda_{T_{c}}(v_{1})$ will not depend on $a$ when the
logarithmic factors do not depend on $a$. Indeed, we consider the general
expression inside the logarithmic functions:

$\displaystyle\begin{split}\frac{T_{c}^{2}(v_{1})}{m^{2}(v_{1})}&=\frac{T_{0}^{2}(v_{1})}{\left[1-\theta^{2}(v_{1})/\gamma(v_{1})\lambda_{T_{c}}(v_{1})\right]m^{2}(v_{2})}\\\
&=\frac{(a-1)T^{2}_{0}(v_{2})}{\left(1-\theta^{2}(v_{2})/\gamma(v_{2})\lambda_{T_{c}}(v_{1})\right)(a-1)m^{2}(v_{2})}\\\
&=\frac{T^{2}_{0}(v_{2})}{\left(1-\theta^{2}(v_{2})/\gamma(v_{2})\lambda_{T_{c}}(v_{1})\right)m^{2}(v_{2})},\end{split}$
(63)

in which the functions $\theta(v_{1}),\gamma(v_{1})$ are all independent of
$a$ as proven earlier. The expression inside the logarithmic functions
actually depends on $\lambda_{T_{c}}(v_{1})$. Substitute the expression Eq.
(63) into Eq. (62), we can finally realize that $a$ no longer appears in the
expression Eq. (62). This proves that $\lambda_{T_{c}}(v_{1})$ does not depend
on $a$. So, the functions
$\lambda_{T_{c}}(v_{1}),\theta(v_{1}),\gamma(v_{1}),T_{c}(v_{1})$ is truly
independent of $a$ when the effective potential without daisy loops.

When the effective potential with daisy loops which is Eq. (54), the critical
temperature $T_{C}(v_{1})$ is not Eq. (51) but $T_{C}(v_{1})\sim
T_{0}(v_{1})$. So $\lambda_{T_{C}}(v_{1})$ still does not depend on $a$.
Furthermore, by the similar proof, $V^{d}(v_{1})$ does not depend on $a$.
Because $V^{d}(v_{1})$ just depends on the ratio
$m^{2}(\mathcal{V})/\mathcal{V}^{2}$. So finally at the critical temperature,
the effective potential with daisy loops $[V^{daisy-S_{3}}_{eff}(v_{1})]$
remains independent of $a$ and deduced that $S$ in the same regardless of
being calculated with $v_{1}$ or $v_{2}$.

This result agrees with the conclusions in Ref. davidson , which concludes
that $\tan\beta$ is not a meaningful parameter in the 2HDM.

Also commented in Sec. II.5, in the Fig. 1 of Ref. dori , in the 2HDM model,
since the ratio between the two phase transition strengths may not depend much
on $\tan\beta$ . In other words, the strength of the phase transition can be
independent of $\tan\beta$. We have clearly demonstrated this in the
2HDM-$S_{3}$ model, thanks to the $S_{3}$ symmetry that separates the two
phase transitions.

###  IV.5 The true critical temperatures

To indicate critical temperatures in the model, the effective potential
without daisy loops is only used. The estimation of daisy loop contributions
will be done in Sec. IV.6. As analyzed in Sec. IV.2, the 2HDM-$S_{3}$ model
will have two critical temperatures which correspond to the two stages.

Since the coefficients are independent of $a$ as shown above, the parameters
of the second phase transition can be expressed in terms of the parameters of
the first phase transition,

$\displaystyle\lambda_{T,v_{1}}$
$\displaystyle=\lambda_{T,v_{2}}=\lambda_{T,v},$ (64)
$\displaystyle\theta_{v_{1}}$ $\displaystyle=\theta_{v_{2}}=\theta_{v},$ (65)
$\displaystyle\gamma_{v_{1}}$ $\displaystyle=\gamma_{v_{2}}=\gamma_{v},.$ (66)
$\displaystyle T^{2}_{0,v_{2}}$
$\displaystyle=T^{2}_{0,v_{1}}/(a-1)=T^{2}_{0,v}/a.$ (67)

With these equalities, the effective potentials of the second stage and the
combined stages can be expressed in terms of the effective potential of the
first phase transition stage

$\displaystyle V_{eff}^{S_{3}}(v_{2})$
$\displaystyle=\dfrac{\lambda_{T}}{4}v_{2}^{4}-\theta
T_{v_{2}}v^{3}_{1}+\gamma(T^{2}_{v_{1}}-T_{0,v_{1}}^{2})v_{1}^{2},$ (68)
$\displaystyle V_{eff}^{S_{3}}(v_{1})$
$\displaystyle=\dfrac{\lambda_{T}}{4}v_{1}^{4}-\theta
T_{v_{1}}v^{3}_{1}+\gamma(T^{2}_{v_{1}}-T_{0,v_{1}}^{2})v_{1}^{2}$ (69)
$\displaystyle=(a-1)^{2}\dfrac{\lambda_{T}}{4}v_{2}^{4}-(a-1)^{2}\theta
T_{v_{2}}v^{3}_{2}+(a-1)^{2}\gamma(T^{2}_{v_{2}}-T_{0,v_{2}}^{2})v_{2}^{2}$
(70) $\displaystyle=(a-1)^{2}V_{eff}(v_{2}),$ (71) $\displaystyle V_{eff}(v)$
$\displaystyle=\dfrac{\lambda_{T}}{4}v^{4}-\theta
T_{v}v^{3}+\gamma(T^{2}_{v}-T_{0,v}^{2})v^{2}$ (72)
$\displaystyle=a^{2}\dfrac{\lambda_{T}}{4}v_{2}^{4}-a^{2}\theta
T_{v_{2}}v^{3}_{2}+a^{2}\gamma(T^{2}_{v_{2}}-T_{0,v_{2}}^{2})v_{2}^{2}=a^{2}V_{eff}(v_{2}).$
(73)

Hence, $V_{eff}(v)\neq
V_{eff}^{S_{3}}(v_{1})+V_{eff}^{S_{3}}(v_{2})=V^{S_{3}}_{eff}$. From these
equations, $V_{eff}(v)$ can be deduced

$V_{eff}(v)=\frac{a^{2}}{(a-1)^{2}+1}V^{S_{3}}_{eff}=f(a)V^{S_{3}}_{eff}.$
(74)

Here when writing down the effective potential of the system in terms of
$V_{eff}(v)$, the correct effective potential of the system has been
multiplied by $\frac{a^{2}}{(a-1)^{2}+1}$ times.

Figure 5: The function $f(a)$, the ratio between $V^{S_{3}}_{eff}$ and
$V_{eff}(v)$

According to Fig. 5, it follows $1<a<2$, and function $f(a)$ has value of $1$
at $a=1$, that is $V_{eff}(v)=V^{S_{3}}_{eff}$, or $v_{1}=0$, but this cannot
be true. The maximum value of function $f(a)$ is $2$ when $a=2$, which is
$v_{1}=v_{2}$. Therefore when the two VEVs are equal to each other, the
difference between the two effective potentials is at maximum. When $a>2$ or
$v_{1}>v_{2}$, it can exchange $v_{1}\longleftrightarrow v_{2}$, so that this
case is similar to the case of $a<2$.

When the effective potential is rewritten as $V_{eff}(v)$ when calculating the
phase transition’s strength, our $S$ is correct. However, now the temperature
of the phase transition $T_{C}$ of the system turns out not to be correct. If
$v_{1}\neq v_{2}$, our system has two phase transition temperatures
$T_{C1}<T_{C}$ and $T_{C2}<T_{C}$, respectively. Hence, writing down the
effective potential for the system as $V_{eff}(v)$ is just a way to compare it
to that of SM (or putting this model in the context of SM, we call it “the SM-
like effective potential”), and now $T_{C}$ is not the true temperature for
the phase transition of the system, and it should be called “the SM-like
critical temperature.”

As commented earlier in Sec. II.5, when studying phase transitions in 2HDM,
studies rarely mention phase transition temperature. Because the analysis of
the phase transition temperature of the two VEVs model would be very difficult
due to their mixing. As analyzed in this section, the $S_{3}$ symmetry
separates the two phase transitions, so it makes the determination of the
transition temperature more obvious.

(a) $V_{eff}^{S_{3}}(v_{2})$ and $V_{eff}^{S_{3}}(v_{1})$ are the first and
second EWWPT when $a=3/2$ and $S=1$.

(b) $V_{eff}^{S_{3}}(v_{2})$ at several temperatures when $a=3/2$ and $S=1$.

Figure 6: The effective potential with $a=3/2$.

Fig. 6(b) shows us the potential $V_{eff}^{S_{3}}(v_{2})$ at different
temperatures. These potentials all have the second non zero minima, with
$m_{H}(v_{2})=150~{}\textrm{GeV}$,
$m_{H^{\pm}}(v_{2})=m_{A}(v_{2})=302.087~{}\textrm{GeV}$. The solid line that
corresponds to $T_{C1}=139.739~{}\textrm{GeV}$, shows that there exists a
potential well between the two minima. This is proof of the existence of the
first-order in the phase transitions.

In Fig. 6(a), the solid line is the shape of $V_{eff}^{S_{3}}(v_{2})$ at
$m_{H}(v_{2})=150~{}\textrm{GeV}$,
$m_{H^{\pm}}(v_{2})=m_{A}(v_{2})=302.087~{}\textrm{GeV}$, and
$T_{C1}=139.739~{}\textrm{GeV}$. The dash-dotted line is
$V_{eff}^{S_{3}}(v_{1})$ when $m_{H}(v_{1})=150/\sqrt{2}~{}\textrm{GeV}$,
$m_{H^{\pm}}(v_{1})=m_{A}(v_{1})=302.087/\sqrt{2}~{}\textrm{GeV}$, and
$T_{C2}=98.761~{}\textrm{GeV}$. The nonzero minimum and maximum of the dash-
dotted are smaller than that of the solid line and $T_{C1}>T_{C2}$, which
shows that the phase transition must occur in two stages. The distance between
the two stages is $\Delta T_{C}=T_{C1}-T_{C2}=40.978~{}\textrm{GeV}$.

To be more intuitively in the comments, we plot the effective potential
$V^{S_{3}}_{eff}(v_{2})$ and $V_{eff}(v)$ in case of $a=2$ as in Fig. 7. The
blue solid line is the potential $V_{eff}(v_{2})$ when
$m_{H}(v_{2})=140~{}\textrm{GeV}$,
$m_{H^{\pm}}(v_{2})=m_{A}(v_{2})=272.647~{}\textrm{GeV}$ and the phase
transition’s temperature $T_{C1}=124.283~{}\textrm{GeV}V$. The dashed line is
for $m_{A}(v)=198~{}\textrm{GeV}$,
$m_{H^{\pm}}(v)=m_{A}(v)=362.284~{}\textrm{GeV}$ and the critical temperature
is $T_{C}=169.288~{}\textrm{GeV}$.

Figure 7: $V_{eff}^{S_{3}}(v_{2})$ and $V_{eff}(v)$ when $v_{1}=v_{2}$ or
$a=2,\tan\beta=1$

Through Fig. 7, when $a=2$, that is when $v_{1}=v_{2}$, the two stages of the
EWPT occur at the same time. Each stage of the phase transition is described
by the blue solid line. The correct temperature for the phase transition of
the system must be $T_{C1}$.

However, also through Fig. 7, if we describe our system using the composite
potential $V_{eff}(v)$, that is if we study the phase transition of the system
in just one stage as in SM, the effective potential of the system at the phase
transition temperature is as the dashed line. Here the temperatures for the
system are $T_{C}>T_{C1}$. By describing the EWPT in only one stage as in the
SM, the two effective potentials were turned in the solid line into the dashed
line. Accidentally, this did not change the strength of the phase transition
but instead increased the critical phase transition’s temperature and VEV,
making them different compared to the correct ones.

###  IV.6 Searching the first-order EWPT and the role of $\tan\beta$

In order to meet a first-order phase transition, the transition strength must
have its value bigger or at least equal to 1. However, there are three unknown
variables $m_{A},m_{H^{\pm}}$ and $m_{H}$ in our problem. Therefore we can
assume that $m_{A}\equiv m_{H^{\pm}}$. This assumption is only intended to
reduce the number of variables in the problem and find the domain for the
masses of the particles in the first-order phase transition. This assumption
should not be applied to the parameters in the Higgs potential. Also assuming
that $m_{A}=m_{H^{\pm}}$ or $m_{H}=m_{H^{\pm}}$, but the results are all the
same. Choosing to use $m_{A}=m_{H^{\pm}}$ which is consistent with the
previous studies and the data for the parameter $\rho$ michela1 ; michelb1 ;
michelb2 ; michelb3 ; michelb4 .

According to the comments in Sec. II.5, especially in Ref. dori , there can be
many suggestions between the three quantities $m_{H^{\pm}},m_{H},m_{A}$ but it
is possible for a first order EWPT, so from the suggestions, we can choose the
scenario $m_{A}=m_{H^{\pm}}$ in the 2HDM-$S_{3}$ model. Because the symmetry
$S_{3}$ does not lose or add any of the three additional Higgs bosons.

Furthermore, the effective potential without daisy loop is first used to
calculate $S$. Then daisy loops will be additionally calculated later.

From there, the domain for the value of the masses cannot be too broad, in
fact, these domains must be closed. Indeed, to have $S>1$, and if we apply the
following conditions altogether: $T^{2}_{0}>0$ and $T_{C}$ must be real or
according to Eq. (51), $1-\theta^{2}/[\gamma\lambda_{T_{c}}]>0$, mass domains
must be closed. The numerical solution for the case of $a=3/2$ for the first
stage of the phase transition that corresponds to $v_{2}$ is given as Fig. 8.

Figure 8: The strength of first stage EWPT with $a=3/2$

The contours of $S=1$ are plotted onto the mass axes and then gradually
increase the value of $S$. With different values of $a$, a range of values of
the phase transition strength was find, $1\leq S<2.8$.

Moreover, according to Fig. 8, when $S$ increases, the domain for the masses
of the particles must narrow down. Therefore, in order to find these domains
for different values of $a$, only plotting $S=1$ with different $a$’s, which
is shown in Fig. 9.

Figure 9: The mass domain of the particles depends on $a$ with $S=1$ in the
first phase transition

Carefully looking at Fig. 9, when $a$ increases, the mass domains of the
particles narrow down. Since $1<a\leq 2$, the mass domain for each value of
$a$ as follows:

$\displaystyle 200\quad\textrm{GeV}<$ $\displaystyle
m_{H^{\pm}}(v_{2})<800\quad\textrm{GeV},$ (75) $\displaystyle 0<$
$\displaystyle m_{H}(v_{2})<600\quad\textrm{GeV}.$ (76)

The second stage of the phase transition that corresponds to $v_{1}$ is
similar to the first stage. And because the strength of the phase transition
does not depend on $a$, if the first stage of the phase transition has the
transition’s strength larger than 1 then so does stage 2. Therefore, the
domains for the masses of the particles at stage 2 can be indicated,

$\displaystyle 200\times\sqrt{(a-1)}\quad\textrm{GeV}<$ $\displaystyle
m_{H^{\pm}}(v_{1})<800\times\sqrt{(a-1)}\quad\textrm{GeV},$ (77)
$\displaystyle 0<$ $\displaystyle
m_{H}(v_{1})<600\times\sqrt{(a-1)}\quad\textrm{GeV}.$ (78)

Thus, according to Eq. (76) and Eq. (78), although $a$ (or $\tan\beta$) does
not affect the phase transition’s strength, it affects the domains of the
masses of the particles to have the first-order phase transition.

Combining Eqs. (75), (76), (77), (78) together, it follows that

$\displaystyle 200\times\sqrt{a}\quad\textrm{GeV}<$ $\displaystyle
m_{H^{\pm}}<800\times\sqrt{a}\quad\textrm{GeV},$ (79) $\displaystyle 0<$
$\displaystyle m_{H}<600\times\sqrt{a}\quad\textrm{GeV}.$ (80)

Notice that the maximum value of $a$ is 2, it follows from Eqs. (79), (80)
that for $a=\sqrt{2}$. So, the maximum value of masses are only increased by
about $1.41$ times. Therefore, the effect of $a$ on the mass domain of the
particles is not too large. This is also easy to see when observing the lines
in Fig. 9, they are very close together.

One more thing, if $a$ is closer to 2, the two stages are also closer to each
other. Hence, $a$ can be used to define the distance between the two-phase
transition stages in this model.

With the analysis of EWPT in 2HDM, the value of $\tan\beta$ is quite wide
(from 1 to 17 in all scenarios as indicated in Sec. II.5), it is almost a free
parameter in the EWPT problem. However, in the 2HDM-$S_{3}$ model, due to the
$S_{3}$ symmetry, the two EWPT stages are separated, thereby highlighting the
role of $\tan\beta$ which determines the width of the mass domain as well as
the gap between the two stages.

The last important part in this section will be estimating the contribution
from daisy loops. Based on Eq. (78), the masses of the exotic Higgs bosons are
usually chosen to be larger than the mass of the top quark. This is in keeping
with the difficulty of detecting these particles today. Because it only
consider the temperature region where EWPT occurs or
$\frac{\mathcal{V}_{c}}{T_{c}}>1$. So
$\frac{m_{A,H,H^{\pm}}(\mathcal{V}_{c})}{T_{c}}\sim 1$, the contribution of
daisy loops of exotic Higgs bosons will be small carrington . Thus when adding
daisy loops, only the daisy loops of gauge boson and SM-like Higgs boson,
i.e., Eqs. (56) are taken into account and neglecting $\Lambda_{h-exotic}$.
But notice that, in the temperature region $T\gg T_{C}$, the contributions
from daisy loops of the exotic Higgs bosons cannot be ignored.

Take a look at the graphs in Fig. 10, the red zone indicates the difference
between $V_{eff}$ and $V_{eff}^{daisy}$, the higher the temperature, the
larger the difference, and the larger the area of the red region. The phase
transition temperatures in this model are similar to SM, they are in the range
of $100$ to $150~{}\textrm{GeV}$. Because
$v=\sqrt{v^{2}_{1}+v^{2}_{2}}=246~{}\textrm{GeV}$ and
$v_{1},v_{2}<246~{}\textrm{GeV}$. So also from Fig. 10, when
$T<150~{}\textrm{GeV}$, the area of the red region is very small.

Thus, also from Fig. 10, daisy loops will increase the phase transition
temperature. Indeed, the second pair of lines in Fig. 10, when
$T=139.739~{}\textrm{GeV}$ is also the phase transition temperature
corresponding to $V_{eff}$ (the upper line). But $V_{eff}^{daisy}$ (bottom
line) has a second minimum that is below the VEV axis, so
$T=139.739~{}\textrm{GeV}$ is not yet the phase transition temperature
corresponding to $V_{eff}^{daisy}$, this phase transition temperature must be
greater than $139.739~{}\textrm{GeV}$. Finally, a sure result is that as the
phase transition temperature increases, the phase transition strength will
decrease.

Figure 10: Difference between $V^{S_{3}}_{eff}$ and $V_{eff}^{daisy}$ with
$m_{H}(v_{2})=150~{}\textrm{GeV},m_{H^{\pm}}(v_{2})=m_{A}(v_{2})=302.087~{}\textrm{GeV}.$

Next to see the effect from daisy loops about $S$ as the above comments. The
masses of $H^{\pm}$ and $H$ are randomly selected, then recalculate $S$ with
and without daisy loops, as shown in Table 5.

$m_{H}(v_{2})$ [GeV] |  $m_{H^{\pm}}(v_{2})$ [GeV] | $S_{nodaisy}$ | $S_{daisy}$ | $S_{daisy}/S_{nodaisy}$
---|---|---|---|---
300 | 500 | 1.116 | 1.0938 | 0.979
400 | 600 | 1.048 | 1.030 | 0.982
250 | 470 | 1.253 | 1.224 | 0.976
220 | 430 | 1.254 | 1.222 | 0.974
180 | 410 | 1.452 | 1.410 | 0.971
Table 5: Strength $S$ with and without daisy loops

The ring-loops reduce $S$ by $2/3$ times r23 ; 1101.4665 . But in this model,
this ratio is about $0.97$ times or the maximum value of strength is about
$2.71$. However, this still ensures the first-order EWPT.

Although the daisy loops do not significantly change the strength of the phase
transition in this model, but the role of reducing divergences in effective
potential calculations cannot be ignored.

In addition, it should be noted that in this model the triggers for the first
order phase transitions are heavy Higgs particles. So the daisy loops of heavy
Higgs have no effect in EWPT. However, other models exist light particles
besides the SM, surely the daisy loops of these particles will have a great
effect on EWPT.

##  V CONCLUSION AND OUTLOOKS

The 2HDM-$S_{3}$ was chosen to analyze the EWPT stages, not only because it is
close to the SM but also because the structure of the EWPT process in this
model is interesting enough for us to find new things. Moreover, although
these models have some features that are new to the SM, they are far more
complicated as they have more new fields and VEVs.

The symmetry breaking process that corresponds to the EWPT in the
2HDM-$S_{3}$, when compared with the SM, is depicted in the below diagram.

2HDM-$S_{3}$ $\displaystyle\Downarrow\text{breaking $v_{2}$}$
$\displaystyle\Downarrow\text{breaking $v_{1}$}$ $\displaystyle U(1)_{Q}$

We first summarize the structural analysis methods of the EWPT from previous
studies along with the method in this article as follows:

  * •

The first method as in Refs. thdm ; Fuyutob : The EWPT process is considered
as in the SM. The effective potential $V_{eff}(v)$ is used and
$v^{2}=v^{2}_{1}+v_{2}^{2}$.

  * •

The second method as in Refs. dori ; mayumi : Assuming the equation Eq. (95),
examining the multistep EWPT process, and the effective potential is a
function of two variables $v_{1}$ and $v_{2}$. But when calculating the
strength of phase transition
$\xi=\frac{v_{c}}{T_{c}},v^{2}=v^{2}_{1}+v^{2}_{2}$, it is still referred to
as one phase as the standard model.

  * •

The third method is in this article: Analyzing the division of the effective
potential into two separate parts $[V_{eff}=V_{eff}(v_{1})+V_{eff}(v_{2})]$
and examining each stage separately.

The EWPT has been intensively considered in the 2HDM. Some remarks are in
order.

In Refs. 124 ; 127 , the EWPT has been considered in the 2HDM type I and II.
In these papers, the masses of heavy particles are not larger than 1 TeV and
their mass difference is not bigger 400 GeV. In our results, the mass region
of heavy particles ranges from 200 to 800 GeV. Therefore, the maximum mass
difference between the heavy particles is about 600 GeV. The lattice
simulation with one-loop effective potential for crystal has been considered
in Ref. 127 , the first order EWPT happens for a scenario $m_{H^{\pm}}=m_{A}$,
and this agrees with our assumption. Therefore, the results in Refs. 124 ; 127
and ours are compatible.

It was shown that the EWPT is related to a significant uplifting of the Higgs
vacuum 125 . The first order phase transition leads to the following condition

$m_{A}>m_{H}+m_{Z}\,.$ (81)

The mass domain of particles in our calculation is also compatible with this.

The first order EWPT is possible if
$580\,\textrm{GeV}<m_{H^{\pm}}<1\,\textrm{TeV}$ 126 . This constraint agrees
with our result $m_{H^{\pm}}<800\times\sqrt{a}\,\textrm{GeV}\sim 1$ TeV as in
Eqs. (79) and (80).

In Ref. 128 , the first order EWPT happens in the 2HDM type I and II with the
mass difference between $H,H^{\pm}$ and $A$ in the range (100, 300) GeV. This
coincides with the mass region in Eqs. (75) and (76).

From Fig. 1 in Ref. 129  it follows that to have the first order phase
transition, the maximal mass difference among $m_{H}$ and $m_{A}$ is about 500
GeV. In our study this value is about 600 GeV. Hence, both results are
consistent.

In Ref. 130 , Fig. 3, the effective potential is plotted in the region of
masses lower than 600 GeV and the value of $\tan\beta$ runs from 1 to 20. The
lines of the effective potential for different $\tan\beta$ values are very
close to each other. The result shows that $S$ is almost independent on
$\tan\beta$, and this supports our conclusion.

All three methods are acceptable approximations. The effectiveness of the
third method in this article is clearly stated in Sec. IV. Mathematical
techniques for analyzing VEVs in the first and second method have been
successfully used to analyze decay channels in multi-VEV models hue ; alves .

Exploring the EWPT process into several stages has been analyzed in many other
models besides 2HDM, as shown in Refs. ptl ; borges . Therefore, the results
of this paper aim to emphasize the feature of the multiphase in 2HDM-$S_{3}$
model and investigate the factors that affect the division of stages as well
as the influence of $\tan\beta$.

The strength of the phase transition does not depend on $a$ or $\tan\beta$,
and to study the full structure of the two stages of the phase transition, the
effective potential must not be written as $V_{eff}(v)$. The greater
$\tan\beta$ is, the narrower the mass domain of the particles in the first
order of the phase transition becomes. The $S_{3}$ symmetry has shown that
there are two subsequential stages in the phase transition process, which has
not been shown clearly in the 2HDM where there are many mixing terms of the
VEVs in the mass domain of the particles.

The $S_{3}$ symmetry could explain the mixing of quarks or this symmetry could
be related to the Yukawa couplings that can affect EWPT processes bpal1 ;
bpal2 ; bpal3 . For example in Ref. braconi , changing the Yukawa coupling
constants results in a first order EWPT. Therefore, the $S_{3}$ symmetry
associated with the quark mixing has an effect on EWPT processes that need to
be further elucidated after these works.

By analyzing the effect of daisy loops, a way of assessing the contribution
from daisy loops: first the EWPT was calculated by using the effective
potential without daisy loops, to estimate the mass domain of particles; then
based on that mass domain, estimating the ratio $m(v)/T$ to consider the
contribution from daisy loops.

In this paper, in order to reduce the number of variables in our problem,
$m_{A}=m_{H^{\pm}}$ has been assumed. As said earlier, this assumption was
only made to find the mass domain of the particles, and must not be applied to
the parameters in the Higgs potential, since the real values of these 2
particles can be different from each other, even though their domains of mass
can be identical. However, from Eq. (40) and this assumption, $l_{3}\approx
0$. This assumption was made by the authors in order to satisfy the data of
the parameter $\rho$ in the 2HDM michela1 ; michelb1 ; michelb2 ; michelb3 ;
michelb4  so that if the 2HDM-$S_{3}$ also satisfies the data of the parameter
$\rho$, this can lead to $l_{3}$ being very small. This is one of the results
of this paper that can lead to research on the parameter $\rho$ in the
2HDM-$S_{3}$.

As stated in remark 4, the assessment of the impact the first stage of EWPT
has on the second stage of EWPT is made through the mixing terms of the VEVs.
However, the investigation of the effects of $\kappa$ is still unclear. All
these effects of $\kappa$ have been renormalized under the minimum conditions
of the later EWPT stage. Therefore, from remark 4, to assess the effects of
$\kappa$, we must rebuild the whole effective potential that contains the
mixing terms of VEVs all over. Assessing this direction is a new incoming, and
interesting job after this paper.

Notice the comments in the Sec. II.5, we focus on scalar decays into heavy
fermions ($A,H,H^{\pm}\longrightarrow tt,tb$), which are the most promising
channels for demonstrating the first-order EWPT through confirmation of
additional bosons shinya . Also, there is another way to check, we can measure
the gravitational waves generated by the EWPT process in future experiments by
LISA shinya .

Finally, through the comments on effective potentials, writing down the
effective potentials in the form of $V_{eff}(v)$ is imprecise, but still, it
is concise and gives accurate predictions for the strength of the phase
transition. However, there will be some small errors in calculating the
corresponding sphalaron energy, but these errors would not be large. Since the
contributions of the effective potential term in the sphaleron energy are
quite small, about 5.5% phong2022 , when writing down the effective potential
of the system in the form of
$V_{eff}(v)=\frac{a^{2}}{(a-1)^{2}+1}V^{S_{3}}_{eff}$, this will make the
sphaleron energy deviate by about 5.5%. From this, to make the calculation of
the sphaleron become more accurate, we only need to replace $V_{eff}(v)$ by
$\frac{(a-1)^{2}+1}{a^{2}}V_{eff}(v)$, and the methods are still the same as
in Refs. Fuyutob ; sphagt1 ; sphagt2 ; Fuyutoc . In addition, if the two
stages of the phase transition in this model occur at the same time or very
close to each other, the bubbles of each of the phase transition stages can
collide with each other or collide with the bubbles from other stages, thereby
causing some big gravitational waves. For that matter, the full estimation of
the contributing terms, as well as the impact of the ratio
$\frac{(a-1)^{2}+1}{a^{2}}$ on the sphaleron energy for the gravitational wave
calculation, will be the extension of this paper.

The method of high-temperature dimensional reduction to the 2HDM to obtain
three-dimensional effective theories that can be used for nonperturbative
simulations 131 . These results can be used to recalculate EWPT in 2HDM, and
to check our results. This is part of the upcoming work.

## ACKNOWLEDGMENTS

V.Q.P. would like to thank Pham Quang Khanh for reading and editing the
article in English and running my small code. H. N. L. is thankful to Van Lang
University.

##  Appendix A Effective action of multiscalar field models

In the $\phi^{4}$ theory and the single field case, calculating the effective
potential from summing the diagrams is shown in Ref. quiros . Let us consider
a toy model described by two neutral scalar fields ($\phi_{i},i=1,2$) with an
action

$S(\phi_{1},\phi_{2})=\int d^{4}x{\cal L}(\phi_{1}(x),\phi_{2}(x)),$ (82)

where

${\cal
L}(\phi_{1}(x),\phi_{2}(x))=\sum_{i=1}^{2}\left(\partial_{\mu}\phi^{*}_{i}\partial^{\mu}\phi_{i}-m^{2}_{\phi_{i}}\phi_{i}^{*}\phi_{i}\right)-\frac{\lambda_{\phi}}{4}\left(\phi^{*}_{1}\phi_{1}+\phi^{*}_{2}\phi_{2}\right)^{2}$
(83)

In the path-integral representation, the generating functional is as the
following:

$Z[J]=\langle 0_{\rm out}\mid 0_{\rm
in}\rangle_{J}\equiv\int\mathcal{D}\phi_{1}\mathcal{D}\phi_{2}\exp\\{i(S(\phi_{1},\phi_{2})+\phi_{1}J_{1}+\phi_{2}J_{2})\\},$
(84)

and

$Z[J]\equiv\exp\\{iW[J]\\},$ (85)

in which

$\phi_{i}J_{i}\equiv\int d^{4}x\phi_{1}(x)J_{1}(x);\quad i=1,2.$ (86)

The effective action $\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]$ is the
Legendre transform of Eq. (85)

$\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]=W[J]-\sum_{i=1}^{2}\int
d^{4}x\frac{\Delta W[J]}{\Delta J_{i}(x)}J_{i}(x)$ (87)

where the VEVs of $\phi_{i}$ are

$\overline{\phi_{i}}(x)=\frac{\delta W[J]}{\delta J_{i}(x)}$ (88)

From Eq. (87) and Eq. (88), the generating functionals can be obtained

$\frac{\delta\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]}{\delta\overline{\phi_{i}}}=\frac{\delta
W[J]}{\delta J_{i}}\frac{\delta
J_{i}}{\delta\overline{\phi_{i}}}-J_{i}-\overline{\phi_{i}}\frac{\delta
J_{i}}{\delta\overline{\phi_{i}}}=-J_{i}.$ (89)

where using of the notation from Eq. (86). Eq. (89) implies in particular that

$\left.\frac{\delta\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]}{\delta\overline{\phi_{i}}}\right|_{J_{i}=0}=0.$
(90)

$Z[J]$ can be expanded in a power series of $J$ and in terms of Green
functions $G_{(n)}$ as

$Z[J]=\sum^{\infty}_{n=0}\frac{i^{n}}{n!}\int d^{4}x_{1}\ldots
d^{4}x_{n}J(x_{1})\ldots J(x_{n})G_{(n)}(x_{1},\ldots,x_{n})$ (91)

and

$iW[J]=\sum^{\infty}_{n=0}\frac{i^{n}}{n!}\int d^{4}x_{1}\ldots
d^{4}x_{n}J(x_{1})\ldots J(x_{n})G^{\ c}_{(n)}(x_{1},\ldots,x_{n}).$ (92)

However, in the next step, expanding the effective action in terms of the one-
particle irreducible Green functions $(\Gamma^{(n)})$,

$\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]=\sum^{\infty}_{n=0}\frac{1}{n!}\int
d^{4}x_{1}\ldots
d^{4}x_{n}\overline{\phi_{1}}(x_{1})\ldots\overline{\phi_{1}}(x_{n_{1}})\overline{\phi_{2}}(x_{n_{1}+1})\ldots\overline{\phi_{2}}(x_{n_{2}})\Gamma^{(n)}(x_{1},\ldots,x_{n}).$
(93)

The number of vertices $n_{1}$ and $n_{2}$ are arbitrary but $n_{1}+n_{2}=2n$.
We compute $\Gamma^{(n)}(p_{i}=0)$ which are the diagrams with $2n$ external
lines. Analyzing this in detail with the case of one-loop as shown Fig. 11
which is represented by the formula:

$\displaystyle\Gamma_{n}(\phi_{1c},\phi_{2c})=$ $\displaystyle
i\frac{1}{2n}\int\frac{d^{4}p}{(2\pi)^{4}}\left(\frac{i}{p^{2}-m^{2}(\phi_{1c},\phi_{2c})+i\varepsilon}\right)^{n}(-i\lambda_{\phi_{1}})^{n_{1}/2}(\phi_{1c})^{n_{1}}(-i\lambda_{\phi_{2}})^{n_{2}/2}(\phi_{2c})^{n_{2}}.$
(94)

$\phi_{1}$$\phi_{1}$$\phi_{1}$$\phi_{1}$$\phi_{2}$$\phi_{2}$$\phi_{2}$$\phi_{2}$
Figure 11: The 1-loop diagram with the $\phi_{1}$ and $\phi_{2}$ external
lines

In Eq. (94), each vertex is a factor $-i\lambda_{\phi_{i}}$, the external line
is the factor $\overline{\phi}_{i}=\phi_{ic}=const$. The above integral is
easy to calculate if the external lines are of the same type (i.e., $n_{1}=0$
or $n_{2}=0$). However, when both $n_{1}$ and $n_{2}$ are nonzero, i.e., the
external lines have both fields, calculating the above integral is not simple.
Also then summing with $n=0$ to infinity, it is unlikely that this infinite
sum converges. So in general, 1-loop contributions cannot be represented as

$V_{1}(\phi_{1c},\phi_{2c})=\frac{1}{2}\int\frac{d^{4}p}{(2\pi)^{4}}\log\left[p^{2}+m^{2}(\phi_{1c},\phi_{2c})\right].$
(95)

The representation of one-loop contribution like the above result is only a
stereotype application from the calculation results of the single-field case
(i.e., from the result,
$V_{1}(u_{c})=\frac{1}{2}\int\frac{d^{4}p}{(2\pi)^{4}}\log\left[p^{2}+m^{2}(u_{c})\right]$
in Ref. quiros , $u_{c}$ is VEV of single scalar field). But the computation
of the diagrams has been ignored. Although this is imprecise, if we
unconditionally accept the one-loop contributions as Eq. (95), it is also a
fairly general estimate of one-loop contributions when considering two fields
at once. This is also a possible method today in the context of calculation
Eq. (94). Nevertheless, it is very difficult.

However, this representation is also true, if we interpret the above result as
implying $\phi_{2c}=\kappa\phi_{1c}$, that is, the above result represents
only one field $\phi_{1}$ or $\phi_{2}$; or ignoring all the diagrams where
$\phi_{1}$ and $\phi_{2}$ are present at the same time. Therefore
$V_{1}(\phi_{1c},\phi_{2c})\equiv V_{1}(\phi_{c})$ with
$\phi_{c}^{2}=\phi^{2}_{1c}+\phi^{2}_{2c}$. This is a very good approximation
that eliminates the difficulty of summing diagrams as shown in Fig. 5. This
approximation has been used in calculating decay channels or diagrams as in
Refs. tadashi ; zhen ; abdesslam . Also, there are basis-independent methods
for the two-Higgs-doublet model shinya ; davidson , it is possible to rewrite
2HDM model under one VEV. But in the sections of article, this approximation
is only imprecise in the EWPT analysis.

The next interesting thing here is what if we could rewrite ${\cal
L}\left(\phi_{1}(x),\phi_{2}(x)\right)={\cal L}\left(\phi_{1}(x)\right)+{\cal
L}\left(\phi_{2}(x)\right)$, when expanding $\phi_{i}$ in terms of
$\phi_{ic}$. At that point, the generating functional is rewritten as

$Z[J]=Z[J_{1}].Z[J_{2}],$ (96)

so that,

$W[J]=W[J_{1}]+W[J_{2}].$ (97)

In other words, the effective potential can be separated into two separate
parts:

$\displaystyle\Gamma[\overline{\phi}_{1},\overline{\phi}_{2}]$
$\displaystyle=W[J]-\int d^{4}x\frac{\delta W[J]}{\delta J_{i}(x)}J_{i}(x)$
(98) $\displaystyle=W[J_{1}]-\int d^{4}x\frac{\delta W[J_{1}]}{\delta
J_{1}(x)}J_{1}(x)+W[J_{2}]-\int d^{4}x\frac{\delta W[J_{2}]}{\delta
J_{2}(x)}J_{2}(x)$ (99)
$\displaystyle=\Gamma[\overline{\phi}_{1}]+\Gamma[\overline{\phi}_{2}]$ (100)

This will make summing the diagrams easier, but they are not always separated
like that. The cases of more than two fields are similarly constructed.

## References

  * (1) A. D. Sakharov, JETP Lett. 5, 24 (1967).
  * (2) M. I. Dine, R. G. Leigh, P. Huet, A. Linde, and D. Linde, Phys. Rev. D 46, 550 (1992).
  * (3) K. Kajantie, M. Laine, K. Rummukainen, and M. Shaposhnikov, Phys. Rev. Lett. 77, 2887 (1996).
  * (4) F. Csikor, Z. Fodor, and J. Heitger, Phys. Rev. Lett. 82, 21 (1999).
  * (5) J. Grant and M. Hindmarsh, Phys. Rev. D 64, 016002 (2001).
  * (6) M. D’Onofrio, K. Rummukainen and A. Tranberg, J. High Energy Phys. 08 (2012) 123.
  * (7) M. D’Onofrio, K. Rummukainen, and A. Tranberg, Phys. Rev. Lett. 113, 141602 (2014).
  * (8) Vo Quoc Phong, Vo Thanh Van, and Hoang Ngoc Long, Phys. Rev. D 88, 096009 (2013).
  * (9) A. Menon, D. E. Morrissey, and C. E. M. Wagner, Phys. Rev. D 70, 035005 (2004).
  * (10) S. W. Ham, S. K. Oh, C. M. Kim, E. J. Yoo, and D. Son, Phys. Rev. D 70, 075001 (2004).
  * (11) M. Bastero-Gil, C. Hugonie, S. F. King, D. P. Roy, and S. Vempati, Phys. Lett. B 489, 359 (2000).
  * (12) A. Menon, D. E. Morrissey, and C. E. M. Wagner, Phys. Rev. D 70, 035005 (2004).
  * (13) J. M. Cline, G. Laporte, H. Yamashita and S. Kraml, J. High Energy Phys. 0907 (2009) 040.
  * (14) A. Azatov and M. Vanvlasselaer, J. High Energy Phys. 09 (2020) 085.
  * (15) S. Kanemura, Y. Okada, and E. Senaha, Phys. Lett. B 606, 361 (2005).
  * (16) G. C. Dorsch, S. J. Huber, and J. M. No, J. High Energy Phys. 10 (2013) 029.
  * (17) S. W. Ham, S-A Shim, and S. K. Oh, Phys. Rev. D 81, 055015 (2010).
  * (18) D. Borah and J. M. Cline, Phys. Rev. D 86, 055001 (2012).
  * (19) A. Ahriche and S. Nasri, Phys. Rev. D 85, 093007 (2012).
  * (20) S. Das, P. J. Fox, A. Kumar, and N. Weiner, J. High Energy Phys. 1011 (2010) 108.
  * (21) D. Chung and A. J. Long, Phys. Rev. D 84, 103513 (2011).
  * (22) M. Carena, N. R. Shaha, and C. E. M. Wagner, Phys. Rev. D 85, 036003 (2012).
  * (23) V. Q. Phong, H. N. Long, V. T. Van, and N. C. Thanh, Phys. Rev. D 90, 085019 (2014).
  * (24) J. Sá Borges and R. O.Ramos, Eur. Phys. J. C 76, 344 (2016).
  * (25) V. Q. Phong, H. N. Long, V. T. Van, and L. H. Minh, Eur. Phys. J. C 75, 342 (2015).
  * (26) J. R. Espinosa, T. Konstandin, and F. Riva, Nucl. Phys. B854, 592 (2012).
  * (27) S. Kanemura, E. Senaha, T. Shindou, and T. Yamada, J. High Energy Phys. 1305 (2013) 066.
  * (28) D. J. H. Chung and A. J. Long, Phys. Rev. D 81, 123531 (2010).
  * (29) G. Barenboim and N. Rius, Phys. Rev. D 58, 065010 (1998).
  * (30) S. Profumo, M. J. Ramsey-Musolf, and G. Shaughnessy, J. High Energy Phys. 0708 (2007) 010.
  * (31) S. Profumo, M. J. Ramsey-Musolf, C. L. Wainwright, and P. Winslow, Phys. Rev. D 91, 035018 (2015).
  * (32) D. Curtin, P. Meade, and C-T. Yu, J. High Energy Phys. 11 (2014) 127.
  * (33) M. Jiang, L. Bian, W. Huang, and J. Shu, Phys. Rev. D 93, 065032 (2016).
  * (34) M. Carena, G. Nardini, M. Quiros, and C. E.M. Wagner, Nucl. Phys. B812, 243 (2009).
  * (35) A. Katz, M. Perelstein, M. J. Ramsey-Musolf, and P. Winslow, Phys. Rev. D 92, 095019 (2015).
  * (36) J. Kozaczuk, S. Profumo, L. S. Haskins, C. L. Wainwright, J. High Energy Phys. 1501 (2015) 144.
  * (37) H. H. Patel and M. J. Ramsey-Musolf, Phys. Rev. D. 88, 035013 (2013).
  * (38) N. Blinov, J. Kozaczuk, D. E. Morrissey, and C. Tamarit, Phys. Rev. D 92, 035012 (2015).
  * (39) S. Inoue, G. Ovanesyan, and M. J. Ramsey-Musolf, Phys. Rev. D 93, 015013 (2016).
  * (40) H. H. Patel and M. J. Ramsey-Musolf, J. High Energy Phys. 1107 (2011) 029.
  * (41) G. W. Anderson and L. J. Hall, Phys. Rev. D 45, 2685 (1992).
  * (42) H. H. Patel, M.J. Ramsey-Musolf, M. Garny and T. Konstandin, J. High Energy Phys. 1207 (2012) 189.
  * (43) J. De Vries, M. Postma, and J. van de Vis, J. High Energy Phys. 1904 (2019) 024.
  * (44) J. de Vries, M. Postma, J. van de Vis , G. White, J. High Energy Phys. 1801 (2018) 089.
  * (45) C. Balazs, G. White, and J. Yue, J. High Energy Phys. 1703 (2017) 030.
  * (46) A. Ahriche, Phys. Rev. D 75, 083522 (2007).
  * (47) A. Ahriche, Eur. Phys. J. C 66, 333 (2010).
  * (48) T. A. Chowdhury and S. Nasri, J. High Energy Phys. 1411 (2014) 096.
  * (49) A. Ahriche and S. Nasri, J. Cosmol. Astropart. Phys. 07 (2013) 035.
  * (50) A. Ahriche, G. Faisel, S. Y. Ho, S. Nasri and J. Tandean, Phys. Rev. D 92, 035020 (2015).
  * (51) A. Ahriche, K. L. McDonald and S. Nasri, Phys. Rev. D 92, 095020 (2015).
  * (52) A. Ahriche, S. M. Boucenna and S. Nasri, Phys. Rev. D 93, 075036 (2016).
  * (53) A. Ahriche, K. Hashino, S. Kanemura and S. Nasri, Phys. Lett. B 789, 119 (2019).
  * (54) K. Fuyuto and E. Senaha, Phys. Rev. D 90, 015015 (2014).
  * (55) K. Fuyuto and E. Senaha, Phys. Lett. B 747, 152 (2015).
  * (56) K. Funakubo and E. Senaha, Phys. Rev. D 79, 115024 (2009).
  * (57) M. Spannowsky and C. Tamarit, Phys. Rev. D 95, 015006 (2017).
  * (58) C. Grojean, G. Servant, and J. D. Wells, Phys. Rev. D 71 036001 (2005).
  * (59) C. Delaunay, C. Grojean, J. D. Wells, J. High Energy Phys. 0804 (2008) 029.
  * (60) A. Kusenko, L. Pearce, and L. Yang, Phys. Rev. Lett. 114, 061302 (2015).
  * (61) V. Q. Phong, H. N. Long, V. T. Van, and L. H. Minh, Eur. Phys. J. C 75, 342 (2015).
  * (62) C. W. Chiang and T. Yamada, Phys. Lett. B 735, 295 (2014).
  * (63) V. Q. Phong, N. C. Thao, and H. N. Long, Phys. Rev. D 97, 115008 (2018).
  * (64) J. R. Espinosa, T. Konstandin, J. M. No and M. Quiros, Phys. Rev. D 78, 123528 (2008).
  * (65) D. Comelli and J.R. Espinosa, Phys.Rev. D 55, 6253 (1997).
  * (66) M. Joyce, Phys. Rev. D 55, 1875 (1997).
  * (67) K. Funakubo and E. Senaha, Phys. Rev. D 79, 115024 (2009).
  * (68) M. Dvornikov and V. B. Semikoz, Phys. Rev. D 87, 025023 (2013).
  * (69) T. M. Gould and I. Z. Rothstein, Phys. Rev. D 48, 5917 (1993).
  * (70) S. Braibant, Y. Brihaye and J. Kunz, Int. J. Mod. Phys. A 08, 5563 (1993).
  * (71) Vo Quoc Phong, Phan Hong Khiem, Ngo Phuc Duc Loc, and Hoang Ngoc Long, Phys. Rev. D 101, 116010 (2020).
  * (72) A. Braconi, Mu-Chun Chen, and G. Gaswint, Phys. Rev. D 100, 015032 (2019).
  * (73) I. Baldes, T. Konstandin, and G. Servant, Phys. Lett. B 786, 373 (2018).
  * (74) Vo Quoc Phong, N. T. Tuong, N. C. Thao, and H. N. Long, Phys. Rev. D 99, 015035 (2019).
  * (75) Vo Quoc Phong and Nguyen Minh Anh, Int. J. Mod. Phys. A 34, 1950073 (2019).
  * (76) D. Das, U. K. Dey, and Palash B. Pal, Phys. Rev. D 96, 031701 (2017).
  * (77) D. Das and U. K. Dey, Phys. Rev. D 89, 095025 (2014); 91, 039905(E) (2015).
  * (78) D. Das, U. K. Dey and Palash B. Pal, Phys. Lett. B 753, 315 (2016).
  * (79) R. N. Mohapatra and S. Nussinov, Phys. Rev. D 60, 013002 (1999).
  * (80) T. Kitabayashi and M. Yasue, Phys. Rev. D 67, 015006 (2003).
  * (81) W. Grimus and L. Lavoura, Phys. Lett. B 572, 189 (2003).
  * (82) Y. Koide, Phys. Rev. D 69, 093001 (2004).
  * (83) P. F. Harrison and W.G. Scott, Phys. Lett. B 557, 76 (2003).
  * (84) R. N. Mohapatra, S. Nasri, and Yu. Hai-Bo, Phys. Lett. B 639, 318 (2006).
  * (85) J. Barranco and F. Gonzalez Canales, A. Mondragon, Phys. Rev. D 82, 073010 (2010).
  * (86) F. González Canales, A. Mondragón, M. Mondragón, U.J. Saldaña Salazar, L. Velasco-Sevilla, Phys. Rev. D 88, 096004 (2013).
  * (87) D. Cogollo and J. P. Silva, Phys. Rev. D 93, 095024 (2016).
  * (88) S. Davidson and H. E. Haber, Phys. Rev. D 72, 035004 (2005); 72, 099902(E) (2005).
  * (89) J. Herrero-García, T. Ohlsson, S. Riad and J. Wirén, J. High Energy Phys. 04 (2017) 130.
  * (90) G. C. Branco, P. M. Ferreira, L. Lavoura, M. N. Rebelo, Marc Sher, Joao P. Silva, Phys. Rep. 516, 1 (2012).
  * (91) M. Quiros, Proc. ICTP Summer School in High-Energy Physics and Cosmology (1999), p. 187, arXiv:hep-ph/9901312.
  * (92) A. Arbey, F. Mahmoudi, O. Stal, and T. Stefaniak, Eur. Phys. J. C 78, 182 (2018).
  * (93) W. Mader, Jae-hyeon Park, G. M. Pruna, D. Stöckinger, and A. Straessner, J. High Energy Phys. 1209 (2012) 125.
  * (94) D. Das and Palash B. Pal, Phys. Rev. D 98, 115001 (2018).
  * (95) M. Aoki, T. Komatsu, and H. Shibuya, Prog. Theor. Exp. Phys. 2022, 063B05.
  * (96) V. Khachatryan et al, Phys. Lett. B 758, 296 (2016).
  * (97) D. Gonçalves, A. Kaladharan, and Y. Wu, Phys. Rev. D 105, 095041 (2022).
  * (98) G. C. Dorsch, S. J. Huber, and J. M. No, J. High Energy Phys. 10 (2013) 029.
  * (99) N. Blinov, S. Profumo, and T. Stefaniak, J. Cosmol. Astropart. Phys. 07 (2015) 028.
  * (100) J. O. Andersen, T. Gorda, A. Helset, L. Niemi, T. V. I. Tenkanen, A. Tranberg, A. Vuorinen, and D. J. Weir, Phys. Rev. Lett. 121, 191802 (2018).
  * (101) K. Enomoto, S. Kanemura, and Y. Mura, J. High Energy Phys. 01 (2022) 104.
  * (102) Y. Kajiyama, H. Okada, and K. Yagyu, Nucl. Phys. B887, 358 (2014).
  * (103) E. Ma, Non-Abelian discrete family symmetries for leptons and quarks, arXiv:hep-ph/0409075.
  * (104) E. Ma and B. Melic, Phys. Lett. B 725, 402 (2013).
  * (105) M. E. Carrington, Phys. Rev. D 45, 2933 (1992).
  * (106) D. Curtin and P. Meade, H. Ramani, Eur. Phys. J. C 78, 787 (2018).
  * (107) A. Katz and M. Perelstein, J. High Energy Phys. 07 (2014) 108.
  * (108) H. E. Haber and A. Pomarol, Phys. Lett. B 302, 435 (1993).
  * (109) P. H. Chankowski, M. Krawczyk, and J. Zochowski, Eur. Phys. J. C 11, 661 (1999).
  * (110) J.-M. Gerard and M. Herquet, Phys. Rev. Lett. 98, 251802 (2007).
  * (111) W. Grimus, L. Lavoura, O.M. Ogreid, and P. Osland, Nucl. Phys. B801, 81 (2008).
  * (112) S. de Visscher, J.M. Gerard, M. Herquet, V. Lemaitre, and F. Maltoni, J. High Energy Phys. 0908 (2009) 042.
  * (113) P. Basler, M. Krause, M. M¨uhlleitner, J. Wittbrodta, and A. Wlotzkaa, J. High Energy Phys. 12 (2017) 086.
  * (114) J. Bernon, L. Biana, and Y. Jiang, J. High Energy Phys. 05 (2018) 151.
  * (115) G. C. Dorsch, S. J. Huber, K. Mimasub, and J.M. Nob, J. High Energy Phys. 02 (2017) 121.
  * (116) P. Basler, M. M¨uhlleitnera, and J. Wittbrodt, J. High Energy Phys. 03 (2018) 061.
  * (117) K. Kainulainen, V. Keus, L. Niemi, K. Rummukainen, T. V.I. Tenkanenc, and V. Vaskonend, J. High Energy Phys. 06 (2019) 075.
  * (118) Wei Su, A. G. Williams, and M. Zhang, J. High Energy Phys. 04 (2021) 219.
  * (119) R. Zhou and L. Bian, Phys. Lett. B 829, 137105 (2022).
  * (120) H. T. Hung, T. T. Hong, H. H. Phuong, H. L. T. Mai, and L. T. Hue, Phys. Rev. D 100, 075014 (2019).
  * (121) A. Alves, E. Ramirez Barreto, A. G. Dias, C. A. de S. Pires, F. S. Queiroz, and P. S. Rodrigues da Silva, Phys. Rev. D 84, 115004 (2011).
  * (122) J. Sa Borges and Rudnei O. Ramos, Eur. Phys. J. C 76, 344 (2016).
  * (123) A. Braconi, Mu-Chun Chen, Geoffrey Gaswint, Phys. Rev. D 100, 015032 (2019).
  * (124) S. Kanemura, M. Takeuchi, and K. Yagyu, Phys. Rev. D 105, 115001 (2022).
  * (125) Vo Quoc Phong, Nguyen Chi Thao, and Hoang Ngoc Long, Eur. Phys. J. C 82, 1005 (2022).
  * (126) N. S. Manton, Phys. Rev. D 28, 2019 (1983).
  * (127) F. R. Klinkhamer and N.S. Manton, Phys. Rev. D 30, 2212 (1984).
  * (128) T. Gorda, A. Helset, L. Niemi, T. V. I. Tenkanena, and D. J. Weira, J. High Energy Phys. 02 (2019) 081.
  * (129) T. Kon, T. Nagura, and T. Ueda, K. Yagyu, Phys. Rev. D 99, 095027 (2019).
  * (130) Z. Liu and Pei-Hong Gu, Nucl. Phys. B915, 206 (2017).
  * (131) A. Arhrib, R. Benbrik, Chuan-Hung Chen, R. Guedes, and R. Santos, J. High Energy Phys. 08 (2009) 035.