# Cosmic Censorship near FLRW spacetimes with negative spatial curvature

David Fajman, Liam Urban  David Fajman Liam Urban Faculty of Physics Faculty
of Mathematics University of Vienna University of Vienna Boltzmanngasse 5
Oskar-Morgenstern-Platz 1 1090 Vienna, Austria 1090 Vienna, Austria
david.fajman@ univie.ac.at liam.urban@ univie.ac.at

###### Abstract.

We consider general initial data for the Einstein scalar-field system on a
closed $3$-manifold $(M,\gamma)$ which is close to data for a Friedman-
Lemaître-Robertson-Walker solution with homogeneous scalar field matter and a
negative Einstein metric $\gamma$ as spatial geometry. We prove that the
maximal globally hyperbolic development of such initial data in the Einstein
scalar-field system is past incomplete in the contracting direction and
exhibits stable collapse into a Big Bang curvature singularity. Under an
additional condition on the first positive eigenvalue of $-\Delta_{\gamma}$
satisfied, for example, by closed hyperbolic 3-manifolds of small diameter, we
prove that the data evolves to a future complete spacetime in the expanding
direction which asymptotes to a vacuum Friedman solution with $(M,\gamma)$ as
the expansion normalized spatial geometry. In particular, the Strong Cosmic
Censorship conjecture holds for this class of solutions in the $C^{2}$-sense.

###### Key words and phrases:

Einstein scalar-field system, stability, blow-up profile, cosmic censorship,
Big Bang singularity

##  1\. Introduction

###  1.1. Setting and main results

We consider the Einstein scalar-field system

(1.1a)
$\displaystyle\text{{Ric}}[\overline{g}]_{\mu\nu}-\frac{1}{2}R[\overline{g}]\overline{g}_{\mu\nu}=$
$\displaystyle\,8\pi T_{\mu\nu}[\overline{g},\phi]$ (1.1b) $\displaystyle
T_{\mu\nu}=$
$\displaystyle\,\color[rgb]{0,0,0}\overline{\nabla}_{\mu}\phi\color[rgb]{0,0,0}\overline{\nabla}_{\nu}\phi-\frac{1}{2}\overline{g}_{\mu\nu}\overline{\nabla}^{\alpha}\phi\overline{\nabla}_{\alpha}\phi$
(1.1c) $\displaystyle\square_{\overline{g}}\phi=$ $\displaystyle\,0\,$

with initial data $(g_{0},k_{0},\pi_{0},\psi_{0})$ on a closed 3-manifold $M$
that admits a negative Einstein metric $\gamma$.111Here and throughout,
$\pi_{0}$ and $\psi_{0}$ prescribe data for $\nabla\phi|_{\Sigma_{t_{0}}}$ and
$\partial_{0}\phi|_{\Sigma_{t_{0}}}$ respectively. In this paper, we determine
the maximal globally hyperbolic development emanating from such initial data
given that it is sufficiently close to the initial data of a homogeneous
solution with a non-trivial scalar field.
In the collapsing direction, we prove a stable Big Bang formation and
curvature blow-up result, which requires the presence of a non-trivial scalar
field. The results complement those in [RS18b, Spe18], which cover flat and
spherical spatial geometry. In the expanding direction, we prove a nonlinear
future stability result of the corresponding vacuum background solution, which
is the Milne model, under a mild condition for the first positive eigenvalue
of $\Delta_{\gamma}$ (see Definition 9.2). As discussed in more details in
Remark 9.3, numerical studies (see [CS99, Ino01]) show that this condition
holds for an analogue of Weeks space, and suggest that this may hold for all
closed hyperbolic $3$-manifolds with sectional curvature $-\frac{1}{9}$ .

Connecting the two regions, we prove the global stability (i.e., past and
future stability) of the spacetime

(1.2a) $\left([0,\infty)\times M,-dt^{2}+a(t)^{2}\gamma\right),$ given a
negative Einstein manifold $(M,\gamma)$ obeying the aforementioned spectral
condition, with (1.2b) $a(0)=0,\
\dot{a}=\sqrt{\frac{1}{9}+\frac{4\pi}{3}C^{2}a^{-4}}$ for some given constant
$C>0$, and the scalar field given by (1.2c) $\partial_{t}\phi=Ca^{-3},\
\nabla\phi=0\,.$ The scale factor consequently exhibits the following
asymptotic behaviour: (1.2d) $a(t)\simeq t^{\frac{1}{3}}\mbox{ as }t\searrow
0\mbox{ and }a(t)\simeq t\mbox{ as }t\nearrow\infty$

The main result can be split into two parts:

######  Theorem 1.1 (Big Bang stability – rough version).

Let $(M,g_{0},k_{0},\pi_{0},\psi_{0})$ be initial data for the Einstein
scalar-field system that is sufficiently close to
$(M,a(t_{0})^{2}\gamma,-\dot{a}(t_{0})a(t_{0})\gamma,0,Ca(t_{0})^{-3})$, where
$C>0$ and $(M,\gamma)$ is a closed Riemannian 3-manifold with
$\text{{Ric}}[\gamma]=-\frac{2}{9}\gamma$ (i.e., a closed negative Einstein
manifold with scalar curvature $-\frac{2}{3}$).

Then, the past maximal globally hyperbolic development $((0,t_{0}]\times
M,\overline{g},\phi)$ of the initial data within the Einstein scalar-field
system (1.1a)-(1.1c) admits a foliation by CMC hypersurfaces
$\Sigma_{s}=t^{-1}(\\{s\\})$ with zero shift. This development remains close
to the FLRW solution described in (1.2a)-(1.2c) in the past of the initial
data slice $\Sigma_{t_{0}}$. In particular, the solution exhibits curvature
blow-up of order $t^{-4}$ and every causal geodesic becomes incomplete as $t$
approaches $0$.

######  Theorem 1.2 (Global stability).

Let $(M,g_{0},k_{0},\pi_{0},\psi_{0})$ be initial data as in Theorem 1.1. In
addition, we suppose that the smallest positive eigenvalue of
$-\Delta_{\gamma}$ acting on scalar functions is strictly greater than
$\frac{1}{9}$.
Then, the initial data admits a maximal globally hyperbolic development
$((0,\infty)\times M,\overline{g},\phi)$ solving the Einstein scalar-field
system that, in addition to the results of Theorem 1.1, is future (causally)
complete. As $t\nearrow\infty$, the solution is attracted by Milne spacetime
in the sense that the expansion normalized variables
$(\bm{g},\bm{k},\nabla\phi,\phi^{\prime})$ (see Definition 9.4) converge
toward $(\gamma,\frac{1}{3}\gamma,0,0)$.

A more detailed statement of Theorem 1.1 is provided in Theorem 8.2. The
additional spectral condition in Theorem 1.2 is discussed at the end of
Subsection 1.3, and the statement itself is proven in Section 10 to be an
extension of the Milne future stability result in Theorem 9.1.

###  1.2. Background material

We now provide context for the previously discussed setting and the results in
Theorems 1.1-1.2:

####  1.2.1. Initial data to the Einstein scalar-field equations

It is well known that the Einstein equations can, via the the 3+1
decomposition, be viewed as an elliptic-hyperbolic system of PDEs (see, for
example, [AM03]). This reduces solving the Einstein equations to two problems:
finding admissible Einstein initial data in physical space, and then solving
the corresponding initial value problem. Regarding the former, initial data to
the Einstein scalar-field system takes the form

$(M,\mathring{g},\mathring{k},\mathring{\pi},\mathring{\psi}),$

where $\mathring{g}$ and $\mathring{k}$ are symmetric $(0,2)$-tensors on $M$,
$\mathring{\pi}$ is a $(0,1)$-tensor (corresponding to $\nabla\phi$) and
$\mathring{\psi}$ is a scalar function (corresponding to the future directed
normal derivative $\partial_{0}\phi$ of the scalar field). The initial data
must satisfy the Hamiltonian and momentum constraints

(1.3a)
$\displaystyle\text{R}[\mathring{g}]+\left({\mathring{k}^{a}}_{a}\right)^{2}-\left({\mathring{k}^{a}}_{b}{\mathring{k}^{b}}_{a}\right)=$
$\displaystyle\,8\pi\left[\lvert\mathring{\psi}\rvert^{2}+\lvert\mathring{\pi}\rvert_{\mathring{g}}^{2}\right]\,,$
(1.3b) $\displaystyle\text{{div}}_{\mathring{g}}\mathring{k}=$
$\displaystyle\,-8\pi\cdot\mathring{\pi}\cdot\mathring{\psi}\,$

(see (2.16a) and (2.16b)), where the indices of $\mathring{k}$ in the first
line are raised with respect to $\mathring{g}$.
We note that, in our argument, we will additionally assume that our initial
data has constant mean curvature so that our gauges can be satisfied initially
– this is enforced on the level of initial data by requiring

$tr_{\mathring{g}}\mathring{k}=-1\,$

which, after embedding, will lead to
$\text{tr}_{\overline{g}|_{\Sigma_{t}}}k=\tau=-3\frac{\dot{a}}{a}$ (see
(2.10)). We will argue in Remark 8.1 why the initial data being near-FLRW
allows us to assume the initial hypersurface to be CMC without loss of
generality.

The results of [FB52, CBG69] show that there exists an embedding222We usually
ignore the embedding in notation.
$\iota:M\hookrightarrow\iota(M)\subset\overline{M}$ and a solution
$(\overline{M},\overline{g},\nabla\phi,\partial_{0}\phi)$ to the Einstein
scalar-field equations such that $\iota(M)=\Sigma_{t_{0}}$ is a Cauchy
hypersurface and such that

$\iota^{\ast}\overline{g}=\mathring{g},\,\iota^{\ast}k=\mathring{k},\,\iota^{\ast}\pi=\mathring{\pi},\,\text{and
}\iota^{\ast}\partial_{0}\phi=\mathring{\psi_{0}}\,.$

We will perturb around initial data $(M,\gamma,-\gamma,0,C)$ which, after
embedding into FLRW spacetime, corresponds to data

$(M\cong\Sigma_{t_{0}},a(t_{0})^{2}\gamma,-\dot{a}(t_{0})a(t_{0})\gamma,0,Ca(t_{0})^{-3})$

at initial time $t_{0}$. Furthermore, the maximal globally hyperbolic
development (MGHD) is unique (up to diffeomorphism), and thus we can assume
$(\overline{M},\overline{g},\nabla\phi,\partial_{0}\phi)$ to be globally
hyperbolic. However, these statements provide little information on the
properties of the MGHD in the future and past of the initial data slice.

####  1.2.2. Strong Cosmic Censorship

In their groundbreaking papers [Haw67, Pen65] on singularity theorems, Hawking
and Penrose established very general criteria for the MGHD of spacetimes to
become causally geodesically incomplete. Many spacetimes of physical relevance
satisfy these criteria, including the spacetimes considered in this article.
While giving us more information on the MGHD than the existence and uniqueness
results mentioned above, a key issue in the application of this mathematical
result to General Relativity is that no statement is made on how precisely the
singularity comes about: In particular, such incompleteness (within a given
regularity class) could either mean that the geodesic is inextendible – which
must be caused by the blow-up of some geometric quantity – or that there exist
multiple inequivalent extensions. While the latter behaviour is exhibited even
for some cosmological spacetimes (see, for example, the Taub solutions
discussed in [CI93]), such behaviour is usually considered to be unphysical
since it would imply a breakdown of determinism. The Strong Cosmic Censorship
Conjecture (SCCC) posits in its most general form that, for generic solutions
to the Einstein equations, this incompleteness instead manifests as
inextendibility at a given level of regularity (e.g.,
$C^{0},C^{2},C^{\infty},\dots$).

In certain frameworks in the homogeneous cosmological setting – i.e. for
homogeneous initial data on a closed spatial hypersurface –, it was shown in
fundamental works by Chrusciel-Rendall [CR95] and Ringström [Rin09] that the
so called Kretschmann scalar
$R_{\alpha\beta\gamma\delta}R^{\alpha\beta\gamma\delta}$ is unbounded where
incompleteness manifests. Thus, it is the driving force behind geodesic
incompleteness in these cases, forcing $C^{2}$-inextendibility of the MGHD.
For the purposes of analyzing cosmologically relevant spacetimes, the SCCC is
hence often rephrased as follows:

######  Conjecture 1.3 (Cosmological SCCC).

(See e.g. [Rin09, Chapter 17]) For generic initial data, the Kretschmann
scalar is unbounded where causal geodesics become incomplete.

Theorem 1.1, in short, shows that this conjecture is rigorously supported in
the case of FLRW spacetimes with negative spatial curvature. More precisely,
the past asymptotics of such spacetimes, determined by initial data on
$\Sigma_{t_{0}}$ as discussed above, are generic in the following sense: There
exists an open neighbourhood of said FLRW data within the set of Einstein
scalar-field initial data such that the solutions past direct causal geodesics
become incomplete, and the incompleteness is driven by blow-up of Kretschmann
scalar with the same asymptotics as the FLRW solution. The global result in
Theorem 1.2 portrays the other side of Cosmic Censorship – as with the past
evolution, near-FLRW data fully determines the future of the spacetime in the
sense that the MGHD is future complete, again showing that this feature of
FLRW spacetimes with negative spatial sectional curvature is generic.

####  1.2.3. FLRW and generalized Kasner spacetimes with scalar fields

On a large scale, the universe is often viewed as spatially homogenous and
isotropic, i.e., no point in space and no direction are distinguishable from
any other point and direction. In 1935, it was shown by Robertson and Walker
that, under a few very natural additional assumptions, this restricts the
class of potential spacetimes to the FLRW class

$\left(I\times\tilde{M},\,\tilde{g}_{FLRW}=-dt^{2}+a(t)^{2}\tilde{\gamma}\right)\,,$

where $(\tilde{M},\tilde{\gamma})$ is a manifold of constant sectional
curvature $\kappa$ and where the scale factor $a$ depends smoothly on $t$.
This holds before taking the Einstein equations into consideration – when
doing so, the matter model determines how space expands within the
cosmological model via $a$. We refer to Lemma 2.3 for the scalar-field
solution for $\kappa=-\frac{1}{9}$, but note that the scale factor behaves
like $t^{\frac{1}{3}}$ for scalar-field matter, regardless of spatial
geometry, and that the Kretschmann scalar blows up at order
$\mathcal{O}\left(t^{-4}\right)$ toward the Big Bang ($t\downarrow 0$).
Spatially flat FLRW spacetimes are a subclass of the closely related
generalized Kasner spacetimes, which are still spatially homogeneous but
anisotropic in general. For scalar field matter, the spacetime metric is given
by

$\displaystyle\overline{g}_{Kasner}=-dt^{2}+\sum_{i=1}^{D}t^{2p_{i}}dx^{i}\otimes
dx^{i},\quad\sum_{i=1}^{D}{p_{i}}=1,\,\sum_{i=1}^{D}p_{i}^{2}=1-A^{2},\quad\overline{\phi}_{Kasner}(t)=A\log(t)\,.$

The standard Kasner family is obtained by considering the vacuum case ($A=0$),
and the spatially flat FLRW spacetime by setting
$D=3,p_{i}=\frac{1}{3},A=\sqrt{\frac{2}{3}}$. If more than one of the Kasner
exponents is non-zero, the generalized Kasner family satisfies the SCCC, also
by exhibiting Kretschmann scalar blow-up of order $t^{-4}$as $t\downarrow 0$
(see [RS18b, (1.8)]).

Kasner spacetimes are of particular relevance to cosmology due to their
relationship with the BKL conjecture: Heuristically, this conjecture states
that the dynamics of cosmological spacetimes generically exhibit chaotic and
highly oscillatory behaviour, often referred to as “Mixmaster” behaviour. This
behaviour is driven by velocity terms within the Einstein equations and is
locally comparable to that of (vacuum) Kasner solutions. However, even if the
BKL picture is to be believed in general, scalar-field (or, more generally,
stiff-fluid) solutions seem to form an exception to it: They have a dampening
effect on said oscillations, thus generating Big Bang stability as shown
rigorously in [RS18b, FRS23] for Kasner spacetimes (for more details, see
Section 1.3). This scenario, often referred to as quiescient cosmology, was
studied in, for example, in [BK73, Bar78, AR01]. With this in mind, both the
aforementioned Kasner results and the results within this article, along with
the prior FLRW results [RS18b, Spe18], confirm this quiescent effect of scalar
fields in cosmology, in contrast to the BKL picture.

We note that one can view this as a scalar field ensuring a specific scenario
in the very early universe given a class of initial data, namely matching the
asymptotic behaviour of the Big Bang singularity. This fits into the recent
use of nonlinear scalar fields in string cosmology, where specific choices of
field are made to specific behaviours (e.g., inflation) in the early universe.
For a recent review, we refer to [CCM+23].

###  1.3. Relation to previous work

Theorem 1.2 is the first theorem about the full global structure of FLRW
spacetimes with negatively curved spatial geometry. For such solutions,
existing results exclusively concern future stability, which we further
discuss below. Besides [Spe18] covering the $\mathbb{S}^{3}$-case, it is the
only open set of initial data for cosmological spacetimes (i.e., without
symmetry assumptions) with $\Lambda=0$ and in absence of accelerated expansion
for which the global (future and past) dynamics are now fully
understood.333For a related future stability result in accelerated expansion,
see [Rin08] which considers scalar fields with a non-trivial potential.

Scalar field matter (and, more generally, semilinear wave equations or fluid
matter) and their asymptotic behaviour on fixed cosmological backgrounds have
been studied extensively, for example in [AR10, AFF19, Bac19, Rin19, BO20,
Rin20, Rin21, Wan21]. While many of the results, in particular [Rin20], manage
to analyze very general classes of equations and spacetime geometries,
including the wave equation on the FLRW backgrounds studied in [AFF19, FU22],
the methods used are often difficult to apply to the full Einstein scalar-
field system. In [FU22], we extended the approach of [AFF19] to be able to
deal with various warped product spacetimes, and in particular FLRW spacetimes
with negatively curved spatial geometry, by using the spatial Laplace operator
to control high order derivatives. The perturbation-adapted analogue of this
strategy is at the basis of the energy method in this paper.

We also note that, by the results of [GaNS19], there are non-trivial waves on
fixed FLRW backgrounds that converge toward the Big Bang singularity, even if,
as demonstrated in [AFF19, FU22], this behaviour is non-generic. Such waves
can give rise to convergent asymptotics on cosmological backgrounds as studied
in [Rin20]. Thus, it will likely be difficult to replace (1.2c) with an
arbitrary non-trivial reference wave while keeping past stability intact.
However, by restricting to an open neighbourhood near the solution described
in (1.2a)-(1.2c), we exclude this non-generic class of solutions.

Theorem 1.1 forms the counterpart to the pioneering works by Rodnianski-Speck
[RS18a, RS18b] and Speck [Spe18], which cover nonlinear Big Bang stability for
FLRW spacetimes with spatial geometry $\mathbb{T}^{3}$ and $\mathbb{S}^{3}$
respectively. These results were extended to Kasner spacetimes in [RS22] with
$\lvert q_{i}\rvert<\frac{1}{6}$, and to the full subcritical regime in
[FRS23], i.e., (generalized) Kasner spacetimes as discussed in Section 1.2.3
with $\max_{i,j,k=1,\dots,D}p_{i}+p_{j}-p_{k}<1$. The former necessitates
considering $1+D$-dimensional Kasner spacetimes with $D\geq 38$, while the
latter result also can be satisfied in $D=3$ for generalized Kasner
spacetimes. Recall that this means, in contrast to our setting, that the
reference spacetime can be anisotropic, even if it the conditions on Kasner
exponents rule out extremely anisotropic regimes. As a result, the analysis
therein becomes significantly more involved, especially at top order, since
approximately monotonic energy identities as used in our work as well as in
[RS18b, Spe18] have not been found in these anisotropic settings.
We note that the argument in [FRS23] relies on identifying an almost-diagonal
structure for the asymptotics of (combined) connection coefficients for an
adapted frame that is carried along by Fermi-Walker transport; this is
precisely where subcriticality enters. Given that these no longer can vanish
in a reference frame adapted to near-hyperbolic spatial geometry, it is
unclear whether this structure is sufficiently maintained.

