Electroweak phase transition and gravitational waves in a two-component dark matter model

We investigate an extension of the Standard Model (SM) with two candidates for dark matter (DM). One of them is a real scalar field and the other is an Abelian gauge field. Except for these two, there is another beyond SM field which has unit charge under a dark $ U_{D}(1) $ gauge symmetry. The model is classically scale invariant and the electroweak symmetry breaks because of the loop effects. Although SM is extended with a new dark symmetry and three fields, because of scale invariance, the parameter space is strictly restricted compared to other two-component DM models. We study both DM phenomenology and electroweak phase transition and show that there are some points in the parameter space of the model consistent with DM relic density and direct detection constraints, while at the same time can lead to first order electroweak phase transition. The gravitational waves produced during the phase transition could be probed by future space-based interferometers such as Laser Interferometer Space Antenna (LISA) and Big Bang Observer (BBO).


Introduction
The standard model (SM) is the most successful model in describing the observed phenomena, however, there are still intriguing questions, such as the nature and the origin of dark matter (DM), waiting to be answered. There is a wide variety of astrophysical and cosmological observations as well as theoretical arguments that led the scientific community to adopt DM as an essential part of the standard cosmological model (for a fascinating review on history of DM see [1]).
The evidence of DM was strong enough that many strategies have been pursued to reveal its particle nature. For decades, the leading theory of particle DM was a single-component thermal relic with weak size couplings and mass, known as a weakly interacting massive particle (WIMP). These one-component scenarios are increasingly constrained by experimental measurements. Therefore, workers in the field compelled to examine more complex models of dark sector including the multicomponent ones where the total relic abundance of DM is due to the existence of multiple DM species . After all, it should not be surprising if the dark sector has multiple species like SM itself.
On the other hand, SM has a crossover rather than a true phase transition [46], while some extensions of the SM, e.g., with DM candidates  , lead to first order phase transitions with gravitational wave (GW) signals. In the case of first order phase transition, just below the critical temperature, the Universe goes from a metastable quasiequilibrium state into a stable equilibrium state, through a process of bubble nucleation, growth, and merger which generates GWs [73][74][75][76][77]. GW signatures are therefore a new window towards new physics, complementary to that provided by the Large Hadron Collider. Another motivation for studying electroweak phase transitions is the requirements to explain the matter-antimatter asymmetry in the universe [78], which one of them, departure from thermal equilibrium, is inevitable in a first-order phase transition.
The first direct detection of GWs was performed in 2015 [79]. The signal came from the strongest astrophysical sources of GWs, i.e., compact binary systems, in frequency from 35 to 250 Hz. Here, we are interested in production of GWs by first-order phase transitions (for a recent review see [80]). For GWs sourced by cosmological phase transitions, the relevant mission is Laser Interferometer Space Antenna (LISA) [81] which is a space-based interferometric gravitational wave detector working with three satellites orbiting the Earth. LISA is most sensitive at frequencies in the range 10 −3 − 10 −2 Hz and its planned launch year is 2034 with a mission life-time of 4 years. The Big Bang Observer (BBO) [82] is a proposed follow-up experiment consisting of four LISA-like detectors.
In this paper, a model with three motives will be presented. We come up with a beyond SM model to provide a solution for DM problem, hierarchy problem, and vacuum instability. As a solution for hierarchy problem [83], the model has classical conformal symmetry. On the other hand, SM suffers from vacuum instability and beyond SM models with bosonic degrees of freedom can probably solve this issue. As mentioned, due to the strict constraints of direct detection on one-component DM models, two component DM models are more appropriate. Therefore, our two-component DM model consists (bosonic) scalar and vector DM. We assume that the dark sector interacts with the SM particles only through the Higgs portal. The dark sector consists of three new fields, a real scalar, a complex singlet and a real vector field and it introduces a dark U D (1) gauge symmetry. Our model is a two-component DM model which both DM particles are bosons: one DM particle is spin zero and the other is spin one. As a potential solution to the hierarchy problem, we constrain our model to be a classically scale-invariant extension of the SM. Within this framework, all the particle masses are generated dynamically, by means of the Coleman-Weinberg mecha.nism [84]. After constructing the model, we study DM phenomenology including relic density and DM-Nucleon cross section. DM relic density is reported by Planck Collaboration [85], and DM-Nucleon cross section is constrained by direct detection experiments such as the LUX [86], PandaX-II [87] and XENON1T [88]. These experiments are gradually approaching the so-called neutrino floor which is the ultimate sensitivity of future direct detection experiments [89]. We also concentrate on investigating the possibility of achieving a strongly first-order electroweak phase transition within the parameter space of the model. To study electroweak phase transition, we present the complete expression of the finite-temperature 1-loop effective potential, including the contributions of the resummed thermal bosonic daisy diagrams, and show that the finite-temperature corrections induce a first-order electroweak phase transition. We identify regions of parameter space of the model which is consistent with DM relic density and direct detection constraints, while simultaneously realizing a strongly first-order electroweak phase transition. The GW signal from the phase transition is sufficiently strong to be detectable by LISA and BBO.
The rest of the paper is organized as follows. In section 2, we introduce our model. Section 3 is dedicated to DM phenomenological constraints. In section 4, we study the electroweak phase transition and GWs spectrum. Our result is given in section 5 where we simultaneously consider DM phenomenolgy and GW spectrum using two benchmark points (BMs) as shown in table 1. Finally, we conclude in section 6.

