Pulsar Timing Arrays and Primordial Black Holes from a Supercooled Phase Transition

An explicit realistic model featuring a supercooled phase transition, which allows us to explain the background of gravitational waves recently detected by pulsar timing arrays, is constructed. In this model the phase transition corresponds to radiative symmetry breaking (and mass generation) in a dark sector featuring a dark photon associated with the broken symmetry. The completion of the transition is ensured by a non-minimal coupling between gravity and the order parameter and fast reheating occurs thanks to a preheating phase. Finally, it is also shown that the model leads to primordial black hole production.


Introduction
Gravitational wave (GW) astronomy has become an extremely active and exciting field of physics after the discovery of GWs from binary black-hole and neutron-star mergers [1][2][3].More recently, the interest in this field has been further increased by the detection of a background of GWs by pulsar timing arrays (PTAs), including 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].
Among several possibilities to interpret the PTA background, an interesting option is first-order phase transitions (PTs).Indeed, these phenomena are not present in the Standard Model (SM) and, therefore, the detection of a GW background generated in this way (see Ref. [8] for a textbook introduction) would be clear evidence for new physics.A PT interpretation of the signals detected by PTAs performed by the NANOgrav collaboration and relevant for the present study was provided in Ref. [9] and other independent discussions of PT interpretations were given e.g. in Refs. .
First-order PTs generally occur when the symmetries are predominantly broken (and masses are then generated) through radiative corrections [33][34][35][36] (see Ref. [37] for a review).This radiative symmetry breaking (RSB) is interesting also because it has the unique property of allowing a perturbative (and thus calculable) description of PTs [36].Another distinctive property of such PTs is the presence of a period of supercooling when the temperature dropped much below the critical temperature [35,36].
The GW background detected by PTAs could be due to a supercooled PT associated with an RSB at a scale around the GeV [28,29].However, so far a concrete model with these characteristics has not been found and, therefore, its existence has not been guaranteed.Indeed, as pointed out in [20], identifying such a model is very challenging.
As a proof of concept, the purpose of this paper is to solve this issue by constructing an explicit realistic model that can explain the PTA background through a supercooled PT associated with RSB.
Besides generating a GW background, supercooled PT associated with an RSB also naturally lead to primordial black hole (PBH) production through the "late-blooming mechanism" (see e.g.Refs.[38][39][40][41][42] for a quantitative introduction): since the false vacuum decay of PTs is statistical for quantum and thermal reasons, in distinct patches the transition can occur at different times.Patches where the transition occurs the latest undergo the longest vacuum-dominated stage and, therefore, develop large over-densities, which collapse into PBHs.These PBHs can then account for a fraction f PBH of the dark matter density.
Another purpose of the present work is to identify regions of the parameter space of the above-mentioned model where the PBH production is significant and the PTA detection is explained at the same time.In the following sections we show that all this is possible and in the last section we highlight the main results achieved.

A minimal model
Let us start by considering a classically scale-invariant dark sector with RSB at the GeV scale.The simplest option is the Coleman-Weinberg construction [33], which features an Abelian gauge symmetry, U(1) D , and a charged scalar field ϕ.The corresponding no-scale Lagrangian is given by where the parameter e d is the U(1) D gauge coupling).The symmetry U(1) D is the one that undergoes RSB.This happens through a radiatively induced VEV of χ at the required energy scale, χ 0 ∼ GeV.As well known, RSB takes place because quantum corrections produce the following one-loop effective potential for χ where and μ is a value of the renormalization group scale µ at which λ χ = 0.In words, β is the β function of the quartic coupling at the energy scale where such coupling vanishes.However, the dark sector here interacts with the SM.Besides gravity, which is very weak, there are two possible operators that can mediate such interactions.One is the Higgs portal where H is the SM Higgs doublet and λ χh is the corresponding quartic coupling.The other one is a kinetic mixing between the SM photon, F µ , and A µ : where F µν ≡ ∂ µ F ν − ∂ ν F µ and η is the mixing parameter.
In general we can also introduce a non-minimal coupling ξ χ appearing in the Lagrangian density through the term where R is the Ricci scalar.The pure-gravity Lagrangian is assumed to be the standard Einstein-Hilbert one at the energies that are relevant for this work.Note that the term in (6) does not lead to an appreciable modification of gravity in our case because χ 0 ∼ GeV ≪ M P , unless one takes ξ χ extremely large, which we do not.However, as will be shown in Sec. 3, that term can ensure that the PT completes.The mixing in (5) can be removed by performing the field redefinition {A µ , F µ } → {A ′ µ , A µ } given by However, after this field redefinition, η still appears in the Lagrangian because with , and because of the interaction term eJ µ F µ between F µ and the SM electromagnetic current J µ , i.e.
The second term in the right-hand side of (8) leads to a background-dependent squared mass of the "dark photon" After RSB χ acquires a VEV χ 0 and this leads to the dark photon mass m d (χ 0 ).Since in our current setup χ 0 is much smaller than the EW scale and we demand the validity of perturbation theory (e d must not be large), m d (χ 0 ) should be much smaller than the EW scale too.
The second term in the right-hand side of ( 9) is instead an interaction between J µ and A ′ µ .Therefore, in order to satisfy the experimental limits η should be small (see [43,44] for reviews).So and the dark photon is coupled to the electromagnetic current through a millicharge eη.The Higgs portal coupling λ χh should also be very small because the mass m χ of δχ ≡ χ − χ 0 is also much smaller than the EW scale: m χ = β χ χ 0 ≪ χ 0 ∼ GeV.Note that the smallness of λ χh and η (and of gravity) ensures that, barring quantum corrections, the scale invariance of the dark sector is only mildly broken by its interactions.

