Anatomy of the Electroweak Phase Transition for Dark Sector induced Baryogenesis

We investigate the electroweak phase transition patterns for a recently proposed baryogenesis model with CP violation originated in the dark sector. The model includes a complex scalar singlet-Higgs boson portal, a $U(1)_l$ gauge lepton symmetry with a $Z^\prime$ gauge boson portal and a fermionic dark matter particle. We find a novel thermal history of the scalar sector, featuring a $Z_2$ breaking singlet vacuum in the early Universe driven by a dark Yukawa coupling, that induces a one-step strongly first order electroweak phase transition. We explore the parameter space that generates the observed matter-antimatter asymmetry and dark matter relic abundance, while being consistent with constraints from electric dipole moment and collider searches, and dark matter direct detection bounds. The complex singlet can be produced via the Higgs portal and decays into Standard Model particles after traveling a certain distance. We explore the reach for long-lived singlet scalars at the {\unit[13]{TeV}} Large Hadron Collider with $\mathcal{L}=\unit[139]{fb^{-1}}$ and show its impact on the parameter space of the model. Setting aside currently unresolved theoretical uncertainties, we estimate the gravitational wave signatures detectable at future observatories.


Introduction
Electroweak baryogenesis (EWBG) is an elegant mechanism that generates the observed baryon asymmetry in the Universe (BAU) at the electroweak phase transition (EWPT) [1,2].The Standard Model (SM), however, can not account for successful EWBG: the Higgs mass is too heavy to allow for a strong first order phase transition and the CP violation source is not sufficient, hence new physics is required.However, new sources of CP violation from particles charged under the SM are strongly constrained by the remarkable results from electric dipole moment (EDM) experiments [3,4].Recently, a new EWBG mechanism in which CP violation occurs in the dark sector and is transmitted to the visible sector via a vector gauge boson Z from a U (1) l gauge lepton symmetry has been proposed [5,6].In such a scenario, a complex scalar singlet S, provides the source of CP violation through its Yukawa coupling to a dark fermion χ charged under U (1) l and leads to a strongly first order electroweak phase transition (SFOEWPT) via its coupling with the Higgs boson field.In this way, the contribution to EDM is suppressed to beyond two-loop level and compatible with current experiments.The dark fermion χ can also serve as an ideal dark matter candidate.
In this work, we study the pattern of the EWPT in the proposed new mechanism [5,6], and investigate the viable parameter space compatible with a successful EWBG, the observed dark matter relic abundance, and phenomenological bounds from dark matter direct detection and collider searches.We also explore potentially observable gravitational wave (GW) signatures.Singlet extensions of the SM addressing the BAU generation and the EWPT have been extensively studied with focus on scalars heavier than half of the Higgs mass [7][8][9][10][11][12].There are also studies on EWPT for light scalars but leaving out the discussion of a complete EWBG model [13,14].Here we compute, both analytically and numerically, the possibility of EWBG for a broad range of scalar masses, from O(1) GeV to a few hundred GeV.More specifically, while a phase transition pattern for EWBG was assumed in [5,6], we now perform a detailed study of the patterns of EWPT for the underlying model.This includes implementing the complete effective finite temperature potential at one loop order with appropriate thermal resummation and performing the nucleation calculation to assure the completeness of the phase transition.
In a broader context of the SM singlet extension, the presence of a dark fermion sector has a particular impact on the scalar potential.The Yukawa term between the singlet S and a dark fermion breaks the Z 2 symmetry (S → −S) explicitly at tree level, and, with a non-zero bare mass of the dark fermion, contributes a tadpole term to the singlet potential at one loop level.To avoid the mixing between the S and the SM Higgs, we introduce counterterms to impose the expectation value of S at the electroweak vacuum at zero temperature to be zero up to one loop order.At finite temperature, thermal corrections from the dark Yukawa coupling will drive the vacuum of S away from zero, leading to a distinct thermal history of the scalar sector: the Universe would go through a one-step first order phase transition from this Z 2 breaking singlet vacuum to the electroweak symmetry breaking vacuum at a lower temperature.The first order phase transition can readily be strong with a tree level barrier between the two vacua, yielding a successful EWBG and detectable GW signals.
Furthermore, the additional singlet scalars have a rich phenomenology, that is to be updated and discussed in this work.If the EW vacuum were to break the Z 2 symmetry along the singlet direction, the singlet scalars would mix with the SM Higgs and thus could be probed by Higgs boson and electroweak precision measurements [15], and also by Higgs exotic decays when the singlet is light [13,14,16].In our study, the scalar potential has a Z 2 symmetry at zero temperature, and, in principle, the singlet would be stable [9] and when below the Higgs decay threshold, it could be probed by Higgs invisible decay searches [13].In our model, however, as the singlet carries U (1) l charge, the heavier singlet can decay via the Z portal to SM leptons promptly.The presence of the dark Yukawa coupling allows the lighter singlet to decay into SM particles through the dark matter χ loop and the Z portal.Thus the decay width of the lighter singlet scalars is generically suppressed by the heavy fermion loop and the small U (1) l gauge coupling constrained by LEP bounds.This leads to distinct signatures for the lighter singlet to decay either in the tracker, muon chamber or outside the detector.We will investigate the reach for long-lived singlet scalars at the 13 TeV LHC with L = 139 fb −1 to probe the parameter space with successful EWBG.
This paper is organized as follows: in Sec. 2, we introduce the scalar potential and review the basic elements of the model, including the specific mechanism of EWBG and its implications for dark matter phenomenology.In Sec. 3, we study the thermal history as well as the critical and nucleation behavior of the EWPT with analytical calculations.We also perform numerical scans for two benchmark scenarios, where viable parameter space for successful SFOEWPT consistent with analytical evaluations is found.In Sec. 4, we discuss the mechanism of baryogenesis and evaluate the parameter space for successful EWBG through numerical scan.In Sec. 5, we re-evaluate the phenomenological discussions for dark matter direct detection, and present an opportunity for long-lived particle searches for the scalar singlets in the newly revealed parameter space.We also compute the GW signature of the model and show potential compatibility with EWBG.We reserve Sec.6 for our conclusions.We collect various technical aspects in the appendices.

The EWBG model
In this section, we briefly introduce the model for EWBG as schematically shown in Fig. (1), with CP violation sourced from the dark sector, and an SFOEWPT where the barrier is provided by the tree-level coupling between the complex singlet and the SM Higgs.In the following, we would present the Lagrangian terms that define the Higgs and gauge boson portals between the SM sector and the dark sector of the theory.Readers can find more specific details of the model in [5,6].In this work, we focus on the EW-scale low SM sector < l a t e x i t s h a 1 _ b a s e 6 4 = " x 8 a b Z 8 b j D t w q u J + i Y P a 2 3 o J Q o k g = " > A A A C N X i c b V B N S w M x E M 3 6 W e t X 1 a O X Y F E 8 l V 0 p 6 L H o x Y t Q 0 V a h u 5 R s O l u D 2 W R J s m J Z t n / A X + N J 0 F / i w Z t 4 transfer particle asymmetry first order phase transition Long-lived particle search

Gravitational wave
Dark matter direct detection z r v z 4 X w 5 3 + P W B W c y s 4 e m 4 P z 8 A o x l p l 0 = < / l a t e x i t > dark sector < l a t e x i t s h a 1 _ b a s e 6 4 = " T X n l y r J 7 y / e 8 q e f < / l a t e x i t >