The model
Our model consists three beyond SM fields which all of them are bosons, namely, a complex singlet φ, a real singlet S, and a real Abelian vector field V µ . In our setup, the complex scalar field φ has unit charge under a dark U D (1) gauge symmetry with the vector field V µ . All of these fields are singlet under SM gauge group. However, the dark sector is invariant under the charge conjugation of U D (1) and parity of S: (2.1) Due these two symmetries, the model can have two component DM. In the dark sector the discrete symmetry V µ → −V µ forbids the kinetic mixing between the SM U Y (1) gauge boson B µ and the vector field V µ which makes V µ stable and a DM candidate. The other DM candidate is due to S → −S symmetry which makes the real singlet field stable.
The Lagrangian is given by and L SM is the SM Lagrangian without the Higgs potential term. We constrain V (H, φ, S) by: • gauge symmetry, • Z 2 symmetry, • scale invariance, and • renormalizablity.
Regarding these constraints V (H, φ, S) reads In this model, dark sector interacts with SM via Higgs portal. Because of Z 2 symmetry, the real singlet field, S, does not get vacuum expectation value (VEV), but electroweak symmetry, as well as Abelian dark symmetry spontaneously break after H and φ develop VEVs. In the unitary gauge, H † = (0, h1 √ 2 ) and φ = h2 √ 2 , the tree level potential is given by The local minimum of the tree level potential defines the VEVs of the fields which we take as ( h 1 , h 2 , S ) = (ν 1 , ν 2 , 0), leading to the following condition The last relation defines the flat direction in field space where the tree level potential is minimum, along this direction V tree (ν 1 , ν 2 , 0) = 0. Now we substitute h 1 → ν 1 + h 1 and h 2 → ν 2 + h 2 which mixes h 1 and h 2 . The mass eigenstates, h and ϕ, can be obtained from the following rotation where tan α = ν 1 /ν 2 ( h = 0 and ϕ = ν = ν 2 1 + ν 2 2 ). We identify h, which is perpendicular to the flat direction, as the SM-like Higgs observed at the LHC with M h = 125 GeV [91,92]. On the other hand, we know from SM that ν 1 = 246 GeV. The field ϕ is along the flat direction, thus its tree level mass is zero. However, 1-loop correction leads to a specific value along flat direction as the minimum of the potential which gives the following mass to ϕ (see subsection 4.1): where M S,V,W,Z,t being the masses for scalar DM, vector DM, W and Z gauge bosons, and top quark, respectively, after symmetry breaking. We can substitute the parameters of the Lagrangian using According to these relations there are five free parameters which we choose them as M S , M V , g, λ φS , and λ S . On the other hand, λ S is irrelevant in DM phenomenology and phase transition studied in this paper, therefore we left with only four parameters.

