Numerical inside view of hypermassive remnant models for GW170817

The first multimessenger observation attributed to a merging neutron star binary provided an enormous amount of observational data. Unlocking the full potential of this data requires a better understanding of the merger process and the early postmerger phase, which are crucial for the later evolution that eventually leads to observable counterparts. In this work, we perform standard hydrodynamical numerical simulations of a system compatible with GW170817. We focus on a single equation of state and two mass ratios, while neglecting magnetic fields and neutrino radiation. We then apply newly developed postprocessing and visualization techniques to the results obtained for this basic setting. The focus lies on understanding the three-dimensional structure of the remnant, most notably the fluid flow pattern, and its evolution until collapse. We investigate the evolution of mass and angular momentum distribution up to collapse, as well as the differential rotation along and perpendicular to the equatorial plane. For the cases that we studied, the remnant cannot be adequately modeled as a differentially rotating axisymetric neutron star. Further, the dominant aspect leading to collapse is the gravitational wave radiation and not internal redistribution of angular momentum. We relate features of the gravitational wave signal to the evolution of the merger remnant and make the waveforms publicly available. Finally, we find that the three-dimensional vorticity field inside the disk is dominated by medium-scale disturbances and not the orbital velocity, with potential consequences for magnetic field amplification effects.

Figure 1

Overview of the evolution phases for the $q=0.9$ case. Time runs from lower right to upper left, while the other two dimensions correspond to the orbital plane. The time coordinate was compressed by a factor 0.05 with respect to the spatial coordinates in geometric units, such that light cones would appear almost orthogonal to the world tube of the remnant. The transparent green surface corresponds to a fixed density of $5\times {10}^{13}\mathrm{g}/{\mathrm{cm}}^3$, highlighting the evolution of the merged NS and the coalescing NS shortly before merger. The solid red surface corresponds to a density of ${10}^{11}\mathrm{g}/{\mathrm{cm}}^3$, as a proxy for the denser parts of the disk. To avoid occlusion, one half-plane was cut away. The blue surface marks the apparent horizon extracted during the simulation. The glowing region serves as a rough sketch of the shock heating during merger. It is the result of a volume rendering of an optically transparent black body light source with temperature proportional to the remnant temperature (shifted into the visible spectrum) and intensity proportional to mass density.

Figure 2

Postmerger evolution of mass density in the remnant for the equal mass case (top panel) and unequal mass case (bottom panel). The solid curve shows the maximum baryonic mass density in the fluid frame. The dotted curve shows an average density given by bulk mass per bulk volume (see Sec. 2d). The time refers to coordinate time. For comparison, the horizontal lines show the central density (dashed) and average bulk density (dash dotted) of the maximum-mass spherical NS solution (we note that the initial agreement with the maximum baryonic mass density is a pure coincidence). The vertical lines mark the retarded time of merger and the formation of an apparent horizon.

Figure 3

Mass distribution in the meridional plane 2 ms before apparent horizon formation, averaged over a time window $\pm 1\mathrm{ms}$. The color scale shows the baryonic mass density. The contours mark densities for which the corresponding isodensity surfaces contain selected mass fractions. We show the contours for 99%, 98%, the mass swallowed by the BH within 1 ms after formation, the maximum mass of uniformly rotating (“Max Unif.”) and nonrotating NS (“Max TOV”). The contour labeled “core” refers to the mass of the nonrotating NS best approximating the remnant core, which is defined in Sec. 3c and almost identical to the maximum nonrotating NS mass.

Figure 4

Like Fig. 3, but showing the time 2 ms after apparent horizon formation. The baryonic mass within the plotted density contours (excluding the BH interior) is given in the label as fractions of the total baryonic mass of the initial data.

Figure 5

Radial distribution of unbound matter versus time after merger. The color scale corresponds to unbound mass per radial distance, where matter is considered unbound according to the geodesic criterion. The horizontal lines mark the time of BH formation.

Figure 6

The GW signals extracted from our simulations (blue) as observed face on at a distance of 40.7 Mpc. We extend the inspiral with the BNS model IMRPhenomD_NRTidalv2 [75] (orange), aligned with the numerical data over the grey shaded region. For comparison, we also include the BBH model IMRPhenomD [76, 77] (gray dashed lines).

Figure 7

