Supercooling in radiative symmetry breaking: theory extensions, gravitational wave detection and primordial black holes

First-order phase transitions, which take place when the symmetries are predominantly broken (and masses are then generated) through radiative corrections, produce observable gravitational waves and primordial black holes. We provide a model-independent approach that is valid for large-enough supercooling to quantitatively describe these phenomena in terms of few parameters, which are computable once the model is specified. The validity of a previously-proposed approach of this sort is extended here to a larger class of theories. Among other things, we identify regions of the parameter space that correspond to the background of gravitational waves recently detected by pulsar timing arrays (NANOGrav, CPTA, EPTA, PPTA) and others that are either excluded by the observing runs of LIGO and Virgo or within the reach of future gravitational wave detectors. Furthermore, we find regions of the parameter space where primordial black holes produced by large over-densities due to such phase transitions can account for dark matter. Finally, it is shown how this model-independent approach can be applied to specific cases, including a phenomenological completion of the Standard Model with right-handed neutrinos and gauged B - L undergoing radiative symmetry breaking.


Introduction
First-order phase transitions (PTs) leave potentially observable footprints that can be evidence for new physics because the Standard Model (SM) does not feature this type of transitions.
One example of such footprints is the spectrum of gravitational waves (GWs) produced by firstorder PTs.GW astronomy has become an extremely active and exciting branch of physics after the discovery of the GWs from binary black hole and neutron star mergers [1][2][3].Very recently, the interest in this field has been further boosted by the detection of a background of GWs by pulsar timing arrays: these include the North American Nanohertz Observatory for Gravitational Waves (NANOGrav), the Chinese Pulsar Timing Array (CPTA), the European Pulsar Timing Array (EPTA) and the Parkes Pulsar Timing Array (PPTA) [4][5][6][7].
First-order PTs always take place when the corresponding symmetry breaking is mostly induced (and masses are generated) radiatively, i.e. through perturbative loop corrections [20].The seminal work on this radiative symmetry breaking (RSB) is Ref. [21] by Coleman and E. Weinberg, which studied a simple toy model (see also Ref. [22] for a recent analysis).The Coleman-Weinberg work was then extended to a more general field theory by Gildener and S. Weinberg [23].Furthermore, thanks to an approximate scale invariance, in the RSB scenario the first-order PTs feature a period of supercooling, when the temperature dropped by several orders of magnitude below the critical temperature [20,24].
Supercooling in the RSB scenario allows us to be confident about the validity of the one-loop approximation and the derivative expansion [20].Moreover, it also ensures that the gravitational corrections to the false vacuum decay are amply negligible whenever the symmetry breaking scale is small compared to the Planck mass, which is, of course, the most interesting case from the phenomenological point of view.
In [20] it was shown that a model-independent description of PTs and the consequent production of GWs in the RSB scenario is possible in terms of few parameters (which are computable once the model is specified) if enough supercooling occurred.Ref. [20] provided a sufficient condition on the amount of supercooling, which ensures that the model independent description is valid.This led to a "supercool expansion" in terms of a quantity that is small when supercooling is large enough.
In this work we investigate whether this condition can be weakened and, if so, how to systematically perform a corresponding "extended supercool expansion".A weaker condition on supercooling is useful because it allows us to describe a larger class of models through the modelindependent approach, without repeating the study of the PTs every time.
Once one establishes that such extended supercool expansion can be performed, one can use it to describe in a model-independent way not only the spectrum of GWs, but also the production of PBHs due to the first-oder PTs.In particular, the amount of dark matter in the form of these PBHs can be determined in terms of the few parameters, which are computable once the model is specified.Moreover, one can identify the model-independent regions of parameter space corresponding to the GW background recently detected by pulsar timing arrays, as well as those excluded by the runs [40] of the Laser Interferometer Gravitational-Wave Observatory (LIGO) [41,42] and Advanced Virgo [43] and those within the reach of future GW detectors.These include the Laser Interferometer Space Antenna (LISA) [44], Cosmic Explorer (CE) [45,46], Einstein Telescope (ET) [47][48][49], the Big Bang Observer (BBO) [50][51][52], the Deci-hertz Interferometer Gravitational wave Observatory (DECIGO) [53,54], etc.
The paper is structured as follows.
• In Sec. 2 we introduce the general class of theories featuring RSB, where the masses are mostly generated radiatively.These theories may include the SM or may be though of as "dark" sectors weakly coupled to the SM.In the same section the model-independent description of RSB and the corresponding PTs in the supercool expansion is reviewed.This is necessary to render the subsequent original sections understandable and to establish our conventions.
• In Sec. 3 we investigate when and how one can extend the validity of the model-independent description of PTs to a larger class of RSB models by weakening the condition on the amount of supercooling.
• Sec. 4 is devoted to the possible applications of such extended supercool expansion to the production of GWs and PBHs through first-oder PTs.We also include in the discussion the background of GWs recently discovered by pulsar timing arrays.
• Since the usefulness of the model-independent approach studied here is mainly due to the fact that one can avoid repeating the analysis of the PT in each RSB model, in Sec. 5 it is shown how to apply it to specific models by considering a couple of examples: a simple illustrative one and a phenomenological completion of the SM featuring right-handed neutrinos below the EW scale and the gauging of the difference B − L between the baryon and lepton numbers, which undergoes RSB.In these examples the accuracy of the extended supercool expansion is also studied.
• Sec.6 provides a detailed summary of the main original results of this paper and the final conclusions.

