Abstract

This article is a guide to theorems on existence and global dynamics of solutions of the Einstein equations. It draws attention to open questions in the field. The local-in-time Cauchy problem,which is relatively well understood, is surveyed. Global results for solutions with various types of symmetryare discussed. A selection of results from Newtonian theory and special relativity that offer usefulcomparisons is presented. Treatments of global results in the case of small data and results on constructingspacetimes with prescribed singularity structure or late-time asymptotics are given. A conjectural pictureof the asymptotic behaviour of general cosmological solutions of the Einstein equations is built up. Somemiscellaneous topics connected with the main theme are collected in a separate section. 1 Introduction

Systems of partial differential equations are of central importance in physics. Only the simplest of these equations can be solved by explicit formulae. Those that cannot are commonly studied by means of approximations. There is, however, another approach that is complementary. This consists in determining the qualitative behaviour of solutions, without knowing them explicitly.

The first step in doing this is to establish the existence of solutions under appropriate circumstances.

Unfortunately, this is often hard, and obstructs the way to obtaining more interesting information.

When partial differential equations are investigated with a view to applications, proving existence theorems should not become an end in itself. It is important to remember that, from a more general point of view, it is only a starting point.

The basic partial differential equations of general relativity are Einstein's equations. In general, they are coupled to other partial differential equations describing the matter content of spacetime.

The Einstein equations are essentially hyperbolic in nature. In other words, the general properties of solutions are similar to those found for the wave equation. It follows that it is reasonable to try to determine a solution by initial data on a spacelike hypersurface. Thus the Cauchy problem is the natural context for existence theorems for the Einstein equations. The Einstein equations are also nonlinear. This means that there is a big difference between the local and global Cauchy problems. A solution evolving from regular data may develop singularities.

A special feature of the Einstein equations is that they are diffeomorphism invariant. If the equations are written down in an arbitrary coordinate system, then the solutions of these coordinate equations are not uniquely determined by initial data. Applying a diffeomorphism to one solution gives another solution. If this diffeomorphism is the identity on the chosen Cauchy surface up to first order then the data are left unchanged by this transformation. In order to obtain a system for which uniqueness in the Cauchy problem holds in the straightforward sense that it does for the wave equation, some coordinate or gauge fixing must be carried out.

Another special feature of the Einstein equations is that initial data cannot be prescribed freely. They must satisfy constraint equations. To prove the existence of a solution of the Einstein equations, it is first necessary to prove the existence of a solution of the constraints. The usual method of solving the constraints relies on the theory of elliptic equations.

The local existence theory of solutions of the Einstein equations is rather well understood.

Section 2 points out some of the things that are not known. On the other hand, the problem of proving general global existence theorems for the Einstein equations is beyond the reach of the mathematics presently available. To make some progress, it is necessary to concentrate on simplified models. The most common simplifications are to look at solutions with various types of symmetry and solutions for small data. These two approaches are reviewed in Sections 3 and 5 , respectively. A different approach is to prove the existence of solutions with a prescribed singularity structure or late-time asymptotics. This is discussed in Section 6 . Section 9 collects some miscellaneous results that cannot easily be classified. Since insights about the properties of solutions of the Einstein equations can be obtained from the comparison with Newtonian theory and special relativity, relevant results from those areas are presented in Section 4 .

The sections just listed are to some extent catalogues of known results, augmented with some suggestions as to how these could be extended in the future. Sections 7 and 8 complement this by looking ahead to see what the final answer to some interesting general questions might be.

They are necessarily more speculative than the other sections but are rooted in the known results surveyed elsewhere in the article. Section 7 also summarizes various results on cosmological models with accelerated expansion.

The area of research reviewed in the following relies heavily on the theory of differential equations, particularly that of hyperbolic partial differential equations. For the benefit of readers with little background in differential equations, some general references that the author has found to be useful will be listed. A thorough introduction to ordinary differential equations is given in [170] . A lot of intuition for ordinary differential equations can be obtained from [186] . The article by Arnold and Ilyashenko [29] is full of information, in rather compressed form. A classic introductory text on partial differential equations, where hyperbolic equations are well represented, is [197] .

Useful texts on hyperbolic equations, some of which explicitly deal with the Einstein equations, are [335, 202, 269, 239, 327, 198, 137] .

An important aspect of existence theorems in general relativity that one should be aware of is their relation to the cosmic censorship hypothesis. This point of view was introduced in an influential paper by Moncrief and Eardley [248] . An extended discussion of the idea can be found in [114] .

This article is descriptive in nature and equations have been kept to a minimum. A collection of relevant equations together with the background necessary to understand the notation can be found in [301] .

2 Local Existence

In this section basic facts about local existence theorems for the Einstein equations are recalled.

Since the theory is well developed and good accounts exist elsewhere (see for instance [148] ), attention is focussed on remaining open questions known to the author. In particular, the questions of the minimal regularity required to get a well-posed problem and of free boundaries for fluid bodies are discussed.

2.1 The constraints

The unknowns in the constraint equations are the initial data for the Einstein equations. These consist of a three-dimensional manifold
$S$
, a Riemannian metric
${h}_{ab}$
, and a symmetric tensor
${k}_{ab}$
on
$S$
, and initial data for any matter fields present. The equations are:

Here
$R$
is the scalar curvature of the metric
${h}_{ab}$
, and
$\rho $
and
${j}_{a}$
are projections of the energy-momentum tensor. Assuming that the matter fields satisfy the dominant energy condition implies that
$\rho \ge ({j}_{a}{j}^{a}{)}^{1/2}$
. This means that the trivial procedure of making an arbitrary choice of
${h}_{ab}$
and
${k}_{ab}$
and defining
$\rho $
and
${j}_{a}$
by Equations ( 1 ) is of no use for producing physically interesting solutions.

$$\begin{array}{c}\begin{array}{ccc}R-{k}_{ab}{k}^{ab}+({h}^{ab}{k}_{ab}{)}^{2}& =& 16\pi \rho ,\\ {\nabla}^{a}{k}_{ab}-{\nabla}_{b}\left({h}^{ac}{k}_{ac}\right)& =& 8\pi {j}_{b}.\end{array}\end{array}$$ | (1) |

The usual method for solving the Equations ( 1 ) is the conformal method [95] . In this method parts of the data (the so-called free data) are chosen, and the constraints imply four elliptic equations for the remaining parts. The case that has been studied the most is the constant mean curvature (CMC) case, where
$trk={h}^{ab}{k}_{ab}$
is constant. In that case there is an important simplification. Three of the elliptic equations, which form a linear system, decouple from the remaining one. This last equation, which is nonlinear, but scalar, is called the Lichnerowicz equation. The heart of the existence theory for the constraints in the CMC case is the theory of the Lichnerowicz equation.

Solving an elliptic equation is a non-local problem and so boundary conditions or asymptotic conditions are important. For the constraints, the cases most frequently considered in the literature are that where
$S$
is compact (so that no boundary conditions are needed) and that where the free data satisfy some asymptotic flatness conditions. In the CMC case the problem is well understood for both kinds of boundary conditions [76, 110, 187] . (Note that a significant mistake in the literature has been corrected in [244] , section 4.) The other case that has been studied in detail is that of hyperboloidal data [9] . The kind of theorem that is obtained is that sufficiently differentiable free data, in some cases required to satisfy some global restrictions, can be completed in a unique way to a solution of the constraints. It should be noted in passing that in certain cases physically interesting free data may not be “sufficiently differentiable” in the sense it is meant here. One such case is mentioned at the end of Section 2.6 . The usual kinds of differentiability conditions that are required in the study of the constraints involve the free data belonging to suitable Sobolev or Hölder spaces. Sobolev spaces have the advantage that they fit well with the theory of the evolution equations (compare the discussion in Section 2.2 ). The question of the minimal differentiability necessary to apply the conformal method has been studied in [243] , where it was shown that the method works for metrics in the Sobolev space
${H}^{s}$
with
$s>3/2$
. It was also shown that each of these solutions can be approximated by a sequence of smooth solutions.

Usually it is not natural to prescribe the values of solutions of the Einstein equations on a finite boundary. There is, however, one case which naturally occurs in physical problems, namely that of the boundary of a black hole. Existence of solutions of the constraints appropriate for describing black holes has been proved by solving boundary value problems in [127] and [244] .

In the non-CMC case our understanding is much more limited although some results have been obtained in recent years (see [191, 93] and references therein). It is an important open problem to extend these so that an overview is obtained comparable to that available in the CMC case.

Progress on this could also lead to a better understanding of the question of whether a spacetime that admits a compact, or asymptotically flat, Cauchy surface also admits one of constant mean curvature. Up to now there have been only isolated examples that exhibit obstructions to the existence of CMC hypersurfaces [34] . Until very recently it was not known whether there were vacuum spacetimes with a compact Cauchy surface admitting no CMC hypersurfaces. In [117] it was shown using gluing techniques (see below) that spacetimes of this type do exist and this fact restricts the applicability of CMC foliations for defining a preferred time coordinate in cosmological spacetimes. Certain limitations of the conformal method in producing non-CMC initial data sets were exhibited in [193] .

