4.4 What has been done?4 Numerical Issues4.2 Construction of initial data

4.3 Evolution equations 

The initial data constructed in one of the ways described above must be evolved with the propagation equations extracted from the conformal field equations. As we remarked upon above, there are several ways to write the evolution equations, and therefore there are also several ways to implement them numerically. To date there exist two different implementations. One is based on the ADM-like formulation of the conformal field equations, while the other uses a tetrad formalism. In the ADM-like formulation developed by Hübner [82Jump To The Next Citation Point In The Article], one introduces the coefficients of the intrinsic metric on the space-like slices and their extrinsic curvature as variables, while in the other formulation by Frauendiener [38Jump To The Next Citation Point In The Article], a tetrad and connection coefficients with respect to that tetrad are the basic variables. Not much difference concerning efficiency and usefulness between these formulations is to be expected, provided the numerical schemes are the same. One minor advantage of the tetrad formulation might be the additional freedom in the orientation of the tetrad which can be used on tex2html_wrap_inline3905 to make the computational tetrad agree with the preferred tetrad on tex2html_wrap_inline3905 (see below).

So far, the discretization of the equations has been rather straightforward. One of the schemes which have been used is the higher dimensional Lax-Wendroff scheme, also called the rotated Richtmyer scheme, a discretization scheme with second order accuracy. It has been employed alone [38Jump To The Next Citation Point In The Article] or together with Strang splitting [82Jump To The Next Citation Point In The Article] to treat the principal part of the equations differently from the source part. Since a second order scheme requires much more computing resources compared to higher order methods to achieve the same accuracy, Hübner started to use the method of lines [73Jump To The Next Citation Point In The Article] with a fourth order scheme to compute the spatial derivatives and fourth order Runge-Kutta for the evolution in time. He reports [82Jump To The Next Citation Point In The Article] that the fourth order method is very much superior to the second order scheme in terms of efficiency. (The feasibility of the method of lines in relativity has been studied by Hungerbühler [84] using pseudo-spectral methods for the spatial derivatives and a combination of Adams-Bashforth and Adams-Moulton schemes for the time evolution.)

The conformal field equations and the propagation equations derived from them are quasi-linear. This implies that the characteristics of the system depend on the current solution and this, in turn, means that one has to be able to change the time-step tex2html_wrap_inline5276 between successive time-slices in order to keep a stable evolution scheme. This is necessary because the schemes are explicit schemes and, therefore, subject to the Courant-Friedrichs-Lewy condition which states that the numerical domain of dependence of a point should always include the analytical domain of dependence. This requirement already excludes the popular leapfrog scheme which is nevertheless used sometimes also for evolving the Einstein equations. A general criterion for computing the maximal time-step allowed in each iteration in arbitrary dimension has been derived in [38Jump To The Next Citation Point In The Article].

Another important point in the development of evolution codes is the numerical treatment of the boundaries. As explained already above, it is one of the advantages in the conformal approach that the outer boundary is not as influential as it is in the conventional approach using the standard Cauchy problem. It was also pointed out that in this case it is enough to impose a boundary condition which results in a numerically stable code because the outer boundary is located in the unphysical region and, therefore, cannot influence the physical space-time.

The proper way to treat the boundary is to prescribe conditions which are compatible with the full conformal field equations, in particular with their restriction to the boundary manifold. This has not been done so far. Since the outer boundary is not important for the physical effects, other ways of dealing with the boundary have been devised. One way is to forget about the restrictions of the conformal field equations to the boundary and to analyze the possible boundary conditions for the propagation equations. To first order, one can define ingoing and outgoing fields on the boundary. Then a sufficiently general boundary condition will be obtained by specifying the ingoing fields in terms of the outgoing ones. Although this boundary treatment is not compatible with the restriction of the conformal field equations to the boundary, it is compatible with the evolution equations. This means that the evolution can remain stable although the solution will not satisfy the constraints in the domain of influence of the boundary, which, however, is always in the unphysical part of space-time. This method has been used in [38Jump To The Next Citation Point In The Article] with satisfactory results. In particular, the boundary did not give rise to non-physical modes. These findings are in agreement with the analysis of numerical boundary conditions by Trefethen [141].

Another, very clever method for dealing with the boundary has been found by Hübner [81]. He realized that it is sufficient to solve the conformal field equations in the physical space-time only, and not necessary to solve them in the unphysical region as long as the characteristics remain such that the information created in the unphysical part of the computational domain cannot reach the physical part. Consequently, in his treatment the grid is divided into three zones: the inner zone, the outer zone, and a transition zone. The inner zone covers the physical space-time (flagged by a positive conformal factor) and some part of the adjacent unphysical region. On this part of the grid the conformal field equations are solved. In the outer zone, which is located in a neighbourhood of the grid boundary, one solves an advection equation which propagates outwards, off the grid. In the transition zone, a sufficiently smooth interpolation between these two systems of equations is solved. The effect is that the boundary condition which has to be imposed on the grid boundaries is very simple and that the noise which is generated in the transition region is propagated away from the physical region outwards towards the grid boundary.

