Advantages of a modified ADM formulation: Constraint propagation analysis of the Baumgarte-Shapiro-Shibata-Nakamura system

Advantages of a modified ADM formulation: Constraint propagation analysis of the Baumgarte-Shapiro-Shibata-Nakamura system

Gen Yoneda*

Department of Mathematical Sciences, Waseda University, Okubo, Shinjuku, Tokyo 169-8555, Japan Hisa-aki Shinkai^†

Computational Science Division Institute of Physical & Chemical Research (RIKEN), Hirosawa, Wako, Saitama 351-0198, Japan 共Received 1 April 2002; published 17 December 2002兲

Several numerical relativity groups are using a modified Arnowitt-Deser-Misner 共ADM兲 formulation for their simulations, which was developed by Nakamura and co-workers共and widely cited as the Baumgarte- Shapiro-Shibata-Nakamura system兲. This so-called BSSN formulation is shown to be more stable than the standard ADM formulation in many cases, and there have been many attempts to explain why this reformu- lation has such an advantage. We try to explain the background mechanism of the BSSN equations by using an eigenvalue analysis of constraint propagation equations. This analysis has been applied and has succeeded in explaining other systems in our series of works. We derive the full set of the constraint propagation equations, and study it in the flat background space-time. We carefully examine how the replacements and adjustments in the equations change the propagation structure of the constraints, i.e., whether violation of constraints 共if it exists兲will decay or propagate away. We conclude that the better stability of the BSSN system is obtained by their adjustments in the equations, and that the combination of the adjustments is in a good balance, i.e., a lack of their adjustments might fail to obtain the present stability. We further propose other adjustments to the equations, which may offer more stable features than the current BSSN equations.

DOI: 10.1103/PhysRevD.66.124003 PACS number共s兲: 04.25.Dm, 04.20.Fy

I. INTRODUCTION

One of the most important current topics in the field of numerical relativity is to find a formulation of the Einstein equations which gives us stable and accurate long-term evo- lution. We all know that simulating space-time and matter based on general relativity is the essential research direction to go in the future, but we do not have a definite recipe for controlling numerical blow-ups. We concentrate our discus- sion on the free evolution of the Einstein equations based on the 3⫹1 共space⫹time兲 decomposition of space-time, which requires solving the constraints only on the initial hypersur- face and monitors the violation共error兲 of the calculation by checking constraints during the evolution.

Over the decades, the Arnowitt-Deser-Misner共ADM兲 关1兴 formulation has been treated as the default by numerical rela- tivists.共More precisely, the version introduced by Smarr and York关2兴was taken as the default, which we denote the stan- dard ADM formulation hereafter.兲 Although the ADM for- mulation mostly works for gravitational collapse or cosmo- logical models in numerical treatments, it does not satisfy the requirement for long-term evolution, e.g., the studies of gravitational wave sources.

As we mentioned in our previous paper 关3兴, we think we can classify the current efforts of formulating equations for numerical relativity in the following three ways: 共i兲apply a modified ADM 关Baumgarte-Shapiro-Shibata-Nakamura 共BSSN兲兴formulation关4,5兴,共ii兲apply a first-order hyperbolic

formulation共see the references, e.g., in关6 – 8兴兲, or共iii兲apply an asymptotically constrained system关9–12兴.

The first refers to using a modified ADM formulation, originally proposed by Nakamura in the late 1980s, and sub- sequently modified by Nakamura-Oohara and Shibata- Nakamura 关4兴. This introduces conformal decomposition of the ADM variables, a new variable for calculating Ricci cur- vature, and adjusts the equations of motion using constraints.

The advantage of this formulation was reintroduced by Baumgarte and Shapiro 关5兴, and therefore this is often cited as the BSSN formulation, which we follow also. The BSSN equations are now widely used in the large-scale numerical computations, including coalescence of binary neutron stars 关13兴and binary black holes关14兴.

The second and third efforts use similar modifications such as introductions of new variables and/or adjustments of the equations, but differ in their purposes: to construct a hyperbolic formulation or to construct a formulation which constraints will decay or propagate away. The latter is in- tended to control numerical evolution so that the constrained manifold is its attractor. While the hyperbolic formulations have been extensively studied in this direction, we think the worrisome point in the discussion is the treatment of the nonprincipal part which is ignored in the hyperbolic formu- lation. As Kidder, Scheel, and Teukolsky 关8兴 reported re- cently, unless we reduce the effect of the nonprincipal part of the equations we may not gain any advantages from the hy- perbolic formulation for numerical results关6,15兴.

Through the series of studies 关3,6,12,16兴, we propose a systematic treatment for constructing a robust evolution sys- tem against perturbative error. We call it an asymptotically constrained共or asymptotically stable兲system if the error de- cays itself. The idea is to adjust evolution equations using

*Electronic address: [email protected]

†Electronic address: [email protected]

constraints 共we call this an adjusted system兲 and to decide the coefficients 共multipliers兲 by analyzing constraint propa- gation equations. We propose to apply an eigenvalue analysis of the propagation equations of the constraints, especially in its Fourier components, so as to include the nonprincipal part in the analysis. The characters of eigenvalues will be changed according to the adjustments to the original evolu- tion equations. We conjectured that the constraint violation that occurred during the evolution will decay共if negative real eigenvalues兲or propagate away共if pure imaginary eigenval- ues兲.