It would be interesting to know whether there is a useful concept of the most general physically reasonable solutions of the constraints representing regular initial configurations. Data of this kind should not themselves contain singularities. Thus it seems reasonable to suppose at least that the metric
${h}_{ab}$
is complete and that the length of
${k}_{ab}$
, as measured using
${h}_{ab}$
, is bounded. Does the existence of solutions of the constraints imply a restriction on the topology of
$S$
or on the asymptotic geometry of the data? This question is largely open, and it seems that information is available only in the compact and asymptotically flat cases. In the case of compact
$S$
, where there is no asymptotic regime, there is known to be no topological restriction. In the asymptotically flat case there is also no topological restriction implied by the constraints beyond that implied by the condition of asymptotic flatness itself [349, 189] . This shows in particular that any manifold that is obtained by deleting a point from a compact manifold admits a solution of the constraints satisfying the minimal conditions demanded above. A starting point for going beyond this could be the study of data that are asymptotically homogeneous. For instance, the Schwarzschild solution contains interesting CMC hypersurfaces that are asymptotic to the metric product of a round 2-sphere with the real line. More general data of this kind could be useful for the study of the dynamics of black hole interiors [292] .

Recently techniques have been developed for gluing together solutions of the constraints (see [117] and references therein). Given two solutions of the constraints it is possible, under very general conditions, to cut a hole in each and connect the resulting pieces by a wormhole to get a new solution of the constraints. Depending on the variant of the method used, the geometry on the original pieces is changed by an arbitrarily small amount, or not at all. This gives a new flexibility in constructing solutions of the constraints with interesting properties.

To sum up, the conformal approach to solving the constraints, which has been the standard one up to now, is well understood in the compact, asymptotically flat and hyperboloidal cases under the constant mean curvature assumption, and only in these cases. For some other approaches see [35, 36, 353] . New techniques have been applied by Corvino [121] to prove the existence of regular solutions of the vacuum constraints on
${\mathbf{R}}^{3}$
that are Schwarzschild outside a compact set. The latter ideas have also flowed into the gluing constructions mentioned above.

2.2 The vacuum evolution equations

The main aspects of the local-in-time existence theory for the Einstein equations can be illustrated by restricting to smooth (i.e. infinitely differentiable) data for the vacuum Einstein equations. The generalizations to less smooth data and matter fields are discussed in Sections 2.3 and 2.5 , respectively. In the vacuum case, the data are
${h}_{ab}$
and
${k}_{ab}$
on a three-dimensional manifold
$S$
, as discussed in Section 2.1 . A solution corresponding to these data is given by a four-dimensional manifold
$M$
, a Lorentz metric
${g}_{\alpha \beta}$
on
$M$
, and an embedding of
$S$
in
$M$
. Here,
${g}_{\alpha \beta}$
is supposed to be a solution of the vacuum Einstein equations, while
${h}_{ab}$
and
${k}_{ab}$
are the induced metric and second fundamental form of the embedding, respectively.

The basic local existence theorem says that, given smooth data for the vacuum Einstein equations, there exists a smooth solution of the equations which gives rise to these data [95] .

Moreover, it can be assumed that the image of
$S$
under the given embedding is a Cauchy surface for the metric
${g}_{\alpha \beta}$
. The latter fact may be expressed loosely, identifying
$S$
with its image, by the statement that
$S$
is a Cauchy surface. A solution of the Einstein equations with given initial data having
$S$
as a Cauchy surface is called a Cauchy development of those data. The existence theorem is local because it says nothing about the size of the solution obtained. A Cauchy development of given data has many open subsets that are also Cauchy developments of that data.

It is intuitively clear what it means for one Cauchy development to be an extension of another.

The extension is called proper if it is strictly larger than the other development. A Cauchy development that has no proper extension is called maximal. The standard global uniqueness theorem for the Einstein equations uses the notion of the maximal development. It is due to Choquet-Bruhat and Geroch [91] . It says that the maximal development of any Cauchy data is unique up to a diffeomorphism that fixes the initial hypersurface. It is also possible to make a statement of Cauchy stability that says that, in an appropriate sense, the solution depends continuously on the initial data. Details on this can be found in [95] .

A somewhat stronger form of the local existence theorem is to say that the solution exists on a uniform time interval in all of space. The meaning of this is not a priori clear, due to the lack of a preferred time coordinate in general relativity. The following is a formulation that is independent of coordinates. Let
$p$
be a point of
$S$
. The temporal extent
$T\left(p\right)$
of a development of data on
$S$
is the supremum of the length of all causal curves in the development passing through
$p$
. In this way, a development defines a function
$T$
on
$S$
. The development can be regarded as a solution that exists on a uniform time interval if
$T$
is bounded below by a strictly positive constant. For compact
$S$
this is a straightforward consequence of Cauchy stability. In the case of asymptotically flat data it is less trivial. In the case of the vacuum Einstein equations it is true, and in fact the function
$T$
grows at least linearly as a function of spatial distance at infinity [110] . It should follow from the results of [211] that the constant of proportionality in the linear lower bound for
$T$
can be chosen to be unity, but this does not seem to have been worked out explicitly.

When proving the above local existence and global uniqueness theorems it is necessary to use some coordinate or gauge conditions. At least no explicitly diffeomorphism-invariant proofs have been found up to now. Introducing these extra elements leads to a system of reduced equations, whose solutions are determined uniquely by initial data in the strict sense, and not just uniquely up to diffeomorphisms. When a solution of the reduced equations has been obtained, it must be checked that it is a solution of the original equations. This means checking that the constraints and gauge conditions propagate. There are many methods for reducing the equations. An overview of the possibilities may be found in [144] . See also [148] .

2.3 Questions of differentiability

Solving the Cauchy problem for a system of partial differential equations involves specifying a set of initial data to be considered, and determining the differentiability properties of solutions.

Thus, two regularity properties are involved – the differentiability of the allowed data, and that of the corresponding solutions. Normally, it is stated that for all data with a given regularity, solutions with a certain type of regularity are obtained. For instance, in Section 2.2 we chose both types of regularity to be “infinitely differentiable”. The correspondence between the regularity of data and that of solutions is not a matter of free choice. It is determined by the equations themselves, and in general the possibilities are severely limited. A similar issue arises in the context of the Einstein constraints, where there is a correspondence between the regularity of free data and full data.

The kinds of regularity properties that can be dealt with in the Cauchy problem depend, of course, on the mathematical techniques available. When solving the Cauchy problem for the Einstein equations, it is necessary to deal at least with nonlinear systems of hyperbolic equations.

(There may be other types of equations involved, but they will be ignored here.) For general nonlinear systems of hyperbolic equations the standard technique is the method of energy estimates.

This method is closely connected with Sobolev spaces, which will now be discussed briefly.

Let
$u$
be a real-valued function on
${\mathbf{R}}^{n}$
. Let

The space of functions for which this quantity is finite is the Sobolev space
${H}^{s}\left({\mathbf{R}}^{n}\right)$
. Here,
$|{D}^{i}u{|}^{2}$
denotes the sum of the squares of all partial derivatives of
$u$
of order
$i$
. Thus, the Sobolev space
${H}^{s}$
is the space of functions, all of whose partial derivatives up to order
$s$
are square integrable.

$$\begin{array}{c}\parallel u{\parallel}_{s}={\left({\sum}_{i=0}^{s}\int \left|{D}^{i}u{|}^{2}\right(x)dx\right)}^{1/2}.\end{array}$$ | (2) |

Similar spaces can be defined for vector valued functions by taking a sum of contributions from the separate components in the integral. It is also possible to define Sobolev spaces on any Riemannian manifold, using covariant derivatives. General information on this can be found in [30] . Consider now a solution
$u$
of the wave equation in Minkowski space. Let
$u\left(t\right)$
be the restriction of this function to a time slice. Then it is easy to compute that, provided
$u$
is smooth and
$u\left(t\right)$
has compact support for each
$t$
, the quantity
$\parallel Du\left(t\right){\parallel}_{s}^{2}+\parallel {\partial}_{t}u\left(t\right){\parallel}_{s}^{2}$
is time independent for each
$s$
. For
$s=0$
this is just the energy of a solution of the wave equation. For a general nonlinear hyperbolic system, the Sobolev norms are no longer time-independent. The constancy in time is replaced by certain inequalities. Due to the similarity to the energy for the wave equation, these are called energy estimates. They constitute the foundation of the theory of hyperbolic equations. It is because of these estimates that Sobolev spaces are natural spaces of initial data in the Cauchy problem for hyperbolic equations. The energy estimates ensure that a solution evolving from data belonging to a given Sobolev space on one spacelike hypersurface will induce data belonging to the same Sobolev space on later spacelike hypersurfaces. In other words, the property of belonging to a Sobolev space is propagated by the equations. Due to the locality properties of hyperbolic equations (existence of a finite domain of dependence), it is useful to introduce the spaces
${H}_{loc}^{s}$
, which are defined by the condition that whenever the domain of integration is restricted to a compact set, the integral defining the space
${H}^{s}$
is finite.