Furthermore, Beyer and Oliynyk have recently shown in [BO21] that over
$\mathbb{T}^{3}$, the Big Bang formation can be localized in the sense that
data given solely on a ball within the initial hypersurface must also cause
stable blow-up on a (smaller) ball on the Big Bang hypersurface. While this
result further indicates that blow-up behaviour of near-FLRW spacetimes might
be, at least, independent of global geometric properties as it seems to be a
localizeable, we note that proof of localized stability crucially relies on
the flatness of the conformal reference spacetime. To be more precise, the
proof relies on extending the local initial data to global data for a Fuchsian
system of metric and matter quantities as well as, again, connection
coefficients for an adapted, Fermi-Walker transported frame. However, the
derivation of the system for the former explicitly seems to use flat spatial
geometry to obtain the necessary Fuchsian form. This forms seems to similarly
be broken as soon as the connection coefficients are not perturbed around $0$,
since this would lead to inhomogeneous error terms of order $t^{-1}$ for the
rescaled variables which are stronger than what the method, so far, accounts
for.

By contrast, in [RS18b, Spe18], the reference frame itself is used in the
commutator method to obtain the necessary energy identities at high orders. In
all of these works, it hence is a priori unclear how one could extend these
methods to the negative spatial Einstein geometry of $(M,\gamma)$. We provide
an alternative approach that, besides establishing the complementary stability
result to [RS18b, Spe18], does not rely on any information on the spatial
geometry of the reference manifold in its methodology (although it is of
course relevant in determing the FLRW reference solution that we are
studying). Instead, we rely on differential operators adapted to the evolved
spatial metric. Hence, we believe that our approach may also prove useful for
stability problems in spatially inhomogeneous (and hence also anisotropic)
settings. In light of [RS22, FRS23] in particular, the main challenge in
achieving this would either be to find approximately monotonic energy
identities with our Bel-Robinson approach that have not been observed
previously, or to also find ways to circumvent the lack thereof.

To obtain Theorem 1.1, we use the Laplace-Beltrami-operator (acting,
respectively, on scalar functions and tensor fields) with respect to the
(rescaled) evolved metric as our commutating operator instead of a fixed
reference frame. This, in turn, leads us to replacing the wave-like system for
metric and second fundamental form exploited in [RS18b, Spe18] by an
evolutionary system in the second fundamental form and Bel-Robinson variables.
The latter technique to show stability of near-vacuum solutions dates back to
[AM04] and even to the pioneering proof of global stability of Minkowski space
in [CK93], which both cover the Einstein vacuum system, as well as to [CK90],
where Bel-Robinson variables were introduced to analyse field equations on
Minkowski space. It has recently also been applied to the future stability of
Milne spacetimes in the massive Einstein Klein-Gordon system by Wang in
[Wan19]. As far as we are aware, this method has not yet been applied to
solutions that are not near-vacuum or in the context of Big Bang singularity
formation.

Toward the Big Bang, the solutions exhibit asymptotically velocity dominated
(AVTD) behaviour in the sense that they behave, to leading order, like
solutions to the Einstein scalar-field equations in CMC gauge with zero shift
with all terms involving spatial derivatives set to zero (the “velocity term
dominated” (VTD) equations). This behaviour also matches results obtained by
studying high regularity solutions (e.g., [AR01]), or related works using
Fuchsian methods that prescribe a behaviour at the singularity and then
develop it locally, often under additional symmetry assumptions (e.g.,
[DHRW02, CBIM04, IM02, FL20]). In particular, this asymptotic behaviour leads
to the same types of “Kasner footprint states” as in [RS18a, RS18b]: As one
approaches the Big Bang, the rescaled variables converge toward tensor fields
on the Big Bang hypersurface that precisely solve the truncated VTD equations.
Further, the distance between the footprints of the FLRW and the perturbed
solution are controlled by the initial data. For example, the rescaled
Weingarten map $a^{3}{k^{a}}_{b}$ converges to ${(K_{Bang})^{a}}_{b}$ on the
Big Bang hypersurface, which is close to
$\frac{\sqrt{4\pi}3}{C}\mathbb{I}^{a}_{b}$, the rescaled FLRW footprint (see
(8.3e) and (8.5c)).

What remains to be considered to obtain Theorem 1.2 is future stability, which
we can reduce to future stability of the vacuum solution in the Einstein
scalar-field system. This solution, called the Milne spacetime, has been shown
to be stable within the set of vacuum solutions – see [AM11] – and a range of
other Einstein systems – see, for example , [Wan19, AF20, FW21, FOW21, BF20,
BFK19] and related work in lower dimensions, e.g. [AMT97, Mon08, Faj17, Faj20,
Mon20]. As such, our contribution to the study of future stability of Milne
spacetimes is that we deal with the massless scalar field matter via corrected
energy estimates which are inspired by work of Choquet-Bruhat and Moncrief in
[CBM01] for vacuum Einstein equations with $U(1)$-symmetry. Out of the works
listed above, only [Wan19, FW21] deal with scalar field matter at all, namely
the massive case. These fields exhibit stronger decay toward the future,
making the matter components easier to deal with than in our analysis.
The additional spectral condition is needed to ensure coercivity of the
corrected scalar field energy. Numerical work, e.g. [CS99, Ino01], does not
suggest that this condition is violated by any closed 3-manifold with constant
sectional curvature $\kappa=-\frac{1}{9}$, and verifies that is is satisfied,
for example, by an analogue of Weeks space in which the metric is
appropriately scaled to have the required sectional curvature. We refer to
Remark 9.3 where this discussed in more detail.

###  1.4. Challenges in the proof

The contracting and expanding regimes of near-FLRW spacetime are analyzed in
two separate and methodologically independent parts. Before providing an
overview of both arguments, we summarize the challenges that arise:

####  1.4.1. Big Bang stability

The main difficulties in establishing Big Bang stability are three-fold:

Firstly, we have to expect that the solutions are asymptotically velocity term
dominated (as argued in Remark 8.3, we end up proving that this is the case),
and thus that rescaled variables at best exhibit the same asymptotic behaviour
as their counterparts in FLRW spacetime, up to a small perturbation in the
asymptotic footprint. For example, note that, in the reference FLRW spacetime,
one has

${(k_{FLRW})^{i}}_{j}=-3\frac{\dot{a}}{a}\mathbb{I}^{i}_{j}\approx-\frac{1}{t}\mathbb{I}^{i}_{j}\,.$

At best, the shear ${\hat{k}^{i}_{j}}$ of the perturbed solution then behaves
like $\frac{\varepsilon}{t}$. In fact, we show that this is the case in
(4.2b). This implies that the contraction rescaled metric
$G_{ij}=a^{-2}g_{ij}$ can only be controlled up to
$\mathcal{O}\left(t^{-c\sqrt{\varepsilon}}\right)$ (see (4.4c)), since one has
$\partial_{t}g_{ij}\approx-2g_{il}{k^{l}}_{j}$ and thus

$\partial_{t}G_{ij}\approx
G_{il}{\hat{k}^{i}}_{j}\approx\frac{\varepsilon}{t}\ast G\,.$

However, to be able to use the structure of the evolution equations to cancel
terms in our energy arguments, we have to work with adapted quantities. For
example, we need to use integration by parts with respect to
$(\Sigma_{t},G_{t})$ to cancel high order scalar field terms with help of the
(rescaled) wave equation that contains $\Delta_{G}$, or to obtain elliptic
estimates from the lapse equation via the operator $\Delta_{G}$ or from the
adapted div-curl-system for $\Sigma$ arising from the constraint equations.

As a result, even the rescaled solution variables will diverge at order
$\mathcal{O}\left(t^{-c\sqrt{\varepsilon}}\right)$ toward the singularity, so
we need to track and control their rate of divergence within the bootstrap
argument. This significantly complicates dealing with nonlinear terms, where
the bootstrap assumptions often cannot be inserted naively. This in turn makes
coercivity of the energies more involved to establish (see Lemma 4.5 and
Remark 4.6), since this only holds up to curvature errors that also diverge
and thus need to be carefully tracked.

Secondly, and in contrast to [RS18b, Spe18], replacing the wave structure of
the geometric evolution in the Einstein equations with our less geometry
dependent Bel-Robinson framework seems to lose regularity at first glance: The
evolution system for the scalar field energy and the geometric energies can be
caricatured as follows, where the superscript for each energy term denotes the
number of derivatives:

$\displaystyle-\frac{d}{dt}\mathcal{E}^{(L)}(\phi,\cdot)\lesssim$
$\displaystyle\,\frac{\varepsilon^{\frac{1}{8}}}{t}\left[\mathcal{E}^{(L)}(\phi,\cdot)+\mathcal{E}^{(L)}(\Sigma,\cdot)\right]+\dots$
$\displaystyle-\frac{d}{dt}\left[\mathcal{E}^{(L)}(\Sigma,\cdot)+\mathcal{E}^{(L)}(W,\cdot)\right]+\dots\lesssim$
$\displaystyle\,\frac{\varepsilon^{\frac{1}{8}}}{t}\left[\mathcal{E}^{(L)}(\Sigma,\cdot)+\mathcal{E}^{(L)}(W,\cdot)\right]+\frac{\varepsilon^{-\frac{1}{8}}}{t}\cdot
a^{4}\mathcal{E}^{(L+1)}(\phi,\cdot)+\dots$

Thus, it seems that we lose derivatives in the scalar field and aren’t able to
close the argument. This is remedied using the div-curl-system in $\Sigma$,
see (2.36a) and (2.36b), which yields a weak estimate of the form

$a^{4}\mathcal{E}^{(L+1)}(\Sigma,\cdot)\lesssim\mathcal{E}^{(L)}{(\phi,\cdot)}+\mathcal{E}^{(L)}(W,\cdot)+\mathcal{E}^{(L)}(\Sigma,\cdot)+\dots\,.$

Combining these estimates to improve the bootstrap assumptions then
necessitates an intricately constructed total energy to balance these
different types of estimates against one another.

Finally, given (1.2c), the rescaled time derivative of the scalar field is not
small and does not become so toward the Big Bang. This leads to various terms
within the core linearized evolutionary system of both matter and geometry
that, if estimated naively, could lead to exponential blow-up toward the
singularity. When such terms occur in the scalar field energy evolution, this
can be dealt with along similar lines as in [RS18b, Spe18], but we incur
additional large terms in our geometric evolution that only cancel using the
explicit form of the Friedman equations, which we highlight in Lemma 7.1 and
its proof.

####  1.4.2. Future and global stability

For Milne stability, the canonical Sobolev energies for the scalar field
variables, i.e.,

$\int_{M}\lvert\phi^{\prime}\rvert_{\bm{g}}^{2}+\lvert\nabla\phi\rvert_{\bm{g}}^{2}\,{\text{{vol}}}_{\bm{g}}$

and higher order analogues, do not obey useful energy estimates. This can be
overcome by adding an indefinite correction term of the type

$\int_{M}\phi^{\prime}(\phi-\overline{\phi}){\text{{vol}}}_{\bm{g}}$

to the canonical energy, see Definition 9.6. This is similar to what was done
in [CBM01] in a $2+1$-dimensional setting, as well as similar to the
indefinite terms we introduce in our geometric energy to control the wave
system in the metric variables, as in previous work on Milne stability in
different matter models, including [AF20, FW21]. That this corrected energy
controls Sobolev norms relies on the aforementioned spectral condition. As a
result, and unlike for past stability, the specific spatial geometry is
crucial in generating decay from energy estimates, even before considering the
geometric evolution.

Moreover, we need to transition from the near-FLRW data used to analyze the
contracting regime to data in the expanding regime on a distant enough future
hypersurface such that it is near-Milne and the future stability result
applies. This requires as a gauge switch from CMC gauge with zero shift to
CMCSH gauge, as well as careful control of the solution variables over a
finite time interval using continuous dependence on initial data. For the
former, close inspection of [FK20] gives us a diffeomorphism close to the
identity that maps the initial data for the metric to new data satisfying the
spatially harmonic gauge condition, thus allowing us to switch gauges without
losing proximity to the reference solution. This is discussed in detail in
Section 10.

###  1.5. Proof outline

####  1.5.1. Big Bang stability

m

The big picture. The key argument in our Big Bang stability proof is a
hierarchized series of energy estimates that establishes the asymptotic
behaviour of solution variables toward the singularity. We rely on a bootstrap
argument which establishes that energies $\mathcal{E}^{(L)}$ (see Definition
3.9) for the scalar field, the rescaled shear, the Bel-Robinson variables, the
lapse and the curvature at worst only diverge slightly. Here, $0\leq
L\leq\color[rgb]{0,0,0}18\color[rgb]{0,0,0}$ denotes the order of derivatives
considered. To this end, we make a bootstrap assumption on the solution norms
$\mathcal{C}$ (see Definition 3.6) which controls the distance of these
rescaled variables as well as the metric itself to their FLRW counterparts in
terms of supremum norms with respect to $G$, where $G=a^{-2}g$ is the rescaled
adapted spatial metric (see Definition 2.9). We refer to Assumption 3.16 and
Remark 3.19 for the detailed bootstrap assumptions and improvements, as well
as to Lemma 3.14 for the underpinning local well-posedness result. That this
bootstrap argument implies Theorem 1.1 follows from a straightfoward
adaptation of the arguments in [RS18b, Theorem 15.1].

We work with evolution-adapted norms even though $G(t,x)$ degenerates toward
the Big Bang singularity. Indeed, since we need to exploit the structure of
the evolutionary equations, it is more convenient to have these adapted
quantities controlled by the solution norms $\mathcal{H}$ and $\mathcal{C}$
directly instead of having to perform changes of metric at that point. Once
the improved energy estimates are shown, a (time-scaled) coercivity notion
(see Lemma 4.5 and the proof of Corollary 7.3) and Sobolev embeddings with
respect to the reference metric $\gamma$ then ensure that these improved
estimates translate to $\mathcal{H}$ and $\mathcal{C}$. This then closes the
bootstrap. To actually achieve this improved energy behaviour, we derive
elliptic energy estimates or integral-type estimates that, once suitably
combined and scaled, yield the desired improvements by straightforwardly
applying the Gronwall lemma. Additionally, note that we assume that the
initial data is close to FLRW data not just in $\mathcal{H}$, which contains
precisely the norms needed to control $\mathcal{C}$ by Sobolev embedding, but
also scaled smallness assumptions at one order higher, contained in the top
order semi-norm $\mathcal{H}_{top}$ (see Assumption 3.10). This is needed to
ensure that the top order energy is small initially, and thus to close the
bootstrap.

Scale factor a(t). The precise structure of the Friedman equations (2.3)-(2.4)
is crucial not only to control time integral quantities up to the Big Bang
hypersurface (see Lemma 2.4), but also to ensure that certain terms in the
evolution that would otherwise cause large divergences contribute with
favourable sign (see the arguments in Lemma 6.2 as well as Lemma 7.1). It
turns out that the sectional curvature entering the Friedman equations
actually is not of key importance to large parts of the Big Bang stability
analysis given that the leading order of the scale factor toward the Big Bang
singularity is determined via the Friedman equation (1.2b) by the matter term,
not the sectional curvature. This indicates that our method might extend to
different settings.

Gauge choice, commutation method and Bel-Robinson variables. We commute the
resulting elliptic-hyperbolic Einstein system with the Laplace-Beltrami
operator $\Delta_{G}$ with respect to the rescaled evolved spatial metric
$G(t,x)$ to obtain higher order energy control. Commuting with this operator
has the advantage of leaving many integration-by-parts identities intact.
These are needed to provide specific cancellations, e.g., to cancel
$\Delta^{\frac{L}{2}+1}\phi$-terms arising from the wave equation when
computing $\partial_{t}\mathcal{E}^{(L)}(\phi,\cdot)$. We also note that the
only feature of the adapted metric we use is that it is close to $\gamma$, and
don’t use any further information on the geometry, e.g., by choosing a
specific reference frame in our commutation method. Further, we employ CMC
gauge with zero shift to avoid badly behaved shift terms (see Remark 1.4).

We still, however, need to deal with the Ricci term in the evolution equation
for the second fundamental form. To this end, we consider the Bel-Robinson
variables $E$ and $B$ which are $\Sigma_{t}$-tangent symmetric tracefree
$(0,2)$-tensors and contain all information of the spacetime Weyl tensor
$W[\overline{g}]$ (see Subsection 2.4). Suitably projecting the Gauss-Codazzi
equations admits additional constraint equations in terms of $E$ and $B$ that
allow us to replace the Ricci tensor at the “cost” of introducing Bel-Robinson
energies into the formalism, see (2.24a) and the rescaled version (2.29c).
Further, $E$ and $B$ satisfy a Maxwell-type system (see Lemma 2.7) that can be
exploited to obtain energy estimates and, as with the other evolution
equations, is well adapted to commutation with $\Delta_{G}$.