This conjecture was confirmed to explain the following numerical behaviors: wave propagation in the Maxwell equa- tions 关12兴, in the Ashtekar version of the Einstein equations 关12兴, and in the ADM formulation 共flat space-time back- ground兲 关16兴. The advantage of this construction scheme is that it can be applied to a formulation which is not a first- order hyperbolic form, such as to the ADM formulation 关3,16兴. We think, therefore, that our proposal is an alternative way to control or predict the violation of constraints. 共We believe that the idea of the constraint propagation analysis first appeared in Frittelli 关17兴, where she derived a hyperbo- licity classification for the standard ADM formulation.兲

The purpose of this paper is to apply this constraint propa- gation analysis to the BSSN formulation, and understand how each improvement contributes to more stable numerical evolution. Together with numerical comparisons with the standard ADM case 关18,19兴, this topic has been studied by many groups with different approaches. Using numerical test evolutions, Alcubierre et al.关20兴found that the essential im- provement is in the process of replacing terms by constraints, and that the eigenvalues of the BSSN evolution equations have fewer ‘‘zero eigenvalues’’ than those of ADM, and they conjectured that the instability can be caused by ‘‘zero eigen- values’’ that violate the ‘‘gauge mode.’’ Miller 关21兴applied von Neumann’s stability analysis to the plane-wave propaga- tion, and reported that BSSN has a wider range of parameters that give us stable evolution. These studies provide some support regarding the advantage of BSSN, while it also showed an example of an ill-posed solution in BSSN共as well as in ADM兲 关22兴.共Inspired by BSSN’s conformal decompo- sition, several related hyperbolic formulations have also been proposed关23–25兴.兲

We think our analysis will offer a new vantage point on the topic, and contribute an alternative understanding of its background. Consequently, we propose a more effective im- provement of the BSSN system that has not yet been tried in numerical simulations.

The construction of this paper is as follows. We review the BSSN system in Sec. II, and also we discuss where the adjustments are applied. In Sec. III, we apply our constraint propagation analysis to show how each improvement works in the BSSN equations, and in Sec. IV we extend our study to seek a better formulation which might be obtained by small steps. We only consider the vacuum space-time throughout the paper, but the inclusion of matter is straight- forward.

II. BSSN EQUATIONS AND THEIR CONSTRAINT PROPAGATION EQUATIONS

A. BSSN equations

We start by presenting the standard ADM formulation, which expresses the space-time with a pair of 3-metric ␥i j

and extrinsic curvature K_{i j}. The evolution equations become

⳵t

A␥i j⫽⫺2␣^Ki j⫹D_i␤j⫹D_j␤i, 共2.1兲

⳵t

AK_{i j}⫽␣^Ri j

ADM⫹␣^KKi j⫺2␣^KikK^k_j⫺D_iD_j␣

⫹共D_i␤^k兲K_{k j}⫹共D_j␤^k兲K_ki⫹␤^kDkK_{i j}, 共2.2兲 where ␣^,␤i are the lapse and shift function and D_i is the covariant derivative on 3-space. The symbol ⳵t

A means the time derivative defined by these equations, and we distin- guish them from those of the BSSN equations⳵t

B, which will be defined in Eqs. 共2.15兲–共2.19兲. The associated constraints are the Hamiltonian constraint H and the momentum con- straints Mi:

H^ADM⫽R^ADM⫹K²⫺K_{i j}K^{i j}, 共2.3兲 Mi

ADM⫽D_jK^j_i⫺D_iK. 共2.4兲 The widely used notation 关4,5兴 is to introduce the vari- ables (␸^,˜␥i j, K, A˜_{i j}, ⌫˜ⁱ) instead of (␥i j,K_{i j}), where

␸⫽共1/12兲log共det␥i j兲, 共2.5兲

␥

˜i j⫽e^⫺4␸␥i j, 共2.6兲 K⫽␥^{i jK}i j, 共2.7兲 A˜

i j⫽e^⫺4␸关K_{i j}⫺共1/3兲␥i jK兴, 共2.8兲

⌫

˜ⁱ⫽⌫˜

jk

i ^˜␥^jk. 共2.9兲

The new variable ⌫˜ⁱ was introduced in order to calculate Ricci curvature more accurately. ⌫˜ⁱalso contributes to mak- ing the system reproduce wave equations in its linear limit.

In the BSSN formulation, Ricci curvature is not calculated as R_{i j}^ADM⫽⳵k⌫i j

k⫺⳵i⌫k j k⫹⌫i j

l⌫lk k⫺⌫k j

l ⌫li

k, 共2.10兲

but

R_{i j}^BSSN⫽R_{i j}^␸⫹R˜_{i j}, 共2.11兲 R_{i j}^␸⫽⫺2D˜_iD˜_j␸⫺2␥^˜i jD˜^kD˜_k␸⫹4共D˜_i␸兲共D˜_j␸兲

⫺4␥^˜i j共D˜^k␸兲共D˜

k␸兲, 共2.12兲

R˜

i j⫽⫺共1/2兲␥^˜lk⳵l⳵k^˜␥i j⫹␥^˜k(i⳵j)⌫˜^k⫹⌫˜^k⌫˜

(i j)k

⫹2␥^˜lm⌫˜

l(i k ⌫˜

j)km⫹^˜␥^lm⌫˜

im k ⌫˜

kl j, 共2.13兲

where D˜_iis a covariant derivative associated with␥^˜i j. These are weakly equivalent, but R_{i j}^BSSNdoes have a wave operator