In the end, the solution of the Cauchy problem should be a function that is differentiable enough so that all derivatives that occur in the equation exist in the usual (pointwise) sense. A square integrable function is in general defined only almost everywhere and the derivatives in the above formula must be interpreted as distributional derivatives. For this reason, a connection between Sobolev spaces and functions whose derivatives exist pointwise is required. This is provided by the Sobolev embedding theorem. This says that if a function
$u$
on
${\mathbf{R}}^{n}$
belongs to the Sobolev space
${H}_{loc}^{s}$
and if
$k<s-n/2$
, then there is a
$k$
times continuously differentiable function that agrees with
$u$
except on a set of measure zero.

In the existence and uniqueness theorems stated in Section 2.2 , the assumptions on the initial data for the vacuum Einstein equations can be weakened to say that
${h}_{ab}$
should belong to
${H}_{loc}^{s}$
and
${k}_{ab}$
to
${H}_{loc}^{s-1}$
. Then, provided
$s$
is large enough, a solution is obtained that belongs to
${H}_{loc}^{s}$
.

In fact, its restriction to any spacelike hypersurface also belongs to
${H}_{loc}^{s}$
, a property that is a priori stronger. The details of how large
$s$
must be would be out of place here, since they involve examining the detailed structure of the energy estimates. However, there is a simple rule for computing the required value of
$s$
. The value of
$s$
needed to obtain an existence theorem for the Einstein equations using energy estimates is that for which the Sobolev embedding theorem, applied to spatial slices, just ensures that the metric is continuously differentiable. Thus the requirement is that
$s>n/2+1=5/2$
, since
$n=3$
. It follows that the smallest possible integer
$s$
is three. Strangely enough, the standard methods only give uniqueness up to diffeomorphisms for
$s\ge 4$
. The reason is that in proving the uniqueness theorem a diffeomorphism must be carried out, which need not be smooth. This apparently leads to a loss of one derivative. In [10] local existence and uniqueness for the vacuum Einstein equations was proved using a gauge condition defined by elliptic equations for which this loss does not occur. In that case the gap of one derivative is eliminated. On the other hand, the occurrence of elliptic equations as part of the reduced Einstein equations with this gauge makes the result intrinsically global, and it is not clear whether it can be localized in space. Another interesting aspect of the main theorem of [10] is that it includes a continuation criterion for solutions. There exists a definition of Sobolev spaces for an arbitrary real number
$s$
, and hyperbolic equations can also be solved in the spaces with
$s$
not an integer [334] .

Presumably these techniques could be applied to prove local existence for the Einstein equations with
$s$
any real number greater than
$5/2$
. In any case, the condition for local existence has been weakened to
$s>2$
using other techniques, as discussed in Section 2.4 .

Consider now
${C}^{\infty}$
initial data. Corresponding to these data there is a development of class
${H}^{s}$
for each
$s$
. It could conceivably be the case that the size of these developments shrinks with increasing
$s$
. In that case, their intersection might contain no open neighbourhood of the initial hypersurface, and no smooth development would be obtained. Fortunately, it is known that the
${H}^{s}$
developments cannot shrink with increasing
$s$
[87] , and so the existence of a
${C}^{\infty}$
solution is obtained for
${C}^{\infty}$
data. It appears that the
${H}^{s}$
spaces with
$s$
sufficiently large are the only spaces containing the space of smooth functions for which it has been proved that the Einstein equations are locally solvable.

What is the motivation for considering regularity conditions other than the apparently very natural
${C}^{\infty}$
condition? One motivation concerns matter fields and will be discussed in Section 2.5 .

Another is the idea that assuming the existence of many derivatives that have no direct physical significance seems like an admission that the problem has not been fully understood. A further reason for considering low regularity solutions is connected to the possibility of extending a local existence result to a global one. If the proof of a local existence theorem is examined closely it is generally possible to give a continuation criterion. This is a statement that if a solution on a finite time interval is such that a certain quantity constructed from the solution is bounded on that interval, then the solution can be extended to a longer time interval. (In applying this to the Einstein equations we need to worry about introducing an appropriate time coordinate.) If it can be shown that the relevant quantity is bounded on any finite time interval where a solution exists, then global existence follows. It suffices to consider the maximal interval on which a solution is defined, and obtain a contradiction if that interval is finite. This description is a little vague, but contains the essence of a type of argument that is often used in global existence proofs. The problem in putting it into practice is that often the quantity whose boundedness has to be checked contains many derivatives, and is therefore difficult to control. If the continuation criterion can be improved by reducing the number of derivatives required, then this can be a significant step toward a global result. Reducing the number of derivatives in the continuation criterion is closely related to reducing the number of derivatives of the data required for a local existence proof.

A striking example is provided by the work of Klainerman and Machedon [210] on the Yang–Mills equations in Minkowski space. Global existence in this case was first proved by Eardley and Moncrief [135] , assuming initial data of sufficiently high differentiability. Klainerman and Machedon gave a new proof of this, which, though technically complicated, is based on a conceptually simple idea. They prove a local existence theorem for data of finite energy. Since energy is conserved this immediately proves global existence. In this case finite energy corresponds to the Sobolev space
${H}^{1}$
for the gauge potential. Of course, a result of this kind cannot be expected for the Einstein equations, since spacetime singularities do sometimes develop from regular initial data. However, some weaker analogue of the result could exist.

2.4 New techniques for rough solutions

Recently, new mathematical techniques have been developed to lower the threshold of differentiability required to obtain local existence for quasilinear wave equations in general and the Einstein equations in particular. Some aspects of this development will now be discussed following [209, 214] .

A central aspect is that of Strichartz inequalities. These allow one to go beyond the theory based on
${L}^{2}$
spaces and use Sobolev spaces based on the Lebesgue
${L}^{p}$
spaces for
$p\ne 2$
. The classical approach to deriving Strichartz estimates is based on the Fourier transform and applies to flat space. The new ideas allow the use of the Fourier transform to be limited to that of Littlewood–Paley theory and facilitate generalizations to curved space.

The idea of Littlewood–Paley theory is as follows (see [1] for a good exposition of this). Suppose that we want to describe the regularity of a function (or, more generally, a tempered distribution)
$u$
on
${\mathbf{R}}^{n}$
. Differentiability properties of
$u$
correspond, roughly speaking, to fall-off properties of its Fourier transform
$\hat{u}$
. This is because the Fourier transform converts differentiation into multiplication. The Fourier transform is decomposed as
$\hat{u}=\sum {\phi}_{i}\hat{u}$
, where
${\phi}_{i}$
is a dyadic partition of unity. The statement that it is dyadic means that all the
${\phi}_{i}$
except one are obtained from each other by scaling the argument by a factor which is a power of two. Transforming back we get the decomposition
$u=\sum {u}_{i}$
, where
${u}_{i}$
is the inverse Fourier transform of
${\phi}_{i}\hat{u}$
. The component
${u}_{i}$
of
$u$
contains only frequencies of the order
${2}^{i}$
. In studying rough solutions of the Einstein equations, the Littlewood–Paley decomposition is applied to the metric itself. The high frequencies are discarded to obtain a smoothed metric which plays an important role in the arguments.

Another important element of the proofs is to rescale the solution by a factor depending on the cut-off
$\lambda $
applied in the Littlewood–Paley decomposition. Proving the desired estimates then comes down to proving the existence of the rescaled solutions on a time interval depending on
$\lambda $
in a particular way. The rescaled data are small in some sense and so a connection is established to the question of long-time existence of solutions of the Einstein equations for small initial data. In this way, techniques from the work of Christodoulou and Klainerman on the stability of Minkowski space (see Section 5.2 ) are brought in.

What is finally proved? In general, there is a close connection between proving local existence for data in a certain space and showing that the time of existence of smooth solutions depends only on the norm of the data in the given space. Klainerman and Rodnianski [214] demonstrate that the time of existence of solutions of the reduced Einstein equations in harmonic coordinates depends only on the
${H}^{2+\epsilon}$
norm of the initial data for any
$\epsilon >0$
. Combining this with the results of [243] gives an existence result in the same space. It is of interest to try to push the existence theorem to the limiting case
$\epsilon =0$
or even to the slightly weaker assumption on the data that the curvature is square integrable. This
${L}^{2}$
curvature conjecture (local existence in this setting) is extremely difficult but interesting progress has been made by Klainerman and Rodnianski in [215] where the structure of null hypersurfaces was analysed under very weak hypotheses. For this purpose the authors developed an invariant form of Littlewood–Paley theory [216] . This uses the asymptotics of solutions of the heat equation on a manifold and is coordinate-independent.

The techniques discussed in this section, which have been stimulated by the desire to understand the Einstein equations, are also helpful in understanding other nonlinear wave equations. Thus, this is an example where information can flow from general relativity to the theory of partial differential equations.

It may be that the technique of using parabolic equations as a tool to better understand hyperbolic equations can be carried much further. In [333] Tao presents ideas how the harmonic map heat flow could be used to define a high quality gauge for the study of wave maps.

2.5 Matter fields