Dark matter phenomenology
In the following we will focus on the DM relic density constraint reported by Planck collaboration [85] and the available data on direct DM detection.

Relic density
It is considered in general that in WIMP scenarios, the DM relic density is inversely proportional to the thermally averaged DM annihilation cross section into SM particles. In the case of two-component DM, the situation is more interesting since there are additional important processes such as conversion of one DM component into another which complicates the analysis. We use the public numerical code micrOMEGAs [93] for solving the two Boltzmann equations governing the cosmological evolution of our DM candidates, the scalar DM (S) and the Vector DM (V ). The Boltzmann equations are determined by three types of processes: • DM Annihilation: pair annihilation of both DM components into SM particles, • DM Conversion: which converts one DM component into another, and • DM Semi-(co)annihilation: which could change the abundances of both DM components The relevant diagrams for the annihilation processes of the two DM components are presented in figure 1. There are no (subleading) semi-annihilation and co-annihilation processes in our scenario. As we see in figure 1, Feynman diagrams are the same for both DM components. In these figure, SM and SM stand for massive SM particles and anti-particles, respectively. Besides DM annihilation into SM particles, the two DM candidate can also annihilate into each other (DM conversion: S S ←→ V V ) which are shown in figure 2. The coupled Boltzmann equations for scalar S and vector V DM are given by: where j runs over SM massive particles and h, ϕ. By changing the variable, x = M/T and Y = n/s, where T is the photon temperature and s is the entropy density, one can rewrite the Boltzmann equations in terms of Y = n/s: where g 1/2 * is the degrees of freedom parameter and M pl is the Planck mass. The last terms in summations are new terms in Boltzmann equations which describe the conversion of DM particles into each other. Since these two cross sections are described by the same matrix element, we expect σ V V →SS v and σ SS→V V v are not independent and their relation is: The interactions between the two DM candidates take place by exchanging two scalar mass eigenstates h and ϕ where the coupling of V to h is suppressed by sin α. Therefore, it is usually the ϕ-mediated diagram that gives the dominant contribution. Notice that the conversion of the heavier particle into the lighter one is relevant. The relic density for each DM candidate is related to Y at the present temperature through Ω S, where h is the Hubble expansion rate at present times in units of 100 (km/s)/M pc, and the total relic density of DM according to the data by the Planck collaboration should be [85], Finally, we define the fraction of the DM density of each component by, Figure 3 depicts the relic density of both DM components as a function of free parameters of the model. According to these figure, λ φS , and M S are not relevant in Ω V (Note that, in general, Ω V depends on M S , because M ϕ is a function of M S ). However, Ω S depends on all four parameters. In this case, besides λ φS and M S which are relevant, λ HS can also change Ω S . On the other hand, considering eq. (2.8), λ HS is dependent in all free parameters: So it is not a surprise that all free parameters can effect the relic density of scalar DM. The minimum in figure 3 (c) at M S M h 2 is due to resonance, and maximum in all four diagrams is due to λ HS 0.