apparently in the flat background limit, so that we can expect more natural wave propagation behavior.

Additionally, BSSN requires us to impose the conformal factor as

␥

˜共ªdet^˜␥i j兲⫽1 共2.14兲 during the evolution. This is a kind of definition, but can also be thought of as a constraint. We will return to this point shortly.

BSSN’s improvements are not only the introductions of new variables, but also the replacement of terms in the evo- lution equations using the constraints. The purpose of this paper is to understand and to identify which improvement works for the stability. Before doing that, we first show the standard set of the BSSN evolution equations:

⳵t

B␸⫽⫺共1/6兲␣^K⫹共1/6兲␤ⁱ共⳵i␸兲⫹共⳵i␤ⁱ兲, 共2.15兲

⳵t

B␥^˜i j⫽⫺2␣^A˜i j⫹^˜␥ik共⳵j␤^k兲⫹^˜␥jk共⳵i␤^k兲

⫺共2/3兲^˜␥i j共⳵k␤^k兲⫹␤^k共⳵k␥^˜i j兲, 共2.16兲

⳵t

BK⫽⫺DⁱD_i␣⫹␣^A˜i jA˜^{i j}⫹共1/3兲␣^K2⫹␤ⁱ共⳵iK兲, 共2.17兲

⳵t BA˜

i j⫽⫺e^⫺4␸共D_iD_j␣兲^TF⫹e^⫺4␸␣共R_{i j}^BSSN兲^TF⫹␣^KA˜i j

⫺2␣^A˜ikA˜^k_j⫹共⳵i␤^k兲A˜_{k j}⫹共⳵j␤^k兲A˜_ki

⫺共2/3兲共⳵k␤^k兲A˜

i j⫹␤^k共⳵kA˜

i j兲, 共2.18兲

⳵t

B⌫˜ⁱ⫽⫺2共⳵j␣兲A˜^{i j}

⫹2␣关⌫˜

jk

i A˜^{k j}⫺共2/3兲␥^{˜i j}共⳵jK兲⫹6A˜^{i j}共⳵j␸兲兴

⫺⳵j关␤^k共⳵k␥^{˜i j}兲⫺^˜␥^{k j}共⳵k␤ⁱ兲⫺^˜␥^ki共⳵k␤^j兲

⫹共2/3兲^˜␥^{i j}共⳵k␤^k兲兴. 共2.19兲 We next summarize the constraints in this system. The normal Hamiltonian and momentum constraints 共the ‘‘kine- matic’’ constraints兲are naturally written as

H^BSSN⫽R^BSSN⫹K²⫺K_{i j}K^{i j}, 共2.20兲 Mi

BSSN⫽Mi

ADM, 共2.21兲

where we use the Ricci scalar defined by Eq. 共2.11兲. Addi- tionally, we regard the following three as the constraints共the

‘‘algebraic’’ constraints兲:

Gⁱ⫽⌫˜ⁱ⫺␥^˜jk⌫˜

jk

i , 共2.22兲

A⫽A˜

i j␥^{˜i j,} 共2.23兲

S⫽␥^˜⫺1, 共2.24兲 where the first two are from the algebraic definition of the variables共2.8兲and共2.9兲, and Eq.共2.24兲is from the require- ment of Eq.共2.14兲. Hereafter, we writeH^BSSNandM^BSSN simply as HandM, respectively.

Taking careful account of these constraints, Eqs. 共2.20兲 and共2.21兲can be expressed directly as

H⫽e^⫺4␸R˜⫺8e^⫺4␸D˜^jD˜

j␸⫺8e^⫺4␸共D˜^j␸兲共D˜

j␸兲

⫹共2/3兲K²⫺A˜_{i j}A˜^{i j}⫺共2/3兲AK, 共2.25兲 Mi⫽6A˜^j_i共D˜_j␸兲⫺2A共D˜_i␸兲⫺共2/3兲共D˜_iK兲

⫹^˜␥^{k j}共D˜_jA˜_ki兲. 共2.26兲 In summary, the fundamental dynamical variables in BSSN are (␸^,˜␥i j,K,A˜

i j,⌫˜ⁱ), a total of 17. The gauge quan- tities are (␣^,␤ⁱ), which is four, and the constraints are (H,Mi,Gⁱ,A,S), i.e., nine components. As a result, four共2 by 2兲components are left which correspond to two gravita- tional polarization modes.

B. Adjustments in evolution equations

Next, we show the BSSN evolution equations 共2.15兲– 共2.19兲again, identifying where the terms are replaced using the constraints 共2.20兲–共2.24兲.

By a straightforward calculation, we get

⳵t B␸⫽⳵t

A␸⫹共1/6兲␣A⫺共1/12兲␥^˜⫺1共⳵jS兲␤^j, 共2.27兲

⳵t B␥^˜i j⫽⳵t

A^˜␥i j⫺共2/3兲␣␥^˜i jA⫹共1/3兲^˜␥^⫺1共⳵kS兲␤^k˜␥i j, 共2.28兲

⳵t BK⫽⳵t

AK⫺共5/3兲␣^KA⫺␣H⫹␣^e⫺4␸共D˜

jG^j兲, 共2.29兲

⳵t BA˜

i j⫽⳵t AA˜

i j⫹关共1/3兲␣␥^˜i jK⫺共2/3兲␣^A˜i j兴A

