4.1 The harmonic formulation

We start by discussing the harmonic formulation of Einstein’s field equations. Like in the potential formulation of electromagnetism, where the Lorentz gauge allows one to cast Maxwell’s equations into a system of wave equations, it was observed early in [134, 269] that Einstein’s equations reduce to a system of wave equations when harmonic coordinates,
are used. There are many straightforward generalizations of these gauge conditions; one of them is to replace the right-hand side of Eq. (4.1) by given source functions  [178, 182, 202].

In order to keep general covariance, we follow [232] and choose a fixed smooth background metric with corresponding Levi-Civita connection , Christoffel symbols , and curvature tensor . Then, the generalized harmonic gauge condition can be rewritten as

In the particular case where and where the background metric is Minkowski in standard Cartesian coordinates, vanishes, and the condition reduces to the harmonic coordinate expression (4.1). However, unlike condition (4.1), Eq. (4.3) yields a coordinate-independent condition for any given vector field on spacetime since the difference between two connections is a tensor field. In terms of the difference, , between the dynamical and background metric, this tensor field can be expressed as
Of course, the coordinate-independence is now traded for the introduction of a background metric , and the question remains of how to choose and the vector field . A simple possibility is to choose and equal to the initial data for the metric, such that initially.

Einstein’s field equations in the gauge are equivalent to the wave system

where is the stress-energy tensor and Newton’s constant. This system is subject to the harmonic constraint

4.1.1 Hyperbolicity

For any given smooth stress-energy tensor , the equations (4.5) constitute a quasilinear system of ten coupled wave equations for the ten coefficients of the difference metric (or equivalently, for the ten components of the dynamical metric ) and, therefore, we can apply the results of Section 3 to formulate a (local in time) well-posed Cauchy problem for the wave system (4.5) with initial conditions

where and are two sufficiently-smooth symmetric tensor fields defined on the initial slice satisfying the requirement that has Lorentz signature such that and the induced metric , , on is positive definite . For detailed well-posed Cauchy formulations we refer the reader to the original work in [169]; see also [85], [164], and [246], which presents an improvement on the results in the previous references due to weaker smoothness assumptions on the initial data.

An alternative way of establishing the hyperbolicity of the system (4.5) is to cast it into first-order symmetric hyperbolic form [164, 18, 286]. There are several ways of constructing such a system; the simplest one is obtained [164] by introducing the first partial derivatives of as new variables,

Then, the second-order wave system (4.5) can be rewritten in the form
where l.o. are lower-order terms not depending on any derivatives of the state vector . The system of equations (4.9, 4.10, 4.11) constitutes a quasilinear first-order symmetric hyperbolic system for with symmetrizer given by the quadratic form
However, it should be noted that the symmetrizer is only positive definite if is; that is, only if the time evolution vector field is time-like. In many situations, this requirement might be too restrictive. Inside a Schwarzschild black hole, for example, the asymptotically time-like Killing field is space-like.

However, as indicated above, the first-order symmetric hyperbolic reduction (4.9, 4.10, 4.11) is not unique. A different reduction is based on the variables , where is the derivative of in the direction of the future-directed unit normal to the time-slices , and . This yields a first-order system, which is symmetric hyperbolic as long as the slices are space-like, independent of whether or not is time-like [18, 286].

4.1.2 Constraint propagation and damping

The hyperbolicity results described above guarantee that unique solutions of the nonlinear wave system (4.5) exist, at least for short times, and that they depend continuously on the initial data , . However, in order to obtain a solution of Einstein’s field equations one has to ensure that the harmonic constraint (4.3) is identically satisfied.

The system (4.5) is equivalent to the modified Einstein equations