Supercool expansion: a recap
In this section the important properties of RSB and the supercool expansion are summarized.This is necessary to explain in a clear way the original results of the subsequent sections.The reader can find in Ref. [20] the proof of any non-trivial statement that is present in this section.We will also define our basic conventions here.
In the RSB scenario the sector responsible for the symmetry breaking is (at least approximately) classically scale invariant and it is thus described in general by a Lagrangian of the form while gravity is assumed to be Einstein's gravity at the energies that are relevant for this work 1 .
Here we consider generic numbers of real scalars ϕ a , Weyl fermions ψ j and vectors V A µ (with field strength F A µν ), respectively.The V A µ are gauge fields and allow us to construct the covariant derivatives D µ ϕ a and D µ ψ j .Of course, in (2.1) all terms are gauge-invariant.Also, the Y a ij are the Yukawa couplings and V ns (ϕ) has the general form where λ abcd are the quartic couplings.
In the RSB mechanism mass scales emerge radiatively from loops because there may be some specific energy at which the potential in (2.2) develops a flat direction, ϕ a = ν a χ, where ν a are the components of a unit vector ν, i.e. ν a ν a = 1, and χ is a single scalar field.So, the RG-improved potential V along ν reads Having a flat direction along ν for the RG energy µ equal to some specific value μ means λ χ (μ) = 0.
Including the one-loop correction the quantum effective potential can be written where and χ 0 is related to μ through a renormalization-scheme-dependent formula.The field value χ 0 is a stationary point of V q and, when β > 0, is also a point of minimum.Therefore, when the conditions are satisfied one has a minimum of V q at a non-vanishing value χ 0 of χ and the fluctuations of χ around χ 0 have mass This non-trivial minimum can generically break global and/or local symmetries and generate the particle masses, with χ 0 being the symmetry breaking scale.EW symmetry breaking can also be induced when there is a term in L of the form where H is the SM Higgs doublet and λ χh is some coupling.Indeed, by evaluating this term at χ = χ 0 we obtain the Higgs squared mass parameter So when λ χh (μ) > 0 the masses of the SM elementary particles are generated.Recalling the wellknown formula that relates µ2 h and the Higgs mass, it is clear that we cannot use this mechanism to generate µ 2 h when χ 0 is much below the EW scale and demand the validity of perturbation theory at the same time.Of course, it is still possible that the SM with an explicit scale-symmetry breaking parameter is weakly coupled to a scale-invariant sector that features RSB.In this case perturbation theory can be compatible with a χ 0 much smaller than the EW scale.
Including now thermal corrections, the general expression of the effective potential V eff at finite temperature T is (in the Landau gauge and at one-loop level 2 ) where the m b and m f are the background-dependent bosonic and fermionic masses, respectively, the sum over b runs over all bosonic degrees of freedom and n b = 1 for a scalar (we work with real scalars) and n b = 3 for a vector degree of freedom.In (2.10) the sum over f , which runs over the fermion degrees of freedom, is multiplied by 2 because we work with Weyl spinors.Also, the thermal functions J B and J F are where and γ E is the Euler-Mascheroni constant (the derivation of the expansions above are given in [63]).In Eq. (2.10) we have included a constant term Λ 0 to account for the observed value of the cosmological constant when χ is set to the point of minimum of V eff .
The PT associated with a radiative symmetry breaking always turns out to be of first order.The absolute minimum of the effective potential is at χ = 0 for T larger than the critical temperature T c , while, for T < T c , is at a non-zero temperature-dependent value.In the latter case the decay rate per unit of spacetime volume, Γ, of the false vacuum into the true vacuum can be computed with the formalism of [64][65][66][67]: where S is the action evaluated at the bounce, i.e. the solution of the differential problem [68] χ A dot and a prime denote a derivative with respect to the Euclidean time t E and the spatial radius r ≡ √ ⃗ x 2 , respectively.A particular solution of (2.15)-(2.16) is the time-independent bounce, for which If the time-independent bounce dominates, the decay rate is [66,67] and S 3 evaluated at the time-independent bounce can be written as follows: As long as perturbation theory holds, in a generic theory with RSB, Eq. (2.1), when T goes below T c the scalar field χ is always trapped in the false vacuum until T is much below T c , i.e. the universe always features a phase of supercooling.If this process is strong enough, in a generic theory of the form (2.1) the full effective action for relevant values of χ can be described by three and only three parameters: χ 0 , β and a real, non-negative and χ-independent quantity g, which plays the role of a "collective coupling" of χ with all fields of the theory.This is possible because the field value χ b around the barrier, which can be defined by Veff (χ b , T ) = 0, turns out to be small compared to T for large-enough supercooling: such that the small-field expansions in (2.11) and (2.12) can be truncated as and the logarithmic term in V q can be written as follows: where m and λ are real and positive functions of T defined by For this effective potential the tunneling process is dominated by the time-independent bounce.The bounce action S 3 computed with the effective potential in (2.27) turns out to be (see also [69,70] for previous calculations).
In general the nucleation temperature T n can be defined as the temperature for which Γ/H 4 I = 1, so, using the fact that the decay is dominated by the time-independent bounce, at where is the Hubble rate associated with the exponential expansion of space that takes place during supercooling.By using the expression of S 3 in (2.29) and the definitions in (2.28) one finds the following solution for T n with and MP is the reduced Planck mass.
In general, the strength of the PT is measured by the parameter α defined as [71,72] where g * (T ) is the effective number of relativistic species at temperature T , in the presence of supercooling and ⟨χ⟩ is the point of absolute minimum of the full effective potential.For an RSB PT α ≫ 1.
Another important parameter to analyse the production of GWs and PBHs is the inverse duration β of the PT that, in models with supercooling, is [34,73,74] where t n is the value of the time t when T = T n .Recalling that the tunneling process is dominated by the time-independent bounce, where H n ≈ H I is the Hubble rate when T = T n .
Note that here we are relying on a small ϵ expansion (a "supercool expansion") and what we have done so far is the analysis at leading order (LO), that is modulo terms of relative order √ ϵ.Including these terms and treating them perturbatively would mean working in the supercool expansion at next-to-leading order (NLO).This can be done by including the term of order x 3/2 in the expansion of J B (x), Eq. (2.11), and is justified if ϵ is small.In Sec. 3 we will explain how to extend the supercool expansion to order-one values of ϵ.
The effective potential at NLO, therefore, includes a cubic-in-χ term and reads where m 2 and λ are defined in (2.28), and g is a real, non-negative and χ-independent parameter defined by This is an extra parameter that is needed for a model-independent description of this scenario at NLO.In general we have g ≤ g. (2.41) To understand why the term cubic in χ in (2.38) can be considered as a small correction in the supercool expansion, one can rescale χ → χ/ √ λ in the bounce action, Eq. (2.14), to obtain Since we eventually need to set T = T n , the term proportional to k has relative order at most √ ϵ times a small number ≈ 1/( √ 2π) (where the LO result S 3 ≈ 4πgT /( √ 12λ) has been used).Working with the supercool expansion at NLO (i.e.treating the cubic term in (2.38) perturbatively at first order) one can then find corrected analytical expressions for T n , S 3 and β, which depend on the extra parameter g.For example, where and χ LO is the LO bounce configuration.Of course, one can then go ahead and compute smaller and smaller corrections.

