Gravitational Waves from a Pati-Salam Phase Transition

We analyse the gravitational wave and low energy signatures of a Pati-Salam phase transition. For a Pati-Salam scale of $M_{PS} \sim 10^5$ GeV, we find a stochastic power spectrum within reach of the next generation of ground-based interferometer experiments such as the Einstein Telescope, in parts of the parameter space. We study the lifetime of the proton in this model, as well as complementarity with low energy constraints including electroweak precision data, neutrino mass measurements, lepton flavour violation, and collider constraints.

Gauge coupling unification in non-supersymmetric theories can therefore occur through multiple steps. One of the most well known examples of an intermediate scale model is the Pati-Salam (PS) model [2,21]. As opposed to fully unified models such as SU (5) or SO(10), PS models can survive at relatively low energies, because they do not induce rapid proton decay [22]. Primordial monopoles created at such low scales have been shown to be inflated away in models with a low scale strong phase transition [23,24]. A low PS scale is phenomenologically attractive, as it implies that unification has low energy consequences which can be constrained by precision data and collider experiments. Furthermore, PS symmetry breaking can in principle lead to a first order phase transition [25], so long as the post-inflation reheating temperature is larger than the breaking scale, which is favoured by the parameters of the theory.
Inhomogeneous cosmic phase transitions are associated with stochastic gravitational wave (GW) spectra. Such gravitational radiation is an important mechanism through which energy is dissipated when bubbles of the new vaccuum collide. The associated gravitational power spectrum is therefore a function of the set of parameters which govern the thermal evolution of the phase transition: the latent heat normalized to the radiation density, α, the speed of the transition β/H, the transition temperature T n and the velocity of the bubble wall upon collision v w . In this work, we compute the thermal parameters of a PS phase transition. We motivate an effective model with three free parameters, which describes the broken direction in the scalar potential and its most important thermal contributions. We study this parameter space, and show that the PS transition may lead to a stochastic spectrum which is observable in the next generation of ground-based interferometer experiments, such as the Einstein Telescope [26][27][28][29][30], and the Cosmic Explorer [31].
Ground-based interferometer experiments are sensitive to GW spectra with relatively low PS transition scales M P S ∼ O(10 5 ) GeV. As such, there are several low-energy experimental directions which may probe the PS-GW parameter space. Firstly, collider searches for righthanded neutrinos and gauge bosons become significant for low values of the scale M P S . Moreover, further low-energy constraints may come from the neutrino sector, as the SM neutrino masses are determined by v L and v R and they have strong contributions to lepton flavour violating processes. Lastly, a GW result may inform future experimental efforts in determining the lifetime of the proton, which we have currently used to set a lower bound on M P S .
The structure of this paper is as follows. In Section 2 we describe the model, its field content, finite temperature potential and the conditions for gauge coupling unification. Section 3 contains the study of the phase transition and the spectrum of gravitational waves it produces, highlighting the optimal scenario for visibility of the GW spectrum in the nextgeneration of experiments. We then run this scenario through some low energy probes in Section 4, including neutrino masses, lepton flavour violation, collider searches and proton decay. We conclude in Section 5 with a summary of the findings and a brief discussion on selected topics for expansion.

A first order Pati-Salam phase transition
The Pati-Salam model [2] is a good candidate for a first order phase transition which peaks within the frequency windows of the next generation of ground-based interferometer experiments. It can admit a low energy symmetry breaking scale, M P S < 10 7 GeV, which implies that the stochastic spectrum peaks within the experimental reach, f peak 10 3 Hz [49]. The model also has a fairly large gauge coupling constant, g 4 0.8, which increases the strength of the phase transition. Moreover, the rank of the broken group, SU (4), is larger than e.g. the electroweak phase transition, such that more latent heat is released [35].
The symmetry group of the PS model is G P S = SU (4) c × SU (2) L × SU (2) R . This unifies the quarks and leptons of a given chirality for each generation into merely two representations of the group. The matter content in this model is therefore embedded in the representations This gauge group and matter content are manifestly left-right symmetric, invariant under exchanges of SU (2) L ↔ SU (2) R . Manifest left-right symmetry (also known as D-parity) forces g L = g R , which makes gauge coupling unification an impossible task. We will thus explicitly break this symmetry by adding purely right-handed fields at the PS scale that ensure g L = g R and gauge coupling unification can occur.
We therefore construct a model with a PS and a left-right (LR) symmetric group, [50][51][52], as intermediate scales from an unified UV model, e.g. SO (10), with the breaking chain, Upon the construction of this model we aim to achieve gauge coupling unification and at the same time keep the PS breaking scale low M P S < 10 7 GeV.