Analogues of the results for the vacuum Einstein equations given in Section 2.2 are known for the Einstein equations coupled to many types of matter model. These include perfect fluids, elasticity theory [43] , kinetic theory, scalar fields, Maxwell fields, Yang–Mills fields, and combinations of these.

An important restriction is that the general results for perfect fluids and elasticity apply only to situations where the energy density is uniformly bounded away from zero on the region of interest.

In particular, they do not apply to cases representing material bodies surrounded by vacuum. In cases where the energy density, while everywhere positive, tends to zero at infinity, a local solution is known to exist, but it is not clear whether a local existence theorem can be obtained that is uniform in time. In cases where the fluid has a sharp boundary, ignoring the boundary leads to solutions of the Einstein–Euler equations with low differentiability (cf. Section 2.3 ), while taking it into account explicitly leads to a free boundary problem. This will be discussed in more detail in Section 2.6 . In the case of kinetic or field theoretic matter models it makes no difference whether the energy density vanishes somewhere or not.

2.6 Free boundary problems

In applying general relativity one would like to have solutions of the Einstein–matter equations modelling material bodies. As will be discussed in Section 3.1 there are solutions available for describing equilibrium situations. However, dynamical situations require solving a free boundary problem if the body is to be made of fluid or an elastic solid. We will now discuss the few results which are known on this subject. For a spherically symmetric self-gravitating fluid body in general relativity, a local-in-time existence theorem was proved in [206] . This concerned the case in which the density of the fluid at the boundary is non-zero. In [285] a local existence theorem was proved for certain equations of state with vanishing boundary density. These solutions need not have any symmetry but they are very special in other ways. In particular, they do not include small perturbations of the stationary solutions discussed in Section 3.1 . There is no general result on this problem up to now.

Remarkably, the free boundary problem for a fluid body is also poorly understood in classical physics. There is a result for a viscous fluid [319] , but in the case of a perfect fluid the problem was wide open until recently. A major step forward was taken by Wu [352] , who obtained a result for a fluid that is incompressible and irrotational. There is a good physical reason why local existence for a fluid with a free boundary might fail. This is the Rayleigh–Taylor instability which involves perturbations of fluid interfaces that grow with unbounded exponential rates (cf. the discussion in [41] ). It turns out that in the case considered by Wu this instability does not cause problems, and there is no reason to expect that a self-gravitating compressible fluid with rotation in general relativity with a free boundary cannot also be described by a well-posed free boundary value problem. For the generalization of the problem considered by Wu to the case of a fluid with rotation, Christodoulou and Lindblad [109] have obtained estimates that look as if they should be enough to obtain an existence theorem. Strangely, it proved very difficult to complete the argument.

This point deserves some further comment. In many problems the heart of an existence proof is obtaining suitable estimates. Then more or less standard approximation techniques can be used to obtain the desired conclusion (for a discussion of this see [148] , Section 3.1). In the problem studied in [109] it is an appropriate approximation method that is missing. More recently Lindblad was able to obtain an existence result using a different approach involving the Nash–Moser theorem and a detailed analysis of the linearized system about a given solution. He treated the incompressible case in [231] while in the case of a compressible fluid with non-vanishing boundary density the linearized analysis has been carried out [230] .

One of the problems in tackling the initial value problem for a dynamical fluid body is that the boundary is moving. It would be very convenient to use Lagrangian coordinates, since in those coordinates the boundary is fixed. Unfortunately, it is not at all obvious that the Euler equations in Lagrangian coordinates have a well-posed initial value problem, even in the absence of a boundary.

It was, however, recently shown by Friedrich [145] that it is possible to treat the Cauchy problem for fluids in general relativity in Lagrangian coordinates.

In the case of a fluid with non-vanishing boundary density it is not only the evolution equations that cause problems. It is already difficult to construct suitable solutions of the constraints. A theorem on this has recently been obtained by Dain and Nagy [128] . There remains an undesirable technical restriction, but the theorem nevertheless provides a very general class of physically interesting initial data for a self-gravitating fluid body in general relativity.

3 Global Symmetric Solutions

An obvious procedure to obtain special cases of the general global existence problem for the Einstein equations that are amenable to attack is to make symmetry assumptions. In this section, we discuss the results obtained for various symmetry classes defined by different choices of number and character of Killing vectors.

The existence of a timelike Killing vector is a strong restriction on solutions of the Einstein equations. It corresponds physically to an equilibrium situation. This case is considered in Section 3.1 . A null Killing vector also represents a strong restriction. It seems to fit badly with the Cauchy problem and is thus far from the discussion in this paper. This case will not be discussed further here. There remains the case where all Killing vectors are spacelike. The classes of solutions can be distinguished by the dimension of the group orbits. The case where the orbits are three-dimensional is that of spatially homogeneous spacetimes, discussed in Section 3.2 . If the group orbits are two-dimensional there are important classes of solutions where the group is threeor two-dimensional. The case of a three-dimensional group includes asymptotically flat spherically symmetric solutions, discussed in Sections 3.3 and 3.4 . It also includes spatially compact spacetimes with spherical, plane and hyperbolic symmetry, for which results are included in Section 3.6 . That section also covers spatially compact solutions with a two-dimensional isometry group. Most results concern the situation where the Killing vectors have no zeroes. This means that when the spacetime is quotiented by the symmetry group in order to get an effective lower-dimensional system, no singularities occur in the reduced equations of motion. This is a big advantage for analytical work. The assumption of exactly two spacelike Killing vectors is not compatible with asymptotic flatness. It does allow asymptotic flatness in all directions except one and this is the subject of Section 3.5 .

A single Killing vector without fixed points gives rise to an interesting class of cosmological models. Because of the perceived closeness to solutions without any Killing vectors, the results on these spacetimes have been placed in Section 5.4 . Spacetimes with one fixed-point free Killing vector which are partly asymptotically flat are mentioned briefly in Section 3.5 . Axially symmetric asymptotically flat spacetimes have one spacelike Killing vector with zeroes on the axis. It seems that no analytical advantage over the general case has been found for these spacetimes and so there are no results to report.

3.1 Stationary solutions

Many of the results on global solutions of the Einstein equations involve considering classes of spacetimes with Killing vectors. A particularly simple case is that of a timelike Killing vector, i.e. the case of stationary spacetimes. In the vacuum case there are very few solutions satisfying physically reasonable boundary conditions. This is related to no hair theorems for black holes and lies outside the scope of this review. More information on the topic can be found in the book of Heusler [181] and in his Living Review [182] (see also [55] where the stability of the Kerr metric is discussed). Anderson [3, 4] has proved uniqueness theorems for stationary vacuum spacetimes under very weak assumptions. The case of phenomenological matter models has been reviewed in [298] .

The account given there will be updated in the following.

The area of stationary solutions of the Einstein equations coupled to field theoretic matter models has been active in recent years as a consequence of the discovery by Bartnik and McKinnon [37] of a discrete family of regular, static, spherically symmetric solutions of the Einstein–Yang–Mills equations with gauge group
$SU\left(2\right)$
. The equations to be solved are ordinary differential equations, and in [37] they were solved numerically by a shooting method. The first existence proof for a solution of this kind is due to Smoller, Wasserman, Yau, and McLeod [324] and involves an arduous qualitative analysis of the differential equations. The work on the Bartnik–McKinnon solutions, including the existence theorems, has been extended in many directions. Recently, static solutions of the Einstein–Yang–Mills equations that are not spherically symmetric were discovered numerically [217] . It is a challenge to prove the existence of solutions of this kind. Now the ordinary differential equations of the previously known case are replaced by elliptic equations. Moreover, the solutions appear to still be discrete, so that a simple perturbation argument starting from the spherical case does not seem feasible. In another development, it was shown that a linearized analysis indicates the existence of stationary non-static solutions [71] . It would be desirable to study the question of linearization stability in this case, which, if the answer were favourable, would give an existence proof for solutions of this kind. It has, however, been argued that solutions of this kind should not exist [343] .

Now we return to phenomenological matter models, starting with the case of spherically symmetric static solutions. Basic existence theorems for this case have been proved for perfect fluids [307] , collisionless matter [279, 271] , and elastic bodies [265] . All these theorems demonstrate the existence of solutions that are everywhere smooth and exist globally as functions of area radius for a general class of constitutive relations. The physically significant question of the finiteness of the mass of these configurations was only answered in these papers under restricted circumstances.

For instance, in the case of perfect fluids and collisionless matter, solutions were constructed by perturbing about the Newtonian case. Solutions for an elastic body were obtained by perturbing about the case of isotropic pressure, which is equivalent to a fluid. Further progress on the question of the finiteness of the mass of the solutions was made in the case of a fluid by Makino [240] , who gave a rather general criterion on the equation of state ensuring the finiteness of the radius.