C, B
< l a t e x i t s h a 1 _ b a s e 6 4 = " q K j e j 0 k R v D q 9 J / g T 3 y B Z / A K 3 q w X 6 8 P 6 t L 5 m r T k r m z k E c 7 C + f w H n 6 6 + + < / l a t e x i t > energy effective theory of the UV complete model for baryogenesis presented in [5], where the U (1) l symmetry, the gauge symmetry promoted from the lepton number, has been spontaneously broken by a heavy singlet field Φ. Integrating out various anomalon fields, the low energy sector of the model contains the SM fields plus the new fields Z µ , S, χ L , χ R . (2.1) The complex singlet S couples to the SM Higgs to provide a barrier for the EWPT.At zero temperature, the tree level effective scalar potential reads where H is the Higgs doublet and Observe that apart from the δV the tree level potential only depends on |S|.Hence, to accommodate a CP violating effect, additional terms depending on S are needed .Such terms can arise from renormalizable, U (1) l invariant terms, involving the heavy singlet field Φ.For this work, we study the case where S has a zero vacuum expectation value (vev) at the EW vacuum, that implies no tadpole at zero temperature.This can be achieved by adjusting the coefficient for the tadpole term to be zero at one-loop order, as will be discussed later.For simplicity, we will turn off the cubic term.The coefficient of the quadratic term, κ 2 S , is a free parameter for scanning and is in general complex.The complex SM SU (2) Higgs doublet H and the complex SM scalar singlet S can be expressed in terms of real fields as follows Using unitary gauge to eliminate the Goldstone fields G i 's, we can rewrite the potential in Eq. (2.2) in terms of (h, s, a) as where . The Higgs quartic coupling and the Higgs vev at T = 0 are fixed by Higgs boson measurements to be λ H = 0.129 and v H = 246 GeV.For convenience, we rewrite the singlet quartic coupling λ S and the coefficient of the CP violating term κ 2 S in Eq. (2.5) in terms of physical parameters, i.e., the masses of the scalars a, s at T = 0, m a,s : As will be shown in Sec.3.2 and Sec.3.4, v S is the value of the singlet field in the a-direction at a zero temperature EW symmetric stationary point of the scalar potential.The existence of this stationary point facilities a SFOEWPT.
The fermion fields χ L and χ R are SM singlets that couple to the singlet composite scalar S. The Yukawa coupling between S and χ is given by χL (m 0 + λ c S)χ R + h.c., (2.7) and the mass of χ reads where λ = |λ c | and θ = arg(λ c )1 .Note that m 0 is induced by the coupling of χ to the heavy singlet scalar field Φ, when it acquires a non-zero vev at the UV.The corresponding Yukawa coupling is assumed to be small such that m 0 is much lighter than the vev Φ and the dark fermions remain dynamical at low energies.We use the freedom of field redefinition to make κ 2 S and m 0 real and positive, leaving λ c complex in general.During the EWPT, the phase of M χ varying in the direction of the expanding bubble wall, can be derived from m 0 , the phases of λ c and the S vev.This is the physical source of CP violation, that will then induce a chiral asymmetry in the χ particles.
As has been mentioned above, the lepton number is promoted to a U (1) l gauge symmetry with an associated Z gauge boson, and the dark fermion χ and the singlet S are assigned certain lepton number charges.Possible anomaly-free UV completions can be found in [17][18][19][20][21][22].The new interactions introduced at low energy are: (2.9) We assume U (1) l with l = e, µ, τ and thus N g = 3, and the charge q ∼ O(1) throughout this study.Given that χ L and χ R carry different U (1) l charges, the chiral asymmetry in the χ sector will give a net U (1) l charge density near the bubble wall, that generates a background for the Z 0 component.Because of the coupling of the SM leptons with the Z , this background further generates a chemical potential for the SM leptons and consequently a net chiral asymmetry for them.As will be discussed in detail later on, solving the corresponding Boltzmann equation, considering the EW sphaleron rate suppressed inside the bubble wall, one obtains a lepton asymmetry.Ultimately, the sphaleron process, which preserves B − L, will generate equal asymmetry in the lepton and baryon sectors to source the observed BAU.
Protected by a Z 2 symmetry in the Lagrangian, χ → −χ, the dark fermion χ is stable and could be a dark matter candidate.The annihilation channels for χ at tree level include annihilating to Z Z , SM lepton pairs, ss, aa and sa.The corresponding Feynman diagrams are shown in Fig. (2).Given the LEP constraints on Z search, the dominant annihilation channel to achieve the correct relic density is χ χ → ss, aa, sa, which requires [6] The range for λ is due to the θ dependence of the annihilation cross section.We will implement this relation in the numerical scanning to generate the observed dark matter relic abundance.To conclude this section, we categorize our model parameters through the following groups: The stability of the scalar potential and the EW thermal history would be affected by the fixed parameters and the scalar potential parameters at tree level, as well as by the dark fermion parameters at loop level.The Z parameters do not enter the scalar potential but would be crucial for the transmission of the baryon asymmetry and thus are relevant for the BAU.They are also of great importance for the phenomenology associated to the Z searches.The scalar singlet phenomenology and the dark matter direct detection constraints involve parameters across the model and will be discussed in detail in later sections.

Anatomy of the electroweak phase transition
In this section, we introduce the one loop order scalar potential, based on which, we analyze the thermal history, the phase transition patterns, and the nucleation requirement.With various boundary conditions, derived at zero and finite temperatures, we identify the viable parameter space that can be compatible with the desired thermal history for EWBG.Towards the end of this section, we show numerical results from parameter scannings to support our analytical calculations for the phase transitions.

The one loop order scalar potential
The one-loop order corrections to the scalar potential in Eq. (2.5) can be calculated via the Coleman-Weinberg (CW) potential in the MS scheme using dimensional regularization as [23] V where Q is the renormalization scale fixed to be the scale of the top quark mass, m B and m F are the field-dependent masses of bosons and fermions [24], n B and n F are the degrees of freedom, c B = 3/2 (5/6) for scalar (vector) bosons, and c F = 3/2.We consider the physical EW vacuum at zero temperature to be: To satisfy the extremum condition at the vacuum at one loop order, we require, with the Higgs mass GeV, and m a,s fixed to be the tree-level relations in Eq. (2.6).We introduce counterterms of the following form: with coefficients fixed by the first and second order derivatives of the CW potential.Note that the Hessian matrix evaluated at the EW vacuum contains an off-diagonal term between s and a due to the dark fermion loop.We would leave out this loop-suppressed off-diagonal term when fixing the counterterms, and thus the small mixing effect in the singlet sector on singlet phenomenology.At a finite temperature T , the one-loop thermal correction to the effective potential is given by [25]: where J B and J F are bosonic and fermionic thermal functions defined as with α = |m| 2 /T 2 .Resummation of higher loop daisy diagrams, that ensures the validity of perturbative expansion near the critical temperature of the phase transition, needs to be included in the full one-loop potential.We employ the Parwani scheme [26] for daisy resummation, replacing the tree-level bosonic squared masses is the self-energy calculated from [25].In our model, the Π i (T ) functions for the fields involved are as follows where g 2 and g 1 are the SM SU (2) and U (1) gauge couplings, and y t is the Yukawa coupling of the top quark with the SM Higgs boson.

Zero temperature boundary conditions
At zero temperature, the Universe should arrive at the physical vacuum, ( h , s , a ) T =0 = (v H , 0, 0), implicating that such a minimum should be the global minimum of the boundedfrom-below (BFB) scalar potential.Such requirements imply several boundary conditions on the zero temperature scalar potential, that we summarize here: 1.The stationary point (v H , 0, 0) is non-tachyonic; On top of requiring a physical Higgs mass m h , this additionally requires non-tachyonic scalar masses m 2 s > 0 and m 2 a > 0. According to Eq. (2.6), this at tree level requires 1 2 while at one loop level, this is guaranteed by the counterterms; 2. The potential is bounded from below; At tree level, this is satisfied as we consider the region where all quartic couplings are positive.At one loop level, the quartic coupling for the singlet receives negative correction from the dark fermion loop.A necessary condition for the one loop potential to be BFB, can be derived by requiring the tree level coupling being larger than the one loop contribution from the dark Yukawa coupling: which can be converted to a lower bound on λ SH via Eq.(2.6): In this work we mainly investigate the low energy part of the model, that is related to the EWPT, around the EW scale.Therefore we check numerically that the potential is BFB up to 10 TeV before the UV sector of the complete model factors in.
3. The physical vacuum at (v H , 0, 0) is the global minimum of the potential; Analysing the stationary points structure of the potential, we start with the treelevel potential given in Eq. (2.5), focusing on the relations among m a , v S and λ SH , neglecting the one loop corrections and hence any constraints on the dark fermion parameters.There are four stationary points (h, s, a) = p i for i = 1, • • • , 4 (considering only the positive solutions due to Z 2 symmetry of the scalar potential at zero temperature), which read For p 1 to be the global minimum, the potential values at p 2 -p 4 must be greater than that at p 1 .Requiring V (p 1 ) < V (p 2 ), we get the following condition: (3.12) H < 0 which is automatically satisfied.One can check analytically that, the necessary condition for the existence of p 4 : The analytical conditions C1 and C2 are derived either at tree level or estimated at loop level.In the numerical calculations we impose the above requirements at one loop level.The benchmarks satisfying all conditions will be the ones used for phase transition and EWBG studies.The analytical conditions C1 and C2 are shown as blue and gray lines, respectively.We see that the boundaries set by the analytical conditions agree well with the numerical results, with a few exceptions caused by loop effects.