Left panels: evolution of GW frequency (with respect to retarded time) $f_{\mathrm{GW}}$ in comparison to twice the maximum rotation rate $f_{\mathrm{rot}}^{\max }$. We also show the radial oscillation frequency, $f_{\text{radial}}$, and the rate of the modulation of the GW frequency, $f_{\text{vibrato}}$ (see Sec. 3c). Vertical lines mark time of merger and apparent horizon formation. For comparison, we show the inspiral GW frequency $f_{\text{model}}$ according to the IMRPhenomD_NRTidalv2 waveform model. Further, $f_{\mathrm{BH}}$ denotes the $l=m=2$, $n=0$ QNM frequency of Kerr BHs using mass and angular momentum extracted from the BH in the simulation, as function of time (see Sec. 2d). Horizontal lines mark the frequency at merger (maximum GW amplitude) and the main peak of the spectrum. Right panels: Power spectrum of the $l=m=2$ component of the GW signal, at distance 40.7 Mpc, in terms of $F\tilde{h}(F)$, where ${\tilde{h}}^2(F)={\tilde{h}}_{+}^2(F)+{\tilde{h}}_{\times }^2(F)$. For comparison we show the design sensitivity curves for various detectors, taken from [78].

Figure 8

Total baryonic mass contained inside surfaces of constant density (blue curve) versus proper volume contained within the same surfaces. The top panel shows the remnant for the equal mass case, the bottom panel for the unequal mass system, both at a time 1 ms before apparent horizon formation. The dot marks bulk mass and bulk volume of the remnant (see Sec. 2d). For comparison, we show bulk mass versus bulk volume (green line) for the sequence of nonrotating NS following the same EOS as the BNS initial data, starting at mass $0.9M_{\odot }$ up to the maximum bulk mass. The intersection with the remnant mass-volume curve defines the core equivalent TOV model, for which we show the mass-volume relation as well (dashed red curve).

Figure 9

Evolution of masses for the equal mass case (top panel) and the unequal mass case (bottom panel). The dashed red curve shows the bulk baryonic mass, the solid blue curve the bulk baryonic mass of the core TOV equivalent (see Sec. 3c), and the horizontal dash-dotted line the maximum bulk baryonic mass of nonrotating NS following the same EOS as the BNS initial data (SFHO). The solid black curve shows the gravitational mass of the BH (Sec. 2d) and the vertical line, the time of first apparent horizon detection. The solid olive curve shows the total baryonic mass present in the computational domain up to BH formation and a constant afterwards. The dashed black curve shows the baryonic mass swallowed by the BH. The latter is computed by first adding the mass inside the domain, excluding the interior of apparent horizons, to the cumulative mass loss by ejecta, and then computing the difference to the time directly before BH formation. The green curve shows the ADM mass, where the shaded area is the contribution attributed to the energy of GW radiation within the computational domain (see Sec. 2d).

Figure 10

Evolution of radial mass distribution in the merger remnant for the equal mass case (top panel) and unequal mass case (bottom panel). The solid black curves show the time evolution of the proper volume within isodensity surfaces containing fractions 0.1,0.5,0.9,0.97, and 0.98 of the total baryonic mass. The volume is given in terms of volumetric radius $R_V$, the radius of an Euclidean sphere of equal volume. We further compare the remnant to the spherical NS solution of maximum bulk mass. For both, we compute the density as function of proper volume within isodensity surfaces. The color scale represents the difference of density at given volumetric radius, normalized to the maximum density of the spherical NS. The volumetric surface radius of the latter is shown as dashed vertical line.

Figure 11

Quasistationary part of the remnant structure, in the frame corotating with $m=2$ density perturbation of remnant. The visualization represents the average over a time window $7\pm 1\mathrm{ms}$ after merger (2 ms before collapse) for the $q=0.9$ case. The inner and outer surfaces (cut open) mark mass densities 0.3 and 0.01 of the maximum density, respectively. The wires represent integral lines of the averaged velocity field, shown inside the dense region enclosed within the outer surface. The top left rendering shows the remnant from a perspective along the rotation axis, the bottom left one from the side, along the longer axes of the deformed core, looking onto a meridional plane that cuts through the secondary vortices.

Figure 12

Remnant structure visualized as in Fig. 11, but for three different time windows. From left to right: $3\pm 1\mathrm{ms}$ after merger, $5\pm 1\mathrm{ms}$ after merger, and $7\pm 1\mathrm{ms}$ after merger (2 ms before collapse). The camera distance remains constant to allow size comparison.

