Gravitational Waves From a Dark (Twin) Phase Transition

In this work, we show that a large class of models with a composite dark sector undergo a strong first order phase transition in the early universe, which could lead to a detectable gravitational wave signal. We summarise the basic conditions for a strong first order phase transition for SU(N) dark sectors with n_f flavours, calculate the gravitational wave spectrum and show that, depending on the dark confinement scale, it can be detected at eLISA or in pulsar timing array experiments. The gravitational wave signal provides a unique test of the gravitational interactions of a dark sector, and we discuss the complementarity with conventional searches for new dark sectors. The discussion includes Twin Higgs and SIMP models as well as symmetric and asymmetric composite dark matter scenarios.


Introduction
Violent phenomena in the early universe can lead to large anisotropic fluctuations in the energy momentum tensor, which can act as sources for gravitational waves (GW). Strong first order phase transitions (PT) are an example of such a phenomenon, and it is well known that they can produce GWs [1][2][3][4]. Once produced, GWs propagate through space almost undisturbed, and can therefore serve as a unique probe of phenomena in the early universe.
Phase transitions in particle physics are usually associated with symmetry breaking, i.e. with the transition of the universe from a symmetric phase to a phase of broken symmetry, as the temperature drops below a critical temperature T * . Within the standard model (SM) of particle physics, at least two phase transitions should take place, associated with the breaking of electroweak symmetry around T * ∼ 100 GeV and with the breaking of chiral symmetry at the time of the QCD phase transition, T * ∼ 0.1 GeV.
Today we know that both the QCD and the electroweak phase transition are not first order, but proceed through a smooth cross-over [5][6][7][8], and can therefore not produce a strong GW signal through the usual mechanism (see however [9]). This can however be changed in models beyond the SM. Extensions of the SM which lead to a strong first order electroweak phase transition are particularly attractive since they can provide one missing ingredient for generating the observed baryon asymmetry of the universe. GW signals from such models have for example been studied in [10][11][12][13][14][15][16][17][18][19][20][21][22][23][24]. It is more difficult to modify the QCD phase transition, although a large neutrino chemical potential could be sufficient to provide a strong first order PT [25]. The resulting signal was studied in [26]. The aim of this work is to point out that gravitational waves could also be produced by a strong PT in a dark or hidden sector. The particular scenario we have in mind is a dark sector with a new SU(N d ) gauge interaction which confines at some scale Λ d . Such models have recently received renewed interest either as models of dark matter [27][28][29][30][31][32][33][34][35][36][37][38][39][40][41] or as part of the low energy sector of so called Twin Higgs models [42][43][44][45][46][47][48]. Different from generic hidden sectors [49], these models provide a preferred mass range and some restrictions on the particle content, such that the frequency range of the potential GW signal can be predicted.
Given that the SM QCD transition is not first order, we will review the known results on the order of the PT in strongly coupled gauge theories in the next section, followed by a discussion of models that fall into this category. In Sec. 3 we calculate the GW spectra that can be produced in these models, and compare them to the sensitivity of current and planned GW detection experiments in Sec. 4. We discuss the complementarity of GW experiments with other searches for dark sectors in Sec. 5, before presenting our conclusions.

