<html><head><meta name="color-scheme" content="light dark"></head><body><pre style="word-wrap: break-word; white-space: pre-wrap;">\makeatletter
\@ifundefined{HCode}                  
{\PassOptionsToPackage{russian,english}{babel}
\documentclass[screen,CRMECA,Unicode,biblatex]{cedram}
\addbibresource{crmeca20240892.bib}
\newenvironment{noXML}{}{}%
\def\tsup#1{$^{{#1}}$}
\def\tsub#1{$_{{#1}}$}
\def\sfrac#1#2{{#1}/{#2}}
\def\thead{\noalign{\relax}\hline}
\def\endthead{\noalign{\relax}\hline}
\def\tbody{\noalign{\relax}\hline}
\def\xtbody{\noalign{\relax}}
\def\tabnote#1{\vskip4pt\parbox{.51\linewidth}{#1}}
\def\xxtabnote#1{\vskip4pt\parbox{.65\linewidth}{#1}}
\def\mathbi#1{\boldsymbol{#1}}
\def\xmathop#1{\mathop{\text{#1}}}
\def\rmDelta{\Delta}
\usepackage{etoolbox}
\usepackage{upgreek}
\usepackage{graphicx}
\usepackage{amsfonts}
\usepackage{amssymb} 
\usepackage{bm} 
\usepackage{slashed}
\csdef{Seqnsplit}{\\}
\def\botline{\\\hline}
\newcounter{runlevel}
\usepackage{multirow} 
\def\nrow#1{\@tempcnta #1\relax%
\advance\@tempcnta by 1\relax%
\xdef\lenrow{\the\@tempcnta}}
\def\morerows#1#2{\nrow{#1}\multirow{\lenrow}{*}{#2}}
\let\MakeYrStrItalic\relax
\def\jobid{crmeca20240892}
%\graphicspath{{/tmp/\jobid_figs/web/}} 
\graphicspath{{./figures/}} 
\def\0{\phantom{0}}
\def\rmmu{\upmu}
\def\mn{\phantom{$-$}}
\def\refinput#1{}
\def\back#1{}
\newcommand{\bnabla}{\boldsymbol{\nabla}}
\newcommand{\bgamma}{\boldsymbol{\gamma}}
\newcommand{\bdelta}{\boldsymbol{\delta}}
\newenvironment{inftab}{}{}
\def\Custref#1#2{\hyperlink{#2}{#1}}
\def\og{\guillemotleft}
\def\fg{\guillemotright}
\makeatletter
\g@addto@macro{\UrlBreaks}{\UrlOrds}
\gappto{\UrlBreaks}{\UrlOrds}
\makeatother
\DOI{10.5802/crmeca.287}
\datereceived{2024-10-30}
\daterevised{2025-01-31}
\dateaccepted{2025-02-06}
\ItHasTeXPublished
\newtheorem{lemmaa}{Lemma}[section]
\newtheorem{theor}[lemmaa]{Theorem}
\newtheorem{rem}[lemmaa]{Remark}
}
{
\documentclass[crmeca,eqnobysec]{article}
\usepackage[T1]{fontenc} 
\usepackage[T2A,T1]{fontenc}
\usepackage{upgreek}
\usepackage{slashed}
\def\bibliographystyle#1{}
\def\bibliography#1{}
\let\refinput\input
\def\CDRdoi{10.5802/crmeca.287}
\def\newline{\unskip\break}
\def\sfrac#1#2{{#1}/{#2}}
\def\xxtabnote#1{\tabnote{#1}}
\def\xmathop#1{\mathop{{#1}}}
\def\xtbody{\tbody}
\def\hypertarget#1#2{\label{#1}}
\newtheorem{defi}{Definition}
\newtheorem{prop}[defi]{Proposition}
\newtheorem{lemmaa}{Lemma}[section]
\newtheorem{theor}[lemmaa]{Theorem}
\newtheorem{rem}[lemmaa]{Remark}
}

\usepackage{moreverb}

\begin{document}

\begin{noXML}

\CDRsetmeta{articletype}{review}

\title{The stability of Minkowski space and its influence on the
mathematical analysis of General Relativity}