where denotes the Ricci tensor, and where if the harmonic constraint holds. From the twice contracted Bianchi identities one obtains the following equation for the constraint variable ,
This system describes the propagation of constraint violations, which are present if is nonzero. For this reason, we call it the constraint propagation system, or subsidiary system. Provided the stress-energy tensor is divergence free, , this is a linear, second-order hyperbolic equation for . Therefore, it follows from the uniqueness properties of such hyperbolic problems that provided the initial data , satisfies the initial constraints
This turns out to be equivalent to solving plus the usual Hamiltonian and momentum constraints; see [169, 286]. Summarizing, specifying initial data , satisfying the constraints (4.15), the corresponding unique solution to the nonlinear wave system (4.5) yields a solution to the Einstein equations.

However, in numerical calculations, one cannot assume that the initial constraints (4.15) are satisfied exactly, due to truncation and roundoff errors. The propagation of these errors is described by the constraint propagation system (4.14), and hyperbolicity guarantees that for each fixed time of existence, these errors converge to zero if the initial constraint violation converges to zero, which is usually the case when resolution is increased. On the other hand, due to limited computer resources, one cannot reach the limit of infinite resolution, and from a practical point of view one does not want the constraint errors to grow rapidly in time for fixed resolution. Therefore, one would like to design an evolution scheme in which the constraint violations are damped in time, such that the constraint hypersurface is an attractor set in phase space. A general method for damping constraints violations in the context of first-order symmetric hyperbolic formulations of Einstein’s field equations was given in [74]. This method was then adapted to the harmonic formulation in [224]. The procedure proposed in  [224] consists in adding lower-order friction terms in Eq. (4.13), which damp constraint violations. Explicitly, the modified system reads

with the future-directed unit normal to the surfaces, and and real constants, where determines the timescale on which the constraint violations are damped.

With this modification the constraint propagation system reads

A mode analysis for linear vacuum perturbations of the Minkowski metric reveals [224] that for and all modes, except those, which are constant in space, are damped. Numerical codes based on the modified system (4.16) or similar systems have been used in the context of binary black-hole evolutions [335, 336, 286, 384, 36, 403, 320], the head-on collision of boson stars [323] and the evolution of black strings in five-dimensional gravity [279], among other references.

For a discussion on possible effects due to nonlinearities in the constraint propagation system; see [185].

4.1.3 Geometric issues

The results described so far guarantee the local-in-time unique existence of solutions to Einstein’s equations in harmonic coordinates, given a sufficiently-smooth initial data set . However, since general relativity is a diffeomorphism invariant theory, some questions remain. The first issue is whether or not the harmonic gauge is sufficiently general such that any solution of the field equations can be obtained by this method, at least for short enough time. The answer is affirmative [169, 164]. Namely, let , , be a smooth spacetime satisfying Einstein’s field equations such that the initial surface is spacelike with respect to . Then, we can find a diffeomorphism in a neighborhood of the initial surface, which leaves it invariant and casts the metric into the harmonic gauge. For this, one solves the harmonic wave map equation (4.2) with initial data

Since equation (4.2) is a second-order hyperbolic one, a unique solution exists, at least on some sufficiently-small time interval . Furthermore, choosing small enough, describes a diffeomorphism when restricted to its image. By construction, satisfies the harmonic gauge condition (4.3).

The next issue is the question of geometric uniqueness. Let and be two solutions of Einstein’s equations with the same initial data on , i.e., , . Are these solutions related, at least for small time, by a diffeomorphism? Again, the answer is affirmative [169, 164] because one can transform both solutions to harmonic coordinates using the above diffeomorphism without changing their initial data. It then follows by the uniqueness property of the nonlinear wave system (4.5) that the transformed solutions must be identical, at least on some sufficiently- small time interval. Note that this geometric uniqueness property also implies that the solutions are, at least locally, independent of the background metric. For further results on geometric uniqueness involving only the first and second fundamental forms of the initial surface; see [127], where it is shown that every such initial-data set satisfying the Hamiltonian and momentum constraints possesses a unique maximal Cauchy development.

Finally, we mention that results about the nonlinear stability of Minkowski spacetime with respect to vacuum and vacuum-scalar perturbations have been established based on the harmonic system [283, 284], offering an alternative proof to the one of [129].