Further information on this issue was obtained in [174] using dynamical systems methods. Makino's criterion was generalized to kinetic theory in [281] . This resulted in existence proofs for various models that have been considered in galactic dynamics and which had previously been constructed numerically (cf. [58, 320] for an account of these models in the non-relativistic and relativistic cases, respectively). In the non-relativistic case dynamical systems methods were applied to the case of collisionless matter in [173] . Most of the work quoted up to now refers to solutions where the support of the density is a ball. For matter with anisotropic pressure the support may also be a shell, i.e. the region bounded by two concentric spheres. The existence of static shells in the case of the Einstein–Vlasov equations was proved in [274] .

In the case of self-gravitating Newtonian spherically symmetric configurations of collisionless matter, it can be proved that the phase space density of particles depends only on the energy of the particle and the modulus of its angular momentum [38] . This is known as Jeans' theorem. It was already shown in [271] that the naive generalization of this to the general relativistic case does not hold if a black hole is present. Recently, counterexamples to the generalization of Jeans' theorem to the relativistic case, which are not dependent on a black hole, were constructed by Schaeffer [318] .

It remains to be seen whether there might be a natural modification of the formulation that would lead to a true statement.

For a perfect fluid there are results stating that a static solution is necessarily spherically symmetric [234] . They still require a restriction on the equation of state, which it would be desirable to remove. A similar result is not to be expected in the case of other matter models, although as yet no examples of non-spherical static solutions are available. In the Newtonian case examples have been constructed by Rein [275] . (In that case static solutions are defined to be those in which the particle current vanishes.) For a fluid there is an existence theorem for solutions that are stationary but not static (models for rotating stars) [172] . At present there are no corresponding theorems for collisionless matter or elastic bodies. In [275] , stationary, non-static configurations of collisionless matter were constructed in the Newtonian case.

Two obvious characteristics of a spherically symmetric static solution of the Einstein–Euler equations that has a non-zero density only in a bounded spatial region are its radius
$R$
and its total mass
$M$
. For a given equation of state there is a one-parameter family of solutions. These trace out a curve in the
$(M,R)$
plane. In the physics literature, pictures of this curve indicate that it spirals in on a certain point in the limit of large density. The occurrence of such a spiral and its precise asymptotic form have been proved rigorously by Makino [241] for a particular choice of equation of state. An approach to these spirals which leads to a better conceptual understanding can be found in [174] .

The existence of cylindrically symmetric static solutions of the Einstein–Euler system has been proved in [56] . For some remarks on the question of stability of spherically symmetric solutions see Section 4.1 .

3.2 Spatially homogeneous solutions

A solution of the Einstein equations is called spatially homogeneous if there exists a group of symmetries with three-dimensional spacelike orbits. In this case there are at least three linearly independent spacelike Killing vector fields. For most matter models the field equations reduce to ordinary differential equations. (Kinetic matter leads to an integro-differential equation.) The most important results in this area have been reviewed in the book [344] edited by Wainwright and Ellis (see, in particular, Part Two of the book). There remain a host of interesting and accessible open questions. The spatially homogeneous solutions have the advantage that it is not necessary to stop at just existence theorems; information on the global qualitative behaviour of solutions can also be obtained.

An important question that was open for a long time concerns the mixmaster model, as discussed in [296] . This is a class of spatially homogeneous solutions of the vacuum Einstein equations, which are invariant under the group
$SU\left(2\right)$
. A special subclass of these
$SU\left(2\right)$
-invariant solutions, the (parameter-dependent) Taub–NUT solution, is known explicitly in terms of elementary functions. The Taub–NUT solution has a simple initial singularity which is in fact a Cauchy horizon. All other vacuum solutions admitting a transitive action of
$SU\left(2\right)$
on spacelike hypersurfaces (Bianchi type IX solutions) will be called generic in the present discussion. These generic Bianchi IX solutions (which might be said to constitute the mixmaster solution proper) have been believed for a long time to have singularities that are oscillatory in nature where some curvature invariant blows up. This belief was based on a combination of heuristic considerations and numerical calculations. Although these together do make a persuasive case for the accepted picture, until recently there were no mathematical proofs of these features of the mixmaster model available. This has now changed. First, a proof of curvature blow-up and oscillatory behaviour for a simpler model (a solution of the Einstein–Maxwell equations) which shares many qualitative features with the mixmaster model, was obtained by Weaver [347] . In the much more difficult case of the mixmaster model itself, corresponding results were obtained by Ringström [312] . Later he extended this in several directions in [311] . In that paper more detailed information was obtained concerning the asymptotics and an attractor for the evolution was identified. It was shown that generic solutions of Bianchi type IX with a perfect fluid whose equation of state is
$p=(\gamma -1)\rho $
with
$1\le \gamma <2$
are approximated near the singularity by vacuum solutions. The case of a stiff fluid (
$\gamma =2$
) which has a different asymptotic behaviour was analysed completely for all models of Bianchi class A, a class which includes Bianchi type IX. Ringström's analysis of the mixmaster model is potentially of great significance for the mathematical understanding of singularities of the Einstein equations in general. Thus, its significance goes far beyond the spatially homogeneous case. According to extensive investigations of Belinskii, Khalatnikov, and Lifshitz (see [226, 47, 48] and references therein), the mixmaster model should provide an approximate description for the general behaviour of solutions of the Einstein equations near singularities. This should apply to many matter models as well as to the vacuum equations. The work of Belinskii, Khalatnikov, and Lifshitz (BKL) is hard to understand and it is particularly difficult to find a precise mathematical formulation of their conclusions. This has caused many people to remain sceptical about the validity of the BKL picture. Nevertheless, it seems that nothing has ever been found to indicate any significant flaws in the final version. As long as the mixmaster model itself was not understood, this represented a fundamental obstacle to progress on understanding the BKL picture mathematically. The removal of this barrier opens up an avenue to progress on this issue. The BKL picture is discussed in more detail in Section 8 .

Some recent and qualitatively new results concerning the asymptotic behaviour of spatially homogeneous solutions of the Einstein–matter equations, both close to the initial singularity and in a phase of unlimited expansion (and with various matter models), can be found in [308, 309, 302, 345, 259, 184, 180] .

These show in particular that the dynamics can depend sensitively on the form of matter chosen.

(Note that these results are consistent with the BKL picture.) The dynamics of indefinitely expanding cosmological models is discussed further in Section 7 .

3.3 Spherically symmetric asymptotically flat solutions

The most extensive results on global inhomogeneous solutions of the Einstein equations obtained up to now concern spherically symmetric solutions of the Einstein equations coupled to a massless scalar field with asymptotically flat initial data. In a series of papers, Christodoulou [97, 96, 99, 98, 100, 101, 102, 106] has proved a variety of deep results on the global structure of these solutions. Particularly notable are his proofs that naked singularities can develop from regular initial data [102] and that this phenomenon is unstable with respect to perturbations of the data [106] . In related work, Christodoulou [103, 104, 105] has studied global spherically symmetric solutions of the Einstein equations coupled to a fluid with a special equation of state (the so-called two-phase model).

Generalization of the results of [97] to the case of a nonlinear scalar field and to the Maxwell–Higgs system have been given by Chae [81, 82] .

The rigorous investigation of the spherically symmetric collapse of collisionless matter in general relativity was initiated by Rein and the author [278] , who showed that the evolution of small initial data leads to geodesically complete spacetimes where the density and curvature fall off at large times. Later, it was shown [282] that independent of the size of the initial data the first singularity, if there is one at all, must occur at the centre of symmetry. This result uses a time coordinate of Schwarzschild type; an analogous result for a maximal time coordinate was proved in [297] . The generalization of these results to the case of charged matter has been investigated in [261] and [260] .

The question of what happens in the collapse of uncharged collisionless matter for general large initial data could not yet be answered by analytical techniques. In [283] , numerical methods were applied to try to make some progress in this direction. The results are discussed below.

Some of the results of Christodoulou have been extended to much more general spacetimes by Dafermos [124] . In this work there are two basic assumptions. The first is the existence of at least one trapped surface in the spacetime under consideration. The second is that the matter content is well behaved in a certain sense which means intuitively that it does not form singularities outside black hole regions. Under these circumstances conclusions can be drawn on the global structure of the spacetime. It contains a black hole with a complete null infinity. It has been shown that collisionless matter has the desired property [125] . It also holds for certain nonlinear scalar fields and this has led to valuable insights in the discussion of the formation of naked singularities in a class of models motivated by string theory [179, 123] .

Despite the range and diversity of the results obtained by Christodoulou on the spherical collapse of a scalar field, they do not encompass some of the most interesting phenomena that have been observed numerically. These are related to the issue of critical collapse. For sufficiently small data the field disperses. For sufficiently large data a black hole is formed. The question is what happens in between. This can be investigated by examining a one-parameter family of initial data interpolating between the two cases. It was found by Choptuik [86] that there is a critical value of the parameter below which dispersion takes place and above which a black hole is formed, and that the mass of the black hole approaches zero as the critical parameter value is approached. This gave rise to a large literature in which the spherical collapse of different kinds of matter was computed numerically and various qualitative features were determined. For reviews of this see [162, 163] . In the calculations of [283] for collisionless matter, it was found that in the situations considered the black hole mass tended to a strictly positive limit as the critical parameter was approached from above. These results were confirmed and extended by Olabarrieta and Choptuik [264] . There are no rigorous mathematical results available on the issue of a mass gap for either a scalar field or collisionless matter, and it is an outstanding challenge for mathematical relativists to change this situation.

