4.3 The BSSN formulation

The BSSN formulation is based on the 3+1 decomposition of Einstein’s field equations. Unlike the harmonic formulation, which has been motivated by the mathematical structure of the equations and the understanding of the Cauchy formulation in general relativity, this system has been mainly developed and improved based on its capability of numerically evolving spacetimes containing compact objects in a stable way. Interestingly, it turns out that in spite of the fact that the BSSN formulation is based on an entirely different motivation, mathematical questions like the well-posedness of its Cauchy problem can be answered, at least for most gauge conditions.

In the BSSN formulation, the three metric γij and the extrinsic curvature Kij are decomposed according to

γij = e4ϕ&tidle;γij, (4.49 ) ( 1 ) Kij = e4ϕ A&tidle;ij + -&tidle;γijK . (4.50 ) 3
Here, K = γijK ij and A&tidle; ij are the trace and the trace-less part, respectively, of the conformally-rescaled extrinsic curvature. The conformal factor 2ϕ e is determined by the requirement for the conformal metric &tidle;γij to have unit determinant. Aside from these variables one also evolves the lapse (α), the shift (i β) and its time derivative (Bi), and the variable
&tidle;Γ i := − ∂ &tidle;γij. (4.51 ) j
In terms of the operator ˆ∂0 = ∂t − βj ∂j the BSSN evolution equations are
2 μ μ ˆ∂0α = − α f (α,ϕ,x )(K − K0 (x )), ( ) (4.52 ) −4ϕ [ i i ] ij 1 2 ˆ∂0K = − e &tidle;D D&tidle;i α + 2∂iϕ ⋅ &tidle;D α + α A&tidle; A&tidle;ij + --K − αS, (4.53 ) 3 ˆ∂0βi = α2G (α, ϕ,xμ)Bi, (4.54 ) ˆ i − 4ϕ μ ˆ &tidle;i i i μ ∂0B = e H (α, ϕ,x )∂0Γ − η (B ,α, x ), (4.55 ) ˆ∂ ϕ = − α-K + 1∂ βk, (4.56 ) 0 6 6 k k 2 k ˆ∂0&tidle;γij = − 2α &tidle;Aij + 2 &tidle;γk(i∂j)β − -&tidle;γij∂kβ , (4.57 ) [ 3 ]TF ˆ∂0A&tidle;ij = e− 4ϕ αR&tidle;ij + αR ϕij − &tidle;DiD&tidle;j α + 4∂(iϕ ⋅D&tidle;j )α + αK &tidle;A − 2 αA&tidle; A&tidle;k + 2 &tidle;A ∂ βk − 2A&tidle; ∂ βk − αe −4ϕˆS , (4.58 ) ij ik j k(i j) 3 ij k ij i kl i 1 ij k kj i 2 ki j ˆ∂0&tidle;Γ = &tidle;γ ∂k∂lβ + -&tidle;γ ∂j∂kβ + ∂k &tidle;γ ⋅ ∂jβ − --∂k&tidle;γ ⋅ ∂jβ 3 [ 3 ] − 2A&tidle;ij∂jα + 2 α (m − 1)∂kA&tidle;ki − 2m-&tidle;DiK + m (&tidle;Γ iA&tidle;kl + 6A&tidle;ij∂jϕ) − Si. (4.59 ) 3 kl
Here, quantities with a tilde refer to the conformal three metric &tidle;γij, which is also used in order to raise and lower indices. In particular, D&tidle;i and k &tidle;Γ ij denote the covariant derivative and the Christoffel symbols, respectively, with respect to &tidle;γij. Expressions with a superscript T F refer to their trace-less part with respect to the conformal metric. Next, the sum &tidle;Rij + Rϕ ij represents the Ricci tensor associated with the physical three metric γ ij, where
1 ( ) &tidle;Rij = − -&tidle;γkl∂k∂l&tidle;γij + &tidle;γk(i∂j)&tidle;Γ k − &tidle;Γ (ij)k∂l&tidle;γlk + &tidle;γls 2&tidle;Γ kl(i&tidle;Γ j)ks + Γ&tidle;kis&tidle;Γ klj , (4.60 ) 2 Rϕij = − 2 &tidle;DiD&tidle;j ϕ − 2&tidle;γij &tidle;Dk &tidle;Dkϕ + 4 &tidle;DiϕD&tidle;j ϕ − 4&tidle;γij &tidle;Dk ϕD&tidle;k ϕ. (4.61 )
The term ˆ &tidle; i ∂0Γ in Eq. (4.55View Equation) is set equal to the right-hand side of Eq. (4.59View Equation). The parameter m in the latter equation modifies the evolution flow off the constraint surface by adding the momentum constraint to the evolution equation for the variable &tidle;Γ i. This parameter was first introduced in [10] in order to compare the stability properties of the BSSN evolution equations with those of the ADM formulation.