Extending the validity of the supercool expansion
In this section we study when and how one can extend the validity of the supercool expansion to cases in which ϵ ∼ 1. (3.1)

Several degrees of freedom
The expansion developed in [20], which we have reviewed in Sec. 2, generally works for ϵ small.However, it also holds for values of ϵ of order one if there are several degrees of freedom, say N , with dominant couplings (all of the same order of magnitude, say τ ) to the flat-direction field χ.Indeed, in this case g defined in (2.21) scales as g ∼ √ N τ , while g defined in (2.40) scales as g ≲3 √ N τ , and so g3 /g 3 ≲ 1/ √ N : the inequality here is due to the fact that g receives contributions only from bosons, while both fermions and bosons contribute to g.As a result, the extra cubic term in the bounce action of Eq. (2.42) gets a further suppression factor (see Eq. (2.39)), which is at least as small as 1/ √ N .On the other hand, • since 1/X = 6 βϵ/g 4 , for order one ϵ the quantity 1/X is still small because β is loop suppressed and so the approximation in (2.25) is still good, • truncating the small-x expansions in (2.11) and (2.12) up to the x 3/2 term is still justified because the higher-order terms involve smaller and smaller coefficients 3 , with the coefficient of the O(x 2 ) term being already quite small, ∼ 1/32.

Improved supercool expansion
On the other hand, if the number of degrees of freedom with a dominant coupling to χ is too small, one instead has g ≈ g and, in this case, the expansion of Sec. 2 breaks down for order 1 values of ϵ (although it still holds for small ϵ).