Thermal history
As a starting point, to analyse the thermal history, we use the high-temperature expansion of the thermal functions in Eq. (3.6) to the lowest order of α given by Such a leading order approximation is sufficient for analytical purposes, as will be shown in the following: a first order phase transition is guaranteed by a tree-level barrier, and the thermal barrier arising from the next leading order in the high-temperature expansion is negligible for achieving a SFOEWPT.With this expansion, the finite temperature potential can be written as The first derivatives of the effective potential at finite temperature, considering the high temperature expansion in Eq. (3.14), read ) The above equations vanish at the stationary points.For the analytical understanding of the thermal history that follows, we will fix θ = π/2 so that the tadpole coefficients A s = 0 and A a = √ 2m 0 λ/6.With this simplification, one can prove that all local minima have s = 0 at arbitrary temperatures, see calculations in Appendix.B2 .Solving for the minima, the thermal history usually goes through several stages as shown by Fig. ( 4) and described below: 1.At sufficiently high temperatures -well above the EW scale and below the U (1) l symmetry breaking scale, the thermal corrections dominate, and the symmetry is restored for h and s, while a has a non-zero vev that minimizes the effective potential.Hence, the thermal history starts at the vacuum as shown in Fig. (4a).We introduce the notation a S to denote the vev of a at the EW symmetry preserving global minimum.The non-restoration of Z 2 symmetry for a is induced by the non-zero A a from the Yukawa interaction in the dark sector and is an important signature of our model.
2. As the Universe cools down but still with the vev along the h direction being zero, the equation for a can be written in the form of a depressed cubic equation with Here we have replaced (−µ 2 S − 2κ 2 S ) with (m 2 a − 1 2 λ SH v 2 H ) using Eq.(2.6).The number of real solutions to Eq. (3.20) is determined by the sign of its discriminant D(T ) with: If D(T ) is positive, there will be only one real solution, while, if it is negative, there will be three real solutions to Eq. (3.20).
For the model parameter space we investigate, all the quartic couplings are positive, rendering c S > 0. Hence, in the case where λ SH < 2m 2 a /v 2 H , p(T ) is always positive, and so is the discriminant D(T ), implying there will only be one real solution to T T EW : The dark Yukawa coupling breaks the Z 2 symmetry and induces a non-zero vev for the a field at temperatures much higher than T EW ≈ 140 GeV.(b) Below T a , the temperature where two additional solutions for the vev of the a field arise, there are three branches along the direction h = 0, two of which are minima and one is a saddle point, with the global minimum being (0, 0, a S ).When the temperature drops to T h , a new minimum (h EW ,0,a EW ) starts to develop from (0, 0, a 2 ).(c) At the critical temperature, the new minimum (h EW ,0,a EW ) becomes degenerate with (0, 0, a S ). a 2 has transformed into a maximum at this temperature.a 1 could be a minimum or saddle point depending on the relation between the temperature and the model parameters.(d) At T = 0, (v H , 0, 0) becomes the global minimum and the origin becomes a maximum.There are other potential minima or maxima given by Eq. (3.11): a 2 , a 1 and a S in the plot are the values for the singlet vev at the three stationary points p 3 and ∓p 2 respectively.Note that we use −p i to denote the negative solutions.a throughout the thermal history.Because of the Z 2 symmetry of the potential in the Higgs direction, i.e. h → −h, once an EW broken stationary point develops, the stationary point (0, 0, a S ) necessarily becomes a maximum in the h direction, and a roll-over to the EW broken stationary point would happen.Note that the inclusion of the thermal cubic term beyond Eq.(3.14) introduces a thermal barrier and the roll-over is promoted to a first-order phase transition, which however is generically weak given the loop-suppression nature of the barrier.We do not investigate this case further.Thus, is a necessary condition for a SFOEWPT to be induced by a sizable tree-level barrier.This is a weaker constraint than C1, however can be used to down-select the parameter space for the numerical scan.
With Eq. (3.23), p(T ) starts out to be positive at high temperatures, and turns negative at sufficiently low temperatures, rendering a negative discriminant D(T ) and more than one solution to Eq. (3.20).We define T a as the temperature at which T a can be solved analytically.
To forbid a roll-over to the EW broken direction until T a , the stationary point (0, 0, a S ) needs to be a minimum along the h direction, leading to the necessary condition with a S (T a ) = 3 4q(T a ).This condition provides a lower bound on λ SH as a function of v S .Notice that C3 also guarantees that the term inside the bracket of Eq. (3.16) is positive before or at T a and thus no EW broken stationary point is developed.
3. Below T a but above the temperature T h , at which h starts to develop a non-zero vev, two more solutions for a arise, as shown in Fig. (4b).The second derivative vanish at a ± = ± |p| 3 which are the stationary points of the first derivative ∂V ∂a .Thus the three real solutions to Eq. (3.20) where ∂V ∂a = 0, denoted as a 1,2,3 with a 3 > a 2 > a 1 should satisfy a 3 > a + > a 2 > a − > a 1 .This also implies that ∂ 2 V /∂a 2 is positive for both a 1 and a 3 and negative for a 2 .To determine the global minimum, we need to compare the potentials at (0, 0, a 1 ) and (0, 0, a 3 ).For the solutions given above, we find that Thus right below T a , the Universe would be at the (0, 0, a 3 ), identified as (0, 0, a S ), i.e. the global minimum where the EW symmetry is preserved.
4. From Eq. (3.16), in order to allow for a non-zero h, the following condition should be satisfied as temperature drops: We define T h as the temperature at which the left-handed side above equals to zero, below which a non-zero h vev starts to develop.The stationary point (0, 0, a 2 ) is the first one to satisfy the above condition as the Universe cools down and it smoothly transforms into an EW broken stationary point.As the left-hand side of Eq. (3.28) also equals to ∂ 2 V /∂h 2 | (0,0,a 2 ) , T h is also the temperature at which ∂ 2 V /∂h 2 | (0,0,a 2 ) starts to be negative.Thus the stationary point (0, 0, a 2 ) transforms from a minimum along the h direction to a maximum as temperature drops below T h .The new stationary point with non-vanishing h will be denoted as (h EW , 0, a EW ).The second derivative in the h direction is guaranteed to be positive for any EW breaking stationary points given that: We expect (h EW , 0, a EW ) to be the minimum where the Universe settles in after the phase transition.The critical temperature T c is defined as the temperature at which V (h EW , 0, a EW ) = V (0, 0, a S ).The SFOEWPT from (0, 0, a S ) to (h EW , 0, a EW ) (Fig. (4c)) proceeds through bubble nucleation at a slightly lower temperature T n .
Since the start and end local minima are separated by a region with |H| = 0 and |S| = 0, a barrier is guaranteed by the tree level coupling between H and S proportional to λ SH .More details on the nucleation condition are to be discussed in Sec.3.5.
5. Finally, the Universe smoothly deforms into the EW vacuum (v H , 0, 0) from (h EW , 0, a EW ) as the temperature drops to zero, with a EW and a 2 go to zero, shown in Fig. (4d).
Analysing the thermal history, in the case of a SFOEWPT, there should be more than one branch of the stationary point solutions in the singlet direction above the phase transition temperature, as shown in Fig. (4).At the nucleation temperature T n , the false singlet vacuum resides at one branch and tunnels to the true EW vacuum, which develops from another branch.Thus, between the false and the true vacuum, a barrier forms in the region where both |S| and |H| are non-zero.Accordingly, we identify the most relevant condition C3 for a SFOEWPT, i.e.Eq. (3.25).Numerical study to be shown later agrees very well with such a condition, together with zero temperature boundary conditions, C1 and C2 we have derived.Notice that due to the difficulties in solving the cubic equations for the stationary points, we do not present an analytical solution for the critical temperature T c -the temperature at which to impose the SFOEWPT condition.However, for light scalar a with mass m a < m h /2, the coupling λ SH is constrained by experiments to be λ SH 10 −2 .In this scenario, various approximations can be reliably applied on top of the high temperature expansion, which results in an approximate analytic form for T c .With T c estimated this way, we could provide a lower bound for λ SH that agrees better with the numerical results for m a < m h /2.We describe such an approach and the derived conditions in the following.