The gauge conditions, which are imposed on the lapse and shift in Eqs. (4.52View Equation, 4.54View Equation, 4.55View Equation), were introduced in [52Jump To The Next Citation Point] and generalize the Bona–Massó condition [62Jump To The Next Citation Point] and the hyperbolic Gamma driver condition [11]. It is assumed that the functions f (α, ϕ,xμ), G (α,ϕ,x μ) and H (α,ϕ, xμ) are strictly positive and smooth in their arguments, and that K (xμ) 0 and ηi(Bj, α,x μ) are smooth functions of their arguments. The choice

μ 2- μ m = 1, f (α,ϕ,x ) = α, K0 (x ) = 0, (4.62 ) μ 3 μ 4ϕ i j μ i G (α, ϕ,x ) = --2, H (α,ϕ,x ) = e , η (B ,α, x ) = ηB , (4.63 ) 4α
with η a positive constant, corresponds to the evolution system used in many black-hole simulations based on 1 + log slicing and the moving puncture technique (see, for instance, [423] and references therein). Finally, the source terms S, Sˆ ij and Si are defined in the following way: denoting by R (3) ij and (4) Rij the Ricci tensors belonging to the three-metric γij and the spacetime metric, respectively, and introducing the constraint variables
( ) H := 1- γij R (3ij)+ 2K2 − &tidle;Aij &tidle;Aij , (4.64 ) 2 3 j 2 j Mi := &tidle;D A&tidle;ij − -&tidle;DiK + 6A&tidle;ijD&tidle; ϕ, (4.65 ) i i ij3 C := &tidle;Γ + ∂j&tidle;γ , (4.66 )
the source terms are defined as
[ ]T F S := γijR (4)− 2H, ˆSij := R (4)+ &tidle;γ ∂ Ck , Si := 2 αm &tidle;γijMj − ∂ˆ0Ci. (4.67 ) ij ij k(i j)
For vacuum evolutions one sets S = 0, ˆSij = 0 and Si = 0. When matter fields are present, the Einstein field equations are equivalent to the evolution equations (4.52View Equation4.59View Equation) setting S = − 4πG (ρ + σ) N, ˆ TF Sij = 8πGN σij, i ik S = 16 πGN m α &tidle;γ jk and the constraints H = 8πGN ρ, Mi = 8 πGN ji and Ci = 0.

When comparing Cauchy evolutions in different spatial coordinates, it is very convenient to reformulate the BSSN system such that it is covariant with respect to spatial coordinate transformations. This is indeed possible; see [77, 82Jump To The Next Citation Point]. One way of achieving this is to fix a smooth background three-metric ˚ γ ij, similarly as in Section 4.1, and to replace the fields ϕ and Γ&tidle;i by the scalar and vector fields

( ) 1 γ i ij ϕ := 12-log ˚- , &tidle;Γ := − ˚Dj &tidle;γ , (4.68 ) γ
where γ and γ˚ denote the determinants of γij and γ˚ij, and ˚Dj is the covariant derivative associated to the latter. If γ˚ ij is flat16 and time-independent, the corresponding BSSN equations are obtained by replacing ˚ ∂k ↦→ Dk and &tidle;Γ kij ↦→ &tidle;Γ kij − Γ˚kij in Eqs. (4.52View Equation4.59View Equation, 4.60View Equation, 4.61View Equation, 4.64View Equation4.66View Equation).

4.3.1 The hyperbolicity of the BSSN evolution equations

In fact, the ADM formulation in the spatial harmonic gauge described in Section 4.2.3 and the BSSN formulation are based on some common ideas. In the covariant reformulation of BSSN just mentioned, the variable &tidle;Γ i is just the quantity V i defined in Eq. (4.40View Equation), where γij is replaced by the conformal metric γ˚ ij. Instead of requiring &tidle;Γ i to vanish, which would convert the operator on the right-hand side of Eq. (4.60View Equation) into a quasilinear elliptic operator, one promotes this quantity to an independent field satisfying the evolution equation (4.59View Equation) (see also the discussion below Equation (2.18) in [390]). In this way, the &tidle;γij − A&tidle;ij-block of the evolution equations forms a wave system. However, this system is coupled through its principal terms to the evolution equations of the remaining variables, and so one needs to analyze the complete system. As follows from the discussion below, it is crucial to add the momentum constraint to Eq. (4.59View Equation) with an appropriate factor m in order to obtain a hyperbolic system.

The hyperbolicity of the BSSN evolution equations was first analyzed in a systematic way in [373Jump To The Next Citation Point], where it was established that for fixed shift and densitized lapse,