⫹关共1/2兲␣^e⫺4␸共⳵k␥^˜i j兲⫺共1/6兲␣^e⫺4␸˜␥i j␥^˜⫺1共⳵kS兲兴G^k⫹␣^e⫺4␸˜␥k(i共⳵j)G^k兲⫺共1/3兲␣^e⫺4␸˜␥i j共⳵kG^k兲, 共2.30兲

⳵t B⌫˜ⁱ⫽⳵t

A⌫˜ⁱ⫹关⫺共2/3兲共⳵j␣兲^˜␥^ji⫺共2/3兲␣共⳵j^˜␥^ji兲⫺共1/3兲␣␥^˜ji˜␥^⫺1共⳵jS兲⫹4␣␥^{˜i j}共⳵j␸兲兴A⫺共2/3兲␣␥^˜ji共⳵jA兲⫹2␣␥^{˜i j}Mj

⫺共1/2兲共⳵k␤ⁱ兲␥^{˜k j}␥^˜⫺1共⳵jS兲⫹共1/6兲共⳵j␤^k兲^˜␥^{i j}␥^˜⫺1共⳵kS兲⫹共1/3兲共⳵k␤^k兲␥^{˜i j˜}␥^⫺1共⳵jS兲⫹共5/6兲␤^k˜␥^⫺2˜␥^{i j}共⳵kS兲共⳵jS兲

⫹共1/2兲␤^k˜␥^⫺1共⳵k␥^{˜i j}兲共⳵jS兲⫹共1/3兲␤^k˜␥^⫺1共⳵j^˜␥^ji兲共⳵kS兲, 共2.31兲

where⳵t

A denotes the part of no replacements, i.e., the terms only use the standard ADM evolution equations in its time derivatives.

From Eqs. 共2.27兲–共2.31兲, we understand that all the BSSN evolution equations are adjusted using constraints.

This fact will give us the importance of the scaling constraint S⫽0 and the trace-free operationA⫽0 during the evolution.

As we have pointed out in the case of adjusted ADM systems 关16,3兴, certain combinations of adjustments 共re- placements兲in the evolution equations change the eigenval- ues of constraint propagation equations drastically. For ex- ample, all negative eigenvalues can be negative real by applying Detweiler’s adjustment 关26兴 or its simplified ver- sion. One common fact we found is that such a case has an adjustment which breaks time-reversal parity of the original equation. That is, with a change of time integration direction

⳵t→⫺⳵t, an adjusted term might become effective if it breaks time-reversal symmetry. 共This time asymmetric fea- ture was first implemented as a ‘‘lambda-system’’ in 关9兴.兲 Unfortunately, for the case of the BSSN equations, Eqs.

共2.27兲–共2.31兲, all the above adjustments keep the time- reversal symmetry, so that we cannot expect direct decays of constraint violation in the present form. We will give the details on this point later.

III. CONSTRAINT PROPAGATION ANALYSIS IN FLAT SPACE-TIME

A. Procedures

We start this section overviewing the procedures and our goals. In our series of previous works 关3,12,16兴, we have concluded that eigenvalue analysis of the constraint propaga- tion equations is quite useful for explaining or predicting how the constraint violation grows.

Suppose we have a set of dynamical variables u^a(xⁱ,t) and their evolution equations

⳵tu^a⫽f共u^a,⳵iu^a, . . .兲, 共3.1兲 and the共first class兲constraints

C^␣共u^a,⳵iu^a, . . .兲⬇0. 共3.2兲 For monitoring the violation of constraints, we propose to investigate the evolution equations of C^␣ 共constraint propa- gation兲,

⳵tC^␣⫽g共C^␣,⳵iC^␣, . . .兲. 共3.3兲 关We do not mean to integrate Eq. 共3.3兲 numerically, but rather to evaluate it analytically in advance.兴In order to ana- lyze the contributions of all right-hand-side terms in Eq.

共3.3兲, we propose to reduce Eq.共3.3兲in ordinary differential equations by Fourier transformation,

⳵tCˆ^␣⫽gˆ共Cˆ^␣兲⫽M^␣_␤Cˆ^␤, 共3.4兲 where C(x,t)^␳⫽兰Cˆ (k,t)^␳exp(ik•x)d³k, and then to analyze the set of eigenvalues, say ⌳_␣, of the coefficient matrix, M^␣_␤, in Eq.共3.4兲. We call⌳’s and M^␣_␤ the constraint am- plification factors共CAFs兲of Eq.共3.3兲and constraint propa-

gation matrix, respectively. Our guidelines to have ‘‘better stability’’ are that共A兲if the CAFs have a negative real part 共the constraints are forced to be diminished兲, then we see more stable evolution than a system which has a positive real part, and共B兲if the CAFs have a nonzero imaginary part共the constraints are propagating away兲, then we see more stable evolution than a system which has zero CAFs. We found heuristically that the system becomes more stable when more

⌳’s satisfy the above criteria 关6,12兴. We note that these guidelines are confirmed numerically for wave propagation in the Maxwell system and in the Ashtekar version of the Einstein system 关12兴, and also for error propagation in Minkowskii space-time using adjusted ADM systems 关16兴. Supporting theorems for guideline 共A兲 were recently dis- cussed关31兴.

The above features of the constraint propagation, Eq.

共3.3兲, will differ when we modify the original evolution equations. If we add 共adjust兲 the evolution equations using constraints

⳵tu^a⫽f共u^a,⳵iu^a, . . .兲⫹F共C^␣,⳵iC^␣, . . .兲, 共3.5兲 then Eq. 共3.3兲will also be modified as