A priori low order $C_{G}$-control. By applying the bootstrap assumptions on
$\mathcal{C}$ to the evolution equations, we can immediately deduce improved
low order estimates in $C_{G}^{l}$ for $l\geq 10$ for the solution variables
by inserting them into the respective evolution equations (see Lemma 4.3), as
well as via the maximum principle for the lapse (see Lemma 4.1). These usually
still diverge slightly, mostly due to the asymptotic behaviour of $G$. However
and crucially to our argument, at order $0$, the renormalized time derivative
$\Psi$ of the wave, the rescaled tracefree part $\Sigma$ of the second
fundamental form and the rescaled Bel-Robinson variable $\bm{E}$ are in fact
$K\varepsilon$-small in $C^{0}_{G}$ on the bootstrap interval (see Lemma 4.2).
If these estimates would not hold, this would lead to terms that diverge at
order $\mathcal{O}\left(a^{-3-c\sqrt{\varepsilon}}\right)$ in the differential
inequalities, and thus causing exponential energy blow-up of order
$\mathcal{O}\left(e^{a^{-c\sqrt{\varepsilon}}}\right)$ that we could no longer
control. This behaviour is closely related to the fact that $\Psi$ and
$\Sigma$ converge toward footprint states on the Big Bang hypersurface that
remain $K\varepsilon$-small (see (8.3c) and (8.3e)), and then pass this
convergence on to $\lvert\bm{E}\rvert_{G}$ (see (8.8a).

Energy estimates and hierarchy. The main part of the analysis is establishing
various energy estimates.

  * •

For the lapse (see Section 5), the relevant estimates are direct results of
the elliptic lapse equations (2.30a)-(2.30b). The non-lapse terms on the right
hand side only diverge slightly toward the Big Bang, in contrast to the
divergence at oder $a^{-4}$ in (2.30a), and thus allows one to show that, at
lower derivative order, the lapse converges to $1$. However, since the right
hand side of (2.30b) contains the scalar curvature of $G$, this estimate loses
derivatives. On the other hand, (2.30a) does not lose derivatives, and the
elliptic nature in fact allows one to lapse energies of order $L+2$ by
energies in $\Sigma$ and the scalar field of oder $L$. This makes it possible
to control the higher order lapse term occuring, for example, in (2.28c),
without losing regularity. Conversely, both of these gains in regularity are
at the cost of losing powers of $a$. In short, (2.30b) is needed to establish
the asymptotic behaviour of the lapse, and (2.30a) to obtain improved energy
bounds as a whole.

  * •

The core matter energy estimate (see Lemma 6.2) relies on delicate
cancellations when computing the time derivative of
$\mathcal{E}^{(L)}(\phi,\cdot)$. While we derive this in a fashion that
differs from the energy flux method used in [Spe18], the necessary
cancellations to arrive at Lemma 6.2 are similar.

  * •

The (rescaled) tracefree component of the second fundamental form $\Sigma$
(see Lemma 6.8) and the (rescaled) Bel-Robinson variables $\bm{E}$ and
$\bm{B}$ (see Lemma 6.6) need to be treated simultaneously to deal with the
leading curvature term in the evolution of the former by inserting a
constraint equation in which $\bm{E}$ occurs as the leading term (see
(2.36d)). However, the matter terms within the evolution of $\bm{E}$ and
$\bm{B}$ contain, firstly, terms where we again need very precise estimates to
show that they do not contribute large $a^{-3}$-divergences, and, secondly,
matter terms that lose one order of derivatives.
This order of regularity can be regained using the momentum constraint
equation (2.36a) and its Bel-Robinson counterpart (2.36b) containing $\bm{B}$,
which leads to a div-curl-system for $\Sigma$ (see Lemma 6.10). This is,
again, at the cost of losing powers of $a$.

  * •

As a result, the core Gronwall argument performed in Proposition 7.2 combines
energies for the matter variables, $\Sigma$ and the Bel-Robinson variables, as
well as energies for $\text{{Ric}}[G]$ and additional scalar field quantities
to control borderline error terms. In particular, the curvature energies are
necessary to handle commutation errors within the energy estimates, and
improved bounds on them need to be obtained to apply the coercivity results in
Lemma 4.5 – else, none of energy improvements would extend to improved Sobolev
norm bounds and the bootstrap argument would not close.
As many of the a priori $C_{G}$-norm estimates add small additional
divergences, it is necessary to perform an induction over derivative orders
within this mechanism to deal with lower order error terms. Since $\Delta_{G}$
is elliptic, it is sufficient to perform this for even orders. Along with
energies at order $L\in 2\mathbb{N}_{0}$, the total energy also includes the
energy controlling $\Sigma$ as well as the scalar field and curvature energies
at order $L+1$, appropriately scaled to account for the degenerate elliptic
estimate for $\Sigma$ from Lemma 6.10. This remedies the derivative loss in
the Bel-Robinson energy and allows one to improve the total energy at each
order until reaching $L=18$, at which point the bootstrap argument can be
closed.

  * •

Note that the metric itself does not enter the core energy mechanism. In fact,
trying to replace control of the Ricci tensor by control of $G$ is likely too
imprecise in dealing with high order curvature errors. Instead, control of
$G-\gamma$ and $\Gamma[G]-\hat{\Gamma}[\gamma]$ is a consequence of a simple
integral energy inequality and the improvements achieved for $\Sigma$ and
matter variables (see Lemma 6.14 and Corollary 7.3). Since we cannot utilize
any additional structure in dealing with the metric, we have to construct our
argument carefully to allow for the metric control to be weaker than what one
gets for the core variables, while still being sufficiently strong to
constitute an improvement and allowing to switch between $H_{G}$ and
$H_{\gamma}$ (and, respectively, $C_{G}$ and $C_{\gamma}$) norms.

We also point to Remark 6.1 for a more detailed sketch of how the integral
inequalities for the core Gronwall argument are structured and how this leads
to the bootstrap improvement for the energies.

####  1.5.2. Future stability and connecting the regions

We follow similar lines as in [AF20, FW21] to prove that near-FLRW spacetimes
in negative spatial geometry are future stable. Since $\partial_{t}\phi$
decays like $a^{-3}\simeq t^{-3}$ in the reference spacetime, the sectional
curvature becomes dominant in the Friedman equations and the scale factor
approaches that of Milne spacetime as $t$ approaches $\infty$. Hence, if one
moves sufficiently far to future, choosing near-FLRW data with a homogeneous
scalar field is equivalent to choosing near-vacuum data. Thus, what we prove
first in Section 9 is future stability of near-Milne spacetimes under the
Einstein scalar-field system. Once this is established, we argue in Section 10
how early near-FLRW initial data evolves to data that is sufficiently close to
Milne for large enough times, which is essentially a consequence of the scale
factor and the (physical) mean curvature approaching that of Milne, up to a
multiplicative constant.

In terms of dealing with geometric and elliptic estimates, we can essentially
carry over the results of [AF20], as was also done in [FW21], by working in
CMCSH gauge and verifying that the matter components are indeed only
perturbative terms within the geometric evolution.

This leaves only the scalar field to be examined. Here, we introduce
corrective terms to the energies (see Definition 9.6) which yield decay
estimates for the corrected scalar field energy (see Lemma 9.16 and Lemma
9.17). That these energies are coercive (see Lemma 9.12 and Lemma 9.13)
requires the aforementioned lower bound for the first positive eigenvalue of
$\Delta_{\gamma}$.

######  Remark 1.4 (Why not use CMCSH gauge to prove Big Bang stability?).

One might consider applying this gauge to Big Bang stability as well since
this is precisely the choice of gauge turning the geometric evolution into a
wave-like system in $(g,k)$, which seems simpler than our chosen approach in
CMC gauge with zero shift. In particular, this would also not rely on any
choice of reference frame, and keep the wave structure of the geometric
evolution intact, unlike when using Bel-Robinson variables. However, the issue
with this approach lies in the shift equation, which would take the following
form for the rescaled shift vector $X=a^{3}\tilde{X}$:

(1.4) $\displaystyle\Delta_{G}X^{l}+\text{{Ric}}[G]^{l}_{m}X^{m}=$
$\displaystyle-2(N+1)(G^{-1})^{im}(G^{-1})^{jn}\Sigma_{ij}\left(\Gamma_{mn}^{l}-\hat{\Gamma}_{mn}^{l}\right)$
$\displaystyle+2(G^{-1})^{im}\nabla_{i}X^{n}\left(\Gamma_{mn}^{l}-\hat{\Gamma}_{mn}^{l}\right)$
$\displaystyle\,+\langle\text{error terms in lapse and matter}\rangle$

As a result, the first term has to be expected to diverge at the same rate as
the metric, i.e., we expect even low order norms of $\tilde{X}$ to behave like
$a^{-3-c\sqrt{\varepsilon}}$ at best up to small prefactors. However,
computing the time derivative of an integral over $\lvert
G-\gamma\rvert_{G}^{2}$ (or derivatives thereof) becomes the integral over the
$(\partial_{t}-\mathcal{L}_{\tilde{X}})$-derivative of this quantity, and
hence we get explicit terms of the form $\mathcal{L}_{\tilde{X}}\gamma$ which
always exist at highest order and diverge worse than $t^{-1}$. In short, the
fact that the metric cannot be expected to converge to a footprint state leads
to leading order terms in the differential energy estimates to carry strongly
divergent pre-factors in CMCSH gauge. This obstructs improvements in a
tentative bootstrap argument.

###  1.6. Paper outline

  * •

Sections 2-8 cover the proof of Big Bang stability:

    * –

In Section 2, we introduce notation and provide the necessary information on
the FLRW background solution as well as the equations relevant to the
subsequent analysis.

    * –

Then, in Section 3, we discuss the solution norms and energies and state the
initial data and bootstrap assumptions.

    * –

In Section 4, improved low order $C_{G}$-norm estimates that follow directly
from the bootstrap assumptions are established, along with additional formulas
and a priori estimates.

    * –

Section 5 concerns the elliptic estimates for the lapse.

    * –

In Section 6, we discuss the energy and Sobolev norm estimates for all other
variables, all of which are integral estimates except for the aforementioned
elliptic estimate for $\Sigma$, as well as a norm bound for $\nabla\phi$ that
is not needed for the energy improvement.

    * –

These are all combined in Section 7 to improve the bootstrap assumptions –
first for the energies, then for $\mathcal{H}$ and finally $\mathcal{C}$.

    * –

In Section 8, we show how this bootstrap argument implies the main Big Bang
stability result (see Theorem 8.2, which is the formal version of Theorem
1.1).

  * •

Section 9 contains the proof of near-Milne future stability.

  * •

In Section 10, we show that this is sufficient for future stability of near-
FLRW spacetimes, proving Theorem 1.2.

  * •

The appendices (Sections 11-12) collect various basic formulas and commutator
expressions as well as error terms and how these can be estimated.

Acknowledgements. D.F. acknowledges support by the Austrian Science Fund
through the project Relativistic Fluids in cosmology (project number P34313).
L.U. acknowledges the support by the START-Project Isoperimetric study of
initial data for the Einstein equations (project number Y963-N35), also by the
Austrian Science Fund. L.U. is a recipient of a DOC Fellowship of the Austrian
Academy of Sciences at the Faculty of Mathematics at the University of Vienna.
L.U. also thanks the German Academic Scholarship Foundation (Studienstiftung
des deutschen Volkes) for their scholarship. The authors thank Klaus Kröncke
and Roman Prosanov for their help in seeking out numerical evidence for the
spectral condition used in Section 9, and Michael Eichmair and the anonymous
referees for their extensive feedback on a previous version of this
manuscript.

##  2\. Big Bang stability: Preliminaries

###  2.1. Notation

####  2.1.1. Foliations

On a spacetime manifold $(\overline{M},\overline{g})$, we assume the existence
of a spacelike Cauchy hypersurface $\Sigma_{t_{0}}$ that is diffeomorphic to
$M$. As argued in Remark 8.1, we can assume without loss of generality that it
has constant mean curvature. We will ultimately show that there exists a time
function $t$ such that the past of $\Sigma_{t_{0}}=t^{-1}(t_{0})$ can be
foliated by $\Sigma_{s}=t^{-1}(s)$ for $s\in(0,t_{0})$, and that where the
solution exists, this is at least possible up to some $T\in(0,t_{0})$. These
constant time surfaces are then also spacelike Cauchy hypersurfaces
diffeomorphic to $M$ and CMC. We will use this notation throughout with little
comment and often simply view $\Sigma_{s}$ as $\\{s\\}\times M$.

####  2.1.2. Metrics

The spacetime metric $\overline{g}$ on $\overline{M}$ takes the general form

$\overline{g}=-n^{2}dt^{2}+g_{ab}dx^{a}dx^{b}$

where $n\equiv n(t,x)$ is the lapse function and $g|_{\Sigma_{t}}\equiv
g|_{\Sigma_{t}}(t,x)$ is a Riemannian metric on $\Sigma_{t}$. We will often
simply denote the spatial metric by $g$. Furthermore, we denote the rescaled
spatial metric $byG_{ij}=a^{-2}g_{ij}$ (see Definition 2.9) and the tensor-
field induced by the matrix inverse of $(G_{ij})$ by $G^{-1}$. Similarly,
$\det g$ and $\det G$ are also meant as the determinants in the matrix sense.
Finally, we define ${\text{{vol}}}_{g}$ and $\mu_{g}$ as the volume form and
volume element with regards to $g$, and the same for $\gamma$ and $G$.

####  2.1.3. Indices and coordinates

Greek indices $\alpha,\beta,\dots,\mu,\nu,\dots$ run from $0$ to $3$,
lowercase latin indices $a,b,\dots$, $i,j,\dots$ from $1$ to $3$. The spatial
indices on some coordinate neighbourhood $V\subseteq\overline{M}$ are always
with regards to the local frame induced by coordinates $(x^{1},x^{2},x^{3})$
on $M$, applied to each $V\cap\Sigma_{t}$ by the standard embedding where this
intersection is non-empty. The index $0$ always denotes components relative to
$\partial_{0}=n^{-1}\partial_{t}$, where $\partial_{t}$ is the derivative
associate to the time function $t$. The Levi-Civita connections associated to
$\overline{g}$, respectively $g$ and $G$, are denoted by $\overline{\nabla}$,
respectively $\nabla$.444Note that $g$ and $G$ have the same Levi-Civita-
connection since, on every hypersurface $\Sigma_{t}$, they are related by a
scalar multiple. Additionally, for the hyperbolic spatial reference metric
$\gamma$ on $M$ (see Definition 2.1), we write the Levi-Civita connection as
$\hat{\nabla}$.

Whenever we raise or lower Greek (resp. Latin) indices without additional
notation, it is with regards to $\overline{g}$ (resp. $g$). When we raise
indices of a tensor $\mathfrak{T}$ with regards to the rescaled spatial metric
$G$, we flag this by writing $\mathfrak{T}^{\sharp}$. We never raise or lower
with respect to $\gamma$. Refer to Section 2.1.9 as to how we distinguish
taking multiple covariant derivatives from index raising.

####  2.1.4. $\Sigma_{t}$-tangent tensors

For any $\Sigma_{t}$-tangent tensor
${\xi^{\alpha_{1}\dots\alpha_{r}}}_{\beta_{1}\dots\beta_{s}}$, we write
${{\xi(t)}^{a_{1}\dots a_{r}}}_{b_{1}\dots b_{r}}$ for the
$\overline{g}$-orthogonal projection of $\xi$ onto the hypersurface
$\Sigma_{t}$. When clear from context, we will drop the time dependency in
notation.

####  2.1.5. Sign conventions

Within this paper, the second fundamental form with regards to $\Sigma_{t}$ is
defined as the $(0,2)$-tensor $k$ given by

$k(X,Y)=-\overline{g}(\overline{\nabla}_{X}\partial_{0},Y)\,,$

where $X$ and $Y$ are $\Sigma_{t}$-tangent vectors. The Riemann curvature
tensor of $\overline{g}$ is taken to be

$\overline{\nabla}_{\alpha}\overline{\nabla}_{\beta}Z_{\gamma}-\overline{\nabla}_{\beta}\overline{\nabla}_{\alpha}Z_{\gamma}={\text{{Riem}}[\overline{g}]_{\alpha\beta\gamma}}^{\delta}Z_{\delta}$

for the covariant vector field $(Z_{\mu})$, and the analogous convention holds
for all other Riemann curvature tensors that appear.

####  2.1.6. Constants

For two nonnegative scalar functions $\zeta_{1},\zeta_{2}$, we write
$\zeta_{1}\lesssim\zeta_{2}$ if and only if there exists a constant $K>0$ such
that $\zeta_{1}\leq K\zeta_{2}$. This implicit constant may depend on
information from the FLRW reference solution at the starting point of the
evolution (in particular on $\gamma$ and $a(t_{0})$, see Definition 2.1) and
combinatorial quantities. We extend this notation to a real function
$\zeta_{1}^{\prime}$ by

$\zeta_{1}^{\prime}\lesssim\zeta_{2}:\Leftrightarrow\max\left(\zeta_{1}^{\prime},0\right)\lesssim\zeta_{2}\,.$

Additionally, we write $\zeta_{1}\simeq\zeta_{2}$ if and only if
$\zeta_{1}\lesssim\zeta_{2}\lesssim\zeta_{1}$ is satisfied.

####  2.1.7. Tensor contractions

We denote by $\bm{\varepsilon}_{\alpha\beta\gamma\delta}$ the Levi Civita
tensor with regards to $\overline{g}$ and define Levi Civita tensor on spatial
hypersurfaces $\Sigma_{t}$ by
$\bm{\varepsilon}[g]_{ijk}=\bm{\varepsilon}_{0ijk}$. Notice that this
corresponds to the Levi Civita tensor associated to $g$. Further,
$\bm{\varepsilon}[G]_{ijk}=a^{-3}\bm{\varepsilon}[g]_{ijk}$ is the Levi Civita
tensor with respect to the rescaled metric $G$ (see (2.27a)).

For $\Sigma_{t}$-tangent $(0,2)$-tensors $A,\tilde{A}$ and vector field $v$,
we define the following objects as in [AM04, Section A.2]:

$\displaystyle A\cdot\tilde{A}$
$\displaystyle=\color[rgb]{0,0,0}{A^{a}}_{b}{\tilde{A}_{a}}^{b}\color[rgb]{0,0,0}=\langle
A,\tilde{A}\rangle_{g}$ $\displaystyle(A\odot_{g}\tilde{A})_{ij}$
$\displaystyle=\color[rgb]{0,0,0}A_{ik}{\tilde{A}^{k}}_{\
j}\color[rgb]{0,0,0}$ $\displaystyle(A\wedge\tilde{A})_{i}$
$\displaystyle={\bm{\varepsilon}_{i}}^{jp}{A_{j}}^{q}\tilde{A}_{qp}$
$\displaystyle(v\wedge A)_{ab}$
$\displaystyle={\bm{\varepsilon}_{a}}^{cd}v_{c}A_{db}+{\bm{\varepsilon}_{b}}^{cd}v_{c}A_{ad}$
$\displaystyle(A\times\tilde{A})_{ij}$
$\displaystyle={\bm{\varepsilon}_{i}}^{ab}{\bm{\varepsilon}_{j}}^{pq}A_{ap}\tilde{A}_{bq}+\frac{1}{3}(A\cdot\tilde{A})g_{ij}-\frac{1}{3}(tr_{g}A\cdot
tr_{g}\tilde{A})g_{ij}$
$\displaystyle(\text{{curl}}A)_{ij}=(\text{{curl}}_{g}A)_{ij}$
$\displaystyle=\frac{1}{2}\left[{\bm{\varepsilon}_{i}}^{cd}\nabla_{d}A_{cj}+{\bm{\varepsilon}_{j}}^{cd}\nabla_{d}A_{ci}\right]$
$\displaystyle(\text{{div}}_{g}A)_{i}$ $\displaystyle=\nabla^{b}A_{ib}$

The operations $\odot_{G}$, $\langle\cdot,\cdot\rangle_{G}$ and
$\text{{div}}_{G}$ are defined analogously, with all indices raised and
lowered by $G$ instead of $g$. Finally, for two $(0,1)$-tensors
$\pi,\tilde{\pi}$, we denote their symmetrized product by

$(\pi\otimes\tilde{\pi})_{ij}=\frac{1}{2}\left(\pi_{i}\tilde{\pi}_{j}+\pi_{j}\tilde{\pi}_{i}\right)\,.$

For pointwise estimates of these quantities, refer to Lemma 11.3.

####  2.1.8. Schematic term notation

We will denote as $\mathfrak{T}_{1}\ast\dots\ast\mathfrak{T}_{l}$, where
$\mathfrak{T}_{i}$ are $\Sigma_{t}$-tangent tensors, any multiple of
$(\mathfrak{T}_{i})$, with regards to the rescaled adapted spatial metric $G$
or as standard multiplication if no summation over indices occurs between
factors. Constant prefactors and contractions with regards to $G$ are also
suppressed in this notation.

When working with terms where such notation is used, we will estimate these
inner products by $\lesssim\prod_{i=1}^{l}\lvert\mathfrak{T}_{i}\rvert_{G}$,
making any constant in front irrelevant, and further we can view any
contraction with regards to $G$ as a product of the non-contracted tensor
$\mathfrak{T}$ with $G$ or $G^{-1}$, and estimate that up to constant by
$\lvert G\rvert_{G}\lvert\mathfrak{T}\rvert_{G}$, where the first factor is
simply $\sqrt{3}$.
For similar products with respect to $\gamma$, we denote them by
$\ast_{\gamma}$.

####  2.1.9. On multiple derivatives of variables

For a scalar function $\zeta$, an $(r,s)$-tensor field $\mathfrak{T}$ and
capitalized integers
$I,J,\color[rgb]{0,0,0}\ldots\color[rgb]{0,0,0}\in\mathbb{N}_{0}$, we denote
by $\nabla^{I}\zeta$ and $\nabla^{I}\mathfrak{T}$ the tensors
$\nabla_{l_{1}}\dots\nabla_{l_{I}}\zeta$
and$\nabla_{l_{1}}\dots\nabla_{l_{I}}{\mathfrak{T}^{i_{1}\dots
i_{r}}}_{j_{1}\dots j_{s}}$. We extend this notation to other covariant
derivatives analogously. To avoid potential ambiguity with an index raised by
$g$, we will apply the following convention:

  * •

If a covariant derivative carries an uppercase letter, a formula with more
than one symbol or a positive integer in its superscript, this refers to
taking a derivative of that order.

  * •

If a covariant derivative carries a lowercase letter or $0$ in its
superscript, this refers to an index.

Further, we will only apply this notation where the precise distribution of
indices is not important (e.g. in schematic notation, see Section 2.1.8).

###  2.2. FLRW spacetimes and the Friedman equations

Herein, we collect the properties of the reference FLRW solution to the
Einstein scalar-field system in CMC-transported coordinates. Our main focus
will lie on the behaviour of the scale factor as determined by the Friedman
equations. Before moving on to that, we collect the information on the spatial
geometry we will need:

######  Definition 2.1 (Hyperbolic reference geometry).

$(M,\gamma)$ is a three-dimensional, connected, closed, orientable Riemannian
manifold with constant sectional curvature $-\frac{1}{9}$, hence Ricci tensor
$\text{{Ric}}[\gamma]=-\frac{2}{9}\gamma_{ij}$ and scalar curvature
$R[\gamma]=-\frac{2}{3}$.

######  Remark 2.2 (Orientability is not a restriction).

We assume that $M$ is orientable for the sake of simplicity. If $M$ should be
non-orientable, we may pass the initial data to the oriented double cover and
solve the problem there. Since the result is equivariant with respect to the
double covering map, this then solves the original problem.

With this in hand, we can express our classical family of solutions to the
Einstein scalar-field system as follows:

######  Lemma 2.3 (FLRW solutions and Friedman equations).

Consider FLRW spacetimes $(\overline{M},{\overline{g}}_{FLRW})$ with
$\overline{M}=(0,\infty)\times M$, where $(M,\gamma)$ is as in Definition 2.1
and where

(2.1) ${\overline{g}}_{FLRW}=-dt^{2}+a(t)^{2}\gamma$

holds for some $a\in C^{\infty}((0,\infty))$, with the conventions $a(0)=0$
and $\dot{a}(T)>0$ for some arbitrary $T>0$. Further, choose a (smooth) scalar
function ${\phi}_{FLRW}$ such that the following holds:

(2.2) $\partial_{t}{\phi}_{FLRW}=C\cdot
a(t)^{-3},\quad\nabla{\phi}_{FLRW}=0,\quad\square_{\overline{g}_{FLRW}}{\phi}_{FLRW}=0$

Such a pair $(\overline{g}_{FLRW},{\phi}_{FLRW})$ solves the Einstein scalar-
field system (1.1a)-(1.1b) on $\overline{M}$ if and only if $a$ satisfies the
Friedman equation

(2.3) $\dot{a}=\sqrt{\frac{1}{9}+\frac{4\pi}{3}C^{2}a^{-4}}\,.$

In particular, one has

(2.4) $\ddot{a}=-\frac{8\pi}{3}C^{2}a^{-5}\,.$

###### Proof.

The first statement follows from explicitly computing
$\text{{Ric}}[\overline{g}]$ as in [O’N83, p.345]. That (2.3) implies (2.4)
follows simply by computing the derivative of $\dot{a}^{2}$. ∎

In the subsequent analysis, the following properties of $a$ that follow from
(2.3) will be crucial for our analysis:

######  Lemma 2.4 (Scale factor analysis).

Let $a$ solve (2.3) with $a(0)=0$. Then $a$ is analytic on $(0,\infty)$ and
extends to a continuous function on $[0,\infty)$ with $a(t)\simeq
t^{\frac{1}{3}}$ being satisfied near $t=0$. Further, for any $p>0$, there
exist constants $c>0$ and $K_{p}>0$, where $c$ is independent of $p$ and
$K_{p}$ depends analytically on $p$, such that, for any $t\in(t,t_{0}]$, one
has

(2.5) $\exp\left(p\int_{t}^{t_{0}}a(s)^{-3}\,ds\right)\leq K_{p}a(t)^{-cp}$

and

(2.6)
$\int_{t}^{t_{0}}a(s)^{-3-p}\,ds\leq\frac{1}{p}a(t)^{-p},\quad\int_{t}^{t_{0}}a(s)^{-3+p}\,ds\lesssim\frac{1}{p}\,.$

Moreover, for any $t\in(0,t_{0}]$ and any $q>0$, there exist constants $c>0$
and $K>0$ which both are independent of $q$ such that one has

(2.7)
$\int_{t}^{t_{0}}a(s)^{-3}\,ds\leq\frac{K}{q}a(t)^{-cq}\,.\color[rgb]{0,0,0}$

Finally, (2.3) also implies

(2.8) $\sqrt{\frac{4\pi}{3}}Ca^{-2}\leq\dot{a}$

######  Remark 2.5.

We will use the estimates in Lemma 2.4 where $p$ is a positive power of
$\varepsilon$ up to algebraic constants. Then, we can simply replace $K_{p}$
in (2.5) by a uniform constant.

###### Proof.

For the first statement, we refer to [FU22, Lemma 2.1] with $\gamma=2$. We
also collect from there555In [FU22, Lemma 2.1], one at first only has
$\int_{t}^{t_{0}}a(s)^{-3}\lesssim\log(t_{0})-\log(t)$ for $t_{0}$ small
enough that we can estimate $a(t)$ by $t^{\frac{1}{3}}$ up to constant.
However, assuming this inequality were to hold up to $\overline{t}>0$ and one
had $t_{0}>\overline{t}$, the contribution
$\int_{\overline{t}}^{t_{0}}a(s)^{-3}\,ds$ only adds a constant that we can
absorb into our notation. Similarly, $a(t)\simeq t^{\frac{1}{3}}$ can be
assumed to hold on $(0,t_{0}]$ for any fixed $t_{0}>0$, and we can ignore this
technicality in proving the integral formulas in Lemma 2.4. that, for
$t<t_{0}$,

$\int_{t}^{t_{0}}a(s)^{-3}\,ds\lesssim 1+\lvert\log(t)\rvert$

is satisfied. Hence, there exists some $c^{\prime}>0$ such that the following
holds (with $t<1$ without loss of generality, see Footnote 5):

$\exp\left(p\int_{t}^{t_{0}}a(s)^{-3}\,ds\right)\leq\exp(c^{\prime}\cdot
p\log(t_{0}))\cdot\exp(-c^{\prime}\cdot p)\leq K_{p}t^{-c^{\prime}p}$

(2.5) follows by applying $a(t)\simeq t^{\frac{1}{3}}$. Noting that
$a^{-3}\simeq\nicefrac{{\dot{a}}}{{a}}$ holds, one further has

(2.9)
$\int_{t}^{t_{0}}a(s)^{-3-p}\,ds\lesssim\int_{a(t)}^{a(t_{0})}y^{-1-p}dy=\frac{1}{p}\left(a(t)^{-p}-a(t_{0})^{-p}\right)\leq\frac{1}{p}a(t)^{-p}\,,$

and the other inequality in (2.6) follows analogously. Finally, (2.7) follows
directly from (2.6) when assuming without loss of generality that
$a\lvert_{(t,t_{0})}$ only takes values in $(0,1)$.

∎

###  2.3. Solutions to the Einstein scalar-field equations in CMC gauge

From here on out, we impose the CMC condition666Recall that $k$ is negative in
our sign convention, see Section 2.1.5.

(2.10) ${k^{l}}_{l}(t,\cdot)=\tau(t)=-3\frac{\dot{a}(t)}{a(t)}\,.$

We use (2.3) and (2.4) to collect the following formulas for the mean
curvature:

(2.11) $\displaystyle\partial_{t}\tau$ $\displaystyle=12\pi
C^{2}a^{-6}+\frac{1}{3}a^{-2}$ (2.12) $\displaystyle\tau^{2}$
$\displaystyle=9\frac{\dot{a}^{2}}{a^{2}}=12\pi C^{2}a^{-6}+a^{-2}\,.$

We consequently define the trace-free component $\hat{k}$ of $k$ as

(2.13) $\hat{k}_{ij}=k_{ij}-\frac{\tau}{3}g_{ij}$

and recall that the future directed unit normal to our foliation is written as

(2.14) $\partial_{0}=n^{-1}\partial_{t}\,.$

With this, we can express the Einstein scalar-field equations in our gauge as
follows:

######  Proposition 2.6 (The Einstein scalar-field system in CMC gauge).

A pair $(\overline{g},\phi)$ solves the Einstein scalar-field equations
(1.1a)-(1.1c) on $I\times M$ in CMC gauge (2.10) for some interval
$I\subseteq(0,t_{0}]$, where the scale factor satisfies (2.3), if and only if
the following equations are satisfied on $I\times M$:
The metric evolution equations

(2.15a) $\displaystyle\partial_{t}g_{ij}=$
$\displaystyle\,-2nk_{ij}=-2n\hat{k}_{ij}+2n\frac{\dot{a}}{a}g_{ij}\,,$
(2.15b) $\displaystyle\partial_{t}\hat{k}_{ij}=$
$\displaystyle\,-\nabla_{i}\nabla_{j}n+n\left[\text{{Ric}}[g]_{ij}-\frac{\dot{a}}{a}\hat{k}_{ij}-2\hat{k}_{il}\hat{k}^{l}_{j}-8\pi\nabla_{i}\phi\nabla_{j}\phi\right]$
$\displaystyle\quad+4\pi
C^{2}a^{-6}(n-1)g_{ij}+\frac{1}{9}\left(3n-1\right)a^{-2}g_{ij}\,,$

the Hamiltonian and momentum constraint equations

(2.16a) $\displaystyle
R[g]+\frac{2}{3}\tau^{2}-\langle\hat{k},\hat{k}\rangle_{g}=$ $\displaystyle
8\pi\left[\lvert\partial_{0}\phi\rvert^{2}+\lvert\nabla\phi\rvert_{g}^{2}\right]\,,$
(2.16b) $\displaystyle\nabla^{l}\hat{k}_{lj}=$
$\displaystyle-8\pi\nabla_{j}\phi\cdot\partial_{0}\phi\,,$

the lapse equation

(2.17a) $\displaystyle\Delta_{g}n=$ $\displaystyle-12\pi
C^{2}a^{-6}-\frac{1}{3}a^{-2}+n\left[\frac{1}{3}a^{-2}+4\pi
C^{2}a^{-6}+\langle\hat{k},\hat{k}\rangle_{g}+8\pi\lvert\partial_{0}\phi\rvert^{2}\right]\,,$
or equivalently by (2.16a) (2.17b) $\Delta_{g}n=-12\pi
C^{2}a^{-6}-\frac{1}{3}a^{-2}+n\left[R[g]-8\pi\lvert\nabla\phi\rvert_{g}^{2}+12\pi
C^{2}a^{-6}+a^{-2}\right]\,,$

and the wave equation

(2.18)
$\square_{\overline{g}}\phi=-\partial_{0}^{2}\phi+n^{-1}g^{ij}\nabla_{i}n\nabla_{j}\phi+\Delta_{g}\phi+\tau\partial_{0}\phi=0\,.$

###### Proof.

These are standard equations that follow from [Ren08, Chapter 2.3] and
applying (2.10)-(2.12). ∎

###  2.4. Bel-Robinson variables

In this subsection, we briefly (re-)establish Bel-Robinson variables and how
they behave within the Einstein scalar-field system.

Recall that the Weyl tensor $W\equiv W[\overline{g}]$ is the trace-free
component of the spacetime curvature and, in the Einstein-scalar-field system,
takes the form

$\displaystyle W_{\alpha\beta\gamma\delta}=$
$\displaystyle\,\text{{Riem}}[\overline{g}]_{\alpha\beta\gamma\delta}-P[\overline{g}]_{\alpha\beta\gamma\delta}\,,$
$\displaystyle P[\overline{g}]_{\alpha\beta\gamma\delta}=$
$\displaystyle\,\frac{1}{2}\left(\overline{g}_{\alpha\gamma}\text{{Ric}}[\overline{g}]_{\beta\delta}-\overline{g}_{\beta\gamma}\text{{Ric}}[\overline{g}]_{\alpha\delta}-\overline{g}_{\alpha\delta}\text{{Ric}}[\overline{g}]_{\gamma\beta}+\overline{g}_{\beta\delta}\text{{Ric}}[\overline{g}]_{\alpha\gamma}\right)-\frac{1}{6}R[\overline{g}]\left(\overline{g}_{\alpha\gamma}\overline{g}_{\beta\delta}-\overline{g}_{\alpha\delta}\overline{g}_{\beta\gamma}\right)$
$\displaystyle=$
$\displaystyle\,4\pi\left(\overline{g}_{\alpha\gamma}\overline{\nabla}_{\beta}\phi\overline{\nabla}_{\delta}\phi-\overline{g}_{\beta\gamma}\overline{\nabla}_{\alpha}\phi\overline{\nabla}_{\delta}\phi-\overline{g}_{\alpha\delta}\overline{\nabla}_{\beta}\phi\overline{\nabla}_{\gamma}\phi+\overline{g}_{\beta\delta}\overline{\nabla}_{\alpha}\phi\overline{\nabla}_{\gamma}\phi\right)$
$\displaystyle\,-\frac{4\pi}{3}\left(\overline{\nabla}^{\rho}\phi\overline{\nabla}_{\rho}\phi\right)\left(\overline{g}_{\alpha\gamma}\overline{g}_{\beta\delta}-\overline{g}_{\alpha\delta}\overline{g}_{\beta\gamma}\right)\,.$

We define the dual $W^{\ast}$ of the Weyl tensor as

$W^{\ast}_{\alpha\beta\gamma\delta}=\frac{1}{2}\bm{\varepsilon}_{\alpha\beta\mu\nu}{W^{\mu\nu}}_{\gamma\delta}\,.$

The electric and magnetic components of the Weyl tensor, referred to as the
Bel-Robinson variables from here on, are now defined as

$E(W)_{\alpha\beta}=W_{\alpha\mu\beta\nu}\partial_{0}^{\mu}\partial_{0}^{\nu}=W_{\alpha
0\beta 0},\
B(W)_{\alpha\beta}=W^{\ast}_{\alpha\mu\beta\nu}\partial_{0}^{\mu}\partial_{0}^{\nu}=W^{\ast}_{\alpha
0\beta 0}\,.$

We note that, conversely, the Weyl tensor can be fully reconstructed from $E$
and $B$ since the following identities hold:

(2.19) $W_{a0c0}=E_{ac},\quad
W_{abc0}=-{\bm{\varepsilon}_{ab}}^{m}B_{mc},\quad
W_{abcd}=-\bm{\varepsilon}_{abi}\bm{\varepsilon}_{cdj}E^{ij}$

By the symmetries of the Weyl tensor as a whole, $E$ and $B$ are symmetric and
one has $E_{0\beta}=0=B_{0\beta}$. Hence, $E$ and $B$ are symmetric, tracefree
$\Sigma_{t}$-tangent $(0,2)$-tensors which we shall simply denote as $E_{ij}$
and $B_{ij}$.
Further, we define

(2.20)
$J_{\beta\gamma\delta}:=\overline{\nabla}^{\alpha}W_{\alpha\beta\gamma\delta},\quad
J^{\ast}_{\beta\gamma\delta}:=\overline{\nabla}^{\alpha}W^{\ast}_{\alpha\beta\gamma\delta}\,.$

By applying the Bianchi identity for $\text{{Riem}}[\overline{g}]$, we gain
the explicit expression

(2.21)
$J_{\beta\gamma\delta}=\frac{1}{2}\left(\overline{\nabla}_{\gamma}\text{{Ric}}[\overline{g}]_{\delta\beta}-\overline{\nabla}_{\delta}\text{{Ric}}[\overline{g}]_{\gamma\beta}\right)-\frac{1}{12}\left(\overline{g}_{\beta\delta}\overline{\nabla}_{\gamma}R[\overline{g}]-\overline{g}_{\beta\gamma}\overline{\nabla}_{\delta}R[\overline{g}]\right)\,.$

Using (1.1a), we collect:

(2.22) $\displaystyle\color[rgb]{0,0,0}J_{i0j}=\color[rgb]{0,0,0}$
$\displaystyle\color[rgb]{0,0,0}\,4\pi\Bigr{[}\nabla_{i}(\partial_{0}\phi)\nabla_{j}\phi+{k^{l}}_{i}\nabla_{l}\phi\nabla_{j}\phi-\partial_{0}\phi\nabla_{i}\nabla_{j}\phi$
$\displaystyle\color[rgb]{0,0,0}\,\phantom{4\pi}-k_{ij}(\partial_{0}\phi)^{2}-n^{-1}\nabla_{i}n\cdot\nabla_{j}\phi\cdot\partial_{0}\phi\Bigr{]}-\frac{2\pi}{3}\Bigr{[}\partial_{0}(\overline{\nabla}^{\alpha}\phi\overline{\nabla}_{\alpha}\phi)\Bigr{]}g_{ij}$
(2.23) $\displaystyle\color[rgb]{0,0,0}J^{\ast}_{i0j}=$
$\displaystyle\,4\pi\bm{\varepsilon}_{lmj}\left(\nabla^{l}\nabla_{i}\phi+{k^{l}}_{i}\partial_{0}\phi\right)\nabla^{m}\phi+\frac{2\pi}{3}\bm{\varepsilon}_{imj}\nabla^{m}\left(\overline{\nabla}^{\alpha}\phi\overline{\nabla}_{\alpha}\phi\right)$

Note that expressions containing
$\overline{\nabla}^{\alpha}\phi\overline{\nabla}_{\alpha}\phi$ can be ignored
throughout our analysis since they are either pure trace or antisymmetric and
thus will cancel in inner products with $E$, $B$, $\hat{k}$ and their rescaled
analogues.
The Bel-Robinson variables then behave as follows:

######  Lemma 2.7 (Constraint and evolution equations for Bel-Robinson
variables).