Figure 13

Rotation profile of the merger remnant 1 ms before collapse, for mass ratios $q=1$ and $q=0.9$, in the equatorial plane, versus circumferential radius. The filled curves show the angular velocity as seen from infinity and the frame dragging contribution (see Sec. 2d). The data have been averaged in $\varphi$ direction and over a time window $\pm 0.5\mathrm{ms}$. For comparison, we show the angular velocity of test particles in prograde circular orbit (dashed curve), and the same for particles orbiting the spinning BH formed after collapse (solid green curve). The solid black curves show the rotation rate and surface radius of uniformly rotating supramassive NS, either at mass shedding limit (right curve) or at minimal angular momentum (left curve). The vertical line marks the radius where the $\varphi$-averaged mass density in the equatorial plane falls below 5% of the central value. The horizontal line marks half of the GW angular frequency.

Figure 14

Two-dimensional differential rotation profile of the remnant for mass ratio $q=0.9$, averaged over the time window $7\pm 1\mathrm{ms}$ after merger. The color scale denotes the rotational frequency in the coordinate system corotating with the perturbation pattern. Negative values signify that (in the inertial frame) the fluid is rotating more slowly than the pattern. The left panel shows the average in $\varphi$ direction versus cylindrical radius and z coordinate. The right panel shows a cut in $(y,z)$ plane, the $y$ axis being roughly aligned with the secondary vortices shown in Fig. 11. The middle panel shows the $(x,z)$ plane instead. The solid curves mark isodensity contours (in the left panel with respect to $\varphi$-averaged density). Radial distances in the equatorial plane and distances along the z axis are both proper distances, allowing us to asses the oblateness of the remnant.

Figure 15

Vorticity lines inside the disk 2 ms before (left) and 5 ms after BH formation (right), for the $q=0.9$ case. The black surface in the right panel marks the apparent horizon. The cut-open solid surfaces in the left panel mark densities $4.0\times {10}^{10}$, $2.2\times {10}^{11}$, $1.7\times {10}^{13}$, and $5.0\times {10}^{14}\mathrm{g}/{\mathrm{cm}}^3$, in the right panel $1.5\times {10}^{10}$ and $8.1\times {10}^{10}\mathrm{g}/{\mathrm{cm}}^3$. Vortex lines are cut off outside the outermost surface and beyond a length cutoff in order to limit cluttering. The color indicates the vorticity magnitude, where lighter color corresponds to larger values. The camera distance is the same in both panels.

Figure 16

Evolution of energy and angular momentum for the equal mass system (top panel) and the $q=0.9$ case (bottom panel). The solid blue curve shows total ADM energy versus angular momentum, and the dashed black line a post-Newtonian approximation (see Sec. 3g). The blue circles represent energy and angular momentum excluding the amount attributed to GW radiation outside $r=100\mathrm{km}$, at times 1 ms after merger and 1 ms before collapse. The black and red dots mark energy and angular momentum of the BH at the end of the simulation and 1 ms after apparent horizon formation. The diagonally oriented red curves show the contributions within spheres of constant coordinate radius to the ADM volume integrals, at regular intervals from 1 ms after merger to 1 ms before apparent horizon formation. The horizontally oriented red curves show the time evolution of the contributions within spheres containing different fixed amounts of baryonic mass, ranging up to the baryonic mass swallowed by the BH within 1 ms after formation. Time increases monotonically from right to left. The shaded region is bounded by the mass shedding limit and smallest possible angular momentum of stable uniformly rotating NSs. The green horizontal line marks the maximum mass of nonrotating NS. The dotted line shows the angular momentum of extremal Kerr BHs.

Figure 17

Distribution of angular momentum at time 1 ms before collapse for the equal mass system (top panel) and the $q=0.9$ case (bottom panel). The solid red curve shows the contribution to the ADM angular momentum within isodensity surfaces as function of the baryonic mass within the same surfaces. The green dot marks the position of the remnant bulk (see Sec. 2d). The solid green line shows the same for spherical surfaces (spherical with respect to simulation coordinates). It is hidden behind the dashed yellow curve, which shows an estimate of the Komar angular momentum (see Sec. 2d). The shaded area shows the values possible for uniformly rotating NS.