Another aspect of Choptuik's results is the occurrence of a discretely self-similar solution. It would seem hard to prove the existence of a solution of this kind analytically. For other types of matter, such as a perfect fluid with linear equation of state, the critical solution is continuously self-similar and this looks more tractable. The problem reduces to solving a system of singular ordinary differential equations subject to certain boundary conditions. This problem was solved in [102] for the case where the matter model is given by a massless scalar field, but the solutions produced there, which are continuously self-similar, cannot include the Choptuik critical solution.

Bizoń and Wasserman [62] studied the corresponding problem for the Einstein equations coupled to a wave map with target
$SU\left(2\right)$
. They proved the existence of continuously self-similar solutions including one which, according the results of numerical calculations, appears to play the role of critical solution in collapse. Another case where the question of the existence of the critical solution seems to be a problem that could possibly be solved in the near future is that of a perfect fluid. A good starting point for this is the work of Goliath, Nilsson, and Uggla [158, 159] . These authors gave a formulation of the problem in terms of dynamical systems and were able to determine certain qualitative features of the solutions (see also [77, 78] ).

A possible strategy for learning more about critical collapse, pursued by Bizoń and collaborators, is to study model problems in flat space that exhibit features similar to those observed numerically in the case of the Einstein equations. Until now, only models showing continuous self-similarity have been found. These include wave maps in various dimensions and the Yang–Mills equations in spacetimes of dimension greater than four. As mentioned in Section 2.3 , it is known that in four dimensions there exist global smooth solutions of the Yang–Mills equations corresponding to rather general initial data [135, 210] . In dimensions greater than five it is known that there exist solutions that develop singularities in finite time. This follows from the existence of continuously self-similar solutions [61] . Numerical evidence indicates that this type of blow-up is stable, i.e. occurs for an open set of initial data. The numerical work also indicates that there is a critical self-similar solution separating this kind of blow-up from dispersion. The spacetime dimension five is critical for Yang–Mills theory. Apparently singularities form, but in a different way from what happens in dimension six. There is as yet no rigorous proof of blow-up in five dimensions.

The various features of Yang–Mills theory just mentioned are mirrored in two dimensions less by wave maps with values in spheres [60] . In four dimensions, blow-up is known while in three dimensions there appears (numerically) to be a kind of blow-up similar to that found for Yang–Mills in dimension five. There is no rigorous proof of blow-up. What is seen numerically is that the collapse takes place by scaling within a one-parameter family of static solutions. The case of wave maps is the most favourable known model problem for proving theorems about critical phenomena associated to singularity formation. The existence of a solution having the properties expected of the critical solution for wave maps in four dimensions has been proved in [59] . Some rigorous support for the numerical findings in three dimensions has been given by work of Struwe [328] . He showed, among other things, that if there is blow-up in finite time it must take place in a way resembling that observed in the numerical calculations.

Self-similar solutions are characteristic of what is called Type II critical collapse. In Type I collapse an analogous role is played by static solutions and quite a bit is known about the existence of these. For instance, in the case of the Einstein–Yang–Mills equations, it is one of the Bartnik–McKinnon solutions mentioned in Section 3.1 which does this. In the case of collisionless matter the results of [264] show that at least in some cases critical collapse is mediated by a static solution in the form of a shell. There are existence results for shells of this kind [274] although no connection has yet been made between those shells whose existence has been proved and those which have been observed numerically in critical collapse calculations. Note that Martín-García and Gundlach [242] have presented a (partially numerical) construction of self-similar solutions of the Einstein–Vlasov system.

3.4 Weak null singularities and Price's law

The results of this section concern spherically symmetric solutions but in order to explain their significance they need to be presented in context. A non-rotating uncharged black hole is represented by the Schwarzschild solution, which contains a singularity. At this singularity the Kretschmann scalar
${R}_{\alpha \beta \gamma \delta}{R}^{\alpha \beta \gamma \delta}$
blows up uniformly and this represents an obstruction to extending the spacetime through the singularity, at least in a
${C}^{2}$
manner.

A rotating uncharged black hole is represented by the Kerr solution in which the Schwarzschild singularity is replaced by a Cauchy horizon. This horizon marks a pathology of the global causal structure of the solution but locally the geometry can be extended smoothly through it. A similar situation is found in the non-rotating charged black hole which is represented by the Reissner–Nordström solution. These facts are worrying since they suggest that black holes may generally lead to causal pathologies. The rotating case is the more physically interesting one, but the charged case is a valuable model problem for the rotating case. Spherical symmetry leads to immense technical simplifications and so only that case will be discussed here. It is the only one where theorems on global existence and qualitative behaviour relevant to this problem are available.

It was early suggested that the Cauchy horizon of the Reissner–Nordström solution should be unstable and that a generic perturbation of the initial data would lead to its being replaced by a Schwarzschild-like singularity. This scenario turned out to be oversimplified. An alternative was suggested by Poisson and Israel [267] . In their picture a generic perturbation of the Reissner–Nordström data leads to the Cauchy horizon being replaced by what they call a weak null singularity. At this singularity the curvature blows up, but the metric can be extended through the singularity in a way which is continuous and non-degenerate. In this situation it is possible to make sense of the causal character of the singularity which turns out to be null. Furthermore, an important invariant, the Hawking mass, blows up at the singularity, a phenomenon known as mass inflation.

All these conclusions were based on heuristic arguments which were later backed up by numerical results [185] .

A mathematical understanding of these effects came with the work of Dafermos [122] . He showed how, starting from a characteristic initial value problem with data given on two null hypersurfaces, one of which is the event horizon, it is possible to prove that a weak null singularity forms and that there is mass inflation. He uses a model with an uncharged scalar field and a static charge and works entirely inside the black hole region.

Ideally one would wish to start with regular data on a standard Cauchy surface and control both formation of the black hole and the evolution in its interior. This requires using some kind of charged matter, e.g., a charged scalar field. This is what was done numerically in [185] . Analytically it remains out of reach at the moment.

In the original heuristic arguments it is important to make statements about the behaviour of the solution outside the black hole and what behaviour on the horizon results. Here there are classical heuristic results of Price [268] for a scalar field on a black hole background. He states that the scalar field falls off in a certain way along the horizon. Let us call this Price's law.

Now a form of Price's law and its analogue for the coupled spherically symmetric Einstein–scalar field system have been proved by Dafermos and Rodnianski [126] . Thus we have come a long way towards an understanding of the problem discussed here. This has required the development of new mathematical techniques and these may one day turn out to be of importance in understanding the nonlinear stability of black holes.

3.5 Cylindrically symmetric solutions

Solutions of the Einstein equations with cylindrical symmetry that are asymptotically flat in all directions allowed by the symmetry represent an interesting variation on asymptotic flatness.

There are two Killing vectors, one translational (without fixed points) and one rotational (with fixed points on the axis). Since black holes are apparently incompatible with this symmetry, one may hope to prove geodesic completeness of solutions under appropriate assumptions. (It would be interesting to have a theorem making the statement about black holes precise.) A proof of geodesic completeness has been achieved for the Einstein vacuum equations and for the source-free Einstein–Maxwell equations in [51] , building on global existence theorems for wave maps [112, 111] .

For a quite different point of view on this question involving integrable systems see [351] . A recent paper of Hauser and Ernst [171] also appears to be related to this question. However, due to the great length of this text and its reliance on many concepts unfamiliar to this author, no further useful comments on the subject can be made here.

Solutions of the Einstein–Vlasov system with cylindrical symmetry have been studied by Fjällborg [140] . He shows global existence provided certain conditions are satisfied near the axis.

Cylindrical symmetry can be generalized by abandoning the rotational Killing vector while maintaining the translational one. This sitation does not seem to have been studied in the literature. It may be that results on solutions with approximate cylindrical symmetry may be obtained using the work of Krieger [218] on wave maps.

3.6 Spatially compact solutions

In the context of spatially compact spacetimes it is first necessary to ask what kind of global statements are to be expected. In a situation where the model expands indefinitely it is natural to pose the question whether the spacetime is causally geodesically complete towards the future.

In a situation where the model develops a singularity either in the past or in the future one can ask what the qualitative nature of the singularity is. It is very difficult to prove results of this kind. As a first step one may prove a global existence theorem in a well-chosen time coordinate.

In other words, a time coordinate is chosen that is geometrically defined and that, under ideal circumstances, will take all values in a certain interval
$({t}_{-},{t}_{+})$
. The aim is then to show that, in the maximal Cauchy development of data belonging to a certain class, a time coordinate of the given type exists and exhausts the expected interval. The first result of this kind for inhomogeneous spacetimes was proved by Moncrief in [246] . This result concerned Gowdy spacetimes. These are vacuum spacetimes with a two-dimensional Abelian group of isometries acting on compact orbits.