Bounce solution and action
In order to extend the class of theories that can be described by the supercool expansion one is, therefore, interested in including the cubic term in (2.38) in the non-perturbative computation of the bounce action and treating the other corrections as perturbations (indeed, they are still small as long as ϵ is at most of order one, as we have seen in Sec.3.1).We will refer to this improvement as the "improved supercool expansion".Let us explain how to construct it.The expression of Veff in (2.38), together with the form of the bounce problem in (2.15)-(2.16),tells us that the characteristic bounce size R b is of order R b ∼ 1/m(T ) ≳ 1/T , where in the second estimate we have used the perturbativity condition that g is not too large.Indeed, the bounce size can be read from the large-r limit of the bounce solution and in this limit the last condition in (2.16) tells us that only the quadratic term in (2.38) matters.Therefore, the bounce solutions are approximately time-independent even including the cubic term in (2.38).
Looking then at (2.18) and redefining [22] r ≡ lρ and χ ≡ ξφ one obtains the bounce action for the new radial variable ρ and the new field φ where By evaluating at the bounce solution one then obtains, like in (2.20), a simplified bounce action with m and k defined in (2.28) and (2.39), respectively, gives and λ defined in (2.28).The quantity λ can also be rewritten by using (2.28) and (2.39) as follows which depends on T only through log(χ 0 /T Ṽeff , which appears in the integrand of the bounce action in (3.4).
We are not able to find the analytic dependence of the bounce action S 3 on κ.However, one can compute the bounce and the corresponding S 3 for several values of κ and then perform a fit [22,76,77].Doing so we find that reproduces the numerical calculations at the ∼ 1% level for the values of λ we are interested in.The result in (3.12) was found by [22] in a specific setup.Here its validity has been established in a model-independent way within the improved supercool expansion.

Nucleation temperature
Inserting the expression in (3.12) into the equation for the nucleation temperature T n in (2.30) leads to a complicated non-polynomial equation in λ.This equation can be partially simplified by dropping the second term on the left-hand side of Eq. (2.30), which is always negligible with respect to the first one because the semiclassical approximation requires S 3 /T large.Within this approximation the equation for λ reads where the value of c 3 is given in (2.29) and a and c are defined in Eq. (2.33) (the term 3 2 log a 2π in c can be dropped as it comes from the second term on the left-hand side of Eq. (2.30)).Here we are interested in the smallest real and positive solution λn ≡ λ(T n ) of Eq. (3.13) for which the straight line a 1 − a 2 λ reaches F ( λ) from below in increasing λ (that corresponds to Γ reaching H 4 I from below).Clearly, such a solution does not always exist for any a 1 and a 2 .First, one must have a 1 ≤ F (0) = 9/2; second, for each given a 1 the parameter a 2 must me smaller than a certain critical value ā2 (a 1 ), which is given in the inset of the right plot of Fig. 2. Fig. 2 also shows as a function of a 1 and a 2 the solution λn (when it exists), which has been obtained numerically.Tables containing the numerical determination of ā2 as a function of a 1 and of λn as a function of a 1 and a 2 can be found at [78].Once we fix the parameters g, β, χ 0 and g the quantities a 1 and a 2 as well as λn and thus the nucleation temperature T n are fixed.Using the obtained solution λn we checked that the PT strength parameter α is large in an RSB PT for realistic and perturbative values of the parameters, even in the improved supercool expansion discussed in this paper.Thus, the plasma effects (such as those studied in Refs.[79,80]) can be neglected in this particular scenario.
One might wonder whether the effect of the spacetime curvature due to H I ̸ = 0 can alter the decay rate.In standard Einstein gravity, this may happen if T n is so small to be comparable with H I .We checked that, whenever a solution for λn exists, this never happens, at least for realistic and perturbative values of the parameters.On the other hand, if a solution for λn does not exist, the effect of the spacetime curvature, along with quantum fluctuations, can eventually become important in the decay rate [33,[81][82][83] and lead to the completion of the transition.

Duration of the phase transition
Using the expression of S 3 in (3.12) and dropping the last term in (2.37), which is negligible in the semiclassical approximation as we have pointed out in Sec.3.2.2,we obtain a formula for the inverse duration of the PT: where F ′ is the derivative of F defined in Eq. (3.13) with respect to λ; note that F is a monotonic decreasing function of λ so −F ′ > 0.
Figs. 3 and 4 show β/H n (computed with the improved supercool expansion) as a function of g and β for fixed values of χ 0 .There g has been set equal to g: when g is significantly lower than g the expansion developed in [20] works well as discussed in Sec.3.1 and there is no need to resort to the improved supercool expansion.Moreover, in Figs. 3 and 4 only values of g and β with ϵ < 3 are displayed 4 : indeed, for large values of ϵ one needs to take into account the higher-order corrections for a good accuracy.In Figs. 3 and 4 β/H n never vanishes although there are values of g and β for which β/H n ∼ 1.The relevant solution of the nucleation temperature equation in (3.13), i.e. λn , ceases to exist before β/H n vanishes.As commented in the last paragraph of Sec.3.2.2,when the solution λn does not exist the effect of the spacetime curvature, as well as quantum fluctuations, can eventually become important in the decay rate.

Applications
Let us now apply the improved approximations developed in Sec. 3 to the production of GWs and PBHs.

Gravitational Waves
In the RSB scenario the dominant source of GWs are vacuum bubble collisions: the energy density of the space where the bubbles move is dominated by the vacuum energy density due χ, which leads to an exponential growth of the corresponding cosmological scale factor.This inflationary behavior dilutes preexisting matter and radiation and, thus, we neglect the GW production due to turbulence and sound waves in the cosmic fluid [73,84,85].
In the presence of supercooling and for α ≫ 1 one finds the following GW spectrum due to vacuum bubble collisions 5 [73] where T r is the reheating temperature after supercooling, H r is the corresponding Hubble rate and f peak is the red-shifted frequency peak today, given by [73] f peak ≈ 3.79 β H r g * (T r ) 100 T r 10 8 GeV Hz.
Ref. [73] used the results of [86] based on the envelope approximation.This is an approximation where all the energy is assumed to be stored in the bubble walls, that are taken to be thin, and at bubble collision one uses as a source for GW production the energy-momentum tensor of the uncollided part of the bubble walls.Studying the collision of two bubbles in a scalar field model with symmetry breaking entirely due to the standard Higgs mechanism, Ref. [87,88] found that this has about 5% accuracy.For ϵ ∼ 1 this is comparable with the precision of the improved supercool expansion when one uses the approximation in (2.25) and neglects the terms 4 In the present improved approximation ϵ is computed by using Eq.(3.10) with λ = λn 5 The spectral density Ω GW is defined as usual as where ρ cr ≡ 3H 2 0 M 2 P is the critical energy density, H 0 is the present value of the Hubble rate and ρ GW is the energy density of the stochastic background. in the small-x expansions of (2.11) and (2.12) of order higher than O(x 3/2 ).In our situation the envelope approximation is expected to capture the dominant source 6 of GWs [94] because, during the exponential growth of the universe, the bubbles expand considerably and in this process the energy gained in the transition from the false to the true vacuum is transferred to the bubble But otherwise H r and T r can depend on the details of the specific model.Reheating can occur e.g.thanks to the Higgs portal coupling in (2.8) or other portal interactions such as a kinetic mixing between the photon and a dark photon (see [95] for a review) that become massive through RSB.
Note also that the dependence of Ω GW and f peak on g * (T r ) is quite weak.
Ref. [20] computed f peak and provided regions where Ω GW (f peak ) is above the sensitivities of several current and proposed GW detectors (including LIGO, Virgo, LISA, ET, CE, BBO and DECIGO); moreover, Ref. [20] found corresponding regions in the space of g, β, χ 0 and g using the supercool expansion at LO and NLO.Then here we focus on the improved supercool expansion.
Combining with the information in Fig. 3, one finds that it is possible to account for the signals detected by pulsar timing arrays for χ 0 ∼ 10 GeV and choosing the basic parameters g, g and β appropriately, as illustrated in Fig. 6.Fig. 7 shows instead the regions where Ω GW (f peak ) is above the sensitivities of Advanced LIGO's and Advanced Virgo's third observing (O3) run (left plot) and LISA with power law sensitivity [113] (right plot) for two non-vanishing values of g and using the improved supercool approximation.The regions in the left plot are thus, remarkably, ruled out.In Fig. 7 we again considered only values of g and β such that ϵ < 3. The parameter χ 0 has been chosen around 10 9 GeV in the left plot of Fig. 7 and 10 4 GeV in the right one because the corresponding f peak is then around the frequency range of LIGO-VIRGO O3 [40] and LISA [113], respectively (see Fig. 5).A χ 0 around 10 9 GeV is relevant e.g. for axion models, while a χ 0 around 10 or 100 TeV could be associated with observable physics at colliders and is relevant e.g. for supersymmetric models or low-scale unified theories such as the Pati-Salam model [114] or Trinification [115].

Primordial black holes
As shown in [20] the PT associated with an RSB is always of first order.Besides having the potential of leaving observable GW footprints, first-order PTs can also naturally generate PBHs because generically lead to large over-densities [8][9][10][11][12][13][14][15][16][17][18][19].One of the main motivations for studying PBHs is the fact that they can account for a fraction f PBH of (or even the whole) dark matter density.

Late-blooming mechanism
One of the mechanism to generate PBHs from first-order PTs is based on the presence of strong supercooling, which generically takes place in the RSB scenario and is a key property for the validity of the supercool expansion.Since the bubble formation process is statistical for both quantum and thermal reasons, distinct causal patches percolate at different times.Patches that percolate the latest undergo the longest vacuum-dominated stage and, therefore, develop large over-densities triggering their collapse into PBHs.This late-blooming mechanism has been studied in a number of papers (see e.g.Refs.[15][16][17][18][19]) and we refer the reader to these works for an introduction to this mechanism.A key feature is that the longer the supercooling period lasts (the smaller β/H n is) the more effective this mechanism is.Following the method illustrated in Ref. [19] and using the improved supercool expansion we have identified regions of the parameter space (shown in Fig. 8) for which PBHs produced by the late-blooming mechanism can account for a significant fraction of the dark matter density in a model-independent way.This was possible because, as discussed in Sec.3.2.3 and shown in Figs. 3 and 4, we can compute β/H n only in terms of the parameters g, g, β and χ 0 for largeenough supercooling: this hypothesis allows us to use the supercool expansion (in Fig. 8 the improved version is used).For all values of these parameters in Fig. 8 the PT is very strong (α > 100 for g * ∼ 10 2 ) and the improved supercool expansion gives a good approximation for the

Other mechanisms?
Several other mechanisms to produce PBHs have been proposed in the literature.Some of these are unrelated to the RSB and strong supercooled PTs and thus we do not discuss them, although, of course, they could contribute to the PBH abundance in specific models.
Another mechanism that can be a priori related to the RSB and strong supercooled PTs in a model-independent way is the one based on bubble collisions [8,9,13].However, in Ref. [13] it was pointed out that bubble collisions during a first-order PT can produce PBHs only if the bubble radii become near-horizon-sized and the bubble walls have a non-negligible thickness when they collide.In RSB PTs this PBH production mechanism is, therefore, suppressed because the bubble walls become very thin after a long period of supercooling (as discussed in Sec.4.1) and also we checked that the bubble radii never become near-horizon-sized for values of χ 0 up to 10 16 GeV.Larger values of χ 0 are not considered here as they require a UV completion of gravity.

Examples of specific models
What we have done so far is a model-independent study of phase transitions and corresponding production of GWs and PBHs in the RSB scenario, which is valid in the supercool expansion or, more generally, in the improved supercool expansion.This formalism can be applied to any RSB model featuring a large-enough amount of supercooling (ϵ at most of order 1).To illustrate the usefulness of these results here we apply them to some concrete models.

A simple model
We start with a simple toy model that can illustrate all essential features of the RSB scenario.The basic requirements of RSB is the existence of a flat direction that is radiatively broken to generate a minimum ( β > 0).This positivity condition can be satisfied by introducing a gauge group, which can generically give positive contributions to the scalar beta functions.Here we take this group to be SU (2) (an Abelian case will be discussed in Sec.5.2).The scalar fields will be organized here in a complex adjoint field A, whose no-scale potential is7 where λ 1 and λ 2 are real couplings.Therefore, there exist a non-trivial flat direction, for which A = A † , at a scale μ where λ 1 + λ 2 = 0.When A = A † the three components A k of A along the Pauli matrices, A = A k σ k /2, can always be transformed through an element of the gauge group SU(2) in a way that only one of these components is not vanishing and positive.We identify this non-zero component with χ.
where g a is the gauge coupling of SU( 2).Evaluating at the scale μ, at which λ 2 = −λ 1 , gives In order to simplify the following discussion we also assume that λ 1 ≪ g a such that we have a single coupling to deal with.In this case the massive background-dependent spectrum only features two spin-1 particles with equal mass, M V = |g a |χ.All the other masses either vanish or are negligibly smaller.So the collective coupling g defined in Eq. (2.21) turns out to be and so β = g 4 3(4π) 2 .
(5.5) Also, g defined in Eq. (2.40) reads Having determined β and g in terms of g one can now use the model-independent analysis based on the improved supercool expansion of Sec.3.2 with only two free parameters: g and χ 0 .At this point it is interesting to quantify the error that one is making in analysing this model with the standard supercool expansion at NLO of Sec. 2 rather than with the improved supercool expansion of Sec.3.2, namely treating the cubic term in (2.38) perturbatively.Fig. 9 (upper plots) shows ϵ and e nlo ≡ max(e 1 , e 2 , e 3 ), (e 1 , e 2 , e 3 ) ≡ c2 (computed for simplicity with the LO formula for T n in (2.32)).The parameter e nlo in (5.7) quantifies the above-mentioned error: the first entry e 1 in the max function is the square of the size of the second term in (2.43) relative to the first one (which is an estimate of the next-tonext-to-leading correction in treating the cubic term in (2.38) perturbatively); the second and third entries, e 2 and e 3 , are instead estimates of the error due to the approximation in (2.25) and to truncating the small-x expansions in (2.11) and (2.12) up to the x 3/2 term8 , respectively.
As one can see, although ϵ is above 1 the quantity e nlo is small, especially for smaller values of g.The reason why this happens is because here one has two massive vector fields for a total of six degrees of freedom and there is, therefore, an extra suppression of the neglected terms as explained in Sec.3.1.Looking at Fig. 9 one also sees that e 1 is larger than e 2 and e 3 , meaning that the improved supercool expansion is a better approximation than the standard supercool expansion at NLO in this case.This is not surprising because the number of degrees of freedom with dominant couplings to the flat-direction field χ is not very large (it is six) and ϵ is not smaller than one in this case.

Radiative electroweak and lepton symmetry breaking
Let us study now an example that is phenomenologically well-motivated.The SM is a very successful model but it clearly has to be extended: neutrino oscillations, dark matter and baryon asymmetry must be accounted for in a phenomenologically complete model.One of the most economical way to achieve this goal is to add three right-handed neutrinos N i featuring Majorana masses below the EW scale (see e.g.[116,117]).
The corresponding Majorana mass terms can be promoted to scale-invariant Yukawa interactions 1 2 y ij AN i N j + h.c.(5.8) in L ns matter of Eq. (2.1) by introducing a charged scalar A with a non-vanishing lepton number (here the y ij are the corresponding Yukawa couplings).Coupling A and the N i to the classically scale-invariant part 9 of the SM through renormalizable dimension-four interactions allows us to build a classically scale-invariant model of the type described in Sec. 2. In order to generate the N i Majorana masses one can then try to realize an RSB of the lepton number along the field direction 10 |A|, which, as we have seen, requires the quartic coupling λ a of the field A to vanish Figure 9: Comparison between ϵ and the error e nlo that one is making in using the standard supercool expansion at NLO of Sec. 2. The upper plots refer to the simple toy model of Sec.5.1 (e nlo is defined in (5.7)), while the lower plots regards Sec.5.2 (e nlo is defined in (5.17)).
at an energy scale where its beta function is positive (see Eq. (2.6)).However, it turns out that in this simple model the Yukawa interactions in (5.8) drives this beta function to negative values when λ a and λ ah are negligibly small.This problem can be elegantly solved by gauging the Abelian U(1) symmetry acting on A. As well known, in order to avoid any gauge anomalies, such new gauge symmetry must correspond to B − L and so we call it U (1) B−L .Therefore all leptons (including the N i ) and quarks as well as the scalar A are charged under this Abelian symmetry.A radiatively-induced vacuum expectation value of A can then generate the N i Majorana masses and induce the tachyonic Higgs mass parameter in Eq. (2.9).This classically scale-invariant model has been previously considered in Ref. [118], but without accounting for dark matter.
The Lagrangian is given by where L ns SM represents the classically scale-invariant SM Lagrangian and the L i are the three families of SM lepton doublets.Here D µ is the covariant derivative with respect to the full gauge group SU (3) C × SU (2) L × U (1) Y × U (1) B−L , i.e. the SM group times the B − L one: which involve the gluons G α µ , the triplet of W bosons W a µ as well as the gauge fields B µ and Here g m takes into account the Abelian mixing between U (1) Y and U (1) B−L .We do not propose this model as UV completion of the SM but just as an effective field theory valid up to the symmetry breaking scale χ 0 .
As discussed around Eq. (2.4), to realize RSB we need the beta function of the quartic coupling of the flat-direction field χ, in this case λ a .Using the general formalism of [119][120][121] we find the following one-loop expression Evaluating now this beta function at a scale where λ a = 0 to compute β defined in (2.5) and neglecting λ ah and y for the reasons explained above (all right-handed neutrino Majorana masses are taken below the EW scale in [116,117]) one finds (5.12) One the other hand, in this setup the background-dependent mass of the new gauge boson, Z ′ , is where we used Eq.(5.10) and the fact that |B − L| = 2 for the new scalar field A, the collective coupling g defined in Eq. (2.21) is and g defined in Eq. (2.40) is . (5.15) Therefore, β = 2g 4 3(4π) 2 . (5.16) Like in the previous section, one can now use the model-independent analysis based on the improved supercool expansion of Sec.3.2 with only two free parameters: g and χ 0 .
Let us quantify the error that one is making in analysing this model with the standard supercool expansion at NLO of Sec. 2. Fig. 9 (lower plots) shows ϵ and e nlo ≡ max(e 1 , e 2 , e 3 ), (e 1 , e 2 , e 3 ) ≡ c2 (5.17) (computed for simplicity with the LO formula for T n in (2.32)).The estimate of the error in (5.17) has been obtained like11 in (5.7).As one can see, although ϵ is above 1 the quantity e nlo is small, especially for smaller values of g.The reason why this happens is again because we have more than one massive degrees of freedom.In this case, however, we have a single massive vector field, Z ′ , rather than two like in Sec.5.1 and so the suppression of the neglected terms is slightly weaker as one can see in Fig. 9. Fig. 9 also shows that the improved supercool expansion is a better approximation than the standard supercool expansion at NLO as e 1 > e 2 and e 1 > e 3 in this case too.Again we attribute this to the fact that the number of degrees of freedom with dominant couplings with the flat-direction field χ is not very large (it is three here) and ϵ is not smaller than one in this case.Since this model is phenomenologically very well motivated, it is also interesting to study the reheating after the supercooling period.In the following discussion we focus on the decay of the flat-direction field χ coming from A into two physical Higgs bosons of mass M h ≈ 125 GeV.This process is induced by the portal interaction λ ah |A| 2 |H| 2 .Expanding around χ = χ 0 one obtains the effective interaction λ ah χ 0 δχ|H| 2 , where δχ ≡ χ − χ 0 .On the other hand, using (2.7) and (5.16) one finds g 2 4π χ 0 . (5.18) The radiative symmetry breaking of B − L also induces EW symmetry breaking according the discussion around (2.9) and a physical Higgs mass M h = √ λ ah χ 0 .So the decay rate of χ in a pair of Higgs particles χ → HH is (when M h is negligible compared to m χ ) The reheating temperature due to this channel may be computed through

.20)
But this formula is only valid if the radiation energy density ρ R does not exceed the vacuum energy density ρ V due to χ (because ρ V represents the full energy budget of the system).If this condition is not satisfied we determine T r as the maximal temperature compatible with ρ R ≤ ρ V , leading to the formula for T r in (4.4).For g of order one, g * ∼ 10 2 and χ 0 ≲ 10 5 GeV the reheating temperature is well above the EW scale and (4.4) holds such that the reheating effectively is fast.Increasing χ 0 lowers T r in (5.20) and ρ R ≤ ρ V can be satisfied.

Summary and conclusions
Let us conclude by providing a summary of the main original results obtained.
• In Sec. 3 we have significantly extended the applicability of the model-independent approach to study PTs proposed in [20], which now works for a larger class of RSB models: the amount of supercooling required for the model-independent approach to work has been extended from ϵ small up to values of ϵ of order 1, where ϵ is defined in Eq. (2.26).First, in Sec.3.1 it was pointed out that the supercool expansion proposed in [20] already gives an accurate model-independent description even for ϵ ∼ 1 if there are several degrees of freedom with dominant couplings to the flat-direction field χ (the one responsible for RSB).In Sec.3.2 it was then explained how to improve the supercool expansion to obtain a good model-independent description for ϵ ∼ 1 and an arbitrary number of (even few) degrees of freedom with dominant couplings to χ.This has been achieved by including, unlike in [20], the cubic term of the effective potential in (2.38) in the non-perturbative computation of the bounce action and treating the other corrections as perturbations (indeed, those are still small as long as ϵ is at most of order one, as we discussed in Sec.3.1).Such "improved supercool expansion" has been used to compute to good accuracy the nucleation temperature T n (and thus the strength α of the PT) as well as the inverse duration β of the PT in terms of few parameters that are fixed once the model is specified: χ 0 : the symmetry breaking scale β: the beta function of the quartic coupling λ χ of χ, evaluated at the scale where λ χ vanishes.
g: a sort of collective coupling of χ to all fields of the theory, which is precisely defined in Eq. (2.21).It is the square root of the sum of the squares of the couplings of χ to all fields.
g: an extra parameter that characterizes the size of the cubic term in (2.38).It is the cube root of the sum of the cubes of the couplings of χ to all bosonic fields, so g ≤ g.
Analytical calculation can be performed to a greater extent in the approach of Sec.3.1, but the improved supercool expansion of Sec.3.2 works for a larger class of models.
• In Sec. 4 such improved supercool expansion has then been applied to study in a modelindependent way the spectrum of GWs and the production of PBHs due to the first-order PT associated with the RSB.We have explained how to determine the GW spectrum (its amplitude and f peak ) in terms of the above mentioned parameters in the hypothesis of fast reheating after supercooling.Among other things, we have found values of f peak and regions of the parameter space that correspond to the GW background recently detected by pulsar timing arrays.Moreover, we have also found regions of the parameter space where the GW spectrum is above the sensitivity of LIGO-VIRGO O3 (which are then ruled out) and others that are within the reach of LISA.Furthermore, we have studied the generic validity of PBH production mechanisms in the RSB scenario for large supercooling.Also, we identified regions of the parameter space where PBHs produced by large over-densities due to an RSB PT can account for a significant fraction of the dark matter density.Other mechanisms for PBH production can be active, however, in specific models.
• In Sec. 5 we have applied the developed model-independent approach to study the PT in two RSB models: a simple illustrative one and a gauged B − L phenomenological completion of the SM featuring right-handed neutrinos below the EW scale.In both these models there are more than one degrees of freedom with dominant couplings to χ, but the number of these degrees of freedom is not much larger than one (it is six in the first model and three in the second one).Then we find that the improved supercool expansion of Sec.3.2 works better than the method of Sec.3.1, which however already allows us to obtain a reasonably-good semi-analytical estimate of the PT properties.Given the phenomenological interest of the B − L model, in the same section we have also studied reheating after supercooling, finding values of the parameters for which the reheating is fast and others for which it is not.

Figure 2 :
Figure 2: The solution λn of Eq. (3.13) as a function of a 1 and a 2 defined in (3.14).The inset in the right plot gives the maximal value of a 2 for a given a 1 such that the solution λn exists.Using the definitions of λ and λ in (3.8) and (2.28) one can extract the nucleation temperature.

Figure 3 :
Figure 3: Inverse duration β of the phase transition in units of the Hubble rate H n as a function of g and β for various values of the symmetry breaking scale χ 0 .Here g = g and ϵ < 3 has been imposed to guarantee the validity of the improved supercool expansion.

Figure 5 :
Figure 5: The peak frequency as a function of g and β in the case of fast reheating and fixing g * (T r ) = 110.Also, g = g and ϵ < 3 has been imposed.

Figure 8 :
Figure 8: Density plots giving the values of β/H n varying g and β.On the lower dashed line the whole dark matter is due to PBHs generated through the late-blooming mechanism (f PBH = 1); the upper dashed line corresponds instead to f PBH = 10 −10 .Here g = g and ϵ < 3 has been imposed.