Misplaced Pages

Revision as of 17:04, 4 July 2008

Numerical relativity is one of the branches of general relativity that uses numerical methods and algorithms to solve and analyze problems. To this end, supercomputers are often employed to study black holes, gravitational waves, neutron stars and many other phenomena governed by Einstein's Theory of General Relativity. A current active research in numerical relativity is the simulation of relativistic binaries and their associated gravitational waves, other branches are also quite active.

Overview

A primary goal of numerical relativity is to study spacetimes whose exact form is not known. The spacetimes so found computational can either be fully dynamical, stationary or static and may contain matter fields or vacuum. In the case of stationary and static solutions, numerical methods may also be used to study the stability of the equilibrium spacetimes. In the case of dynamical spacetimes, the problem may be divided into the initial value problem and the evolution, each requiring different methods.

Numerical relativity is applied to many areas, such as cosmological models, critical phenomena, perturbed black holes and neutron stars, and the coalescence of black holes and neutron stars, for example. In any of these cases, Einstein's equations can be formulated in several ways that allow us to evolve the dynamics. While Cauchy methods have received a majority of the attention, characteristic and Regge calculus based methods have also been used. All of these methods begin with a snapshot of the gravitational fields on some hypersurface, the initial data, and evolve these data to neighboring hypersurfaces.

Like all problems in numerical analysis, careful attention is paid to the stability and convergence of the numerical solutions. In this line, much attention is paid to the gauge conditions, coordinates, and various formulations of the Einstein equations and the effect they have on the ability to produce accurate numerical solutions.

Numerical relativity research is distinct from work on classical field theories as many techniques implemented in these areas are inapplicable in relativity. Many facets are however shared with large scale problems in other computational sciences like computational fluid dynamics, electromagnetics, and solid mechanics. Numerical relativists often work with applied mathematicians and draw insight from numerical analysis, scientific computation, partial differential equations, and geometry among other mathematical areas of specialization.

The "3+1" Approach to Numerical Relativity

Most active research in NR has involved a 3+1, or "space${\displaystyle \times }$time" decomposition of Einstein's equations. Here the spacetime of interest is sliced up into a stack (or "foliation") of spacelike three-dimensional hypersurfaces. That is, cross-sections are taken of the spacetime, so that each point in the spacetime lies on exactly one slice. Successive slices then represent--loosely speaking--different instants in time. More specifically, the hypersurfaces are parametrized by a coordinate usually denoted by t, which can be interpreted as a universal coordinate time. (Note that this time does not necessarily coincide with the proper time of any observer.)

If we now restrict attention to a single hypersurface, we see that it can be treated just like any other geometric space. For example, the space has some Riemannian metric. This three-dimensional metric is usually written as ${\displaystyle \gamma _{ij}}$, so that the line element is given by
${\displaystyle dl^{2}=\gamma _{ij}dx^{i}dx^{j}\ .}$
Here, the coordinates ${\displaystyle x^{i}}$ are just some coordinates placed on the hypersurface, in some smooth way, and the indices i and j run over the three dimensions of the hypersurface.

Next, we need some way to relate the successive hypersurfaces. We imagine a set of observers moving through time along the direction normal to the hypersurface. If we choose one such observer, and denote the proper time by ${\displaystyle \tau }$, we can relate it to the time coordinate t by way of the "lapse function" ${\displaystyle \alpha }$:
${\displaystyle d\tau =\alpha (t,x^{k})dt\ .}$
Similarly, we need to relate the position of the observer on one slice to its position on the next slice, using coordinates. If the change in time coordinate is dt, and the change in spatial coordinates is ${\displaystyle dx^{i}}$, we can relate the two as:
${\displaystyle dx^{i}=-\beta ^{i}(t,x^{k})\,dt\ .}$
The functions ${\displaystyle \beta ^{i}}$ are known as the "shift vector".

The way in which we foliate a spacetime is certainly not unique, as there is clearly no need (other than smoothness) for the spatial coordinates on neighboring slices to be related to each other uniquely. This means that the functions ${\displaystyle \alpha (t,x^{k})}$ and ${\displaystyle \beta ^{i}(t,x^{k})}$ can be freely chosen. These two functions determine our choice of coordinate system, and are known as "gauge functions" for this reason.

The metric of the full spacetime can now be written in the form
${\displaystyle ds^{2}=\left(-\alpha ^{2}+\beta ^{i}\beta ^{j}\gamma _{ij}\right)dt^{2}+2\beta ^{i}\gamma _{ij}dtdx^{j}+\gamma _{ij}dx^{i}dx^{j}\ .}$
In computations for General Relativity, the desired "solution" is just the metric, given as a function of the coordinates. From the above equation, then, we see that objective is to find ${\displaystyle \gamma _{ij}(t,x^{k})}$.