Critical conditions with light scalars
For the light singlet scalar with m a < m h /2, the minimal requirement from Higgs invisible decays imposes that λ SH 10 −2 .Considering the lower bound for λ SH imposed by Eq. (3.23), m a is bounded to be m a 17 GeV.The orders of λ S and v S consistent with the conditions C1-C3 are estimated to be λ S 10 −2 and v S ∼ 10 3 GeV.Given that we are mainly interested in SFOEWPT, where h EW (T c )/T c > 1, the temperature involved in the following calculation is bounded to be T < v H .We can solve for T c analytically with reliable approximations, by evolving the field vevs from their T = 0 values to higher temperatures.We consider the potential of the form given by Eq. (3.14) to get analytical insights.Based on our analysis above, the SFOEWPT takes place between the two minima (0, 0, a S ) and (h EW , 0, a EW ), which are smooth deformations of the two stationary points at T = 0: p 2 = (0, 0, v S ) and p 1 = (v H , 0, 0) as in Eq. (3.11).The first derivatives of the potential at finite temperatures can be written as Taylor expansions around the stationary points at T = 0. Expanding around the EW vacuum (v H , 0, 0), we have where m 2 h and m 2 a both take their T = 0 values at the EW vacuum, and the last term in both equations sums over all non-negative integers satisfying m + n + p = 3. Solving Eq. (3.30)=0 and keeping only the leading-order terms (see the validity argument in Appendix.C), we get the finite temperature vev of the h field at the EW vacuum Substituting this into Eq.(3.31), the leading-order contribution to the finite temperature vev of a at the EW vacuum is found to be Expanding around the singlet vacuum (0, 0, v S ), the first-order derivative along the a direction with h = s = 0 reads with the last term being a summation over all non-negative integers satisfying n + p = 3. Solving for the roots and keeping the leading-order terms, we get the finite temperature vev of a at the singlet vacuum (See details in Appendix.C) (3.35) The critical temperature T c is when the two stationary points are degenerate which can be written as an expansion around T = 0 to the second order of the field vevs and temperature: (3.37) Substituting the field vevs derived above, we get a quadratic equation of T 2 c : with (3.39) T c can be solved as Similar to condition C3, a necessary condition for SFOEWPT with light singlet scalars at the critical temperature given by Eq. (3.40) reads where a S at T c can be estimated by Eq. (3.35).However, as the temperature dependent term in Eq. (3.35) is found to be negligible for light singlets, we will use a S (T c ) ≈ v S for the following calculations.Both C3 and C3' would give a lower bound for λ SH as a function of v S as shown in Fig. (5), while C3' is expected to give a more restrictive bound since T c is supposed to be lower than T a for a SFOEWPT.

Nucleation
In this section, we derive a semi-analytical condition for nucleation to proceed.The potential must be away from the "thin wall limit" for the tunneling to happen.The thin wall limit refers to the case where the energy difference between the minima of the potential is significantly smaller than the height of the barrier separating the two phases.The tunneling rate is known to be highly suppressed in such cases.Denoting the singlet phase as φ s , the EW phase as φ h and the barrier as φ b , this requirement could be formulated as where O(0.1) < ∆ < 1 can be chosen empirically to be consistent with the numerical results [13].This condition should be checked at the temperature where the phase transition completes.Since this is only an approximate condition, we would check this condition at T = 0 to gain analytical insights.In this case, φ s , φ h and φ b would be zero temperature stationary points p 2 , p 1 and p 4 given in Sec.3.2.The ratio ∆ is chosen to be 0.7 to be consistent with the numerical results.
To summarize, we have derived four necessary conditions (five for m a < m h /2) analytically from requiring T = 0 boundary conditions, first order phase transtion and nucleation.These requirements, which should be checked by numerical calculations at one loop level, will be discussed in the following section.

Numerical scan
We use CosmoTransitions [27] for numerical studies of the phase transitions.This tool starts from the local minima of the scalar potential at zero temperature, tracing the evolution of these minima as temperature increases until it hits some saddle point, which indicates the possible existence of other local minima.At such temperatures, the tool will check the existence of other local minima in the field space close to the saddle point.If it finds one, it will trace back in temperature the evolution of this new minimum until its appearance.With this procedure, the tool can locate all local minima in field space at different temperatures, and thus determine all critical temperatures at which phase transitions may happen.However, the phase transitions identified this way by CosmoTransitions are not physically allowed to happen unless they satisfy the following conditions: • First order phase transition should proceed via bubble nucleation.By requiring that the expectation value for one bubble to nucleate per Hubble volume is ∼ O(1), one can define the bubble nucleation temperature T n at which S 3 (T )/T ∼ 140, with S 3 being the 3-dimensional Euclidean action of the instanton integrated over the bubble wall [27].Around T n , when temperature drops, S 3 (T )/T usually drops, and the nucleation will be more efficient.A first order phase transition can happen only if the nucleation condition given above is satisfied at some point within the phase transition temperature window.
• For a specific starting phase (the global minimum at higher temperature), if there are more than one eligible ending phases for it to transit to, only the transition occurring at the highest temperature can actually happen, since it occurs prior to all others in the cooling-down history of the Universe.
For a successful EWBG, the first order EWPT needs to be strong enough such that the sphaleron rate inside the bubble is suppressed to preserve the produced baryon asymmetry.The strength of the phase transition could be measured by the order parameter h EW (T )/T at the phase transition temperature.We use the following criterion for a SFOEWPT in our numerical scan: where h EW (T ) and h S (T ) are the Higgs vevs inside and outside the bubble at T c or T n .Note that the requirement h S (T )/T < 0.5 is to ensure the sphaleron process is efficient outside the bubble.
We use CosmoTransitions to numerically study the thermal history of the parameter space satisfying zero temperature boundary conditions for the two benchmarks presented in Fig. (3).Fig. (5) shows the parameter space that gives SFOEWPT ( ) as defined by h EW (T c )/T c > 1 as well as the parameter space that satisfies the nucleation condition S 3 /T 140 and SFOEWPT defined by h EW (T n )/T n > 1 (•).We also show the analytical conditions C1-C4 (C3' included for light singlet) in Fig. (5).For the nucleation condition C4, the empirical parameter ∆ is chosen to be 0.7.Our numerical results agree well with the analytical conditions with only a few exceptions on the boundary.This is due to the difference between the scalar potentials used for the analytical and numerical calculations (C2 and C3) as explained in the following, or the approximate nature of the analytical conditions (C3' and C4).For analytical calculations, we use the high temperature expanded thermal potential up to the leading order plus tree-level zero-temperature potential, while for the numerical calculation, we use the full expression for the thermal potential plus one-loop order zero-temperature potential and daisy resummation in the Parwani scheme.