The (improved) supercool expansion
As recalled in the introduction, in a PT associated with RSB supercooling always occur.The field χ is trapped in the false minimum, χ = 0, for a long time during which the spacetime expands exponentially with Hubble rate where M P is the reduced Planck mass.Eq. ( 11) can be obtained by using the trace of the Einstein equations and noting that the potential energy density difference between the false and the true vacuum is, according to Eq. ( 2), βχ 4 0 /16.After this exponential expansion, χ goes towards the true minimum, χ 0 , and reheating occurs up to a temperature T r at most given by where g * (T ) is the effective number of relativistic species at temperature T .The amount of supercooling that actually takes place depends on the model and can be parameterized by [36] where g plays the role of a "collective coupling" of χ with all fields of the theory: in the model of Sec. 2, The larger the supercooling, T ≪ χ 0 , the smaller ϵ.
As shown in [36] and [29], if ϵ ≪ 1, one can describe the PT, as well as the GW spectrum and the PBH production, with good accuracy by means of a small-ϵ expansion (called "supercool expansion").If T r is given by (12), at leading order (LO) one can parameterize the relevant physics in terms of χ 0 , g, β and g * only, while at next-to-leading order (NLO) one must introduce an extra parameter g, which in the model of Sec. 2 is g = g/ 6 √ 3.At LO the nucleation temperature T n (defined as the temperature at which the false-vacuum decay rate per unit volume Γ v equals H 4 I ), can be computed with the simple formula with If one wants to describe not only models with ϵ ≪ 1, but also those with ϵ of order 1, one can again describe the relevant physics perturbatively in the amount of supercooling, but through an improved version of the supercool expansion [29], where one has to exactly take g into account at any order [29].In this "improved supercool expansion" T n can be computed with a more complex procedure explained in [29].
Let us now quantify the error that one is making in analysing the model of Sec. 2 with the standard supercool expansion at NLO. Fig. 1  The inset shows the three sources of errors defined in (16).Varying χ 0 around the GeV scale does not show significant differences: a plot corresponding to another value of χ 0 is shown in the appendix.
max(e 1 , e 2 , e 3 ) and (e 1 , e 2 , e 3 ) defined by (e 1 , e 2 , e 3 ) ≡ c2 (computed for simplicity with the LO formula for T n in ( 15)) with c3 = 31.6915... .As discussed in [29], e 1 quantifies the above-mentioned error, while e 2 and e 3 are estimates of the errors present also in the improved supercool expansion at leading order.Fig. 1 tells us that such improved supercool expansion is valid and a better approximation than the standard supercool expansion as ϵ ∼ 1 and {e 1 > e 2 , e 1 > e 3 } in this model.Therefore, in the rest of the paper we will use the improved supercool expansion.
In general the existence of T n restricts the parameter space [29].One might wonder if the PT always complete, especially in the region of the parameter space where T n does not exist.If Γ v is not large enough to ensure the completion of the PT the temperature keeps decreasing and eventually the effect of the spacetime curvature, as well as quantum fluctuations, can become important in the decay rate [45][46][47][48].This can lead to the completion of the transition.Indeed, during the exponential expansion, the non-minimal coupling in (6), with ξ χ positive, gives a tachyonic mass contribution, −6ξ χ H 2 I χ 2 , to the effective potential of χ, which makes the transition complete when, eventually, the spacetime curvature dominates over thermal fluctuations.This contribution allows the field to reach the absolute minimum of the effective potential and ensures that percolation (i.e. when most of the vacuum energy is converted into radiation, see Sec. 4) occurs.We, therefore, assume that ξ χ is positive at least when this mechanism is needed to complete the PT.