Direct detection
The direct detection experiments aim to study DM-Nucleon interactions. These events induced by DM particles from the Milky Way's halo. The Standard Halo Model assumes that the DM particles are distributed in an isotropic isothermal sphere with a Maxwellian velocity distribution. The local DM density ρ 0 adopted for the interpretation of direct detection experiments is ρ 0 = 0.3 GeV /c 2 /cm 3 . The possibility of DM direct detection Planck Planck Figure 3: Variation of DM Relic density respect to parameter space in the form of WIMPs was first discussed in [94]. The idea is simple: since in most scenarios the WIMP carries no electric charge, therefore it will not interact with the atomic electrons, however, DM particles can elastically scatter off the atomic nucleus and the momentum transfer gives rise to a nuclear recoil which might be detectable. In our model both DM candidates interact with quarks via Higgs portal, see figure 4, which results in a spin independent DM-Nucleon cross section.  The relevant DM-quark interaction terms in Lagrangian are: where q stands for quarks and DM scattering on nuclei have a characteristic energy scale of the order of 1 GeV and very low momentum exchange between the DM and the nucleon. On the other hand, the low DM velocity allows to consider this process in the nonrelativistic limit. Therefore, the scattering of DM with nucleons can be described by effective four field interactions between the DM and the SM quarks. After integrating the scalar mediators out, the low-energy 5-dimensional effective interaction of the DM with quarks will be From this, it is possible to obtain effective interactions between the DM particle and a nucleon which gives DM-Nucleon cross section of scalar DM and vector DM [95] where M N is the nucleon mass and f N = 0.3 parametrizes the Higgs-Nucleon coupling.
The Higgs portal scattering that we discussed here led to spin independent interactions of the DM with nuclei. Some of the present constraints on DM-Nucleon spin independent interactions come from the world leader experiments, such as LUX [86], PandaX-II [87] and XENON1T [88]. Finally, the DARWIN experiment [96], with sensitivity close to the irreducible background coming from scattering of SM neutrinos on nucleons (the so-called neutrino floor [89]), would be the ultimate DM detector. In Section 5, we constrain the model with the results of the PandaX-II experiment [87] which set an upper limit on the spin-independent WIMP-Nucleon cross section with the lowest exclusion at M DM = 40 GeV: PandaX-II : σ DM-N 8.6 × 10 −47 cm 3 . (3.14) Note that in above constraint, it is assumed that the local DM density is only provided by one DM specie. However, in our scenario both scalar and vector DM contribute to the local DM density. Assuming that the contribution of each DM candidates to the local DM density is the same as their contribution to the relic density, many authors constrain the rescaled DM-Nucleon cross sections, i.e., ξ S σ s and ξ V σ v , with experimental results. However, both DM candidates contribute the DM signatures and one must combine both signatures. For the large DM mass case, where DM energy is much larger than detector threshold energy, the statistical combination is easy and the direct detection constraint reads [90] ξ for M 40 GeV.
We have depicted DM-Nucleon cross section as a function of free parameters of the model in figure 5. Again, λ φS is irrelevant in vector DM phenomenology, e.g., σ v . In all diagrams there are some dips where α s or α v vanishes. According to eq. (3.12), If λ HS = 0, then we expect that dips occur at the same place where M h M ϕ . But, in general , the dips occur at different places as we see in figure 5. Another interesting feature is the double dips of σ v in (d) diagram of figure 5. If we solve α v = 0 (⇒ M h = M ϕ ) to find M V , we will find two solutions corresponding these double dips.
In our model, indirect detection limits are not competitive with the ones from direct detection and we will not explicitly discuss them here.

Electroweak phase transition
In order to study the electroweak phase transition in our model, we need to construct the effective potential. In our scenario, the effective potential is a function of the scalar field ϕ and temperature T . As the Universe cools down, the VEV varies from ϕ = 0 to ϕ = ν = 0. In the following subsections we derive 1-loop effective potential and study GWs during the cosmological phase transition.