Baryon asymmetry
The baryon asymmetry can be generated during a SFOEWPT when the Universe tunnels from the electroweak symmetric vacuum to the broken one via bubble nucleation.Both vev of the SM Higgs field and the singlet field S change during the phase transition.The real and imaginary parts of the complex singlet across the bubble wall during the phase transition can be modeled by ) where the coefficients are chosen to match the field vevs.For the a(s) field, with a S (s S ) and a EW (s EW ) being the vevs in the EW symmetric and broken phases at the EWPT respectively, we have and analogously for the s field.For general θ values, a S (s S ) and a EW (s EW ) are both non-zero.L w is the characteristic scale of the bubble wall width, and we will set it to be a random number in the range (1/T n , 10/T n ).The dark fermion mass can thus be written with explicit spatial coordinate dependence in the rest frame of the bubble wall as The Yukawa interaction between χ L,R particle and the S background contributes to the CP violating (CPV) source calculated as [28] Non-vanishing S CPV requires non-vanishing arg(M χ (z)) .As discussed in Appendix.B, θ = π/2 leads to s S = s EW = 0 and thus M χ (z) is real according to Eq. (4.4), rendering S CPV = 0.With m 0 = 0, both tadpole coefficients A a and A s vanish, which also gives s S = s EW = 0, rendering arg(M χ (z)) = θ + π/2 be a constant and thus the S CPV vanishes.Therefore, a non-vanishing S CPV requires θ = π/2, non-vanishing m 0 and λ.Non-vanishing S CPV also implicitly requires non-vanishing κ 2 S as discussed in Sec. 2. This can also be seen from Eq. (3.17) and Eq.(3.18), which show that κ 2 S = 0 implies arg(S) = −θ and thus M χ (z) is a real number.
We comment on the bubble wall velocity we use, which may introduce theoretical uncertainties to the calculation [29][30][31][32][33].The bubble wall velocity, crucial for the calculation of BAU, is in principle determined by the phase transition dynamics.However, the calculation for the bubble wall velocity is convoluted due to the non-equilibrium nature of first order phase transitions.In our calculation for BAU, we take a simplified approach that is commonly used in phenomenological studies: treat v w as a free parameter and truncate the Boltzmann equation to the leading moments [28].This approach only applies to subsonic bubble wall velocities, v w < 1/ √ 3, up to theoretical uncertainties.In the following, we shall choose v w to be a random number in the range (0, 0.5).
Following [5], we remind the readers on the generation of the BAU.The CPV source leads to nonzero particle chiral asymmetries in the dark sector defined as: ) with n i being the number density of the corresponding particle.As the sum ξ χ L (z) + ξ χ R (z) vanishes due to the conservation of the global U (1) χ symmetry, we only need to consider the evolution of ξ χ L (z) according to the diffusion equation with D being the diffusion constant and Γ m the transport rate.The solution to this diffusion equation is given by where G(z) is the Green's function.The chiral asymmetries imply a net U (1) l charge density near the bubble wall as which will yield a Coulomb background of the Z potential As the singlet vev along the bubble wall is small compared to the U (1) l symmetry breaking sacle, the change of M Z along the bubble wall is negligible.This Z background effectively acts as a chemical potential µ L L (z) = µ l R (z) = g Z 0 (z) for the SM leptons and sources the net chiral asymmetry in the SM lepton sector with its thermal equilibrium value given by In the presence of sphaleron, which would change lepton and baryon numbers while preserving the SM B − L, equal asymmetries would be generated for the SM lepton and baryon numbers, ∆n B = ∆n L L .These asymmetries would evolve towards their equilibrium values following the rate equation given by The sphaleron rate Γ sph at nucleation temperature is considered to be unsuppressed outside the bubble, and exponentially suppressed inside the bubble: with Γ 0 10 −6 T n , M sph = 4πh EW (T n )B/g 2 being the sphaleron mass inside the bubble, and B 1.96 a fudge factor depending on the Higgs mass and g 2 [34].The solution to the rate equation is given by The observed baryon asymmetry is quantified by the baryon-to-entropy ratio and is measured to be Having reviewed the dependence of the baryogenesis mechanism on the model parameters, we discuss the large-scale numerical scan we performed to identify the parameter space that produces the observed baryon asymmetry in Eq. (4.16).We require zero temperature boundary conditions, and use CosmoTransitions to identify the parameter space compatible with nucleation and a SFOEWPT (v n /T n > 1).In addition, we require the χ dark matter candidate to yield the observed relic density.The model parameters are highly correlated and restricted by the above requirements and the Z search bounds as summarized below: 1. m a : The upper bound for m a is around 600 GeV, which is set by Eq. (3.23) in addition with the perturbativity requirement λ SH < 4π.
2. m 0 : The dark fermion has to be at least heavier than 2m s to open up the annihilation channel to a pair of singlet scalars.The upper bound is chosen to be 1 TeV to avoid the decoupling of the dark fermion from the thermal history.
3. λ would be chosen to be in the range given by Eq. (2.10) to give the observed dark matter relic density.
4. m s : Its difference to m a is chosen to be in the range (5,100) GeV to enable the annihilation channel of dark matter to both of the two singlet scalars.

v S :
The upper bound for v S comes from the requirement that the lower bound for λ SH from C1 must be lower than the upper bound for λ SH from C2, which is found to be (4.17) The upper bound for v S is smaller for larger value of λ S,BFB , which occurs for larger dark fermion mass m 0 and stronger dark Yukawa coupling λ.With heavier singlets, larger m 0 and thus larger λ are required to obtain the right relic density, hence v S is constrained to be smaller.In the numerical scan, the range for v S is chosen to be (1, 2.5) TeV for light singlets with m a < m h /2, and (10, 700) GeV for heavy singlets with m a,s > m h /2.
6. λ SH : The BFB condition C1 sets a lower bound on λ SH , while the global minimum condition C2 sets an upper bound on it.Note that the upper bound is a tree level bound and would be less likely to be fulfilled if the dark fermion loop becomes sizable.
7. g : This is chosen to give the observed baryon asymmetry, that is parametrically proportional to ∼ g 2 /M 2 Z .
8. M Z : The parameter region for M Z < 10 GeV is highly constrained by Kaon and B meson decay [35,36], which we will not consider here.On the other hand, for light m 0 , say m 0 < M Z /2, larger g is required to achieve the observed BAU, which would be excluded entirely by the Z searches [5].Thus we focus on the case with large m 0 , especially m 0 > M Z /2.
The points giving successful nucleation, observed baryon asymmetry and relic abundance are shown together with the experimental constraints set by LEP [37,38] in Fig. (6).Gray points in this plot are excluded by LEP constraints, while red and cyan points are not.Cyan points are excluded by the dark matter direct detection bounds as will be discussed in the next section.This plot shows a linear correlation between the logs of the two parameters as expected from the expression for baryon asymmetry, which is proportional to ∼ g 2 /M 2 Z .

Phenomenology
In this section, for the parameter space compatible with the EWBG and the dark matter relic abundance, we update the phenomenology on the dark matter direct detection bound on the dark fermion, discuss the search for the singlet scalars at the collider, and study the gravitational wave signatures of the SFOEWPT .We refer the readers to [5,6] for a detailed discussion on the Z search and the EDM, while relevant bounds have been applied for the numerical scan.

DM direct detection
The most stringent bound from dark matter direct detections for the dark matter mass in our model comes from nuclear recoil experiments.The nuclear recoil of dark matter occurs through Z or scalar exchange at one loop order as shown by Fig. (7) (left).The cross section of the scattering process through Z exchange is given by where α = g 2 /(4π), µ p = m 0 m p /(m 0 + m p ).The f function is the loop factor given by (5.2) where q 2 = −4µ 2 v 2 is the square momentum transfer, with v 10 −3 being the typical halo dark matter velocity, µ = m 0 m Xe /(m 0 + m Xe ), and Λ being the renormalization scale chosen to be 1 TeV.In the limit q 2 → 0, which is generically true for our parameter region, The scattering process through h exchange shown in Fig. (7) (right) could also be important as it involves dark matter scattering off both protons and neutrons, and λ, λ SH can be sizable compared to the g coupling involved in the Z exchange channel.The cross section of the scattering process through h exchange in the q → 0 limit is given by with the subscript n in µ n representing neutron or proton, and (5.4) The form factor y n = 0.3m n /v H for the Higgs coupling to nucleon is taken from [39].In the numerical calculations, we also take into account the interference term between the two channels in Fig. (7).We show benchmark points (red) on the m 0 -g /M Z plane in Fig. (8) that produce the observed baryon asymmetry, and satisfy LEP constraints, and pass the current DM direct detection bound from XENON1T [40].We also show the upper bound that passes DM direct detection when considering only the Z channel (blue line) in Fig. (8).According to Eq. (5.1), the direct detection cross section depends on the model parameters g and M Z and is proportional to g 4 /M 4 Z in the low q limit.As dark matter mass m 0 increases from 100 GeV to 1 TeV, its number density decreases by a factor of 10, and the experimental constraints on the direct detection cross section is weakened by one order of magnitude, corresponding to a weaker limit on g /M Z by about a factor of 2. For cyan points below the blue line, the contribution to DM direct detection is dominated by the h exchange channel, where λ SH can be sizable.This happens at the larger m 0 region, where χ can annihilate to singlet scalars with large m s,a and sizable λ SH .Future dark matter direct detection from XENONnT will probe the regions with cross section roughly two orders smaller [41].This will probe most of the points in Fig. (8) where g /M Z is smaller and the dominant contribution to dark matter direct detection comes from the h exchange channel.