Our next point is concerned with the extraction of the radiative information from the numerically generated data. This is the part of the entire numerical process where the superiority of the conformal approach becomes apparent. How does one determine the radiative field? First of all, one needs to find tex2html_wrap_inline3905 on the current time-slice. Since tex2html_wrap_inline3905 is the surface on which the conformal factor tex2html_wrap_inline3721 vanishes and since tex2html_wrap_inline3721 is explicitly known during the evolution the location of tex2html_wrap_inline3905 is a simple task. The next problem is concerned with the orientation of the tetrad on tex2html_wrap_inline3905 . The asymptotic quantities are defined with respect to a specific geometrically characterized tetrad, a Bondi frame. But, in general, this tetrad is completely unrelated to the ``computational'' tetrad used for the evolution. Therefore, one needs to find the transformation from one to the other at each point of tex2html_wrap_inline3905 . Without going into too much details (see [38Jump To The Next Citation Point In The Article, 42Jump To The Next Citation Point In The Article, 41Jump To The Next Citation Point In The Article]) we remark that most asymptotic quantities, in particular the radiation field, are of a local character so they can be read off without constructing a Bondi frame. This is rather fortunate because there are global issues involved in the transformation from the computational tetrad to the Bondi frame. These have implications for the determination of global quantities like the Bondi energy-momentum four-vector, but they have no effect on the radiation field, which is defined as that (complex) component of tex2html_wrap_inline3877 which is entirely intrinsic to tex2html_wrap_inline3905 :


Here, tex2html_wrap_inline5298 is a null-vector tangent to the generators of tex2html_wrap_inline3905, i.e., tex2html_wrap_inline5302, and tex2html_wrap_inline5304 is any complex space-like null-vector which is orthogonal to tex2html_wrap_inline5298 . It is useful to require the space-like vector to be tangent to the intersection of tex2html_wrap_inline3905 with the current time-slice. Augmenting these two vectors by a further real null-vector tex2html_wrap_inline5310 yields a null-tetrad which is fixed up to rotations in the (two-dimensional) tangent space of that intersection and boosts in the plane orthogonal to it. The behaviour of tex2html_wrap_inline3697 under these transformations is that of a GHP-weighted quantity [67, 117] with boost weight -2 and spin weight -2. This corresponds to the quadrupole-like character of the gravitational radiation field. However, tex2html_wrap_inline3697 really depends only on the null-vector tex2html_wrap_inline5298 . For, suppose we perform a null-rotation around tex2html_wrap_inline5298, then tex2html_wrap_inline5304 transforms into tex2html_wrap_inline5326 for some complex valued function tex2html_wrap_inline5184 on tex2html_wrap_inline3905 . But tex2html_wrap_inline3697 is invariant under this transformation. So in order to find tex2html_wrap_inline3697 it is only necessary to transform from the given computational tetrad to the tetrad specified above which is rather straightforward. In fact, the computation of tex2html_wrap_inline3697 involves only the combination of certain components of the gravitational field with powers of tex2html_wrap_inline5338 .

The final step in the correct determination of the radiation is to find the correct time parameter. Suppose we follow a specific null generator of tex2html_wrap_inline3905 crossing through successive time-slices. On each slice we compute tex2html_wrap_inline3697 on that generator. Then we obtain the radiation emitted by the source into the direction specified by the generator as a function of our computational time parameter. Since the time coordinate is rather arbitrary, this means that the wave form determined so far has no physical meaning. The problem is already present in Maxwell's theory: Suppose we have an emitter which sends out a pure sine wave. A detector far away from the source cannot determine the absolute frequency of the signal because the relative velocity of emitter and receiver might be non-zero but the detector should also find a pure sine signal. However, this will be true only if the detector records the signal as a function of proper time. Any other time parameter along the detector's world-line will not produce a pure sine.

What one needs to do in the general case is to select among all parameters along the generator a specific, geometrically distinguished one, namely a Bondi parameter. A generator and such a parameter along it can be understood as a certain limit of freely falling observers with proper time clocks as they move towards infinity [42]. Bondi parameters are obtained as solutions of an ordinary, linear, second order differential equation, which is conformally invariant.

The computation of the Bondi energy-momentum is a global procedure, i.e., it depends on properties of the entire cut of tex2html_wrap_inline3905 with the current time-slice. There are two steps involved in this procedure. First, one needs to obtain the asymptotic translation group (see e.g. [118]) on each cut. This provides four functions on the cut which are then, in a second step, integrated against the ``mass aspect'' which is another function obtained from the ``Coulomb'' part tex2html_wrap_inline5346 of the gravitational field, and the ``news function'' which is a combination of components of the Ricci tensor and connection coefficients. The first step, the determination of the translation group, is the global step because it involves solving a second order elliptic equation on the cut. These issues are discussed in more detail in [41].

All these procedures for finding the relevant data on tex2html_wrap_inline3905 have been worked out analytically and they have also been tested (at least in part) numerically [38Jump To The Next Citation Point In The Article]. The tests have been performed under the assumption that null-infinity admits toroidal cuts, which has the advantage that one can actually compare the numerical results with analytical expressions because a whole class of exact solutions [135Jump To The Next Citation Point In The Article, 80Jump To The Next Citation Point In The Article] is known to exist. Admittedly, such space-times are rather unphysical, but since most of the extraction procedures are local there is no doubt that they will also work in more realistic cases.

4.4 What has been done?4 Numerical Issues4.2 Construction of initial data

image Conformal Infinity
Jörg Frauendiener
© Max-Planck-Gesellschaft. ISSN 1433-8351
Problems/Comments to livrev@aei-potsdam.mpg.de