Models with First Order Phase Transition
Near the QCD confinement scale Λ QCD , the dynamics of QCD is governed by three flavours, two of which are almost massless, while the strange quark mass is of order Λ QCD . Lattice studies [5,6,50] have shown that for these values of the quark masses, the QCD PT is a weak cross-over.
However this is not a generic result for QCD and similar theories, but more a consequence of the precise values of m u ≈ m d and m s in the SM. The QCD phase diagram for arbitrary m u,d and m s can be summarised in the so called Columbia plot, which is reproduced in Fig. 1, based on [51]. The pure Yang-Mills limit m u,d , m s → ∞ is known to have a strong first order PT [52] from the restoration of a global Z 3 center symmetry at low temperatures. The opposite m u,d , m s → 0 limit, i.e. the theory with three exactly massless quarks, also features a strong first order transition, related to the breakdown of the SU (3) × SU (3) chiral symmetry [53].
Here we are interested in SU (N d ) theories with N d ≥ 3 and n f massless 1 Dirac fermions in the fundamental representation, and with a confinement scale Λ d . To guarantee the existence of a confining phase we further impose n f < 4N d , in order to stay outside of the conformal window.
For pure Yang-Mills theories, n f = 0, the confinement phase transition is related to the restoration of a global Z N d ⊂ SU (N d ) center symmetry, which is broken in the high temperature phase. Lattice simulations have shown that this PT is strong first order for N d ≥ 3 [54]. For the N d = 3 case there is also the classic argument of [52]: In the case of a second order PT, the critical behaviour of the theory should be described by a Z 3 symmetric effective theory at an infrared fixed point. Since no fixed points are known for theories with Z 3 symmetry, the PT must be first order.
The case of non-zero n f for N d = 3 and N d → ∞ is discussed in [53]. For n f = 1 there is no symmetry breaking, and therefore no phase transition. For n f ≥ 3 the PT is first order 2 for any value of N d ≥ 3. In the large N d limit, the PT is also first order for n f = 2. The only case that is currently not understood is n f = 2 for finite N d , where no analytic arguments can be applied and lattice results are difficult to obtain.
To summarise, SU (N d ≥ 3) theories with n f massless flavours have a first order PT if either n f = 0 or 3 ≤ n f < 4N d .
The above discussion hopefully convinced the reader that strong PTs are generic in extensions of the SM which feature new confining gauge symmetries. In the following we will discuss a few examples which are well motivated either from a dark matter perspective or by naturalness arguments. These physically motivated scenarios furthermore provide a preferred range for the confinement scale Λ d , which allows us to estimate the temperature of the PT, T * ∼ Λ d , and therefore make a prediction for the GW spectrum.
Composite Dark Matter I (CDM1) In this class of models, the dark matter candidate is the lightest baryonic bound state of a SU (N d ) dark sector with n f dark quarks which are neutral with respect to all SM interactions. This allows the DM to be light, since it will only communicate with the SM through heavy mediators. The most natural realisation of these models is in the context of asymmetric dark matter, where the DM number density is related to the baryon asymmetry (see e.g. [55][56][57] for reviews). The measured DM density then implies a mass range of O(5) GeV for the DM particle, and therefore motivates T * ≈ Λ d ∼ (1 − 10) GeV.
As a concrete benchmark we will consider the dark QCD model [32,40], which, as the name suggests, consists of a SU (N d = 3) theory with n f ≥ 3 dark quarks, such that a strong PT is guaranteed.
Composite Dark Matter II (CDM2) Similar to the previous case, the DM is a baryonic bound state of a new SU (N d ) theory. The important difference now is that the dark quarks carry electroweak quantum numbers. Models of these type were considered for example in [21,22], where mass spectra and form factors are calculated from first principles on the lattice. Constraints on the invisible width of the Z-boson immediately suggest a confinement scale at or above the 100 GeV scale, while unitarity arguments on DM freeze-out place a limit of Λ d O(100) TeV. For fermionic DM (i.e. for odd N d ), direct detection furthermore imposes a bound of M DM ∼ Λ d 10 TeV [21]. 1 More precisely, n f is the number of fermions with m Λ d . 2 The behaviour can be understood by studying the effective theory near the confinement scale, which consists of n 2 f massless Goldstone bosons from the breaking of the SU (n f ) × SU (n f ) × U (1)A chiral symmetry. For N d → ∞, the anomaly is negligible, and the PT is first order due to the absence of infrared fixed points in the effective theory, provided that n f ≥ √ 3. For finite N d and n f = 3 the anomaly contribution to the effective lagrangian (i.e. the η mass) is cubic in the Goldstone fields, and therefore alone sufficient to render the PT first order. For n f > 3 the anomaly does not affect the PT behaviour, so the large N d result applies and the PT is first order. What remains unresolved is the n f = 2 case for finite N d . For more details the reader is referred to [53].
For scalar DM (even N d ) lower DM masses are allowed [22]. The benchmark model considered in [22] has N d = 4 with n f = 4 flavours and should therefore feature a strong first order PT. 3 Twin Higgs (TH) These models attempt to solve the hierarchy problem without introducing coloured partners for the top quark [42,43], and have recently received renewed attention [45][46][47]. This is achieved by adding a twin sector to the SM with an approximate Z 2 parity symmetry, with the minimal requirement that there should be fermionic partners for the top and bottom quarks which are charged under a new SU (N d = 3) interaction. The approximate Z 2 symmetry constrains Λ d to lie in the (1 − 10) GeV range. Furthermore the twin top and bottom partner masses are necessarily much larger than Λ d , such that the theory at the PT scale is SU (3) with n f = 0, which again predicts a strong first order PT.