The mathematics of the 3+1 split of Einstein's equations was laid out in a series of papers by Arnowitt, Deser and Misner (ADM) in the early 1960s; the form actually used for most simulation efforts until the mid-1990s was described by York (1979). ADM showed that data could be specified on an initial hypersurface, and Einstein's equations then show how that initial data evolves onto later hypersurfaces. More specifically, ADM showed that Einstein's equations separate into four "constraint" equations, and a set of "evolution" equations.

The four constraint equations are equations which only refer to quantities defined on each individual hypersurface, and must be satisfied on each hypersurface--including the initial slice. The evolution equations then describe how those quantities change in time. (In technical terms, ADM showed that Einstein's equations form a constrained hyperbolic system.)

Computational Difficulties

It can be shown that if the constraint equations are satisfied initially, and the variables are evolved exactly according to the evolution equations, then the constraints will remain satisfied. In theory, this means that all we need to do is specify initial data which satisfy the constraints, and evolve according to the evolution equations, and we will have a correct solution to Einstein's equations. Unfortunately, discretising the problem produces complications.

Methodology

Elliptic Equations

In the case of producing stationary spacetimes, static spacetimes or producing initial data, three dimensional elliptic equations must be solved. In some cases, symmetry may effectively reduce the dimensionality of the problem. The well posedness of these problems depend crucially on the specification of the boundary values. It is common to suppose that isolated systems in general relativity are asymptotically flat, which gives dirichlet boundary conditions if compact coordinates are used or a Robin boundary condition which suitably approximates asymptotic flatness if the domain is truncated for computational simplicity. In the case of spacetimes containing a black hole, an inner surface is specified as an apparent horizon which yields a Robin boundary condition.

The resulting partial differential equations are then solved using finite difference methods, finite element methods, spectral methods or any other numerical technique for solving elliptic equations. Other specific problems in general relativity also require the solution of elliptic equations.

Evolution

In general, an evolution equation expresses the time derivative of a variable in terms of spatial derivatives of that variable and other variables. Implementing these equations on a computer, then, requires methods to take spatial derivatives of the functions (spatial differentiation), and use the resulting time derivative to evolve the variables (time stepping).

Spatial Differentiation

For D spatial dimensions in a simulation where physical phenomena are reduced to a set of difference equations (restricting attention to finite-difference methods) with a corresponding D-dimensional array of N "mesh points", the CPU time is of order ${\displaystyle {\mathcal {O}}(N^{D+1})}$ and the memory required is of order ${\displaystyle {\mathcal {O}}(N^{D})}$. Typically, the accuracy of finite-difference methods improves with the number of points as ${\displaystyle N^{2}}$ (second-order), ${\displaystyle N^{4}}$ (fourth-order), or some times ${\displaystyle N^{6}}$ (sixth-order).

An alternative technique is to approximate each variable in a spectral series. Examples of spectral series include Fourier series, expansion in terms of Chebyshev polynomials, and spherical harmonics. These approximations can then be differentiated exactly. The advantage of this technique is the speed with which its accuracy improves (its convergence rate); the accuracy of spectral methods typically improves as ${\displaystyle e^{N}}$. Compared to ${\displaystyle N^{2}}$ or ${\displaystyle N^{4}}$ for finite differencing, this is a vast improvement.

The disadvantage of spectral methods is their sensitivity to smoothness of the functions they are used with. Especially in situations involving matter, shocks and discontinuities can arise easily. In those cases, spectral methods are not worth the trouble, and finite differencing methods are preferred. Nonetheless, many simulations of solutions of the vacuum Einstein equations can be significantly improved by using spectral methods. The great majority of Numerical Relativity research groups currently use finite differencing methods, for reasons of simplicity.

Time Stepping

The most common techniques for time stepping are Iterated Crank-Nicolson and Runge-Kutta. Iterated Crank-Nicolson is a second-order convergent method, which means that its errors decrease with time step ${\displaystyle \Delta t}$ as ${\displaystyle 1/(\Delta t)^{2}}$. Runge-Kutta methods are somewhat more general. The "classic" Runge-Kutta algorithm is fourth-order convergent, but related algorithms can be made arbitrarily convergent.

There is a penalty, however, for moving to higher-order methods. These methods require more computations, and therefore run more slowly. The method used in a given simulation depends on balancing the accuracy with the time required. Many simulations which use second-order finite differencing, for example, use Iterated Crank-Nicolson time stepping because higher accuracy is not always needed.

See also

• Mathematics of general relativity
• Post-Newtonian expansion

Notes

1. York, J W "Kinematics and Dynamics of General Relativity" from "Sources of Gravitational Radiation: Proceedings of the Battelle Seattle Workshop" ed. Larry Smarr, pp. 83-126 (Cambridge University Press 1979)

External links

• Initial Data for Numerical Relativity — A review article which includes a nice (technical) discussion of numerical relativity.
• Rotating Stars in Relativity — A (technical) review article about rotating stars, with a section on numerical relativity applications.
• A Relativity Tutorial at Caltech — A basic introduction to concepts of Numerical Relativity.
• Numerical relativity on arxiv.org

• Computational physics
• Mathematical methods in general relativity