⳵tC^␣⫽g共C^␣,⳵iC^␣, . . .兲⫹G共C^␣,⳵iC^␣, . . .兲. 共3.6兲 Therefore, the problem is how to adjust the evolution equa- tions so that their constraint propagation satisfies the above criteria as much as possible.

B. BSSN constraint propagation equations

Our purpose in this section is to apply the above proce- dure to the BSSN system. The set of the constraint propaga- tion equations, ⳵t(H,Mi,Gⁱ,A,S)^T, turns to be quite long and not elegant共it is not a first-order hyperbolic and includes many nonlinear terms兲, and we put them in the Appendix. In order to understand the fundamental structure, we hereby show an analysis on the flat space-time background.

For the flat background metric g_␮␯⫽␩␮␯, the first-order perturbation equations of Eqs. 共2.27兲–共2.31兲 can be written as

⳵t(1)␾⫽⫺共1/6兲⁽¹⁾K⫹共1/6兲共␬_␾⫺1兲⁽¹⁾A, 共3.7兲

⳵t(1)␥^˜i j⫽⫺2⁽¹⁾A˜_{i j}⫺共2/3兲共␬␥˜⫺1兲␦i j(1)A, 共3.8兲

⳵t(1)K⫽⫺共⳵j⳵j(1)␣兲⫹共␬K1⫺1兲⳵j(1)G^j⫺共␬K2⫺1兲⁽¹⁾H, 共3.9兲

⳵t(1)A˜_{i j}⫽⁽¹⁾共R_{i j}^BSSN兲^TF⫺⁽¹⁾共D˜_iD˜_j␣兲^TF

⫹共␬A1⫺1兲␦k(i共⳵j) (1)G^k兲

⫺共1/3兲共␬A2⫺1兲␦i j共⳵k(1)G^k兲, 共3.10兲

⳵t(1)⌫˜ⁱ⫽⫺共4/3兲共⳵i(1)K兲⫺共2/3兲共␬_⌫˜ 1⫺1兲共⳵i(1)A兲

⫹2共␬_⌫˜ 2⫺1兲⁽¹⁾Mi, 共3.11兲

where we introduced parameters␬^{’s, all}␬⫽0 reproduce the no-adjustment case from the standard ADM equations, and all ␬⫽1 correspond to the BSSN equations. We express them as

␬adjª共␬␸,␬␥˜,␬K1,␬K2,␬A1,␬A2,␬_⌫˜ 1,␬_⌫˜ 2兲. 共3.12兲 Constraint propagation equations at the first order in the flat space-time, then, become

⳵t(1)H⫽关␬_␥^˜⫺共2/3兲␬_⌫˜ 1⫺共4/3兲␬_␸⫹2兴⳵j⳵j(1)A

⫹2共␬_⌫˜ 2⫺1兲共⳵j(1)Mj兲, 共3.13兲

⳵t(1)Mi⫽关⫺共2/3兲␬K1⫹共1/2兲␬A1

⫺共1/3兲␬A2⫹共1/2兲兴⳵i⳵j(1)G^j⫹共1/2兲␬A1⳵j⳵j(1)Gⁱ

⫹关共2/3兲␬K2⫺共1/2兲兴⳵i(1)H, 共3.14兲

⳵t(1)Gⁱ⫽2␬_⌫˜ 2(1)Mi⫹关⫺共2/3兲␬_⌫˜ 1⫺共1/3兲␬^˜_␥兴共⳵i(1)A兲, 共3.15兲

⳵t

(1)S⫽⫺2␬␥˜(1)A, 共3.16兲

⳵t(1)A⫽共␬A1⫺␬A2兲共⳵j(1)G^j兲. 共3.17兲 We will discuss the CAFs of Eqs.共3.13兲–共3.17兲.

C. Effect of adjustments

We check the CAFs of the BSSN equations in detail. The list of examples is shown also in Table I. Hereafter, we let k²⫽k_x²⫹k_y²⫹k_z² for Fourier wave numbers.

共i兲No-adjustment case,␬adj⫽(all zeros). This is the start- ing point of the discussion. In this case,

X_CAFs⫽„0共⫻7兲,⫾

冑

⫺k²…,

i.e., „0(⫻7),⫾pure imaginary(one pair)…. In the standard ADM formulation, which uses (␥i j,K_{i j}), CAFs are (0,0,

⫾pure imaginary)关16兴. Therefore, if we do not apply adjust-

ments in the BSSN equations, the constraint propagation structure is quite similar to that of the standard ADM formal- ism.

共ii兲For the BSSN equations,␬adj⫽(all 1s), X_CAFs⫽„0共⫻3兲,⫾

冑

⫺k²共three pairs兲…,

i.e., „0(⫻3),⫾pure imaginary(three pairs)…. The number of pure imaginary CAFs is increased over that of case 共i兲, and we conclude this is the advantage of adjustments used in the BSSN equations.

共iii兲NoS-adjustment case. All the numerical experiments so far apply the scaling conditionSfor the conformal factor

␸^{. The}S-originated terms appear many places in the BSSN equations 共2.15兲–共2.19兲, so that we suspect nonzero S is a kind of source of the constraint violation. However, since all S-originated terms do not appear in the flat space-time back- ground analysis关no adjusted terms in Eqs.共3.7兲–共3.11兲兴, our analysis is independent of theSconstraint.共Note that we do not deny the effect ofSadjustment in other situations.兲