One-loop effective potential
The effective potential were initially studied at 1-loop level by Coleman and Weinberg [84]. A few years later, Gildener and Weinberg, presented their formulation for a scale invariant theory with many scalar fields [97]. As we discussed in Sec. 2, along the flat direction, the tree-level potential is zero, therefore, 1-loop corrections are dominated. The 1-loop effective potential at zero temperature is given by and Λ is the renormalization group (RG) scale. The other parameters are: C k = 3/2 (5/6) for scalars/spinors (vectors), M k for the measured mass of particles, and g k for the number of degrees of freedom of the particle k (it is positive for bosons and negative for fermions). In order to have a non-zero VEV, the potential (4.1) should have a minimum at ϕ = 0: Considering both eq. (4.1) and eq. (4.4) one can substitute RG scale Λ and find a final expression for the 1-loop potential in terms of b coefficient and the true vacuum expectation value ν: and for the mass of ϕ we have M 2 ϕ = 8bν 2 which gives eq. (2.7). Vacuum stability of a model depends on the behavior of the effective potential. If the vacuum of the effective potential is a global minimum, then the vacuum is absolutely stable. Vacuum stability up to Planck scale puts constraint on the parameters of models. After the discovery of 125 GeV Higgs boson at the LHC, we know that for SM the vacuum is not stable if no new physics is assumed (for a brief review of vacuum stability in SM see [98]). In conformal models tree-level potential is zero along flat direction and 1-loop contribution determines the behavior of the effective potential. According to conditions (4.3), the effective potential (4.5) has global minimum if b > 0. To fulfill this condition, considering eq. (4.2), we need new bosonic degrees of freedom with constrained masses. Albeit, for full treatment of vacuum stability, one should obtain one-loop β-functions and solve renormalization group equations (RGEs) in order to derive running coupling constants (see, e.g., [99,100] where classically scale-invariant non-Abelian extensions of the SM are constructed satisfying perturbativity and stability up to the Planck scale.) Apart from 1-loop zero-temperature potential (4.5) , the 1-loop corrections at finite temperature also contribute to the effective potential which is given by [101] where J B,F are thermal functions: We approximate thermal functions J B (x) and J F (x) in terms of modified Bessel functions of the second kind, K 2 (x), which is a good approximation both in high and low temperature regimes [61]. We also consider resummed daisy graphs contribution given by [102] V daisy (ϕ, T ) = n k=1 g k T 4 12π where the sum runs only over scalar bosons and longitudinal degrees of freedom of the gauge bosons 1 . The thermal masses in (4.9), Π k (T ), are given by Finally, our effective potential contains Gildener-Weinberg term (4.5) and finitetemperature contributions (4.6) and (4.9): (4.11) In our calculations, in order to get V ef f (0, T ) = 0 at all temperatures, we subtract a constant term from potential: Having potential (4.11), now we are ready to study phase transition.