Scalar field content
The minimum set of scalar fields for a valid PS model needs to be sufficient to trigger spontaneous symmetry breaking (SSB) of every step in the breaking chain. This requires the following set of fields which trigger the symmetry breaking of G SM , G LR and G P S respectively. In addition we add a few more scalar fields. A left handed triplet field ∆ L = {10, 3, 1}, which gives masses to the neutrinos via type II seesaw [53,54]. A right-handed coloured triplet, Ω R = {15, 1, 3}, to explicitly break manifest LR symmetry. And two adjoint coloured fields, Ξ 2,3 = {15, 1, 1} to help with gauge coupling unification. At the scale at which the PS group is broken, when Ξ 1 acquires a vacuum expectation value (vev), Ξ 1 = v, the off-diagonal gauge bosons G, the scalar fields Ω R , Ξ 2,3 and coloured components of ∆ LR get integrated out, with the following masses 1 where ρ 1 is a portal coupling between Ω R and Ξ 1 . The component of Ξ 1 that acquires the vev (Ξ v 1 ) gets a mass M Ξ v 1 = √ 2λ 1 v ≡ M P S , with λ 1 its quartic coupling. After PS symmetry breaking, the remaining scalar fields decompose into representations of the LR group as The field δ R is now responsible for the breaking of the LR group into G SM , and the vev of φ triggers electroweak symmetry breaking. Additionally, the field δ L acquires a vev at the same time, which has consequences for neutrino masses, as we will see later . The vevs of these fields can be expressed as where the SM vev is the combination v 2 SM = v 2 u + v 2 d . Lastly, after LR symmetry breaking, the gauge bosons associated with SU (2) R and U (1) B−L , as well as the δ R and δ L fields get masses that look like with the approximation that M W R M W L and M Z R M Z L and their mixing is negligible [55].

Gauge coupling unification
Much of the motivation for Pati-Salam models comes from their ultraviolet completion in SO(10) or E 6 , where all fermions are unified into a single representation of the group [3][4][5]. Although we will not worry about the details of the UV completion beyond the GUT scale, we enforce the unification of the gauge couplings as it implies a relation between the different energy scales which is the source of complementarity between low energy and gravitational wave searches.
The one-loop gauge Renormalization Group Equations (RGEs) for the gauge couplings are µ dg a dµ = b a 16π 2 g 3 a . (2.7) for a each element in a direct product of Lie groups, and b a a parameter that controls the slope of the RGE flow for g a and depends on group properties and the field content as [56] b a = 2 3 where C 2 (G a ) is the Casimir of the group G a , R f and R s the representations of fermion and spinor fields, respectively, S(R i ) is the Dynkin index of the representation R i and d ⊥ (R i ) its dimension in the groups orthogonal to G a . Equation (2.7) can be solved analytically for each step of the breaking chain, and iterated from M GU T to M Z by using matching conditions at each scale. In fact for α a = g 2 a 4π and t = 1 2π log µ the solution becomes a linear system of equations of the form for m steps in the breaking chain, i = 1, 2, 3 labels the SM gauge couplings and ∆t j = t j −t j−1 . The mechanism to achieve gauge coupling unification described above relies on the assumption that after each symmetry breaking, there is an Effective Field Theory remaining where many of the fields of the full theory have been integrated out, and these fields have all masses equal to the symmetry breaking scale. This is not generally the case and thus the matching conditions at each energy scale depend on the masses of these fields through the threshold corrections [57,58] with α j and α i the couplings of the theory before and after SSB, respectively. The threshold corrections at each scale can be computed as [59] λ ij (µ) = 1 12π where g, f and s label the vector, fermion and scalar fields integrated out at µ, and M i are their masses. The masses of the fields in the intermediate scales, eqs.(2.4)-(2.6), are mostly fixed by the symmetry breaking conditions and the structure of the potential, and these contribute towards the threshold corrections whenever their masses stray from the respective scales at  which they are integrated out. The mass of the field Ω R , integrated out at the PS scale, causes considerable threshold corrections since its mass is dominated by its portal coupling to Ξ 1 , as it will be seen later when we discuss gravitational waves. Lastly most the masses of the fields at the GUT scale are unconstrained, since they depend on the field content and SSB mechanism at the GUT scale, which we do not consider here. Hence we take the liberty of setting these masses to values that assist in achieve gauge coupling unification within the desired ranges of relevant mass scales. After adding threshold corrections we find that in our scenario with gauge coupling unification can be achieved for any value of M P S ∈ (2.9 × 10 3 , 2.25 × 10 7 ) GeV, and the remaining scales and couplings can be obtained in terms of M P S . The ranges for other relevant quantities are Fig. 1 shows the one-loop RGE evolution of the gauge couplings in this model, for the choice of M P S = 10 5 GeV, which we will later motivate as the optimal choice for the detection of gravitational waves. As can be noticed in the figure, around M P S and M GU T there are changes in slopes and discontinuities in the matching on the gauge couplings. These are a consequence of the threshold corrections described above, where the strong effect of the SO(10) corrections can be readily spotted.

