††institutetext: Harvard University, 17 Oxford Street, Cambridge, MA, 02139,
USA

# Consequences of a Stabilizing Field’s Self-Interactions for RS Cosmology

Rashmish K. Mishra and  Lisa Randall

###### Abstract

It has been argued that the Randall-Sundrum (RS) phase transition rate is
suppressed when the holographic theory corresponds to a large $N$ Yang-Mills
and when the stabilizing field has a small mass. Here we argue that self-
interactions can alleviate the latter suppression. We consider a cubic term in
the bulk potential for the Goldberger-Wise (GW) scalar that is responsible for
stabilizing the RS geometry. Adding a cubic term suffices to separate the two
roles of the GW stabilization: generating a large hierarchy and triggering
confinement. We study the resulting radion potential and the dynamics of the
early universe phase transition. For a negative coefficient of the cubic term,
the effect of the cubic becomes important in the infra-red, and the resulting
radion potential is deeper, thereby increasing the radion mass while
maintaining a large hierarchy. Staying within the radion effective field
theory, we calculate the rate of bubble nucleation from the hot phase to the
confined RS phase, both in thin and thick wall limits. The cubic term enhances
the rate and allows relaxing the condition on the maximum number of colors
$N_{\text{max}}$ of the dual theory for which the phase transition can be
completed. Importantly, this reduces the amount of supercooling that the false
vacuum undergoes, increases the peak frequency of the gravitational waves (GW)
produced from bubble collisions, and reduces the strength of the GW signal.
The reduced GW signal is however still within the reach of proposed space-
based GW detectors.

##  1 Introduction

The Randall-Sundrum (RS) framework based on warped extra-dimensional geometry
is an elegant way to understand various hierarchies in the Standard Model
(SM), and provides a rich framework for exploring several directions in both
formal and phenomenological research Agrawal:2022rqd . The early cosmological
history of these models is a confluence of many interesting phenomena. At low
temperatures, the RS phase is the thermodynamically stable phase. At high
temperatures, the RS phase is only metastable, while the stable phase is given
by a black-brane geometry, which in the dual field theory corresponds to the
deconfined phase. In minimal constructions, these two phases are separated by
a barrier Creminelli:2001th , and a transition between the two phases proceeds
by bubble nucleation.

Starting in the deconfined phase at high temperatures, the rate of transition
to the confined phase is very suppressed in the minimal models, and they
generically supercool past the critical temperature. The field configuration
that allows tunneling from one phase to another is a gravitational instanton,
the grammar for which is an active area of investigation. Early cosmological
history of RS models is therefore tied to the physics of confinement,
supercooling, and gravitational instantons. All these effects have direct
phenomenological relevance—the confinement/deconfinement phase transition and
its order is relevant for potential gravitational wave (GW) signals.
Supercooling is important for estimating the peak frequency and abundance for
the GW signal Randall:2006py ; Konstandin:2011dr , the present abundance of
relics from that era, as well as when and whether the phase transition
completes.111See ref. Levi:2022bzt  for a general discussion of supercooling
at both weak and strong coupling. The characteristic peak frequency and the
frequency dependence of the GW abundance have an imprint of the energy scales
involved and can be probed in present and proposed GW detectors LISA:2017pwj ;
Baker:2019nia ; Seto:2001qf ; Kawamura:2011zz ; Yagi:2011wg ; Isoyama:2018rjb
; Crowder:2005nr ; Corbin:2005ny ; Harry:2006fi . Finally, the exploration of
gravitational instantons is relevant for all of these since it affects the
rate of transition. These considerations make an exploration of the
theoretical and phenomenological aspects of phase transitions in RS framework
a well-motivated direction to pursue. If the phase transition is first-order,
the resulting GWs can provide access to the cosmological history of the
universe Kosowsky:1991ua ; Kosowsky:1992vn ; Kosowsky:1992rz ;
Kamionkowski:1993fg  and point towards yet to be discovered beyond the
Standard Model (BSM) physics.

In a realistic UV complete warped scenario, we generically expect significant
IR modifications, since in the dual picture, the theory is close to
confinement, and is very far away from a conformal field theory (CFT). With
this motivation, in the present work we allow for an IR modification in the RS
framework stabilized by a Goldberger-Wise (GW) scalar. We focus on how such a
modification affects various cosmological features.

More concretely, we consider a quadratic and a cubic term in the bulk
potential for the stabilizing GW field, with both their coefficients negative,
so that the GW profile grows in the IR, and the effect of the cubic becomes
important in the IR as well. Generically one expects even higher-order terms
in the potential, but for the present purposes, a cubic suffices to model the
IR modification. A non-zero cubic allows splitting the two roles played by the
GW mechanism—a logarithmic running which gives a large hierarchy, and the
triggering of the IR brane. As we will see, a non-zero cubic term can change
the shape of the radion potential while maintaining a large hierarchy. This
translates to a modification of the free energy and the bounce action for the
phase transition. In the present work, we will focus on the modifications on
the RS phase only and stay in the regime where the backreaction is not
important.

A suppressed rate of phase transition in the minimal models can be tracked
down to two parametric reasons: large $N_{c}$ and small $\delta$. Here $N_{c}$
is the number of colors in the dual theory, and $\delta$ characterizes CFT
breaking in the IR where the phase transition takes place (and is a function
of the parameters of the stabilization mechanism). For theoretical control, we
need $N_{c}\gg 1$. The minimal RS models, stabilized with a quadratic bulk GW
potential, also have $\delta\ll 1$. The present work can be understood as a
way to enhance the rate by increasing $\delta$, and is similar in spirit to
other such attempts in the literature. Addressing the $N_{c}$ suppression
directly will require a better modeling of the IR dynamics, and including
backreaction, which we will present in a future work.

Cosmological aspects of RS models have been studied both in the minimal
scenario Creminelli:2001th , and with interesting variations Nardini:2007me ;
Konstandin:2010cd ; Konstandin:2011dr ; vonHarling:2017yew ; Dillon:2017ctw ;
Bunk:2017fic ; Bruggisser:2018mrt ; Baratella:2018pxi ; Megias:2018sxv ;
Pomarol:2019aae ; Agashe:2019lhy ; Fujikura:2019oyi ; Megias:2020vek ;
Bigazzi:2020phm ; Agashe:2020lfz ; Agrawal:2021alq ; Csaki:2023pwy ;
Girmohanta:2023sjv ; Eroncel:2023uqf  that were aimed towards addressing the
issues in the minimal scenario. The authors in ref. Creminelli:2001th  pointed
out that without a stabilization mechanism, if the RS model is in the
deconfined phase at high temperatures, it never completes a transition to the
confined phase. With a stabilization mechanism from a GW scalar, the rate for
the phase transition is extremely suppressed, and is very constraining on what
constitute viable physical parameters. Ref. Konstandin:2010cd  pointed out
that by including back-reaction, the situation is less constraining than
originally claimed. In subsequent work, ref. vonHarling:2017yew  used
corrections to the radion potential coming from the QCD condensate, which
reduced the barrier depth, and hence enhanced the rate for the phase
transition. The authors in Baratella:2018pxi  further required the 4D CFT
degrees of freedom to flow to another fixed point after QCD confinement,
thereby gaining more theoretical control, and gave a geometrical picture for
this situation where the mass of the GW scalar is tachyonic Pomarol:2019aae .
The authors in Agashe:2019lhy ; Agashe:2020lfz  used similar ingredients where
they modeled the IR flow to another fixed point by a special form of bulk
potential for the GW field. In an orthogonal direction, the authors of ref.
Agrawal:2021alq  used finite temperature corrections to the combined potential
of radion and other light fields to address the issue. Ref. Csaki:2023pwy
used a relevant UV deformation to model the CFT breaking in the IR, which
resulted in a deeper radion potential. See also Girmohanta:2023sjv ;
Eroncel:2023uqf  for other interesting alternatives.

Compared to the previous approaches, the present work is different in the
following ways: the modified bulk potential we consider is generic, and is
expected in a realistic UV completion. We stay within the radion EFT parameter
space. We use a combination of analytical and numerical methods to handle the
computations. We consider both thin and thick wall limits, for completeness.
Our approach is a first step towards systematically including strongly coupled
IR effects. In some respects our results mimic those of Agashe:2019lhy ;
Agashe:2020lfz , who explicitly separated the two roles of the radion by
having one CFT to establish the hierarchy and the second to establish the
mass. Though perhaps under less control, our model is more generic and
reproduces the second feature organically. Once the GW field is sufficiently
big, playing a role in triggering the IR phase transition, it also contributes
to the radion mass and eases the phase transition completion.

The outline of the paper is the following. In sec. 2 we set up the notation
and obtain the GW scalar profile for both phases. In sec. 3 we derive the
radion potential, and in sec. 4 we derive the free energy of the two phases as
a function of temperature. In sec. 5 we calculate the rate of phase
transition, first in thin-wall approximation and then away from it, focusing
on the role of the IR effects. Sec. 6 contains the results, and a conclusion
follows in sec. 7. Technical details of the calculations are presented in app.
A, B, C.

##  2 5D action and scalar profiles

To set the stage we first fix some notation. We consider 5D spacetime with
locally constant negative cosmological constant. The general solution can be
parameterized as

$\displaystyle\differential s^{2}$
$\displaystyle=-e^{-2r}\left(1-e^{4(r-r_{h})}\right)\,\differential
t^{2}+e^{-2r}\differential\vec{x}^{2}+\frac{\differential
r^{2}}{1-e^{4(r-r_{h})}}\>,$ (1)

which is the metric for Anti-de Sitter Schwarzschild (AdSS) space, with the
asymptotic boundary of AdS space at $r\to-\infty$ and a horizon $r=r_{h}$,
extended in the transverse directions (a black brane). We have set the AdS
scale $\ell_{\text{AdS}}=1$ here, and in what follows. In the limit of
$r_{h}\to\infty$, this reduces to the Poincare patch of AdS space.

A UV brane with appropriate (positive) tension, whose location can be chosen
to be at $r_{\text{uv}}=0$, cuts off the UV region of the spacetime and
ensures a normalizable 4d graviton. For the $r_{h}=\infty$ case, the spacetime
will be truncated at the IR brane located at $r=r_{\text{ir}}$, again with an
appropriately chosen (negative) tension. We therefore have two spacetime
metrics to consider, which we refer to as the BB (Black Brane) and RS (Randall
Sundrum) spacetimes respectively:

RS $\displaystyle:\>\>\differential s^{2}=-e^{-2r}\,\differential
t^{2}+e^{-2r}\differential\vec{x}^{2}+\differential
r^{2}\>,\qquad\qquad\qquad\qquad\qquad\>\>0\leq r\leq r_{\text{ir}}\>,$ BB
$\displaystyle:\>\>\differential
s^{2}=-e^{-2r}\left(1-e^{4(r-r_{h})}\right)\,\differential
t^{2}+e^{-2r}\differential\vec{x}^{2}+\frac{\differential
r^{2}}{1-e^{4(r-r_{h})}}\>,\>\>0\leq r\leq r_{h}\>.$ (2)

The RS and BB spacetimes are dual to the confined and the deconfined phases
respectively in the field theory. The Hawking radiation from the horizon in
the BB geometry gives a temperature to the black hole, which is a function of
the location of the horizon. The effect of a finite temperature $T$ can be
studied in the Euclidean version of the spacetime, with the Euclidean time
$t_{E}$ identified with a period $\beta=1/T$. For $r_{h}\neq 0$, the Euclidean
continuation of eq. (1) with a periodic $t_{E}$ is smooth only when

$\displaystyle\beta=\pi\exp(r_{h})\>,$ (3)

whereas any $\beta$ is fine for $r_{h}=\infty$ (RS background). Since we need
to consider the dynamics at finite temperature, we will work in Euclidean
compactified time.

The location of the IR brane is a modulus in the RS geometry, and needs a
stabilization mechanism to have a fixed value. Such a stabilization can be
provided by a 5D GW scalar $\chi$ with appropriately chosen parameters such
that it gets a profile along the extra dimension, and generates a potential
for the field $r_{\text{ir}}(x)$. In this work we consider more general
potentials and boundary conditions than the original proposal
Goldberger:1999uk ; Goldberger:1999un , and argue that these modifications are
well-motivated and geared towards modeling the IR dynamics appropriately.
Specifically, we consider a 5D scalar $\chi$, with the action

$\displaystyle S_{\chi}$ $\displaystyle=\int
d^{5}x\sqrt{g}\left(-\frac{1}{2}\left(\partial\,\chi\right)^{2}-V_{B}(\chi)\right)-\sum_{i}\int
d^{4}x\,\sqrt{g_{i}}\,V_{i}(\chi)\>,$ (4)

with a bulk potential $V_{B}(\chi)$ and boundary potential(s) $V_{i}(\chi)$
which set the boundary conditions. Note that $i$ takes the value $(uv,ir)$ for
the RS case, and only $(uv)$ for the BB case. We would like to solve for the
profile of $\chi$ in the two backgrounds, and in the limit of small back-
reaction.