First order phase transition and gravitational waves
The first order phase transitions in the early Universe leave imprints in GWs which could be detected in the future. Many beyond Standard Models predict a first-order phase transition at the electroweak scale. This transition also provides an explanation for the matter-antimatter asymmetry in our Universe. In the first order phase transition, just below the critical temperature, the Universe goes from a metastable false vacuum into a stable true vacuum, through a process of bubble nucleation, growth, and merger. Such a first-order phase transition may occur in the early Universe and naturally produces GWs [73][74][75][76][77]. In the following we will study the dynamics of first-order phase transition and search for the parameter points of our model that can cause such transitions.
The effective potential (4.11), at some critical temperature T C , have two degenerate minima separated by a high barrier: one in ϕ = 0 and the other in ϕ = ν C = 0: By solving these two equations, one can obtain ν C and T C . Although, all independent parameters of the model contribute in the effective potential, we find that daisy term is negligible compared to other terms, therefore λ φS and λ S are irrelevant and dynamic of the phase transition only depends on M S , M V , and g. In figure 6, we have depicted ν C and T C as a function of free parameters of our model. Above the critical temperature ϕ = 0 is the true vacuum and the symmetry is not broken. As the Universe cools down, the temperature drops below the critical one and we expect a phase transition from the false vacuum ϕ = 0 to the true vacuum ϕ = 0 via thermal tunneling at finite temperature. Once this transition has happened, bubbles of the broken phase form in the sea of the symmetric phase and spread throughout the universe converting the false vacuum into the true one.
The bubbles formation starts after the temperature drops below T C , however it goes sufficiently fast to fill the universe with bubbles of the new phase only at some lower temperature, the nucleation temperature T N , where the corresponding euclidean action is S E = S 3 (T N )/T N ∼ 140 2 . The function S 3 (T ) is the three-dimensional Euclidean action for a spherical symmetric bubble given by where ϕ satisfies the differential equation which minimizes S 3 : with the boundary conditions: dϕ dr r=0 = 0, and ϕ(r → ∞) = 0. (4.15) In order to solve eq. (4.14) and find the Euclidean action (4.13), we have used AnyBubble package [106]. In figure 6, we have also depicted ν N and T N as a function of g, M S , and M V . The stochastic GW background produced by strong first-order electroweak phase transitions comes from three contributions: • bubble walls collisions and shocks in the plasma, • sound waves to the stochastic background after bubble collisions but before expansion has dissipated the kinetic energy in the plasma, and • turbulence forming after bubble collisions.
These three processes may coexist, and each one contributes to the stochastic GW background: All of the above contributions are controlled by four thermal parameters: • the nucleation temperature, T N , • the strength parameter which is the ratio of the free energy density difference between the true and false vacuum and the total energy density, α, where ρ * is • the inverse time duration of the phase transition, β, • and the velocity of the bubble wall, v w , which is anticipated to be close to 1 for the strong transitions [107].
In the process of the GW production, first the bubbles of the stable phase collide and merge. This stage is subdominant compared to the subsequent stages of GW production, unless the bubbles grow as large as the Hubble length itself. The bubble collision contribution is given by [108] Ω coll (f )h 2 = 1.67 where S coll parametrises the spectral shape given by After the bubble collision, the shells of fluid kinetic energy continue to expand into the plasma as sound waves. These different waves overlap and source the dominant contribution to the GW signal which is given by 3 [110] Ω The spectral shape of S sw is where Hz. (4.25) Finally, the last stage is the turbulent phase which its contribution to the GW spectrum is given by [111] Ω where Hz. (4.28) In eq. (4.27), the parameter h * is the value of the inverse Hubble time at GW production, redshifted to today, h * = 1.65 × 10 −5 T N 100 In the formulas of GW spectrum we have used [112,113] κ = 1 1 + 0.715 α (0.715 α + 4 27 where the parameters κ, κ v , and κ turb denote the fraction of latent heat that is transformed into gradient energy of the Higgs-like field, bulk motion of the fluid, and MHD turbulence, respectively.

Results
After having studied DM phenomenology and electroweak phase transition of our twocomponent DM model in the previous sections, we will now concentrate on the question of whether it is feasible to correctly reproduce the known features of DM and first order phase transition at the same time. On the DM side, we have relic density constraint as well as the upper bound of DM-Nucleon cross section obtained in direct detection experiments. On the other hand, we are looking for strong first order phase transition which leads to GW production at the early Universe. Obviously, the above requirements will impose constraints on the parameter space of the model, which is the subject of the present section. In order to obtain the parameter space consistent with DM relic density (see eq. (3.6)) and direct detection constraint (see eq. (3.15)), we should scan over four independent parameters of the model, i.e., M S , M V , g, and λ φS . To do so, regarding the strong constraint on DM-Nucleon cross section, we first obtain g and λ using the following equations Now, looking for the correct value of DM relic density (3.6) and regarding perturbativity constraints (all couplings < 4π)), we scan over random values of M S and M V , while we choose 0.9 g < g < 1.1 g and 0.9 λ < λ φS < 1.1 λ. In this way, according to eq. (3.13), we restrict ourselves to the small values of σ s and σ v which can possibly evade the direct detection constraint (3.15). According to this strategy, we obtain DM relic density for 40 GeV M DM 2000 GeV.
In figure 8 (a) and (b), the parameter space consistent with DM relic density is obtained based on our strategy. In (c) and (d), for these parameters, rescaled DM-Nucleon cross sections, i.e., ξ S σ s and ξ V σ v , are also depicted. As we see, there are some points between the PandaX-II direct detection bound and the neutrino floor which can be probed in the future direct detection experiments. Although DM relic density and direct detection experiments restrict the model, there are some parts of the parameter space which is not excluded yet. In figure 8 (c) and (d), we have also depicted the fraction of vector DM, ξ V . As we see, for some scatter points, vector DM is dominated, while for the other points DM mostly consists of scalar DM. Note that, although our strategy strongly restricts the parameter space in order to get unconstrained DM-Nucleon cross section, still there are some points which violate direct detection constraint (3.15). This points are depicted with hollow circles and they should be excluded even when they are below the upper bound of Pandax-II. Therefore, considering ξ S σ s and ξ V σ v separately is not enough and in order to constrain two component DM models by direct detection experiments, one should combine both signatures. Now looking for the first order electroweak phase transition and GW, we scan over the parameter space once more. In accordance with the DM results, we choose the same range of the parameters as figure 8. The result is shown in figure 9. In study of phase transition and GW, the relevant parameters are M S , M V , and g. We found that for around 14 percent of the parameter space scanned here, first order phase transition can occur, generating a stochastic background of GWs. Most peaks of the GW spectrum are detectable by LISA and BBO detectors 4 .
So far, we have studied DM phenomenology and phase transition separately. Now we consider both aspects simultaneously which gives us figure 10. In this figure, we have used the scatter points already obtained in the study of DM phenomenology, and saw if they also produce first order phase transition and GWs. In general, the DM freeze-out temperature T F ∼ M DM /20 may be greater than nucleation temperature T N in some parameter points. In this case, as the phase transition is not completed, the freeze-out can be affected. For phenomenology of a late phase transition see [117]. However, for the parameter points in figure 10, we have compared freeze-out temperature with nucleation temperature and find out T F < T N where T F ∼ max(M S , M V )/20. Therefore, this issue does not affect our result and the DM properties would not be modified between T F and the present day, at least for the parameter space considered in figure 10 (where 0.31 T F /T N 0.47).
As figure 10 implies, for some points in the parameter space, the model is consistent with DM constraints, while at the same time generates first order electroweak phase transition and GWs detectable by LISA and BBO. To be more explicit, we have chosen two benchmark points given in table 1. In this table all relevant quantities, including independent parameters of the model, DM properties, and phase transition parameters, are given. #  For these benchmark points, the GW spectrum is depicted in figure 11. For both benchmark points, around the peak, the dominated contribution of GW signal is sound wave. For the first benchmark point, scalar DM is dominant, while for the second one, vector DM makes most of DM relic density. The peak of the GW spectrum for both benchmark points falls within the observational window of LISA.
Finally, we should mention that sound wave contribution dominates for all scatter points in figure 10 and it is almost indistinguishable from the sum of the three sources, at least around the GW peak. However, far from the peak, this is not necessarily the case. For example, after the fracture of the first benchmark curve around 0.1 HZ, the bubble collision contribution will be dominated.