Thermal potential
As described in the previous subsections, the Pati-Salam symmetry is broken when Ξ 1 acquires a vacuum expectation value. In the absense of new fermions charged under SU (4) C , the resulting phase transition can be described by the Pati-Salam gauge coupling, the portal couplings between Ξ 1 and other scalars, and the scalar field potential in the Ξ 1 direction. The only portal coupling that can be large without undesirable low energy consequences is the mixed quartic between Ξ 1 and Ω R . Here we will consider this term to set the effective mass of the Ω R field. As described in the previous subsection, the SU (4) C gauge coupling constant is fixed as a function of the PS mass. Then, our parameter space is limited to the two parameters in the Ξ 1 potential, and a single portal coupling. 2 With these considerations in mind we can approximate the scalar potential as follows It is convenient to reparametrize the zero-temperature scalar potential in terms of the overall scale and vev, We will see that the strength of the phase transition can be effectively determined based on the zero-temperature ratio v/Λ and the size of the portal coupling ρ 1 [35].
At one loop, we consider Coleman-Weinberg contributions and thermal corrections to the potential (2.14), In the above equation, the sum is over all bosons (n i denotes the multiplicity factors), GB refers to gauge bosons, µ is the renormalization scale (in our analysis, we will assume µ ∼ T ) 3 , and the mass terms are field dependent and given by 4 for the physical field φ, Goldstone modes G, gauge bosons GB and the scalar Ω R respectively. We also include Debye masses, given by [60], to delay the breakdown of perturbation theory at high temperature. The Debye masses can be approximated in the high temperature limit as 5 Note that we have assumed the combination of the Debye mass and the mass parameter, Π Ω R − µ 2 , is negligible compared to the field dependent mass m Ω R . In the next sections, we use the full 1-loop thermal potential in the φ direction, to find the thermal parameters of the phase transition.
3 Gravitational Wave spectrum from a PS phase transition

Strength of the phase transition
At high temperature T v, the scalar potential (2.22) has a single minimum at Ξ = 0. As the sector cools, a second minimum develops with Ξ = 0. The minima are degenerate at the critical temperature, Some intuition for the strength of the gravitational wave signal can be developed from the ratio φ C /T C for different parameter choices. In Fig. 2 this ratio is shown for different values of the portal coupling ρ 1 the ratio of zero temperature variables v/Λ. Here we have fixed the gauge coupling according to the relation in the previous subsection (2.12) and set M P S = 10 5 GeV. It is seen that the ratio peaks at around φ C /T C = 5 for large v/Λ, and portal 4 Radiative corrections to the tree level masses due to CW contributions have no appreciable effect at the target accuracy. 5 Going beyond the high temperature limit for the Debye masses requires solving a self consistency condition, which is outlined in ref [61]. coupling strength of around ρ 1 ∼ 10 −1 . This can be understood in the following way: as ρ 1 increases, two competing effects occur. The thermal mass term increases, such that the critical temperature is lower. At the same time, the Coleman-Weinberg potential contributes an effective interaction term, which drives the value of φ C smaller. The CW potential depends on ρ 2 1 , while the thermal potential depends on ρ 1 in the high temperature limit. Ultimately, the behavior of φ C /T C results from the balance of the two effects.