12σϕ α = e (4.69 )
the evolution system (4.53View Equation, 4.56View Equation4.59View Equation) is strongly hyperbolic for σ > 0 and m > 1∕4 and symmetric hyperbolic for m > 1 and 6σ = 4m − 1. This was shown by introducing new variables and enlarging the system to a strongly or symmetric hyperbolic first-order one. In fact, similar first-order reductions were already obtained in [196, 188Jump To The Next Citation Point]. However, in [373Jump To The Next Citation Point] it was shown that the first-order enlargements are equivalent to the original system if the extra constraints associated to the definition of the new variables are satisfied, and that these extra constraints propagate independently of the BSSN constraints H = 0, Mi = 0 and Ci = 0. This establishes the well-posedness of the Cauchy problem for the system (4.69View Equation, 4.53View Equation, 4.56View Equation4.59View Equation) under the aforementioned conditions on σ and m. Based on the same method, a symmetric hyperbolic first-order enlargement of the evolution equations (4.52View Equation, 4.53View Equation, 4.56View Equation4.59View Equation) and fixed shift was obtained in [52Jump To The Next Citation Point] under the conditions f > 0 and 4m = 3f + 1 and used to construct boundary conditions for BSSN. First-order strongly-hyperbolic reductions for the full system (4.52View Equation4.59View Equation) have also been recently analyzed in [82Jump To The Next Citation Point].

An alternative and efficient method for analyzing the system consists in reducing it to a first-order pseudodifferential system, as described in Section 3.1.5. This method has been applied in [308Jump To The Next Citation Point] to derive a strongly hyperbolic system very similar to BSSN with fixed, densitized lapse and fixed shift. This system is then shown to yield a well-posed Cauchy problem. In [52Jump To The Next Citation Point] the same method was applied to the evolution system (4.52View Equation4.59View Equation). Linearizing and localizing, one obtains a first-order system of the form ˚s ˚ Uˆt = P (ik)ˆU = iβ ksUˆ + α Q (ik )Uˆ. The eigenvalues of Q(ik) are 0, ±i, √ -- ±i m, √ -- ±i μ, √ -- ± f, √ ---- ± GH, √ -- ± κ, where we have defined μ := (4m − 1)∕3 and κ := 4GH ∕3. The system is weakly hyperbolic provided that

f > 0, μ > 0, κ > 0, (4.70 )
and it is strongly hyperbolic if, in addition, the parameter m and the functions f, G, and H can be chosen such that the functions
κ m − 1 6(m − 1)κ f-−--κ, μ-−-κ-, 4m--−-3κ--- (4.71 )
are bounded and smooth. In particular, this requires that the nominators converge to zero at least as fast as the denominators when f → κ, μ → κ or 3 κ → 4m, respectively. Since κ > 0, the boundedness of κ ∕(f − κ) requires that f ⁄= κ. For the standard choice m = 1, the conditions on the gauge parameters leading to strong hyperbolicity are, therefore, f > 0, κ > 0 and f ⁄= κ. Unfortunately, for the choice (4.62View Equation, 4.63View Equation) used in binary black-hole simulations these conditions reduce to
e4ϕ ⁄= 2α, (4.72 )
which is typically violated at some two-surface, since asymptotically, α → 1 and ϕ → 0 while near black holes α is small and ϕ positive. It is currently not known whether or not the Cauchy problem is well posed if the system is strongly hyperbolic everywhere except at points belonging to a set of zero measure, such as a two-surface. Although numerical simulations based on finite-difference discretizations with the standard choice (4.62View Equation, 4.63View Equation) show no apparent sign of instabilities near such surfaces, the well-posedness for the Cauchy problem for the BSSN system (4.52View Equation4.59View Equation) with the choice (4.62View Equation, 4.63View Equation) for the gauge source functions remains an open problem when the condition (4.72View Equation) is violated. However, a well-posed problem could be formulated by modifying the choice for the functions G and H such that f ⁄= κ and f,κ > 0 are guaranteed to hold everywhere.

Yet a different approach to analyzing the hyperbolicity of BSSN has been given in [219Jump To The Next Citation Point, 220Jump To The Next Citation Point] based on a new definition of strongly and symmetric hyperbolicity for evolution systems, which are first order in time and second order in space. Based on this definition, it has been verified that the BSSN system (4.69View Equation, 4.53View Equation, 4.56View Equation4.59View Equation) is strongly hyperbolic for σ > 0 and m > 1∕4 and symmetric hyperbolic for 6σ = 4m − 1 > 0. (Note that this generalizes the original result in [373] where, in addition, m > 1 was required.) The results in [220Jump To The Next Citation Point] also discuss more general 3+1 formulations, including the one in [308] and construct constraint-preserving boundary conditions. The relation between the different approaches to analyzing hyperbolicity of evolution systems, which are first order in time and second order in space, has been analyzed in [221].