If $(\overline{g},\phi)$ is a classical solution to the Einstein scalar-field
system (1.1a)-(1.1b) in CMC gauge (see (2.10)), $E$ and $B$ satisfy the
following constraint equations:

(2.24a) $\displaystyle E$
$\displaystyle\,=\text{{Ric}}[g]+\frac{2}{9}\tau^{2}g+\frac{\tau}{3}\hat{k}-\hat{k}\odot_{g}\hat{k}-4\pi(\nabla\phi\otimes\nabla\phi)-\left(\frac{8\pi}{3}\left\lvert\partial_{0}\phi\right\rvert^{2}+\frac{4\pi}{3}\lvert\nabla\phi\rvert_{g}^{2}\right)g$
(2.24b) $\displaystyle B$ $\displaystyle\,=-\text{{curl}}\hat{k}$

Further, they satisfy the following evolution equations:

(2.25a) $\displaystyle\partial_{t}E_{ij}$
$\displaystyle=n\text{{curl}}{B}_{ij}-(\nabla n\wedge
B)_{ij}-\frac{5}{2}n\left(E\times k\right)_{ij}-\frac{2}{3}n(E\cdot
k)g_{ij}-\frac{\tau}{2}n\cdot
E_{ij}-\color[rgb]{0,0,0}\frac{n}{2}\left(J_{i0j}+J_{j0i}\right)\color[rgb]{0,0,0}$
(2.25b) $\displaystyle\partial_{t}B_{ij}$
$\displaystyle=-n\text{{curl}}{E}_{ij}+(\nabla n\wedge
E)_{ij}-\frac{5}{2}n\left(B\times k\right)_{ij}-\frac{2}{3}n(B\cdot
k)g_{ij}-\frac{\tau}{2}n\cdot
B_{ij}-\color[rgb]{0,0,0}\frac{n}{2}\left(J_{i0j}^{\ast}+J_{j0i}^{\ast}\right)\color[rgb]{0,0,0}$

###### Proof.

For (2.25a)-(2.25b), we refer to [AM04, (3.11a)-(3.11b)].777Note that there is
a minor error in the statement in [AM04], where the authors forget to
symmetrize the $J$-tensors when applying the symmetrization to [AM04, (3.14)].
This error seems to also occur in [CK93, (7.2.2jk)]. (2.24a)-(2.24b) follow as
in [Wan19, (3.63a)-(3.63b)] from contracting the Gauss-Codazzi constraints. ∎

######  Remark 2.8 (Inital data for Bel-Robinson variables).

Since the Weyl tensor vanishes over FLRW spacetimes, so do
$E(W[\overline{g}_{FLRW}])$ and $B(W[\overline{g}_{FLRW}])$. Furthermore, note
that given initial data
$(M,\mathring{g},\mathring{k},\mathring{\pi},\mathring{\psi})$ on
$\Sigma_{t_{0}}$ in the sense discussed in Section 1.2.1, and defining
$\hat{\mathring{k}}=\mathring{k}-\frac{\tau}{3}\mathring{g}$, we can use
(2.24a) and (2.24b) to define the following $(0,2)$-tensors:

(2.26a) $\displaystyle\mathring{E}=$
$\displaystyle\,\text{{Ric}}[\mathring{g}]+\frac{2}{9}\tau^{2}\mathring{g}+\frac{\tau}{3}\hat{\mathring{k}}-\hat{\mathring{k}}\odot_{\mathring{g}}\hat{\mathring{k}}-4\pi\left(\mathring{\pi}\otimes\mathring{\pi}\right)-\left(\frac{8\pi}{3}\left\lvert\mathring{\psi}\right\rvert_{\mathring{g}}^{2}+\frac{4\pi}{3}\left\lvert\mathring{\pi}\right\rvert_{\mathring{g}}^{2}\right)\mathring{g}$
(2.26b) $\displaystyle\mathring{B}=$
$\displaystyle\,-\text{{curl}}_{\mathring{g}}\hat{\mathring{k}}$