Singlet scalar searches at LHC
In this section, we explore the collider searches for the new singlet S in our model.The physical states from the S field are its real part s and imaginary part a, with a lighter than s.The new scalars s and a can be produced via the Higgs portal.As S has zero vev at zero temperature, both s and a have to be pair-produced via an on-shell or off-shell Higgs boson.In the low energy effective field theory, the coupling between S and Z can be UV model dependent.In the simple case where S carries U (1) l charge, the real singlet s can decay to aZ with Z further decays to SM leptons.Z can be on-shell or off-shell depending on the mass difference between s and a.We find that s has a decay length mostly smaller than 0.01 cm for our parameter space.The signature of s singlets via Higgs production thus contains four leptons and two a singlets in the final states, which we probe with four prompt leptons final state.
The light singlet a can decay via the dark matter χ loop to on-shell/off-shell Z and then further decay to multiple SM leptons, as shown in Fig. (9) (right).As the decay width of this channel is suppressed by both the heavy χ loop and the small g coupling, a can Figure 8: Benchmark points producing the observed BAU, relic abundance on the plane of the combined dark gauge parameter ratio 10 4 g /M Z vs the dark fermion mass parameter m 0 .Red points satisfy dark matter direct detection constraints, while cyan points do not.The dark gauge parameters g and M Z are calculated to give the observed baryon asymmetry, which is proportional to ∼ g 2 /M 2 Z .The Z exchange channel dominates the direct detection cross section in most of the parameter space except for the region where m 0 is heavy and g /M Z is small.The bound for g /M Z (blue line) is derived from the Z exchange channel cross section, which is proportional to g 4 /M Z at the low energy limit.This bound weakens as the dark matter mass increases, since the direct detection constraint is weaker for heavier dark matter.When m s,a < m h /2, the multi-lepton production cross section at LHC via the singlet scalars is given by σ(gg → ss(aa) → l + l − + X) ∼ σ(gg → h) × Br(h → ss(aa)), while when cτ[cm] 1 100 10 000 100 000 m s,a ≥ m h /2, the corresponding cross sectin σ(gg → ss(aa) → l + l − + X) ∼ σ(gg → ss(aa)) [5].We neglect the branching ratio Br(Z → l + l − ) ∼ O(1) for any leptons.In Fig. (10) (right), we present the theoretical estimation of σ(gg → aa → l + l − + X) for all the points that pass all our considerations discussed in the previous sections, as a function of the singlet scalar mass m a .
For the prompt-decay regime (cτ < 0.1 mm) with m s,a < m h /2, the bounds on cross section are obtained from the limits at 13 TeV on Br(h → 4l) = 2 × 10 −5 for singlet s with 139 fb −1 data [16], while for the singlet a, searches for 5l final states constrains Br(h → 5l) 10 −5 with 35.9 fb −1 data [42].We further apply these bounds to the parameter space where m s,a > m h /2, i.e. singlet s or a is produced via off-shell SM Higgs as an approximation to obtain the points passing collider searches.
In the case of long-lived singlet a and s, we utilize existing searches [42][43][44][45] for long-lived particles decaying into a pair of electrons or muons, where benchmarks with on-shell Higgs production at different Higgs mass values have been considered.For the lifetime between 0.1 − 10 4 mm, a and s can be constrained by the searches for lepton pairs with displaced vertices, while for longer lifetimes they can be constrained by searches for muon pairs reconstructed in the muon chamber (see [46] and references therein).Though bounds on the cross section σ(gg → ss(aa) → l + l − + X) depend on m s (m a ) and the Higgs boson being on-shell or off-shell, we take the bounds for our parameter regions, as listed in Tab.(1), to be the ones that are satisfied by all benchmarks studied in [42][43][44][45].This treatment would result in more stringent bounds for our case due to the fact that our final state leptons can be less energetic than those considered in the benchmark searches.For scalars decaying outside the detector region (cτ > km), the bounds are obtained from Br(h → invisible) = 0.145 with 139 fb −1 data [47].
We will project these limits to evaluate the reach of the 13 TeV LHC with luminosity L = 139 fb −1 and probe our parameter space.The projection is performed as follows: For the parameter region with 0.1 mm < cτ < km where the SM background is negligible, the limit on cross sections roughly scale as 1/L.Also with negligible SM backgrounds, the limit on cross section barely depends on the energy of the collisions, hence we can take the results from 8 TeV to be approximately valid for 13 TeV and only scale the luminosity.For the parameter region with singlet scalars decaying promptly (cτ < 0.1 mm) or outside the detector (cτ > km), the limits on cross section should scale as 1/ √ L in addition to the energy-dependence from the background cross section.In this case, we consider directly limits from 13 TeV and scale them with the luminosity.The cross section values that can be probed at 95% confidence level at the 13 TeV LHC with luminosity L = 139 fb −1 are summarized in Tab.(1), for various values of scalar masses and lifetimes, which are comparable to the bounds recently reviewed in [16] for the scalar singlet with its mass below m h /2.With more dedicated searches with displaced multi-lepton final states, these bounds can be improved.We leave this study for the future.
Taking both a and s scalar singlet into consideration, the points that can be probed at 13  • m s,a < m h /2, the scalar s can be pair-produced from an on-shell Higgs boson and then promptly decays to 4l + X.The parameter space will be constrained by the Higgs exotic decay to multi-leptons with its branching ratio constrained to be smaller than 2 × 10 −5 [16].The small open diamonds to the left of the dashed line in Fig. (10)(right) fall in this case.Given that: with Γ h = 4.07 MeV being the SM Higgs total width, λ SH is thus constrained to be λ SH 10 −4 by Higgs exotic decays.
• m s > m h /2 and m a < m h /2, the singlet a is produced via on-shell Higgs while s via off-shell Higgs.In this scenario, a usually decays outside the detector and the parameter space will be constrained by Higgs invisible decay, which requires λ SH 0.01 (σ < 6380 fb).The solid circle to the left of the dashed line in Fig. (10)(right) fall in this case.Searches for singlet s prompt decay give weaker bounds, e.g.λ SH 0.07 for m s ∼ 70 GeV, with larger upper bound on λ SH for larger m s .
• m s,a > m h /2, both singlet a and s are produced via off-shell Higgs.Searching for s can be efficient in probing the parameter regions due to the fact that constraints are stronger for particles with shorter lifetimes as observed from Tab. (1).However, as a is lighter than s, its production cross section can be much larger, which makes the  [42][43][44][45].The first number in the first row of each column is the bound obtained by searches for 5 leptons [42] which could be applicable to singlet a searches if a decays promptly, while the second number is the bound from Higgs exotic decay [16] to four leptons, applicable to singlet s searches if s decays promptly.
a The bounds are approximated as the ones obtained for on-shell SM Higgs (ma,s < m h /2).This is a conservative way to obtain all points passing collider searches.
searches for long-lived a particles more efficient.We present the parameter space that can only be probed by the singlet a search with larger open diamonds.
We can improve on probes of our parameter space through long-lived particle searches by using the detector in a more creative way, e.g.including the timing information [48,49] or considering new LHC auxiliary detectors [50].With an expected luminosity increase by roughly a factor of 20, HL-LHC will probe cross sections a factor of 4-20 smaller, allowing us to probe a broader range of parameter space.For example, consider the searches for singlet a decaying outside the detector, while LHC requires Br(h → invisible) < 0.145, HL-LHC will constrain Br(h → invisible) < 0.031, probing all the solid circles with m a < m h /2 in Fig. (10)(right).

Gravitational wave signature
The stochastic gravitational wave generated during an SFOEWPT could become detectable in current and future GW detectors [51].There are three main sources for GWs from a cosmological first-order phase transition: Collisions of the bubble walls and shocks in the plasma, sound waves in the plasma after the bubble collisions and before the kinetic energy gets dissipated by bubble expansion, and the magnetohydrodynamic (MHD) turbulence in the plasma after bubble collisions.The strengths of these three sources depend on the dynamics of the phase transition, especially the bubble wall velocity v w .We focus on the non-runaway bubble wall scenario with a subsonic terminal velocity v w < 1/ √ 3 for GW calculation, which is compatible with our BAU calculation.In this scenario, the main contributions to GW signals come from the bulk motion of the fluid, since the energy in the scalar field is negligible.The total power spectrum is a linear combination of contributions from sound waves and MHD turbulence: with h being the reduced Hubble constant satisfying H 0 = h × 100km/s/Mpc.Calculations for the GW spectra from the above two sources are presented in [51] as a function of several key parameters.Here we comment on such parameters, with T * denoting the temperature of the thermal bath when the GW is produced.
1.The fraction β/H * , with β being the inverse time duration of the phase transition and H * being the Hubble constant at temperature T * , can be evaluated as T * is usually taken to be T n for phase transitions without significant supercooling, as is the case for our parameter space.
2. The ratio of the vacuum energy density released during the phase transition to that of the radiation bath, denoted as α, is evaluated as where φ EW and φ S denote the phases inside and outside the bubble wall at the time of the phase transition, and ρ * rad = g * π 2 T 4 * /30, with g * = 106.75being the relativistic degrees of freedom in the plasma at T * ; 3. The fraction of the total vacuum energy released during the phase transition converted into the bulk motion of the fluid, is given by For subsonic bubble walls with v w < 0.5, the second expression applies.κ v is used to calculate the GW signals from sound waves.The same quantity but for MHD turbulence κ turb is given by with typically in the range of 5% − 10% [51,52].We use = 0.05 for a conservative evaluation.
In Fig. (11), we show the GW signals calculated for the benchmark points satisfying all the considerations discussed in previous sections and the comparison to the sensitivities of various proposed GW detectors covering the relevant frequency range [53]: LISA, DECIGO, BBO, Einstein Telescope (ET), MAGIS-100 and MAGIS-Space [54,55], and AEDGE [56].The peaks of our GW signals occur between 10 −3 − 10 Hz, which can be covered by LISA, AEDGE, DECIGO and BBO.GW signals for benchmarks shown in the plot are strong enough to be observed by these detectors based on evaluations presented above.Note that our approach of treating the bubble wall velocity as a free parameter would introduce uncertainties to the GW signature calculations, as it should be determined by the specific phase transition dynamics.There are also alternative calculations [57][58][59] which takes into account the expansion of the Universe and the finite lifetime for the sound waves.Such calculations in general yield weaker signal strengths and lower peak frequencies than our current approach [51].More investigation addressing such theoretical uncertainties is needed to be conclusive, which we will leave for future studies.