Strong hyperbolicity for different versions of the gauge evolution equations (4.52View Equation, 4.54View Equation, 4.55View Equation), where the normal operator ˆ∂ 0 is sometimes replaced by ∂ t, has been analyzed in [222Jump To The Next Citation Point]. See Table I in that reference for a comparison between the different versions and the conditions they are subject to in order to satisfy strong hyperbolicity. It should be noted that when m = 1 and ˆ∂0 is replaced by ∂t, additional conditions restricting the magnitude of the shift appear in addition to f > 0 and f ⁄= κ.

4.3.2 Constraint propagation

As mentioned above, the BSSN evolution equations (4.52View Equation4.59View Equation) are only equivalent to Einstein’s field equation if the constraints

i ℋ := H − 8πGN ρ = 0, ℳi := Mi − 8πGN ji = 0, C = 0 (4.73 )
are satisfied. Using the twice contracted Bianchi identities in their 3+1 decomposed form, Eqs. (4.45View Equation, 4.46View Equation), and assuming that the stress-energy tensor is divergence free, it is not difficult to show that the equations (4.52View Equation4.59View Equation) imply the following evolution system for the constraint fields [52Jump To The Next Citation Point, 220Jump To The Next Citation Point]:
ˆ 1- j 2 −4ϕ &tidle;ij k 2α- ∂0ℋ = − α D (α ℳj ) − αe A &tidle;γki∂jC + 3 K ℋ, (4.74 ) α3 ( [ ] ) ˆ∂0ℳj = ---Dj (α−2ℋ ) + αK ℳj + ℳi ∂jβi + Di α &tidle;γk(i∂j)Ck TF , (4.75 ) 3 ˆ∂0Ci = 2αm &tidle;γijℳj. (4.76 )
This is the constraint propagation system for BSSN, which describes the propagation of constraint violations, which are usually present in numerical simulations due to truncation and roundoff errors. There are at least three reasons for establishing the well-posedness of its Cauchy problem. The first reason is to show that the unique solution of the system (4.74View Equation4.76View Equation) with zero initial data is the trivial solution. This implies that it is sufficient to solve the constraints at the initial time t = 0. Then, any smooth enough solution of the BSSN evolution equations with such data satisfies the constraint propagation system with ℋ = 0, ℳj = 0 and Ci = 0, and it follows from the uniqueness property of this system that the constraints must hold everywhere and at each time. In this way, one obtains a solution to Einstein’s field equations. However, in numerical calculations, the initial constraints are not exactly satisfied due to numerical errors. This brings us to the second reason for having a well-posed problem at the level of the constraint propagation system; namely, the continuous dependence on the initial data. Indeed, the initial constraint violations give rise to constraint violating solutions; but, if these violations are governed by a well-posed evolution system, the norm of the constraint violations is controlled by those of the initial violations for each fixed time t > 0. In particular, the constraint violations must converge to zero if the initial constraint violations do. Since the initial constraint errors go to zero when resolution is increased (provided a stable numerical scheme is used to solve the constraints), this guarantees convergence to a constraint-satisfying solution.17 Finally, the third reason for establishing well-posedness for the constraint propagation system is the construction of constraint-preserving boundary conditions, which will be explained in detail in Section 6.

The hyperbolicity of the constraint propagation system (4.74View Equation4.76View Equation) has been analyzed in [220Jump To The Next Citation Point, 52Jump To The Next Citation Point, 81Jump To The Next Citation Point, 80Jump To The Next Citation Point], and [315Jump To The Next Citation Point] and shown to be reducible to a symmetric hyperbolic first-order system for m > 1∕4. Furthermore, there are no superluminal characteristic fields if 1∕4 < m ≤ 1. Because of finite speed of propagation, this means that BSSN with 1∕4 < m ≤ 1 (which includes the standard choice m = 1) does not possess superluminal constraint-violating modes. This is an important property, for it shows that constraint violations that originate inside black hole regions (which usually dominate the constraint errors due to high gradients at the punctures or stuffing of the black-hole singularities in the turducken approach [156, 81, 80]) cannot propagate to the exterior region.

In [353Jump To The Next Citation Point] a general result is derived, showing that under a mild assumption on the form of the constraints, strong hyperbolicity of the main evolution system implies strong hyperbolicity of the constraint propagation system, with the characteristic speeds of the latter being a subset of those of the former. The result does not hold in general if “strong” is replaced by “symmetric”, since there are known examples for which the main evolution system is symmetric hyperbolic, while the constraint propagation system is only strongly hyperbolic [108Jump To The Next Citation Point].

  Go to previous page Go up Go to next page