These are easily seen to be symmetric, and the constraints (1.3a) and (1.3b)
on the initial data ensure that they are also tracefree. Hence, any choice of
initial data for the Cauchy problem immediately also contains a unique choice
of initial data for the Bel-Robinson variables that is consistent with
solutions to the Einstein scalar-field equations in CMC gauge.

###  2.5. Rescaled variables and equations

It will be more convenient to work the rescaled and shifted solution variables
to measure their distance from the FLRW reference solution. In this
subsection, we introduce the renormalized solution variables and restate the
Einstein scalar-field system in terms of these variables.

######  Definition 2.9 (Rescaled variables for Big Bang stability).

We will consider the rescaled variables

(2.27a) $\displaystyle
G_{ij}=a^{-2}g_{ij},\quad(G^{-1})^{ij}=a^{2}g^{ij},\quad\Sigma_{ij}=a\hat{k}_{ij}$
(2.27b) $\displaystyle N=n-1$ (2.27c) $\displaystyle\bm{E}_{ij}=a^{4}\cdot
E_{ij},\quad\bm{B}_{ij}=a^{4}\cdot B_{ij}$ (2.27d)
$\displaystyle\Psi=a^{3}\partial_{0}\phi-C,$

We note that the scaling of $\bm{B}$ in (2.27c) is not the asymptotic
rescaling of $B$ – in fact, we expect $B$ to have (approximate) leading order
$a^{-2}$ as one can see in (4.4g). However, keeping this scaling parallel to
that of $\bm{E}$ makes the structurally very similar evolution equations
significantly easier to deal with. We also don’t rescale $N$ asymptotically,
unlike [RS18b, Spe18], but note that $N$ converges to $0$ at an order slightly
below $a^{4}$ at low orders (see (3.17h) and (8.3a)).

######  Proposition 2.10 (The rescaled Einstein scalar-field system).

The Einstein scalar-field system in CMC gauge as in Propopsition 2.6 are
solved by $(g,\hat{k},n,\nabla\phi,\partial_{0}\phi)$ if the rescaled
variables $(G,\Sigma,N,\nabla\phi,\Psi,\bm{E},\bm{B})$ as in Definition 2.9
solve the following set of equations:888We refer to Lemma 11.3 for the
scalings that occur when switching between tensor field operations like
$\wedge$ and $\wedge_{G}$. The rescaled metric evolution equations

(2.28a) $\displaystyle\partial_{t}G_{ij}=$
$\displaystyle-2(N+1)a^{-3}\Sigma_{ij}+2N\frac{\dot{a}}{a}G_{ij}$ (2.28b)
$\displaystyle\partial_{t}(G^{-1})^{ij}=$
$\displaystyle\,2(N+1)a^{-3}(\Sigma^{\sharp})^{ij}-2N\frac{\dot{a}}{a}(G^{-1})^{ij}$
(2.28c) $\displaystyle\partial_{t}\Sigma_{ij}=$
$\displaystyle\,-a\nabla_{i}\nabla_{j}N+(N+1)\left[a\text{{Ric}}[G]_{ij}-2a^{-3}(\Sigma\odot_{G}\Sigma)_{ij}-8\pi
a\nabla_{i}\phi\nabla_{j}\phi\right]$ $\displaystyle+4\pi C^{2}a^{-3}\cdot
NG_{ij}+\frac{1}{9}\left(3N+2\right)aG_{ij}+N\frac{\dot{a}}{a}\Sigma_{ij}$
(2.28d) $\displaystyle\partial_{t}(\Sigma^{\sharp})^{a}_{b}=$
$\displaystyle\,\tau N(\Sigma^{\sharp})^{a}_{b}-a\nabla^{\sharp
a}\nabla_{b}N+(N+1)a\left[\left(\text{{Ric}}[G]^{\sharp}\right)^{a}_{b}+\frac{2}{9}\mathbb{I}^{a}_{b}\right]\,$
$\displaystyle-8\pi(N+1)a\nabla^{\sharp a}\phi\nabla_{b}\phi+N\left(4\pi
C^{2}a^{-3}+\frac{1}{9}a\right)\mathbb{I}^{a}_{b}\,,$

the rescaled Hamiltonian and momentum constraints

(2.29a) $\displaystyle R[G]+\frac{2}{3}-a^{-4}\langle\Sigma,\Sigma\rangle_{G}$
$\displaystyle=8\pi\left[a^{-4}\Psi^{2}+2Ca^{-4}\Psi+\lvert\nabla\phi\rvert_{G}^{2}\right]$
(2.29b) $\displaystyle\nabla^{\sharp m}\Sigma_{ml}$
$\displaystyle=-8\pi\nabla_{l}\phi(\Psi+C)$ with their Bel-Robinson analogues
(2.29c) $\displaystyle\bm{E}_{ij}=$
$\displaystyle\,a^{4}\left(\text{{Ric}}[G]_{ij}+\frac{2}{9}G_{ij}\right)+\frac{\tau}{3}a^{3}\Sigma_{ij}-(\Sigma\odot_{G}\Sigma)_{ij}-4\pi
a^{4}\nabla_{i}\phi\nabla_{j}\phi$
$\displaystyle\,-\left[\frac{4\pi}{3}a^{4}\lvert\nabla\phi\rvert_{G}^{2}+\frac{8\pi}{3}\Psi^{2}+\frac{16\pi}{3}C\Psi\right]G_{ij}$
(2.29d) $\displaystyle\bm{B}_{ij}=$
$\displaystyle\,-\color[rgb]{0,0,0}a^{2}\text{{curl}}_{G}\Sigma_{ij}\color[rgb]{0,0,0}\,,$

the rescaled lapse equation

(2.30a) $\Delta N=\left(12\pi
C^{2}a^{-4}+\color[rgb]{0,0,0}\frac{1}{3}\color[rgb]{0,0,0}\right)N+(N+1)a^{-4}\left[\langle\Sigma,\Sigma\rangle_{G}+8\pi\Psi^{2}+16\pi
C\Psi\right]$ or equivalently (2.30b) $\Delta N=\left(12\pi
C^{2}a^{-4}+\color[rgb]{0,0,0}\frac{1}{3}\color[rgb]{0,0,0}\right)N+(N+1)\left[R[G]+\frac{2}{3}-8\pi\lvert\nabla\phi\rvert_{G}^{2}\right]\,,$

the rescaled evolution equations for the Bel-Robinson variables999

(2.31a) $\displaystyle\partial_{t}\bm{E}_{ij}=$
$\displaystyle\,\left(3-N\right)\frac{\dot{a}}{a}\bm{E}_{ij}\color[rgb]{0,0,0}-a^{-1}(\nabla
N\wedge_{G}\bm{B})_{ij}+(N+1)a^{-1}\text{{curl}}_{G}\bm{B}_{ij}\color[rgb]{0,0,0}$
$\displaystyle\,-(N+1)\left[\color[rgb]{0,0,0}\frac{5}{2}a^{-3}(\bm{E}\times_{G}\Sigma)_{ij}\color[rgb]{0,0,0}+\frac{2}{3}a^{-3}\langle\bm{E},\Sigma\rangle_{G}G_{ij}\right]$
$\displaystyle\,+4\pi(N+1)a^{-3}(\Psi+C)^{2}\Sigma_{ij}-4\pi(N+1)\dot{a}a^{3}\nabla_{i}\phi\nabla_{j}\phi\color[rgb]{0,0,0}+4\pi
a\nabla_{(i}N\nabla_{j)}\phi(\Psi+C)\color[rgb]{0,0,0}$ $\displaystyle\,-4\pi
a(N+1)\left[\nabla_{i}\Psi\nabla_{j}\phi+\nabla_{j}\Psi\nabla_{i}\phi+\color[rgb]{0,0,0}(\Sigma^{\sharp})^{l}_{(i}\nabla_{j)}\phi\color[rgb]{0,0,0}\nabla_{l}\phi-(\Psi+C)\nabla_{i}\nabla_{j}\phi\right]$
$\displaystyle\,+(N+1)\left[\frac{2\pi}{3}a^{6}\partial_{0}\left(a^{-6}(\Psi+C)^{2}+a^{-2}\lvert\nabla\phi\rvert_{G}^{2}\right)+4\pi\frac{\dot{a}}{a}(\Psi+C)^{2}\right]G_{ij}$
(2.31b) $\displaystyle\partial_{t}\bm{B}_{ij}=$
$\displaystyle\,\frac{\dot{a}}{a}\left(3-N\right)\bm{B}_{ij}\color[rgb]{0,0,0}+a^{-1}(\nabla
N\wedge_{G}\bm{E})_{ij}-(N+1)a^{-1}\text{{curl}}_{G}\bm{E}_{ij}\color[rgb]{0,0,0}$
$\displaystyle\,-(N+1)\left[\color[rgb]{0,0,0}\frac{5}{2}a^{-3}(\bm{B}\times_{G}\Sigma)_{ij}+\frac{2}{3}a^{-3}\langle\bm{B},\Sigma\rangle_{G}G_{ij}\color[rgb]{0,0,0}\right]\,$
$\displaystyle\,-4\pi(N+1)\bm{\varepsilon}[G]_{lmj}\left(a^{3}\nabla^{\sharp
l}\nabla_{\color[rgb]{0,0,0}j\color[rgb]{0,0,0}}\phi\nabla^{\sharp
m}\phi+a^{-1}{(\Sigma^{\sharp})^{l}}_{i}\nabla^{\sharp m}\phi(\Psi+C)\right)$

and the rescaled wave equation

(2.32a) $\partial_{t}\Psi=a\langle\nabla
N,\nabla\phi\rangle_{G}+a(N+1)\Delta\phi-3\frac{\dot{a}}{a}N(\Psi+C)$ along
with the evolution equation (2.32b)
$\partial_{t}\nabla\phi=a^{-3}(N+1)\nabla\Psi+a^{-3}(\Psi+C)\nabla N\,.$

Finally, we collect the rescaled Ricci evolution equation

(2.33) $\displaystyle\partial_{t}\text{{Ric}}[G]_{ab}=$
$\displaystyle\,a^{-3}(N+1)(\Delta_{G}\Sigma_{ab}-\nabla^{\sharp
d}\nabla_{a}\Sigma_{db}-\nabla^{\sharp d}\nabla_{b}\Sigma_{da})$
$\displaystyle\,+a^{-3}\nabla^{\sharp
d}N(2\nabla_{d}\Sigma_{ab}-\nabla_{a}\Sigma_{db}-\nabla_{b}\Sigma_{da})$
$\displaystyle\,-a^{-3}\left(\nabla_{a}N(\text{{div}}_{G}\Sigma)_{b}+\nabla_{b}(\text{{div}}_{G}\Sigma)_{a}\right)+\Delta_{G}N(a^{-3}\Sigma_{ab}+\frac{\tau}{3}G_{ab})$
$\displaystyle\,\color[rgb]{0,0,0}-\color[rgb]{0,0,0}a^{-3}\Bigr{(}\nabla^{\sharp
d}\nabla_{a}N\cdot\Sigma_{db}+\nabla^{\sharp
d}\nabla_{b}N\cdot\Sigma_{da}\Bigr{)}+\frac{\tau}{3}\nabla_{a}\nabla_{b}N$

as well as (in a coordinate neighbourhood) the Christoffel evolution equation

(2.34) $\displaystyle\partial_{t}\Gamma_{ij}^{k}[G]=$
$\displaystyle\,\frac{1}{2}(G^{-1})^{kl}\left(\nabla_{i}(\partial_{t}G_{jl})+\nabla_{j}(\partial_{t}G_{il})-\nabla_{l}(\partial_{t}G_{ij})\right)$
$\displaystyle=$
$\displaystyle\,-(N+1)a^{-3}\left[\nabla_{i}{(\Sigma^{\sharp})^{k}}_{j}+\nabla_{j}{(\Sigma^{\sharp})^{k}}_{i}-\nabla^{\sharp
k}\Sigma_{ij}\right]$
$\displaystyle\,-a^{-3}\left[\nabla_{i}N{(\Sigma^{\sharp})^{k}}_{j}+\nabla_{j}N{(\Sigma^{\sharp})^{k}}_{i}-\nabla^{\sharp
k}N\Sigma_{ij}\right]$
$\displaystyle\,+\frac{\dot{a}}{a}\left[\nabla_{i}N\cdot\mathbb{I}^{k}_{j}+\nabla_{j}N\cdot\mathbb{I}^{k}_{i}-\nabla^{\sharp
k}N\cdot G_{ij}\right]$

###### Proof.

For the first identity in (2.34), we refer to [CBCLN06, Lemma 2.27] and insert
the evolution equation (2.28a). Otherwise, all equations simply follow by
computing the effects of rescaling on the equations from Proposition 2.6
(respectively the Ricci evolution equation as in [Ren08, Chapter 2.3, (2.32)])
as well as the Bel-Robinson evolution equations (2.25a)-(2.25b) and constraint
equations (2.24a)-(2.24b). Notice that one already finds a solution to the
system in Proposition (2.6) with the rescaled variables excluding (2.31a),
(2.31b), (2.29c) and (2.29d). Conversely, all of the rescaled equations are
satisfied by a solution to Proposition 2.6 at sufficiently high regularity.
Hence, solving the full rescaled system is always sufficient to solve the
Einstein system in Proposition 2.6 and they are equivalent if the initial data
is regular enough to ensure sufficiently high regularity of solutions.
∎

###  2.6. Commuted equations

We collect Laplace-commuted versions of the equations for the rescaled
variables in Proposition 2.10 in this subsection. For the sake of brevity, we
won’t state all possible commutations for every equation, but restrict
ourselves to the ones we actually need within the bootstrap argument. We also
refer to Appendix 11.2 for expressions for commutators of spatial differential
operators with each other and with $\partial_{t}$.

The terms written down explicity in Lemma 2.11 are ones that dominate the
evolution behaviour or that are the largest higher order terms, both of which
require careful treatment within the bootstrap argument. The error terms are
broadly categorized into three groups:

  * •

“Borderline” terms are terms that critically contribute to the fact that the
energies diverge toward the Big Bang singularity. This almost always takes the
form of adding energy terms at the same order as the evolved variable scaled
by factors of the type $\varepsilon a^{-3}$ or $\varepsilon
a^{-3-c\sqrt{\varepsilon}}$, which causes the energies to slightly diverge
since $a^{-3}$ is barely not integrable (see (2.6)).

  * •

“Junk” terms are terms that are subcritical in the sense that they lead to
integrable error terms, or terms that only contain lower order derivatives of
the solution variables.

  * •

“Top” order terms (which only appear in (2.38a) and (2.38b)) are terms that
are terms that are junk terms for low order energies, but become borderline
terms at top order.

All of these error terms are tracked schematically in Section 11.3. Since we
will only need $L^{2}_{G}$-bounds on these error terms, which are given in
Section 11.4, we will treat them as notational “black boxes” outside of the
appendix.

######  Lemma 2.11 (Laplace-commuted rescaled equations).

Let $L\in 2\mathbb{N},\,L\geq 2$. With error terms as defined in Appendix
11.3, the system in Proposition 2.10 leads to the following Laplace-commuted
equations:
The Laplace-commuted rescaled evolution equations for the second fundamental
form

(2.35)
$\partial_{t}\Delta^{\frac{L}{2}}\Sigma=-a\nabla^{2}\Delta^{\frac{L}{2}}N+a(N+1)\Delta^{\frac{L}{2}}\text{{Ric}}[G]+\mathfrak{S}_{L,Border}+\mathfrak{S}_{L,Junk}\,,$

the Laplace-commuted rescaled momentum constraint equations

(2.36a)
$\text{{div}}_{G}\Delta^{\frac{L}{2}}\Sigma=-8\pi(\Psi+C)\left[\nabla\Delta^{\frac{L}{2}}\phi+\Delta^{\frac{L}{2}-1}\text{{Ric}}[G]\ast\nabla\phi\right]+\nabla\Delta^{\frac{L}{2}-1}\text{{Ric}}[G]\ast\Sigma+\mathfrak{M}_{L,Junk}$
and (2.36b)
$\color[rgb]{0,0,0}\text{{curl}}_{G}\Delta^{\frac{L}{2}}\Sigma=-a^{-2}\Delta^{\frac{L}{2}}\bm{B}+\bm{\varepsilon}[G]\ast\nabla\Delta^{\frac{L}{2}-1}\text{{Ric}}[G]\ast\Sigma+\tilde{\mathfrak{M}}_{L,Junk}\color[rgb]{0,0,0}\,,$
the Laplace-commuted rescaled Hamiltonian constraint equations (2.36c)
$\displaystyle\,\Delta^{\frac{L}{2}}R[G]+a^{-4}\sum_{I_{1}+I_{2}=L}\nabla^{I_{1}}\Sigma\ast\nabla^{I_{2}}\Sigma$
$\displaystyle=$ $\displaystyle\,16\pi
Ca^{-4}\Delta^{\frac{L}{2}}\Psi+a^{-4}\sum_{I_{1}+I_{2}=L}\left[\nabla^{I_{1}}\Psi\ast\nabla^{I_{2}}\Psi+\nabla^{I_{1}+1}\phi\ast\nabla^{I_{2}+1}\phi\right]$
and (2.36d)
$\Delta^{\frac{L}{2}}\text{{Ric}}[G]=a^{-4}\Delta^{\frac{L}{2}}\bm{E}-\frac{\tau}{3}a^{-1}\Delta^{\frac{L}{2}}\Sigma+\mathfrak{H}_{L,Border}+\mathfrak{H}_{L,Junk}\,,$
the Laplace-commuted rescaled lapse equations (2.37a)
$\displaystyle\Delta^{\frac{L}{2}+1}N=$ $\displaystyle\,\left(12\pi
C^{2}a^{-4}+\frac{1}{3}\right)\Delta^{\frac{L}{2}}N+16\pi
Ca^{-4}\cdot\Delta^{\frac{L}{2}}\Psi+\mathfrak{N}_{L,Border}+\mathfrak{N}_{L,Junk}\,,$
(2.37b) $\displaystyle\nabla\Delta^{\frac{L}{2}+1}N=$
$\displaystyle\,\left(12\pi
C^{2}a^{-4}+\frac{1}{3}\right)\nabla\Delta^{\frac{L}{2}}N+16\pi
Ca^{-4}\cdot\nabla\Delta^{\frac{L}{2}}\Psi+\mathfrak{N}_{L+1,Border}+\mathfrak{N}_{L+1,Junk}\,,$

the Laplace-commuted rescaled Bel-Robinson evolution equations

(2.38a) $\displaystyle\partial_{t}\Delta^{\frac{L}{2}}\bm{E}=$
$\displaystyle\,\frac{\dot{a}}{a}\left(3-N\right)\Delta^{\frac{L}{2}}\bm{E}\color[rgb]{0,0,0}+(N+1)a^{-1}\text{{curl}}_{G}\Delta^{\frac{L}{2}}\bm{B}-a^{-1}\nabla\Delta^{\frac{L}{2}}N\wedge_{G}\bm{B}\color[rgb]{0,0,0}$
$\displaystyle\,+4\pi
C^{2}a^{-3}(N+1)\Delta^{\frac{L}{2}}\Sigma\color[rgb]{0,0,0}+4\pi
a(\Psi+C)\nabla\Delta^{\frac{L}{2}}N\otimes\nabla\phi\color[rgb]{0,0,0}$
$\displaystyle\,+4\pi a(\Psi+C)(N+1)\nabla^{2}\Delta^{\frac{L}{2}}\phi-8\pi
a(N+1)\left(\nabla\phi\otimes\nabla\Delta^{\frac{L}{2}}\Psi\right)$
$\displaystyle\,+\mathfrak{E}_{L,Border}+\mathfrak{E}_{L,top}+\mathfrak{E}_{L,Junk}$
(2.38b) $\displaystyle\partial_{t}\Delta^{\frac{L}{2}}\bm{B}=$
$\displaystyle\,\frac{\dot{a}}{a}(3-N)\Delta^{\frac{L}{2}}\bm{B}\color[rgb]{0,0,0}-(N+1)a^{-1}\text{{curl}}_{G}\Delta^{\frac{L}{2}}\bm{E}+a^{-1}\nabla\Delta^{\frac{L}{2}}N\wedge_{G}\bm{E}\color[rgb]{0,0,0}$
$\displaystyle\,+a^{3}\bm{\varepsilon}[G]\ast\nabla^{2}\Delta^{\frac{L}{2}}\phi\ast\nabla\phi+\mathfrak{B}_{L,Border}+\mathfrak{B}_{L,top}+\mathfrak{B}_{L,Junk}\,,$

