Evolving black hole-neutron star binaries in general relativity using pseudospectral and finite difference methods

We present a code for solving the coupled Einstein-hydrodynamics equations to evolve relativistic, self-gravitating fluids. The Einstein field equations are solved in generalized harmonic coordinates on one grid using pseudospectral methods, while the fluids are evolved on another grid using shock-capturing finite difference or finite volume techniques. We show that the code accurately evolves equilibrium stars and accretion flows. Then we simulate an equal-mass nonspinning black hole-neutron star binary, evolving through the final four orbits of inspiral, through the merger, to the final stationary black hole. The gravitational waveform can be reliably extracted from the simulation.

Figure 1

A snapshot at $t=0.35$ of the density, pressure, and 3-velocity at each gridpoint on a line perpendicular to the shock interface. The lines show the exact solution. Circles, squares, and triangles show the numerical pressure, density, and velocity, respectively.

Figure 2

The deviation of the density $\rho$ from its exact equilibrium ${\rho }_{\mathrm{ex}}$. $L_{\mathrm{dev}}(f)$ is a relative L2 norm: $L_{\mathrm{dev}}(f)=L2(f-f_{\mathrm{ex}})/L2(f_{\mathrm{ex}})$, where the norms sum only over points outside the horizon. The top panel is a relativistic Bondi flow, the middle panel is a perturbed Bondi flow, and the bottom panel is a moving black hole flow. The deviations in the top and bottom panels are scaled for second-order convergence.

Figure 3

The deviations from known equilibrium values of the density $\rho$ and 4-metric ${\psi }_{\alpha \beta }$ for various grids. As in Fig. 2, $L_{\mathrm{dev}}(f)=L2(f-f_{\mathrm{ex}})/L2(f_{\mathrm{ex}})$. In the top panel, the metric is fixed and the fluid is evolved. In the other panels, the fluid is fixed and the metric is evolved. In the middle panel, the computational domain contains the stellar surface. In the bottom panel, the domain is entirely inside the star.

Figure 4

A neutron star with a realistic equation of state with both matter and metric evolved. Two runs are shown, corresponding to lower resolution PS and fluid grids and higher resolution grids. We plot the normalized L2 norm of the change in density and the normalized L2 norm of the violation in the generalized harmonic constraints.

Figure 5

The central density (relative to its initial value) for simulations of the four isolated stars. At the initial time, each star is subjected to a perturbation in the form of a 1% pressure depletion throughout the star. The inset shows the same thing but with a larger density range, so that later stages of the collapse of star B can be seen.

Figure 6

A representation of the choices of grid resolutions used in our reported runs. $N_F$ is the effective total number of grid points on the fluid finite difference grid, which includes “mirror” points given by the reflection symmetry about the equatorial plane; $N_{\mathrm{PS}}$ is the total number of collocation points on the spectral grid.

Figure 7

The normalized constraint energy $\parallel \mathcal{C}\parallel$ defined by Eq. (71) in 39. In this plot and all plots below, we use units in which the initial irreducible mass of the black hole is unity.

Figure 8

The baryonic density and polytropic constant of the neutron star, normalized to their initial values.

Figure 9

The irreducible mass of the black hole.

Figure 10

The proper separation between the black hole horizon and the neutron star center.

Figure 11

The mass and angular momentum of the binary. Solid lines are ADM surface integrals. Dashed lines show the expected changes due to losses by the observed gravitational waves. $M$, $J$, and $r{\psi }_4$ are extracted at radius $r=300$.

Figure 12

Snapshots of the density and horizon shape on the equatorial plane at times $t=1730$, $t=1850$, and $t=1970$. Shades of gray represent $\rho /{\rho }_{\max }$, the density relative to its current maximum. The black object is the apparent horizon. The density and horizon are shown in inertial coordinates. The distorted rectangle is the outer edge of the fluid grid, which is a fixed rectangle in moving coordinates but moves in the inertial frame.

Figure 13

The evolution of the neutron star and black hole masses from the beginning of the star’s disruption to the final stationary state. Here, $M_0$ is the rest mass of the neutron star matter outside the hole, $M_{\mathrm{irr}}$ is the black hole’s irreducible mass, and $M_{\mathrm{c}}$ is the Christodoulou mass of the black hole, defined as $\sqrt{M_{\mathrm{irr}}^2+J^2/(4M_{\mathrm{irr}}^2)}$, where $J$ is the black hole’s spin.

Figure 14

The trajectory in inertial coordinates of the neutron star center from $t=0$ to $t=t_{\mathrm{disr}}$, and that of the black hole horizon center from $t=0$ to $t=3000$. The black hole line switches from dashed to solid at $t=t_{\mathrm{disr}}$. Both the black hole and the neutron star remain on the equatorial plane throughout the evolution.

Figure 15

The gravitational wave signal for the entire inspiral and merger (excluding the initial burst of junk radiation). The quadrupole contribution dominates, so we plot only the $l=m=2$ amplitudes.

Figure 16

The constraint energy $\parallel \mathcal{C}\parallel$, black hole mass $M_{\mathrm{c}}$, and baryonic rest mass $M_0$ for the merger phase in evolutions carried out at two different resolutions.