Strongly Interacting Massive Particle (SIMP)
The key idea behind SIMP DM is that the relic density is determined through freeze-out of 3 → 2 instead of 2 → 2 annihilations [58]. The correct relic abundance is then obtained for DM masses in the 100 MeV range. One simple realisation of this mechanism is a SU (N d ) dark sector with dark pions as DM candidates. The 3 → 2 processes can then be induced by the Wess-Zumino-Witten term [59,60]. The condition for the existence of the WZW term is n f ≥ 3, which again allows for a strong first order PT, this time in the 100 MeV − GeV range.
All these models should give rise to a GW signal from the strong PT, and the resulting spectra and detection possibilities will be discussed in the following sections. Before moving on, we should mention a few other scenarios which also contain new strongly coupled sectors. The first group consists of Hidden Valley [49] and vector-like confinement models [61], which propose the existence of new confining sectors, communicating with the SM either through heavy mediators or directly via SM gauge interactions. If the particle content is such that the models fall into one of the strong PT regions, then also these models will give rise to a GW signal. Furthermore there are so called mirror-world models [62], where the dark sector consists of an exact copy of the SM. Here one does not expect a GW signal for two reasons. First, if the quark content is an exact copy of QCD, there will not be a first order PT. If for some reason the quark masses are modified, such that a strong PT happens, cosmological constraints on the mirror photon require the mirror sector to have a lower temperature. Since the energy density is proportional to the fourth power of the temperature, and the intensity of a GW signal is proportional to the energy density, this will strongly suppress any GW signal from a mirror world.

Gravitational Wave Spectra
Gravitational waves produced at a time t * (or equivalently at a temperature T * ) will propagate undisturbed in the expanding universe, therefore their frequency f and their fraction of the critical energy density Ω GW will decrease as a −1 and a −4 , respectively, where a(t) is the scale factor [10,63]. Denoting by a * and a 0 the scale factors at time of production and today, entropy conservation (sa 3 = const) implies Here T 0 = 2.725 K = 2.348 × 10 −13 GeV is the temperature of the CMB, g * ,s (g 0,s = 3.91) is the effective number of relativistic degrees of freedom contributing to the total entropy at the time of production (today), and the entropy density at temperature T is given by and g i counts the internal degrees of freedom of the i-th particle. It follows that the frequency today can be expressed as where we have used the Hubble rate at time of production, H * = 4π 3 g *