The choice of the bulk and brane localized potentials is governed by the
dynamics one wants to model. In the dual interpretation, the bulk GW potential
can be mapped to the renormalization flow of a deformation. For the effect to
become more important in the IR, the deformation should be relevant, and
higher-order terms in the beta function should become important as the
deformation grows. A constant mass for the GW scalar corresponds to a constant
beta function. Higher order terms in the beta function correspond to higher
order terms in the GW potential Rattazzi:2000hs ; Chacko:2013dra . We add a
cubic term in the bulk GW potential, in addition to the mass term, to model
the higher-order terms in the beta function. With this motivation, we choose
the bulk potential to be

$\displaystyle V_{B}(\chi)$
$\displaystyle=2\epsilon_{2}\chi^{2}+\frac{4}{3}\epsilon_{3}\chi^{3}\>,$ (5)

with $\epsilon_{2}<0,\epsilon_{3}<0$. The first sign, $\epsilon_{2}<0$,
corresponds to a relevant deformation of the dual theory, and ensures a
logarithmic running for small $|\epsilon_{2}|$. The second sign,
$\epsilon_{3}<0$, ensures the deformation gets larger in the IR. We have kept
to a cubic term in the potential, which is sufficient to capture the effect of
strong coupling as we show later, although there can in general be quartic and
higher-order terms as well.222Ref. Agashe:2020lfz  used a GW potential with up
to quartic terms to model a specific flow in the dual theory, from a UV fixed
point to an IR fixed point. The dynamics we want to model in the present work
is different.

The mass term and the self-interaction term in the GW potential model
different aspects of the dynamics. A small mass allows a large running,
whereas the self-interaction term, which for small $\epsilon_{3}$ is important
only in the IR, models the existence of a more complicated radion potential in
the IR after confinement. It is important to separate these two effects, and
study the resulting effect on radion potential and the phase transition, which
is the motivation for this work.

For simplicity we choose a UV brane potential that fixes the value of $\chi$
at the UV brane. The IR brane potential is chosen to allow the scalar to
adjust its value, and again for simplicity we fix its derivative. These
features can be modeled by the potentials

$\displaystyle V_{\text{uv}}(\chi)$ $\displaystyle=\beta_{\text{uv}}(\chi-
v_{\text{uv}})^{2}\>,\;\>\beta_{\text{uv}}\to\infty\>,$ $\displaystyle
V_{\text{ir}}(\chi)$ $\displaystyle=2\alpha_{\text{ir}}\chi\>.$ (6)

Note that these are simplifications, and a more complete analysis should allow
for mixed boundary conditions at both the UV and the IR. Also note that the IR
boundary condition is not relevant for the BB case—the GW profile is required
to have a “regular” behavior at the horizon. For a truly holographic
interpretation, we would not necessarily have a boundary condition imposed
directly in the IR but we follow the original analysis and do this for
simplicity.

Given the bulk and boundary potentials, the GW field develops a profile
$\chi(r)$ along the radial direction, which is different for the RS and BB
backgrounds because the IR boundary condition is different for the two cases.
While an exact solution can be obtained numerically, it is possible to get
approximate analytical solutions by dividing the bulk into different parts
where different terms in the equation dominate. Leaving the details to
appendix A, the leading order solutions are given as

$\displaystyle\chi_{\text{RS}}(r)$
$\displaystyle=-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}+\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}\left(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}}\right)}\>,\>\>0\leq
r\leq r_{\text{ir}}\>,$ (7) $\displaystyle\chi_{\text{BB}}(r)$
$\displaystyle=\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}\left(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}}\right)}\>,\qquad\qquad\qquad\>\>\>\>\>\>0\leq
r\leq r_{h}\>.$ (8)

The solutions are obtained in the limit of $|\epsilon_{2}|\ll
1,v_{\text{uv}}\ll 1,r_{\text{ir}}\gg 1,|\epsilon_{2}|r_{\text{ir}}\lesssim
1,|\epsilon_{3}|r_{\text{ir}}\lesssim 1$ ($r_{\text{ir}}\to r_{h}$ in the
corresponding conditions for the BB solution). In the case of
$\epsilon_{2}<0$, $\chi$ becomes singular at
$r_{s}=-(1/\epsilon_{2})\log(1+\epsilon_{2}/v_{\text{uv}}\epsilon_{3})$ and is
an artifact of the analytical solution failing to be a good approximation. By
appropriately choosing parameters $r_{s}$ can be made large and as long as
$r_{s}\gg r_{\text{ir}},r_{h}$, we can trust the solutions.

A few comments are in order about these solutions. First, for the RS solution,
the term proportional to $\alpha_{\text{ir}}$ can be important only very close
to $r=r_{\text{ir}}$. Second, the requirement of regularity of the solution in
the BB background is the same as dropping the term proportional to
$\alpha_{\text{ir}}$ in the RS solution, at least to leading order in the
approximation. Finally, these approximate solutions clarify the effect of
$\epsilon_{3}$. For small $r$, the term proportional to $v_{\text{uv}}$ can be
written as $v_{\text{uv}}\exp(-(\epsilon_{2}+v_{\text{uv}}\epsilon_{3})r)$.
The self-interaction term therefore effectively acts as an additive effect to
$\epsilon_{2}$. As $r$ increases, this effect compounds.

A good diagnostic to quantify the effect of a non-zero $\epsilon_{3}$ is to
look at the ratio of the mass of the radion to the IR scale. The model
requires a large hierarchy, which is the result of logarithmic running. A non-
zero $\epsilon_{3}$ allows higher order terms in the radion potential to
balance in the IR. We will have more to say about this in the next section.

Figure 1 shows a comparison between the approximate solutions in eqs. (7), (8)
with numerical solutions, for some choices of parameters. The procedure to
numerically obtain regular solutions in the BB background are discussed in
appendix A.

Figure 1:  GW scalar profile with and without self-interactions, for the RS
background (left) and the BB background (right), for other parameters fixed.
Solid lines (numerical) show good agreement with the approximate solutions (in
dashed). A negative (positive) $\epsilon_{3}$ leads to more (less) growth of
the profile in the IR.

##  3 Radion potential

Given the profile for $\chi(r)$, we can calculate the resultant radion
potential. The basic idea is that as one varies $r_{\text{ir}}$, the energy
contained in the potential for $\chi(r)$ changes and for a choice of
parameters, a suitable minimum can be obtained. The purely gravitational part
of the action gives a kinetic term for $r_{\text{ir}}$. Evaluating the $\chi$
action on the approximate solution for $\chi(r)$ in eq. (7), one obtains the
potential for $r_{\text{ir}}$. In terms of the field
$\varphi=e^{-r_{\text{ir}}}$, the 4D action is given as

$\displaystyle S$ $\displaystyle=\int
d^{4}x\sqrt{g}\left(-12M_{5}^{3}\left(\partial\varphi\right)^{2}-V(\varphi)\right)\>,$
$\displaystyle V(\varphi)$
$\displaystyle=24M_{5}^{3}\kappa^{4}\,\varphi^{4}\left(1+\frac{a_{2}}{24M_{5}^{3}\kappa^{4}}\frac{\lambda\varphi^{\epsilon_{2}}}{1-\lambda\varphi^{\epsilon_{2}}}-\frac{a_{3}}{24M_{5}^{3}\kappa^{4}}\,\log(1-\lambda\varphi^{\epsilon_{2}})\right),\>$
$\displaystyle\lambda=\frac{v_{\text{uv}}\epsilon_{32}}{1+v_{\text{uv}}\epsilon_{32}},\>\epsilon_{32}=\frac{\epsilon_{3}}{\epsilon_{2}},\>a_{2}=-\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}-\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}+2\alpha_{\text{ir}},\>a_{3}=\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\>.$
(9)

The details of the computation for $V(\varphi)$ are given in app. B. Note that
$\varphi$ is not canonically normalized in our notation. We have pulled an
overall factor of $M_{5}^{3}$ outside from the potential, and the parameter
$\kappa\lesssim 1$ for small back-reaction.333We are working in the glueball
normalization where the quartic scales as $N_{c}^{2}\sim M_{5}^{3}$. In the
limit of $\lambda\varphi^{\epsilon_{2}}\ll 1$ (note that for this we need
$\lambda\sim\mathcal{O}(v_{\text{uv}})\ll 1$, since $\varphi^{\epsilon_{2}}$
can be $\mathcal{O}(1)$, for $\epsilon_{2}<0$), the potential can be expanded
in a power series:

$\displaystyle V(\varphi)$
$\displaystyle=\varphi^{4}\left(b_{0}+b_{1}\lambda\varphi^{\epsilon_{2}}+b_{2}\lambda^{2}\varphi^{2\epsilon_{2}}+b_{3}\lambda^{3}\varphi^{3\epsilon_{2}}+\cdots\right)\>.$
(10)

The coefficients $b_{i}$ are readily calculable given the explicit form of the
potential, and are functions of
$\alpha_{\text{ir}},v_{\text{uv}},\epsilon_{2}$ and $\epsilon_{3}$. In the
limit of $\epsilon_{3}\to 0$, only $b_{0}$ and $b_{1}$ are non-zero, and the
potential simplifies to the familiar racetrack form. Note that in this limit,
$\lambda\to 0$ but it is balanced by the $\epsilon_{3}$ in the denominator of
$a_{2}$ and $a_{3}$ in eq. (9).

$\displaystyle V(\varphi)$ $\displaystyle\underset{\epsilon_{3}\to
0}{=}\>\>24M_{5}^{3}\kappa^{4}\,\varphi^{4}\left(1-\frac{1}{48M_{5}^{3}\kappa^{4}}\alpha_{\text{ir}}v_{\text{uv}}\varphi^{\epsilon_{2}}\right)$
$\displaystyle=24M_{5}^{3}\kappa^{4}\,\varphi^{4}\left(1-\frac{1}{1+\epsilon/4}\left(\frac{\varphi}{\varphi_{\text{min}}}\right)^{\epsilon_{2}}\right)\>.$
(11)

The $b_{2}$ and higher order terms in eq. (10) mimic the effects of strong
dynamics—as the running coupling grows to be big, the higher order terms
become important Rattazzi:2000hs . In fact, we will see that for a relevant
choice of parameters, we have to include several terms in the expansion. The
effect of higher order terms is that for a small $\varphi_{\text{min}}$ (large
hierarchy), various combinations of terms can balance each other and give a
minimum in the radion potential in which case the second derivative of the
potential at the minimum can be enhanced. This is the usual expectation that
for a strong breaking in the IR, the radion is not parametrically light
anymore. To understand this enhancement, let’s first consider the generic form
of the radion potential

$\displaystyle V(\varphi)$
$\displaystyle=b_{0}\,\varphi^{4}\,P(\varphi^{\epsilon_{2}})\>,$ (12)

where $P(x)$ is a polynomial in $x$ of some given order, with the first term
$1$ (since we have factored out an overall constant in eq. (12)). Expanding
eq. (9), some of the terms in $P(x)$ are explicitly given as

$\displaystyle b_{0}\,P(x)$ $\displaystyle=b_{0}+\sum_{i\geq
1}b_{i}x^{i}\>,\qquad b_{0}=24M_{5}^{3}\kappa^{4}\>,$ $\displaystyle b_{1}$
$\displaystyle=-v_{\text{uv}}\left(\frac{1}{2}\alpha_{\text{ir}}-2\alpha_{\text{ir}}\frac{\epsilon_{3}}{\epsilon_{2}}+\frac{1}{32}\alpha_{\text{ir}}^{2}\epsilon_{3}\right)\left(1+\frac{v_{\text{uv}}\epsilon_{3}}{\epsilon_{2}}\right)^{-1}\>,$
$\displaystyle b_{2}$
$\displaystyle=-v_{\text{uv}}\left(\frac{v_{\text{uv}}\epsilon_{3}}{\epsilon_{2}}\right)\left(\frac{3}{4}\alpha_{\text{ir}}-2\alpha_{\text{ir}}\frac{\epsilon_{3}}{\epsilon_{2}}+\frac{1}{32}\alpha_{\text{ir}}^{2}\epsilon_{3}\right)\left(1+\frac{v_{\text{uv}}\epsilon_{3}}{\epsilon_{2}}\right)^{-2}\>,$
$\displaystyle b_{3}$
$\displaystyle=-v_{\text{uv}}\left(\frac{v_{\text{uv}}\epsilon_{3}}{\epsilon_{2}}\right)^{2}\left(\frac{3}{4}\alpha_{\text{ir}}-2\alpha_{\text{ir}}\frac{\epsilon_{3}}{\epsilon_{2}}+\frac{1}{32}\alpha_{\text{ir}}^{2}\epsilon_{3}\right)\left(1+\frac{v_{\text{uv}}\epsilon_{3}}{\epsilon_{2}}\right)^{-2}\>,$
$\displaystyle\vdots\>\>\>$ (13)

For $v_{\text{uv}}\lesssim 1$ and $\epsilon_{3}/\epsilon_{2}\lesssim 1$,
higher order terms are successively smaller. However, near the minimum,
$x=\varphi^{\epsilon_{2}}\gtrsim 1$ so that higher powers of $x$ are bigger.
Therefore, near the minimum, terms in $P(x)$ are products of successively
increasing and successively decreasing factors. For certain choice of
parameters, a combination of terms balance each other. Since $x\gtrsim 1$, we
still get a large hierarchy. This discussion also makes it clear that for a
small $\epsilon_{3}$, the higher order terms can contribute without changing
the hierarchy too much when $\alpha_{\text{ir}}$ and $v_{\text{uv}}$ take
larger values. At such values, $P$ is enhanced and so is the second derivative
of the potential.