Conclusions
In this work we focus on studying the anatomy of the electroweak phase transition for a recently proposed novel mechanism for electroweak baryogenesis, where the new required source for CP violation resides in a dark sector.Introducing dark CP violation for a successful EWBG evades the stringent constraints imposed by measurements of electron and neutron electric dipole moments.This new EWBG mechanism involves a new fermionic particle χ and a complex scalar S in the dark sector that together can generate a non-vanishing CPV source.The global lepton number symmetry U (1) l , with l = e + µ + τ , needs to be promoted to gauge symmetry, with the associated Z gauge boson that acts as a portal to transfer the CPV source from the dark sector to the SM sector.More in detail, as the Universe undergoes a strong first order phase transition from an EW preserving vacuum to the EW breaking one, the singlet vev is changing along the bubble wall and generates a varying non-vanishing imaginary component for the dark fermion mass, that creates a CPV induced non-zero chiral asymmetry for the dark fermion χ .This chiral asymmetry is transferred to a chiral asymmetry in the SM lepton sector via the Z portal.The EW sphaleron process in the SM can convert the chiral asymmetry in the SM lepton sector to baryon asymmetry which will be preserved once entering the bubble of EW breaking vacuum.The Higgs portal, connecting the SM sector with the dark sector via its interactions with the singlet complex scalar, allows for a SFOEWPT.Interestingly enough, the fermionic particle χ introduced can also serve as a dark matter candidate.In [6], the generated baryon asymmetry at the EWPT was calculated, assuming a sufficiently strong EWPT and that nucleation takes place.Here we perform a meticulous analysis of the phase transition dynamics.
Exploring the phase transitions, we find distinct thermal histories for the scalar sector: The Yukawa coupling between χ and S, together with the bare mass for χ, breaks the Z 2 symmetry (S → −S) explicitly at tree level and introduces a tadpole term for S at finite temperature.This generically leads to stationary points with non-vanishing vev of S at high temperatures.Later on, the Universe goes through a one-step SFOEWPT to converge to the EW symmetry breaking vacuum.The EWPT proceeding in this way is generically strong due to the barrier between the two vacua induced by the quartic coupling between the complex singlet and the SM Higgs boson.We derive two zero temperature boundary conditions to secure that the potential is bounded from below and the EW breaking vacuum is the global minimum.Furthermore, we derive two necessary conditions to achieve SFOEWPT and successful nucleation.The parameter space constrained analytically agrees well with that found by numerical calculations using CosmoTransitions.
Once we identify the parameter space with successful SFOEWPT, nucleation and dark matter relic abundance for the dark fermion, we evaluate the parameter space that can generate the observed baryon asymmetry while satisfying the current LEP bound from Z searches.We obtain values of g , the U (1) l gauge coupling, compatible with the current LEP bounds for Z masses between 10 GeV and 1 TeV.We then examine this scenario via dark matter direct detection searches.Current limit from XENON1T constrains g /M Z 10 −4 GeV −1 for dark matter masses between 10 GeV to 1TeV.While most of the parameter space is constrained through dark matter scattering with nucleons via the Z channel, there are regions, with sufficiently small g , for which the h channel becomes important.Smaller scattering cross sections governed by h channel will be probed by future direct detection experiments, e.g., XENONnT.In this analysis, we assume χ to be the total dark matter relic density, while assuming sub-component dark matter can relax the direct detection bounds opening up more regions of parameter space.
Another interesting probe of this scenario is via searches for the singlet scalars at LHC.The singlet scalars s and a (real and imaginary components of the complex singlet S) can be produced via the Higgs portal and subsequently decay to SM leptons through the Z portal at tree level (s) or via a dark fermion loop (a).For the singlet a, with the decay width suppressed by both the heavy χ loop and the small g coupling, a can be long-lived particles.We calculate the decay lifetime of the singlet scalars for the parameter space leading to successful EWBG and other phenomenological considerations.In the absence of dedicated displaced multi-lepton searches, we utilize previous constraints from displaced di-lepton searches to estimate the regions of parameter space probed by current data at the 13 TeV LHC with luminosity L = 139 fb −1 .We find that the searches for s and a can be complementary to constrain our parameter space: while s can decay relatively faster than a, the production cross section for a can be much larger than that for s.Thus our parameter space with s not so heavy are mostly constrained by searches for s, while with heavy s, the parameter space can instead be constrained by searches for a singlet.The projection to HL-LHC can be performed straightforwardly, which can probe a broader parameter space for cross sections a factor of 4 to 20 smaller.We expect that powerful probes for long lived scalars will come from dedicated LHC searches for displaced multi-lepton signals and we encourage the LHC collaborations to look into these topologies.
We further evaluate the gravitational wave spectra for the parameter space compatible with the observed BAU, dark matter relic abundance and all other phenomenological constraints.The peak of the GW signal for benchmark points arises between 10 −3 − 10 Hz, and the power spectrum is strong enough to be observed by LISA, AEDGE, DECIGO and BBO.We note that there are various theoretical uncertainties associated with the evaluation of the GW signal and the BAU, which need to be further addressed to be conclusive on the simultaneous realization of the EWBG mechanism and strong detectable GW signals.
Summarizing, we have performed a complete analysis, including the anatomy of a sufficiently strong first order EWPT, for a novel mechanism of EWBG with CPV in the dark sector.We show its capability to explain the observed BAU and the dark matter relic abundance.The model can be tested by future dark matter direct detection experiments, dedicated search for long-lived particles at the LHC and HL-LHC, as well as future GW laboratory probes.
For the case with θ = 0 and thus vanishing A a , s = 0 and a = 0 can be a local stationary point.Following a similar calculation, at such a point, The determinant is Det = 8λ S a 2 κ 2 S , (B.5) which is positive for positive λ S and κ 2 S .Together with positive diagonal terms, this indicates that a solution with s = 0 and a = 0 could be a valid minimum in this case.At finite tempreatures, the Universe would reside at this CP breaking vacuum if such a stationary point develops to be a global minimum of the potential.One observes that S CPV = 0 according to Eq. (4.5) and successful EWBG can be realized based on our numerical study.With θ = 0, thus all model parameters real, such finite temperature CP violation is spontaneous and has no zero temperature residue.
In the case with m 0 = 0, both tadpole coefficients vanish.Following the discussion for vanishing A s , we know the local minima of the potential should have s = 0.