45
T 2 * M Pl , and assumed that all species are in thermal equilibrium at T = T * , i.e. g * = g * ,s . For the fraction of energy density in gravitational waves today we similarly obtain where we used that ρ crit /ρ * crit = H 2 0 /H 2 * and H 0 = 2.13 × h × 10 −42 GeV. It is noteworthy that the intensity of the GW signal is independent of the production temperature T * (except for the implicit dependence of g * on T * ).
The most sensitive frequency regions of pulsar timing arrays and satellite based experiments are in the nano-Hz and milli-Hz range, respectively. To get an idea about the detectability of GWs from a strong dark PT we will therefore need to understand the spectrum of the produced GWs. For this, we will closely follow the discussion of [26].
Gravitational Waves are sourced by tensor fluctuations of the energy momentum tensor of the primordial plasma. During first order phase transitions both bubble collisions [64,65] and magnetohydrodynamical (MHD) turbulence [66][67][68][69] provide sources of gravitational waves. As functions of the conformal wave number k = 2πaf , the GW spectra produced by either source can be approximated by [26] dΩ Eqn. (5) is based on [70,71] while Eqn. (6) is a fit to the numerical results of [72]. Here H * is the conformal Hubble parameter H = Ha at T = T * , and Ω r0 is the radiation energy density today. The quantities that determine the GW spectrum are the bubble nucleation rate β (the duration of the PT is β −1 ), the bubble velocity v and the relative energy density in the source, Ω S * = ρ S * /ρ * ,crit = Ω * GW . The temperature of the PT enters through the dependence of H * on T * . The duration of the PT is usually taken as (1 − 100)% of a Hubble time, and therefore β = (1 − 100)H [2]. To understand the relation with the physical frequency, remember that the conformal frequency is related to the conformal wave number via af = k/(2π). Furthermore using H = Ha we see that f * /H * = af * /H * = (k/H * )/(2π), which together with Eqn. (3) allows us to translate the GW spectra into physical frequencies.
In a given theory, the dynamics of the phase transition, and therefore the parameters β, v and Ω S * , are in principle calculable. For the strongly coupled models considered here they are however not known, and can only be estimated using lattice simulations. We will therefore take β, v and Ω S * as additional input parameters, with values motivated by results of analyses in weakly coupled models. Following [26], we will use Ω S * = 0.1 and β = 10H * , but v = 1.0 instead of 0.7. We are now in a position to study the location of the peaks of the GW signals from bubble collisions and MHD turbulence. The bubble collision signal is triangular shaped with a maximum at k/β = 4 √ 3 ≈ 1.3, while the MHD turbulence peaks at somewhat larger wave numbers k/β ≈ π 2 /v. To obtain physical frequencies, we use Eqn.
In Fig. 2 we show the location of the frequency peaks as function of the PT temperature T * and β. As expected from Eqn. (3), the peak frequencies increase linearly with the transition temperature T * and with β/H * .
The source term Ω S * can be different for bubble collision and turbulence. Here we will assume that equal amounts of energy act as source for Ω Recent simulations of first order PTs suggest that sound waves generated by the expansion of bubbles could be the dominant source of GWs from these transitions [73][74][75]. Sound waves continue propagating through the early universe after the PT is finished, and decay on a timescale H * . Compared to the above discussed spectra, they will therefore not be suppressed as much by the velocity of the transition β, and the signal could be increased by a factor (β/H * ) compared to the bubble collision signal, but with a spectrum decaying as k −3 . This could potentially boost the signal, in particular for cases where the PT is fast, i.e. β/H 1.

Detectability
In the previous section, we have seen that the peak frequencies of GW signals from GeV-TeV scale PTs are of order (10 −6 − 10 −3 ) Hz. Furthermore it is important to note that a broad spectral region around the peak is populated by GWs, from (10 −10 − 1) Hz.
GWs with frequencies down to 10 −5 Hz can be probed by satellite based experiments like eLISA [76], however the sensitivity quickly degrades below 10 −3 Hz. On the other end of the spectrum, pulsar timing arrays (PTA) can probe frequencies in the (10 −9 − 10 −7 ) Hz range. In Fig. 3 we overlay the expected GW signal for different model parameters with the expected sensitivities of current and planned GW detection experiments (based on [77]).
Clearly the most promising signals are those from models with a PT temperature in the TeV range, where the peak region of the GW spectrum falls right into the most sensitive frequency range of satellite experiments. Here the signal should even be detectable for choices of the parameters that are less optimistic than those used for Fig. 3. Models of the CDM2 type naturally fall into this region, but also the CDM1 models can be viable with a confinement scale in that region.
The fact that TeV dark sectors predict an observable GW signal is not surprising, since here we are in the energy range of the electroweak PT or beyond, and the observability in particular of TeV scale strong PTs has been noted before [79,80]. The novelty here is that the dynamics leading to this strong PT does not have to be connected to the electroweak sector of the SM, and is therefore not in tension with the non-observation of new physics at the LHC.
For models with T * ∼ (1 − 10) GeV the situation is a bit more difficult, since the signal peak ends up in a frequency region where neither PTA nor (e)LISA are sensitive. Looking at the T * = 3 GeV curve in Fig. 3 more closely, we see that in the best case scenario, for β = H, both PTAs and LISA would be able to detect parts of the GW spectrum. For larger β the signal quickly drops out of the PTA sensitivity region, however LISA remains sensitive. This is due to the increase of the observed frequency with (β/H), which partially compensates the overall (H/β) 2 drop of the signal in the LISA sensitivity region. Therefore there is a chance to detect a GW signal from the CDM1 and TH models, even if the predicted frequency range is not optimal.
Recently the EPTA experiment has reported the first limit on a stochastic GW background [78], which probes values of Ω GW h 2 of order 10 −8 in the (10 −9 − 10 −7 ) Hz range. Their limit can not directly be displayed in Fig. 3 since it depends on the assumed spectrum, but it is nevertheless interesting to note that GWs from a 100 MeV scale dark sector could already be detectable.
Another way of interpreting these results is that GW searches have the potential to observe signals of completely unknown physics in the early universe, while the absence of a signal in a certain frequency range would place (very weak) bounds on the dynamics of new physics at the corresponding energy scale.