In summary, a non-zero $\epsilon_{3}$, in conjunction with other parameters of
the GW sector, can increase the radion mass while generating a similar
hierarchy. To illustrate this, we will work with four benchmark parameters
A,B,C,D in table 1, for which $\varphi_{\text{min}}\sim 10^{-16}$.444In table
1, $\varphi_{\text{min}}$ is calculated for
$M_{5}^{3}=N_{c}^{2}/16\pi^{2},N_{c}=1$. For a different $N_{c}$,
$\alpha_{\text{ir}}$ and $v_{\text{uv}}$ have to be adjusted to keep
$\varphi_{\text{min}}$ fixed. The mass of the physical radion and $T_{c}$
(defined in eq. (20)) do not change with $N_{c}$. These parameters are chosen
with certain self-consistency conditions in mind. Requiring to stay in the
radion EFT, we need to ensure that the radion mass is at most or slightly
smaller than the Kaluza-Klein (KK) scale which sets the mass of other KK
modes. This also ensures that the back-reaction on the geometry from the GW
scalar can be ignored. Another requirement is to have
$T_{c}/\varphi_{\text{min}}\lesssim 1$ so that temperature corrections to the
potential can be ignored, at least in the vicinity of the minimum. $T_{c}$ is
set by the value of the potential at the minimum $V(\varphi_{\text{min}})$ (as
discussed in the next section).

| $\kappa$ | $\epsilon_{2}$ | $\epsilon_{3}$ | $\alpha_{\text{ir}}$ | $v_{\text{uv}}$ | $\varphi_{\text{min}}\times 10^{16}$ | $V^{\prime\prime}(\varphi_{\text{min}})/\varphi_{\text{min}}^{2}$ | $-V(\varphi_{\text{min}})/\varphi_{\text{min}}^{4}$
---|---|---|---|---|---|---|---|---
A | $10^{-1/4}$ | -1/25 | 0 | 1/10 | 1/14 | 1.47 | 0.002 | $10^{-4}$
B | $10^{-1/4}$ | -1/25 | -1/100 | 5/2 | 1/14 | 1.09 | 0.005 | $3\times 10^{-4}$
C | $10^{-1/4}$ | -1/25 | -1/90 | 5/2 | 1/5 | 0.86 | 0.032 | $2\times 10^{-3}$
D | $10^{-1/4}$ | -1/25 | -1/81 | 5/2 | 1/3 | 0.59 | 0.135 | $8\times 10^{-3}$
Table 1: Benchmark choice of parameters to show the effect of self-
interaction.

Fig. 2(a) shows the radion potential for these parameters. For B,C,D, with
$\epsilon_{3}\neq 0$, a deeper potential at the minimum can be clearly seen.
Fig. 2(b) shows the value of various terms in the series expansion of the
derivative of the potential near the minimum. It is clear that for A, the
first and the second terms balance, for B, the first and the third terms
balance, for C, the second and the third terms balance, and for D, the second
term balances the third and fourth terms together. For C,D, $v_{\text{uv}}$
was increased to make sure the higher order terms are enhanced. Further, the
magnitude of the dominant term is largest in D, which is correlated with the
largest value of the second derivative at the minimum. One can estimate the
second derivative at the minimum once we know which terms balance which.
Starting with the $m^{\text{th}}$ term in the potential,
$V_{m}=b_{m}\,\varphi^{4+m\epsilon_{2}}$, the derivative and the second
derivative are

$\displaystyle
V^{\prime}_{m}=b_{m}\,(4+m\epsilon_{2})\,\varphi^{3+m\epsilon_{2}}\>,\qquad
V^{\prime\prime}_{m}=b_{m}\,(4+m\epsilon_{2})(3+m\epsilon_{2})\,\varphi^{2+m\epsilon_{2}}\>.$
(14)

The $m^{\text{th}}$ and $n^{\text{th}}$ term in $V^{\prime}$ can balance near
the minimum if

$\displaystyle
b_{m}\,(4+m\epsilon_{2})\,\varphi_{\text{min}}^{3+m\epsilon_{2}}\sim-
b_{n}\,(4+n\epsilon_{2})\,\varphi_{\text{min}}^{3+n\epsilon_{2}}\>.$ (15)

Considering these two terms, and using the above, the second derivative at the
minimum is

$\displaystyle V^{\prime\prime}/\varphi_{\text{min}}^{2}\sim
b_{m}\,(4+m\epsilon_{2})(m-n)\epsilon_{2}\,\varphi^{m\epsilon_{2}}\sim
b_{n}\,(4+n\epsilon_{2})(n-m)\epsilon_{2}\,\varphi^{n\epsilon_{2}}\>.$ (16)

Using $b_{0}=24M_{5}^{3}\kappa^{4}$, $b_{m>0}\sim
v_{\text{uv}}\left(v_{\text{uv}}\epsilon_{3}/\epsilon_{2}\right)^{m-1}\alpha_{\text{ir}}$
and $\varphi^{\epsilon_{2}}\sim(10^{-16})^{-1/25}\sim 4.3$ near the minimum,
the above estimate for $V^{\prime\prime}/\varphi_{\text{min}}^{2}$ matches the
numbers in table 1 obtained numerically. For parameters C,D, the dominant term
at the minimum scales as $1/\epsilon_{2}$ (i.e. $m=2$), and this
parametrically cancels the $\epsilon_{2}$ factor in the numerator in eq. (16).

(a) Radion potentials near the minimum, for the parameter choices in table 1.
A deeper potential results from a non-zero $\epsilon_{3}$.

(b) Various terms in the series expansion of the derivative of the radion
potential, for the parameter choices in table 1. The terms balancing near the
minimum are, top left: first and second, top right: first and third, bottom
left: second and third, and bottom right: second, third and fourth.

Figure 2: Radion potentials, for the parameter choices in table 1.

A word of caution for the obtained radion potential: for
$\epsilon_{2}<0,\epsilon_{3}<0$, which is the case at hand, the potential has
a singularity at $\varphi_{s}=(1/\lambda)^{1/\epsilon_{2}}$. The singularity
is coming from the analytical solution for $\chi(r)$ breaking down, and since
we used this analytical solution to obtain the radion potential we see it in
the potential too. In a more complete but numerically tedious calculation,
this singularity would not appear. The potential in (9) therefore cannot be
trusted for large $r_{\text{ir}}$, or equivalently for $\varphi\to 0$. This is
to be kept in mind when we discuss the thick wall results in sec. 5 (and the
computational details in app. C), which probe small values of $\varphi$ as the
temperature becomes small. At a finite temperature $T$, the potential is
expected to be modified for $\varphi\lesssim T$. As we will explain later, we
use the potential only in the $\varphi\gtrsim T$ region. The expression for
the radion potential is therefore useful as long as $\varphi_{s}\ll T$. This
also means we cannot extend the analysis to arbitrarily low temperatures. As
we will show, before running into this issue a significant reduction in bounce
action can be achieved.

##  4 Free energy of the phases

For both RS and BB phases, the free energy gets contributions from the
gravitational and the GW sectors. Since we are not including back-reaction,
the gravitational calculation is the same as reported in the literature
Creminelli:2001th . The free energy is UV sensitive, so it is more sensible to
talk about the difference in the free energies between the phases. At the
minimum, the difference in free energies coming purely from the gravitational
part of the action is given by:

$\displaystyle\left(F_{\text{GR, BB}}-F_{\text{GR, RS}}\right)_{\text{min}}$
$\displaystyle=-2\pi^{4}M_{5}^{3}T^{4}+C=-\frac{\pi^{2}}{8}\,N_{c}^{2}\,T^{4}+C\>,$
(17)

where $C$ is a finite constant to be determined later, and in the last
equality we have used the relation between $M_{5}$ and the number of colors
$N_{c}$ of the dual theory: $16\pi^{2}M_{5}^{3}=N_{c}^{2}$. As discussed in
ref. Creminelli:2001th , in this computation the location of the horizon
$r_{h}$, parameterized by a temperature $T_{h}=\pi\exp(r_{h})$, is made a
dynamical variable. For generic values of $T_{h}$ (in the Euclidean picture),
there is a conical singularity at the horizon. Regulating the singularity one
gets a contribution to the free energy for $T_{h}\neq T$. At the minimum,
$T_{h}=T$, the conical singularity disappears, and $T$ is related to $r_{h}$
by eq. (3). In deriving eq. (17), $\beta$ in the RS phase is adjusted to match
the geometry at the UV cutoff (see refs. Witten:1998zw ; Creminelli:2001th
for details).

We can effectively define eq. (17) to be the free energy of the BB phase at
the minimum. Taking $C\equiv 2\pi^{4}M_{5}^{3}T_{c}^{4}$, the free energy of
the BB phase at the minimum can be written as

$\displaystyle F_{\text{BB,
min}}=2\pi^{4}M_{5}^{3}\left(T_{c}^{4}-T^{4}\right)\>.$ (18)

At the moment, $T_{c}$ is a parameter and the free energy of the BB phase is
positive or negative depending on whether $T<T_{c}$ or $T>T_{c}$.

In the presence of a GW scalar, the free energy of the BB phase gets an
additional contribution. This can be computed by first solving for the scalar
profile in the BB background, which is then used to calculate the free energy.
The scalar contribution is subleading compared to the purely gravitational
contribution, and we will not include it here. This is self-consistent with
not including back-reaction from the scalar on the BB geometry.

Similarly, the free energy of the RS phase gets an additional contribution in
the presence of a GW scalar. Staying within the radion EFT, and at low enough
temperatures, the free energy is given by the radion potential itself. The
radion potential is normalized so that it vanishes at the minimum. This
amounts to tuning the cosmological constant in the 4D EFT to zero. We
therefore have

$\displaystyle F_{\text{GW,
RS}}=V(\varphi)-V(\varphi_{\text{min}})\>,\>\>F_{\text{GW, RS, min}}=0\>.$
(19)

We further require that the free energies of the two phases are equal in the
$\varphi\to 0$ and $T\to 0$ limits, which gives

$\displaystyle-V(\varphi_{\text{min}})=2\pi^{4}M_{5}^{3}T_{c}^{4}\>.$ (20)

This fixes $T_{c}$ in terms of the radion potential. Putting everything
together, we have

$\displaystyle F_{\text{BB, min}}$
$\displaystyle=2\pi^{4}M_{5}^{3}\left(T_{c}^{4}-T^{4}\right)\>,$ (21)
$\displaystyle F_{\text{RS}}$
$\displaystyle=V(\varphi)-V(\varphi_{\text{min}})\>,$ (22) $\displaystyle
F_{\text{RS, min}}$ $\displaystyle=0\>,$ (23) $\displaystyle F_{\text{BB,
min}}-F_{\text{RS, min}}$
$\displaystyle=2\pi^{4}M_{5}^{3}\left(T_{c}^{4}-T^{4}\right)\>.$ (24)

The meaning of $T_{c}$ is made clear by eq. (24): At the critical temperature
$T=T_{c}$, the difference in the free energy of the two phases is zero. For
$T>T_{c}$, the BB phase has a lower free energy than the RS phase and is the
preferred phase thermodynamically. As the temperature drops below $T_{c}$, the
RS phase has a lower free energy and becomes the thermodynamically preferred
phase.

##  5 Dynamics of the phase transition

We assume that after the end of inflation and reheating, the RS model is in
the BB phase, which is the thermodynamically stable phase at very high
temperatures. As the temperature drops, the RS phase becomes thermodynamically
favorable. Since both these phases are local minima of free energy, they are
separated by some kind of barrier in the field space, and the phase transition
is first order, proceeding by bubbles of true vacuum nucleating inside the
false vacuum. The standard prescription to calculate the rate of such a phase
transition Coleman:1977py ; Callan:1977pt  can be applied, but with important
modifications pointed out in ref. Creminelli:2001th , which we mention next.

To calculate the tunneling rate across a barrier, one has to identify the
relevant fields that interpolate between the true and the false vacuum. In the
present case, staying within the radion EFT, the relevant field for the
tunneling rate calculation on the RS side is the radion itself. The radion is
a composite, and a different field has to be identified in the BB phase. Ref
Creminelli:2001th  assumed the relevant degree of freedom in the BB phase to
be $T_{h}$ (temperature of the horizon), for which a potential can be written,
but the kinetic term is not known (see however ref. Bigazzi:2020phm  for a
discussion of this aspect). The way out was provided by noticing that the path
in the field space that interpolates between the false and the true vacua has
some features in the phenomenological cases of interest. One can split it into
three regions: $i)$ the BB region, $ii)$ the $\varphi\lesssim T$ region, and
$iii)$ the $T\lesssim\varphi\lesssim\varphi_{\text{min}}$ region. For
$T\lesssim T_{c}$ (i.e. in the thin wall limit), the calculable contribution
from region $iii)$ is parametrically larger than the combined contribution
from $i)$ and $ii)$, so that it provides a useful estimate for the full
action. As $T/T_{c}$ reduces, the contribution from region $iii)$ also
reduces. This means that if at a given $T/T_{c}$, the contribution from region
$iii)$ is too small, one cannot ignore the other two contributions.
Generically, we can take the other regions to contribute $\mathcal{O}(1)$
amount to the action, so that the actually computable action, from region
$iii)$, can only be trusted for values of $T/T_{c}$ such that it is at least
$\mathcal{O}(1)$. All in all, this means the above approach cannot be extended
to very small values of $T/T_{c}$. Also, note that the radion potential has a
singularity at
$\varphi_{s}=(1/\lambda)^{1/\epsilon_{2}}\ll\varphi_{\text{min}}$ coming from
a breakdown of the approximations that were used to derive it. Since
$\varphi\gtrsim T$ in region $iii)$, we do not need to worry about this as
long as $T$ is not too small. With these considerations in mind, we now
calculate the bounce action in both the thin and thick wall limits