共iv兲NoA-adjustment case. The trace共or traceout兲condi- tion for the variables is also considered necessary共e.g.,关27兴兲. This can be checked with ␬adj⫽(␬^,␬^,1,1,1,1,␬,1), and we get

X_CAFs⫽„0共⫻3兲,⫾

冑

⫺k²共three pairs兲…,

independent of ␬. Therefore, the effect of A adjustment is unimportant according to this analysis, i.e., on flat space- time background.共Note that we do not deny the effect ofA adjustment in other situations.兲

共v兲No Gⁱ-adjustment case. The introduction of ⌫ⁱ is the key in the BSSN system. If we do not apply adjustments by Gⁱ关␬adj⫽(1,1,0,1,0,0,1,1)兴, then we get

X_CAFs⫽„0共⫻7兲,⫾

冑

⫺k²…,

which is the same as case共i兲. That is, adjustments due toGⁱ terms are effective to make progress from the ADM method.

TABLE I. Summary of Sec. III C: contributions of adjustment terms and effects of introductions of new constraints in the BSSN system.

The center column indicates whether each constraint is taken as a component of constraints in each constraint propagation analysis共‘‘use’’兲, and whether each adjustment is on共‘‘adj’’兲. The column ‘‘diag?’’ indicates diagonalizability of the constraint propagation matrix. The right column shows CAFs, where Im and Re mean pure imaginary and real eigenvalue, respectively. Case共0兲 共standard ADM兲is shown in关16兴.

No. Constraints共number of components兲 diag? CAFs

in text H共1兲 Mi共3兲 Gⁱ共3兲 A共1兲 S共1兲 in Minkowskii background

共0兲 standard ADM use use yes 共0,0,Im,Im兲

共i兲 BSSN no adjustment

use use use use use yes 共0,0,0,0,0,0,0,Im,Im兲

共ii兲 the BSSN use⫹adj use⫹adj use⫹adj use⫹adj use⫹adj no 共0,0,0,Im,Im,Im,Im,Im,Im兲 共iii兲 noSadjustment use⫹adj use⫹adj use⫹adj use⫹adj use no no difference in flat background 共iv兲 noAadjustment use⫹adj use⫹adj use⫹adj use use⫹adj no 共0,0,0,Im,Im,Im,Im,Im,Im兲

共v兲 noGⁱ adjustment

use⫹adj use⫹adj use use⫹adj use⫹adj no 共0,0,0,0,0,0,0,Im,Im兲 共vi兲 noMiadjustment use⫹adj use use⫹adj use⫹adj use⫹adj no 共0,0,0,0,0,0,0,Re,Re兲 共vii兲 noHadjustment use use⫹adj use⫹adj use⫹adj use⫹adj no (0,0,0,Im,Im,Im,Im,Im,Im)

共vi兲 No Mi-adjustment case. This can be checked with

␬adj⫽(1,1,1,1,1,1,1,␬), and we get X_CAFs⫽„0,⫾

冑

⫺␬^k2共two pairs…,

⫾

冑

⫺k²共⫺1⫹4␬⫹兩1⫺4␬兩兲/6,

⫾

冑

⫺k²共⫺1⫹4␬⫺兩1⫺4␬兩兲/6兲.

If ␬⫽0, then „0(⫻7),⫾

冑

^k2/3…, which is „0(⫻7),

⫾real value…. Interestingly, these real values indicate the ex- istence of the error-growing mode together with the decaying mode. Alcubierre et al.关20兴found that the adjustment due to the momentum constraint is crucial for obtaining stability.

We think that they picked up this error-growing mode. For- tunately at the BSSN limit (␬⫽1), this error-growing mode disappears and turns into a propagation mode.

共vii兲 No H-adjustment case. The set ␬adj

⫽(1,1,1,␬,1,1,1,1) gives

X_CAFs⫽„0共⫻3兲,⫾

冑

⫺k²共3 pairs兲…,

independently of ␬. Therefore the effect ofHadjustment is unimportant according to this analysis, i.e., on a flat space- time background.共Note again that we do not deny the effect of Hadjustment in other situations.兲

These tests are on the effects of adjustments. We will consider whether much better adjustments are possible in the next section.

We list the above results in Table I. 共Table I includes a column of diagonalizability of constraint propagation matrix M, the importance of which was pointed out in 关31兴.兲 The most characteristic points of the above are 共v兲 and 共vi兲, which denote the contribution of the momentum constraint adjustment and the importance of the new variable ⌫˜ⁱ. It is quite interesting that the unadjusted BSSN equations 关case 共i兲兴 does not have apparent advantages from the ADM sys- tem. As we showed in共v兲and共vi兲, if we missed a particular adjustment, then the expected stability behavior occasionally gets worse than the starting ADM system. Therefore, we conclude that the better stability of the BSSN formulations is obtained by their adjustments in the equations, and the com- bination of the adjustments is in a good balance. That is, a lack of their adjustments might fail to bring about the stabil- ity of their system.

IV. PROPOSALS OF IMPROVED BSSN SYSTEMS In this section, we consider the possibility of whether we can obtain a system which has much better properties, whether more pure imaginary CAFs or negative real CAFs.

A. Heuristic examples

1. A system which has eight pure imaginary CAFs One direction is to seek a set of equations which make fewer zero CAFs than the standard BSSN case关point共ii兲in the previous section兴. Using the same set of adjustments in Eqs. 共3.7兲–共3.11兲, CAFs are written in general as