Gravitational Wave spectrum
To find the thermal parameters governing the phase transition, we find classical solutions to the Euclidean equations of motion (the scalar bounce solution), which describe the nucleation of an O(3) bubble of the true vacuum in a medium of the false vacuum. We solve the Euclidean equations of motion by varying the initial conditions via a simple bisection method. We give more details on our calculation of the thermal parameters in appendix B. The thermal parameters can be used to predict the stochastic gravitational wave spectrum using a combination of analytic and lattice studies, as reviewed in appendix C.
Informed by the results in Fig. 2, we vary the portal coupling between 0.1 < ρ 1 < 0.5 and the ratio of scales 2 < v/Λ < 6 and calculate the nucleation temperature T N as well as the thermal parameters α and β/H. The results of our parameter scan are smoothed using a local quadratic regression with tri-cube weights [62]. We show the behaviour of the thermal parameters with v/Λ and ρ 1 in Fig. 3. It is seen that the latent heat α asymptotes for large v/Λ and reaches a maximum for ρ 1 ∼ O(10 −1 ) which is expected from the behaviour of the  order parameter shown in Fig. 2. Likewise, the speed of the transition β/H also reaches a minimum value for v/Λ ∼ 5 and coupling ρ 1 ∼ 3 × 10 −1 . The nucleation temperature is mostly determined by the scale M P S , but shows a small dependence on the coupling strength ρ 1 .
The thermal parameters can be used to find the stochastic gravitational wave spectra, a calculation we review in appendix C. We plot contours of peak amplitude of the gravitational wave spectra in the (β/H, α) plane in Fig. 4, with a selection of benchmark points from our study. The benchmarks shown here represent an evenly spaced grid with ρ 1 = {1, 2, 3, 4, 5} × 10 −1 and v/Λ = {2, 3, 4, 5, 6}. The colour scaling in this plot gives the frequency at the peak of the sound wave spectrum, which is the dominant contribution.
For reference, the thicker dashed line in Fig. 4 gives the anticipated peak sensitivity of the Einstein telescope [26]. We expect the sound wave peak to be visible at the Einstein telescope [26] for β/H ∼ O(10 3 ), T n 10 5 GeV and α 0.07 [47] with the Cosmic Explorer [31] allowing slightly higher values of (T n × β/H).  Figure 4. The sound wave (in black) and turbulence (in blue) spectra of a Pati-Salam phase transition for the benchmark points described in the text. Here we have assumed v w = 1 as motivated in appendix C. To find the peak of the turbulence spectrum, which depends on two scales, we used the fiducial value T N = 10 5 GeV. The colour scaling denotes the peak frequency for the sound wave spectrum; the peak frequency of the turbulence spectrum is expected to be smaller but of the same order of magnitude for these benchmarks. The thicker line shows the peak sensitivity of the Einstein telescope [26].

Complementarity with low energy probes
GeV are Using these values in this section we will study some low-energy signatures of this model, such as neutrino masses, lepton flavour violation, collider searches and proton decay.