We leave the technical details of the bounce action calculation for thin and
thick wall cases to app. C, and give the final expressions here. For the thin
wall case, the bounce action is

$\displaystyle S_{b}$
$\displaystyle=\frac{4}{3\pi^{7}M_{5}^{6}}\left(\frac{S_{1}}{T_{c}^{3}}\right)^{3}\frac{T_{c}/T}{\left(1-T^{4}/T_{c}^{4}\right)^{2}}\>,$
$\displaystyle S_{1}$
$\displaystyle=\sqrt{48M_{5}^{3}}\int_{T}^{\varphi_{\text{min}}}\differential\varphi\sqrt{V(\varphi)-V(\varphi_{\text{min}})}\>.$
(25)

For the thick wall case, the approach is numerical. We have to minimize the
action

$\displaystyle S_{b}$
$\displaystyle=\int\differential^{4}x\left(12M_{5}^{2}(\partial\varphi)^{2}+V(\varphi)+2\pi^{4}M_{5}^{3}T^{4}\right)\>,$
(26)

subject to the boundary conditions $\varphi^{\prime}(0)=0$ and
$\varphi^{\prime}(\varphi\approx 0)=-(\pi^{2}/\sqrt{6})T^{2}$. The second
condition comes from equating the energy across the bubble boundary
Agashe:2020lfz . Working with rescaled quantities

$\displaystyle\frac{1}{24M_{5}^{3}}\widetilde{V}(\varphi)\equiv\kappa^{4}\,\varphi^{4}\,v(\varphi)\>,\>y=\kappa\,\varphi\,T^{-1}\>,\>\>x=\kappa\,r\,T\>,$
(27)

the thick wall bounce action looks like

$\displaystyle S_{b}$ $\displaystyle=\frac{96\pi
M_{5}^{3}}{\kappa^{3}}\int\differential
x\,x^{2}\left(\frac{1}{2}\left(\frac{\differential y}{\differential
x}\right)^{2}+y^{4}\,v(Ty/\kappa)+\frac{\pi^{4}}{12}\right)\>.$ (28)

The overall factor of $M_{5}^{3}\sim N_{c}^{2}$ makes the $N_{c}^{2}$
dependence of the bounce action manifest.

For the $O(3)$ symmetric bubbles with bounce action $S_{b}$, the tunnelling
probability per unit time per unit 3-volume, at temperature $T$ is given by

$\displaystyle\Gamma(T)=T^{4}e^{-S_{b}(T)}\>,$ (29)

where we have ignored constant $\mathcal{O}(1)$ multiplicative factors in the
above. For the phase transition to complete, we need the probability in a
Hubble volume to be at least of order 1, which translates to

$\displaystyle\Gamma(T)\gtrsim H^{4}\>,$ (30)

where $H$ is the Hubble constant in the BB phase, and is fixed by the
Friedmann equations

$\displaystyle H^{2}=\frac{\rho_{\text{total}}}{3M_{\text{pl}}^{2}}\>.$ (31)

The energy density $\rho_{\text{total}}$ gets contributions from the vacuum
energy and from the radiation (since the false vacuum is at a temperature
$T$). Recall that the vacuum energy is tuned to zero at the RS minimum, so
that when BB phase is meta-stable, it has higher free energy than the RS
phase, and therefore has a positive vacuum energy. This is given as

$\displaystyle\rho_{\text{vac, BB}}$
$\displaystyle=2\pi^{4}M_{5}^{3}T_{c}^{4}\>.$ (32)

The energy density from radiation scales as $T^{4}$ and quickly becomes
subdominant as $T\lesssim T_{c}$, which is necessary for the RS phase to
become stable. The condition for the phase transition to complete becomes

$\displaystyle
T^{4}e^{-S_{b}(T)}>\frac{4\pi^{8}}{9}\frac{M_{5}^{6}T_{c}^{8}}{M_{\text{Pl}}^{4}}\>,$
(33)

or equivalently

$\displaystyle
S_{b}(T)<-\log\left(\frac{4\pi^{8}}{9}\frac{M_{5}^{6}T_{c}^{8}}{T^{4}M_{\text{Pl}}^{4}}\right)\equiv
S_{b}^{\text{max}}(T/T_{c})\>.$ (34)

Defining the nucleation temperature $T_{n}$ as the temperature at which
$\Gamma(T_{n})=H^{4}$, we get

$\displaystyle S_{b}(T_{n})=S_{b}^{\text{max}}(T_{n}/T_{c})\>.$ (35)

Ignoring order 1 factors, for $T_{c}$ around the TeV scale,
$S_{b}^{\text{max}}\sim 140$ and changes slowly as a function of $T/T_{c}$.
Using $M_{\text{pl}}^{2}=M_{5}^{3}$ ($\ell_{\text{AdS}}=1$),
$S_{b}^{\text{max}}$ is independent of the number of colors $N_{c}$. The
bounce action $S_{b}(T)$ on the other hand, scales as $N_{c}^{2}$. Using this
we can effectively obtain a bound on the maximum $N_{c}$ that allows a phase
transition to complete, at a given $T_{n}/T_{c}$ (i.e. the amount of
supercooling):

$\displaystyle
N_{\text{max}}(T_{n}/T_{c})=\sqrt{\frac{S_{b}^{\text{max}}(T_{n}/T_{c})}{S_{b,N_{c}=1}(T_{n})}}\>.$
(36)

In the next section, for the parameter choices in table 1, we present the
results for $S_{b}$ and $S_{b}^{\text{max}}$, both in the thin wall limit and
away from it. We find that the bounce action is reduced for some of the
parameter choices, and is correlated with increasing the radion mass.

##  6 Results

###  6.1 Bounce action and maximum $N$

The left panel of fig. 3(a) shows the bounce action (normalized by
$16\pi^{2}M_{5}^{3}=N_{c}^{2}$) in the thin wall limit (solid lines) and away
from it (dashed lines), for the parameter choices in table 1. Also shown is
the maximum bounce action (dotted) beyond which the rate is too small to
compete with the Hubble expansion. We note a couple of things in fig. 3(a).
The thin wall curves are only applicable for $T/T_{c}\lesssim 1$, and for all
the parameter choices, they are above the $S_{b}^{\text{max}}$ line, so that
the phase transition rate is too small even for $N_{c}=1$, in the thin-wall
limit. The phase transition completes in the thick-wall case, at a temperature
$T_{n}$ where the $S_{b}$ and $S_{b}^{\text{max}}$ curves intersect (which in
turn depends on $N_{c}$ because $S_{b}$ scales as $N_{c}^{2}$). For a non-zero
$\epsilon_{3}$, the thick wall bounce action is smaller than the
$\epsilon_{3}=0$ case (e.g. red vs blue curves in fig. 3(a)), and is a slowly
varying function of $T/T_{c}$. Further, the variation in the thick wall bounce
action as a function of $T/T_{c}$ is larger in the presence of a non-zero
$\epsilon_{3}$, which has consequences for GW signal, as we explain later. For
convenience, we also show the corresponding radion potential in the right
panel of fig. 3(a). We can notice a correlation between the second derivative
of the radion potential at the minimum and the bounce action. When the second
derivative is larger, i.e. the physical radion is heavier, the potential is
deeper and the bounce action is lower. For a given amount of supercooling, say
$T/T_{c}=10^{-3}$, the bounce action can be lowered approximately by a factor
of 30 due to a non-zero $\epsilon_{3}$ (e.g. comparing the red and the blue
curves in fig. 3(a)).

(a) Left: The bounce action normalized by $16\pi^{2}M_{5}^{3}=N_{c}^{2}$, for
the parameter choices in table 1. Both thin wall (solid lines) and thick wall
(dashed lines) estimates for the bounce action are shown. Also shown in dotted
is the maximum value of bounce action $S_{b}^{\text{max}}$, beyond which the
phase transition rate to too small to compete with the Hubble expansion.
Right: The corresponding radion potential for the parameter choices in table
1.

(b) Maximum number of colors for the phase transition to complete, as a
function of $T_{n}/T_{c}$, for the parameter choice in table 1.

Figure 3: Bounce action, potential, and the maximum number of colors for the
phase transition to complete, for the parameters in table 1.

Using eq. (36), we can calculate the maximum color $N_{\text{max}}$ in the
dual theory for which the phase transition can complete for a given amount of
supercooling. Figure 3(b) shows $N_{\text{max}}$ as a function of
$T_{n}/T_{c}$. For some amount of supercooling, one can have
$N_{\text{max}}\sim 10$ (e.g. the orange and red curves in fig. 3(b)). Such
values of $N_{\text{max}}$ are obtained at $T_{n}/T_{c}\sim 10^{-6}$, where
the bounce action is becoming smaller than $\mathcal{O}(1)$ (e.g. the orange
and red curves in Fig. 3(a)). $N_{\text{max}}\sim 10$ is therefore the most
one can hope in the present analysis, because we have only computed a part of
the bounce action that was supposed to be dominant. Once this calculable part
is reduced beyond $\mathcal{O}(1)$ values, one cannot use this as a complete
answer. Taking the lowest bounce action to be $\mathcal{O}(1)$, eq. (36) gives
$N_{\text{max}}\sim\sqrt{S_{b}^{\text{max}}}\sim 12$ for
$S_{b}^{\text{max}}\sim 140$.

The amount of supercooling that the metastable phase experiences has important
phenomenological implications. During this period of supercooling, the
universe inflates, and this dilutes any matter abundances generated before the
phase transition. If the framework is to address dark matter abundance and
baryon asymmetry with sufficient supercooling, they must be generated after
the phase transition completes. Supercooling also has implications for
primordial black hole generation mechanisms, cosmological aspects of axion
model building, topological defects and so on (e.g. see Pomarol:2019aae ;
Gouttenoire:2023naa  and references therein). A key more generic observation
is that bubble collision at the end of the phase transition can generate
gravitational waves, whose frequency and abundance are dependent on the amount
of supercooling. For all these applications but especially with an eye to the
latter, we estimate the peak frequency and the GW abundance from bubble
collisions in the next subsection.

###  6.2 Gravitational wave signal

First-order phase transitions can give rise to stochastic gravitational wave
signals, which can potentially be detected in ground and space-based GW
detectors, according to their characteristic frequency and abundance (see
refs. Caprini:2015zlo ; Caprini:2019egz  for a review). There are three
sources of GW production in a first-order phase transition—bubble collision,
sound waves, and turbulence in the plasma. The last two sources require
detailed numerical analysis; here we focus on the bubble collision as a source
of stochastic GW production.

The signal in GW can be characterized by the peak frequency $f_{p}$, and the
frequency dependence of the fractional abundance $\Omega_{\text{GW}}\,h^{2}$.
For the GWs generated by bubble collision, these two quantities are given by
Caprini:2015zlo

$\displaystyle f_{p}$ $\displaystyle=0.037\text{ mHz
}\>\Bigg{(}\frac{\beta}{H}\Bigg{)}\>\Bigg{(}\frac{T_{*}}{\text{TeV}}\Bigg{)}\>\Bigg{(}\frac{g_{*}}{100}\Bigg{)}^{1/6}\>,$
$\displaystyle\Omega_{\text{GW}}\,h^{2}(f)$ $\displaystyle=1.3\times
10^{-6}\>\>\Bigg{(}\frac{H}{\beta}\Bigg{)}^{2}\>\Bigg{(}\frac{100}{g_{*}}\Bigg{)}^{1/3}\frac{3.8(f/f_{p})^{2.8}}{1+2.8(f/f_{p})^{3.8}}\>.$
(37)

Here we have assumed that when the transferred latent heat is large compared
to the energy of the surrounding plasma, all the latent heat is transferred to
the bubble wall, and the bubble wall velocity is ultra-relativistic. $H$ here
is the Hubble scale during the phase transition, $T_{*}$ is the temperature of
the radiation bath right after the phase transition, and $g_{*}$ is the number
of relativistic degrees of freedom in the plasma during the phase transition.
The parameter $\beta/H$ is related to the duration of the phase transition
maggiore2018gravitational , and can be calculated from the bounce action as:

$\displaystyle\frac{\beta}{H}$
$\displaystyle=\left.-\frac{\differential\log\Gamma}{\differential\log
T}\right|_{T=T_{n}}\approx-4+\left.\frac{\differential
S_{b}}{\differential\log T}\right|_{T=T_{n}}\>,$ (38)

