Entropy-limited higher-order central scheme for neutron star merger simulations

Numerical relativity simulations are the only way to calculate exact gravitational waveforms from binary neutron star mergers and to design templates for gravitational-wave astronomy. The accuracy of these numerical calculations is critical in quantifying tidal effects near merger that are currently one of the main sources of uncertainty in merger waveforms. In this work, we explore the use of an entropy-based flux-limiting scheme for high-order, convergent simulations of neutron star spacetimes. The scheme effectively tracks the stellar surface and physical shocks using the residual of the entropy equation thus allowing the use of unlimited central flux schemes in regions of smooth flow. We perform the first neutron star merger simulations with such a method and demonstrate up to fourth-order convergence in the gravitational waveform phase. The scheme reduces the phase error up to a factor 5 when compared to state-of-the-art high-order characteristic schemes and can be employed for producing faithful tidal waveforms for gravitational-wave modeling.

Figure 1

Simple wave. Numerical solution of the rest-mass density with $n=800$ grid points and CFL factor 0.125 for the reconstruction schemes CENO3 (using purple cross) and WENOZ (using green circle). The initial profile and the exact solution are also included using dashed blue and solid red lines, respectively.

Figure 2

Profiles of the rest-mass density (green filled circle), velocity (blue filled circle) and pressure (red filled circle) for the special-relativistic Sod test at $t=0.6$. Top: WENO5 reconstruction. Bottom: GODUNOV reconstruction. The solution is computed on a grid of 1600 points with resolution $\mathrm{\Delta }x=1.25\times {10}^{-3}$. Solid black lines are the exact solutions.

Figure 3

Profiles of the rest-mass density (green filled circle), velocity (blue filled circle) and pressure (red filled circle) for the special-relativistic blast wave 1 test [71] at $t=0.4$. Top: WENO5 reconstruction. Bottom: GODUNOV reconstruction. The solution is computed on a grid of 800 points with resolution $\mathrm{\Delta }x=1.25\times {10}^{-3}$. Solid black lines are the exact solutions.

Figure 4

Profiles of the rest-mass density (green filled circle), velocity (blue filled circle) and pressure (red filled circle) for the special-relativistic blast wave 2 test [71] at $t=0.4$. Top: WENO5 reconstruction. Bottom: GODUNOV reconstruction. The solution is computed on a grid of 800 points with resolution $\mathrm{\Delta }x=1.25\times {10}^{-3}$. Solid black lines are the exact solutions.

Figure 5

Two-dimensional profiles of the rest-mass density $\rho$ (right half of the plots) and entropy production function $\nu$ (left half of the plots) for a static Tolmann-Oppenheimer-Volkoff (TOV) star in a dynamical spacetime with $\mathrm{\Gamma }$-law EOS at $t=1000$. Top (left to right): CENO3, WENO5 and WENOZ reconstruction schemes with $n=128$. Bottom (left to right): The WENO5 scheme with increasing resolution of 64, 96 and 128 grid points, respectively. The other two schemes show similar behavior. Dashed lines denote the surface of the NS at $t=0$.

Figure 6

The TOV star in a dynamical spacetime with $\mathrm{\Gamma }$-law EOS. Top: central rest-mass density evolution for simulations of a single spherical star with $n=96$ for different EFL reconstruction schemes. Bottom: one-dimensional rest-mass density profiles in the $x$ direction at time $t=1000$ with $n=96$. The profiles in the y and z direction are identical.

Figure 7

Evolution of the $L_1$ distance $\parallel \rho (t)-\rho (0){\parallel }_1$ for a static TOV star in a dynamically evolved spacetime with a $\mathrm{\Gamma }$-law EOS. All three EFL reconstruction schemes are presented. Dashed lines show results scaled to second order.

Figure 8

Central rest-mass density evolution for simulations of a stationary RNS in a dynamical spacetime with $\mathrm{\Gamma }$-law EOS and $n=96$. Different EFL reconstruction schemes are shown.

Figure 9

Velocity profile of a stationary rotating neutron star in a dynamical spacetime with $\mathrm{\Gamma }$-law EOS. Top: one-dimensional profile of the velocity component ${\upsilon }_y$ along the $x$ direction at time $t=1000$ (four periods) with $n=128$. Bottom: the ${\upsilon }_y$ profile of the WENO5 scheme with increasing resolution.

Figure 10

Evolution of the $L_1$ distance $\parallel \rho (t)-\rho (0){\parallel }_1$ for a rotating neutron star in a dynamically evolved spacetime with $\mathrm{\Gamma }$-law EOS. Dashed lines show results scaled to second order.

Figure 11

Time evolution of the proper distance between the NSs for the ten-orbit bam:97 simulation. Different grid resolutions are presented for the EFL-WENOZ method (solid lines) and HO-LLF-WENOZ method (dashed lines).

Figure 12

Two-dimensional hybrid plots depicting the entropy production and rest-mass density profiles across different stages of the bam:27 simulation. Notice that bam:27 uses a $\mathrm{\Gamma }$-law EOS. Left to right: inspiral, merger, postmerger.

Figure 13

Two-dimensional hybrid plots depicting the entropy production and rest-mass density profiles across different stages of the bam:100 simulation. Notice that bam:100 uses a SLy EOS. Left to right: inspiral, merger, postmerger.

Figure 14

Conservation and convergence of the total rest mass on the refinement level $l=2$. Top: three-orbit bam:100 simulation. All three reconstruction schemes are presented for different grid resolutions. Dotted lines show results scaled to third order. Bottom: ten-orbit bam:97 simulation. The WENOZ reconstruction scheme is presented for different grid resolutions. Dotted lines show results scaled to fourth order.

Figure 15

Evolution and convergence of the central rest-mass density on the refinement level $l=1$ for the ten-orbit bam:97 simulation. Dotted lines show results scaled to fifth order.

Figure 16

$L_2$ norm of the Hamiltonian constraint on the refinement level $l=1$ for the BNS simulations. Top: three-orbit bam:100 simulation. All three reconstruction schemes are shown. Bottom: ten-orbit bam:97 simulation. Dotted lines show results scaled to second order.

Figure 17

Amplitude (blue) and real part (green) of strain $R{h}_{22}$ as well as instantaneous frequency $M{\omega }_{22}$ (red) of the GW signal obtained from the ten-orbit bam:97 simulation using the EFL-WENOZ scheme. The time of merger is defined as the first peak in the amplitude $A_{22}$ and is indicated by a black vertical line at $u_{\mathrm{mrg}}\approx 2400M$.

Figure 18

GW phase difference convergence rate study for the three-orbit bam:100 simulation: Solid lines represent phase differences between runs with consecutively increasing resolutions; dashed and dotted lines correspond to rescaled differences of HIG-FIN and MIG-HIG differences, respectively, where the scaling factor is computed from (24) using an integer convergence rate $p$; the latter was determined experimentally by looking for the best visual overlap between dashed, dotted and solid lines of the same color (which we failed to do with the CENO3 and WENO5 reconstructions); the green shaded regions indicate the obtained error budget as discussed in the text; the gray shaded regions mark the differences in merger times between runs with consecutively increasing resolutions.

Figure 19

GW phase difference convergence rate study for ten-orbit bam:97 simulation. Left panel: EFL-WENOZ method. Right panel: HO-LLF-WENOZ method.

Figure 20

GW phase difference convergence rate study for higher modes for ten-orbit bam:97 simulation using the entropy flux limiter with WENOZ reconstruction. Left panel: (3, 2) mode. Right panel: (4, 4) mode.

Figure 21

Faithfulness as a function of the resolution for the bam:97 simulation.