Complementarity
Signals of the new physics models introduced in Sec. 2 are being searched for at collider experiments and in direct and indirect dark matter searches. The possibility to observe such a signal in those experiments always relies on sufficiently strong non-gravitational interactions of the dark sector with the SM. Instead the GW signal is unique in the sense that it probes the gravitational effects of a model at very early times, and, while the detectability of a signal eventually still depends on model parameters, those are relatively independent of whether the model is detectable at colliders or in dark matter searches.
In CDM1 type models, or more general, hidden valley models with a new strong interaction, the dark sector is neutral under all SM interactions, and communicates only through a heavy mediator. If that mediator is in the TeV range, both direct detection experiments and the LHC could discover these models, but for larger masses it becomes increasingly difficult. From the DM perspective on the other hand, mediator masses of 10s of TeV are acceptable, such that part of the parameter space will remain unexplored in the near future. The GW signal instead is completely independent of the mediator mass, such that it could be detected even if the model was not discovered before, and thanks to the straightforward connection of mass scales and GW frequencies, could even motivate future experimental efforts.
The situation is similar for TH models, which mainly communicate with the SM through the Higgs portal. Here the parameter that controls detectability is not a mass scale but the smallness of the mixing angle of the Higgs with its twin partner. First hints for such a model could come from deviations in Higgs couplings and from exotic Higgs decays, but even then it would be difficult to uncover the whole structure. A GW signal in the right energy range could provide much information about the dynamics of these models near their confinement scale.
Finally the CDM2 type models have many detectable features at hadron colliders. However the upper limit on their mass scale is of order 100 TeV, which comes from the requirement that dark matter is not overproduced. Such high scales are not in reach of the current or next generation of collider experiments, but seem very accessible by GW searches, since the signal moves into the most sensitive frequency region of LISA and other satellite based experiments.
Overall we see that GW experiments provide a unique window to explore the dynamics of these models in the early universe, even if they are not discovered in the near future. Of course it would be even more interesting if such a dark sector was discovered at the LHC. The GW signal would then provide an important and unique probe of the physics of those models in the early universe.
Moving beyond the concrete models discussed above, perturbative unitarity constrains the mass of thermal DM to be below 110 TeV [81,82]. For composite non-perturbative DM this limit does not apply directly 4 , instead a lower bound on the radius of the extended object can be obtained, R (100 TeV) −1 . It is reasonable to expect the radius R to be of order of the inverse mass, which again implies an upper bound on the DM mass of order 100 TeV. GW signals could therefore be a unique probe of the thermal DM paradigm.

Conclusions
Models beyond the standard model with a confining dark sector can lead to unexpected phenomenological signatures. Here we have explored the possibility to detect gravitational waves due to a first order phase transition at the confinement scale Λ d . The main messages from this paper are as follows: • Different from QCD, dark sectors with QCD-like interactions can undergo strong first order phase transitions, with only mild constraints on the particle content of the theories.
• Several classes of new physics models that are currently being explored fulfil the criteria for first order PTs. The physics problems these models are trying to address, either dark matter or naturalness, constrain the confinement scales and therefore the temperature range of the phase transition.
• The GW signals originating from these dark phase transitions are in the detectable frequency range of future GW experiments, either at (E)LISA and BBO for high scale models, or in PTA experiments for the lower end of the spectrum.
Depending on other aspects of the model, GW signals will either provide complementary information about the models in question, or might even be the the best option to find evidence for these models of new physics. Different from the electroweak PT, a null result at LHC will not strongly disfavour a strong PT in a dark sector, although of course a confirmation of their existence would be more exciting.
A shortcoming of the present work is a lack of precise quantitative predictions. The bubble velocity v as well as the time scale of the phase transition β −1 and the energy fraction Ω S * are currently unknown, and are set to optimistic (but not unrealistic) values. Two approaches seem possible to improve upon this situation: On one side, lattice simulations could be used to measure quantities like the latent heat and the surface tension, which are related to the above parameters and can be used to obtain a more quantitative prediction for the GW spectra. Alternatively, one could attempt to construct a holographic dual for some of these theories, and analyse the PT in that setup.