where we have used eq. (29) in the last equality above. Note that $\beta/H$ is
a function of $T_{n}$, the temperature at which the phase transition can
proceed, which in turn depends on $N_{c}$. A small $\beta/H$ decreases the
peak frequency, but crucially increases the abundance. Figure 4 shows
$\beta/H$ as a function of $T_{n}/T_{c}$ for the parameter choices in table 1.
On each curve, $N_{c}$ varies, which changes $T_{n}$. Values of $N_{c}=2,5,10$
are shown on the individual curves. The behavior of $\beta/H$ with and without
$\epsilon_{3}$ is very different, as seen for example by the blue and the red
curves in fig. 4. When $\epsilon_{3}=0$, $\beta/H$ is small and decreases as
$T_{n}/T_{c}$ decreases, unlike the $\epsilon_{3}\neq 0$ case. When the
deformation in the IR is very small so that one is close to a CFT, the bounce
action is weakly dependent on temperature Agashe:2020lfz , and $\beta/H$ is
close to zero. This is different than the behavior seen for $\epsilon_{3}\neq
0$.

Figure 4: The parameter $\beta/H$ as a function of $T_{n}/T_{c}$ for the
parameter choices in table 1. On each curve, $N$ varies as a parameter. Points
with $N_{c}=2,5,10$ are shown by markers.

The increased peak frequency and reduced GW signal from a non-zero
$\epsilon_{3}$ is still within the reach of proposed space-based GW
interferometers such as LISA LISA:2017pwj ; Baker:2019nia , DECIGO Seto:2001qf
; Kawamura:2011zz ; Yagi:2011wg ; Isoyama:2018rjb , and BBO Crowder:2005nr ;
Corbin:2005ny ; Harry:2006fi . In fig. 5 we show the expected GW abundance as
a function of frequency, for some values of $\beta/H$, and compare to the
projected reach of the experiments (taking $g_{*}=100$ and $T_{*}=\text{TeV}$
in eq. (37)). For a given value of $N_{c}$, the values of $\beta/H$ are
different for the four parameter choices in table 1, and are at different
values of $T_{n}/T_{c}$. The left panel of fig. 5 shows the GW signal for
$\beta/H=20,50,200,350$ (corresponding to $N_{c}=2$, shown by “x” in fig. 4).
The value of $\beta/H$ changes significantly for a non-zero $\epsilon_{3}$:
from $20$ ($\epsilon_{3}=0$, $T_{n}/T_{c}=0.002$, parameter A) to $350$
($\epsilon_{3}=-1/81$, $T_{n}/T_{c}=0.7$, parameter D). Such different values
of $\beta/H$ correspond to very different peak frequencies and fractional
abundances (e.g. the red and blue curves in the left panel of fig. 5). For
comparison, we also show the GW signal for $\beta/H=10,30,120$ in the right
panel of fig. 5 (corresponding to $N_{c}=5$, shown by “o” in fig. 4). For both
choices of $N_{c}$ we see that while the signal strength is reduced and peak
frequency is increased due to a non-zero $\epsilon_{3}$, there is still a
possibility of discovery.

Figure 5: GW abundance as a function of frequency for some values of
$\beta/H$: values in the left panel correspond to $N_{c}=2$ (indicated in fig.
4 by a “x”), values in the right panel correspond to $N_{c}=5$ (indicated in
fig. 4 by a “o”). Also shown is the projected reach from LISA LISA:2017pwj ;
Baker:2019nia , DECIGO Seto:2001qf ; Kawamura:2011zz ; Yagi:2011wg ;
Isoyama:2018rjb , and BBO Crowder:2005nr ; Corbin:2005ny ; Harry:2006fi . We
have taken $g_{*}=100$ and $T_{*}=\text{TeV}$.

##  7 General comments

In this work we have argued for including self-interaction terms in the bulk
potential of the stabilizing scalar. We considered a bulk potential with a
quadratic and a cubic term, with signs chosen such that both terms grow in the
IR. We stayed in the limit of small back-reaction and within the radion EFT.
In the presence of a cubic term, and for a large hierarchy, the defining
features of the radion potential are changed, and the radion mass is
increased. The same effect also reduced the bounce action and thereby
increased the rate for a transition from the hot phase to the RS phase.
Equivalently, this reduced the amount of supercooling needed before completing
the phase transition and increased the maximum number of colors
$N_{\text{max}}$ for which the phase transition can complete, for a given
amount of supercooling. For a choice of parameters, we were able to have
$N_{\text{max}}\sim 10$. We also discussed the resulting GW signals from
bubble collisions, and showed that in the presence of self-interactions, the
parameter $\beta/H$ that characterizes the frequency and abundance of the GWs
does not get too small, unlike the case when the CFT breaking in the IR is
small and $\beta/H$ is close to zero.

In fig. 6 we show $N_{\text{max}}$ and $\beta/H$ for a moderate amount of
supercooling $T_{n}/T_{c}=10^{-4}$, as a function of the mass squared of the
physical radion. The figure summarizes the main point of the paper—the
presence of a self-interaction term increases the mass of the physical radion
and the same effect also reduces the bounce action to increase
$N_{\text{max}}$, while disfavoring a small $\beta/H$. As the radion mass gets
close to the KK scale, we cannot trust the result entirely, and a complete 5D
calculation would be necessary. Our conclusion is similar to Csaki:2023pwy ,
even though the underlying dynamics driving the radion mass up is different.

Figure 6: The effect on $N_{\text{max}}$ (left) and $\beta/H$ (right) as the
radion mass squared increases due to a non-zero $\epsilon_{3}$. The curves are
shown for a fixed $T_{n}/T_{c}=10^{-4}$.

Our work is a first step towards systematically including the effects of
strong coupling in the IR, in the dynamics of the phase transition. In this
work we have chosen parameters such that the radion is lighter than the other
KK masses, allowing us to use the radion EFT. A more generic situation is when
the radion mass is of the same order as other KK masses, which would drive the
calculation out of the radion EFT and a full 5D calculation would be needed.
The back-reaction would be important in the IR and it would change the free
energy of both the RS and the BB phase. The stability of the phases can change
the order of the phase transition, and a full 5D gravitational instanton
computation would be needed to address the question properly. Rather than
modeling choices, these are general effects. We will address some of these
issues in a more general model in a future work.

## Acknowledgements

We would like to thank P. Creminelli, P. Du and A. Pomarol for comments and
discussions at various stages. We would also like to thank N. DePorzio, J.
Lodman and W. L. Xu for collaborations at an early stage. The work of RKM and
LR is supported by the National Science Foundation under Grant Nos.
PHY-1620806, PHY-1748958 and PHY-1915071, the Chau Foundation HS Chau postdoc
award, the Kavli Foundation grant “Kavli Dream Team,” and the Moore Foundation
Award 8342. Part of this work was completed at Aspen Center for Physics, which
is supported by National Science Foundation grant PHY-2210452.

##  Appendix A Numerical and approximately analytical method

In this appendix we briefly outline the procedure to obtain the approximate
analytical solutions for the scalar $\chi$, for both BB and RS backgrounds. We
also describe the numerical method to obtain the profile for $\chi$ in the BB
background.

The equation of motion for $\chi$ in the background of eq. (1) is

$\displaystyle\left(1-e^{4(r-r_{h})}\right)\frac{\differential^{2}\chi}{\differential
r^{2}}-4\frac{\differential\chi}{\differential r}-\frac{\differential
V_{B}(\chi)}{\differential\chi}-\sum_{i}\delta(r-r_{i})\frac{\differential
V_{i}(\chi)}{\differential\chi}=0\>,$ (39)

where $i$ runs over $(uv,ir)$ for the RS case and only over $(uv)$ for the BB
case. The range of $r$ is $0\leq r\leq r_{\text{ir}}$ for the RS case and
$0\leq r\leq r_{h}$ for the BB case. We will work with a rescaled coordinate
$y=r/r_{\text{ir}}$ for the RS case, and $y=r/r_{h}$ for the BB case, so that
$0\leq y\leq 1$ for both backgrounds. Given the bulk and brane potentials in
eq. (5), (6), the equation and the boundary conditions for the RS and BB
backgrounds are

RS:
$\displaystyle\frac{1}{r_{\text{ir}}^{2}}\frac{\differential^{2}\chi}{\differential
y^{2}}-\frac{4}{r_{\text{ir}}}\frac{\differential\chi}{\differential
y}-4\epsilon_{2}\chi-4\epsilon_{3}\chi^{2}=0,\>\>\chi(0)=v_{\text{uv}},\>\>\chi^{\prime}(1)=-r_{\text{ir}}\alpha_{\text{ir}}\>,$
(40) BB:
$\displaystyle\left(\frac{1-e^{4r_{h}(y-1)}}{r_{h}^{2}}\right)\frac{\differential^{2}\chi}{\differential
y^{2}}-\frac{4}{r_{h}}\frac{\differential\chi}{\differential
y}-4\epsilon_{2}\chi-4\epsilon_{3}\chi^{2}=0,\>\>\chi(0)=v_{\text{uv}}\>.$
(41)

###  A.1 RS background

In the limit $r_{\text{ir}}\gg 1,|\epsilon_{2}|r_{\text{ir}}\lesssim
1,|\epsilon_{3}|r_{\text{ir}}\lesssim 1$, eq (40) can be solved approximately
by singular perturbation theory and boundary layer analysis bender1999advanced
. Close to $r=y=0$, $\chi$ varies slowly and the derivatives are small. To
leading order in $1/r_{\text{ir}}$, one can drop the second derivative in eq.
(40), which results in a first-order differential equation, readily solved.
Applying the boundary condition $\chi(0)=v_{\text{uv}}$ to the first-order
differential equation resulting from dropping the second derivative term in
eq. (40), we get

$\displaystyle\chi_{\text{left}}(y)=\frac{v_{\text{uv}}e^{-\epsilon_{2}r_{\text{ir}}y}}{1+v_{\text{uv}}\epsilon_{3}\left(\frac{1-e^{-\epsilon_{2}r_{\text{ir}}y}}{\epsilon_{2}}\right)}\>,$
(42)

where the subscript “left” refers to the fact that the solution is valid only
on the left side of the interval. As $y\to 1$, $\chi$ has to change fast to
match the boundary condition at $y=1$. In a small region near $y=1$, of width
$\mathcal{O}(1/r_{\text{ir}})$, $\chi$ itself is small, but changes fast so
that the derivatives are large. In this region (referred to as “right”) one
only keeps the first and the second derivative terms in eq. (40), which gives

$\displaystyle\chi_{\text{right}}(y)=-\frac{\alpha_{\text{ir}}}{4}e^{4r_{\text{ir}}(y-1)}+C\>,$
(43)

where we have applied the boundary condition
$\chi^{\prime}(1)=-r_{\text{ir}}\alpha_{\text{ir}}$, and $C$ is an
undetermined constant at the moment. Using asymptotic matching in a region
where both the solutions are valid, and requiring the same functional form,
the constant $C$ can be fixed. The final solution is given as (now switching
to $r=r_{\text{ir}}\,y$)

$\displaystyle\chi_{\text{RS}}(r)=-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}+\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}})}\>.$
(44)

For consistency, we can check if the equations of motion and the boundary
conditions are satisfied. In the limit of $r_{\text{ir}}\gg 1,v_{\text{uv}}\ll
1,|\epsilon_{2}|\ll 1,|\epsilon_{3}|\ll 1$, the errors on the boundary and in
the bulk are small, and under control. Note that for $\epsilon_{2}<0$, the
denominator of eq. (42) can vanish at some $r$. Clearly that is outside the
validity of the approximation since the second derivative would not be small
there anymore. To be consistent, we therefore need
$-\epsilon_{2}r_{\text{ir}}<\log(1+\epsilon_{2}/v_{\text{uv}}\epsilon_{3})$.
For $\epsilon_{2}<0,\epsilon_{3}<0$, in the limit of small $\epsilon_{2}$,
this condition simplifies to $r_{\text{ir}}<1/v_{\text{uv}}|\epsilon_{3}|$,
which can be satisfied for small $v_{\text{uv}}$.

We can understand intuitively the role of $\epsilon_{3}$ as follows. For small
enough $r$ one can expand the exponential in the denominator in the second
term of eq. (44), and approximate it as follows:

$\displaystyle\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}})}\approx\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}r}\approx
v_{\text{uv}}e^{-(\epsilon_{2}+v_{\text{uv}}\epsilon_{3})r}\>.$ (45)

This makes it clear that for small $r$, $\epsilon_{3}$ effectively increases
$\epsilon_{2}$. This effect compounds as $r$ increases.

###  A.2 BB background

Similar to the RS case, eq. (41) can be solved in the limit of large $r_{h}$.
For $y=r/r_{h}\ll 1$, the exponential term in the coefficient of
$\chi^{\prime\prime}$ can be dropped, and we have the same solution as eq.
(42). Note that the UV boundary condition is the same for both RS and BB
cases. As $y\to 1$, unlike the RS case, the coefficient of
$\chi^{\prime\prime}$ becomes small, so that one can again drop the second
derivative term. Since the only condition on $\chi(r)$ is to be regular at
$y=1$, the leading order solution is

$\displaystyle\chi_{\text{BB}}(r)=\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}})}\>.$
(46)

To be more precise, we can do a Taylor series expansion around $y=1$ to solve
the equation near $y=1$ and match it to $\chi_{\text{left}}(r)$ at some
intermediate value. This procedure gives the same result as above, to leading
order.

