Advances in Astronomy

Volume 2017, Article ID 6127031, 12 pages

https://doi.org/10.1155/2017/6127031

## Analysis of the Conformally Flat Approximation for Binary Neutron Star Initial Conditions

^{1}Center for Astrophysics, Department of Physics and Center for Research Computing, University of Notre Dame, Notre Dame, IN 46556, USA^{2}Center for Astrophysics, Department of Physics, University of Notre Dame, Notre Dame, IN 46556, USA^{3}Hanoi National University of Education, 136 Xuan Thuy, Hanoi, Vietnam^{4}Joint Institute for Nuclear Astrophysics (JINA), University of Notre Dame, Notre Dame, IN 46556, USA

Correspondence should be addressed to Grant J. Mathews; ude.dn@swehtamg

Received 19 September 2016; Accepted 30 November 2016; Published 9 January 2017

Academic Editor: Ignazio Licata

Copyright © 2017 In-Saeng Suh et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The spatially conformally flat approximation (CFA) is a viable method to deduce initial conditions for the subsequent evolution of binary neutron stars employing the full Einstein equations. Here we analyze the viability of the CFA for the general relativistic hydrodynamic initial conditions of binary neutron stars. We illustrate the stability of the conformally flat condition on the hydrodynamics by numerically evolving ~100 quasicircular orbits. We illustrate the use of this approximation for orbiting neutron stars in the quasicircular orbit approximation to demonstrate the equation of state dependence of these initial conditions and how they might affect the emergent gravitational wave frequency as the stars approach the innermost stable circular orbit.

#### 1. Introduction

The epoch of gravitational wave astronomy has now begun with the first detection [1, 2] of the merger of binary black holes by Advanced LIGO [3]. Now that the first ground based gravitational wave detection has been achieved, observations of binary neutron star mergers should soon be forthcoming. This is particularly true as other second generation observatories such as Advanced VIRGO [4] and KAGRA [5] will soon be online. In addition to binary black holes, neutron star binaries are thought to be among the best candidate sources gravitational radiation [6, 7]. The number of such systems detectable by Advanced LIGO is estimated [7–14] to be of order several events per year based upon observed close binary-pulsar systems [15, 16]. There is a difference between neutron star mergers and black hole mergers; however, in that neutron star mergers involve the complex evolution of the matter hydrodynamic equations in addition to the strong gravitational field equations. Hence, one must carefully consider both the hydrodynamic and field evolution of these systems.

To date there have been numerous attempts to calculate theoretical templates for gravitational waves from compact binaries based upon numerical and/or analytic approaches (see, e.g., [17–26]). However, most approaches utilize a combination of post-Newtonian (PN) techniques supplemented with quasicircular orbit calculations and then applying full GR for only the last few orbits before disruption. In this paper we analyze the hydrodynamic evolution in the spatially conformally flat metric approximation (CFA) as a means to compute stable initial conditions beyond the range of validity of the PN regime, that is, near the last stable orbits. We establish the numerical stability of this approach based upon many orbit simulations of quasicircular orbits. We show that one must follow the stars for several orbits before a stable quasicircular orbit can be achieved. We also illustrate the equation of state (EoS) dependence of the initial conditions and associated gravitational wave emission.

When binary neutron stars are well separated, the post-Newtonian (PN) approximation is sufficiently accurate [27]. In the PN scheme, the stars are often treated as point masses, either with or without spin. At third order, for example, it has been estimated [28–30] that the error due to assuming the stars are point masses is less than one orbital rotation [28] over the ~16,000 cycles that pass through the LIGO detector frequency band [7]. Nevertheless, it has been noted in many works [25, 31–42] that relativistic hydrodynamic effects might be evident even at the separation (~10–100km) relevant to the LIGO window.

Indeed, the templates generated by PN approximations, unless carried out to fifth and sixth order [28, 29], may not be accurate unless the finite size and proper fluid motion of the stars are taken into account. In essence, the signal emitted during the last phases of inspiral depends upon the finite size and the equation of state (EoS) through the tidal deformation of the neutron stars and the cut-off frequency when tidal disruption occurs.

Numeric and analytic simulations [43–51] of binary neutron stars have explored the approach to the innermost stable circular orbit (ISCO). While these simulations represent some of the most accurate to date, simulations generally follow the evolution for a handful of orbits and are based upon initial conditions of quasicircular orbits obtained in the conformally flat approximation. Accurate templates of gravitational radiation require the ability to stably and reliably calculate the orbit initial conditions. The CFA provides a means to obtain accurate initial conditions near the ISCO.