The area of the orbits defines a natural time coordinate (areal time coordinate). Moncrief showed that in the maximal Cauchy development of data given on a hypersurface of constant time, this time coordinate takes on the maximal possible range, namely
$(0,\infty ).$
This result was extended to more general vacuum spacetimes with two Killing vectors in [50] . Andréasson [16] extended it in another direction to the case of collisionless matter in a spacetime with Gowdy symmetry. This development was completed in [20] where general cosmological solutions of the Einstein–Vlasov system with two commuting spacelike Killing vectors were treated. Corresponding results for spacetimes with hyperbolic symmetry were obtained in [19] .

In all of these cases other than Gowdy the areal time coordinate was proved to cover the maximal globally hyperbolic development, but the range of the coordinate was only shown to be
$({R}_{0},\infty )$
for an undetermined constant
${R}_{0}>0$
. It was not known whether
${R}_{0}$
was necessarily zero except in the Gowdy case. This issue was settled in [195] for the vacuum case with two commuting Killing vectors and this was extended to include Vlasov matter in [348] . It turns out that in vacuum
${R}_{0}=0$
apart from the exceptional case of the flat Kasner solution and an unconventional choice of the two Killing vectors. With Vlasov matter and a distribution function which does not vanish identically,
${R}_{0}=0$
without exception. The corresponding result in cosmological models with spherical symmetry was proved in [336] where the case of a negative cosmological constant was also included. For solutions of the Einstein–Vlasov system with hyperbolic symmetry the question is still open, although the homogeneous case was treated in [336] .

Another attractive time coordinate is constant mean curvature (CMC) time. For a general discussion of this see [292] . A global existence theorem in this time for spacetimes with two Killing vectors and certain matter models (collisionless matter, wave maps) was proved in [295] . That the choice of matter model is important for this result was demonstrated by a global non-existence result for dust in [294] . As shown in [194] , this leads to examples of spacetimes that are not covered by a CMC slicing. Results on global existence of CMC foliations have also been obtained for spherical and hyperbolic symmetry [289, 72] .

A drawback of the results on the existence of CMC foliations just cited is that they require as a hypothesis the existence of one CMC Cauchy surface in the given spacetime. More recently, this restriction has been removed in certain cases by Henkel using a generalization of CMC foliations called prescribed mean curvature (PMC) foliations. A PMC foliation can be built that includes any given Cauchy surface [178] and global existence of PMC foliations can be proved in a way analogous to that previously done for CMC foliations [177, 176] . These global foliations provide barriers that imply the existence of a CMC hypersurface. Thus, in the end it turns out that the unwanted condition in the previous theorems on CMC foliations is in fact automatically satisfied. Connections between areal, CMC, and PMC time coordinates were further explored in [19] . One important observation there is that hypersurfaces of constant areal time in spacetimes with symmetry often have mean curvature of a definite sign. Related problems for the Einstein equations coupled to fields motivated by string theory have been studied by Narita [252, 253, 254, 255] .

Once global existence has been proved for a preferred time coordinate, the next step is to investigate the asymptotic behaviour of the solution as
$t\to {t}_{\pm}$
. There are few cases in which this has been done successfully. Notable examples are Gowdy spacetimes [113, 190, 116] and solutions of the Einstein–Vlasov system with spherical and plane symmetry [272] . These last results have been extended to allow a non-zero cosmological constant in [336] . Progress in constructing spacetimes with prescribed singularities will be described in Section 6 . In the future this could lead in some cases to the determination of the asymptotic behaviour of large classes of spacetimes as the singularity is approached. Detailed information has been obtained on the late-time behaviour of a class of inhomogeneous solutions of the Einstein–Vlasov system with positive cosmological constant in [337, 338] (see Section 7.6 ).

In the case of polarized Gowdy spacetimes a description of the late-time asymptotics was given in [116] . A proof of the validity of the asymptotic expansions can be found in [199] . The central object in the analysis of these spacetimes is a function
$P$
that satisfies the equation
${P}_{tt}+{t}^{-1}{P}_{t}={P}_{\theta \theta}$
. The picture that emerges is that the leading asymptotics are given by
$P=Alogt+B$
for constants
$A$
and
$B$
, this being the form taken by this function in a general Kasner model, while the next order correction consists of waves whose amplitude decays like
${t}^{-1/2}$
, where
$t$
is the usual Gowdy time coordinate. The entire spacetime can be reconstructed from
$P$
by integration. It turns out that the generalized Kasner exponents converge to
$(1,0,0)$
for inhomogeneous models. This shows that if it is stated that these models are approximated by Kasner models at late times it is necessary to be careful in what sense the approximation is supposed to hold.

General (non-polarized) Gowdy models, which are technically much more difficult to handle, have been analysed in [316] . Interesting and new qualitative behaviour was found. This is one of the rare examples where a rigorous mathematical approach has discovered phenomena which had not previously been suspected on the basis of heuristic and numerical work. In the general Gowdy model the function
$P$
is joined by a function
$Q$
and these two functions satisfy a coupled system of nonlinear wave equations. Assuming periodic boundary conditions the solution at a fixed time
$t$
defines a closed loop in the
$(P,Q)$
plane. (In fact it is natural to interpret it as the hyperbolic plane.) Thus the solution as a whole can be represented by a loop which moves in the hyperbolic plane. On the basis of what happens in the polarized case it might be expected that the following would happen at late times. The diameter of the loop shrinks like
${t}^{-1/2}$
while the centre of the loop, defined in a suitable way, moves along a geodesic. In [316] Ringström shows that there are solutions which behave in the way described but there are also just as many solutions which behave in a quite different way. The shrinking of the diameter is always valid but the way the resulting small loop moves is different. There are solutions where it converges to a circle in the hyperbolic plane which is not a geodesic and it continues to move around this circle forever. A physical interpretation of this behaviour does not seem to be known.

Ringström has also obtained important new results on the structure of singularities in Gowdy spacetimes. They are discussed in Section 8.4 .

4 Newtonian Theory and Special Relativity

To put the global results discussed in this article into context it is helpful to compare with Newtonian theory and special relativity. Some of the theorems that have been proved in those contexts and that can offer insight into questions in general relativity will now be reviewed. It should be noted that even in these simpler contexts open questions abound.

4.1 Hydrodynamics

Solutions of the classical (compressible) Euler equations typically develop singularities, i.e.

discontinuities of the basic fluid variables, in finite time [322] . Some of the results of [322] were recently generalized to the case of a relativistic fluid [167] . The proofs of the development of singularities are by contradiction and so do not give information about what happens when the smooth solution breaks down. One of the things that can happen is the formation of shock waves and it is known that, at least in certain cases, solutions can be extended in a physically meaningful way beyond the time of shock formation. The extended solutions only satisfy the equations in the weak sense. For the classical Euler equations there is a well-known theorem on global existence of weak solutions in one space dimension which goes back to [157] . This has been generalized to the relativistic case. Smoller and Temple treated the case of an isentropic fluid with linear equation of state [323] while Chen analysed the cases of polytropic equations of state [84] and flows with variable entropy [85] . This means that there is now an understanding of this question in the relativistic case similar to that available in the classical case.

In space dimensions higher than one there are no general global existence theorems. For a long time there were also no uniqueness theorems for weak solutions even in one dimension. It should be emphasized that weak solutions can easily be shown to be non-unique unless they are required to satisfy additional restrictions such as entropy conditions. A reasonable aim is to find a class of weak solutions in which both existence and uniqueness hold. In the one-dimensional case this has recently been achieved by Bressan and collaborators (see [68, 70, 69] and references therein).

It would be desirable to know more about which quantities must blow up when a singularity forms in higher dimensions. A partial answer was obtained for classical hydrodynamics by Chemin [83] .

The possibility of generalizing this to relativistic and self-gravitating fluids was studied by Brauer [66] .

There is one situation in which a smooth solution of the classical Euler equations is known to exist for all time. This is when the initial data are small and the fluid initially is flowing uniformly outwards. A theorem of this type has been proved by Grassin [161] . There is also a global existence result due to Guo [164] for an irrotational charged fluid in Newtonian physics, where the repulsive effect of the charge can suppress the formation of singularities.

A question of great practical interest for physics is that of the stability of equilibrium stellar models. The linear stability of a large class of static spherically symmetric solutions of the Einstein–Euler equations within the class of spherically symmetric perturbations has been proved by Makino [240] (cf. also [228] for the Newtonian problem). A nonlinear stability result for solutions of the Euler-Poisson system was proved in [276] under the assumption of global existence. The spectral properties of the linearized operator for general (i.e. non-spherically symmetric) perturbations in the Newtonian problem have been studied by Beyer [54] . This could perhaps provide a basis for a stability analysis, but this has not been done.

4.2 Kinetic theory

Collisionless matter is known to admit a global singularity-free evolution in many cases. For self-gravitating collisionless matter, which is described by the Vlasov–Poisson system, there is a general global existence theorem [266, 236] . There is also a version of this which applies to Newtonian cosmology [280] . A more difficult case is that of the Vlasov–Maxwell system, which describes charged collisionless matter. Global existence is not known for general data in three space dimensions, but has been shown in two space dimensions [154, 155] and in three dimensions with one symmetry [153] or with almost spherically symmetric data [270] .