We would like to check the approximate solution (46) with the numerical
solution for the profile. Without an explicit boundary condition to be applied
at $y=1$, it is not clear a priori how to numerically solve for $\chi$ in the
BB geometry. For this we apply a method that is akin to matching in an
intermediate region. We solve the full equation numerically with given
$\chi(0)$, and for some value of $\chi^{\prime}(0)$, varying
$\chi^{\prime}(0)$. We also expand the equation in a Taylor series expansion
around $y=1$ which analytically fixes $\chi^{\prime}(1)$ in terms of
$\chi(1)$. As we vary $\chi^{\prime}(0)$ in the numerical solution, we check
whether the analytical relation between $\chi(1)$ and $\chi^{\prime}(1)$ is
satisfied, which uniquely determines $\chi^{\prime}(0)$, and hence the
numerically regular solution.

In the special case of $\epsilon_{3}=0$, the solution is a linear combination
of hypergeometric functions Creminelli:2001th , and the correct linear
combination regular at the horizon is easily identified. Figure 1 shows a
comparison between the solutions: for $\epsilon_{3}\neq 0$, between the
approximate solutions obtained in this appendix and the numerical solution
obtained by the method discussed in the previous paragraph, and for
$\epsilon_{3}=0$, between the exact solution and the approximate solution.

##  Appendix B Radion potential

Starting with the purely gravitational 5D action

$\displaystyle S=\int
d^{5}x\sqrt{g}\left(-2M_{5}^{3}\mathcal{R}[g]-\Lambda_{5}\right)-\sum_{i=uv,ir}\int
d^{4}x\sqrt{g_{i}}\>T_{i}\>,$ (47)

where $-\Lambda_{5}=T_{\text{uv}}=-T_{\text{ir}}=24M_{5}^{3}$ (setting
$\ell_{\text{AdS}}=1$), and plugging back the metric

$\displaystyle ds^{2}=-e^{-2r}dx^{2}+dr^{2}\>,\qquad 0\leq r\leq
r_{\text{ir}}\>,$ (48)

but making $r_{\text{ir}}$ a 4D field $r_{\text{ir}}(x)$, one generates a
kinetic term for $r_{\text{ir}}(x)$. In terms of
$\varphi=\exp(-r_{\text{ir}})$, the 4D action looks like

$\displaystyle S=-12M_{5}^{3}\int d^{4}x(\partial\varphi)^{2}\>.$ (49)

At this point, there is no potential for $\varphi$ and it is a modulus. To
generate a potential to stabilize the geometry, we add a GW scalar $\chi$ with
the action defined in eq. (4) and solve for the background value of $\chi$
which is a function of $r$ due to the choice of boundary and bulk potentials.
Evaluating $S_{\chi}$ on $\chi(r)$ gives the potential for $\varphi$, for
which we outline the steps now.

Plugging $\chi(r)$ in $S_{\chi}$, we get

$\displaystyle S_{\chi}$
$\displaystyle=\int\differential^{4}x\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(-\frac{1}{2}\left(\frac{\differential\chi}{\differential
r}\right)^{2}-V_{B}(\chi)\right)-\sum_{i=uv,ir}e^{-4r_{i}}V_{i}(\chi(r_{i}))\>.$
(50)

Using $\chi$ equations of motion
$\chi^{\prime\prime}-4\chi^{\prime}-\partial_{\chi}V_{B}(\chi)=0$ and an
integration by parts, the bulk term of $S_{\chi}$ evaluates to

$\displaystyle\int\differential^{4}x\left(-\frac{1}{2}\left.\left(e^{-4r}\chi\frac{\differential\chi}{\differential
r}\right)\right|_{0}^{r_{\text{ir}}}-\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(V_{B}(\chi)-\frac{1}{2}\chi\frac{\differential
V_{B}}{\differential\chi}\right)\right)\>.$ (51)

Adding the contribution from the $uv,ir$ localized terms (in eq. (50)) to the
bulk contribution in eq. (51), the potential is given as

$\displaystyle V(r_{\text{ir}})$
$\displaystyle=\frac{1}{2}\left.\left(e^{-4r}\chi\frac{\differential\chi}{\differential
r}\right)\right|_{0}^{r_{\text{ir}}}+\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(V_{B}(\chi)-\frac{1}{2}\chi\frac{\differential
V_{B}}{\differential\chi}\right)+V_{\text{uv}}(\chi(0))+e^{-4r_{\text{ir}}}V_{\text{ir}}(\chi(r_{\text{ir}}))\>.$
(52)

For Dirichlet boundary condition in UV, $V_{\text{uv}}(\chi)=0$. Given the
bulk potential in eq. (5),
$V_{B}-(1/2)\chi\partial_{\chi}V_{B}=-(2/3)\epsilon_{3}\chi^{3}$. Using
$\chi(r)$ from eq. (7), we need to evaluate the integral

$\displaystyle-\frac{2}{3}\epsilon_{3}\,\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}+\frac{v_{\text{uv}}e^{-\epsilon_{2}r}}{1+v_{\text{uv}}\epsilon_{3}\left(\frac{1-e^{-\epsilon_{2}r}}{\epsilon_{2}}\right)}\right)^{3}\>.$
(53)

To proceed, we need some identities related to hypergeometric functions. First
note that each term in the above is an integral of a general form, with a
closed-form answer expressed in terms of hypergeometrics:

$\displaystyle\int\differential
r\frac{e^{ar}}{(1+be^{cr})^{n}}=\frac{e^{ar}}{(1+be^{cr})^{n-1}}\,{}_{2}F_{1}(1,\,1-n+\frac{a}{c},\,1+\frac{a}{c},\,be^{cr})\>.$
(54)

Two other useful identities are Cvitkovi__2017 :

$\displaystyle{}_{2}F_{1}(a,b,c;z)$
$\displaystyle=(1-z)^{c-a-b}{}_{2}F_{1}(c-a,c-b,c;z)\>,$
$\displaystyle{}_{2}F_{1}(a+\delta\lambda,b,c+\lambda)$
$\displaystyle\approx\left(1-\delta z\right)^{-b}\>,\qquad|\lambda|\gg
1\>,|\delta|\leq 1\>.$ (55)

To simplify the further expressions, we define

$\displaystyle\epsilon_{32}=\frac{\epsilon_{3}}{\epsilon_{2}}\>,\>\>\lambda=\frac{v_{\text{uv}}\epsilon_{32}}{1+v_{\text{uv}}\epsilon_{32}}\>,\>\>Y(r)=\frac{\lambda
e^{-\epsilon_{2}r}}{1-\lambda e^{-\epsilon_{2}r}}\>,$ (56)

in terms of which $\chi(r)$ can be written as

$\displaystyle\chi(r)=-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}+\frac{1}{\epsilon_{32}}Y(r)\>.$
(57)

We now expand the integrand in eq. (53), use eqs. (54), (55) for
simplification and the notation of eq. (56). Calling the four terms from
expanding eq. (53) as $t_{1},t_{2},t_{3},t_{4}$, we have

$\displaystyle t_{1}$
$\displaystyle=-\frac{2}{3}\epsilon_{3}\,\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}\right)^{3}=\frac{1}{96}\epsilon_{3}\alpha_{\text{ir}}^{3}e^{-12r_{\text{ir}}}\int_{0}^{r_{\text{ir}}}\differential
re^{8r}$
$\displaystyle=\frac{1}{768}\epsilon_{3}\alpha_{\text{ir}}^{3}\left(\varphi^{4}-\varphi^{12}\right)\>.$
(58) $\displaystyle t_{2}$
$\displaystyle=-\frac{2}{3}\epsilon_{3}\,\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\,3\,\left(-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}\right)\left(\frac{1}{\epsilon_{32}}Y(r)\right)^{2}$
$\displaystyle=\frac{1}{2}\frac{\epsilon_{2}^{2}}{\epsilon_{3}}\alpha_{\text{ir}}\,e^{-4r_{\text{ir}}}\,\int_{0}^{r_{\text{ir}}}\differential
r\frac{\lambda^{2}e^{-2\epsilon_{2}r}}{(1-\lambda
e^{-\epsilon_{2}r})^{2}}=\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\,e^{-4r_{\text{ir}}}\,\left.\Big{(}-Y(r)+\log(1+Y(r))\Big{)}\right|_{0}^{r_{\text{ir}}}$
$\displaystyle=\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\varphi^{4}\left(-Y(r_{\text{ir}})+\log(1+Y(r_{\text{ir}}))+\frac{\lambda}{1-\lambda}+\log(1-\lambda)\right)\>.$
(59) $\displaystyle t_{3}$
$\displaystyle=-\frac{2}{3}\epsilon_{3}\,\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\,3\,\left(-\frac{\alpha_{\text{ir}}}{4}e^{4(r-r_{\text{ir}})}\right)^{2}\left(\frac{1}{\epsilon_{32}}Y(r)\right)$
$\displaystyle=-\frac{1}{8}\epsilon_{2}\alpha_{\text{ir}}^{2}e^{-8r_{\text{ir}}}\int_{0}^{r_{\text{ir}}}\differential
r\frac{\lambda e^{(4-\epsilon_{2})r}}{(1-\lambda
e^{-\epsilon_{2}r})}=\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}e^{-8r_{\text{ir}}}\left.\left(e^{4r}{}_{2}F_{1}(1,\frac{4}{\epsilon_{2}},1+\frac{4}{\epsilon_{2}};\lambda^{-1}e^{\epsilon_{2}r})\right)\right|_{0}^{r_{\text{ir}}}$
$\displaystyle=\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}e^{-8r_{\text{ir}}}\left.\Big{(}-e^{4r}Y(r)\Big{)}\right|_{0}^{r_{\text{ir}}}=\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}\left(-\varphi^{4}Y(r_{\text{ir}})+\varphi^{8}\left(\frac{\lambda}{1-\lambda}\right)\right)\>.$
(60) $\displaystyle t_{4}$
$\displaystyle=-\frac{2}{3}\epsilon_{3}\,\int_{0}^{r_{\text{ir}}}\differential
r\,e^{-4r}\left(\frac{1}{\epsilon_{32}}Y(r)\right)^{3}=-\frac{2}{3}\frac{\epsilon_{2}^{3}}{\epsilon_{3}^{2}}\int_{0}^{r_{\text{ir}}}\differential
r\frac{\lambda^{3}e^{-(4+3\epsilon_{2})r}}{(1-\lambda
e^{-\epsilon_{2}r})^{3}}$
$\displaystyle=-\frac{1}{6}\frac{\epsilon_{2}^{3}}{\epsilon_{3}^{2}}\left.\left(e^{-4r}{}_{2}F_{1}(3,-\frac{4}{\epsilon_{2}},1-\frac{4}{\epsilon_{2}};\lambda^{-1}e^{\epsilon_{2}r})\right)\right|_{0}^{r_{\text{ir}}}=-\frac{1}{6}\frac{\epsilon_{2}^{3}}{\epsilon_{3}^{2}}\left.\Big{(}-e^{-4r}Y^{3}(r)\Big{)}\right|_{0}^{r_{\text{ir}}}$
$\displaystyle=\frac{1}{6}\frac{\epsilon_{2}^{3}}{\epsilon_{3}^{2}}\left(\varphi^{4}Y^{3}(r_{\text{ir}})-\left(\frac{\lambda}{1-\lambda}\right)^{3}\right)\>.$
(61)

We have used $|\epsilon_{2}|\ll 1$ in the above to simplify the hypergeometric
functions. Collecting everything together, keeping to leading order in
$\varphi$, and dropping overall constants, the integral in eq. (53) is given
as:

$\displaystyle\varphi^{4}\left(A+BY+CY^{3}+D\log(1+Y)\right)\>,\>\>Y=\frac{\lambda\varphi^{\epsilon_{2}}}{1-\lambda\varphi^{\epsilon_{2}}}$
$\displaystyle
A=\frac{1}{768}\epsilon_{3}\alpha_{\text{ir}}^{2}-\frac{\lambda}{1-\lambda}+\log(1-\lambda)\>,\>\>B=-\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}-\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\>,\>\>C=\frac{1}{6}\frac{\epsilon_{2}^{3}}{\epsilon_{3}^{2}}\>,\>\>D=\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\>.$
(62)

In addition to all the terms in eq. (52), there is a potential generated by
the detuning of the IR brane tension, and is given by

$\displaystyle V_{\text{detune}}(\varphi)$
$\displaystyle=\widetilde{\tau}\,\varphi^{4}\>,$ (63)

where $\widetilde{\tau}$ is a free parameter at this point. Together with
above, including all the contributions in eq. (52), and keeping to linear
order in $\epsilon_{2}$, the radion potential is given as

$\displaystyle
V(\varphi)=\varphi^{4}\left(a_{1}+a_{2}Y+a_{3}\log(1+Y)\right)\>,\>\>Y=\frac{\lambda\varphi^{\epsilon_{2}}}{1-\lambda\varphi^{\epsilon_{2}}}$
$\displaystyle
a_{1}=\tau\>,\>\>a_{2}=-\frac{1}{32}\epsilon_{2}\alpha_{\text{ir}}^{2}-\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}+2\alpha_{\text{ir}}\>,\>\>a_{3}=\frac{1}{2}\frac{\epsilon_{2}}{\epsilon_{3}}\alpha_{\text{ir}}\>,$
(64)