Neutrino masses
As required by the observation of neutrino oscillations [63][64][65][66][67], left-handed neutrinos have nonvanishing masses. However, studies of the CMB by the Planck satellite have imposed a strong upper limit on the sum of the neutrino masses m ν < 0.23 eV [68]. Left-right symmetric models, such as the intermediate step of the model described in Section 2, naturally contain a right-handed neutrino field N and a left-handed triplet δ L , which can make the active neutrinos light via type I and type II seesaw mechanisms [53,54]. The mass matrix of neutrinos in this scenario is where M D is the Dirac-type mass of the neutrinos and M L and M R are the left and righthanded Majorana masses. The former mass arises from the Yukawa coupling of the neutrino field to the SM Higgs after electroweak symmetry breaking, whereas the latter masses are generated dynamically through the vacuum expectation values of δ L/R , denoted by v L and v R respectively. In the limit of small active-sterile mixing, we can write the light left-handed neutrino masses as where y ν is the Yukawa coupling of the neutrinos and ζ is the coupling of the scalar triplets δ L and δ R to the lepton fields, which we have taken to be equal, ζ L = ζ R = ζ as a relic of D parity from the GUT scale. In models with breaking of LR manifest symmetry the vev v L can be obtained as a function of the scale of the scale of D-parity breaking [18,69,70]. In this scenario, D-parity is broken already at the GUT scale, and the only coupling of δ L to the field responsible is through Ω R . Hence we can express v L as [18] v with η the coupling between δ L and Ω R (see Appendix A) and M the dimensionful coupling between Ω R and the D-parity breaking field. Identifying M LR ≡ v R , we can rewrite eq.(4.4) as Electroweak precision data restricts v L 5 GeV in order to keep the electroweak ρ parameter under control [71]. This translates into an upper limit on η M 2 8.2 × 10 14 GeV −1 . The parameter M controls the splitting between the GUT scale and the mass of Ω R , which we want to remain at around M P S (c.f. eq. (2.4)). M is then constrained as M 1.6 × 10 −12 GeV and it thus forces a strong upper bound on η 2.03 × 10 −9 .
The parameters ζ and y ν are unconstrained in this model, but M LR depends on the Pati-Salam scale M P S via gauge coupling unification. For the gravitational wave scenario studied in the previous section, M P S = 10 5 GeV, M LR ∼ 12 TeV and ρ 1 ∼ 0.3. For the maximum values allowed for M and η, the condition that the sum of neutrino masses is below the CMB limit [68] becomes y ν > ζ(−1.54 × 10 −11 + 1.001ζ) (4.7) There is another contribution to neutrino masses that we have not considered here arising from loop corrections involving heavy leptoquarks [72]. These contributions are rather small and do not modify the conclusions of Eq. (4.6) significantly. Therefore we will not discuss them any further.
The mixing of sterile to active neutrinos Θ is the source of many of the contributions from heavy neutrinos to low energy observables, including electroweak precision observables (EWPO), which set a upper limit of |Θ| 2 10 −3 . This mixing is given by Θ = M D M −1 R , c.f. eq. (4.3), which using the parameters ζ and y ν transforms into y ν < 1.562 ζ. (4.8)
In this model we have made the simplifying assumption that M δ L ∼ M δ R ∼ M LR , and also we have that M W R = 1 2 g R M LR and M N = ζM LR . If we take the scenario that maximizes the detection of gravitational waves, M P S = 10 5 GeV, then M LR = 12.15 × 10 3 GeV, g L = 0.627 and g R = 0.376. The active-sterile mixings are not fixed by the scenario, but electroweak precision data has put an upper limit on their values which, as a conservative limit, we can take as |Θ αI | 2 < 10 −3 [92]. This results in the branching ratios which then gives an upper limit for ζ so as to satisfy the limits in eq. 4.9, ζ < 0.113, due to the most constraining of the observables, namely µ → eee. Future experiments measuring µ − e conversion, such as COMET [93] and Mu2e [94] aim to reach the limit of R N (µ − e) < 10 −16 . A positive signal from either of those experiments would fix ζ for this model, which would strengthen our predictions and motivation for the complementarity with gravitational wave detection. If no such signal is found, the new limits would further constrain the value of ζ. Taking the form factor X Al = 0.8 [18] this new upper limit would drop to ζ < 0.093.