Conclusion
In this work we studied a two-component DM model as an extension of the SM with classical scale symmetry. It realizes electroweak symmetry breaking through Gildener-Weinberg mechanism and gives a natural solution to the hierarchy problem. The model consists three new fields: a real scalar, a complex scalar, and a vector field which two of them, the real scalar field and vector field, can play the role of DM. Our two-component DM model is obtained by adding them to the SM via the Higgs-portal. The Boltzmann equations for both DM components were solved numerically and in order to determine a region of parameter space that is consistent with Planck and PandaX-II data a scan over four dimensional parameter space was performed. After introducing the model and investigating DM phenomenology, we focused on the phase transition dynamics. With the aim of exploring the nature and the strength of the electroweak phase transition, the full finite-temperature effective potential of the model at one loop level has been obtained. Despite the absence of a barrier in the zerotemperature potential, it was demonstrated that the finite-temperature effects induce a barrier between the symmetric and the broken phase vacua, and thereby give rise to a first-order electroweak phase transition which can generate GWs.
The spectrum of these GWs can be described in terms of only four properties: the bubble nucleation temperature T N , the strength parameter α, the transition rate parameter β, and the bubble wall speed v w . These are all computable from the underlying particle physics model and therefore are functions of the independent parameters of our model. The space missions such as LISA and BBO could detect these GWs if the phase transition took place at scale of electroweak symmetry-breaking. LISA and BBO are particle physics experiments, as well as astrophysical observatories. Albeit, much should be done to realize the goal of making LISA into a particle physics experiment to complement the Large Hadron Collider (LHC). On the other hand, although continuing efforts at the LHC will be able to examine some of the beyond Standard Models, there remain many that cannot be probed through collider experiments on a timescale as good as LISA and BBO, if at all.
Our results indicate the model can survive DM relic density and direct detection constraints, while at the same time produce GWs during the first order electroweak phase transition. A positive GW signal at LISA and BBO would most likely point toward new physics at the TeV scale such as the classically scale invariant potential studied here.