where we have absorbed terms to define an overall $\tau$ and used exact values
for $\chi(0)$ and $\chi^{\prime}(r_{\text{ir}})$ (these are the specified
boundary conditions and are only approximately satisfied by the approximate
solution). The $\epsilon_{3}\to 0$ limit is finite, and in this limit
$a_{2}Y+a_{3}\log(1+Y)$ reduces to
$-(1/2)\alpha_{\text{ir}}v_{\text{uv}}\varphi^{\epsilon_{2}}$. Note that for
$Y$ to stay finite, we need $r_{\text{ir}}<(1/\epsilon_{2})\log\lambda$.

To make it clear what are the reasonable values for parameters, and to keep
the $N$ dependence of the dual theory clear, we define $\tau\equiv
24M_{5}^{3}\kappa^{4}$. Since $M_{5}^{3}\sim\frac{N^{2}}{16\pi^{2}}$, this
makes it clear that $\tau\sim\frac{N^{2}}{16\pi^{2}}\kappa^{4}$ (which is the
glueball normalization) and we need $\kappa\lesssim 1$ for small back-
reaction. The radion potential can be rewritten as

$\displaystyle V(\varphi)$
$\displaystyle=24M_{5}^{3}\,\kappa^{4}\,\varphi^{4}\,v(\varphi/\varphi_{\text{min}};\epsilon_{2},\epsilon_{3},\alpha_{\text{ir}},v_{\text{uv}})\>,$
(65)

where the function $v$ does not have an analytical expression (since
$\varphi_{\text{min}}$ does not have an analytical expression), but is easily
obtained numerically. The dimensionless function $v$ encodes the information
about the breaking of the CFT, and plays a crucial role when calculating the
bounce action. Putting everything together, the radion action is

$\displaystyle
S=24M_{5}^{3}\int\differential^{4}x\left(-\frac{1}{2}(\partial\varphi)^{2}-\kappa^{4}\,\varphi^{4}\,v(\varphi/\varphi_{\text{min}};\epsilon_{2},\epsilon_{3},\alpha_{\text{ir}},v_{\text{uv}})\right)\>.$
(66)

##  Appendix C Bounce action in the thin and thick wall limits

Starting with the action

$\displaystyle S$
$\displaystyle=\int\differential^{4}x\left(-12M_{5}^{2}(\partial\varphi)^{2}-V(\varphi)\right)\>,$
(67)

the bounce action is obtained by looking for a solution to $\varphi$ that
minimizes the Euclidean action

$\displaystyle S_{b}$
$\displaystyle=\int\differential^{4}x\left(12M_{5}^{2}(\partial\varphi)^{2}+\widetilde{V}(\varphi)\right)\>,\qquad\widetilde{V}(\varphi)=V(\varphi)-C\>,$
(68)

and evaluating the action on the solution. The constant $C$ is chosen to
subtract the contribution from the false vacuum, and is given by
$-2\pi^{4}M_{5}^{3}T_{c}^{4}=V(\varphi_{\text{min}})$ for the thin wall case,
and $-2\pi^{4}M_{5}^{3}T^{4}$ for the thick wall case. At zero temperature,
$\varphi$ is assumed to be a function of the combination
$\rho^{2}=\vec{x}\cdot\vec{x}+t_{E}^{2}$ and this is referred as the $O(4)$
symmetric solution. At finite temperature, there is another saddle that can
dominate. For inverse temperature $\beta=1/T$, the Euclidean time is made
periodic with period $\beta$. The field $\varphi$ is assumed to be a function
of the combination $r^{2}=\vec{x}\cdot\vec{x}$ and this is referred to as the
$O(3)$ solution. In rest of the discussion we focus on the $O(3)$ solution
only.

For $O(3)$ symmetric solutions, the action is given more explicitly as
(ignoring temperature corrections to the potential)

$\displaystyle S_{b}$ $\displaystyle=\int_{0}^{\beta}\differential
t_{E}\int\differential^{3}x\left(12M_{5}^{2}(\partial\varphi)^{2}+\widetilde{V}(\varphi)\right)=\frac{4\pi}{T}\int\differential
r\,r^{2}\left(12M_{5}^{2}\left(\frac{\differential\varphi}{\differential
r}\right)^{2}+\widetilde{V}(\varphi)\right)\>.$ (69)

The equations of motion are

$\displaystyle 24M_{5}^{3}\left(\frac{\differential^{2}\varphi}{\differential
r^{2}}+\frac{2}{r}\frac{\differential\varphi}{\differential
r}\right)=\frac{\differential\widetilde{V}(\varphi,T)}{\differential\varphi}\approx\frac{\differential\widetilde{V}(\varphi)}{\differential\varphi}\>.$
(70)

One of the boundary conditions is $\varphi^{\prime}(r=0)=0$ and we will have
more to say about the second boundary condition later.

###  C.1 Thin wall case

In the limit of the two vacua being very degenerate, the solution for
$\varphi$ is such that it has a large region where it is constant (and equal
to its value in the true vacuum) and then changes quickly to the value in the
false vacuum. In this limit, the wall is small compared to the size of the
bubble, hence this is called the thin-wall solution. In this limit, one can
ignore the first derivative term in the equations of motion. Using the
identity $2\differential^{2}y/\differential x^{2}=\differential/\differential
y(\differential y/\differential x)^{2}$, the equations can be solved to give

$\displaystyle 12M_{5}^{3}\left(\frac{\differential\varphi}{\differential
r}\right)^{2}=\widetilde{V}(\varphi)\>,$ (71)

where we used the boundary condition $\varphi^{\prime}(r=0)=0$ and the fact
that at $r=0$, $\varphi=\varphi_{\text{min}}$ and
$\widetilde{V}(\varphi_{\text{min}})=0$. To evaluate the action on this
solution, a somewhat indirect approach is more intuitive. First note that for
a bubble of size $R$, the field is mostly constant for $r\lesssim R$, changes
quickly in the vicinity of $r=R$, and is again constant afterwards. We can
split the action into these three regions. In the region $r\lesssim R$, we can
drop the derivative, $\widetilde{V}(\varphi)$ is a constant, and the integrand
is proportional to $r^{2}$. For $r\sim R$, the factor of $r^{2}$ can be
approximated to be $R^{2}$, and we have to keep both the derivatives and the
potential inside the integral. For $r\gtrsim R$ there is no contribution. We
therefore have

$\displaystyle
S_{b}=\frac{4\pi}{T}\left(\widetilde{V}(\varphi)\int_{0}^{R}\differential
r\,r^{2}+R^{2}\int_{r\sim R}\differential
r\left(12M_{5}^{3}\left(\frac{\differential\varphi}{\differential
r}\right)^{2}+\widetilde{V}(\varphi)\right)\right)\>.$ (72)

Using $\widetilde{V}(\varphi)=\Delta F=F_{\text{FV}}-F_{\text{TV}}>0$ as the
difference in the free energies between the two vacua, we get

$\displaystyle S_{b}=\frac{4\pi}{T}\left(\frac{1}{3}\Delta
FR^{3}+R^{2}S_{1}\right)\>,\>\>\>S_{1}=\int_{r\sim
R}dr\left(12M_{5}^{3}\left(\frac{\differential\varphi}{\differential
r}\right)^{2}+\widetilde{V}(\varphi)\right)\>.$ (73)

Changing variables from $r$ to $\varphi$, $S_{1}$ can be written as

$\displaystyle S_{1}$
$\displaystyle=\sqrt{48M_{5}^{3}}\int_{0}^{\varphi_{\text{min}}}\differential\varphi\sqrt{\widetilde{V}(\varphi)}\>,$
(74)

and is independent of $R$. Since $\varphi$ changes from $0$ to
$\varphi_{\text{min}}$ in the $r\sim R$ region, this fixes the limits of the
integration in eq. (74).555Technically the lower limit is $\varphi\sim T$,
since we are estimating the action from region $iii)$ (see discussion in main
text). Taking the lower limit to be approximately zero does not change the
estimate. Equation (73) gives $S_{b}$ as a function of $R$, which is minimized
at $R=-2S_{1}/\Delta F$, at which the bounce action is given as

$\displaystyle S_{b}=\frac{16\pi}{3T}\frac{S_{1}^{3}}{(\Delta F)^{2}}\>.$ (75)

The free energy difference $\Delta F$ can be written in terms of the critical
temperature $T_{c}$ as

$\displaystyle\Delta F$
$\displaystyle=F_{\text{FV}}-F_{\text{TV}}=2\pi^{4}M_{5}^{3}\left(T^{4}-T_{c}^{4}\right)\>,$
(76)

using which, the bounce action becomes

$\displaystyle S_{b}$
$\displaystyle=\frac{4}{3\pi^{7}M_{5}^{6}}\left(\frac{S_{1}}{T_{c}^{3}}\right)^{3}\frac{T_{c}/T}{\left(1-T^{4}/T_{c}^{4}\right)^{2}}\>.$
(77)

One can also calculate the $\varphi$ profile for thin-wall case, using eq.
(71) as

$\displaystyle-\sqrt{12M_{5}^{3}}\int_{\varphi_{\text{min}}}^{\varphi}\frac{\differential\varphi}{\sqrt{\widetilde{V}(\varphi)}}=r\>,$
(78)

where we have chosen the negative sign of the square root, and used the
boundary condition $\differential\varphi/\differential r=0$ at $r=0$, which by
eq. (71) is at $\varphi=\varphi_{\text{min}}$.

###  C.2 Thick wall case

In the thick wall case, we minimize the action

$\displaystyle S_{b}$
$\displaystyle=\int\differential^{4}x\left(12M_{5}^{2}(\partial\varphi)^{2}+V(\varphi)+2\pi^{4}M_{5}^{3}T^{4}\right)\>,$
(79)

One of the boundary conditions to be satisfied is the standard one:
$\varphi^{\prime}(0)=0$ ($r$ being a radial coordinate). The second boundary
condition is more subtle. The usual second boundary condition is to require
$\varphi(r\to\infty)=\varphi_{\text{false vacuum}}$. In the case at hand,
$\varphi$ is not a dynamical variable in the other phase. Inside the bubble,
close to the boundary, we have $\varphi\sim T\approx 0$, and the energy is
only in the gradient. Outside the bubble and close to the boundary, the energy
is proportional to $T^{4}$. Requiring the energies to match at the bubble
boundary we get Agashe:2019lhy

$\displaystyle
12M_{5}^{3}\left.\left(\frac{\differential\varphi}{\differential
r}\right)^{2}\right|_{\varphi\,\approx\,0}$
$\displaystyle=2\pi^{4}M_{5}^{3}T^{4}\>.$ (80)

We choose the negative sign in the square root because $\varphi$ starts near
$\varphi_{\text{min}}$ at the center of the bubble and decreases to zero at
the boundary, thereby the derivative is negative when $\varphi\to 0$. Putting
everything together, the equations and the boundary conditions are

$\displaystyle 24M_{5}^{3}\left(\frac{\differential^{2}\varphi}{\differential
r^{2}}+\frac{2}{r}\frac{\differential\varphi}{\differential
r}\right)=\frac{\differential\widetilde{V}(\varphi)}{\differential\varphi}\>,\qquad\left.\frac{\differential\varphi}{\differential
r}\right|_{r=0}=0,\>\>\left.\frac{\differential\varphi}{\differential
r}\right|_{\varphi\approx 0}=-\frac{\pi^{2}}{\sqrt{6}}T^{2}\>.$ (81)

Since the location at which $\varphi\approx 0$ is not known a priori, one has
to approach the problem indirectly. The strategy is to solve the differential
equation numerically with $\varphi(0)=\varphi_{0},\varphi^{\prime}(0)=0$,
calculate $r_{*}$ such that $\varphi(r_{*})\approx 0$, and adjust
$\varphi_{0}$ till $\varphi^{\prime}(r_{*})$ has the appropriate value. For
numerical convenience we define

$\displaystyle\frac{1}{24M_{5}^{3}}\widetilde{V}(\varphi)\equiv\kappa^{4}\,\varphi^{4}\,v(\varphi)\>,$
$\displaystyle y=\kappa\,\varphi\,T^{-1}\>,\>\>x=\kappa\,r\,T\>.$ (82)

With this rescaling, the action in eq. (69) looks like

$\displaystyle S_{b}$ $\displaystyle=\frac{96\pi
M_{5}^{3}}{\kappa^{3}}\int\differential
x\,x^{2}\left(\frac{1}{2}\left(\frac{\differential y}{\differential
x}\right)^{2}+y^{4}\,v(Ty/\kappa)+\frac{\pi^{4}}{12}\right)\>.$ (83)

The differential equation and the boundary conditions in the rescaled
coordinates are

$\displaystyle\frac{\differential^{2}y}{\differential
x^{2}}+\frac{2}{x}\frac{\differential y}{\differential
x}=\frac{\differential}{\differential
y}y^{4}v(Ty/\kappa)\>,\>\>\left.\frac{\differential y}{\differential
x}\right|_{x=0}=0\>,\>\>\left.\frac{\differential y}{\differential
x}\right|_{y\approx 0}=-\frac{\pi^{2}}{\sqrt{6}}\>.$ (84)

