Simulating Binary Neutron Star Mergers with Finite-temperature Equations of State: The influences of the slope of the symmetry energy and artificial heating

TL;DR

This work investigates how the slope of the symmetry energy, $L_{ m sym}$, in finite-temperature nuclear EOSs affects binary neutron star mergers. Using three nucleonic EOSs (DDLS30, DDLS50, DDLS70) with fixed $Q_{ m sat}$ and neutrino-transport-enabled NR simulations, the authors quantify a thermally induced modification to the tidal deformability, obtaining a thermally corrected parameter $ m ilde\Lambda_{sim}$ that can differ notably from the cold value. They show that GW signals are largely governed by macroscopic parameters and that the post-merger peak frequency $f_2$ is only weakly sensitive to $L_{ m sym}$, whereas the ejecta composition and resulting $r$-process yields exhibit stronger $L_{ m sym}$-dependence, suggesting EM counterparts as a more promising probe. A key methodological advance is the construction of $ m ilde\\Lambda_{sim}$ by mapping the dynamical state to a 1D EOS and solving TOV, which improves agreement with GW models for the HH scenario, particularly for DDLS30 where thermal effects are strongest. Overall, the paper highlights the potential of electromagnetic observations, in combination with thermally informed GW modeling, to constrain $L_{ m sym}$, while acknowledging the challenges in establishing a direct GW-to-$L_{ m sym}$ mapping.

Abstract

We present a new set of numerical-relativity simulations of merging binary neutron stars, aiming to identify possible observable signatures of the slope of the symmetry energy $L_{\rm sym}$. To achieve this goal, we employ a set of equations of state based on a parameterization of the covariant density functional theory of nuclear matter that allows controlled variations of $L_{\rm sym}$ and the skewness $Q_{\rm sat}$, holding the latter fixed. For a set of our simulations, we identify a steep energy gradient in the equation of state at subsaturation densities, which acts as a source of heating with subsequent stiffening produced by thermal support. Accounting for related structural modifications in the tidal deformability reconciles our results with theoretical expectations. On the other hand, we show that gravitational waves are unlikely to distinguish the role of $L_{\rm sym}$. In contrast to this, we find that the ejecta composition is significantly altered in our simulations, which employ an M1 moment scheme, when $L_{\rm sym}$ is varied. Based on our extracted dynamical ejecta properties, we compute r-process yields and find that they are distinct for the different $L_{\rm sym}$, especially at lower mass numbers $A \lesssim 120$. This suggests that electromagnetic counterparts are more likely to exhibit signatures; however, a direct connection to $L_{\rm sym}$ remains a challenge, given the complex interplay between details of the ejecta properties and the kilonova signal.

Figures (10)

• Figure 1: Equations-of-state for cold, $\beta$-equilibrated models (left panel), employed for the construction of the initial data. Mass-radius (middle panel) and tidal deformability-mass (right panel) diagrams for the simulated EOSs. Markers correspond to our simulated configurations.
• Figure 2: Specific internal energy per baryon $\epsilon$ as a function of scaled density $n_b/n_{\rm sat}$ and temperature $T$ for the ${\rm DDLS}30$ (upper panel), where a clear discontinuity can be seen, and ${\rm DDLS}50$ (lower panel) EOSs at a typical $Y_e = 0.04$, found in outer layers of a NS. The markers labeled as $\epsilon_i$ represent an initial point at temperature $T_i = 0.1~{\rm MeV}$ (e.g., found in the beginning of a simulation) and $n_i/n_{\rm sat}=0.287$, while markers labeled as $\epsilon_f$ are defined with same internal energy as the initial one, but at the immediate predecessor tabulated density $n_f/n_{\rm sat}=0.262$.
• Figure 3: Compositional profile of the NSs in the initial data along the $x$-axis, where $x=0$ coincides with the center of the NSs, while transition to the artificial atmosphere is seen in the increase of $Y_e$.
• Figure 4: Evolution of the central rest-mass density (upper panel), maximum temperature (middle panel) and ratio of mass-averaged thermal-to-total pressure (lower panel) of matter elements with $n_b/n_{\rm sat} = [0.1,2.0]$, roughly corresponding to the the hot annulus formed in the post-merger. A moving average filter with a $0.2~{\rm ms}$ time window is employed.
• Figure 5: Temperature (left half) and electron fraction (right half) of the remnants for simulations DDLS$30$ (left panel), DDLS$50$ (middle panel) and DDLS$70$ (right panel), on the $x-y$ plane (bottom panels) and $x-z$ plane (top panels) at the end of the simulations $t-t_{\rm mrg}\approx 52~{\rm ms}$. Contour lines mark isodensity surfaces with $n_b/n_{\rm sat} = [10^{-3}, 10^{-2}, 10^{-1}, 10^0, 2]$. The softer DDLS$30$ EOS produces the hottest annular region between one and two saturation densities, followed by DDLS$50$ and DDLS$70$, respectively, with the least electronized remnant (see text). Likewise, the low density, comparatively colder polar region of all remnants exhibits a sustained neutrino-wind with increased $Y_e$ for higher $L_{\rm sym}$.