Collider Searches
Low scale Pati-Salam and left-right symmetric models predict light exotic particles that can be visible at the LHC. In particular, the lightest exotic states produced in our model after LR symmetry breaking are the heavy right-handed neutrinos N j , the right-handed gauge bosons W R and Z R and the left-handed scalar triplet δ L , with masses in (2.6).
Right handed neutrinos can be produced directly on-shell at the LHC from the decay of a W L boson. The primary process for detection of right-handed neutrinos at the LHC is pp → W → N l → W ll → lljj, where the two final state leptons have the same sign [95]. ATLAS and CMS reported exclusion limits on searches for same sign dilepton final states for mass ranges of 100 GeV < M N < 500 GeV [96] and 20 GeV < M N < 1600 GeV [97], respectively. In the heavy mass range, above the Z resonance, M N > 90 GeV, CMS has the strongest exclusion power which is almost linear in the M N − |Θ eN | 2 plane, so its limit can be approximated as M N |Θ eN | 2 1.5 TeV. (4.12) This can be translated to the parameters of our model using the relations M N ζv R and Θ eN ∼ yν v SM ζv R and for the chosen value of M P S = 10 5 GeV, as y 2 ν < 1.98 × 10 4 ζ 3 for ζ > 7.4 × 10 −3 . (4.13) In addition to the same sign dilepton search, CMS reported results on searches for heavy neutrinos in three lepton final states [98]. The limits of the search for M N > 100 GeV are rather similar to the dilepton search. For smaller masses below the Z resonance, the exclusion limits of this search are among the strongest in the literature on par with the results of DELPHI [99], and it effectively excludes all neutrino masses for Θ eN > 10 −5 . So we have the constraint y ν < 0.156 ζ for ζ < 7.4 × 10 −3 . (4.14) In the case that the gauge boson W R can be produced at the LHC, another channel opens for the production of right-handed neutrinos where the W R takes the place of W L in the decay chain. In this channel, the two leptons in the final states can have either the same or opposite signs, depending on the Majorana or Dirac nature of the neutrinos [70]. Both ATLAS and CMS reported strong exclusion limits for W R and N in searches with two same and opposite sign leptons and two jets final states [100,101]. The limits from both experiments reach up to M W R > 4.7 TeV for 500 GeV< M N < 3 TeV, for the simplified model where g L = g R and maximal coupling |Θ| = 1.
However, in the cases where g L = g R , as it is our model, the constraint is slightly relaxed. In order to assess the effect of these searches on our model we make a very rough comparison of the number of events predicted in our model for the same-sign eejj signal region with the measured data by ALTAS and CMS at 36 fb −1 [100,101]. The cross-section of this model for this process can be estimated to be (to leading order in ζ) 6 σ(pp → W R → eejj) ≈ 5.356 y −2 ν ζ 4 fb (4.15) ATLAS and CMS reported a number of observed events of 11 and 4, and predicted background events 11.2 and 2.6, respectively. Using the reported efficiencies for the high W R mass region of 0.54 (ATLAS) and 0.57 (CMS) we find that, at 95% CL y ν < 4.010 ζ 2 − 77.222 ζ 4 (ATLAS), y ν < 4.500 ζ 2 − 86.539 ζ 4 (CMS). (4.16) The other heavy gauge boson in the theory with a mass low enough to be relevant for collider searches is Z R . High mass resonance searches in the dilepton invariant mass from ATLAS and CMS have imposed strong constraints on the mass of Z R [102,103]. Both experiments give a simplified model limit in the range M Z R > (3.5, 4.0) TeV. It has been shown that for models with g L = g R the limits on Z resonances are much weaker [104]. In any case, for M P S = 10 5 GeV, the mass of M Z R in our model is fixed to M Z R ≈ 14.5 TeV and hence it is not affected from the current experimental limits.
Finally let us consider searches for triplet scalar bosons δ L . Most relevant to us are searches for doubly charged scalar bosons, δ ++ L , in same sign diboson final states [105]. The lower limits set by such searches are of the order of a few hundred GeV, depending on the vev v L . Similarly to the Z R case above, for M P S = 10 5 GeV, M δ L ≈ M LR = 12.152 TeV, and thus the limits do not affect the outcome of this model. So far the LHC experiments have reported analyses on 36 fb −1 of collected data. The increased sensitivity of future upgrades of the LHC will impose stronger constraints on the masses of exotic states. The high-luminosity LHC (HL-LHC) is projected to collect up to 1 ab −1 of data at 14 TeV, and the hypothetical upgrade to the Very Large Hadron Collider (VLHC) will push the energy frontier to 100 TeV with a projected luminosity of 10 ab −1 [106]. Discovery of any of N , W R , Z R or δ L in either HL-LHC or VLHC would conclude in strong evidence towards a LR symmetric model at low scales, which motivates a low scale PS phase transition leading to a GW spectra observable in the next generation of experiments.