As discussed before, since the value of $x$ at which $y\approx 0$ is not known
before solving the equation, we trade that boundary condition with
$y(0)=y_{0}$, calculate $x_{*}$ such that $|y(x_{*})|<\delta=10^{-1}$ and
adjust $y_{0}$ till $y^{\prime}(x_{*})+\pi^{2}/\sqrt{6}=0$. Note that eq. (84)
depends on $T/T_{c}$ and one has to solve it for different values of $T/T_{c}$
to get the temperature dependence of $S_{b}$.

Figure 7 compares the results for the thin and thick wall cases for a generic
choice of parameters: left shows the bounce action as a function of $T/T_{c}$
for the thin and the thick wall cases, right shows the scalar profile for
several values of $T/T_{c}$ for the thick wall case, and for the thin wall
case where $T\approx T_{c}$ (the two vacua being almost degenerate). As the
thin-wall limit is approached (i.e. the bubble radius gets bigger), the
results for the thick-wall calculation converge to the thin-wall results.

Figure 7: Thin vs thick wall for
$\tau=2,\alpha_{\text{ir}}=5,v_{\text{uv}}=1/5,\epsilon_{2}=-1/10,\epsilon_{3}=0$.
Left shows the bounce action and right shows the bubble profile.

When computing the bounce action and the bubble profiles for $\epsilon_{3}\neq
0$, one technical difficulty is to be kept in mind. As discussed earlier, the
radion potential $V(\varphi)$ has a singularity at a finite but non-zero
$\varphi_{s}=(1/\lambda)^{1/\epsilon_{2}}$, where
$\lambda=v_{\text{uv}}(\epsilon_{3}/\epsilon_{2})/(1+v_{\text{uv}}(\epsilon_{3}/\epsilon_{2}))$.
Given that $T_{c}\lesssim\varphi_{\text{min}}$, if the thick wall profile is
computed at $T/T_{c}\ll 1$, the relevant $\varphi$ probed is of the order
$(T/T_{c})\varphi_{\text{min}}\ll\varphi_{\text{min}}$, and for small enough
$T$, one can be sensitive to the singular point $\varphi=\varphi_{s}$. Since
$\varphi_{s}\ll\varphi_{\text{min}}$, one can calculate for intermediate
temperatures without running into this issue. In this work we consider only
values of $T/T_{c}$ such that we are not sensitive to the singularity.
Further, note that the boundary condition requires imposing a condition when
$y\approx 0$ or equivalently $\varphi\approx 0$. Numerically we require
$\left|y(x_{*})\right|\leq\delta=10^{-1}$ or equivalently
$\left|\varphi\right|\leq T(\delta/\kappa)\sim 10^{-1}T$. As long as $T$ is
not too small, we do not probe the singular region of the potential.

## References

  * (1) P. Agrawal, C. Cesarotti, A. Karch, R. K. Mishra, L. Randall, and R. Sundrum, “Warped Compactifications in Particle Physics, Cosmology and Quantum Gravity,” in Snowmass 2021.  3, 2022.  arXiv:2203.07533 [hep-th].
  * (2) P. Creminelli, A. Nicolis, and R. Rattazzi, “Holography and the electroweak phase transition,” JHEP 03 (2002) 051, arXiv:hep-th/0107141.
  * (3) L. Randall and G. Servant, “Gravitational waves from warped spacetime,” JHEP 05 (2007) 054, arXiv:hep-ph/0607158.
  * (4) T. Konstandin and G. Servant, “Cosmological Consequences of Nearly Conformal Dynamics at the TeV scale,” JCAP 12 (2011) 009, arXiv:1104.4791 [hep-ph].
  * (5) N. Levi, T. Opferkuch, and D. Redigolo, “The supercooling window at weak and strong coupling,” JHEP 02 (2023) 125, arXiv:2212.08085 [hep-ph].
  * (6) LISA Collaboration, P. Amaro-Seoane et al., “Laser Interferometer Space Antenna,” arXiv:1702.00786 [astro-ph.IM].
  * (7) J. Baker et al., “The Laser Interferometer Space Antenna: Unveiling the Millihertz Gravitational Wave Sky,” arXiv:1907.06482 [astro-ph.IM].
  * (8) N. Seto, S. Kawamura, and T. Nakamura, “Possibility of direct measurement of the acceleration of the universe using 0.1-Hz band laser interferometer gravitational wave antenna in space,” Phys. Rev. Lett. 87 (2001) 221103, arXiv:astro-ph/0108011.
  * (9) S. Kawamura et al., “The Japanese space gravitational wave antenna: DECIGO,” Class. Quant. Grav. 28 (2011) 094011.
  * (10) K. Yagi and N. Seto, “Detector configuration of DECIGO/BBO and identification of cosmological neutron-star binaries,” Phys. Rev. D 83 (2011) 044011, arXiv:1101.3940 [astro-ph.CO]. [Erratum: Phys.Rev.D 95, 109901 (2017)].
  * (11) S. Isoyama, H. Nakano, and T. Nakamura, “Multiband Gravitational-Wave Astronomy: Observing binary inspirals with a decihertz detector, B-DECIGO,” PTEP 2018 no. 7, (2018) 073E01, arXiv:1802.06977 [gr-qc].
  * (12) J. Crowder and N. J. Cornish, “Beyond LISA: Exploring future gravitational wave missions,” Phys. Rev. D 72 (2005) 083005, arXiv:gr-qc/0506015.
  * (13) V. Corbin and N. J. Cornish, “Detecting the cosmic gravitational wave background with the big bang observer,” Class. Quant. Grav. 23 (2006) 2435–2446, arXiv:gr-qc/0512039.
  * (14) G. M. Harry, P. Fritschel, D. A. Shaddock, W. Folkner, and E. S. Phinney, “Laser interferometry for the big bang observer,” Class. Quant. Grav. 23 (2006) 4887–4894. [Erratum: Class.Quant.Grav. 23, 7361 (2006)].
  * (15) A. Kosowsky, M. S. Turner, and R. Watkins, “Gravitational radiation from colliding vacuum bubbles,” Phys. Rev. D 45 (1992) 4514–4535.
  * (16) A. Kosowsky and M. S. Turner, “Gravitational radiation from colliding vacuum bubbles: envelope approximation to many bubble collisions,” Phys. Rev. D 47 (1993) 4372–4391, arXiv:astro-ph/9211004.
  * (17) A. Kosowsky, M. S. Turner, and R. Watkins, “Gravitational waves from first order cosmological phase transitions,” Phys. Rev. Lett. 69 (1992) 2026–2029.
  * (18) M. Kamionkowski, A. Kosowsky, and M. S. Turner, “Gravitational radiation from first order phase transitions,” Phys. Rev. D 49 (1994) 2837–2851, arXiv:astro-ph/9310044.
  * (19) G. Nardini, M. Quiros, and A. Wulzer, “A Confining Strong First-Order Electroweak Phase Transition,” JHEP 09 (2007) 077, arXiv:0706.3388 [hep-ph].
  * (20) T. Konstandin, G. Nardini, and M. Quiros, “Gravitational Backreaction Effects on the Holographic Phase Transition,” Phys. Rev. D 82 (2010) 083513, arXiv:1007.1468 [hep-ph].
  * (21) B. von Harling and G. Servant, “QCD-induced Electroweak Phase Transition,” JHEP 01 (2018) 159, arXiv:1711.11554 [hep-ph].
  * (22) B. M. Dillon, B. K. El-Menoufi, S. J. Huber, and J. P. Manuel, “Rapid holographic phase transition with brane-localized curvature,” Phys. Rev. D 98 no. 8, (2018) 086005, arXiv:1708.02953 [hep-th].
  * (23) D. Bunk, J. Hubisz, and B. Jain, “A Perturbative RS I Cosmological Phase Transition,” Eur. Phys. J. C 78 no. 1, (2018) 78, arXiv:1705.00001 [hep-ph].
  * (24) S. Bruggisser, B. Von Harling, O. Matsedonskyi, and G. Servant, “Electroweak Phase Transition and Baryogenesis in Composite Higgs Models,” JHEP 12 (2018) 099, arXiv:1804.07314 [hep-ph].
  * (25) P. Baratella, A. Pomarol, and F. Rompineve, “The Supercooled Universe,” JHEP 03 (2019) 100, arXiv:1812.06996 [hep-ph].
  * (26) E. Megías, G. Nardini, and M. Quirós, “Cosmological Phase Transitions in Warped Space: Gravitational Waves and Collider Signatures,” JHEP 09 (2018) 095, arXiv:1806.04877 [hep-ph].
  * (27) A. Pomarol, O. Pujolas, and L. Salas, “Holographic conformal transition and light scalars,” JHEP 10 (2019) 202, arXiv:1905.02653 [hep-th].
  * (28) K. Agashe, P. Du, M. Ekhterachian, S. Kumar, and R. Sundrum, “Cosmological Phase Transition of Spontaneous Confinement,” JHEP 05 (2020) 086, arXiv:1910.06238 [hep-ph].
  * (29) K. Fujikura, Y. Nakai, and M. Yamada, “A more attractive scheme for radion stabilization and supercooled phase transition,” JHEP 02 (2020) 111, arXiv:1910.07546 [hep-ph].
  * (30) E. Megias, G. Nardini, and M. Quiros, “Gravitational Imprints from Heavy Kaluza-Klein Resonances,” Phys. Rev. D 102 no. 5, (2020) 055004, arXiv:2005.04127 [hep-ph].
  * (31) F. Bigazzi, A. Caddeo, A. L. Cotrone, and A. Paredes, “Fate of false vacua in holographic first-order phase transitions,” JHEP 12 (2020) 200, arXiv:2008.02579 [hep-th].
  * (32) K. Agashe, P. Du, M. Ekhterachian, S. Kumar, and R. Sundrum, “Phase Transitions from the Fifth Dimension,” JHEP 02 (2021) 051, arXiv:2010.04083 [hep-th].
  * (33) P. Agrawal and M. Nee, “Avoided deconfinement in Randall-Sundrum models,” JHEP 10 (2021) 105, arXiv:2103.05646 [hep-ph].
  * (34) C. Csáki, M. Geller, Z. Heller-Algazi, and A. Ismail, “Relevant Dilaton Stabilization,” arXiv:2301.10247 [hep-ph].
  * (35) S. Girmohanta, S. J. Lee, Y. Nakai, and M. Suzuki, “Multi-brane cosmology,” arXiv:2304.05586 [hep-ph].
  * (36) C. Eröncel, J. Hubisz, S. J. Lee, G. Rigo, and B. Sambasivam, “New Horizons in the Holographic Conformal Phase Transition,” arXiv:2305.03773 [hep-ph].
  * (37) W. D. Goldberger and M. B. Wise, “Modulus stabilization with bulk fields,” Phys. Rev. Lett. 83 (1999) 4922–4925, arXiv:hep-ph/9907447.
  * (38) W. D. Goldberger and M. B. Wise, “Phenomenology of a stabilized modulus,” Phys. Lett. B 475 (2000) 275–279, arXiv:hep-ph/9911457.
  * (39) R. Rattazzi and A. Zaffaroni, “Comments on the holographic picture of the Randall-Sundrum model,” JHEP 04 (2001) 021, arXiv:hep-th/0012248.
  * (40) Z. Chacko, R. K. Mishra, and D. Stolarski, “Dynamics of a Stabilized Radion and Duality,” JHEP 09 (2013) 121, arXiv:1304.1795 [hep-ph].
  * (41) E. Witten, “Anti-de Sitter space, thermal phase transition, and confinement in gauge theories,” Adv. Theor. Math. Phys. 2 (1998) 505–532, arXiv:hep-th/9803131.
  * (42) S. R. Coleman, “The Fate of the False Vacuum. 1. Semiclassical Theory,” Phys. Rev. D 15 (1977) 2929–2936. [Erratum: Phys.Rev.D 16, 1248 (1977)].
  * (43) C. G. Callan, Jr. and S. R. Coleman, “The Fate of the False Vacuum. 2. First Quantum Corrections,” Phys. Rev. D 16 (1977) 1762–1768.
  * (44) Y. Gouttenoire and T. Volansky, “Primordial Black Holes from Supercooled Phase Transitions,” arXiv:2305.04942 [hep-ph].
  * (45) C. Caprini et al., “Science with the space-based interferometer eLISA. II: Gravitational waves from cosmological phase transitions,” JCAP 04 (2016) 001, arXiv:1512.06239 [astro-ph.CO].
  * (46) C. Caprini et al., “Detecting gravitational waves from cosmological phase transitions with LISA: an update,” JCAP 03 (2020) 024, arXiv:1910.13125 [astro-ph.CO].
  * (47) M. Maggiore, Gravitational Waves: Volume 2: Astrophysics and Cosmology.  Oxford University Press, 2018.
  * (48) C. M. Bender and S. A. Orszag, Advanced mathematical methods for scientists and engineers I: Asymptotic methods and perturbation theory, vol. 1.  Springer Science & Business Media, 1999.
  * (49) M. Cvitković , A.-S. Smith, and J. Pande, “Asymptotic expansions of the hypergeometric function with two large parameters—application to the partition function of a lattice gas in a field of traps,” Journal of Physics A: Mathematical and Theoretical 50 no. 26, (Jun, 2017) 265206. https://doi.org/10.1088%2F1751-8121%2Faa7213.