Preheating and reheating
Can one has a sufficiently high reheating temperature T r in this model through direct decays of χ into SM particles?The answer to this question seems negative for the following reason.
The only tree-level channels for the decay of δχ into SM particles are χ → A ′ A ′ (HH) → 4 fermions, where the intermediate A ′ s and Higgs bosons H are virtual because an onshell A ′ or H is much heavier than δχ in our setup.The decays of A ′ and H into fermions are allowed by the interaction with the electromagnetic current in (10) and the Yukawa couplings, respectively.Dropping the external legs, the amplitudes for these processes for light-enough leptons are of order e 2 d χ 0 (eη) 2 /m 4 d (χ 0 ) and λ χh χ 0 y 2 SM /M 4 h , respectively, where the denominators m 4 d (χ 0 ) and M 4 h are due to the low energy limit of the two A ′ or H propagators and y SM is the SM Yukawa coupling of a light SM fermion.Taking the absolute value squared of the amplitude and integrating over the four-body phase space to obtain the decay rate Γ one finds or even smaller because of four-body phase-space suppressions.The factor m 7 χ has been inserted for dimensional reasons.As usual the scattering channels (such as χχ → A ′ A ′ (HH)) are not efficient in converting the energy density of χ into SM radiation density because the inverse scattering processes (A ′ A ′ (HH) → χχ) balance them and, therefore, the net average energy transfer vanishes.Moreover, setting η compatible with the experimental limits, χ 0 ≲ 10 GeV and g to realistic values compatible with the PTA GW signal, one finds that the channels in ( 17)-(18) can reheat the universe only up to a reheating temperature below the MeV scale, too small to realize a sufficiently fast reheating, where T r is rather given by the last equation in (12).
Let us then consider the preheating phase [49][50][51][52], when production of particles interacting with χ occurs as a result of the time dependence of this field through parametric resonance.In this model the only particle with a sizable coupling with χ is A ′ and its field equations are where det ĝ is the determinant of the spacetime metric, µ and we have used the smallness of η.After the false vacuum decay has taken place mostly of the energy stored in χ is thus transferred to A ′ because H 4 I is tiny compared to the potential density in the setup under study.
At some stage χ undergoes small oscillations around χ 0 .In such a regime one can analytically show that most of the energy density stored in χ is efficiently transferred to A ′ as we now do.In that case, the dependence of χ on the cosmic time t is given by with Φ ≪ χ 0 , and the background dependent mass m d (χ) ∼ e d χ 0 is generically extremely large compared to H I (see Eq. ( 11)).So we can neglect the expansion of the universe.
Recalling that η must be very small and setting the unitary gauge, where J ′ µ = 0, Eq. ( 19) then reads Taking the divergence of this equation and inserting the result back into Eq.( 21) gives In the regime (20) with Φ ≪ χ 0 the second term in the righthand side of the equation above is negligible with respect to the first one because suppressed by extra powers of χ 0 in the denominator and of m χ = βχ 0 ≪ χ 0 in the numerator: So the field equations of A ′ simply reduce to that of a scalar with a quartic portal interaction with χ 2 .Working in the Fourier space of momenta ⃗ k, all three polarizations A ′ k of the dark photons are then well described by the equation with k ≡ | ⃗ k|.This equation can now be brought into the wellknown Mathieu equation: where a k ≡ 4(k 2 + e 2 d χ 2 0 )/m 2 χ , q ≡ 4e 2 d χ 0 Φ/m 2 χ and z ≡ m χ t/2 + π/4.Therefore, as explained in [52,53] (see also Ref. [54] for a more recent analysis), there is an exponential growth of the dark photon field, A ′ k ∼ exp(µ (n) k z), which corresponds to an equally explosive growth of dark photon occupation numbers of quantum fluctuations when k belongs to some resonance bands labelled by an arbitrary integer n.These ranges of values of k always exist and the constants of the exponential growth, the µ (n)  k , allow to efficiently transfer all the energy density stored in the inflaton to the dark photon in a negligible cosmic time.This exponential growth cannot last forever because of the backreaction of A ′ on χ: when A ′ grows χ feels a growing effective mass, which makes χ tend to a vanishing value as time passes.However, the preheating phase allows to transfer the vacuum energy density ρ V to A ′ , which can then decay into pairs of SM fermions, thanks to the interaction with J µ in (10).The corresponding tree-level decay rate into lepton pairs is where m l is the lepton mass and the phase space factor Φ s is Unlike for the decays of δχ into four SM fermions, which we have previously discussed, setting η compatible with the experimental limits, χ 0 ≲ 10 GeV and g compatible with the PTA GW signal and all observational constraints, one finds that these processes allow for a very fast reheating.Including the hadronic decay channel of A ′ reinforces this result.Because ρ V represents the full energy budget of the system, the value of T r given by Eq. ( 12) is the correct value of the reheating temperature in this model.Therefore, given g, χ 0 and g * (T r ), the temperature T r can be computed through Eq. ( 12): where the expression of β in terms of g in Eq. ( 14) has been used.