Proton decay
Unified theories typically introduce baryon number violating operators, which render the proton unstable and can lead to rapid proton decay [107][108][109]. This is certainly true for SO(10) models, whose off-diagonal gauge and scalar bosons couple to both quarks and leptons an can mediate nucleon decays [20,110,111]. This is not the case, however, for Pati-Salam models, where the gauge sector preserves baryon and lepton number independently and only selected scalar sectors can mediate the transition, none of which we include in our model [22].
Therefore, the only source of proton decay in our model arises from the leptoquarks at the GUT scale. The half-life of the proton in this scenario, with mass m p , then can be approximated as [20,112] where M X ∼ M GU T is the mass scale of the mediator and λ X its coupling to the quarks and leptons, which corresponds to g GU T for a gauge mediator. Proton decay transitions can occur in a number of different channels, e.g. p → e + π 0 , p → e + K 0 , etc. [113,114]. The most constraining limit was imposed by the Super-Kamiokande collaboration to the process p → e + π 0 , with a half-life lower bound of τ p > 1.29 × 10 34 years [115].
In our particular scenario, with the optimal PS scale for gravitational wave detection, M P S = 10 5 GeV, gauge coupling unification fixes M GU T ∼ 1.9 × 10 16 and g GU T ∼ 0.96, which gives a proton half-life of τ p ∼ 7 × 10 35 years, larger than the experimental limit. This prediction for proton decay is not too far from the SuperK bound and in fact it is fairly close to the projected limit expected to be reached by HyperK [116] of τ > 1.3 × 10 35 years. A positive measurement of proton decay is the smoking gun for unified theories, in particular if the measured decay rate is close to the predicted in our model, it would further motivate the scenario with M P S = 10 5 GeV where the peak amplitude of GW spectra is within sensitivity of the Einstein telescope. Otherwise, if proton decay is not observed, a stronger upper limit of the half-life of the proton would fall within range of our prediction and therefore a more detailed calculation of the decay rate and RGE evolution would need to be performed in order to assess the survivability of the model.

Discussion and conclusion
The Pati-Salam phase transition is a unique candidate for a gravitational wave spectrum which peaks within the frequency windows of ground-based interferometer experiments. If such a signal is observed, complementarity with low-energy experiments can be used to probe the Pati-Salam parameter space.
The strength of the phase transition and the corresponding gravitational wave signal depend most importantly on the degrees of freedom with a large coupling to the broken direction. Therefore, we considered an effective model with four free parameters: the PS gauge coupling g 4 , the PS scale M P S , the portal coupling ρ 1 , and the ratio of parameters in the tree-level scalar potential v/Λ. We found that an observation of a broken power-law spectrum of gravitational waves which peaks for f ∼ [10 − 1000] can be explained by a Pati-Salam model with scale M P S ∼ 10 5 GeV. An argument from gauge coupling unification fixes the Pati-Salam coupling g 4 as a function of this scale. The amplitude of the power spectrum, then, is a function of the portal coupling and the zero-temperature combination (v/Λ). As was demonstrated in Section 3, the peak of the spectrum may be observable if at least one of the portal couplings is sizable (ρ 1 0.1), and for particular zero-temperature parameters in the scalar potential (v/Λ 2).
For a Pati-Salam scale of M P S ∼ 10 5 GeV, many of the low energy constraints described in Section 4 impose limits on the parameters y ν and ζ. 7 We summarize these constraints in Fig. 5, including collider constraints for decays of heavy neutrinos (blue and orange), decays of W R (green), the cosmological limit on the neutrino masses (purple), LFV constraints (red) and the limit from EWPO (brown). 8 As can be seen in the figure, a large part of the parameter space is excluded by several searches. However, there is still a narrow band where this model parameters are allowed. The future projections of several experiments are depicted with dashed and dotted lines, with the projected limit on µ − e conversion from COMET and Mu2e in dashed red, and the limits for W R searches in dashed green (HL-LHC) and dotted green (VLHC). These future searches will be able to explore the parameter space more thoroughly and further constrain the model. The included set of low energy probes is but a subset of the possible relevant phenomenological observables of PS and LR models, which we have chosen to elucidate the complementarity with GW searches. We defer the computation of other relevant observables, such has neutrinoless double beta decay or electric dipole moments to further work.
Finally we note that phase transitions in models of Grand Unified Theories are often associated with the formation of cosmic defects. One-dimensional defects, cosmic strings, may decay into gravitational radiation [117,118]. However, cosmic strings associated with the energy scales studied in this work -M P S ∼ 10 5 GeV -will have a dimensionless string tension of Gµ ∼ 10 −29 and will therefore not lead to any observational signatures. Primordial monopoles created during the PS phase transition can be reduced to acceptable limits during late time inflation [24]. Light monopoles may be produced at colliders [119], however PS monopoles have a mass of the order ∼ 10 6 GeV [120] which is beyond the reach of current experimental searches by ATLAS [121] and MoEDAL [122].