The spatially conformally flat approximation to GR was first developed in detail in [32]. That original formulation, however, contained a mathematical error first pointed out by Flanagan [52] and subsequently corrected in [34]. This error in the solution to the shift vector led to a spurious NS crushing prior to merger. The formalism discussed below is for the corrected equations. Here, we discuss the hydrodynamic solutions as developed in [31–34, 53, 54]. This CFA formalism includes much of the nonlinearity inherent in GR and leads set of coupled, nonlinear, elliptic field equations that can be evolved stably. We also note that an alternative spectral method solution to the CFA configurations was developed by [55, 56], and approaches beyond the CFA have also been proposed [48]. However, our purpose here is to clarify the viability of the hydrodynamic solution without the imposition of a Killing vector or special symmetry. This approach is the most adaptable, for example, to general initial conditions such as that of arbitrarily elliptical orbits and/or arbitrarily spinning neutron stars.

Here, we summarize the original CFA approach and associated general relativistic hydrodynamics formalism developed in [32, 34, 53, 54] and illustrate that it can produce stable initial conditions anywhere between the post-Newtonian to ISCO regimes. We quantify how long this method takes to converge to quasiequilibrium and demonstrate the stability by subsequently integrating up to ~100 orbits for a binary neutron star system. We also analyze the EoS dependence of these quasicircular initial orbits and show how these orbits can be used to make preliminary estimates [57] of the gravitational wave signal for the initial conditions.

This paper is organized as follows. In Section 2 the basic method is summarized and in Section 3 a number of code tests are performed in the quasiequilibrium circular orbit limit to demonstrate the stability of the technique. The EoS dependence of the initial conditions and associated gravitational wave frequency are analyzed in Section 4. Conclusions are presented in Section 5.

#### 2. Method

##### 2.1. Field Equations

The solution of the field equations and hydrodynamic equations of motion were first solved in three spatial dimensions and explained in detail in the 1990s in [31, 32] and subsequently further reviewed in [53, 58]. Here, we present a brief summary to introduce the variables relevant to the present discussion.

One starts with the slicing of space-time into the usual one-parameter family of hypersurfaces separated by differential displacements in a time-like coordinate as defined in the () ADM formalism [59, 60].

In Cartesian , , isotropic coordinates, proper distance is expressed aswhere the lapse function describes the differential lapse of proper time between two hypersurfaces. The quantity is the shift vector denoting the shift in space-like coordinates between hypersurfaces. The curvature of the metric of the 3-geometry is described by a position-dependent conformal factor times a flat-space Kronecker delta (). This conformally flat condition (together with the maximal slicing gauge, ) requires [60]where is the extrinsic curvature tensor and are 3-space covariant derivatives. This conformally flat condition on the metric provides a numerically valid initial solution to the Einstein equations. The vanishing of the Weyl tensor for a stationary system in three spatial dimensions guarantees that a conformally flat solution to the Einstein equations exists.

One consequence of this conformally flat approximation to the three-metric is that the emission of gravitational radiation is not explicitly evolved. Nevertheless, one can extract the gravitational radiation signal and the back reaction via a multipole expansion [32, 61]. An application to the determination of the gravitational wave emission from the quasicircular orbits computed here is given in [57]. The advantage of this approximation is that conformal flatness stabilizes and simplifies the solution to the field equations.

As a third gauge condition, one can choose separate coordinate transformations for the shift vector and the hydrodynamic grid velocity to separately minimize the field and matter motion with respect to the coordinates. This set of gauge conditions is key to the present application. It allows one to stably evolve up to hundreds and even thousands of binary orbits without the numerical error associated with the frequent advocating of fluid through the grid.

Given a distribution of mass and momentum on some manifold, then one first solves the constraint equations of general relativity at each time for a fixed distribution of matter. One then evolves the hydrodynamic equations to the next time step. Thus, at each time slice a solution to the relativistic field equations and information on the hydrodynamic evolution is obtained.

The solutions for the field variables , , and reduce to simple Poisson-like equations in flat space. The Hamiltonian constraint [60] is used to solve for the conformal factor [32, 62]In the Newtonian limit, the RHS is dominated [32] by the proper matter density , but in strong fields and compact neutron stars there are also contributions from the internal energy density , pressure , and extrinsic curvature. The source is also significantly enhanced by the generalized curved-space Lorentz factor where is the time component of the relativistic four velocity and are the covariant spatial components. This factor, , becomes important near the last stable orbit as the specific kinetic energy of the stars rapidly increases.