Gravitational waves and primordial black holes
In the improved supercool expansion one can also obtain accurate expressions for β/H n , where β is the inverse duration of the PT and H n is the Hubble rate at T = T n , as well as the GW spectrum and its red-shifted frequency peak today f peak [29].Moreover, this expansion can be used to describe the PBH late-blooming mechanism in terms of the few parameters χ 0 , g, β and g [29].We now make use of these results in the context of the minimal model of Sec. 2.
Fig. 2 shows the predictions of this model for β/H n , f peak and T n setting χ 0 at the GeV scale.The plot shows that there are values of g such that this model reproduces the GW background detected by PTAs and is not excluded by other observations.As explained in [42], f PBH can be computed once the temperature in (28) and β/H n are known.In our model, these two quantities are functions of g (for given χ 0 and g * ).In particular β/H n decreases very rapidly when g approaches from above a critical value which is around 0.8.This is the manifestation of a general model-independent feature of phase transitions in the presence of RSB: when one approaches a critical line in the {g, β} plane (at fixed g and χ 0 ) β/H n decreases very rapidly (see Figs. 3 and 4 of Ref. [29]).On the other hand, in [42] it was shown that when β/H n decreases the fraction f PBH increases exponentially in this scenario.Therefore, f PBH is a function of g (for fixed χ 0 and g * ) that is very rapidly growing when g approaches the critical value from above.In particular the values of g corresponding to f PBH = 1 and f PBH = 10 −10 and the microlensing constraints [55][56][57][58] (the strongest constraints on PBHs in the setup of the upper plot in Fig. 2) are given in Fig. 2.These PBHs weight about a solar mass [28].
In the appendix we provide the corresponding plots for other values of χ 0 around the GeV scale.These have the same qualitative feature.Also, as expected, they show that when χ 0 goes above the GeV scale the values of g that can reproduce the PTA GW background tend to disappear.Moreover, note that using the values of g * (T r ) in the SM (see e.g.[59]), Eq. ( 28) leads to a reheating temperature of order 10 2 MeV for the realistic values of the parameters, which explain the background of GWs found by PTAs.
Furthermore, although PBHs are produced, one cannot account for the entire dark matter abundance through PBHs, as clear in the upper plot of Fig. 2  g should be of order one, which leads, using χ 0 ∼ GeV, to an almost GeV-scale dark photon, not to a sub-MeV one [43].One can complete the dark matter abundance as well as reproduce the observed neutrino oscillations and matter antimatter asymmetry by introducing, for example, three right-handed neutrinos along the lines of Refs.[60][61][62][63][64][65][66].Since neutrinos are not charged under U(1) D and do not appear in the electromagnetic current, the preheating and reheating analysis of Sec. 4 remains valid.

Conclusions
We have constructed and studied an explicit realistic model featuring RSB (and thus an associated supercooled PT) that explains the background of GWs found by PTAs.The model features a simple dark sector, which involves the field χ responsible for RSB and a hidden photon that can interact with the SM particles through a mixing with the SM photon.The completion of the transition is ensured by a non-minimal coupling between χ and R.Moreover, fast reheating after the supercooled PT occurs thanks to a preheating phase, where the vacuum energy density stored in χ is transferred to the dark photon, which can then decay to SM fermions via the mixing.Finally, we have also identified regions of the parameter space where one has significant PBH production.
The work presented here opens several research directions.As an outlook example, it would be interesting to generalize the dark sector: the model presented here is the minimal RSB one able to reproduce the GW background detected by PTAs, but one can introduce other scalar, vector and fermion fields with sizable couplings to χ.Moreover, it would be valuable to extend the preheating analysis reported here to a general PT associated with RSB in order to identify in a modelindependent way the parameter space corresponding to fast reheating.Finally, the interaction between the theoretical and observational community (PTAs, the LIGO and VIRGO collaborations for GW detection, etc.) will be vital to further probe this scenario.

Figure 2 :
Figure 2: Upper plot.β/H n , and (in the inset) f peak in the minimal model of Sec. 2 and computed with the improved supercool expansion.Here g * (T r ) = 6, although the plot depends weakly on g * .Values of g corresponding to the green solid line reproduce the NANOGrav signal at 68% C.L.; all displayed values of g correspond to ϵ at most of order 1, so the improved supercool expansion can be trusted.The shaded region is constrained by microlensing.Lower plot.T n as a function of χ 0 and g in the improved supercool expansion.
and in Fig. A.4.Also, the dark photon in this model cannot contribute significantly to dark matter: as evident in Fig. 2 (upper plot) and Fig. A.4,
shows ϵ as well as e nlo ≡ Figure 1: Comparison between ϵ and the error e nlo that one is making in using the standard supercool expansion at NLO in the minimal model of Sec. 2.