A Scalar potential
The scalar potential of the Pati-Salam model at zero temperature can be written as where V Ψ refers to the terms in the potential that contain only the field Ψ, and these are Since the fundamental representation of SU (2) is real, the fieldΦ = τ 2 Φ * τ 2 transforms as Φ. Henceforth we call Φ 1 = Φ and Φ 2 =Φ. The self-interaction term V Φ now looks like The term V ΞΞ contains interactions among Ξ (1,2,3) of the form where (i, j, k, l) = (1, 2, 3) and not all i, j, k, l are equal. The term V ΞΨ contains the portal couplings of the fields Ξ i with the rest and they are of the type The term V ΩΨ is fairly similar to V ΞΨ and looks like The last two terms, V ∆∆ and V ∆Φ have the same form as in LR symmetric models This scalar potential contains all possible terms allowed by the gauge symmetries. In the work above we have chosen to remove a few of them setting their couplings to zero, e.g. all portal couplings vanish ρ ij = ρ ij = ρ ij = ρ ijk = ρ ijk = 0 save for the first one, that we have renamed in the text as ρ 11 = ρ 1 = 0.

B Thermal parameters
The nucleation temperature, which approximates the collision temperature very well when the phase transition occurs quickly, is conventionally defined as the temperature for which a volume fraction e −1 is in the true vacuum state. This corresponds approximately to where S E is the Euclidean action evaluated on the bounce which approximates a tanh function and can be solved by bisection or perturbing a tanh ansatz [123,124]. 9 The speed of the phase transition can be calculated from the rate of change of the euclidean action Lastly, the most important parameter governing the amplitude of the relic gravitational waves will be the latent heat released in the transition (normalized to the radiation density) where ρ * = π 2 g * T 4 /30. 9 In our analysis, we assume a radiation dominated universe to relate the nucleation temperature and time.
where is the fraction of the energy in the plasma is expressed as turbulence; in our results, we use = 0.05. The spectral form is given by . (C.7) The Hubble rate at the transition temperature as well as the peak frequency are given by h * = 16.5µHz T N 100GeV g * 100 1/6 (C.8) β H g * 100 1/6 (C. 9) respectively. Recent work [30] has shown that it is difficult to satisfy the criteria that the phase transition completes and the sound waves last longer than a Hubble time. The consequence of this is that the sound waves are likely overestimated and the turbulence is likely underestimated. The suppression of the sound wave peak is naively estimated to be suppressed by a factor [30] HR where we have used the relation for the rms fluid velocity U F ∼ 3 4 κ f α. The last equality holds for the points in our scan, which have 0.07 < α < 0.09. However, the precise suppression factor is subject to future lattice simulations. Similarly, [30] argued that the turbulence factor has been underestimated, however much uncertainty remains about the precise form of the turbulence spectrum in general.