In a similar manner [32, 62], the Hamiltonian constraint, together with the maximal slicing condition, provides an equation for the lapse function,

Finally, the momentum constraints yields [60] an elliptic equation for the shift vector [34, 52],whereHere are the spatial components of the momentum density one-form as defined below.

We note that, in early applications of this approach, the source for the shift vector contained a spurious term due to an incorrect transformation between contravariant and covariant forms of the momentum density as was pointed out in [34, 52]. As illustrated in those papers, this was the main reason why early hydrodynamic calculations induced a controversial additional compression on stars causing them to collapse to black holes prior to inspiral [31]. This problem no longer exists in the formulation summarized here.

##### 2.2. Relativistic Hydrodynamics

To solve for the fluid motion of the system in curved-space time it is convenient to use an Eulerian fluid description [63]. One begins with the perfect fluid stress-energy tensor in the Eulerian observer rest frame, where is the relativistic four velocity one-form.

By introducing the usual set of Lorentz contracted state variables it is possible to write the relativistic hydrodynamic equations in a form which is reminiscent of their Newtonian counterparts [63]. The hydrodynamic state variables are the coordinate baryon mass density the coordinate internal energy density the spatial three velocityand the covariant momentum density

In terms of these state variables, the hydrodynamic equations in the CFA are as follows: the equation for the conservation of baryon number takes the form The equation for internal energy evolution becomes Momentum conservation takes the formwhere the last term in (15) is the contribution from the radiation reaction potential as defined in [32, 57]. In the construction of quasistable orbit initial conditions, this term is set to zero. Including this term would allow for a calculation of the orbital evolution via gravitational wave emission in the CFA. However, there is no guarantee that this is a sufficiently accurate solution to the exact Einstein equations. Hence, the CFA is useful for the construction of initial conditions.

##### 2.3. Angular Momentum and Orbital Frequency

In the quasicircular orbit approximation (neglecting angular momentum in the radiation field), this system has a Killing vector corresponding to rotation in the orbital plane. Hence, for these calculations the angular momentum is well defined and given by an integral over the space-time components of the stress-energy tensor [64]; that is, Aligning the -axis with the angular momentum vector then gives

To find the orbital frequency detected by a distant observer corresponding to a fixed angular momentum we employ a slightly modified derivation of the orbital frequency compared to that of [53]. In asymptotically flat coordinates the angular frequency detected by a distant observer is simply the coordinate angular velocity; that is, one evaluates

In the ADM conformally flat () curved space, our only task is then to deduce from code coordinates. For this we make a simple polar coordinate transformation keeping our conformally flat coordinates, so Now, the code uses covariant four velocities, . This gives . Finally, one must density weight and volume average over the fluid differential volume elements. This givesThis form differs slightly from that of [53] but leads to the similar results.

A key additional ingredient, however, is the implementation of a grid three velocity that minimizes the matter motion with respect to and . Hence, the total angular frequency to a distant observer , and in the limit of rigid corotation, , where .

For the orbit calculations illustrated here we model corotating stars, that is, no spin in the corotating frame. This minimizes matter motion on the grid. However, we note that there is need at the present time of initial conditions for arbitrarily spinning neutron stars and the method described here is equally capable of supplying those initial conditions.

As a practical approach the simulation [32] of initial conditions is best run first with viscous damping in the hydrodynamics for sufficiently long time (a few thousand cycles) to relax the stars to a steady state. Subsequently, one can run with no damping. In the present illustration we examine stars at large separation which are in quasiequilibrium circular orbits and stable hydrodynamic configurations. In the initial relaxation phase the orbits are circularized by damping any radial velocity components. During the evolution, the rigorous conservation of angular momentum is imposed by adjusting the orbital angular velocity in (20) such that (17) remains constant. The simulated orbits described in this work span the time from the last several minutes up to orbit inspiral. Here, we illustrate the stability of the multiple orbit hydrodynamic simulation and examine where the initial conditions for the strong field orbit dynamics computed here deviates from the post-Newtonian regime.

#### 3. Code Validation

##### 3.1. Code Tests

To evolve stars at large separation distance it is best [53] to decompose the grid into a high resolution domain with a fine matter grid around the stars and a coarser domain with an extended grid for the fields. Figure 1 shows a schematic of this decomposition from [54].