C Critical temperature with light singlet scalars
From the Taylor expansion of the first derivative of the potential around the physical minimum at T = 0, we can derive the expression for h EW (T ) given by Eq. (3.32).For the light singlet scalar with mass m a < m h /2, the Higgs and singlet mixing factor is constrained by the singlet searches at colliders, with the minimal requirement λ SH 10 −2 from Higgs invisible decays.The orders of other parameters in Eq. (3.32) are estimated to be: v H ∼ 10 2 GeV, v S ∼ 10 3 GeV, λ H ∼ 10 −1 , λ S ∼ 10 −2 , c H ∼ 10 −1 , c S ∼ 10 −2 and A a ∼ 10 GeV, which would be consistent with the conditions C1-C3.Solving Eq. (3.30)=0, we get the finite temperature vev of the h field at the EW vacuum: where The fourth-order derivatives ∂ 4  i V are determined by the quartic couplings and the thermal coefficients, among which only the following terms do not vanish We can estimate the orders of ∆x i in Eq. (C.1) based on the known information: (a) SFOEWPT condition, h EW (T c )/T c > 1, indicating that the temperature involved in our calculation satisfies T < v H ; (b) m h m a , indicating that the barrier around the EW vacuum is much higher in the h direction than in the a direction, and thus ∆h ∆a ∼ T c .Based on these estimations, we see that the leading order term among the three fractional terms under the square root in Eq. (C.1) is are estimated to be of the same order following similar arguments for h EW (T ), and thus the second fractional term in Eq. (C.3) is ∼ O(1) and taken to be 1 for leading order estimation given by Eq. (3.33).
Similar estimations can be performed to a S (T ), which is obtained by solving Eq. (3.34)=0 as given below ) n T p , and ∂ 4 i V is of order 10 −2 .The last fractional term under the square root is suppressed by the small quartic coupling and (∆x i /v S ) 3  1, and thus is sub-leading compared to the other terms.Eq. (3.35) is obtained by keeping only the leading order terms in the expansion.
9 e z N d N u D V R 8 E H u / N Z G Z e m H C m j e u + O j O z c / M L i 6 W l 8 v L K 6 t p 6 Z W O z r W W q KL S o 5 F J d h 0 Q D Z w J a h h k O 1 4 k C E o c c r s L b k 5 F / d Q d K M y k u z S C B I C Z 9 w S J G i b F S t 7 J b 3 h s O M 7 / 4 K J O K i D 7 k 2 F d x d n E 2 1 E C N V H n u + 9 1 K 1 a 2 5 B f B f 4 k 1 I F U 3 Q 7 F a + / J 6 k a Q z C U E 6 0 7 n h u Y o K M K M M o h 7 z s p x o S Q m 9 J H z q W C h K D D r J i i x z v W q W H I 6 n s E w Y X 6 s + O j M R a D + L Q V s b E 3 O j f 3 k j 8 z + u k J j o K M i a S 1 I C g 4 0 F R y r G R e J Q N 7 j F l T + Y D S w h V z O 6 K 6 Q 1 R h B q b 4 N Q U S g Q F P n V H d j9 e 3 2 b l / U 7 m L 2 k f 1 L x 6 r X 5 e r z a O J 6 m V 0 D b a Q f v I Q 4 e o g U 5 R E 7 U Q R Q / o E T 2 j F + f J e X P e n Y 9 x 6 Y w z 6 d l C U 3 A + v w H M r q 2 7 < / l a t e x i t > Z 0 < l a t e x i t s h a 1 _ b a s e 6 4 = " t K S 4 A f 8 w m R 5 b 1 f 4 8 C p r t Q S + m s y g = " > A A A C r 3 i c b V H b i h N B E O 0 Z b 2 t 0 N a v g i w 8 W B n G F J M x I x A U J L P r i 4 w p m s z g z x p p O T d K k L 0 N 3 j 5 A N + Q H / 0 L / w E + x k 8 2 B 2 L W g 4 f U 5 d u k + V t R T O J 8 n v K L 5 1 + 8 7 d e w f 3 W w 8 e H j 5 6 3 D 5 6 c u 5 M Y z m N u J H G X p T o S A p N I y + 8 p I v a E q p S 0 r h c f N r o 4 5 9 k n T D 6 q 1 r b n V m l j c 2 t 7 p 7 x b 2 d s / O D y q H p 9 0 V B h L T N o 4 Z K H s + U g R R g V p a 6 o Z 6 U W S I O 4 z 0 v U n 9 1 m 9 + 0 y k o q F 4 0 t O I e B y N B B 1 S j L S R B t X a 5 W y W u P m e x G c I T 1 L o S p 4 E S E 5 s D b 5 L b l g J D B h W S 7 O a g r r p H N d d x r 1 x m O j 1 r w r E i q D M 3 A O r o A D b k A T P I A W a A M M X s A b e A c f 1 q v 1 a X 1 Z 3 4 v W k l X M n I I l W P N f K w e k 4 g = = < / l a t e x i t > CP < l a t e x i t s h a 1 _ b a s e 6 4 = " D C s 8 4 W 9 p W O o I Z 3 o 0 s 9 w 5 6 T J w r Y e u H A 4 5 1 7 u v c e P G F X a c b 6 s l d W 1 9 Y 3 N 0 l Z 5 e 2 d 3 b 9 8 + O G y r M J a Y t H D I Q t n 1 k S K M C t L S V D P S j S R B 3 G e k 4 0 8 a m d 9 5 I F

Figure 1 :
Figure 1: Low energy sector of the model and probes for EWBG and a dark matter candidate.

Figure 3 :
Figure 3: Results of the zero temperature numerical scans at one loop level on the v S − λ SH plane with fixed parameters (a) m a = 10 GeV, m 0 = 20 GeV, λ = 0.17 and (b) m a = 150 GeV, m 0 = 200 GeV, λ = 0.4.The benchmarks satisfying zero temperature boundary conditions are shown in red, benchmarks not satisfying BFB are shown in blue and benchmarks where the EW vacuum is not the global minimum are shown in gray.Analytical conditions C1 (loop-level, BFB) and C2 (tree-level, global minimum) are shown by the blue and gray lines.The gray points below the blue line, which should also be colored blue, in the right plot have a global minimum different from the EW minimum, despite satisfying C2, due to the significant one loop effects to the small quartic coupling.

Figure 4 :
Figure 4: Patterns of symmetry breaking at different temperatures on the a − h plane with s = 0. We use × to denote a minimum, ⊗ to denote a saddle point, • to denote a maximum, and × ( ) to denote a stationary point that may be minimum (maximum).(a)

Figure 5 :
Figure 5: Results of phase transitions on the v S − λ SH plane for points satisfying the zero temperature boundary conditions for the same two benchmark scenarios as in Fig. (3): (a) m a = 10 GeV, m 0 = 20 GeV, λ = 0.17 and (b) m a = 150 GeV, m 0 = 200 GeV, λ = 0.4.Points satisfying the nucleation and SFOEWPT conditions, S 3 /T 140 and h EW (T n )/T n > 1, are shown as red solid circles, and points satisfying the condition h EW (T c )/T c > 1 are shown as red open circles.The analytical bounds derived in Sec. 3, C1-C4 are shown by the lines.The shaded regions are excluded by the corresponding conditions.C3', the condition derived at T c , is shown for m a = 10 GeV only, since the derivation for T c only applies to the case with m a < m h /2.The perturbativity requirement λ S < 4π is shown for m a = 150 GeV, while the λ S for m a = 10 GeV is always small enough to satisfy this requirement.The lowest y-axis value for λ SH for these figures is set by Eq. (3.23).

Figure 6 :
Figure 6: Benchmark points producing the observed baryon asymmetry, and the correct dark matter relic density, projected in the dark gauge parameter space defined by g and M Z .Gray points are excluded by LEP constraints -the main collider constraint in this mass region.Red and cyan points survived collider searches, although the latter are excluded by direct dark matter detection searches.

Figure 7 :
Figure 7: Direct detection channels of the dark matter candidate χ.

Figure 9 :
Figure 9: Main production channel (left) and decay channel (right) for the light singlet scalar a.

- 1 Figure 10 :
Figure 10: (left) Lifetime of the real singlet a as a function of λg 2 m a /m χ (x-axis) vs m a /m Z (y-axis); (right) Constraints from singlet scalar searches at 13 TeV LHC with L = 139 fb −1 .Solid circles pass both a and s search constraints.While most of small open diamonds are excluded by prompt decays of s, the larger open diamonds are only excluded by searches for the long-lived particle a singlet.The grey dashed vertical line is the m h /2 threshold.The left and right panels use the same color coding for the lifetime of the singlet a.
TeV LHC with luminosity L = 139 fb −1 are shown in Fig. (10) (right) with open diamonds, while points with solid circles require future collider searches.The constraints on different parameter spaces are summarized as follows:

Figure 11 :
Figure 11: Gravitational wave signals of the SFOEWPT from benchmark points generating the observed baryon asymmetry, the dark matter relic abundance, and satisfying all phenomenological constraints with m a < m h /2 (dark green curves) and m a > m h /2 (light green curves).The power-law integrated sensitivity curves for LISA, MAGIS-100, MAGIS-Space, AEDGE, DECIGO, BBO, and ET are shown for comparison.

Table 1 :
13unds on (h → ss(aa) → l + l − + X) projected to13TeV LHC with luminosity L = 139 fb −1 that are satisfied by all the benchmark scenarios (via on-shell Higgs at different values of mass) with similar kinematic regions studied in 3c H T 2 λ H v 2 H , since the other two terms are either suppressed by the small quartic couplings involving the singlets, or the small change of the h vev.Furthermore, since T < v H , 3c H T 2 λ H v 2 H 1, Eq. (3.32) is obtained by keeping only the leading order terms in the Taylor expansion.Solving Eq. (3.31)=0, we get the finite temperature vev of a at the EW vacuum | (v H ,0,0,T =0) (h − v H ) m a n T p , and ∂ 4 i V is of order 10 −2 since only the terms involving the derivative with respect to a contribute.The terms ∂ 4 i V (∆x i ) 3 AaT 2 and λ SH c H T 2 2λ H m 2