X_CAFs⫽„0,⫾

冑

⫺k²␬A1␬_⌫˜ 2共2 pairs兲,

⫾complicated expression,

⫾complicated expression….

The terms in the first line certainly give four pure imaginary CAFs 共two positive and negative real pairs兲 if ␬A1␬_⌫˜ 2⬎0 (⬍0). Keeping this in mind, by choosing ␬adj

⫽(1,1,1,1,1,␬,1,1), we find

X_CAFs⫽„0,⫾

冑

⫺k²共2 pairs兲,

⫾

冑

⫺k²共2⫹␬⫹兩␬⫺4兩兲/6,

⫾

冑

⫺k²共2⫹␬⫺兩␬⫺4兩兲/6…. Therefore, the adjustment␬adj⫽(1,1,1,1,1,4,1,1) gives

X_CAFs⫽„0,⫾

冑

⫺k²共4 pairs兲…,

which is one step advanced from BSSN’s according our guidelines.

We note that such a system can be obtained in many ways, e.g., ␬adj⫽(0,0,1,0,2,1,0,1/2) also gives four pairs of pure imaginary CAFs.

2. A system which has negative real CAF

One criterion to obtain a decaying constraint mode 共i.e., an asymptotically constrained system兲is to adjust an evolu- tion equation as it breaks time-reversal symmetry关16,3兴. For example, we consider an additional adjustment to the BSSN equation as

⳵t␥^˜i j⫽⳵t

B^˜␥i j⫹␬SD␣␥^˜i jH, 共4.1兲

which is a similar adjustment of the simplified Detweiler type共SD兲 关26兴that was discussed in关3兴. The first-order con- straint propagation equations on the flat background space- time become

⳵t(1)H⫽⳵j⳵j(1)A⫺共3/2兲␬SD⳵j⳵j(1)H,

⳵t

(1)Mi⫽共1/6兲⳵i

(1)H⫹共1/2兲⳵j⳵j (1)Gⁱ,

⳵t

(1)Gⁱ⫽⫺⳵i

(1)A⫹共1/2兲␬SD⳵i

(1)H⫹2⁽¹⁾Mi,

⳵t(1)A⫽⫺共⳵j⳵j(1)␣兲^TF⫹共⁽¹⁾R^BSSN_{j j} 兲^TF,

⳵t(1)S⫽⫺2⁽¹⁾A⫹3␬SD(1)H,

where we wrote only additional terms to Eqs.共3.13兲–共3.17兲. The CAFs become

X_CAFs⫽„0共⫻2兲,⫾

冑

⫺k²共three pairs兲,共3/2兲k²␬SD…, in which the last one becomes negative real if ␬SD⬍0.

3. Combination of Secs. IV A 1 and IV A 2 Naturally we next consider both adjustments,

⳵t␥^˜i j⫽⳵t

B␥^˜i j⫹␬SD␣␥^˜i jH, 共4.2兲

⳵tA˜_{i j}⫽⳵t

BA˜_{i j}⫺␬8␣^e⫺4␸˜␥i j⳵kG^k, 共4.3兲 where the second one produces the eight pure imaginary CAFs. The additional terms in the constraint propagation equations 共3.13兲–共3.17兲are

⳵t(1)H⫽⳵j⳵j(1)A⫺共3/2兲␬SD⳵j⳵j(1)H,

⳵t(1)Mi⫽共1/6兲⳵i(1)H⫹共1/2兲⳵j⳵j(1)Gⁱ⫺␬8⳵i⳵k(1)G^k,

⳵t(1)Gⁱ⫽⫺⳵i(1)A⫹共1/2兲␬SD⳵i(1)H⫹2⁽¹⁾Mi,

⳵t(1)A⫽⫺3␬8⳵k(1)G^k,

⳵t(1)S⫽⫺2⁽¹⁾A⫹3␬SD(1)H. We then obtain

X_CAFs⫽„0,⫾

冑

⫺k²共3 pairs兲, 共3/4兲k²␬SD⫾

冑

^k2关⫺␬8⫹共9/16兲k²␬SD

2 兴…,

which reproduces Sec. IV A 1 when␬SD⫽0,␬8⫽1, and Sec.

IV A 2 when ␬8⫽0. These CAFs can become „0, pure

imaginary 共three pairs兲, complex numbers with a negative real part共one pair兲…, with an appropriate combination of␬8

and␬SD.

B. Possible adjustments

In order to break time-reversal symmetry of the evolution equations关3,9,16兴, the possible simple adjustments are共i兲to add H, S, or Gⁱ terms to the equations of ⳵t␾^, ⳵t^˜␥i j, or

⳵t⌫˜ⁱ, or共ii兲to addMiorAterms to⳵tK or⳵tA˜_{i j}. We write them generally, including the proposal of Sec. IV A 2, as

⳵t␾⫽⳵t

B␾⫹␬_␾H␣H⫹␬_␾G␣^D˜kG^k, 共4.4兲

⳵t^˜␥i j⫽⳵t

B␥^˜i j⫹␬SD␣␥^˜i jH⫹␬␥˜G1␣␥^˜i jD˜_kG^k

⫹␬_␥^˜_G2␣␥^˜k(iD˜

j)G^k⫹␬_␥^˜_S1␣␥^˜i jS⫹␬^˜_␥S2␣^D˜iD˜

jS, 共4.5兲