A model system which has attracted some interest (see [18] ) is the Nordström–Vlasov system where the Vlasov equation is coupled to a scalar field as in Nordström's theory of gravitation. This is not a physically correct model but may be useful for obtaining mathematical insights. A similar procedure was used to look for numerical insights in [321] . At the moment the state of knowledge concerning this system can be summed up by saying that it is roughly equal to that available for the Vlasov–Maxwell system.

The nonlinear stability of static solutions of the Vlasov–Poisson system describing Newtonian self-gravitating collisionless matter has been investigated using the energy–Casimir method. For information on this see [165] and its references. The energy–Casimir method has been applied to the Einstein equations in [350] .

For the classical Boltzmann equation, global existence and uniqueness of smooth solutions has been proved for homogeneous initial data and for data that are small or close to equilibrium. For general data with finite energy and entropy, global existence of weak solutions (without uniqueness) was proved by DiPerna and Lions [132] . For information on these results and on the classical Boltzmann equation in general see [79, 80] . Despite the non-uniqueness it is possible to show that all solutions tend to equilibrium at late times. This was first proved by Arkeryd [27] by non-standard analysis and then by Lions [235] without those techniques. It should be noted that since the usual conservation laws for classical solutions are not known to hold for the DiPerna–Lions solutions, it is not possible to predict which equilibrium solution a given solution will converge to. In the meantime, analogues of several of these results for the classical Boltzmann equation have been proved in the relativistic case. Global existence of weak solutions was proved in [134] . Global existence and convergence to equilibrium for classical solutions starting close to equilibrium was proved in [156] . On the other hand, global existence of classical solutions for small initial data is not known. Convergence to equilibrium for weak solutions with general data was proved by Andréasson [15] . Until recently there was no existence and uniqueness theorem in the literature for general spatially homogeneous solutions of the relativistic Boltzmann equation. A paper claiming to prove existence and uniqueness for solutions of the Einstein–Boltzmann system which are homogeneous and isotropic [250] contains fundamental errors. These problems were corrected in [263] and a global existence theorem for the special relativistic Boltzmann equation was obtained.

In [262] this was generalized to a global existence theorem for LRS Bianchi type I solutions of the Einstein–Boltzmann system.

Further information on kinetic theory and its relation to general relativity can be found in the Living Review of Andréasson [17] .

4.3 Elasticity theory

There is an extensive literature on mathematical elasticity theory but the mathematics of self-gravitating elastic bodies seems to have been largely neglected. An existence theorem for spherically symmetric elastic bodies in general relativity was mentioned in Section 3.1 . More recently, Beig and Schmidt [42] proved an existence theorem for static elastic bodies subject to Newtonian gravity, which need not be spherically symmetric. This was extended to rotating bodies and special relativity in [44] .

5 Global Existence for Small Data

An alternative to symmetry assumptions is provided by “small data” results, where solutions are studied that develop from data close to those for known solutions. This leads to some simplification in comparison to the general problem, but with present techniques it is still very hard to obtain results of this kind.

5.1 Stability of de Sitter space

In [141] , Friedrich proved a result on the stability of de Sitter space. He gives data at infinity but the same type of argument can be applied starting from a Cauchy surface in spacetime to give an analogous result. This concerns the Einstein vacuum equations with positive cosmological constant and is as follows. Consider initial data induced by de Sitter space on a regular Cauchy hypersurface.

Then all initial data (vacuum with positive cosmological constant) near enough to these data in a suitable (Sobolev) topology have maximal Cauchy developments that are geodesically complete.

The result gives much more detail on the asymptotic behaviour than just this and may be thought of as proving a form of the cosmic no hair conjecture in the vacuum case. (This conjecture says roughly that the de Sitter solution is an attractor for expanding cosmological models with positive cosmological constant.) This result is proved using conformal techniques and, in particular, the regular conformal field equations developed by Friedrich. An alternative proof of this result which extends to all higher even dimensions was given in [6] . For some comments on the case of odd dimensions see [304] .

There are results obtained using the regular conformal field equations for negative or vanishing cosmological constant [143, 146] , but a detailed discussion of their nature would be out of place here (cf. however Section 9.1 ).

5.2 Stability of Minkowski space

Another result on global existence for small data is that of Christodoulou and Klainerman on the stability of Minkowski space [108] . The formulation of the result is close to that given in Section 5.1 , but now de Sitter space is replaced by Minkowski space. Suppose then that initial data for the vacuum Einstein equations are prescribed that are asymptotically flat and sufficiently close to those induced by Minkowski space on a hyperplane. Then Christodoulou and Klainerman prove that the maximal Cauchy development of these data is geodesically complete. They also provide a wealth of detail on the asymptotic behaviour of the solutions. The proof is very long and technical. The central tool is the Bel–Robinson tensor, which plays an analogous role for the gravitational field to that played by the energy-momentum tensor for matter fields. Apart from the book of Christodoulou and Klainerman itself, some introductory material on geometric and analytic aspects of the proof can be found in [65, 107] , respectively. The result for the vacuum Einstein equations was generalized to the case of the Einstein–Maxwell system by Zipser [354] .

In the original version of the theorem, initial data had to be prescribed on all of
${\mathbf{R}}^{3}$
. A generalization described in [211] concerns the case where data need only be prescribed on the complement of a compact set in
${\mathbf{R}}^{3}$
. This means that statements can be obtained for any asymptotically flat spacetime where the initial matter distribution has compact support, provided attention is confined to a suitable neighbourhood of infinity. The proof of the new version uses a double null foliation instead of the foliation by spacelike hypersurfaces previously used and leads to certain conceptual simplifications. A detailed treatment of this material can be found in the book of Klainerman and Nicolò [212] .

An aspect of all this work which seemed less than optimal was the following. Well-known heuristic analyses by relativists produced a detailed picture of the fall-off of radiation fields in asymptotically flat solutions of the Einstein equations, known as peeling. It says that certain components of the Weyl tensor decay at certain rates. The analysis of Christodoulou and Klainerman reproduced some of these fall-off rates but not all. More light was shed on this discrepancy by Klainerman and Nicolò [213] who showed that if the fall-off conditions on the initial data assumed in [108] are strengthened somewhat then peeling can be proved.

A much shorter proof of the stability of Minkowski space has been given by Lindblad and Rodnianski [233] . It uses harmonic coordinates and so is closer to the original local existence proof of Choquet-Bruhat. The fact that this approach was not used earlier is related to the fact that the null condition, an important structural condition for nonlinear wave equations which implies global existence for small data, is not satisfied by the Einstein equations written in harmonic coordinates.

Lindblad and Rodnianski formulated a generalization called the weak null condition [232] . This is only one element which goes into the global existence proof but it does play an important role. The result of Lindblad and Rodnianski does not give as much detail about the asymptotic structure as the approach of Christodoulou and Klainerman. On the other hand it seems that the proof generalizes without difficulty to the case of the Einstein equations coupled to a massless scalar field.

5.3 Stability of the (compactified) Milne model

The interior of the light cone in Minkowski space foliated by the spacelike hypersurfaces of constant Lorentzian distance from the origin can be thought of as a vacuum cosmological model, sometimes known as the Milne model. By means of a suitable discrete subgroup of the Lorentz group it can be compactified to give a spatially compact cosmological model. With a slight abuse of terminology the latter spacetime will also be referred to here as the Milne model. The stability of the latter model has been proved by Andersson and Moncrief (see [8, 11] ). The result is that, given data for the Milne model on a manifold obtained by compactifying a hyperboloid in Minkowski space, the maximal Cauchy developments of nearby data are geodesically complete in the future.

Moreover, the Milne model is asymptotically stable in the sense that any other solution in this class converges towards the Milne model in terms of suitable dimensionless variables.

The techniques used by Andersson and Moncrief are similar to those used by Christodoulou and Klainerman. In particular, the Bel–Robinson tensor is crucial. However, their situation is much simpler than that of Christodoulou and Klainerman, so that the complexity of the proof is not so great. This has to do with the fact that the fall-off of the fields towards infinity in the Minkowksi case is different in different directions, while it is uniform in the Milne case. Thus it is enough in the latter case to always contract the Bel–Robinson tensor with the same timelike vector when deriving energy estimates. The fact that the proof is simpler opens up a real possibility of generalizations, for instance by adding different matter models.

5.4 Stability of the Bianchi type III form of flat spacetime

Another vacuum cosmological model whose nonlinear stability has been investigated is the Bianchi III form of flat spacetime. To obtain this model, first do the construction described in Section 5.3 with the difference that the starting solution is three-dimensional Minkowski space.

Then, take the metric product of the resulting three-dimensional Lorentz manifold with a circle.

This defines a flat spacetime that has one Killing vector, which is the generator of rotations of the circle. It has been shown by Choquet-Bruhat and Moncrief [94] that this solution is stable under small vacuum perturbations preserving the one-dimensional symmetry. More precisely, they proved the result only for the polarized case. This restriction was lifted in [88] . As in the case of the Milne model, a natural task is to generalize this result to spacetimes with suitable matter content.