the Laplace-commuted rescaled matter evolution equations

(2.39a) $\displaystyle\partial_{t}\Delta^{\frac{L}{2}}\Psi=$
$\displaystyle\,a\langle\nabla\Delta^{\frac{L}{2}}N,\nabla\phi\rangle_{G}+a(N+1)\Delta^{\frac{L}{2}+1}\phi-3C\frac{\dot{a}}{a}\Delta^{\frac{L}{2}}N+\mathfrak{P}_{L,Border}+\mathfrak{P}_{L,Junk}$
(2.39b) $\displaystyle\partial_{t}\nabla\Delta^{\frac{L}{2}}\phi=$
$\displaystyle\,a^{-3}(N+1)\nabla\Delta^{\frac{L}{2}}\Psi+Ca^{-3}\nabla\Delta^{\frac{L}{2}}N+\mathfrak{Q}_{L,Border}+\mathfrak{Q}_{L,Junk}$
as well as (also allowing $L=0$ for (2.39d)) (2.39c)
$\displaystyle\partial_{t}\nabla_{l}\Delta^{\frac{L}{2}}\Psi=$
$\displaystyle\,a\nabla_{l}\nabla^{\sharp
j}\Delta^{\frac{L}{2}}N\nabla_{j}\phi+a(N+1)\nabla_{l}\Delta^{\frac{L}{2}+1}\phi-3C\frac{\dot{a}}{a}\nabla_{l}\Delta^{\frac{L}{2}}N$
$\displaystyle\,+\left(\mathfrak{P}_{L+1,Border}\right)_{l}+\left(\mathfrak{P}_{L+1,Junk}\right)_{l}$
(2.39d) $\displaystyle\partial_{t}\Delta^{\frac{L}{2}+1}\phi=$
$\displaystyle\,a^{-3}(N+1)\Delta^{\frac{L}{2}+1}\Psi+Ca^{-3}\Delta^{\frac{L}{2}+1}N+\mathfrak{Q}_{L+1,Border}+\mathfrak{Q}_{L+1,Junk}$

and the Laplace-commuted rescaled Ricci evolution equations

(2.40a) $\displaystyle\partial_{t}\Delta^{\frac{L}{2}}\text{{Ric}}[G]_{ij}=$
$\displaystyle\,a^{-3}\left(\Delta^{\frac{L}{2}+1}\Sigma_{ij}-2\nabla^{\sharp
m}\nabla_{(i}\Delta^{\frac{L}{2}}\Sigma_{j)m}\right)$
$\displaystyle\,-\frac{\dot{a}}{a}\left(\nabla_{i}\nabla_{j}\Delta^{\frac{L}{2}}N+\Delta^{\frac{L}{2}+1}N\cdot
G_{ij}\right)+\left(\mathfrak{R}_{L,Border}\right)_{ij}+\left(\mathfrak{R}_{L,Junk}\right)_{ij}$
(2.40b)
$\displaystyle\partial_{t}\nabla_{k}\Delta^{\frac{L}{2}}\text{{Ric}}[G]_{ij}=$
$\displaystyle\,a^{-3}\left(\nabla_{k}\Delta^{\frac{L}{2}+1}\Sigma_{ij}-2\nabla_{k}\nabla^{\sharp
m}\nabla_{(i}\Delta^{\frac{L}{2}}\Sigma_{j)m}\right)$
$\displaystyle\,-\frac{\dot{a}}{a}\left(\nabla_{k}\nabla_{i}\nabla_{j}\Delta^{\frac{L}{2}}N+\nabla_{k}\Delta^{\frac{L}{2}+1}N\cdot
G_{ij}\right)$
$\displaystyle\,+\left(\mathfrak{R}_{L+1,Border}\right)_{ijk}+\left(\mathfrak{R}_{L+1,Junk}\right)_{ijk}\,.$

###### Proof.

The equations (2.36c),(2.36d) and (2.37a) are obtained by simply applying
$\Delta^{\frac{L}{2}}$ on both sides of (2.29a),(2.29c) and (2.30a)
respectively. For (2.36a) and (2.36b), we additionally use the commutator
formulas (11.7e) and (11.7f), while for the evolution equations, we apply the
respective commutator of $\partial_{t}$ and spatial derivatives as collected
in Lemma 11.7 and commute Laplacians past $\nabla$ and curl where needed (see
the commutators in Lemma 11.5). The commutators with $\partial_{t}$ only cause
additional borderline and junk terms that don’t substantially influence the
behaviour, while the spatial commutators often lead to high order curvature
terms, for example the Ricci terms in (2.36a), that need to be more carefully
tracked. ∎

######  Remark 2.12 (Simplified junk term notation).

For junk terms that occur in an inner product with a tracefree symmetric
tensor, any terms that are pure trace will immediately cancel and thus don’t
need to be taken into consideration for the following estimates, even if they
have to be written down in the junk terms. Hence, we will denote with a
superscript “$\parallel$” (for example $\mathfrak{H}_{L,Junk}^{\parallel}$) on
a schematic error term the expressions that arise when dropping all terms of
the form $\zeta\cdot G$ for some scalar function $\zeta$ that occur in this
term’s definition (see, for example, (11.11d)).

##  3\. Big Bang stability: Norms, energies and bootstrap assumptions

Herein, we state the norms and energies we use to control the solution
variables. These will allow us to state our initial data and bootstrap
assumptions, and we then provide which improvement we aim to achieve for the
latter. Note that we don’t provide the coerciveness of our energies
immediately (and actually cannot, at least not in a manner useful to our
analysis), but will establish Sobolev norm control in the proof of Corollary
7.3, the key ingredient being Lemma 4.5. Furthermore, we collect a local well-
posedness statement from previous work in Section 3.4.

###  3.1. Norms

Recall that $\gamma$ is the hyperbolic spatial reference metric on $M$
introduced in Definition 2.1, which we view as a metric on any foliation
hypersurface $\Sigma_{t}$ (see Section 1.2.1), and $G$ is the rescaled spatial
metric arising from the evolution (see Definition 2.9).

######  Definition 3.1 (Pointwise norms and volume forms).

We denote by $\lvert\cdot\rvert_{\gamma}$ (resp.
$\lvert\cdot\rvert_{G(t,\cdot)}$) the pointwise norm with regards to $\gamma$
(resp. $G(t,\cdot)$). For the sake of simplicity, we define
$\lvert\zeta\rvert_{\gamma}=\lvert\zeta\rvert_{G(t,\cdot)}=\lvert\zeta(t,\cdot)\rvert$
for any scalar function $\zeta$ on $\Sigma_{t}$. The volume forms on
$\Sigma_{t}$ with respect to $\gamma$ and $G(t,\cdot)$ are written as
${\text{{vol}}}_{\gamma}$ and ${\text{{vol}}}_{G(t,\cdot)}$ (or just
${\text{{vol}}}_{G}$).

######  Definition 3.2 ($L^{2}$-norms).

Let $\mathfrak{T}$ be a $\Sigma_{t}$-tangent $(r,s)$-tensor field (for
$r,s\geq 0$). Then, we define:

(3.1a)
$\displaystyle\|\mathfrak{T}\|_{L^{2}_{\gamma}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{L^{2}_{\gamma}(\Sigma_{t})}^{2}$
$\displaystyle\,:=\int_{M}\lvert\mathfrak{T}(t,\cdot)\rvert_{\gamma}^{2}\,{\text{{vol}}}_{\gamma},$
(3.1b)
$\displaystyle\|\mathfrak{T}\|_{L^{2}_{G}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{L^{2}_{G}(\Sigma_{t})}^{2}$
$\displaystyle\,:=\int_{M}\lvert\mathfrak{T}(t,\cdot)\rvert_{G(t,\cdot)}^{2}\,{\text{{vol}}}_{G(t,\cdot)}\,$

######  Definition 3.3 (Sobolev norms).

Let $\mathfrak{T}$ be as above and $J\in\mathbb{N}_{0}$. We define:

(3.2a)
$\displaystyle\|\mathfrak{T}\|_{\dot{H}^{J}_{\gamma}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{\dot{H}^{J}_{\gamma}}^{2}=$
$\displaystyle\,\int_{\Sigma_{t}}\left\lvert\hat{\nabla}^{J}\mathfrak{T}(t,\cdot)\right\rvert_{\gamma}^{2}\,{\text{{vol}}}_{\gamma}$
(3.2b)
$\displaystyle\|\mathfrak{T}\|_{\dot{H}^{J}_{G}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{\dot{H}^{J}_{G}}^{2}=$
$\displaystyle\,\int_{\Sigma_{t}}\left\lvert\nabla^{J}\mathfrak{T}(t,\cdot)\right\rvert_{G(t,\cdot)}^{2}\,{\text{{vol}}}_{G(t,\cdot)}$
(3.2c)
$\displaystyle\|\mathfrak{T}\|_{H^{J}_{\gamma}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{H^{J}_{\gamma}}^{2}=$
$\displaystyle\,\sum_{k=0}^{J}\|\mathfrak{T}\|_{\dot{H}^{k}_{\gamma}(\Sigma_{t})}^{2}$
(3.2d)
$\displaystyle\|\mathfrak{T}\|_{H^{J}_{G}(\Sigma_{t})}^{2}=\|\mathfrak{T}(t,\cdot)\|_{H^{J}_{G}}^{2}=$
$\displaystyle\,\sum_{k=0}^{J}\|\mathfrak{T}\|_{\dot{H}^{k}_{G}(\Sigma_{t})}^{2}$

######  Definition 3.4 (Supremum norms).

For $\mathfrak{T}$ as above and $J\in\mathbb{N}_{0}$, we set:

(3.3a)
$\displaystyle\|\mathfrak{T}\|_{\dot{C}^{J}_{\gamma}(\Sigma_{t})}=\|\mathfrak{T}(t,\cdot)\|_{\dot{C}^{J}_{\gamma}}=$
$\displaystyle\,\sup_{p\in\Sigma_{t}}\left\lvert\hat{\nabla}^{J}\mathfrak{T}(t,\cdot)\right\rvert_{\gamma},$
$\displaystyle\|\mathfrak{T}\|_{C^{J}_{\gamma}(\Sigma_{t})}=\sum_{k=0}^{J}\|\mathfrak{T}\|_{\dot{C}^{k}_{\gamma}(\Sigma_{t})}$
(3.3b)
$\displaystyle\|\mathfrak{T}\|_{\dot{C}^{J}_{G}(\Sigma_{t})}=\|\mathfrak{T}(t,\cdot)\|_{\dot{C}^{J}_{G}}=$
$\displaystyle\,\sup_{p\in\Sigma_{t}}\left\lvert\nabla^{J}\mathfrak{T}(t,\cdot)\right\rvert_{G(t,\cdot)},$
$\displaystyle\|\mathfrak{T}\|_{C^{J}_{G}(\Sigma_{t})}=\sum_{k=0}^{J}\|\mathfrak{T}\|_{\dot{C}^{k}_{G}(\Sigma_{t})}$

######  Remark 3.5 (Time dependence is suppressed in notation).

When the choice of $t$ and $\Sigma_{t}$ is clear from context, we will often
drop time dependences of $G$, $\lvert\cdot\rvert_{G}$, ${\text{{vol}}}_{G}$
and $\mathfrak{T}$, surpress the hypersurface $\Sigma_{t}$ in the Sobolev and
supremum norms, and simply write $\int_{M}$ instead of $\int_{\Sigma_{t}}$.
For example, we write

$\|\mathfrak{T}\|_{L^{2}_{G}}^{2}=\int_{M}\lvert\mathfrak{T}\rvert_{G}^{2}\,{\text{{vol}}}_{G}\,.$

######  Definition 3.6 (Solution norms).

We define the following norms to measure the size of near-FLRW solutions:

(3.4a) $\displaystyle\mathcal{H}=$
$\displaystyle\,\|\Psi\|_{\color[rgb]{0,0,0}H^{18}_{G}\color[rgb]{0,0,0}}+\|\nabla\phi\|_{\color[rgb]{0,0,0}H^{17}_{G}\color[rgb]{0,0,0}}+a^{2}\|\nabla\phi\|_{\color[rgb]{0,0,0}\dot{H}^{18}_{G}\color[rgb]{0,0,0}}$
$\displaystyle\,+\|\Sigma\|_{H^{18}_{G}}+\|\bm{E}\|_{H^{18}_{G}}+\|\bm{B}\|_{H^{18}_{G}}$
$\displaystyle\,+\|G-\gamma\|_{H^{18}_{G}}+\|\text{{Ric}}[G]+\frac{2}{9}G\|_{H^{16}_{G}}+a^{-4}\|N\|_{H^{16}_{G}}+a^{-2}\|N\|_{\dot{H}^{17}_{G}}+\|N\|_{\dot{H}^{18}_{G}}$
(3.4b) $\displaystyle\color[rgb]{0,0,0}\mathcal{H}_{top}=\color[rgb]{0,0,0}$
$\displaystyle\,\color[rgb]{0,0,0}a^{2}\|\Psi\|_{{\dot{H}^{19}_{G}}}\color[rgb]{0,0,0}\color[rgb]{0,0,0}+a^{4}\|\nabla\phi\|_{\color[rgb]{0,0,0}\dot{H}^{19}_{G}\color[rgb]{0,0,0}}\color[rgb]{0,0,0}\color[rgb]{0,0,0}+a^{2}\|\Sigma\|_{\dot{H}^{19}_{G}}\color[rgb]{0,0,0}\color[rgb]{0,0,0}+a^{2}\|\text{{Ric}}[G]+\frac{2}{9}G\|_{\dot{H}^{17}_{G}}\color[rgb]{0,0,0}\,$
(3.4c) $\displaystyle\mathcal{C}=$
$\displaystyle\,\|\Psi\|_{\color[rgb]{0,0,0}C^{16}_{G}\color[rgb]{0,0,0}}+\|\nabla\phi\|_{\color[rgb]{0,0,0}C^{15}_{G}\color[rgb]{0,0,0}}+\|\Sigma\|_{C^{16}_{G}}+\|\bm{E}\|_{C^{16}_{G}}+\|\bm{B}\|_{C^{16}_{G}}$
$\displaystyle\,+\|G-\gamma\|_{C^{16}_{G}}+\|\text{{Ric}}[G]+\frac{2}{9}G\|_{C^{14}_{G}}+a^{-4}\|N\|_{C^{14}_{G}}+a^{-2}\|N\|_{\dot{C}^{15}_{G}}+\|N\|_{\dot{C}^{16}_{G}}$
(3.4d) $\displaystyle\mathcal{C}_{\gamma}=$
$\displaystyle\,\|\Psi\|_{\color[rgb]{0,0,0}C^{16}_{\gamma}\color[rgb]{0,0,0}}+\|\nabla\phi\|_{\color[rgb]{0,0,0}C^{15}_{\gamma}\color[rgb]{0,0,0}}+\|\Sigma\|_{C^{16}_{\gamma}}+\|\bm{E}\|_{C^{16}_{\gamma}}+\|\bm{B}\|_{C^{16}_{\gamma}}$
$\displaystyle\,+\|G-\gamma\|_{C^{16}_{\gamma}}+\|\text{{Ric}}[G]+\frac{2}{9}G\|_{C^{14}_{\gamma}}+a^{-4}\|N\|_{C^{14}_{\gamma}}+a^{-2}\|N\|_{\dot{C}^{15}_{\gamma}}+\|N\|_{\dot{C}^{16}_{\gamma}}$

######  Remark 3.7 (Choice of metric and controlling Christoffel symbols).

We could equivalently also phrase $\mathcal{H}$ in terms of $\gamma$-norms, or
predominately use $\mathcal{C}_{\gamma}$ instead of $\mathcal{C}$, since we
include the norms on $G-\gamma$ and
$\text{{Ric}}[G]+\frac{2}{9}G=(\text{{Ric}}[G]+\frac{2}{9}\gamma)+\frac{2}{9}(G-\gamma)$.
We will demonstrate in Lemma 7.4 how $H_{G}$ and $C_{\gamma}$ norms can be
used to control $H_{\gamma}$ and $C_{G}$ norms. We will also indicate how the
initial data and bootstrap assumptions for $\mathcal{C}_{\gamma}$ and
$\mathcal{C}$ are equivalent in Remarks 3.11 and 3.18. The main reason for
this is that, by successively replacing local coordinates in the expressions
of $\Gamma-\hat{\Gamma}$ by $\hat{\nabla}$, one has

(3.5)
$\|\Gamma-\hat{\Gamma}\|_{C^{l-1}_{\gamma}\color[rgb]{0,0,0}(M)\color[rgb]{0,0,0}}\lesssim
P_{l}(\|G-\gamma\|_{C^{l}_{\gamma}(M)},\|G^{-1}-\gamma^{-1}\|_{C^{l}_{\gamma}(M)})\,.$

We choose to work predominately with norms in terms of the rescaled metric
since quantities appearing in the Einstein system are naturally contracted by
$G$, not $\gamma$, and we commute with differential operators associated with
$G$.

######  Remark 3.8 (Redundancies in the solution norms).

The solution norms $\mathcal{H}$, $\mathcal{C}$ and $\mathcal{C}_{\gamma}$
aren’t “optimal” in the sense that controlling the norms of
$\Psi,\,\nabla\phi,\Sigma$ and $G-\gamma$ is entirely sufficient to gain the
claimed control (up to constant) on $N$ via the lapse equation, $\bm{E}$ and
$\bm{B}$ via to the constraint equations and $\text{{Ric}}[G]+\frac{2}{9}G$
via local coordinates. We choose to include all variables in the norms and
subsequent assumptions mainly for the sake of convenience.

###  3.2. Energies

The fundamental objects used to control the solution variables are the
energies that take the following form:

######  Definition 3.9 (Energies).

Let $l\in\mathbb{N}_{0}$. We define:

(3.6a) $\displaystyle\mathcal{E}^{(l)}(\phi,t)=$
$\displaystyle\,(-1)^{l}\int_{M}\Psi\Delta^{l}\Psi-a^{4}\phi\Delta^{l+1}\phi\,{\text{{vol}}}_{G}$
$\displaystyle=\begin{cases}\int_{M}\lvert\Delta^{\frac{l}{2}}\Psi\rvert^{2}+a^{4}\lvert\nabla\Delta^{\frac{l}{2}}\phi\rvert_{G}^{2}\,{\text{{vol}}}_{G}&l\text{
even}\\\
\int_{M}\lvert\nabla\Delta^{\frac{l-1}{2}}\Psi\rvert_{G}^{2}+a^{4}\lvert\Delta^{\frac{l+1}{2}}\phi\rvert_{G}^{2}\,{\text{{vol}}}_{G}&l\text{
odd}\end{cases}$ (3.6b) $\displaystyle\mathcal{E}^{(l)}(W,t)=$
$\displaystyle\,(-1)^{l}\int_{M}\langle\bm{E},\Delta^{l}\bm{E}\rangle_{G}+\langle\bm{B},\Delta^{l}\bm{B}\rangle_{G}\,{\text{{vol}}}_{G}$
(3.6c) $\displaystyle\mathcal{E}^{(l)}(\Sigma,t)=$
$\displaystyle\,(-1)^{l}\int_{M}\langle\Sigma,\Delta^{l}\Sigma\rangle_{G}\,{\text{{vol}}}_{G}$
(3.6d) $\displaystyle\mathcal{E}^{(l)}(\text{{Ric}},\cdot)=$
$\displaystyle\,(-1)^{l}\int_{M}\left\langle\text{{Ric}}[G]+\frac{2}{9}G,\Delta^{l}\left(\text{{Ric}}[G]+\frac{2}{9}G\right)\right\rangle_{G}\,{\text{{vol}}}_{G}$
(3.6e) $\displaystyle\mathcal{E}^{(l)}(N,\cdot)=$
$\displaystyle\,(-1)^{l}\int_{M}\langle
N,\Delta^{l}N\rangle_{G}\,{\text{{vol}}}_{G}$

Usually, we will use integration by parts to distribute derivatives within the
energies as in (3.6). Further, we introduce the notation

(3.7) $\displaystyle\mathcal{E}^{(\leq l)}$
$\displaystyle=\sum_{i=0}^{l}\mathcal{E}^{(i)}$

for any of the energies above.

For any $l\in\mathbb{N}_{0}$ and any smooth functions $f_{1},f_{2}\in
C^{\infty}(\mathbb{R}_{+})$, we have

(3.8) $f_{1}\cdot
f_{2}\cdot\mathcal{E}^{(2l+1)}\leq\frac{f_{1}^{2}}{2}\mathcal{E}^{(2l)}+\frac{f_{2}^{2}}{2}\mathcal{E}^{(2l+2)}\,.$

Performing the calculation for $\Sigma$ as an example, we have:

$\displaystyle\mathcal{E}^{(2l+1)}(\Sigma,\cdot)$
$\displaystyle\,=-\int_{M}\langle\Delta^{l}\Sigma,\Delta^{l+1}\Sigma\rangle_{G}\,{\text{{vol}}}_{G}\leq\int_{M}\lvert\Delta^{l}\Sigma\rvert_{G}\lvert\Delta^{l+1}\Sigma\rvert_{G}\,{\text{{vol}}}_{G}\leq\sqrt{\mathcal{E}^{(2l)}(\Sigma,\cdot)}\sqrt{\mathcal{E}^{(2l+2)}(\Sigma,\cdot)}$

Now, (3.8) then follows from the Young inequality. As a consequence, we also
have

(3.9) $\mathcal{E}^{(\leq
2l)}\lesssim\sum_{m=0}^{l}\mathcal{E}^{(2m)}\,,\,\mathcal{E}^{(\leq
2l+1)}\lesssim\sum_{m=0}^{l+1}\mathcal{E}^{(2m)}\,.$

This allows us to largely restrict our analysis to even order energies,
outside of how we close the bootstrap argument at top order.

###  3.3. Assumptions on the initial data

With the necessary solution norms and energies now defined, we can now state
what we assume near-FLRW initial data to satisfy:

######  Assumption 3.10 (Near-FLRW initial data).

For some small enough $\varepsilon\in(0,1)$ and the solution norms
$\mathcal{H}\color[rgb]{0,0,0},\mathcal{H}_{top}\color[rgb]{0,0,0}$ and
$\mathcal{C}$ as in Definition 3.6, we assume the rescaled inital data to be
close to that of the FLRW solution in Lemma 2.3 in the following sense:

(3.10)
$\mathcal{H}(t_{0})+\mathcal{H}_{top}(t_{0})+\mathcal{C}(t_{0})\lesssim\varepsilon^{2}$

The assumptions on
$\mathcal{H}\color[rgb]{0,0,0}+\mathcal{H}_{top}\color[rgb]{0,0,0}$ also imply
the following:

(3.11)
$\displaystyle\mathcal{E}^{(\leq\color[rgb]{0,0,0}18\color[rgb]{0,0,0})}(\phi,t_{0})+\mathcal{E}^{(\leq
18)}(\Sigma,t_{0})+\mathcal{E}^{(\leq 18)}(W,t_{0})+\mathcal{E}^{(\leq
16)}(\text{{Ric}},t_{0})$
$\displaystyle+\|\nabla\phi\|_{H^{18}_{G}}^{2}+\mathcal{E}^{(\color[rgb]{0,0,0}18\color[rgb]{0,0,0})}(N,t_{0})+a(t_{0})^{-4}\mathcal{E}^{(\color[rgb]{0,0,0}17\color[rgb]{0,0,0})}(N,t_{0})+a(t_{0})^{-8}\mathcal{E}^{(\color[rgb]{0,0,0}\leq
16\color[rgb]{0,0,0})}(N,t_{0})$
$\displaystyle\color[rgb]{0,0,0}+a(t_{0})^{4}\mathcal{E}^{(19)}(\phi,t_{0})+a(t_{0})^{4}\mathcal{E}^{(19)}(\Sigma,t_{0})+a(t_{0})^{4}\mathcal{E}^{(17)}(\text{{Ric}},t_{0})\color[rgb]{0,0,0}$
$\displaystyle\,\lesssim\varepsilon^{4}$

######  Remark 3.11 (Initial data size in $C_{\gamma}(t_{0})$).

Notice that by (3.5), (3.10) implies that

(3.12a)
$\|\Gamma-\hat{\Gamma}\|_{C^{15}_{G}(\color[rgb]{0,0,0}\Sigma_{t_{0}}\color[rgb]{0,0,0})}\lesssim\varepsilon^{4}\,,$
and arguing along similar lines and using $L^{2}-L^{\infty}$-estimates, also
(3.12b)
$\|\Gamma-\hat{\Gamma}\|_{H^{17}_{G}(\color[rgb]{0,0,0}\Sigma_{t_{0}}\color[rgb]{0,0,0})}\lesssim\varepsilon^{4}\,.$

In particular, since moving from $C_{\gamma}^{l}$ to $C_{G}^{l}$ only requires
control on Christoffel symbols to order $l-1$ for general tensors and $l-2$
for scalar functions, as well as zero order control on $G-\gamma$, it follows
from (3.10) that

(3.13) $C_{\gamma}(t_{0})\lesssim\varepsilon^{2}$

We refer to the proof of Lemma 7.4 for a more detailed term analysis and how a
similar argument also applies to the Sobolev norms.

######  Remark 3.12 (Redundancies in the initial data assumptions).

Similar to Remark 3.8, one could also reduce the initial data assumptions in
(3.10), especially at top order. In particular, we highlight that the Bel-
Robinson energy can be entirely controlled by the other terms that occur due
to the additional scaled $\Sigma$-energy at order $19$, or vice versa we could
drop the latter in favour of the former. This will be reflected in Lemma 6.10.

######  Remark 3.13 (The volume form).

Let $\mu_{G}$ and $\mu_{\gamma}$ denote the volume elements of $G$ and
$\gamma$ respectively. Since the determinant is a smooth map on invertible
matrices, the initial data assumptions on $G-\gamma$ also imply

(3.14)
$\|\mu_{G}-\mu_{\gamma}\|_{C^{0}_{G}(\Sigma_{t_{0}})}=\|\mu_{G}-\mu_{\gamma}\|_{C^{0}_{\gamma}(\Sigma_{t_{0}})}\lesssim\varepsilon^{2}\,.$

Consequently, we have

$\|{\text{{vol}}}_{G}-{\text{{vol}}}_{\gamma}\|_{C^{0}_{\gamma}(\Sigma_{t_{0}})}=\mu_{\gamma}^{-1}\|{\text{{vol}}}_{\gamma}\|_{C^{0}_{\gamma}(\Sigma_{t_{0}})}\|\mu_{G(t_{0},\cdot)}-\mu_{\gamma}\|_{C^{0}_{\gamma}(\Sigma_{t_{0}})}\lesssim\varepsilon^{2}$

and, since
$\|G^{-1}-\gamma^{-1}\|_{C^{0}_{\gamma}(\Sigma_{t_{0}})}\lesssim\varepsilon^{2}$
also follows by a von Neumann series argument from the initial data assumption
on $G-\gamma$,

(3.15)
$\|{\text{{vol}}}_{G}-{\text{{vol}}}_{\gamma}\|_{C^{0}_{G}(\Sigma_{t_{0}})}\lesssim\varepsilon^{2}\,.$

###  3.4. Local well-posedness and continuation criteria

For everything that follows, we need to establish that the initial data
assumptions above also ensure local well-posedness. For the core system, we
translate the local well-posedness result for stiff fluids in [RS18b] to the
sub-case of the scalar field system. While statement and proof there are for
vanishing spatial sectional curvature and what corresponds to choosing
$C=\sqrt{\nicefrac{{2}}{{3}}}$, the arguments for our setting are completely
analogous.

######  Lemma 3.14 (Local well-posedness and continuation criteria for the
Einstein scalar-field system (Big Bang version), see [RS18b, Theorem 14.1]).

Let $N\geq 4$ be an integer and
$(M,\mathring{g},\mathring{k},\mathring{\pi},\mathring{\psi})$ be geometric
initial data to the Einstein scalar-field system (see Section 1.2.1) and
assume that one has

$\|\mathring{g}-a(t_{0})^{2}\gamma\|_{H^{N+1}_{\gamma}(M)}+\|\mathring{k}+\frac{\tau(t_{0})}{3}\cdot
a(t_{0})^{2}\gamma\|_{H^{N}_{\gamma}(M)}+\|\mathring{\pi}\|_{H^{N}_{\gamma}(M)}+\|\mathring{\psi}-Ca(t_{0})^{-3}\|_{H^{N}_{\gamma}(M)}<\infty\,$

as well as, for some sufficiently small $\eta^{\prime}>0$,

$\|\mathring{\psi}-Ca(t_{0})^{-3}\|_{C^{0}_{\gamma}(M)}\leq\eta^{\prime}\,.$

Then, the CMC-transported Einstein scalar-field system (respectively the
rescaled system) is locally well-posed in the following sense: The initial
data $(\mathring{g},\mathring{k},\mathring{\pi},\mathring{\psi})$ launches a
unique classical solution $(g,k,n,\nabla\phi,\partial_{t}\phi)$ to
(2.16a)-(2.16b), (2.15a)-(2.15b), (2.18) and (2.17a) on $[t_{1},t_{0}]\times
M$ for some $t_{1}\in(0,t_{0})$ that satisfies
${k^{l}}_{l}=-3\frac{\dot{a}}{a}$ and $n>0$, launches a solution to the
Einstein scalar-field system and such that the variables enjoy the following
regularity:

$\displaystyle g\in$ $\displaystyle
C^{N-1}_{dt^{2}+\gamma}([t_{1},t_{0}]\times M)\cap
C^{0}([t_{1},t_{0}],H_{\gamma}^{N+1}(M))$ $\displaystyle k\in$ $\displaystyle
C^{N-2}_{dt^{2}+\gamma}([t_{1},t_{0}]\times M)\cap
C^{0}([t_{1},t_{0}],H_{\gamma}^{N}(M))$ $\displaystyle\nabla\phi\in$
$\displaystyle C^{N-2}_{dt^{2}+\gamma}([t_{1},t_{0}]\times M)\cap
C^{0}([t_{1},t_{0}],H_{\gamma}^{N}(M))$ $\displaystyle\partial_{t}\phi\in$
$\displaystyle C^{N-2}_{dt^{2}+\gamma}([t_{1},t_{0}]\times M)\cap
C^{0}([t_{1},t_{0}],H_{\gamma}^{N}(M))$ $\displaystyle n\in$
$\displaystyle\,C^{N}_{dt^{2}+\gamma}([t_{1},t_{0}]\times M)\cap
C^{0}([t_{1},t_{0}],H_{\gamma}^{N+2}(M))$

The rescaled variables $(G,\Sigma,N,\nabla\phi,\Psi)$ enjoy the analogous
regularity. If $(\mathfrak{t},t_{0}]$ is the maximal interval on which the
above statements hold, then one either has $\mathfrak{t}=0$ or one of the
following blow-up criteria are satisfied:

  1. (1)

The smallest eigenvalue of $g(t_{m},\cdot)$ converges to 0 for some sequence
$(t_{m},x_{m})\subseteq(\mathfrak{t},t_{0}]\times M$ with
$t_{m}\downarrow\mathfrak{t}$.

  2. (2)

$n(t_{m},x_{m})$ converges to $0$ for some sequence
$(t_{m},x_{m})\subseteq(\mathfrak{t},t_{0}]\times M$ with
$t_{m}\downarrow\mathfrak{t}$.

  3. (3)

$\left(\lvert\partial_{0}\phi\rvert^{2}+\lvert\nabla\phi\rvert_{g}^{2}\right)(t_{m},x_{m})$
converges to $0$ for some sequence
$(t_{m},x_{m})\subseteq(\mathfrak{t},t_{0}]\times M$ with
$t_{m}\downarrow\mathfrak{t}$.

  4. (4)

$s\in(\mathfrak{t},t_{0}]\mapsto\|g\|_{C^{2}_{\gamma}(\Sigma_{s})}+\|k\|_{C^{1}_{\gamma}(\Sigma_{s})}+\|n\|_{C^{2}_{\gamma}(\Sigma_{s})}+\|\partial_{t}\phi\|_{C^{1}_{\gamma}(\Sigma_{s})}+\|\nabla\phi\|_{C^{1}_{\gamma}(\Sigma_{s})}$
is unbounded.

###### A note on the proof.

Note that the additional initial data requirement in the stiff-fluid setting
that the pressure is strictly positive is covered by the smallness assumption
on $\mathring{\psi}-Ca(t_{0})^{2}$, since the pressure corresponds to
$\lvert\mathring{\psi}\rvert^{2}+\lvert\mathring{\pi}\rvert_{\mathring{g}}^{2}$
and the assumptions on $\partial_{t}\phi$ and $\nabla\phi$ ensure that (after
embedding) this quantity behaves like
$C^{2}a(t_{0})^{-6}+\mathcal{O}\left(\eta^{\prime}\right)$ at
$\Sigma_{t_{0}}$. ∎

######  Corollary 3.15 (Local well-posedness for the Bel-Robinson variables).

Under the assumptions of Lemma 3.14, the Bel-Robinson variables $E$ and $B$
corresponding to the Lorentzian metric $\overline{g}=-n^{2}dt^{2}+g$ satisfy
the equations (2.24a)-(2.24b), are the unique classical solutions to the
evolution equations (2.25a)-(2.25b), and satisfy

$E,B\in C^{N-3}([t_{1},t_{0}]\times M)\cap
C([t_{1},t_{0}],H^{N-1}_{\gamma}(M))$

###### Proof.

That $E$ and $B$ satisfy the constraint equations, solve the evolution
equations and have the stated regularity on the interval of existence is a
direct consequence of Lemma 3.14 and the computations in Section 2.4.
Furthermore, with initial data derived from the constraint equations as in
Remark 2.8, the hyperbolic system (2.24a)-(2.24b) launches a unique solution
satisfying the regularity above that must then be $(E,B)$. ∎

For sufficiently regular initial data ($N\geq 21$), it follows that

$\mathcal{E}^{(\leq 19)}(\phi,\cdot),\mathcal{E}^{(\leq
18)}(W,\cdot),\mathcal{E}^{(\leq 19)}(\Sigma,\cdot),\mathcal{E}^{(\leq
17)}(\text{{Ric}},\cdot),\|G-\gamma\|_{H^{18}_{G}}\in C^{1}([t_{1},t_{0}])\,,$

and similarly the square of any supremum norm occuring in $\mathcal{C}$ is
continuously differentiable on $[t_{1},t_{0}]$. Strictly speaking, we would
need to assume this additional regularity on our initial data for the
computations in the following sections (especially Section 6) to hold.
However, since smooth functions are dense in $H^{l}(M)$ for any
$l\in\mathbb{N}_{0}$, any bounds on $\mathcal{H}(t)$ and $\mathcal{C}(t)$ that
we prove assuming sufficient regularity at $\Sigma_{t_{0}}$ then immediately
extend to data only satisfying the regularity implied by (3.10).

Thus, from here on out, we will assume without loss of genearlity that all
energies and squared norms are continuously differentiable on the domain of
existence, and similarly all variables are continuously differentiable for the
lower order $C_{G}$-norm improvements in Section 4.2.

###  3.5. Bootstrap assumption

To keep an overview of the entire bootstrap argument, we state all of the
assumptions and comprehensively list how we intend to improve them.

######  Assumption 3.16 (Bootstrap assumption).

Fix some $t_{Boot}\in[0,t_{0})$. Further, let $c_{0}>0$, let
$\sigma\in(\varepsilon^{\frac{1}{8}},1]$ be suitably small such that
$c_{0}\sigma<1$, and $K_{0}>0$ a suitable constant. For any
$t\in(t_{Boot},t_{0}]$, we assume

(3.16) $\mathcal{C}(t)\leq K_{0}\varepsilon a(t)^{-c_{0}\sigma}\,.$

######  Remark 3.17.

More explicitly, (3.16) means

(3.17a)
$\displaystyle\|\Psi\|_{\color[rgb]{0,0,0}C^{16}_{G}\color[rgb]{0,0,0}}\leq$
$\displaystyle\,K_{0}\varepsilon a^{-c_{0}\sigma}$ (3.17b)
$\displaystyle\|\nabla\phi\|_{\color[rgb]{0,0,0}C^{15}_{G}\color[rgb]{0,0,0}}\leq$
$\displaystyle\,K_{0}\varepsilon a^{-c_{0}\sigma}$ (3.17c)
$\displaystyle\|\Sigma\|_{C^{16}_{G}}\leq$ $\displaystyle\,K_{0}\varepsilon
a^{-c_{0}\sigma}$ (3.17d) $\displaystyle\|\bm{E}\|_{C^{16}_{G}}\leq$
$\displaystyle\,K_{0}\varepsilon a^{-c_{0}\sigma}$ (3.17e)
$\displaystyle\|\bm{B}\|_{C^{16}_{G}}\leq$ $\displaystyle\,K_{0}\varepsilon
a^{-c_{0}\sigma}$ (3.17f)
$\displaystyle\|\text{{Ric}}[G]+\frac{2}{9}G\|_{C^{14}_{G}}\leq$
$\displaystyle\,K_{0}\varepsilon a^{-c_{0}\sigma}$ (3.17g)
$\displaystyle\|G-\gamma\|_{C^{16}_{G}}\leq$ $\displaystyle\,K_{0}\varepsilon
a^{-c_{0}\sigma}$ (3.17h)
$\displaystyle\|N\|_{C^{14}_{G}}+a^{2}\|N\|_{\dot{C}^{15}_{G}}+a^{4}\|N\|_{\dot{C}^{16}_{G}}\leq$
$\displaystyle\,K_{0}\varepsilon a^{4-c_{0}\sigma}$ (3.17i)
$\displaystyle\color[rgb]{0,0,0}\|\Gamma-\hat{\Gamma}\|_{C^{15}_{G}}\leq\color[rgb]{0,0,0}$
$\displaystyle\color[rgb]{0,0,0}\,K_{0}\varepsilon
a^{-c_{0}\sigma}\color[rgb]{0,0,0}$

######  Remark 3.18 (Bootstrap assumptions with respect to $\gamma$).

Note again that we could equivalently make the above bootstrap assumptions
with respect to $H_{\gamma}$\- and $C_{\gamma}$-norms: For example, the
assumptions (3.17i) and (3.17g) imply

$\displaystyle\|\zeta\|_{C^{l}_{\gamma}}\lesssim$
$\displaystyle\,a^{-c\sigma}\|\zeta\|_{C^{l}_{G}}+\|\zeta\|_{C^{\lceil\frac{l-1}{2}\rceil}_{\gamma}}\varepsilon
a^{-c\sigma},\quad\|\mathfrak{T}\|_{C^{l}_{\gamma}}\lesssim
a^{-c\sigma}\|\mathfrak{T}\|_{C^{l}_{G}}+\|\mathfrak{T}\|_{C^{\lceil\frac{l}{2}\rceil}_{\gamma}}\varepsilon
a^{-c\sigma}$

for any smooth function $\zeta\in C^{\infty}(\Sigma_{t})$, any
$\Sigma_{t}$-tangent tensor $\mathfrak{T}$ and a constant $c>0$. This is
essentially a direct consequence of (3.5), and we will prove an improved
version of this rigorously in Lemma 7.4. Applying this to each norm in
$\mathcal{C}$, we get

(3.18) $\mathcal{C}_{\gamma}\lesssim\varepsilon a^{-c\sigma}$

for some updated constant $c\geq c_{0}$.

######  Remark 3.19 (Strategy for the bootstrap improvement).

Our goal is to improve the $C$-norm estimate to

$\mathcal{C}\leq
K_{1}\varepsilon^{\frac{9}{8}}a^{-c_{1}\varepsilon^{\frac{1}{8}}}\,,$

where $c_{1},K_{1}>0$ are positive constants independent of $\sigma$ and
$\varepsilon$ . Notice how this is actually an improvement if we choose
$\sigma$ suitably and then choose $\varepsilon$ sufficiently small: Any update
between $K_{0}$ and $K_{1}$ can be balanced out since we gain at least the
additional prefactor $\varepsilon^{\frac{1}{8}}$ in each estimate, which we
can then choose to have been suitably small. Similarly, we improve the power
of $a$ if we have
$\varepsilon^{\frac{1}{8}}\cdot\sigma^{-1}<\frac{c_{0}}{c_{1}}$. If we then
retroactively choose $\sigma$ large enough compared to $\varepsilon$ but small
overall – for example $\sigma=\varepsilon^{\frac{1}{16}}$ – and then ensure
that $\max\\{c_{0},c_{1}\\}\varepsilon^{\frac{1}{16}}<1$ as well as
$c_{1}\varepsilon^{\frac{1}{16}}<c_{0}$ are satisfied by choosing
$\varepsilon$ to have been small enough, we have strictly improved the
bootstrap assumptions.

######  Remark 3.20 (Conventions within the bootstrap argument).

Throughout the rest of the argument, we tacitly assume $t\in(t_{Boot},t_{0}]$
if not stated otherwise, and we assume $\varepsilon$ and $\sigma$ to be
sufficiently small. In the proof of Theorem 8.2, we will choose
$\sigma=\varepsilon^{\frac{1}{16}}$, but this explicit choice will not be used
or needed up to that point. Finally, we allow $c\geq c_{0}$ be a constant that
we may update from line to line, and will similarly deal with prefactors by
“$\lesssim$”-notation where the constant may change in each line. These
updates will always be independent of $\sigma$ and $\varepsilon$, but may
depend on $t_{0}$, and the quantities arising from the FLRW reference
solution. Hence, we not only assume $c_{0}\sigma<1$, but $c\sigma<1$
throughout the argument.

##  4\. Big Bang stability: A priori estimates

In this section, we collect strong low order $C_{G}$-norm estimates that
follow as an immediate consequence from the bootstrap assumptions, starting
with key estimates at the base level and followed by weaker, but still
improved estimates at higher levels. Finally, we collect a differentiation
formula for integrals with respect to ${\text{{vol}}}_{G}$ as well as a
Sobolev estimate that lays the groundwork for energy coercivity. In
particular, using the strong $C_{G}$-norm estimates, said estimate proves that
moving between energies and norms at most an error involving lower order
energies of the controlled variable and curvature energies, scaled by
$a^{-c\sqrt{\varepsilon}}$.

###  4.1. Strong $C^{0}_{G}$-estimates

First, we establish a pointwise bound on the lapse that actually holds
irrespective of the bootstrap assumptions:

######  Lemma 4.1 (Maximum principle for the lapse).

The lapse remains positive and bounded throughout the evolution:

(4.1) $n=N+1\in(0,3]$

###### Proof.

Let $t\in\mathbb{R}_{+}$ be arbitrary and let $n_{min}$ be the minimum of $n$
over $\Sigma_{t}$ at $(t,x_{min})$. Then, $(\Delta_{g}n)(t,x_{min})>0$ holds.
If $n_{min}$ were nonpositive, (2.17a) would lead to the following
contradiction:

$0\geq-12\pi C^{2}a^{-6}-\frac{1}{3}a^{-2}+n_{min}\left[\frac{1}{3}a^{-2}+4\pi
C^{2}a^{-6}+\langle\hat{k},\hat{k}\rangle_{g}+8\pi\lvert\partial_{0}\phi\rvert^{2}\right]=\Delta_{g}n(t,x_{min})>0$

This shows $n>0$, and the upper bound follows analogously. ∎

The following estimate will be essential in dealing with borderline terms
throughout the bootstrap argument:

######  Lemma 4.2 (Strong $C^{0}_{G}$ estimates).

The following estimates hold:

(4.2a) $\displaystyle\|\Psi\|_{C^{0}_{G}}\lesssim$
$\displaystyle\,\varepsilon$ (4.2b)
$\displaystyle\|\Sigma\|_{C^{0}_{G}}\lesssim$ $\displaystyle\,\varepsilon$
(4.2c) $\displaystyle\|\bm{E}\|_{C^{0}_{G}}\lesssim$
$\displaystyle\,\varepsilon$

###### Proof.

(4.2a): From (2.32a), we obtain the following using Lemma 4.1 for $n$, the
bootstrap assumptions (3.17h) and (3.17b) and that $\dot{a}\simeq a^{2}$ by
(2.3):

$\displaystyle\lvert\partial_{t}\Psi\rvert\lesssim$
$\displaystyle\,\varepsilon a^{5-c\sigma}+\varepsilon
a^{1-c\sigma}+\varepsilon
a^{1-c\sigma}\lvert\Psi\rvert+{\varepsilon}a^{1-c\sigma}$

After integration, we thus obtain using the initial data assumption (3.10):

$\displaystyle\lvert\Psi(t)\rvert\lesssim$
$\displaystyle\,\lvert\Psi(t_{0})\rvert+\int_{t}^{t_{0}}\varepsilon
a(s)^{1-c\sigma}\,ds+\int_{t}^{t_{0}}{\varepsilon}a(s)^{1-c\sigma}\lvert\Psi(s)\rvert\,ds$
$\displaystyle\lesssim$
$\displaystyle\,\varepsilon\left(1+\int_{t}^{t_{0}}a(s)^{1-c\sigma}\,ds\right)+\int_{t}^{t_{0}}{\varepsilon}a(s)^{1-c\sigma}\lvert\Psi(s)\rvert\,ds$

By (2.6), the integral over $a^{1-c\sigma}$ is bounded since $c\sigma<1$, so
the Gronwall lemma now yields (4.2a).

(4.2b): Notice that

(4.3)
$\partial_{t}\lvert\Sigma\rvert_{G}^{2}={(\partial_{t}\Sigma^{\sharp})^{l}}_{m}{(\Sigma^{\sharp})^{m}}_{l}+{(\Sigma^{\sharp})^{l}}_{m}{(\partial_{t}\Sigma^{\sharp})^{m}}_{l}=2{(\partial_{t}\Sigma^{\sharp})^{l}}_{m}{(\Sigma^{\sharp})^{m}}_{l}\leq
2\lvert\partial_{t}\Sigma^{\sharp}\rvert_{G}\lvert\Sigma\rvert_{G}\,.$

Now, we consider (2.28d) and, using the bootstrap assumptions (3.17h), (3.17f)
and (3.17b), get:

$\displaystyle\lvert\partial_{t}\Sigma^{\sharp}\rvert_{G}\lesssim$
$\displaystyle\,\tau\lvert
N\rvert\lvert\Sigma^{\sharp}\rvert_{G}+\lvert\nabla^{\sharp}\nabla
N\rvert_{G}a+\lvert N+1\rvert
a\left\lvert\text{{Ric}}[G]^{\sharp}+\frac{2}{9}G^{\sharp}\right\rvert_{G}+\lvert
N+1\rvert a\lvert\nabla^{\sharp}\phi\nabla\phi\rvert_{G}$
$\displaystyle+\sqrt{3}\lvert N\rvert\cdot\left(4\pi
C^{2}a^{-3}+\frac{1}{9}a\right)$ $\displaystyle\lesssim$
$\displaystyle\,{\varepsilon}a^{1-c\sigma}\lvert\Sigma\rvert_{G}+{\varepsilon}a^{1-c\sigma}$

We can now apply Lemma 11.2 with $f=\lvert\Sigma\rvert_{G}^{2}$, and thus have
along with (3.10) and (4.3):

$\lvert\Sigma\rvert_{G}(t)\leq\lvert\Sigma\rvert_{G}(t_{0})+\int_{t}^{t_{0}}\lvert\partial_{t}\Sigma^{\sharp}\rvert_{G}(s)\,ds\lesssim\varepsilon$

(4.2c): Using the constraint equation (2.29c) and that $\langle
G,\bm{E}\rangle_{G}=\text{tr}_{G}\bm{E}=0$, one sees

$\lvert\bm{E}\rvert_{G}^{2}=\left\langle
a^{4}\left(\text{{Ric}}[G]+\frac{2}{9}G\right)-\dot{a}a^{2}\Sigma-\Sigma\odot_{G}\Sigma-4\pi
a^{4}\nabla\phi\nabla\phi,\bm{E}\right\rangle_{G}\,.$

Then, applying the bootstrap assumptions (3.17f) and (3.17b) shows the Ricci
and matter terms are bounded by $\varepsilon
a^{4-c\sigma}\lvert\bm{E}\rvert_{G}$, and the a priori estimate (4.2b) along
with $\dot{a}a^{2}\simeq 1$ by (2.3) bounds the remaining terms by
$\varepsilon\lvert\bm{E}\rvert_{G}$. The statement then follows by dividing by
$\lvert\bm{E}\rvert_{G}$ and taking the supremum. ∎

Note that, in the proof of (4.2b), it was essential that we used (2.28d)
instead of (2.28c), since using the latter would incur terms of the type
$\lvert\partial_{t}G\rvert_{G}\lvert\Sigma\rvert_{G}^{2}$and
$a^{-3}\lvert\Sigma\rvert_{G}^{3}$ when computing the time derivative of
$\lvert\Sigma\rvert_{G}^{2}$, which, at this point, behave like $\varepsilon
a^{-3-c\sigma}\lvert\Sigma\rvert_{G}^{2}$, and thus not yield the sharp
estimate (or even an improved estimate) that we will need to control
borderline terms.

###  4.2. Strong low order $C_{G}$-norm estimates

Now, we can prove the main supremum norm estimates in this section:

######  Lemma 4.3 (Strong low order $C_{G}$-norm estimates).

The following estimates hold:

(4.4a) $\displaystyle\|\Psi\|_{C^{13}_{G}}\lesssim$
$\displaystyle\,\varepsilon a^{-c\sqrt{\varepsilon}}\,$ (4.4b)
$\displaystyle\|\Sigma\|_{C^{12}_{G}}\lesssim$ $\displaystyle\,\varepsilon
a^{-c\sqrt{\varepsilon}}\,$ (4.4c)
$\displaystyle\|G-\gamma\|_{C^{12}_{G}}\lesssim$
$\displaystyle\sqrt{\varepsilon}a^{-c\sqrt{\varepsilon}}\,$ (4.4d)
$\displaystyle\|G^{-1}-\gamma^{-1}\|_{C^{12}_{G}}\lesssim$
$\displaystyle\sqrt{\varepsilon}a^{-c\sqrt{\varepsilon}}\,$ (4.4e)
$\displaystyle\|\nabla\phi\|_{C^{12}_{G}}\lesssim$
$\displaystyle\,\sqrt{\varepsilon}a^{-c\sqrt{\varepsilon}}$ (4.4f)
$\displaystyle\|\text{{Ric}}[G]+\frac{2}{9}G\|_{C^{10}_{G}}\lesssim$
$\displaystyle\,\sqrt{\varepsilon}a^{-c\sqrt{\varepsilon}}$ (4.4g)
$\displaystyle\|\bm{B}\|_{C^{11}_{G}}\lesssim$
$\displaystyle\,\varepsilon\color[rgb]{0,0,0}a^{2-c\sqrt{\varepsilon}}\color[rgb]{0,0,0}$
(4.4h) $\displaystyle\|\bm{E}\|_{C^{12}_{G}}\lesssim$
$\displaystyle\,\varepsilon a^{-c\sqrt{\varepsilon}}$

###### Proof.

Before going into the individual estimates, we collect the following
commutator term estimates from the expressions in (11.10a)-(11.10b):

(4.5) $\displaystyle\|[\partial_{t},\nabla^{J}]\zeta\|_{C^{0}_{G}}\lesssim$
$\displaystyle\,a^{-3}\|N+1\|_{C^{J-1}_{G}}\|\Sigma\|_{C^{J-1}_{G}}\|\zeta\|_{C^{J-1}_{G}}+\frac{\dot{a}}{a}\|N\|_{C^{J-1}_{G}}\|\zeta\|_{C^{J-1}_{G}}$
(4.6)
$\displaystyle\|[\partial_{t},\nabla^{J}]\mathfrak{T}\|_{C^{0}_{G}}\lesssim$
$\displaystyle\,a^{-3}\|N+1\|_{C^{J}_{G}}\left(\|\nabla^{J}\Sigma\|_{C^{0}_{G}}\|\mathfrak{T}\|_{C^{0}_{G}}+\|\Sigma\|_{C^{J-1}_{G}}\|\mathfrak{T}\|_{C^{J-1}_{G}}\right)+\frac{\dot{a}}{a}\|N\|_{C^{J}_{G}}\|\mathfrak{T}\|_{C^{J-1}_{G}}$

With this in hand, we will prove each estimate by iterating over the
derivative order as long as the bootstrap assumptions can be applied. In each
step, we use the previously obtained estimates at lower order to control the
commutator term (with some additional care for $\mathfrak{T}=\Sigma$ which we
need to consider first), while we can use similar arguments to those at order
$0$ to control the “core” of the evolution equations.

To start out, we apply (4.2b) on $\Sigma$ and the bootstrap assumption (3.17h)
on $N$ to the rescaled evolution equations (2.28a)-(2.28b) and deduce

(4.7) $\lvert\partial_{t}G^{\pm 1}\rvert_{G}=\lvert\partial_{t}(G^{\pm
1}-\gamma^{\pm 1})\rvert_{G}\lesssim\varepsilon a^{-3}+\varepsilon
a^{1-c\sigma}\lesssim\varepsilon a^{-3}\,.$

(4.4b): We assume

(4.8) $\|\Sigma\|_{{C}^{J-1}_{G}}\lesssim\varepsilon a^{-c\sqrt{\varepsilon}}$

to be satisfied for some $J\in\\{1,\dots,12\\}$ (For $J=1$, this is true by
(4.2b)). Observe the following:

$\partial_{t}\lvert\nabla^{J}\Sigma\rvert_{G}^{2}=2\langle\partial_{t}\nabla^{J}\Sigma,\nabla^{J}\Sigma\rangle_{G}+\partial_{t}G^{-1}\ast\nabla^{J}\Sigma\ast\nabla^{J}\Sigma$

Now, we commute (2.28d) with $\nabla^{J}$: As before,
$\nabla^{J}\partial_{t}\Sigma$ is bounded by $\varepsilon a^{-c\sigma}$ for
any admissible $J$. Hence and using (4.7),

$\partial_{t}\lvert\nabla^{J}\Sigma\rvert_{G}^{2}\lesssim\varepsilon
a^{-3}\lvert\nabla^{J}\Sigma\rvert_{G}^{2}+\left(\varepsilon
a^{1-c\sigma}+\|[\partial_{t},\nabla^{J}]\Sigma\|_{C^{0}_{G}}\right)\lvert\nabla^{J}\Sigma\rvert_{G}$

is satisfied. Looking at the commutator term using (4.6), we have with (4.2b)
that

$\|[\partial_{t},\nabla^{J}]\Sigma\|_{C^{0}_{G}}\lesssim\varepsilon
a^{-3}\|\Sigma\|_{\dot{C}^{J}_{G}}+a^{-3}\cdot\|\Sigma\|_{C^{J-1}_{G}}^{2}+\varepsilon
a^{1-c\sigma}\|\Sigma\|_{C^{J-1}_{G}}\,.$

Altogether, we obtain

$\partial_{t}\lvert\nabla^{J}\Sigma\rvert_{G}^{2}\lesssim\left(\varepsilon
a^{-3}\|\Sigma\|_{\dot{C}^{J}_{G}}+\varepsilon
a^{-c\sigma}+\varepsilon^{2}a^{-3-c\sqrt{\varepsilon}}\right)\lvert\nabla^{J}\Sigma\rvert_{G}\,.$

With Lemma 11.2 as well as the initial data assumption (3.10) and the integral
formula (2.6) with $p=c\sqrt{\varepsilon}$, this implies

$\lvert\nabla^{J}\Sigma\rvert_{G}(t)\lesssim\,\int_{t}^{t_{0}}\varepsilon
a^{-3}\|\Sigma\|_{\dot{C}^{J}_{G}(\Sigma_{s})}\,ds+\varepsilon\left(1+\sqrt{\varepsilon}a^{-c\sqrt{\varepsilon}}\right)$

and consequently, after taking the supremum on the left and applying the
Gronwall lemma,

$\|\Sigma\|_{\dot{C}^{J}_{G}(\Sigma_{s})}\lesssim\varepsilon
a^{-c\sqrt{\varepsilon}}\,.$

Combining this with (4.8) proves the statement up to order $J$, and hence
shows (4.4b) by iterating the argument up to $J=12$.

(4.4a): We again assume that

(4.9) $\|\Psi\|_{C^{J-1}_{G}}\lesssim\varepsilon a^{-c\sqrt{\varepsilon}}$

holds for $J\in\\{1,2,\dots,13\\}$. Observe that

$\lvert\partial_{t}\nabla^{J}\Psi\rvert_{G}\lesssim
a\|N+1\|_{C^{J+1}_{G}}\|\nabla\phi\|_{C^{J+1}_{G}}+\frac{\dot{a}}{a}\|\nabla
N\|_{C^{J}_{G}}(1+\|\Psi\|_{C^{J}_{G}})+\|[\partial_{t},\nabla^{J}]\Psi\|_{C^{0}_{G}}$

By (3.17h), (3.17b) and (3.17a), the first two summands can be bounded (up to
constant) by $\varepsilon a^{1-c\sigma}$. By (4.9), (4.4b) and (3.17h) and
using (4.5), the commutator term is bounded (up to constant) by
$\varepsilon^{2}a^{-3-c\sqrt{\varepsilon}}$. Altogether,

$\lvert\partial_{t}\nabla^{J}\Psi\rvert_{G}\lesssim\varepsilon
a^{1-c\sigma}+\varepsilon^{2}a^{-3-c\sqrt{\varepsilon}}$

follows. Inserting this and (4.7) into

$\lvert\partial_{t}\left(\lvert\nabla^{J}\Psi\rvert_{G}^{2}\right)\rvert\leq\lvert\partial_{t}G^{-1}\rvert_{G}\lvert\nabla^{J}\Psi\rvert_{G}^{2}+2\lvert\partial_{t}\nabla^{J}\Psi\rvert_{G}\cdot\lvert\nabla^{J}\Psi\rvert_{G}$

implies, with Lemma 11.2,

$\displaystyle\lvert\nabla^{J}\Psi\rvert_{G}(t)\leq$
$\displaystyle\,\lvert\nabla^{J}\Psi\rvert(t_{0})+\int_{t}^{t_{0}}\left(\frac{1}{2}\lvert\partial_{t}G^{-1}\rvert\lvert\nabla^{J}\Psi\rvert_{G}+\lvert\partial_{t}\nabla^{J}\Psi\rvert_{G}\right)(s)\,ds$
$\displaystyle\lesssim$
$\displaystyle\,\varepsilon^{2}+\int_{t}^{t_{0}}\left(\varepsilon
a(s)^{-3}\lvert\nabla^{J}\Psi(s,\cdot)\rvert_{G}+\varepsilon
a(s)^{1-c\sigma}+\varepsilon^{2}a(s)^{-3-c\sqrt{\varepsilon}}\right)\,ds\,.$

We obtain using (2.6):

$\lvert\nabla^{J}\Psi\rvert_{G}(t)\lesssim\varepsilon
a(t)^{-c\sqrt{\varepsilon}}+\int_{t}^{t_{0}}\varepsilon
a(s)^{-3}\lvert\nabla^{J}\Psi(s,\cdot)\rvert_{G}\,ds$

The Gronwall lemma, applying (2.5) and taking the supremum over $\Sigma_{t}$
then implies $\lvert\nabla^{J}\Psi\rvert_{\dot{C}^{J}_{G}}\lesssim\varepsilon
a^{-c\sqrt{\varepsilon}}$. This proves (4.2a) by iterating over $J$ and adding
up the individual seminorms.

(4.4c)-(4.4d): Note that (4.7) implies (4.4c) at order $0$ since one has

$\lvert\partial_{t}(\lvert
G-\gamma\rvert_{G})^{2}\rvert\lesssim\lvert\partial_{t}G^{-1}\rvert_{G}\lvert
G-\gamma\rvert_{G}^{2}+\lvert\partial_{t}(G-\gamma)\rvert_{G}\lvert
G-\gamma\rvert_{G}\lesssim\varepsilon a^{-3}(1+\lvert
G-\gamma\rvert_{G})\lvert G-\gamma\rvert_{G}$

which we can apply the Gronwall lemma to after integrating, along with (2.7)