\alttitle{La stabilit\'e de l'espace de Minkowski et son influence sur
l'analyse math\'ematique de la Relativit\'e G\'en\'erale}

\author{\firstname{Jacques} \lastname{Smulevici}\CDRorcid{0000-0003-0671-4656}}
\address{Laboratoire Jacques-Louis Lions (LJLL), Sorbonne Universit\'e,
CNRS and Institut Universitaire de France, France}
\email{jacques.smulevici@sorbonne-universite.fr}

\begin{abstract}
The study of the stability of Minkowksi space has had a tremendous
influence on the mathematical analysis of General Relativity and in a
wider context, on geometric hyperbolic partial differential equations. 
This article proposes an overview of the context, ideas and part of the
legacy surrounding this subject, focusing particularly on the original
proof of Christodoulou--Klainerman. 
\end{abstract}

\begin{altabstract}
L'\'etude de la stabilit\'e de l'espace de Minkowski a eu une influence
majeure sur l'analyse math\'ematique de la Relativit\'e G\'en\'erale
et, plus largement, sur les \'equations aux d\'eriv\'ees partielles
hyperboliques g\'eom\'etriques. Cet article propose un aperÃ§u du
contexte, des id\'ees et d'une partie de l'h\'eritage entourant ce
sujet, en mettant particuli\`erement l'accent sur la d\'emonstration
originale de Christodoulou--Klainerman.
{\vspace*{-2pt}}
\end{altabstract}

\keywords{\kwd{Einstein equations}\kwd{General Relativity}\kwd{Asymptotic stability}}

\altkeywords{\kwd{\'Equations d'Einstein}\kwd{Relativit\'e G\'en\'erale}\kwd{Stabilit\'e asymptotique}}

\maketitle

\end{noXML}

\numberwithin{equation}{section}

{\vspace*{-2pt}}

\section{Introduction}

{\vspace*{-2pt}}

Since the work of Minkowski and its unification of space and
time~\cite{Minkowski1909}, Special Relativity is beautifully encoded in
the geometry of the simplest Lorentzian manifold, the Minkowski space
$(\mathbb{R}^{1+3}, \eta)$, where $\eta$ is the Lorentzian metric whose
matrix of components in the standard Cartesian coordinates $(t,x,y,z)$
of $\mathbb{R}^{1+3}$ is given by $\left(\eta_{\alpha
\beta}\right)_{0\le \alpha, \beta \le 3}= \mathrm{diag}(-1, 1, 1, 1)$. 

The notion of simplest Lorentzian manifold can even be formalized,
since by the Killing--Hopf Theorem \cite[Chapter~8]{ONeill1983},
Minkowski space is the unique simply connected complete Lorentzian
space of vanishing constant sectional curvature. For the vacuum
Einstein equations, $\mathrm{Ric}(g)=0$,  it should therefore be
regarded as the trivial solution, in a similar way that the scalar
function $\phi :=0$ is the trivial solution to a non-linear homogeneous
wave equation of the type $\Box \phi :=- \partial_t^2 \phi + \Delta
\phi= F(\phi, \partial \phi, \partial^2 \phi)$, for some non-linearity
$F$ vanishing at $(0,0,0)$. 

The study of the stability properties of the Minkowski space is
natural, both from a physical or a mathematical point of view, once the
Einstein equations are understood as a system of evolution equations;
this is the problem of understanding the dynamics in a neighborhood of
an equilibrium, in fact, of the simplest equilibrium of General
Relativity. The subject was also clearly of interest to Choquet--Bruhat
and, as for many other mathematical problems of General Relativity, she
produced several of the earlier works attempting to address it. In
fact, the well-posedness statement of Choquet--Bruhat on the Einstein
equations\footnote{In particular, the estimates of \cite[Chapter
III]{choquet1952theoreme} leading to a reformulation of the existence
theory into a fix point argument implicitly address an analysis of the
difference of solutions. See also~\cite{MR0291607} for a subsequent
treatment in Sobolev spaces which focuses on Cauchy stability and
states the application to the construction of near Minkowskian
spacetimes for finite time.} already provides Cauchy stability, hence
the local in time stability of the Minkowski space, as long as initial
data close to Minkowski space solving the constraint equations could be
constructed. This topic she addressed in~\cite{ctx13996358890006616}
and~\cite{MR0341358}, where she also mentioned as an open problem the
question of long time evolution of such data: 

\smallskip

\begin{quote}
\fontsize{8.65}{10.2}\selectfont
{The important question whether any asymptotically flat
excitations of the gravitational field remain nonsingular for all time
is still open.}
\end{quote}

\smallskip

This question will of course only be fully resolved through the
monumental work of Christodoulou--Klainerman, completed in 1989 and
published in 1993 in~\cite{CKM} establishing the global-in-time
stability of the Minkowski space. 

\begin{theor}[\cite{CKM}]
Under a global smallness assumption, any strongly asymptotically flat,
maximal, initial data leads to a unique globally hyperbolic regular
solution to the Einstein equations which is geodesically complete and
globally asymptotically flat. 
\end{theor}

The initial value problem for the Einstein equations, including the
relevant notion of initial data are reviewed below in 
Section~\ref{se:ceivp}, while we recall the maximal foliation in 
Section~\ref{se:max}. On top of constructing a global solution,
detailed asymptotics for the solution were obtained in~\cite{CKM}, such
as decay estimates for the curvature or for the second fundamental form
of the maximal foliation used. 

Prior to~\cite{CKM}, a work of Friedrich~\cite{MR0868737} had already
established a form of future stability for hyperboloidal data close to
Minkowski, but at the time of writing, it was not known how to reach
the required hyperboloidal data from asymptotically flat initial data. 
The gluing constructions of Corvino~\cite{MR1794269},
Chrusciel--Delay~\cite{MR1902228} and Corvino--Schoen~\cite{MR2225517}
were eventually able to produce such data, bypassing the issue of
understanding the dynamics of asymptotically flat data in a
neighborhood of spatial infinity and the construction of a small piece
of null infinity emanating from it, the \textit{exterior problem}.  
Nonetheless, since the interior dynamics are heavily influenced by the
asymptotics of the data, it is clear that one should aim to address the
stability problem for the larger possible class of data, and in fact,
several interesting physical phenomena can only be captured by data
which are not strongly localized\footnote{See for
instance~\cite{Kehrberger_2024} for historical details and a
presentation of the physical relevance of slowly decaying asymptotic
flat initial data.}.

The work of Christodoulou--Klainerman~\cite{CKM} not only did provide
the first complete proof of the stability of Minkowski space, but
thanks to the precise asymptotics of the solutions obtained, also did
contain the first rigorous derivation of the laws of gravitational
radiation, such as the \textit{Bondi mass loss
formula}~\cite{Bondi:1960jsa}. The non-linear \textit{memory effect} of
Christodoulou~\cite{PhysRevLett.67.1486} was also obtained essentially
as an outgrowth of the estimates obtained in the proof. 

The stability of Minkowski space also constitutes a landmark in the
field as the first proof in this domain of research which was too long
to be published as a standard journal article; making harder the task
of summarizing its content and highlighting its key ideas\footnote{See
nonetheless~\cite{MR1157847} for an introduction and background
materials.}; and nowadays one could say that it has become a sort of
rite of passage to study or revisit (an excerpt of) the proof or write
some form of extension.

Starting in the late 70's, the study of the long-time behaviour of
non-linear wave equations with small data had seen major
development\footnote{See for instance~\cite{Zuily1993-1994} for a
pertinent review.} and in particular, the standard overall strategy
based on energy and decay estimates, wrapped up in a bootstrap
argument, had already seen many applications for global or long-time
existence. Nonetheless, using this in the context of the Einstein
equations was still challenging a priori. Choquet--Bruhat had in
particular identified in~\cite{MR0315300} a form of instability in wave
coordinates, essentially proving that the components of the metric in
wave coordinates near Minkowski cannot all be uniformly bounded in the
standard linear spaces, such as the energy space. Thus, new important
building blocks were needed, in particular
\begin{itemize}
\item A robust method to derive decay estimates for tensorial wave
equations, including improved decay for certain components or
derivatives of the tensor fields. 
\item A way to bypass the instability highlighted by Choquet--Bruhat.
In~\cite{CKM}, a key component to avoid the instability is a novel
geometric framework which is adapted to the solution: the vector fields
used for commutation and the foliation are not fixed a priori but
constructed throughout the proof. 
\end{itemize}
The first item above was initially developed
in~\cite{cpa3160430202}, following the
development of commutation vector field methods for wave
equations~\cite{37c1f190d8824edab970cf9c04037bb7} as well as the
understanding of the null
condition~\cite{cpa3160390205, MR0837683}. The
failure of the standard null condition is precisely what drove the
instability mechanism identified by Choquet--Bruhat and thus, one of the
aims of the second item was to provide a better framework to capture
the null structural properties of the Einstein equations. Eventually,
Lindblad--Rodnianski~\cite{MR1994592} identified that the Einstein
equations in wave coordinates actually verify a weaker form of the null
condition, that allows global existence and decay estimate, while the
main effects of the instability of Choquet--Bruhat remain only as a
small (logarithmic) growth for the worse component of the
metric~\cite{MR2134337, MR2680391}. 

The aim of the present article is to review some ideas and context
behind the stability of the Minkowski space as well as to present some
recent developments concerning the study of the long time dynamics for
the Einstein equations or, more generally, geometric hyperbolic partial
differential equations,  which uses or extend some important ideas of
the proof. 

\subsection{Overwiew}

In Section~\ref{se:bm}, we briefly review standard material concerning
the Einstein equations. Besides fixing some notations, introducing the
evolution problem in General Relativity and emphasizing the importance
of gauge conditions for the Einstein equations, we also aim here to
explain why the approach used for instance for orbital type stability
results in the case of critical or subcritical non-linear SchrÃ¶dinger
or wave equations does not apply here. In Section~\ref{se:geqw}, we
review the problem of global existence for quasi-linear wave equations,
provide a quick overview of the original Klainerman vector field method
for scalar fields and the null condition in $3$ spatial dimensions and
end the section with a very brief discussion of the weak null condition
and the approach of Lindblad--Rodnianski. Section~\ref{se:ckp} is
devoted to the presentation of the work of Christodolou--Klainerman,
focusing on geometric aspects. Though nowdays, there are several
aspects of the proof that one would treat differently, we have chosen
to focus on the original arguments. In Section~\ref{se:desm}, we review
several extensions and developments on the theme of the stability of
Minkowski space and we conclude in Section~\ref{se:con} with some works
which were influenced by the analysis on the stability of Minkowski
space but address other kinds of problems.

\subsection{Convention}

We use standard Einstein summation convention. Greek indices $\alpha,
\beta, \ldots $ are used for tensors defined on spacetimes. Latin
indices $i, j,k$ are used for tensors on spacelike hypersurfaces or the
manifold of the initial data. The convention for the signature of
Lorentzian metric is $- ++\ldots\,$.  If $g_{\alpha \beta}$ denotes the
components of a metric $g$, $g^{\alpha \beta}$ denotes the components
of the inverse metric and we also used standard raising and lowering of
indices on tensors throughout the article. Partial derivatives will be
denoted either as $\partial_{\alpha} \psi$ or $\psi_{, \alpha}$.
Starting from Section~\ref{se:ckp}, we follow the convention
of~\cite{CKM} that bold quantities such as $\mathbf{R}$, $\mathbf{D}$
are associated to the spacetime geometry and thus denotes the spacetime
Riemann tensor and the Levi-Civita connexion, while $R$, $\nabla$
refers to the geometry of a spacelike slice. 

The article is written assuming standard knowledge of differential and
Riemannian geometry as well as partial differential equations, but does
not assume a priori knowledge concerning General Relativity or the
Einstein equations. 

\section{Background material}\label{se:bm}
\subsection{The Einstein equations and their hyperbolicity}

General Relativity can be understood as the study of
$4$-dimensional\footnote{The theory can easily be written in other
dimensions, but on top of being the physical dimension, dimension $4$
plays a distinguished role for the problem of the stability of
Minkowski space, cf.\ Section~\ref{se:nc} below.} Lorentzian manifolds
$(\mathcal{M}, g)$ solving the Einstein equations 
{\begin{equation} \label{eq:EE}
\mathrm{Ric}(g) - \tfrac{1}{2}\mathrm{S}(g)\cdot g + \Lambda g =
{T}, 
\end{equation}}\unskip
where $\mathrm{Ric}(g)$ and $\mathrm{S}(g)$ denotes respectively the
Ricci and scalar curvature of the Lorentzian metric $g$, $\Lambda \in
\mathbb{R}$ is the so-called cosmological constant and ${T}$ is
a source term, arising from matter fields, the so-called
energy--momentum tensor. This source term, ${T}$, should not be
considered as pure external data since the evolution of matter fields,
whether they are constituted of electromagnetism fields, fluids,
particles,\,\ldots\ are themselves subject to an evolution equation
depending the metric $g$. Thus, unless ${T} = 0$, the Einstein
equations need to be completed by the evolution equations for the
matter to get a closed system. When ${T}=0$, Equation
(\ref{eq:EE}) reduces to $\mathrm{Ric}(g) = \Lambda g$, the
\textit{vacuum Einstein equations}. Minkowski space is the trivial
solution to (\ref{eq:vEE}) with $\Lambda=0$, while the de-Sitter and
Anti-de-Sitter space solve the equation with $\Lambda &gt; 0$ and $\Lambda
&lt; 0$ respectively. In the rest of the article\footnote{Apart from
Section~\ref{se:ems}, when we will review results concerning $T \neq
0$.}, we thus focus on $\Lambda=0$ together the vacuum condition
${T}=0$, in which case the equations reduce to the Ricci flat
condition
{\begin{equation} \label{eq:vEE}
\mathrm{Ric}(g) = 0. 
\end{equation}}\unskip

Perhaps the most important property of the Einstein equations is that
they are geometric: given $\phi : \mathcal{M} \rightarrow \mathcal{M}$ a
diffeomorphism, $g$ solves  (\ref{eq:vEE}) if and only if $\phi^{\star}
g$ solves (\ref{eq:vEE}), where $\phi^{\star} g$ denotes the pullback
of the metric $g$ by $\phi$. 

This already implies that one cannot expect to solve the Einstein
equations based purely on a traditional, direct, application of the
methods from partial differential equations. Indeed, those typically
gives existence and uniqueness at the same time and here uniqueness can
only hold in the geometric sense, that is to say up to the above
diffeomorphism invariance. 

More directly, the Ricci tensor of a metric $g$, viewed as a
(non-linear) partial differential operator on the components of $g$,
takes the form
{\begin{eqnarray} \label{eq:Ricdiffg}
\mathrm{Ric}(g)_{\alpha \beta} = - \tfrac{1}{2} \Box_g g_{\alpha
\beta} + \partial_\alpha \Gamma_\beta + \partial_\beta \Gamma_\alpha +
Q_{\alpha \beta} (\partial g, \partial g ), 
\end{eqnarray}}\unskip
where 
\begin{itemize}
\item $\Box_g$ is the (scalar) wave operator: 
{$$
\Box_g \psi := \dfrac{1}{\sqrt{|\det g |}} \left(\partial_\alpha 
\sqrt{|\det g |} g^{\alpha \beta} \partial_\beta \psi \right).
$$}\unskip
\item $\Gamma_\beta$ are traces of the Christoffel symbols
$\Gamma^{\rho}_{\alpha \gamma}$: 
{$$
\Gamma _\beta := g_{\beta \rho} g^{\alpha \gamma} \Gamma^{\rho}_{\alpha
\gamma} := g_{\beta \rho} g^{\alpha \gamma} \tfrac{1}{2} g^{\rho \delta}
(g_{\delta \gamma, \alpha} + g_{\alpha \delta, \gamma} - g_{\alpha
\gamma, \delta} ).
$$}\unskip
\item $Q_{\alpha \beta} (\partial g, \partial g )$ denotes a quadratic
form in the first derivatives of $g$, with coefficients depending on
$g$.
\end{itemize}

For a fixed metric $g$,  $\Box_g$ is a second-order hyperbolic
operator, but the terms $\partial_\alpha \Gamma_\beta + \partial_\beta
\Gamma_\alpha$ contains second-order derivatives of $g$ and thus a
priori spoil the hyperbolicity of $\mathrm{Ric}(g)$. To retrieve the
hyperbolicity, one needs to address the geometric invariance mentioned
above, for instance by considering the equations in a particular
coordinate system. The standard one, which of course was used by
Choquet--Bruhat in her local well-posedness
proof~\cite{choquet1952theoreme}, is given by the so-callled
\textit{wave coordinates}. 

\begin{defi}
\upshape
A local coordinate system $\left(\mathcal{U} \subset \mathcal{M},
\phi :=(x^{\alpha})_{\alpha=0,1,2,3}: \mathcal{U} \rightarrow
\mathbb{R}^4\right) $ of a Lorentzian manifold $(\mathcal{M}, g)$ is
said to consist of \textit{wave coordinates} if, $\forall
\alpha=0,\ldots, 3$, the function $x{^\alpha}$ verifies the wave equation
{\begin{equation} \label{eq:wcc}
\Box_g x^\alpha=0.
\end{equation}}\unskip
\end{defi}

Since in general, one has $\Box_g x^\alpha= - g^{\gamma \beta}
\Gamma^{\alpha}_{\gamma \beta}=-g^{\beta \alpha} \Gamma_\beta$, we have
as a consequence: 

\begin{lemmaa}
In wave coordinates, the vacuum Einstein equations (\ref{eq:vEE})
reduces to a system of the form 
{\begin{eqnarray} \label{ee:red}
\Box_g g_{\alpha \beta} = Q_{\alpha \beta} (\partial g, \partial g ), 
\end{eqnarray}}\unskip
\end{lemmaa}

The Equations (\ref{ee:red}) are a system of \textit{quasi-linear} wave
equations. For such equations, standard energy methods can be applied
to prove local existence\footnote{See for instance \cite[Chapter
1.4]{MR1715192}.} in high regularity classes. 

Another important property, which again can be seen to arise from the
geometric invariance, is that (\ref{ee:red}) are scale invariant: if
$g_{\alpha \beta}(x^\alpha)$ solves (\ref{ee:red}), restricting
appropriately the domain of definition, so does $^{\lambda}g(x) :=
g_{\alpha \beta}(\lambda\, x^\alpha)$, for any $\lambda &gt; 0$. In three
spatial dimensions, the Sobolev space left invariant by this scaling is
the homogeneous space $\overset{\scriptscriptstyle\circ}{H}{}^{3/2}$.
In particular, this implies that one cannot hope to study the Einstein
equations, even locally and thus a fortiori globally in time, without
controlling more than one derivative of the metric components, for
instance by taking derivatives of the equations. 

\subsection{The constraint equations and the initial value problem}
\label{se:ceivp}

Given any spacelike hypersurface $\Sigma$, the Gauss--Codazzi equations
implies, for any solution of (\ref{eq:vEE}), a system of constraint
equations for the first and second fundamental forms $h,\, k$
associated to $\Sigma$, 
{\begin{eqnarray}
\mathrm{S}(h) + (\mathrm{tr}\,k)^2- \vert k \vert^2 &amp;=&amp; 0, 
\label{eq:const1} \Seqnsplit
\mathrm{div}\, k - \mathrm{d}(\mathrm{tr}\,k) &amp;=&amp; 0, \label{eq:const2}
\end{eqnarray}}\unskip
where $\mathrm{S}(h)$ denotes the scalar curvature of the Riemannian
metric $h$, and operators and traces are taken with respect to $h$.

The initial value problem in General Relativity then consists, given a
Riemannian manifold $(\Sigma,h)$ and a $2$-tensor $k$ verifying
(\ref{eq:const1})--(\ref{eq:const2}), in the construction of a
Lorentzian manifold $(\mathcal{M},g)$ solving (\ref{eq:vEE}) together
with an embedding $i: \Sigma \hookrightarrow \mathcal{M}$ such that
$(h,k)$ coincides with the first and second fundamental forms of $i$.
Such an $(\mathcal{M},g)$ is then called a \textit{development} of the
data $(\Sigma, h, k)$.  

The Choquet--Bruhat theorem~\cite{choquet1952theoreme}, completed by
the Choquet---Bruhat--Geroch construction of the maximal
solution\footnote{See also the modern treatment by
Sbierski~\cite{MR3447847} for an improved construction of the maximal
solution, which contrary to~\cite{cmpx1103841822} does not rely on
Zorn's lemma.}~\cite{cmpx1103841822},  then provides, given any initial
data set $(\Sigma, h, k)$, the existence and uniqueness, up to the
diffeomorphism invariance, of a \textit{maximal globally
hyperbolic}\footnote{Here by \textit{globally hyperbolic}, we mean the
existence of a \textit{Cauchy hypersurface}, that is an achronal
hypersurface such that any inextendible causal curve register once and
only once on this hypersurface, cf.\ \cite[Chapter 14]{ONeill1983}.}
development. 

\subsection{Asymptotic flatness} \label{se:af}

The Euclidean data $(\mathbb{R}^3, \delta, 0)$ form of course a
solution to (\ref{eq:const1})--(\ref{eq:const2}) and its maximal Cauchy
development is given by the Minkowski spacetime. For any concept of
perturbations of the Minkowski spacetime to make sense, one thus also
needs to understand how to construct initial data sets which are
arbitrary close to $(\mathbb{R}^3, \delta, 0)$. In particular, we would
like the pertubations from the Euclidean data to be \textit{localized},
i.e.\ to tend to $0$ near infinity. Physically, this corresponds to
\textit{asymptotically flat} initial data sets and gives rise to
spacetimes describing \textit{isolated systems}. From the point of view
of analysis, this localization is important, as it is needed for the
derivation of decay estimates, which are at the core of the stability
proofs. To introduce this quantitatively, we consider initial data sets
$(\Sigma,h,k)$ such that $\Sigma$ possesses an asymptotic end, that is,
the complement of a compact set of $\Sigma$ is diffeomorphic to the
complement of a ball in $\mathbb{R}^3$. This allows one to introduce a
coordinate system $(x^1, x^2, x^3)$ in a neighborhood of infinity $r:=
\sqrt{ \sum_{i=1}^3 (x^i)^2}\rightarrow+\infty$. The original
proof~\cite{CKM} then applies to \textit{strongly asymptotically flat
initial data}, which by definition verifies the following
expansion\footnote{Recall that $f=o_m(r^{-k})$ if $\partial_l f =
o(r^{-k-|l|})$, for any multi-index $l$ of length $| l | \le
m$.}
{\begin{eqnarray}
h_{ij} &amp;=&amp; (1 + 2M/r) \delta_{ij} + o_4 (r^{-3/2}), 
\label{eq:metrdata} \Seqnsplit
k_{ij} &amp;=&amp; o_3( r^{-5/2}),
\end{eqnarray}}\unskip
for a non-negative number $M\ge 0$. The positive
mass-theorem~\cite{cmpx1103904790, MR0626707} in fact guarantees that
for such data, the mass $M$ is strictly positive unless the initial
data is flat. In particular, one cannot consider compactly supported
perturbations of the flat data. The construction of asymptotically flat
initial data near the Euclidean data was first considered
in~\cite{ctx13996358890006616} and~\cite{MR0341358}. The problem of
producing initial data with arbitrary decay has been more recently
addressed in~\cite{stj:idms} and~\cite{MR4625954}. We refer to
\cite[Chapter 3]{MR2391586} for a more precise introduction to
asymptotic flatness in General Relativity. 

\section{Global existence for quasi-linear wave equations and the basic
vector field method} \label{se:geqw}
\subsection{A strategy for quasi-linear wave equations} 
\label{se:qlns}

In view of the form of the reduced Equations (\ref{ee:red}), it is
natural to first think of the stability of the Minkowski spacetime as a
specific case of global existence for quasi-linear wave equations of
the form
{\begin{eqnarray} 
(- \partial_t^2 + \Delta) u_i + F^{\alpha \beta}(u)
\partial _{\alpha} \partial _\beta(u_i) &amp;=&amp; Q_{i}(\partial u, \partial
u), \label{eq:QLW} \Seqnsplit
u (t=0) &amp;=&amp; u_0, \Seqnsplit
\partial _t u  (t=0) &amp;=&amp; u_1, 
\end{eqnarray}}\unskip
for an unknown $u:(t,x) \rightarrow (u_i)(t,x) \in \mathbb{R}^N$,
non-linearities such that $F^{\alpha\beta}(u)=F^{\beta \alpha}(u)$ are
smooth functions of $u$ vanishing at $0$, $Q_{i}(\partial u, \partial
u)$ which are quadratic forms in the first derivatives of $u$ and
initial data $(u_0,u_1)$ which are small, sufficiently regular and
localized. 

By now, the standard approach for this type of problems is essentially
based on
\begin{itemize}
\item A scale of weighted energy spaces to control the data and the
solutions. We denote the top order energy norm at time $t$ of a
solution $u$ by $E_N[u(t)]$, the index $N$ referring to the order in
the scale, i.e.~$E_N[u]$ controls up to $N$ derivatives of $u$ in a
weighted $L^2$-type space. 
\item Energy estimates allowing to propagate global bounds of the
solutions, in particular at top order of regularity, of the form
{\begin{equation} \label{es:energy}
E_N[u(t)] \le C E_N[u(0)] + \int_0^t \mathrm{NL}(s)\,\mathrm{d}s, 
\end{equation}}\unskip
where $\mathrm{NL}$ denotes error terms arising from the
non-linearities and $C&gt;0$ is independent of $u$. 
\item Decay estimates that allow to control the non-linearities,
leading, for instance, to an estimate of the form
{\begin{equation}\label{es:decay}
|\mathrm{NL}(s)| \lesssim \frac{E_N[u(s)] }{(1+s)^q}, 
\end{equation}}\unskip
with $q &gt; 0 $. 
\item A bootstrap argument: by local existence, for initial data such
that $E_N[u(0)] \le \epsilon$, a solution $u$ exists on some interval
$[0, T]$, for some $T &gt;0$, and such that $E_N[u(t)] \le 2 C \epsilon $,
for all $t \in [0,T]$. Assuming $q&gt;1$ for the time being in
(\ref{es:decay}), a GrÃ¶nwall argument combined with the energy
estimate (\ref{es:energy}) and the decay estimate (\ref{es:decay}) then
allows to improve the bootstrap assumption to $E_N[u(t)] \le  C
\epsilon +O(\epsilon^2) \le 3/2 C \epsilon.$ By a standard continuity
argument, it then follows that the solution exists for all times and
verifies, for all $t \ge 0$,  
{$$
E_N[u(t)] \le C \epsilon +O(\epsilon^2).
$$}\unskip
\end{itemize}
To be more specific and check under which conditions, one can obtain
the estimate (\ref{es:decay}) with $q&gt;1$, we recall below the standard
vector field method, as introduced
in~\cite{37c1f190d8824edab970cf9c04037bb7}, that provides both a first
definition of the spaces $E_N$ and the decay estimate~(\ref{es:decay}).

\subsection{The vector field method of Klainerman} \label{se:sk}

We consider the wave operator $\Box := -\partial_t^2 + \Delta$, where
$\Delta$ is the usual Laplacian in $\mathbb{R}^n$ and recall first the
basic energy estimate. 	

For $\phi$ and $S$ defined on $\mathbb{R}_t \times \mathbb{R}^n_x$,
decaying sufficiently fast in $x$ and solving $\Box \phi = S,$ the
energy 
{$$
E[\phi](t):= \int_{x \in \mathbb{R}^n} (\partial_t \phi(t,x) )^2+
|\nabla_x \phi(t,x) |^2 \,\mathrm{d}x
$$}\unskip
verifies for all $t$, 
{$$
E[\phi](t) = E[\phi](0) + \int_0^t \int_{\mathbb{R}^n}
|\partial\phi\cdot S| \,\mathrm{d}x \,\mathrm{d}s, 
$$}\unskip
which one can (at least formally) derive by multiplying the wave
equation by $\partial_t \phi$ and integration by parts. Moreover,  if
$\phi$ solves instead the perturbed wave equation 
{$$
(-\partial_t^2 + \Delta ) \phi + F^{\alpha \beta}\partial_\alpha
\partial_\beta \phi=S,
$$}\unskip
with $F^{\alpha \beta}=F^{\beta\alpha}$ sufficiently small
(${\Vert F \Vert} _{L^\infty_{t,x}} &lt; 1/4$), then the above is replaced by the
energy estimate: 
{\begin{equation}\label{es:ener}
E[\phi](t) \le 2 E[\phi](0) +2 \int _0^t \int _{\mathbb{R}^n} 
|\partial \phi| \cdot (|S| + |\partial \phi| |\partial F| )
\,\mathrm{d}x \,\mathrm{d}s. 
\end{equation}}\unskip

To construct higher order norms that are still propagated by the
equation, one then considers vector fields that commute with the linear
wave operator: 

\begin{lemmaa}\label{lem:kvf}
Let $Z$ be any of the vector fields 
\begin{itemize}
\item Translations $\partial_t, \partial_{x^i}$, 
\item Rotations $\Omega_{ij}= x^i \partial_{x^j} - x^j \partial_{x^i}$,
\item Hyperbolic rotations $ \Omega_{0i} = t\partial_{x^i} + x^i
\partial_t $. 
\item Scaling $ S=  t \partial_t + \sum_i x^i \partial_{x^i}$. 
\end{itemize}
Then, $[ \Box, Z ] = 0$, apart from the scaling which verifies $[
\Box , S ] = 2 \Box$. 
\end{lemmaa}

Denote by $\mathcal{A}$ the set of all such vector fields and recall
that the $Z$ vector fields form an algebra: $[ Z_1, Z_2]  \in
\textrm{span}\, \mathcal{A}$, for all $Z_1, Z_2 \in \mathcal{A}$. We
introduce an ordering on $\mathcal{A}$ and use standard multi-index
notations below. In particular, for a multi-index $\alpha$, $|\alpha|$
denotes its length and $Z^\alpha$ denotes a differential operator of
order $|\alpha|$ obtained as a composition of $ |\alpha|$ vector fields
within $\mathcal{A}$.

The decay estimate is then essentially incoporated in the following
weighted Sobolev inequality.

\begin{prop}\label{prop:KS}
Let $\psi :=\psi(t,x)$ defined on $[0, T] \times \mathbb{R}_x^n$ be a
sufficiently regular function. Then, $\forall (t,x) \in [0, T] \times
\mathbb{R}_x^n$,
{\begin{equation}
|\psi(t,x)| \lesssim \dfrac{1}{(1+t+|x|)^{(n-1)/2}(1+|t-|x||)^{1/2}}
\sum_{|\alpha| \le (n+2)/2,\, Z^\alpha \in \mathcal{A}^{|\alpha|}}
\Vert Z^\alpha \psi(t,\,.\, ) \Vert _{L^2_x}. \label{ineq:ks}
\end{equation}}\unskip
\end{prop}

Let $\psi$ be a solution to the linear wave equation $\Box \psi=0$. 
Applying (\ref{ineq:ks}) to $\partial \psi$, we obtain
{\begin{eqnarray*}
|\partial \psi(t,x)| 
&amp;\lesssim &amp; \dfrac{1}{(1+t+|x|)^{(n-1)/2}(1+|t-|x||)^{1/2}} \sum_{|\alpha| \le (n+2)/2,\,\, Z \in \mathcal{A}^{|\alpha|}} \Vert Z^\alpha \partial \psi (t,\,.\, ) \Vert _{L^2_x} \\
&amp;\lesssim &amp; \dfrac{1}{(1+t+|x|)^{(n-1)/2}(1+|t-|x||)^{1/2}} \sum_{|\alpha| \le (n+2)/2,\,,\ Z \in \mathcal{A}^{|\alpha|}} \Vert \partial Z^\alpha \psi (t,\,.\, ) \Vert _{L^2_x} \\
&amp;\lesssim &amp; \dfrac{1}{(1+t+|x|)^{(n-1)/2}(1+|t-|x||)^{1/2}} \sum_{|\alpha| \le (n+2)/2,\,,\ Z \in \mathcal{A}^{|\alpha|}} E[ Z^\alpha \psi]^{1/2}(0),
\end{eqnarray*}}\unskip
where we have used that $\Box Z^\alpha \psi=0$ and hence $E[
Z^\alpha \psi ](t)= E[ Z^\alpha \psi ](0)$ from energy conservation. 

Thus, we have obtain a proof of \textit{decay} for solutions to the
homegeneous wave equation. Moreover, we see that a natural scale of
energy norms is given by
{$$
E_N[u(t)] := \sum_{|\alpha| \le N,\,,\,Z \in \mathcal{A}^{|\alpha|} } E
[Z^\alpha u(t)].
$$}\unskip

\begin{rem}\label{rem:mor}
\upshape
We see from (\ref{ineq:ks}) that away from the the light cone region
$t=|x|$, for instance, in the interior region, $|x| &lt; ({1}/{2})t$, the
strong decay $|\partial \phi| \lesssim t^{-3/2}$ holds. Provided the
initial data enjoy enough decay, this estimate can be
improved\footnote{In fact, for compactly supported initial data,
Huygens' principle ensures that the support of the solution is located
in a cylindrical strip around the light-cone.}. One way to capture
extra decay is to consider a stronger weighted energy norm that also
make use of the multiplier $K_0:= 2 t S+ (-t^2+|x|^2)\partial_t$. This
vector field, first used by Morawetz in~\cite{cpa3160150303}, 
also played an
important role in the estimates of~\cite{CKM}, cf.\ Section~\ref{se:brt}
below. 
\end{rem}

\begin{rem} \label{rem:mvsc}
\upshape
In the strategy described above, vector fields play two different kind
of roles. First, the overall energy estimates can be obtained by
\textit{multiplying} the wave equation by $X(\phi)$, for $X=X^\alpha
\partial_{x^\alpha}$ a vector field, and integration over suitable
regions. For instance, the basic energy estimate (\ref{es:ener}) relies
on the vector field multiplier $X=\partial_t$, while we introduce just
above the multiplier $X=K_0$, leading to weighted energy estimates. The
other occurrence of the vector fields is via \textit{commutation}, cf.\ 
Lemma~\ref{lem:kvf} and Proposition~\ref{prop:KS}. The existence of
these commutation vector fields is highly connected with the symmetry
of the equations, hence is less likely to be useful outside from
questions related to the stability of Minkowski. On the other hand, the
use of multiplier vector fields has an extremely wide range of
applications. In fact, it is part of the \textit{general method of
energy currents}, cf.\ \cite{cpa3160150303,
MR1739321}, and many recent developments in General Relativity and
geometric hyperbolic partial equations rely on some applications or
extensions of these principles, cf.\ Section~\ref{newvfm} below.
\end{rem}

\subsection{Global existence with $n \ge 4$}

Returning to the scheme introduced in Section~\ref{se:qlns}, commuting
the equations $N$-times and applying the energy estimate
(\ref{es:ener}), we obtain the top order energy estimate
(\ref{es:energy}) with an error term $NL$ which can be estimated as
{$$
| NL |(s) \lesssim E_N[u(s)] \cdot \Vert \partial^\alpha u(s) 
\Vert _{L^\infty_x}, 
$$}\unskip
with $|\alpha | \le N/2$. If $N/ 2 \ge (n+2)/2$, we can then appeal to
the inequality (\ref{ineq:ks}), to obtain (\ref{es:decay}) with
$q=(n-1)/2$. In particular, with $n \ge 4$, we have $q &gt;1$ and we can
improve the bootstrap assumptions as explained in
Section~\ref{se:qlns}, while we see that we only have $q=1$ if $n=3$,
so that this strategy a priori fails.

\subsection{The case $n=3$: the null condition} \label{se:nc}

With $n=3$, the linear decay and a GrÃ¶nwall argument in
(\ref{es:ener}) a priori gives a logarithmic growth. That this growth
can in fact lead to blow up for certain systems of quasi-linear wave
equations was proved in the celebrated works of John~\cite{MR0600571,
MR0808321}. Nonetheless, global existence can be recovered provided the
non-linearities verify certain structural properties, in particular the
\textit{null condition}~\cite{MR0837683, MR820070}, which guarantees
that the non-linear terms always contain at least a derivative
\textit{tangential} to the light-cone. 

For instance, for the semi-linear terms $Q_i$ in (\ref{eq:QLW}), such a
condition would imply, for $t\ge 0$, that they can be written as a
linear combination with bounded coefficients of products of the form
{$$
\underline{\partial} u \cdot \partial u,
$$}\unskip
where $\partial $ denote all possible derivatives, while
$\underline{\partial} \in \{\partial_t+ \partial_r,\, e_A,\, A =1, 2 \}
$, with $e_A$ derivatives tangential to the spheres of constant
$(t,r)$. The relevance of these tangential derivatives is that they
verify an improved decay estimate. In fact, for any tangential
derivative $\underline{\partial}$, one has a decomposition of the form
{$$
( 1+ |t| + |x| )\underline{\partial} = 
\sum_{Z^k \in \mathcal{A}} c_k Z^k, 
$$}\unskip
for some uniformly bounded coefficients $c_k$, and the $Z$ vector
fields commute with the wave operator from Lemma \ref{lem:kvf},  the
improved decay is clear. Using the extra decay, an estimate similar to
(\ref{es:decay}) and global existence can then be derived with some
minor modifications\footnote{At top order, it is harder to fully
exploit the null condition, so the top order norm $E_N$ is often
allowed to grow slowly in $t$, while lower order norms stay uniformly
bounded. See for instance \cite[inequality (5.17) in Chapter 5, p.
78]{Sogge}. Stronger top order estimates without growth can in fact
also be proven using additional spacetime estimates. Those are related
to the so-called \textit{ghost weight} of Alinhac~\cite{MR1856402}.
Although the equations take on a different form, the top order
estimates in~\cite{CKM} are also free of additional growth.}.  

\subsection{The Einstein equations in wave gauge}

Unfortunately, as noted by Choquet--Bruhat in~\cite{MR0315300}, the
Einstein equations in the wave gauge fail to verify the classical null
condition. In fact, she showed later that there is no generalization of
the classical null condition for the Einstein
equations~\cite{andp200051203508},
independently of the choice of coordinates. 

The failure of the null condition for the $Q_{\alpha \beta}$ quadratic
terms in (\ref{ee:red}) is actually relatively easy to capture. Indeed,
the $Q_{\alpha \beta}$ terms can be decomposed into a sum of null
forms, that verify the null condition, cubic terms and higher order
terms which will decay fast, and extra $P$ quadratic forms of the form
{$$
P_{\alpha \beta} (\partial h, \partial h ) = \eta^{\rho \rho'}
\eta^{\sigma \sigma'}\partial_{\alpha} h_{\rho \sigma} \partial_\beta
h_{\rho'\sigma'} - \tfrac{1}{2}\eta^{\rho \rho'} \eta^{\sigma
\sigma'}\partial_{\alpha} h_{\rho \rho'} \partial_{\beta} h_{\sigma
\sigma'}, 
$$}\unskip
where $h= g-\eta$ is the pertubation of the metric $g$ compared to the
Minkowski metric $\eta$.  The above expression for $P$ does not verify
the classical null condition, for instance $P_{tt}$ will contain
products of $t$-derivattives that do not factorize with a tangential
derivative to the light-cone. Nonetheless, it still contains some
structure, apparent in the various contractions between $\eta$ and the
two copies of $\partial h$. Moreover, from the wave gauge condition,
one derive 
{$$
\eta^{\alpha \beta} \partial_{\alpha} h_{\beta \mu} = \eta^{\alpha
\beta} \tfrac{1}{2} \partial_{\mu} h_{\alpha \beta} +
\mathrm{h.o.t}, 
$$}\unskip
from which one can rewrite non-tangential derivatives of certain
components of $h$ as tangential derivatives. All together, one can show
that while the system does not verify the null condition, and as a
consequence, that some components of $h$ will not behave as in the
linear case and suffer a small (at least logarithmic) growth, this
growth does not propagate through the whole system thanks to a
\textit{nilpotent} type structure. These weaker structural properties,
called the \textit{weak null structure} by
Lindblad--Rodnianski~\cite{MR1994592}, allowed them to eventually
provide an alternative proof of stability of Minkowski spacetime, first
for Schwarzschildean initial data in~\cite{LR05} and then general
asymptotically flat initial data in~\cite{LR10}. 

\begin{rem}
\upshape
This wave gauge approach explained above is less geometrical than the
original proof of Christodoulou--Klainerman.  In particular, geometric
conclusions of the stability of Minkowski spacetime, such as curvatures
estimates, the Bondi Mass law, or even the completeness of null
infinity\footnote{See Chapter 17 in~\cite{CKM}.}, which are rather
direct consequences of the global existence and the estimates obtained
in~\cite{CKM}, require an extra analysis on top of the global existence
of solutions in wave coordinates. On the other hand, the fact that
global existence (from the point of view of pdes) can be obtained
independently of a strong control on the geometry is somehow a
testament of the flexibility of the method and can be advantageous. See
Section~\ref{se:desm} below for many developments based on either of
the two strategies.
\end{rem}

\section{The Christodoulou--Klainerman geometric approach to the
stability of Minkowski space} \label{se:ckp}

We describe in this section some aspects of the proof of~\cite{CKM}.
This proof is arguably too complex to summarize accurately in a few
pages, so we have chosen to focus this section mostly on geometric
notions and the overall strategy.

\subsection{The metric and the time foliation}

As explained in Section~\ref{se:bm}, one is forced to make gauge
choices to study the Einstein equations as an evolution problem and we
start by introducing the one used in~\cite{CKM}. To this end, we
consider $(\mathcal{M}, g)$, the maximal globally hyperbolic
development of evolution of strongly asymptotically flat initial data
posed on $\mathbb{R}^3$ as introduced in Sections \ref{se:ceivp} and
\ref{se:af}. By global hyperbolicity, it follows that
$\mathcal{M}\simeq \mathbb{R} \times \mathbb{R}^3$ and in particular,
one can introduce a time function $t$ and coordinates $(x^1, x^2, x^3)$
relative to which the metric takes the form
{$$
g = -\phi^2 (t,x) \,\mathrm{d}t \otimes \,\mathrm{d}t + g_{ij}(t,x)
\,\mathrm{d}x^i \otimes \,\mathrm{d}x^j. 
$$}\unskip
The function $t$ can be chosen so that its level sets $\Sigma_t :=\{ p
\in \mathcal{M},\, t(p)=t\}$ provide a foliation of $\mathcal{M}$, the
normalized condition $\phi \rightarrow 1$ as $x\rightarrow \infty$ is
verified on each $\Sigma_t$, and such that $\Sigma_0= \left\{
t=0\right\}$ coincides with the\footnote{More precisely, $\Sigma_0$
coincides with the image of the initial data set by the embedding
introduced in the definition of a development, cf.\ 
Section~\ref{se:ceivp}.} initial data set. We denote by  $\mathbf{D}$
the Levi-Civita connexion of $g$ and by $\nabla$ the one associated to
the spacelike metric $g_{ij}$. 

Let $k_{ij}$ denotes the component of the second fundamental form
associated to the $\Sigma_t$ foliation
{\begin{equation} \label{def:k}
k_{ij} := \frac{1}{2 \phi} \partial_t g_{ij}. 
\end{equation}}\unskip
To decompose the tensor fields, it is useful to introduced the
orthonormal frame $(T, e_i)$, defined by $T = ({1}/{\phi}) \partial_t$
while the $e_i$ are tangential to $\Sigma_t$ and Fermi propagated along
$T$, so that 
{$$
{\mathbf{D}}_T\, e_i = \phi^{-1} (\nabla_i \phi ) T. 
$$}\unskip

\subsection{The Weyl tensor}

The Weyl tensor associated to $g$ is the traceless part of the
curvature tensor
{$$
{\mathbf{C}}_{\alpha \beta \gamma \delta} := 
{\mathbf{R}}_{\alpha \beta \gamma \delta} - \tfrac{1}{2} 
(g_{\alpha \gamma}
{\mathbf{R}}_{\beta \delta}+g_{\beta \delta}
{\mathbf{R}}_{\alpha \gamma}-g_{\beta \gamma}
{\mathbf{R}}_{\alpha \delta}-g_{\alpha \delta}
{\mathbf{R}}_{\beta \gamma}) + \tfrac{1}{6}
(g_{\alpha \gamma} g_{\beta \delta}-g_{\alpha \delta}g_{\beta \gamma} ) 
{\mathbf{S}}(g).
$$}\unskip
In particular, for a solution to the Einstein equations, the Weyl
tensor coincides with the curvature tensor. It verifies the Bianchi
equations, 
{\begin{equation}
{\mathbf{D}}_{\mathop{[} \epsilon} {\mathbf{C}}_{\alpha \beta\mathop{]} \gamma \delta} = 0,
\end{equation}}\unskip
and, on top of being traceless by construction, enjoys the same
symmetry properties as the curvature tensor. 

Instead of considering all the components of the Weyl tensor, it is
preferable to introduce an irreducible decomposition, the
\textit{electro-magnetic} decomposition, which takes into account its
symmetries. Given a vector field, here the unit normal $T$, one
considers the tensors $E$ and $H$, defined by 
{\begin{eqnarray*}
E_{\alpha \beta} &amp;=&amp; {\mathbf{C}}_{\mu \alpha \nu \beta} T^{\mu} T^{\nu}, \\
H_{\alpha \beta} &amp;=&amp; {\,}^{\star} 
{\mathbf{C}}_{\mu \alpha \nu \beta} T^{\mu} T^{\nu}, 
\end{eqnarray*}}\unskip
where ${\,}^{\star} C$ is the \textit{Hodge} dual of $C$
{$$
{\,}^{\star} {\mathbf{C}}_{\mu \alpha \nu \beta} = 
\tfrac{1}{2}\varepsilon_{\mu \alpha \rho \sigma}
{\mathbf{C}}^{\rho \sigma}_{\phantom{ab}\nu \beta}, 
$$}\unskip
with $\varepsilon_{\mu \alpha \rho \sigma }$ denoting the components of
the volume form. 

From the symmetric tensors $E$ and $H$, one can reconstruct the whole
Weyl tensor and moreover, from the Bianchi equations, it follows that
$(E,H)$ verifies a symmetric hyperbolic system similar to the Maxwell
equations. For instance, the evolution equation for $E$ takes the form
{\begin{eqnarray*}
-\phi^{-1}\partial_t E_{ab} + \mathrm{curl}\, H_{ab} =  Q_{ab}(\nabla
\phi, E, H, k), 
\end{eqnarray*}}\unskip
where 
$\mathrm{curl}\,H_{ab} := ({1}/{2})
({\varepsilon}_{a}^{cd}\nabla_d H_{cb} +
 {\varepsilon}_{b}^{cd}\nabla_d H_{ca})$ and $Q_{ab}$ is a
quadratic in $\nabla \phi, E, H, k$, which is similar to the
differential Amp\`ere--Maxwell law. 

\subsection{Maximal foliations and elliptic systems} \label{se:max}

One important ingredient of the setup is to consider the Weyl tensor as
one of the the main unknown and only recover a control on the metric
via the definition of $k$, Equation (\ref{def:k}) and the Gauss and
Codazzi equation associated to the hypersurfaces $\Sigma_t$. However,
one still needs to fix the foliation, that is to constrain the
hypersurfaces $\Sigma_t$. Here, one considers \textit{maximal}
hypersurfaces\footnote{In particular,~\cite{CKM} requires the initial
data to already be maximal, i.e.~that $\mathrm{tr}\, k =0$ on the
initial surface. The existence of such complete maximal hypersurfaces
in general, asymptotically flat Lorentzian manifold (under a barrier
assumption) was established by Bartnik in~\cite{bart:max}.}, i.e.\ one
imposes $\mathrm{tr}\,k=0$ on each $\Sigma_t$. With this choice, the
lapse function $\phi$ verifies the elliptic equation
{$$
\Delta \phi = |k| ^2 \phi. 
$$}\unskip
Moreover, $k$ also verifies an elliptic div-curv type system 
{\begin{eqnarray*}
\mathrm{curl}\, k_{ij} &amp;=&amp; H_{ij}, \\
\mathrm{div}\, k &amp;=&amp; 0,
\end{eqnarray*}}\unskip
where $H$ is the magnetic part of the electro-magnetic decomposition of
the Weyl tensor. Thus, one sees that a good control on the Weyl tensor
leads via elliptic estimates to a good control of the second
fundamental form and all the metric quantities. 

\subsection{The Bel--Robinson tensor and the vector field method applied
to the Bianchi equations} \label{se:brt}

At the linearized level, one thus needs to understand the following
problem: given a tensor field $W:=W_{\alpha \beta \gamma \delta}$ that
verifies all the symmetry and traceless properties of the Weyl tensor
as well as the Bianchi equations
{$$
\partial_{\mathop{[} \epsilon} W_{\alpha \beta\mathop{]} \gamma \delta}=0,
$$}\unskip
on Minkowski space, how can one prove estimates on $W$? In fact, to
prove a stability statement, similar to the approach highlighted in
Section~\ref{se:qlns}, one needs both energy type estimates at the top
order and decay estimates. Moreover, having in mind the importance of
the null condition (cf.\ Section~\ref{se:nc}) and thus of structural
properties of the non-linearity,  we will see here, already at the
level of the components of $W$, that different components of $W$ (akin
to different derivatives for scalar waves) enjoy different decay
properties. 

First, the energy estimates are proven using the so-called
\textit{Bel--Robinson} tensor $Q$
{$$
Q_{\alpha \beta \gamma \delta} :=  W_{\alpha \mu \beta \nu} 
W_{\gamma\,\,\delta}^{\mu\,\,\nu}+ {\,}^{\star} W_{\alpha \mu \beta \nu}
{\,}^{\star} W_{\gamma\,\,\delta}^{\mu\,\,\nu}.
$$}\unskip
The Bel--Robinson tensor plays the role of the standard energy--momentum
tensor for matter fields. It is divergence free for $W$ solving the
Bianchi equations and it is coercive whenever contracted with $4$
identically directed timelike vector fields.

On Minkowski space, a vector field method leading to energy estimates
and weighted Sobolev inequalities, which encodes the decay of the
fields, similar to the one for the scalar field presented in
Section~\ref{se:sk} can then be implemented. These ideas were first
presented in~\cite{cpa3160430202}. Compared
with the simplest case of a scalar field of Section~\ref{se:sk}, a
notable difference comes the fact that, since $W$ is a tensor field, it
has different components, and the different components enjoy different
decay properties. These different decay properties are related to the
directions of propagation of the evolution equations, and thus to the
null conditions. To capture them, it is thus necessary to decompose $W$
into components tangential or transversal to the null cones. Let thus
$L= \partial_t + \partial_r$, $\underline{L}= \partial_t - \partial_r$
and complete $L, \underline{L}$ by a pair of orthonormal vector fields
$e_A$, $A=1,2$, tangential to the spheres of constant $(t,r)$. We then
decompose $W$ on the null frame $(L, \underline{L}, e_A)$: 
{\begin{eqnarray*}
\underline{\alpha}_{AB} &amp;:=&amp; W_{A \underline{L} B \underline{L}}, \quad 
\alpha _{AB}:= W_{A L B L}, \\
\underline{\beta}_A &amp;:=&amp; \tfrac{1}{2}W_{A \underline{L} \underline{L} L}, \quad 
\beta _A := \tfrac{1}{2}W_{A L \underline{L} L}, \\
\rho &amp;:=&amp; \tfrac{1}{4}W_{\underline{L} L \underline{L} L}, \quad 
\sigma := \tfrac{1}{4}\epsilon^{\alpha \beta }W_{\alpha \beta \underline{L} L}. 
\end{eqnarray*}}\unskip
Using the Morawetz vector field $K_0$ (cf.~Remark \ref{rem:mor}), which
rewrites in the null frame as $K_0 = ({1}/{2}) ( v^2 L + u^2
\underline{L})$ and the Bel--Robinson tensor $Q$, one can then
compute the energy density
{\begin{eqnarray*}
Q( K_0, T,  T, T) = 
\tfrac{1}{8} u^2 |\underline{\alpha}|^2+ 
\tfrac{1}{8} v^2 |\alpha|^2+ 
\tfrac{1}{4} (3 u^2 +v^2) |\underline{\beta}|^2 + 
\tfrac{1}{4} (3 v^2 +u^2) |\beta|^2 + 
\tfrac{3}{4} (u^2+v^2) (\rho + \sigma),
\end{eqnarray*}}\unskip
where $u=t-r$, $v=t+r$, from which the different weights for the
different components are already apparent. 

\subsection{The Eikonal equation and geometrically defined vector
fields}

In implementing the above estimates on an asymptotically flat spacetime
which is not Minkowski space, one important difficulty arises from the
long-range effect of the $1/r$ term of the initial data,
cf.~(\ref{eq:metrdata}), the \textit{Schwarzschildean} tail of the
metric. Indeed, the geodesic null cones of the spacetime then diverge
logarithmically from the corresponding ones of Minkowski space, and
their shear, which vanishes in Minkowski space, also has bad asymptotic
properties. To avoid these divergences from spoiling the estimates, the
strategy of~\cite{CKM} is to adapt the vector fields which are used in
the construction of the energy norms, using the geometry of the current
spacetime rather than that of the fixed Minkowski space.  For instance,
the vector field $K_0$ is replaced by 
{\begin{equation}\label{def:mor}
K:= \tfrac{1}{2} (v^2 e_+ + u^2 e_-),
\end{equation}}\unskip
where 
\begin{enumerate}[(1)]
\item $u$ is an \textit{optical} function, that is to say, for $r \ge
r_0/2$, it solves\footnote{To avoid caustics, $u$ is actually
constructed as a solution to (\ref{eq:eik}) in an exterior region, $r
\ge r_0/2$, for some $r_0(t)&gt;0$ and matched to an another optical in
the supplementary interior region.} the Eikonal equation, 
{\begin{eqnarray}\label{eq:eik}
g (\mathrm{d}u, \mathrm{d}u) = 0, 
\end{eqnarray}}\unskip
\item $e_\pm$ are null vectors: $e_\pm = T\pm N$, where $T$ is the
forward unit normal to $\Sigma_t$ and $N$ is the outward unit normal,
tangential to $\Sigma_t$, to the spheres $S_{t,u}$ of constant $(t,u)$.
\item $v= 2r -u$, where $r$ is the area radius of the $S_{t,u}$
spheres. 
\end{enumerate}
Similarly, the angular momentum vector fields $\Omega_{ij}$ replacing
the flat ones (cf.~Lemma \ref{lem:kvf}) are also modified, for instance,
we want $\Omega_{ij}(u)=0$, and thus since $u$ is now constructed
dynamically, so must the $\Omega_{ij}$. 

\subsection{Curvature estimates: commutation, error terms and the
principle of signature} 

As explained in Section~\ref{se:brt}, the components of the Weyl
curvature tensor are estimated using the Bel--Robinson tensor. Moreover,
to obtain higher order estimates as well as obtained improved decay,
one also needs to estimate components of derivatives of $W$ such as
$\widehat{\mathcal{L}}_T W$ where $\widehat{\mathcal{L}}_X W$ for a
vector field $X$ is the modified Lie derivative\footnote{The point of
the modified Lie derivative is that it has the property of preserving
the tracefreeness of Weyl fields: $\widehat{\mathcal{L}}_X
W^{\gamma}_{\beta \gamma \delta }=0$.} defined by
{\begin{eqnarray*}
\widehat{\mathcal{L}}_X W_{\alpha \beta \gamma \delta } &amp;:=&amp;
\mathcal{L}_X W_{\alpha \beta \gamma \delta} - \tfrac{1}{2}
\left(
{\,}^{(X)} \uppi^{\mu}_{\alpha}W_{\mu \beta \gamma \delta}+
{\,}^{(X)} \uppi^{\mu}_{\beta} W_{\alpha\mu  \gamma \delta}+
{\,}^{(X)} \uppi^{\mu}_{\gamma}W_{\alpha \beta \mu  \delta}+
{\,}^{(X)}\uppi^{\mu}_{\delta}W_{\alpha\beta  \gamma \mu }
\right) \\
&amp;&amp; -\, \tfrac{3}{8} \mathrm{tr}^{(X)} \uppi W_{\alpha \beta \gamma \delta},
\end{eqnarray*}}\unskip
where ${}^{(X)}\uppi_{\alpha \beta} := ({1}/{2}) \mathcal{L}_X g_{\alpha
\beta}$, with $\mathcal{L}_X$ the standard Lie derivative acting on
tensor fields, is the deformation tensor of the vector field $X$. 

In the case of the $T$ vector field, commuting the Bianchi equations
${\mathbf{D}^\alpha }  W_{\alpha \beta \gamma \delta}=0$, one obtains 
{$$
{\mathbf{D}^\alpha } (\widehat{\mathcal{L}}_T W_{\alpha \beta
\gamma \delta}) = J(T;W)_{\beta \gamma \delta}, 
$$}\unskip
for some error term $J(T;W)$ which is a quadratic form in the
deformation tensor of $T$ and the components of $W$. 

Consider the Bel--Robinson tensor $Q(\widehat{\mathcal{L}}_T W)$
associated to $\widehat{\mathcal{L}}_T W$. Because of the additional
$T$ derivatives, the initial data of $\widehat{\mathcal{L}}_T W$ decay
faster than that of $W$, so that $Q(\widehat{\mathcal{L}}_T W)(K, K, K,
T)$, where $K$ is the Morawetz vector field, has finite integral on the
initial hypersurface, despite the strong weights in the definition of
$K$, cf.~(\ref{def:mor}). The divergence theorem applied between the
hypersurfaces $\Sigma_t$ and $\Sigma_0$ leads to 
{\begin{eqnarray*}
\int_{\Sigma_t} Q(\widehat{\mathcal{L}}_T W)(K, K, K, T) &amp;=&amp;
\int_{\Sigma_0} Q(\widehat{\mathcal{L}}_T W)(K, K, K, T) + 
\int_0^t \,\mathrm{d}t' 
\int_{\Sigma_{t'}} \phi\, {\mathrm{div}}\, Q
(\widehat{\mathcal{L}}_T W)_{\beta \gamma \delta} K^\beta K^\gamma
K^\delta \\
&amp;&amp; +\, \frac{3}{2}
\int_0^t \,\mathrm{d}t' 
\int_{\Sigma_{t'}} \phi\, {Q} (\widehat{\mathcal{L}}_T W)_{\alpha \beta \gamma \delta} 
{\,}^{(K)}\uppi^{\alpha\beta} K^\gamma K^\delta,
\end{eqnarray*}}\unskip
where the first spacetime error term comes from commuting the Bianchi
equations and thus depends on $J(T;W)$ while the second comes from the
fact that the vector field $K$ is not an exact conformal Killing field.
Similar integral formulas hold for other derivatives of $W$ with the
amount of $K$ vector fields contracted with the Bel--Robinson tensor
depending on how fast the data decay. To prove energy estimates for the
components of $W$ and its derivatives, one thus needs to control the
error terms using decay estimates, and to avoid logarithmic
divergences, to exploit the special structure of certain products,
similar to Sections~\ref{se:qlns} and \ref{se:nc}. Each of these error
terms can be expressed as a product of null components of the
deformation tensors of $K$ or $T$ and the null components of
$\widehat{\mathcal{L}}_T W$ and one needs to exploit the different
decay or integrability rates of each term. Since there are many of
those, a convenient, systematic, approach is given by the
\textit{principle of conservation of signature}: given any tensor
tangential to the spheres $S_{t,u}$, obtained by contraction of a
spacetime tensor via contraction with the null vectors $e_3$ or $e_4$
and projections, one assigns a signature defined as the difference
between the total number of contractions with the outgoing null vector
$e_4$ and the total number of contractions with the ingoing vector
$e_3$. For instance, $\alpha_{AB}= W(e_4, e_A, e_4 , e_B)$ is given the
signature $s(\alpha)=+2$ and $\underline{\alpha}_{AB}= W(e_3, e_A, e_3
, e_B)$ the signature $s(\underline{\alpha})=-2$.

By the principle of conservation of signature, given any tensor field
$U$ which can be written as a multilinear form in an arbitrary number
of tensors $U_1,\ldots,U_p$ with coefficients depending only the
metric, the signature of any null term of $U$, expressed in terms of
the null components of $U_1,\ldots,U_p$ is equal to the sum of the
signatures of each term in the decomposition. As an illustration, the
component $Q_{AB34}$ of the Bel--Robinsor tensor has signature $0$. By
definition, $Q$ is a quadratric form in the null curvature components
of a Weyl field, and by the principle of conservation of signature, all
quadratic terms must have signature $0$, leaving as only option
products of the form $\alpha \cdot \underline{\alpha}$, $\beta \cdot
\underline{\beta}$,  $\rho \cdot \rho$, $\rho \cdot \sigma$ and $\sigma
\cdot \sigma$. This can used in a (quasi)-systematic to obtain the
correct form of each error term. 

\subsection{Bootstrap assumptions and the last slice argument}

Similar to Section~\ref{se:qlns}, the overall strategy is based on a
bootstrap argument. Thus, one assume the existence of $t^{\star}$ such
that within $[0,t^{\star}]$, quantitative estimates, both pointwises and
in energy type norms holds. As an illustration, the $L^2$ exterior top
order bootstrap assumptions on $\alpha$ and $\underline{\alpha}$ take
the form
{$$
\Vert (1+u^2)^2 r^q \slashed \nabla^q \underline{\alpha}
\Vert _{L^2( r \ge r_0/2)} + \Vert r^{q+2} \slashed \nabla^q  \alpha
\Vert _{L^2( r \ge r_0/2)} \le \epsilon_0, \quad q=0, 1, 2,
$$}\unskip
and the aim is then to improve these bounds and all the other bootstrap
assumptions.

One important twist on the bootstrap argument concerns the optical
function $u$. The Eikonal equation is a non-linear transport equation
and thus, solutions are determined by their initial data.
In~\cite{CKM}, in the exterior region, $u$ is initialized from the
future, that is to say, from the last slice $\Sigma_{t^{\star}}$. Indeed,
fixing $u$ on a spacelike slice is equivalent to choosing a radial
foliation of the slice. If this is done on the initial slice, since $u$
verify a transport equation, there is no reason why the $S(t,u)$ would
get rounder as $t \rightarrow +\infty$ which would forbid the
improvement of certain bootstrap assumptions. Instead, one initializes
$u$ on the last slice $t^{\star}$ by first solving for a radial foliation
on $\Sigma_{t^{\star}}$ verifying improved bounds. In particular, the
optical function ${}^{(t^{\star})} u$ depends a priori on the time-slab
$[0,t^{\star}]$. After the bootstrap assumptions are improved, the local
existence provides an extended time-slab $[0,t^{\star}+ \delta]$,
$\delta&gt;0$, and thus a new optical function ${}^{(t^{\star}+\delta)}u$
must be constructed with initialization on $\Sigma_{t^{\star}+\delta}$.
In the limit, one then also needs to prove the convergence of the
optical functions ${}^{(t)}u$ as $t \rightarrow + \infty$. 

\section{Developments and extensions of the stability of Minkowski
space} \label{se:desm}

We present here a survey of result which extends or complete the
stability of the Minkowski space. The exposition is thematical rather
than historical. 

\subsection{Stability under weaker assumptions on the initial data:
Bieri and Shen works} \label{se:bs}

The initial data in the original work of Christodoulou--Klainerman were
assumed to be strongly asymptotically flat, cf.\ Section~\ref{se:af}.
In~\cite{MR2531716}, Bieri proved the stability of the Minkowski space
for initial data verifying 
{\begin{eqnarray*}
g_{ij} &amp;=&amp; \delta_{ij} + o_{q+1}(r^{-\frac{s-1}{2}}), \\
k_{ij} &amp;=&amp; o_{q+1}(r^{-\frac{s+1}{2}}),
\end{eqnarray*}}\unskip
with $s=2$. The number of derivatives on the initial data required was
also lower to $q=2$, compared to up to $4$ derivatives in~\cite{CKM}.
Note that under this assumption, one cannot apriori define the ADM
mass, which requires $s &gt; 2$.  This result has been recently improved
in the work of Shen~\cite{shen2023globalstabilityminkowskispacetime},
which requires only $s&gt;1$. A main difference of the analysis
in~\cite{shen2023globalstabilityminkowskispacetime} compared
to~\cite{MR2531716} and~\cite{CKM} comes from the construction of the
optical function, which is no longer initialized on a last slice, but
rather on the initial slice and a timelike curve playing the role of
the axis of symmetry. Another analytical difference comes from
abandoning of the classical vector field method, in particular, the
weighted norms obtained via contraction with the $K$ vector field,
which is replaced by the $r^p$-weighted method of
Dafermos--Rodnianski~\cite{97898143046340032}. 

\subsection{The double null gauge and the exterior problem}

With spatial dimension $3$, due to the slow decay of waves near the
light-cone $t=|x|$, many steps in the analysis of quasi-linear wave
equations, including for the stability of Minkowski space, must treat
differently an interior region of the form $r \le ({r_0}/{2})(t)$, from
the exterior region $r \ge ({r_0}/{2})(t)$, where $r_0(t)$ is the
area-radius of the spheres of intersection between the initial cone
$u=0$ and $\Sigma_t$. In the exterior region, it can be advantageous to
not consider the Einstein equations in the maximal gauge as
in~\cite{CKM}, but in the \textit{double null gauge}, where the
spacetime metric can be put under the form
{\begin{eqnarray*}
g:= - 2 \Omega^2 (\mathrm{d}u\otimes \,\mathrm{d}\underline{u} + 
\,\mathrm{d}\underline{u}\otimes \,\mathrm{d} u) + 
\gamma_{ab} (\mathrm{d}\omega^a -X^a \,\mathrm{d}u) \otimes
(\mathrm{d}\omega^b -X^b \,\mathrm{d}u).
\end{eqnarray*}}\unskip
Here, both $u$ and $\underline{u}$ are optical functions, i.e.~solutions
to the Eikonal equations, while $\omega^a$ are coordinates on the
spheres $S_{u, \underline{u}}$ of constants $(u, \underline{u})$. The
level sets of $u$ and $\underline{u}$, $C(u)$,
$\underline{C}(\underline{u})$ are characteristic cones and thus the
double null gauge is naturally adapted to the equations. In particular,
it captures better the null properties of the equations and it is also
advantageous when studying radiation phenomena taking place on the
asymptotic null cone $C_{\underline{u}}$,\unskip\break $\underline{u}
\rightarrow + \infty$. 

The exterior problem can then be understood as follows: given an
asymptotically flat initial data set $(\Sigma, h, k)$, such that for
some compact set $K$, the initial data on $\Sigma \setminus K$ verifies
a smallness assumption. Let $(\mathcal{M}, g)$ be the maximal Cauchy
development of $(\Sigma, h, k)$, and let $\mathcal{M}_+$ be the domain
of dependence of $\Sigma \setminus K$. Show that $\mathcal{M}_+$ can be
foliated by cones $(C(u)$, $\underline{C}(\underline{u}))$ such that
the outgoing cones $C(u)$ are complete and such that the curvature and
Ricci coefficients are controlled in a quantitative manner. For
instance, for strongly asymptotically flat data, $\alpha$ and
$\underline{\alpha}$ must verify the decay estimates
{\begin{eqnarray*}
r^{7/2} |\alpha | \lesssim C, \quad
r |u|^{5/2} |\underline{\alpha}| \lesssim C,
\end{eqnarray*}}\unskip
where $C$ depends on the initial data and $r$ is the area radius of the
spheres $S_{u, \underline{u}}$. 

The exterior problem was first addressed in the work of
Klainerman--NicolÃ²~\cite{MR1946854} and was recently revisited in the
work of Shen~\cite{shen2023stabilityminkowskispacetimeexterior}, which
in particular generalized the class of initial data considered, as well
as Hintz~\cite{MR4645018},  whose work is based on generalized wave
coordinates instead of the double null gauge. 

The double null foliation is particularly well-adapted to the study of
characteristic initial value problems, cf.~\cite{imrnrnr201}.
Together with the maximal gauge and the wave gauge and their
generalizations, it constitutes now one of the classical ways to write
the Einstein equations, with numerous spectacular applications, such
as~\cite{MR2488976}. 

As a complement to the exterior problem, let us also mention the
\textit{spacelike-characteristic} initial value problem, which one
could coin the interior problem, in which the the initial data is posed
on a compact 3-disk and on a future complete null hypersurface
emanating from its boundary. The global stability of the Minkowski
space in this context was proved by
Graf~\cite{graf2021globalnonlinearstabilityminkowski} based on a gauge
which is partially initialized on a timelike curve, as
in~\cite{shen2023globalstabilityminkowskispacetime}. 

\subsection{The $1+2$-dimensional case}

In~\cite{MR3498988} and~\cite{MR3805301}, Huneau studied the stability
of the Minkowski space with a translational symmetry. Due to the
translational invariance, the data are no longer asymptotically flat as
in Section~\ref{se:af} since they enjoy no decay in the direction of
symmetry. Thanks to the symmetry, the equations reduce to a
quasi-linear system of wave equations in $1+2$ dimensions. Since the
decay of waves improves with the number of spatial dimensions, the
setting is harder than the $1+3$ case. General quasi-linear wave
equations verifying the null condition were first studied by
Alinhac~\cite{MR1856402}. The work of Huneau must of course address all
the extra difficulties coming from the Einstein equations. In the $1+3$
dimensional case, to address the presence of the Schwarzschildean tail
of the metric, the term of order $1/r$ proportional to the mass,
in~\cite{MR2680391}, an ansatz for the perturbation was introduced to
solve only for the remainder terms. In the work of Huneau, another
ansatz must be introduced, since the explicit Schwarzchild solution no
longer captures the asymptotics. Instead, here, the background metrics
are no longer explicit and Huneau constructs them based on the (radial)
Einstein--Rosen waves.  Another important modification compared
to~\cite{MR2680391} concerns the usage of \textit{dynamical,
generalized} wave coordinates, where the coordinates $x^\alpha$ now
verify inhomogeneous wave equations instead of (\ref{eq:wcc}), with an
inhomogeneous term that depends on a non-trivial background metric as
well as the actual metric (but not its derivatives).  This allows to
substract off some of worse source terms in the wave equations for the
metric\unskip\break perturbations.

\subsection{The stability of Minkowski space for Einstein-matter systems} 
\label{se:ems}
\subsubsection{The Einstein equations coupled to wave equations}

One natural way to extend the stability results of the Minkowski space
is to consider the Einstein equations coupled to matter. The matter
models which are also governed by wave or tensorial wave equations are
clear first candidates to be studied, since the additional fields can a
priori be treated by the methods already developed for the Einstein
equations. The Maxwell case was addressed by Zipser in~\cite{MR2531716}
based on the Christodoulou--Klainerman framework and subsequently
Loizelet~\cite{MR2214599} using the wave coordinates approach of
Lindblad--Rodnianski.  The Einstein-scalar field is actually contained
in the analysis of Lindblad--Rodnianski~\cite{MR2680391}, which
interestingly, thanks to the possible matter fields, shows that the
positivity of the ADM mass is logically independent from the stability
properties of the Minkowski space. Indeed, one can produce initial data
with negative mass by considering the Einstein-scalar field equations
with an negative sign in front of the energy--momentum tensor in
(\ref{eq:EE}) and yet the main theorem of~\cite{MR2680391} still
applies. 

In~\cite{MR3254347}, Speck studied the Einstein equations coupled to
non-linear electromagnetism, which contains additional difficulties due
to the intrinsic non-linearities of the matter equations, even on a
fixed background.   Recently, Kauffman and Lindblad have considered the
case of the Einstein equations coupled to a charged scalar
field~\cite{MR4653401}. To treat the case of non-vanishing charge, very
precise estimates on certain metric coefficients are needed, that go
beyond those contained in~\cite{MR2680391}. In the vacuum case, those
were derived by Lindblad~\cite{MR3638312} where a full set of sharp
asymptotics for the metric components is obtained in generalized wave
coordinates. 

\subsubsection{Coupling between massless and massive wave equations:
the Einstein--Klein--Gordon system}

In the Einstein--Klein--Gordon system (EKG) system, the Einstein
equations are coupled to a massive wave equation. At the linear level,
the Klein--Gordon field $\psi$ verifies the Klein--Gordon equation
{\begin{eqnarray} \label{eq:kg}
\Box \psi+m^2 \psi =0, 
\end{eqnarray}}\unskip
with mass $m &gt;0$. Equation (\ref{eq:kg}), while still sharing with the
same principal symbol with the wave equation, does not respect the
conformal symmetries of Minkowski space due to the extra mass term. As
a consequence, the conformal vector fields, such as the scaling $S$, do
not enjoy good commutation properties with the equation. Yet, an
alternative vector field approach for the study of (\ref{eq:kg}) was
developed in~\cite{MR1199196}, exploiting in particular a novel
foliation of the forward light-cone by hyperboloidal spacelike
hypersurfaces. The coupling to the Einstein equations is still
non-trivial: if one were to strictly apply the methods of
either~\cite{CKM} or~\cite{MR2214599} in (\ref{eq:EE}), commuting
vector fields would naturally land over to the energy--momentum tensor
$T$ and force us a priori to also commute the matter equations via the
same vector fields. As a consequence, one is forced to also modify the
analysis of the Einstein equations in a non-trivial way. The stability
of the Minkowski space for the EKG system was first obtained in the
work of LeFloch--Ma in~\cite{MR3535896} for Schwarschildean initial data
and then in~\cite{lm:nssgmf} for data with more general asymptotics,
using an extension of the Klainerman hyperboloidal approach. An
independent proof was obtained by Ionescu--Pausader~\cite{MR4422074},
which includes slowly decaying data, based on a novel framework for the
study of quasi-linear wave equations, incorporating both Fourier
analysis and physical space methods. In~\cite{1586224841},
Wang initiated an approach to this problem based on intrinsic
hyperboloids, instead of using the hyperboloidal foliation of the
background Minkowski space. 

The methods developed for the study of systems of wave-Klein--Gordon
equations have found further applications, such as in the work of
LeFloch--Ma addressing $f(R)$ theories of gravity~\cite{MR4750422},
where the Einstein equations are replaced by a fourth order system of
wave-Klein--Gordon equations, or the analysis of the Kaluza--Klein
spacetime by Huneau, Stingo and Wyatt on the stability of the
Kaluza--Klein
spacetime~\cite{huneau2023globalstabilitykaluzakleinspacetime}. Here
the analysis reduces to a quasi-linear system of wave equations on a
product space $\mathbb{R}^{1+3} \times \mathbb{S}^1$, corresponding to
$1+4$ dimensional Minkowski spacetime with one direction compactified
on a circle. The resulting system can then be viewed as coupling wave
equations to an infinite sequence of Klein--Gordon equations with
different masses. 

\subsubsection{The Einstein--Vlasov system}

Colissionless kinetic matter is a popular matter model widely used in
astropphysics or cosmology. In this case, the Einstein equations are
coupled to a kinetic equation, of unknown $f$ called a \textit{Vlasov}
field, representing a distribution of particles, which is defined on a
hypersurface of the tangent bundle
{\begin{eqnarray}
f: \mathcal{P} \rightarrow \mathbb{R}_+ , \quad 
\mathcal{P} := \{(x,p) \in T \mathcal{M},\,\, g_{x}(p,p)=-m^2, \,
p\text{ future-oriented}\}.
\end{eqnarray}}\unskip
Here, $m$ denotes the mass of the particles and the massless case $m=0$
has different analytical properties from the massive one $m&gt;0$.  The
evolution equation for $f$ is that it must be propagated by geodesics,
so that it verifies the transport equation 
{\begin{eqnarray} \label{eq:ve}
v^\alpha \partial_{x^\alpha} f - p^\alpha p^\gamma
\Gamma^{\beta}_{\alpha \gamma} \partial_{p^\beta} f = 0.
\end{eqnarray}}\unskip
The energy--momentum tensor in (\ref{eq:EE}) is given by the second
moment of $f$. For instance, in the flat case $g=\eta$, $\mathcal{P}$
can be parametrized by 
{\begin{eqnarray}
\mathbb{R}^{1+3}_{x} \times \mathbb{R}^3_p &amp; \rightarrow &amp; \mathcal{P}
\Seqnsplit
( x, p ) &amp; \rightarrow &amp; \left(x, (p^0, p^i):=\left(\sqrt{ m^2+|p|},
p\right) \right)
\end{eqnarray}}\unskip
and the energy--momentum tensor of a Vlasov field can then be written as
{$$
T_{\alpha \beta}[f] := \int_{(p_i) \in \mathbb{R}^3} p_\alpha p_\beta f
\frac{\mathrm{d}p}{p^0}.
$$}\unskip
A novel difficulty comes here from the fact that the Einstein equations
and the kinetic one (\ref{eq:ve}) are of quite different nature a
priori, and thus, it is not clear that they can be treated by a common
approach. For instance, if one tries to commute (\ref{eq:EE}) with
vector fields, it is clear that we will be forced to consider
derivatives of the energy--momentum tensor, and thus, a priori of $f$
and (\ref{eq:ve}).

Within spherical symmetry, the nonlinear stability of the Minkowski
space for the Einstein--Vlasov system was initiated by Rein and
Rendall~\cite{rr-92} in the massive case and Dafermos in the case of
massless particules~\cite{d-06}. Outside from spherical symmetry, the
massless case has been addressed first by Taylor, in which the Einstein
equations are treated as in~\cite{CKM,MR2488976}, in the double null
gauge, while the analysis of the the Vlasov field relied on the
construction of Jacobi fields and assumed compact support on the
initial Vlasov field. Another proof was given in~\cite{MR4302757} by
Bigorgne--Fajman--Joudioux--Smulevici--Thaller, in wave coordinates, in
the case of non-compact support for the Vlasov field. The massive case
has first been addressed in~\cite{MR4241806, MR4062482} by
Fajman--Joudioux--Smulevici and Lindblad--Taylor and more recently by
Wang in~\cite{wang2022globalstabilityminkowskispacetime}, using an
extension of the Fourier analysis techniques of Ionescu--Pausader to
address both the Einstein and the Vlasov part of the equations. 

To address the issue of the coupling between the Einstein and the
Vlasov equations, a systematic, vector field approach was first
developped in~\cite{FJS} and form the core of the strategy
of~\cite{MR4241806} and~\cite{MR4302757}.  Another point of interest
comes from the analysis of the non-linear wave/kinetic field
interactions, where structures similar to the null condition or the
weak null condition have been discovered. For instance,  given $g$ a
solution to the wave equation and $f$ a Vlasov field, the non-linear
terms $\partial_{x^i} g . \partial_{p_i} f$ and $[\partial_{x^i} +
({v^i}/{\sqrt{1+|v|^2}}) \partial_t ](g).  \partial_{p_i} f$ verifies
stronger estimates compared to generic products $\partial_x g .
\partial_p f$, due to identities such as
{\begin{eqnarray}
\partial_{x^i} g  \,\partial_{p_i} = \frac{Z(g)}{t} \, \partial_{v}
- \frac{x}{t}\frac{\partial_t g}{\sqrt{1+|v|^2}}\widehat{Z}+ 
\frac{\partial_t g }{\sqrt{1+|v|^2}} \left( S + \frac{|x|^2-t^2}{t}
\partial_t \right), 
\end{eqnarray}}\unskip
where $Z, S$ (respectively $\widehat{Z}$) are vector fields commuting
with the wave operator (respectively, the relativistic transport
operator on Minkowski space).

\subsection{New and newer vector field methods and other approaches for
stability} \label{newvfm}

After Minkowski space, one would naturally like to understand the
stability properties of the other stationary solutions to the Einstein
equations, which, in $3+1$ dimensions\footnote{Together with the
Minkowski space, it is indeed conjectured that these are the only
stationary solutions of the vacuum Einstein equations~\cite{cch:lr}
with $\Lambda=0$ while the picture in higher dimension or with
non-vanishing $\Lambda$ is considerably more complex.} are given by the
Schwarzschild and Kerr family. To describe all the recent progress on
this fundamental problem goes well beyond the scope of this article,
but we wish only to highlight some developments concerning the approach
of Section~\ref{se:geqw}. 

First, since these geometries have less symmetry than the Minkowski
space, the commutation vector field approach to decay and the strategy
presented in Sections~\ref{se:qlns} and \ref{se:sk} needs to be
considerably altered, see~Remark \ref{rem:mvsc}. 

In the context of the linear scalar wave equation $\Box_{g_{0}}
\psi=0$ on a stationary asymptotically flat geometry, such as the Kerr
or Schwarzshild spacetime, a first important new tool is the derivation
of local integrated decay estimates, which in simple setting can be
stated as
{\begin{eqnarray} \label{eq:ide}
\int_{0}^\infty \int_{ r \le R} \mathcal{E}_{r \le R,\, \chi}
\lesssim \mathcal{E}(0), 
\end{eqnarray}}\unskip
where $\mathcal{E}_{r \le R, \,\chi}(\tau)$ denotes a degenerate
energy flux supported in a compact region $r \le R$ which controls the
first derivative of $\psi$ in $L^2$ multiplied by a function $0 \le
\chi(r)\le 1$ and $\mathcal{E}(0)$ denotes the initial energy of the
solution.  The point of the function $\chi$ is that it degenerates in
the so-called \textit{trapped region}, the physical space location of
null geodesics which never escapes to infinity or fall in the black
hole in the case of Schwarzschild and Kerr. Starting from an estimate
of the form (\ref{eq:ide}), one can then try to derive pointwise
estimates. The $r^p$-weighted vector
field~\cite{97898143046340032}, mentioned already in
Section~\ref{se:bs}, provides now a systematic and robust way to derive
these pointwise estimates. In particular, one advantage compared to
previous methods is that it avoids completely the use of commutator or
multiplier vector fields having $t$ weights. 

In the recent~\cite{dafermos2022quasilinearwaveequationsasymptotically}
by Holzegel--Dafermos--Rodnianski--Taylor, yet a new and more robust
strategy for the stability of stationary states has been developed. The
main feature is that the analysis is essentially localized in dyadic
time slabs and the fact that two main ingredients for stability are a
``black box'' linear inhomogeneous energy estimate on exactly
stationary metrics, together with a novel physical space top order
identity that need not capture the structure of trapping and is robust
to perturbation. On Schwarschild and slowly rotating Kerr, this
estimate relies essentially on the red-shift effect~\cite{MR3098640}.
For a more precise description of the development of this new strategy
(cf.~the work of Mavrogiannis~\cite{0360530220232295035}),
see the introduction of~\cite{dafermos2022quasilinearwaveequationsasymptotically}. 

\subsection{Development of microlocal analysis and other tools for
stability problems}

The black hole geometries have all in common the presence of trapped
rays. In the high frequency approximation, waves can be localized
arbitrarily close to these specific geodesics  cf.~\cite{MR3397001}, a
clear obstacle to decay, which is then linked to the size and stability
properties of the trapped set\footnote{In particular, the trapped set
of Schwarschild or Kerr is \textit{normally hyperbolic},
cf.~\cite{Hintz_2021} for a recent analysis on general asymptotically
flat spacetimes.}. One way to capture this analytically is via
microlocal and semi-classical analysis, and thus, many recent progress
on the stability properties of black holes relies on frequency
decomposition and pseudo-differential calculus\footnote{We note also
that local integrated decay estimates in this setting are typically
replaced by resolvent estimates. See for
instance~\cite{millet2023optimaldecaysolutionsteukolsky} for a recent
application of microlocal analysis to the Teukolsky equation on Kerr
and the analysis of the resolvent.}.

Interestingly, at least in slowly rotating Kerr, physical space
techniques, exploiting the existence of a Killing tensor field
(replacing the role of Killing vector fields), can still be used to
prove decay, cf.\ Andersson--Blue~\cite{MR3418531}, and in fact, the work
of Giorgi--Klainerman--Szeftel~\cite{giorgi2022waveequationsestimatesnonlinear}, o
ne of the pieces in the recent proof of non-linear stability of slowly
rotating Kerr by Klainerman and
Szeftel~\cite{klainerman2021kerrstabilitysmallangular}, relied on the
Andersson--Blue approach. 

Beyond pseudo-differential calculus, the \textit{b-calculus} of Melrose
has found many applications in general relativity. For instance,
in~\cite{hintzvasy:stab}, Hintz and Vasy have given a novel proof of
the stability of Minkowski space, based on generalized wave
coordinates, the weak null structure of Lindblad--Rodnianski and the
b-calculus, obtaining as a by-product of these techniques
polyhomogeneous expansions of the metric components. The proof of the
stability of slowly Kerr--Sitter spacetimes~\cite{MR3816427} by
Hintz--Vasy also relies heavily on such methods. 

\section{Conclusions: beyond stability in relativity and small data
problems} \label{se:con}

The methods introduced in~\cite{CKM} are now the basic blocks of a
large literature concerning evolution problems in hyperbolic geometric
partial differential equations, including many works which go beyond
the problem of stability of stationary states. For instance, the use of
curvature estimates, of the null structure equations and the
construction and analysis of optical functions are essential in the
other two monumental works of Christodoulou, the formation of trapped
surfaces~\cite{MR2488976} and the formation of shocks~\cite{MR2284927}.

It is also interesting to note that the null condition, which is a key
element of the small data theory of quasi-linear waves, as explained in
Section~\ref{se:nc} is also a key structural property for many
low-regularity well-posedness issues, starting from the null form
estimates of Machedon,~\cite{MR1231427}, to the resolution of the
bounded $L^2$ curvature conjecture by
Klainerman--Rodnianski--Szeftel~\cite{krs:bL2}. In fact, the Bel--Robinson
tensor is also used in~\cite{krs:bL2} to estimate the curvature, and
the null structure\footnote{Many more ingredients are of course used
in~\cite{krs:bL2} to exploit the null condition, such as a quasi-linear
Yang--Mills type formulation of the equations.} which allowed
in~\cite{CKM} to control the error terms via decay is now used to prove
bilinear and trilinear estimates with improved regularity. Finally, we
also remark that ~\cite{krs:bL2} also relies on the maximal foliation. 

The Lindblad--Rodnianski approach has also had numerous developements.
On top of those mentioned in the previous section, which concerns
stability problems of stationary solutions, in~\cite{MR4439921}, Luk--Oh
have addresed the stability of large dispersive solutions to the
Einstein--scalar--field, which are very far from the Minkowski metric,
relying in particular on the understanding of the precise asymptotic
behavior of the metric components in generalized wave coordinates, as
in~\cite{l-17}. In the recent work of Touati~\cite{touati:goaev}
concerning the geometric optics approximation of the Einstein
equations,  the nilpotent structure and the the weak null condition of
the Einstein equations in generalized wave coordinates, as identified
by Lindblad--Rodnianski, also plays a key role in the analysis.

In conclusion, as we hope to have convey in this survey article, the
global non-linear analysis of the Minkowski space is sufficiently
complex that a thorough analysis and a full proof required a deep
understanding of geometric quasilinear wave equations, the Einstein
equations and their fine structure. Thanks to this, the main ideas
introduced in its resolution have and will continue to deeply influence
many of the subsequent, important, developments in General Relativity.

\section*{Declaration of interests}

The authors do not work for, advise, own shares in, or receive funds
from any organization that could benefit from this article, and have
declared no affiliations other than their research organizations.

\back{}

\printbibliography
\refinput{crmeca20240892-reference.tex}

\end{document}
</pre></body></html>