⳵tK⫽⳵t

BK⫹␬KM␣␥^˜jk共D˜_jMk兲, 共4.6兲

⳵tA˜_{i j}⫽⳵t

BA˜_{i j}⫹␬AM1␣␥^˜i j共D˜^kMk兲⫹␬AM2␣共D˜_(iMj)兲

⫹␬AA1␣␥^˜i jA⫹␬AA2␣^D˜iD˜

jA, 共4.7兲

⳵t⌫˜ⁱ⫽⳵t

B⌫˜ⁱ⫹␬_⌫^˜_H␣^D˜iH⫹␬_⌫^˜_G1␣Gⁱ

⫹␬_⌫^˜_G2␣^D˜jD˜jGⁱ⫹␬_⌫^˜_G3␣^D˜iD˜jG^j, 共4.8兲 TABLE II. Possible adjustments which make a real-part CAF negative共Sec. IV B兲. The column of adjustments is nonzero multipliers in terms of Eqs.共4.4兲–共4.8兲, which all violate time-reversal symmetry of the equation. The column ‘‘diag?’’ indicates diagonalizability of the constraint propagation matrix. Neg./Pos. means negative/positive, respectively.

Adjustment CAFs diag? Effect of the adjustment

⳵t␾ ␬_␾H␣H „0,0,⫾冑⫺k²(*^3),8␬_␾Hk²… no ␬_␾H⬍0 makes 1 Neg.

⳵t␾ ␬␾G␣D˜

kG^k „0,0,⫾冑⫺k²(*2),long expressions… yes ␬_␾G⬍0 makes 2 Neg. 1 Pos.

⳵t˜␥

i j ␬SD␣␥˜

i jH „0,0,⫾冑⫺k²(*^3),(3/2)␬SDk²… yes ␬SD⬍0 makes 1 Neg. Sec. IV A 2

⳵t˜␥_{i j} ␬_␥^˜_G1␣␥˜_{i j}D˜_kG^k „0,0,⫾冑⫺k²(*2),long expressions… yes ␬_␥^˜G1⬎0 makes 1 Neg.

⳵t˜␥

i j ␬_␥^˜G2␣␥˜

k(iD˜

j)G^k „0,0,(1/4)k²␬_␥^˜G2⫾冑^k2(⫺1⫹k²␬_␥^˜G2/16)

⫻(*2),long expressions…

yes ␬_␥^˜_G2⬍0 makes 6 Neg. 1 Pos. Sec. IV B 2 a

⳵t˜␥

i j ␬_␥^˜S1␣␥˜

i jS „0,0,⫾冑⫺k²(*^3),3␬_␥^˜S1… no ␬_␥^˜_S1⬍0 makes 1 Neg.

⳵t˜␥_{i j} ␬␥˜S2␣D˜_iD˜_jS „0,0,⫾冑⫺k²(*^3),⫺␬_␥^˜_S2k²… no ␬_␥^˜S2⬎0 makes 1 Neg.

⳵tK ␬KM␣␥˜^jk(D˜

jMk) „0,0,0,⫾冑⫺k²(*^{2), (1/3)}␬KMk²

⫾(1/3)

冑

^k2(⫺9⫹k²␬KM2

)…

no ␬KM⬍0 makes 2 Neg.

⳵tA˜

i j ␬AM1␣␥˜

i j(D˜^kMk) „0,0,⫾冑⫺k²(*^3),⫺␬AM1k²… yes ␬AM1⬎0 makes 1 Neg.

⳵tA˜

i j ␬AM2␣(D˜

(iMj)) „0,0,⫺k²␬AM2/4⫾冑^k2(⫺1⫹k²␬AM2/16)

⫻(*2),long expressions…

yes ␬AM2⬎0 makes 7 Neg Sec. IV B 1

⳵tA˜

i j ␬AA1␣␥˜

i jA „0,0,⫾冑⫺k²(*^3),3␬AA1… yes ␬AA1⬍0 makes 1 Neg.

⳵tA˜_{i j} ␬AA2␣D˜_iD˜_jA „0,0,⫾冑⫺k²(*^3),⫺␬AA2k²… yes ␬AA2⬎0 makes 1 Neg.

⳵t⌫˜ⁱ ␬⌫˜H␣D˜ⁱH „0,0,⫾冑⫺k²(*^3),⫺␬AA2k²… no ␬_⌫^˜_H⬎0 makes 1 Neg.

⳵t⌫˜ⁱ ␬⌫˜G1␣Gⁱ „0,0,(1/2)␬⌫˜G1⫾

冑

⫺k²⫹␬_⌫^˜_G1

2 (*2),long exp.… yes ␬⌫˜G1⬍0 makes 6 Neg. 1 Pos. Sec. IV B 2 b

⳵t⌫˜ⁱ ␬⌫˜G2␣D˜^jD˜

jGⁱ „0,0,⫺(1/2)␬⌫˜G2⫾

冑

⫺k²⫹␬_⌫˜G2

2 (*2),long exp.… yes ␬_⌫^˜_G2⬎0 makes 2 Neg. 1 Pos.

⳵t⌫˜ⁱ ␬⌫˜G3␣D˜ⁱD˜_jG^j „0,0,⫺(1/2)␬⌫˜G3⫾

冑

⫺k²⫹␬_⌫^˜_G3

2 (*2),long exp.… yes ␬⌫˜G3⬎0 makes 2 Neg. 1 Pos.