https://orcid.org/0000-0001-6982-1008
ABSTRACT
Massive neutron stars (NSs) with lifetimes of at least several seconds are expected to be the result of a sizable fraction of NS mergers. We study their formation using a large set of numerical relativity simulations. We show that they are initially endowed with angular momentum that significantly exceeds the mass-shedding limit for rigidly rotating equilibria. We find that gravitational wave (GW) emission is not able to remove this excess angular momentum within the time over which solid body rotation should be achieved. Instead, we argue that the excess angular momentum could be carried away by massive winds. Long-lived merger remnants are also formed with larger gravitational masses than those of rigidly rotating NSs having the same number of baryons. The excess mass is likely radiated in the form of neutrinos. The evolution of long-lived remnants on the viscous time-scale is thus determined by the interplay of finite-temperature effects, mass ejection, and neutrinos with potentially dramatic consequences for the remnants’ properties and stability. We also provide an empirical fit for the spin of the remnant at the end of its viscous evolution as a function of its final mass, and we discuss the implications for the magnetar model of short gamma-ray bursts (SGRBs). Finally, we investigate the possible electromagnetic signatures associated with the viscous ejecta. Massive outflows possibly resulting from the formation of long-lived remnants would power unusually bright, blue kilonova counterparts to GW events and SGRBs whose detection would provide smoking gun evidence for the formation of long-lived remnants.
1 INTRODUCTION
The outcome of neutron star (NS) mergers depends on the total mass of the system and on the poorly known equations of state (EOS) of dense nuclear matter (Shibata 2016, and references therein). Binaries with mass larger than ∼1.3–1.7 times the maximum mass for a non-rotating NS result in prompt black hole (BH) formation (Hotokezaka et al. 2011; Bauswein, Baumgarte & Janka 2013). Binaries with lower masses, but above the maximum mass of isolated rigidly rotating NSs, result in the formation of hypermassive neutron stars (HMNSs) temporarily supported against gravitational collapse by the large differential rotation (Baumgarte, Shapiro & Shibata 2000; Rosswog & Davies 2003; Shibata & Taniguchi 2006; Baiotti, Giacomazzo & Rezzolla 2008; Sekiguchi et al. 2011; Palenzuela et al. 2015; Bernuzzi et al. 2016). Even lower mass systems produce NS remnants that are stable on the spin-down time-scale (seconds to hours), called supramassive NSs (SMNSs), or stable massive NSs (MNSs) if their mass is below the maximum mass of a non-rotating NS (Giacomazzo & Perna 2013; Hotokezaka et al. 2013; Foucart et al. 2016a; Kastaun, Ciolfi & Giacomazzo 2016; Ciolfi et al. 2017; Kiuchi et al. 2018).
In the case of the binary NS merger GW170817 (Abbott et al. 2017a,b,c), the most likely outcome was an HMNS (Margalit & Metzger 2017; Shibata et al. 2017; Radice et al. 2018). However, the formation of a long-lived remnant for GW170817 is not completely ruled out (Yu & Dai 2017; Ai et al. 2018; Geng et al. 2018; Li et al. 2018). Indeed, the formation of SMNSs is expected to occur in a sizable fraction of mergers (Lasky et al. 2014; Gao, Zhang & Lü 2016; Piro, Giacomazzo & Perna 2017). This expectation has been recently reinforced by the discovery of two double NS systems with total gravitational masses as low as 2.5 M⊙ (Martinez et al. 2017; Stovall et al. 2018). Long-lived remnants have also been invoked to explain late time X-rays excess seen in some short gamma-ray bursts (SGRBs; Dai & Lu 1998a,b; Zhang & Meszaros 2001; Dai et al. 2006; Metzger, Quataert & Thompson 2008a; Rowlinson et al. 2010; Bucciantini et al. 2012; Rowlinson et al. 2013; Metzger & Piro 2014; Ciolfi & Siegel 2015; Lü et al. 2015; Margalit, Metzger & Beloborodov 2015; Rezzolla & Kumar 2015; Gao et al. 2016; Geng et al. 2016; Siegel & Ciolfi 2016a,b; Murase et al. 2018).
The identification of the outcome of the merger of binary NS systems with different masses would yield a precise measurement of the maximum mass of NSs (e.g. Lasky et al. 2014; Lawrence et al. 2015; Margalit & Metzger 2017; Piro et al. 2017; Drago et al. 2018; Most et al. 2018; Rezzolla, Most & Weih 2018; Ruiz, Shapiro & Tsokaros 2018). This, in turn, would constrain the EOSs of matter at extreme densities (Lattimer 2012). It is therefore important to identify signatures indicative of the formation of long-lived remnants. The presence of temporarily extended X-ray activity immediately after a merger would be one indication that a BH did not form in a dynamical time-scale after the merger (Metzger et al. 2008a; Rowlinson et al. 2013; Zhang 2013; Metzger & Piro 2014; Siegel & Ciolfi 2016a,b; Wang et al. 2016; Gao et al. 2017; Sun, Zhang & Gao 2017; Murase et al. 2018). Another would be the change in the character of the optical counterpart to the merger due to the irradiation of the ejecta by the central object (Metzger & Fernández 2014; Lippuner et al. 2017), the production of magnetized outflows (Metzger, Thompson & Quataert 2018), or the thermalization of the spin-down luminosity of the remnant (Yu, Zhang & Gao 2013; Metzger & Piro 2014; Gao et al. 2015; Kisaka, Ioka & Nakar 2016; Siegel & Ciolfi 2016a,b; Gao et al. 2017). Finally, long-lived remnants could be revealed by the appearance of bright radio flares raising on time-scales of years from the merger (Gao et al. 2013; Metzger & Bower 2014; Gompertz et al. 2015; Hotokezaka & Piran 2015; Fong et al. 2016; Horesh et al. 2016).
In this work, we employ general-relativistic merger simulations with realistic microphysics to study the formation of long-lived remnants and discuss their evolution during the subsequent viscous time-scale. We show that massive and supramassive NSs are born with angular momenta significantly exceeding the mass-shedding limit for uniformly rotating NSs and, as a consequence, they are likely to give rise to massive outflows over the viscous time-scale. These could produce luminous kilonova counterparts that would be smoking gun evidence for the formation of massive or supramassive NSs if detected by future UV/optical/infrared follow ups on gravitational wave (GW) events or SGRBs. We also constrain the spin of the remnants, and we discuss the implication of our results for the magnetar model of SGRBs and the role of thermal effects for the stability of the merger remnant. In our discussion ‘remnant’ is used to indicate all gravitationally bound matter left after the merger. Conversely, where needed, we use the expression ‘NS remnant’ to denote the centrally condensed part of the remnant.
2 MERGER REMNANTS
2.1 Simulation set-up
Our analysis is based on the results of about 35 NS merger simulations performed with the WhiskyTHC code (Radice & Rezzolla 2012; Radice, Rezzolla & Galeazzi 2014a,b). Our simulations span a range of total gravitational masses Mg = M1 + M2 between 2.4 and 3.4 M⊙, and mass ratios q = M2/M1 between 0.85 and 1.0. We adopt four tabulated nuclear EOSs broadly consistent with current laboratory and astrophysical constraints: the DD2 EOS (Hempel & Schaffner-Bielich 2010; Typel et al. 2010), the BHBΛϕ EOS (Banik, Hempel & Bandyopadhyay 2014), the LS220 EOS (Lattimer & Swesty 1991), and the SFHo EOS (Steiner, Hempel & Fischer 2013). We include an approximate treatment of neutrino cooling using the scheme discussed in Radice et al. (2016a). Results from 29 of these simulations have already been presented in Radice et al. (2018) and Zappa et al. (2018). Our data set also contains one simulation modelling the merger of a (1.35 + 1.35) M⊙ binary using the DD2 EOSs and including the effects of neutrino absorption using the M0 scheme presented in Radice et al. (2016a). Neutrino absorption does not significantly affect the outcome of the merger, but its inclusion is necessary for a quantitative prediction of the electromagnetic (EM) counterparts (Perego, Radice & Bernuzzi 2017). Neutrinos determine the properties of the ejecta, and in particular their electron fraction, especially in the polar regions (Sekiguchi et al. 2015; Foucart et al. 2016b; Radice et al. 2016a). The electron fraction, in turn, is the most important parameter determining the nucleosynthetic yields, the nuclear heating rates, the opacities of the outflows from NS mergers, and consequently their optical/infrared signatures (Lippuner & Roberts 2015). We also performed five additional simulations at 30 per cent higher resolution to check for convergence in our results.
2.2 Time-scales
We evolve each configuration for ∼3–4 orbits to merger and for at least |$20\, {\rm ms}$| after merger, or to BH formation, if this occurs earlier. We track the evolution of the total angular momentum J by integrating the flux radiated by the system in GWs following Damour et al. (2012) and Bernuzzi et al. (2012). The integrated JGW is then subtracted from the angular momentum of the binary computed by the initial data solver. We estimate the numerical uncertainty in the determination of JGW to be less than few per cent. Indeed, the discrepancy between standard and high-resolution simulations is below 3 per cent for all of the binaries we have simulated at both resolutions. As in previous studies, we find that gravitational angular momentum losses in the post-merger remnant subside within |${\sim }10{-}20\, {\rm ms}$| after merger (Bernuzzi et al. 2016; Radice, Bernuzzi & Ott 2016b; Zappa et al. 2018). At the end of our simulations the GW radiation time-scale for angular momentum loss |$\tau _{\rm GW} = J/\dot{J}_{\rm GW}$| is typically larger than 0.5 s and rapidly increasing. This is shown in Fig. 1, where we compute τGW averaged over the last millisecond of evolution. We want to stress that, because the GW emission is rapidly decaying with time, the estimate in Fig. 1 represents a conservative lower limit. The GW time-scale should be compared to the time-scale for angular momentum transport due to turbulent viscosity. The latter is expected to be τvisc ≲ 100 ms (Hotokezaka et al. 2013; Kiuchi et al. 2018). Consequently, viscosity is the dominant mechanism determining the evolution of the remnant past the point where we interrupt our simulations. We remark that the effective viscosity due to small-scale turbulence would further reduce the GW luminosity and, hence, increase the GW time-scale (Radice 2017; Shibata & Kiuchi 2017).
Gravitational wave time-scale |$\tau _{\rm GW} = J/\dot{J}_{\rm GW}$| averaged over the last millisecond of evolution for binaries producing massive or supramassive NS remnants. We find |$\tau _{\rm GW} \gtrsim 0.5\, {\rm s}$|, which is longer than the expected viscous time-scale |$\tau _{\rm visc} \lesssim 100\, {\rm ms}$| (see the main text). Note that τGW grows rapidly past the initial |$10{-}15\, {\rm ms}$| after merger, so the values reported here represent a lower limit.
2.3 Remnant angular momentum
We show a summary of the final outcome of our simulations in Fig. 2. As typically done in the literature, we indicate simulations for which BH formation occurs within 1 ms or less after merger as ‘BH’. We categorize the other binaries according to their total baryonic mass Mb: if Mb is larger than the maximum baryonic mass of rigidly rotating NSs, as predicted by the zero-temperature neutrino-less beta-equilibrated EOSs, then the merged object is considered to be an HMNS. Otherwise, we distinguish between MNS and SMNS depending on whether Mb is smaller or larger than the maximum baryonic mass for a non-rotating NS, respectively. Despite the naming convection, it is important to remark that the outcome of mergers with masses close to the demarcation line between SMNS and HMNS is likely to depend on many factors besides the maximum mass for rigidly rotating NSs. As discussed below, mass-loss, angular momentum transport, and finite-temperature effects could all either stabilize low-mass HMNSs or trigger an early collapse for high-mass SMNSs. For these binaries the distinction between SMNSs and HMNSs might not be predictive of the evolution of the remnant over time-scales t ∼ τvisc.
Merger outcome and angular momentum at the end of the simulations. The grey shaded area shows the set of all rigidly rotating equilibrium configurations. The gravitational mass on the right axis corresponds to that of an equal-mass binary having the baryonic mass indicated by the left axis. At the end of the GW radiation time-scale the merger remnant has significantly more angular momentum than the maximum allowed for rigidly rotating configurations.
We use the publicly available code RNS (Stergioulas & Friedman 1995) to construct equilibrium sequences for rigidly rotating NSs. The sequences are constructed assuming zero temperature and neutrino-less beta equilibrium. For brevity, we refer to these equilibria as being ‘cold’. The grey shaded regions in Fig. 2 show the range they span. For a fixed J lower and upper boundaries of the shaded areas are set by the mass shedding and maximum mass limit, respectively. The tip of the shaded region marks the maximum baryonic mass configuration supported by each EOS in the case of rigid rotation. Keeping in mind the caveats we have just discussed, we label binaries with Mb larger than this limit either as HMNSs or as BHs, depending on whether a BH was promptly formed in the simulations or not. Our analysis shows that MNS and SMNS are endowed with significantly more angular momentum than that corresponding to the mass-shedding limit for equilibrium configurations. This can be seen from the fact that the fast GW-drive phase of NS mergers always ends well outside on the right of the shaded areas in Fig. 2.
Our results exclude the possibility that the SMNSs formed in binary mergers could collapse due to the lack of sufficient angular momentum support, as proposed in Ma et al. (2018). These binaries would appear on the left of the grey shaded area in Fig. 2. Moreover, we can also exclude the possibility that the angular momentum of SMNS remnants could be distributed in such a way as to leave to central part of the remnant unstable to gravitational collapse. The reason is that the rotational profiles of NS merger remnants have a minimum at their centre (Shibata & Taniguchi 2006; Kastaun & Galeazzi 2015; Endrizzi et al. 2016; Kastaun et al. 2016; Ciolfi et al. 2017; Hanauske et al. 2017), so the remnant’s core is expected to spin-up during the viscous evolution (Radice 2017). Consequently, the gravitational collapse of a hypothetical low-mass binary, if it occurs, must happen dynamically during the merger and would have been seen in the simulations.
We find that massive or supramassive remnants need to shed excess angular momentum before they can settle into equilibrium configurations. The removal of angular momentum has to occur within the viscous time-scale, which is too rapid for additional GW losses to play a significant role. Consequently, angular momentum losses must be driven by viscous effects and will likely be accompanied by mass-losses. Moreover, because the mass-shedding limit moves to lower J with decreasing Mb, this process could very effectively generate large outflows.
2.4 Viscous-driven ejecta
We estimate an upper limit to the amount of material that could be unbound by viscous processes after merger using 3D data taken at the end of our simulations. We integrate the baryonic mass and the fluid angular momentum densities1 on a series of cylindrical shells. In doing so, we implicitly assume that the space–time is close to stationary and axisymmetric at the end of our simulation. We find that the angular momentum of the system estimated in this way agrees with that measured by integrating the GW flux to within 1 per cent for all models, apart from an outlier, the LS220 binary with (1.4 + 1.2) M⊙, for which the disagreement is 4 per cent. We start from the outer edge of the grid, and we progressively subtract their contribution to the total mass and angular momentum. We proceed in our integration until the resulting Mb and J enter the region spanned by rigidly rotating equilibrium configurations. This estimate is clearly an upper limit to the viscous outflow, because it assumes that the each ejected fluid element only carries away the angular momentum it had at the beginning. In reality, because of the viscous angular momentum transport, the outer edge of the disc will be endowed with some of the angular momentum initially at smaller cylindrical radii. We remark that the main underlying assumptions of our analysis are that the minimum energy state of the system is achieved when a uniformly rotating star is formed (e.g. Hartle & Sharp 1967) and that the dynamics is dominated by the action of viscosity, which drives the system towards this minimum energy state.
Our results are illustrated in Fig. 3 for the DD2 binary (1.35 + 1.35) M⊙ simulated with neutrino reabsorption, which we take as our fiducial binary. The procedure we have just discussed generates the lower edge of the blue band in Fig. 3 representing the range of possible outcomes for the viscous evolution. The starting point for the viscous evolution is the end of the GW-dominated phase of the evolution – and the end of our simulation – marked by the star symbol in Fig. 3. We find that this binary could eject up to |${\sim } 0.17$| M⊙ of material during its viscous evolution. The upper boundary of the blue band in the figure is the unlikely case in which angular momentum is removed without any outflow.
Estimated outcomes for the viscous evolution of a (1.35 + 1.35) M⊙ binary simulated with the DD2 EOS and neutrino cooling/heating. The grey shaded area shows the set of all rigidly rotating equilibrium configurations. The solid line is a conservative estimate of the mass ejection and a possible trajectory for the viscous evolution. The blue shaded region denotes the range of all possible outcomes of the viscous evolution, which we tentatively classify according to the underlying ejection mechanism. The first (disc ejecta) regime corresponds to the ejection of matter due to the nuclear recombination of the accretion disc. The second regime (remnant ejecta) is due to viscous instabilities in the merger remnant. Overall, we find that the merger remnant has enough angular momentum to unbind up to |${\sim }0.17$| M⊙ of material.
We remark that the presence of neutrino-driven winds from the disc might alter the dynamics with respect to the simple viscous spreading model we have considered for our analysis. On the one hand, extant post-merger simulations without viscosity find that the mass entrained by the neutrino-driven wind should only be of few 10−3 M⊙ (Dessart et al. 2009; Perego et al. 2014; Martin et al. 2015; Fujibayashi et al. 2017). So neutrino-driven winds should only amount to a small correction to the viscous outflow. On the other hand, neutrino heating could play an important role, together with nuclear recombination, in unbinding material that has been transported to less gravitationally bound regions by viscosity (Lippuner et al. 2017). High-resolution general-relativistic magnetohydrodynamics (GRMHD) studies of the evolution of post-merger accretion discs with neutrinos will be needed to quantify the relative importance of nuclear recombination and neutrino heating.
Our conservative estimate of the viscous ejecta for our fiducial DD2 (1.35 + 1.35) M⊙ binary decreases by |${\sim }0.01$| M⊙ when neutrino absorption is not included in the simulation (lower panel of Fig. 4). The reason is that the inclusion of neutrino absorption inflates the outer part of the accretion disc in the region |$\varpi \gtrsim 80\, {\rm km}$|. This pushes some of the material to larger radii, where it can be unbound with a smaller expenditure of angular momentum (equation 2). The inner part of the remnant is only weakly affected, so this effect is muted when computing the upper limit on the viscous ejecta.
Upper limit of the viscous ejecta (upper panel) and conservative estimate (lower panel) as a function of the total baryonic mass of the binary. The masses are normalized to the maximum mass for uniformly rotating equilibria MRNS. Supramassive and massive merger remnants are expected to eject up to |${\sim } 0.2$| M⊙ of material.
We point out that the evaporation of the disc due to its nuclear recombination is not specific to binaries forming long-lived remnants. Indeed, it is expected to occur even if the central object is a BH (Beloborodov 2008; Metzger, Piro & Quataert 2008b; Lee et al. 2009; Fernández & Metzger 2013; Metzger & Fernández 2014; Fernández et al. 2015; Siegel & Metzger 2017). However, while BHs formed in NS mergers are well below the Kerr limit (Kiuchi et al. 2009; Kastaun et al. 2013; Bernuzzi et al. 2016; Zappa et al. 2018), long-lived remnants necessarily have to dissipate a significant fraction of their angular momentum within the viscous time (Fig. 2). Consequently, the case of a long-lived remnant is qualitatively and quantitatively different and could result in more massive outflows. For this reason, we distinguish two possible components of the viscous ejecta: the ‘disc’ and the more general ‘remnant’ ejecta. The first component is due to the recombination of the discs, while the second is due to the settling of a long-lived remnant to a uniformly rotation equilibrium. We tentatively identify the disc ejecta component with our conservative estimate of the ejecta and the remnant ejecta component with everything exceeding the conservative estimate.
We repeat the analysis for 14 other binaries producing long-lived remnants. Note that we exclude from this analysis five of our binaries for which the 3D data necessary for this analysis have been lost. Our results are shown in Fig. 4. We find that the formation of massive or supramassive NSs in binary mergers could be accompanied by the ejection of up to |${\sim }0.2$| M⊙ of material within few viscous time-scales. The more conservative estimate using equation (2) yields viscous ejecta mass ∼4 times smaller. Of the five high-resolution binaries we perform to quantify the numerical uncertainty of our simulation three form a long-lived remnant: the (1.35 + 1.35) M⊙ binaries with the BHBΛϕ and DD2 EOSs, and the (1.4 + 1.2) M⊙ binary with the DD2 EOSs. The typical numerical uncertainties in the determination of the ‘disc’ and ‘remnant’ ejecta are less than 25 and 13 per cent, respectively. We conclude that |${\sim }0.05{-}0.2$| M⊙ of material should be generically ejected during the viscous phase of the evolution of long-lived NS-merger remnants.
2.5 Stability of the remnants and neutrino emission
Our simulations indicate that long-lived remnants from binary NS mergers are not only born with excess angular momentum, but also with excess gravitational mass compared to cold rigidly rotating equilibria. This is shown in Fig. 5. We find that long-lived NS merger remnants have gravitational masses |${\sim }0.08$| M⊙ larger than the corresponding equilibrium models having the same baryonic mass, but zero temperature. Given the long GW time-scale and the neutrino luminosities at the end of our simulations, we can infer that most of the excess of gravitational binding energy will be radiated in the form of neutrinos. The cooling time-scale for the massive NS remnant is of ∼2–3 s (Sekiguchi et al. 2011). These conditions are analogous to those found in newly born NSs in core-collapse supernovae (CCSNe; Burrows, Mazurek & Lattimer 1981; Burrows & Lattimer 1986; Pons et al. 1999; Fischer et al. 2010; Hüdepohl et al. 2010; Roberts et al. 2012; Roberts & Reddy 2016).
Difference between the gravitational mass of long-lived merger remnants and that corresponding to rigidly rotating equilibrium configurations having the same number of baryons. Each point represents a simulation. Note that this estimate does not account for the binding energy of the material ejected by viscous-driven wind. However, this should amount to at most a few per cent correction to the reported values. Masses on the x-axis are normalized by the maximum mass for a rigidly rotating NS predicted by the EOS MRNS. We find that long-lived merger remnants need to liberate |${\sim } 0.08$| M⊙ of gravitational binding energy before settling down.
Differently from CCSNe, however, the temperatures reached in mergers are such that the maximum mass for a stable rigidly rotating ‘hot’ NS remnant is actually smaller than that for cold equilibria, as pointed out by Kaplan et al. (2014). They found that uniformly rotating configurations with temperature profiles similar to those found in simulations can support |${\sim }0.1 $| M⊙ less baryonic mass than cold configurations. On the one hand, finite temperature and finite neutrino chemical potential only contribute a modest |${\sim }10{{\ \rm per\ cent}}$| increase of the pressure in the core of the merger remnant, at densities |${\sim }10^{15}\, {\rm g} \, {\rm cm}^{-3}$|, so finite temperature cannot stabilize the NS remnant against gravitational collapse. On the other hand, thermal support inflates the mantle of the NS remnant, i.e. the region with subnuclear densities. Because of the extended envelope, uniformly rotating sequences reach the mass-shedding limit at lower angular frequencies (Kaplan et al. 2014). This implies that a merger NS remnant that is formally supramassive according to the cold EOSs could actually be hypermassive. In other words, it is possible to form supramassive NS remnants with baryonic masses and thermodynamical profiles for which there is no rigidly rotating equilibrium. These NS remnants could either shed their excess mass or collapse to a BH within their viscous evolution.
2.6 Spin of long-lived NS remnants
Rotational periods for rigidly rotating NSs at the mass-shedding limit. This corresponds to the spin period of a long-lived merger remnant after viscosity has erased the differential rotation.
Fitting coefficients a and b (see equation 3) for the spin of long-lived remnants at the end of the viscous evolution and maximum error for |$M_\mathrm{ b} \gt 2$| M⊙ in milliseconds e.
| EOS | a | b | e |
|---|---|---|---|
| 2H | −0.27 | 1.18 | 0.05 |
| ALF2 | −0.23 | 0.85 | 0.04 |
| APR | −0.21 | 0.69 | 0.12 |
| BHBlp | −0.27 | 0.91 | 0.03 |
| DD2 | −0.20 | 0.93 | 0.04 |
| ENG | −0.20 | 0.77 | 0.04 |
| H4 | −0.35 | 0.94 | 0.02 |
| LS220 | −0.34 | 0.82 | 0.06 |
| MPA1 | −0.17 | 0.84 | 0.02 |
| MS1 | −0.21 | 1.10 | 0.02 |
| MS1b | −0.20 | 1.07 | 0.03 |
| NL3 | −0.23 | 1.11 | 0.03 |
| SFHo | −0.27 | 0.74 | 0.03 |
| SLy | −0.25 | 0.72 | 0.06 |
| TM1 | −0.31 | 1.02 | 0.03 |
| TMA | −0.35 | 0.96 | 0.02 |
| EOS | a | b | e |
|---|---|---|---|
| 2H | −0.27 | 1.18 | 0.05 |
| ALF2 | −0.23 | 0.85 | 0.04 |
| APR | −0.21 | 0.69 | 0.12 |
| BHBlp | −0.27 | 0.91 | 0.03 |
| DD2 | −0.20 | 0.93 | 0.04 |
| ENG | −0.20 | 0.77 | 0.04 |
| H4 | −0.35 | 0.94 | 0.02 |
| LS220 | −0.34 | 0.82 | 0.06 |
| MPA1 | −0.17 | 0.84 | 0.02 |
| MS1 | −0.21 | 1.10 | 0.02 |
| MS1b | −0.20 | 1.07 | 0.03 |
| NL3 | −0.23 | 1.11 | 0.03 |
| SFHo | −0.27 | 0.74 | 0.03 |
| SLy | −0.25 | 0.72 | 0.06 |
| TM1 | −0.31 | 1.02 | 0.03 |
| TMA | −0.35 | 0.96 | 0.02 |
Fitting coefficients a and b (see equation 3) for the spin of long-lived remnants at the end of the viscous evolution and maximum error for |$M_\mathrm{ b} \gt 2$| M⊙ in milliseconds e.
| EOS | a | b | e |
|---|---|---|---|
| 2H | −0.27 | 1.18 | 0.05 |
| ALF2 | −0.23 | 0.85 | 0.04 |
| APR | −0.21 | 0.69 | 0.12 |
| BHBlp | −0.27 | 0.91 | 0.03 |
| DD2 | −0.20 | 0.93 | 0.04 |
| ENG | −0.20 | 0.77 | 0.04 |
| H4 | −0.35 | 0.94 | 0.02 |
| LS220 | −0.34 | 0.82 | 0.06 |
| MPA1 | −0.17 | 0.84 | 0.02 |
| MS1 | −0.21 | 1.10 | 0.02 |
| MS1b | −0.20 | 1.07 | 0.03 |
| NL3 | −0.23 | 1.11 | 0.03 |
| SFHo | −0.27 | 0.74 | 0.03 |
| SLy | −0.25 | 0.72 | 0.06 |
| TM1 | −0.31 | 1.02 | 0.03 |
| TMA | −0.35 | 0.96 | 0.02 |
| EOS | a | b | e |
|---|---|---|---|
| 2H | −0.27 | 1.18 | 0.05 |
| ALF2 | −0.23 | 0.85 | 0.04 |
| APR | −0.21 | 0.69 | 0.12 |
| BHBlp | −0.27 | 0.91 | 0.03 |
| DD2 | −0.20 | 0.93 | 0.04 |
| ENG | −0.20 | 0.77 | 0.04 |
| H4 | −0.35 | 0.94 | 0.02 |
| LS220 | −0.34 | 0.82 | 0.06 |
| MPA1 | −0.17 | 0.84 | 0.02 |
| MS1 | −0.21 | 1.10 | 0.02 |
| MS1b | −0.20 | 1.07 | 0.03 |
| NL3 | −0.23 | 1.11 | 0.03 |
| SFHo | −0.27 | 0.74 | 0.03 |
| SLy | −0.25 | 0.72 | 0.06 |
| TM1 | −0.31 | 1.02 | 0.03 |
| TMA | −0.35 | 0.96 | 0.02 |
Our estimated spin periods are significantly smaller than those typically inferred for the progenitors of SGRBs with extended emission in the context of the magnetar model. Those are typically found to be |${\sim } 10\, {\rm ms}$| (Fan, Wu & Wei 2013; Gompertz et al. 2013). A possible way to resolve the tension with the magnetar model would be to assume that GW losses could continue past the viscously driven phase of the evolution and spin-down the remnant over a time-scale of many seconds to minutes (Fan et al. 2013; Gao et al. 2016). GW emission might be supported by secular instabilities in the remnant (Chandrasekhar 1970; Friedman & Schutz 1978; Lai & Shapiro 1995; Stergioulas 2003; Corsi & Meszaros 2009; Doneva, Kokkotas & Pnigouras 2015; Paschalidis et al. 2015; East et al. 2016b; East, Paschalidis & Pretorius 2016a; Lehner et al. 2016; Radice et al. 2016b), or by deformations due to a strong toroidal field (Fan et al. 2013).
We remark that the GW luminosity of the one-armed instability during the first |${\sim }50\, {\rm ms}$| of the post-merger evolution is |${\sim }10^{51}\, {\rm erg}\, {\rm s}^{-1}$| and does not show strong evidence for decay (Radice et al. 2016b). If the one-armed instability were to persist without damping, then it would remove all of the NS remnant rotational energy, which is |${\sim }10^{53}\, {\rm erg}$| (e.g. Margalit & Metzger 2017), over a time-scale of |${\sim }100\, {\rm s}$|. This time-scale is compatible with the spin-down time-scale inferred from the magnetar model (Fan et al. 2013). If so, the GW signal from the one-armed instability would be detectable by Adv. LIGO (Aasi et al. 2015) and Virgo (Acernese et al. 2015) up to a distance of |${\sim }100\, {\rm Mpc}$| for optimally oriented sources (Radice et al. 2016b).
Alternatively, it is possible that SGRBs with extended emission could be the result of the accretion-induced collapse of white dwarfs (Dessart et al. 2009; Abdikamalov et al. 2010; Bucciantini et al. 2012), although the host environments and the offsets from the host galaxy of SGRBs are more consistent with the expectations from NS mergers (Berger 2014; Kumar & Zhang 2014).
3 ELECTROMAGNETIC SIGNATURES
Matter ejected during merger and the subsequent viscous evolution synthesizes heavy elements through the so-called rapid neutron capture process (r-process; Lattimer & Schramm 1974; Eichler et al. 1989; Meyer 1989; Freiburghaus, Rosswog & Thielemann 1999; Korobkin et al. 2012; Wanajo et al. 2014; Just et al. 2015; Martin et al. 2015; Lippuner et al. 2017; Thielemann et al. 2017; Hotokezaka, Beniamini & Piran 2018). The resulting abundances depend sensitively on the neutron richness (i.e. on the electron fraction Ye), entropy, and expansion velocity of the material (e.g. Hoffman, Woosley & Qian 1997; Lippuner & Roberts 2015; Thielemann et al. 2017). Different ejection channels produce outflows with different properties resulting in different nucleosynthetic yields. For the conditions relevant to NS mergers, the nucleosynthesis outcome depends mainly on Ye. For Ye ≳ 0.25, the production of nuclei stops at mass numbers A ∼ 120. The production of lanthanides is possible for Ye ≲ 0.25, while even more neutron-rich material (Ye ≲ 0.15) is necessary to synthesize actinides (Lippuner & Roberts 2015).
The radioactive decay of the freshly synthesized r-process nuclei in the ejecta powers a UV/optical/infrared transient: the kilonova (sometimes also called macronova; Li & Paczynski 1998; Kulkarni 2005; Metzger et al. 2010; Roberts et al. 2011; Barnes & Kasen 2013; Kasen, Badnell & Barnes 2013; Tanaka & Hotokezaka 2013; Grossman et al. 2014; Rosswog et al. 2014; Rosswog et al. 2017). Its properties depend primarily on the rate at which radioactivity deposits heat in the material and on the time-scale over which the expanding matter becomes transparent to thermal photons. The ejecta composition is key to set the photon opacity of the ejecta, κ. In particular, the presence of lanthanide and actinides is expected to significantly increase κ, delaying the kilonova peak and shifting its spectrum to larger wavelengths (Barnes & Kasen 2013; Kasen et al. 2013; Tanaka & Hotokezaka 2013).
The detection of a transient compatible with a kilonova (AT2017gfo; Arcavi et al. 2017; Coulter et al. 2017; Drout et al. 2017; Evans et al. 2017; Kasliwal et al. 2017; Nicholl et al. 2017; Smartt et al. 2017; Soares-Santos et al. 2017; Tanvir et al. 2017) in association to GW170817 confirmed our present understanding of NS mergers and gave, for the first time, the possibility to constrain their ejecta properties and nucleosynthetic yields (Kasen et al. 2017; Rosswog et al. 2018). The analysis of the light curves and of the spectrum revealed the presence of a bright, blue, component peaking at |${\sim } 1\, {\rm d}$| after the merger, which is thought to have been powered by material moving at ∼0.3c. This was followed by a redder component peaking at ∼5 d and originating from more opaque and more slowly expanding material [Chornock et al. 2017; Cowperthwaite et al. 2017; Drout et al. 2017; Nicholl et al. 2017; Perego et al. 2017; Tanaka et al. 2017; Tanvir et al. 2017; Villar et al. 2017; see however Waxman et al. (2017), Yu & Dai (2017), and Li et al. (2018) for alternative interpretations].
We estimate the properties of the kilonova signature associated with the formation of long-lived merger remnants using the semi-analytical model we introduced in Perego et al. (2017). This includes the contribution of ejecta with different physical origin, geometry, and thermodynamical properties (details below). We calibrated the free parameters in this model using AT2017gfo in Perego et al. (2017). For the calculation of the light curves, we assume azimuthal symmetry and discretize the solid angle in 30 slices, equally spaced in cos (θ), θ being the polar angle. We place the observer at a distance of 40 Mpc and at a relative inclination of 45° with respect to the symmetry axis.
We take the DD2 (1.35 + 1.35) M⊙ binary with neutrino heating as our fiducial model. We use simulation data for the dynamical ejecta, i.e. the part of the material unbound at the time of merger, and we vary the amount of the secular ejecta to explore the range of all possible outcomes of the viscous evolution. For the former component, we consider azimuthally averaged profiles of the mass, Ye, and expansion velocity of the ejecta from the merger simulation. We assume low effective photon opacity |$\kappa _{\rm blue} = 1.0 \ {\rm cm^2\, g^{-1}}$| for the ejecta with Ye ≥ 0.25. We assume lanthanide-rich opacity |$\kappa _{\rm red} = 10 \ {\rm cm^2\, g^{-1}}$| if Ye < 0.25.
We also include an ejecta component due to the neutrino ablation of the outer layers of the accretion disc. Note that this is a distinct component of the ejecta from the viscous outflows. Following Perego et al. (2014) and Martin et al. (2015), we assume that 5 per cent of the disc is ejected in the form of a wind. The mass of the disc at the end of our simulation is 0.16 M⊙, so the wind amounts to |$8\times 10^{-3}$| M⊙ of material. Since neutrino-driven winds are only mildly neutron rich, we assume low effective photon opacity for this ejecta component (κej;wind = κblue).
As discussed in the previous section, we subdivide the viscous outflow in two parts: disc and remnant viscous ejecta. The first is assumed to be due to the nuclear recombination of the accretion disc, and would be present also for a short-lived remnant. The second is due to the viscous outflow from the massive NS and is expected only for long-lived remnants. The disc component is expected to display a broad distribution in Ye which would translate in an effective opacity intermediate between the high and low opacities of lanthanide-poor and -rich material, respectively. Ejecta with these properties is sometimes referred to as the purple component (e.g. Tanaka et al. 2017; Villar et al. 2017). For this component we take |$\kappa _{\rm ej; disc} = \kappa _{\rm purple} = 5 \ {\rm cm^2\, g^{-1}}$|, which is consistent with the AT2017gfo photometry after the first few days (Perego et al. 2017). We assume the remnant viscous ejecta to be less neutron rich than either the dynamical ejecta or the disc wind ejecta because of the neutrino irradiation from the remnant (Fujibayashi et al. 2018), and consequently we assume its opacity to be κej;NS = κblue. Our results do not qualitatively change if we assume |$\kappa _{\rm ej; NS} = 5 \ {\rm cm^2\, g^{-1}}$|, but there are quantitative differences, see Appendix A. We assume the disc viscous ejecta to have a sin 2(θ) mass distribution as in Perego et al. (2017) and the remnant viscous ejecta to be isotropic. Expansion velocities for both viscous outflows are taken to be spatially isotropic and with an rms value of 0.06c (Perego et al. 2017). As we explore the range of possibilities, we first switch on the disc viscous ejecta and increase it up to a maximum value of 0.05 M⊙, then we add the remnant viscous ejecta up to the upper limit found in the previous section |$M_{\rm ej}^{\max } = 0.17$| M⊙.
We remark that our model does not include the thermalization of the spin-down luminosity from the merger remnant, which would further enhance the kilonova signal (Yu et al. 2013; Metzger & Piro 2014; Gao et al. 2015; Kisaka et al. 2016; Siegel & Ciolfi 2016a,b; Gao et al. 2017). We will explore this possibility in future works.
In Fig. 7 we show light curves obtained from our kilonova model for three representative photometric bands, namely V, J, and K. The coloured bands correspond to increasing values of the ejecta mass in the viscous components. Light curves generated by varying the amount of the disc viscous ejecta span the light shaded bands. The light curves generated by varying the amount of the remnant viscous ejecta span the dark shaded regions. The most relevant properties of each light curve as a function of the total ejected mass are summarized in Fig. 8. There we present the peak times, magnitudes, and (temporal) widths of the kilonova signal. The latter are defined as the time interval about the peak where the light curve varies by one magnitude.
Kilonova colour light curves for our fiducial binary (see the main text). The coloured bands correspond to the possible outcomes to the viscous evolution shown in Fig. 3. The coloured solid lines correspond to the conservative estimate of the ejecta mass derived in Section 2. The black lines are the prediction for a BHNS merger also ejecting 0.18 M⊙ of material (see the main text for the details). The viscous outflows launched with the formation of long-lived NS merger remnants could power unusually bright kilonova light curves.
Kilonova peak time (upper left panel), peak magnitude (upper central panel), and peak duration (upper right panel) for our fiducial binary as a function of the mass of the viscous ejecta. Kilonova peak magnitude versus peak time (lower left panel) and peak duration versus peak time (lower right panel). We find strong correlation between these key quantities and the ejecta mass. Note the effect of the low-opacity (κej;NS = 1 cm2 g−1) remnant ejecta for Mej ≥ 0.06 M⊙. A bright, slowly evolving kilonova with a blue component at early times would be a clear evidence for the formation of a massive or supramassive NS remnant in a binary NS merger.
Increasing the amount of the viscous ejecta boosts the transient brightness in all bands. However, the V-band peak time and duration are only marginally affected by the presence of a large viscous ejecta. Conversely, a large viscous ejection produces significantly brighter peaks in the J and K infrared bands. The peaks are shifted to later times and have larger temporal widths. Notably, the increase of mass in the remnant wind produces a second peak in the J band at times longer than 1 d. This peak becomes the dominant one when the remnant viscous ejecta is turned on. The K band is the most sensitive to changes in the amount of the viscous ejecta which affect its peak brightness, time, and duration.
The merger of a NS and a BH can also result in the dynamic ejection of up to |${\sim }0.1$| M⊙ of material and in the formation of massive accretion discs (Shibata & Taniguchi 2006; Duez et al. 2008; Etienne et al. 2008; Etienne et al. 2009; Pannarale, Tonita & Rezzolla 2011; Foucart 2012; Foucart et al. 2014; Foucart et al. 2015; Kyutoku et al. 2015). Extreme mass ratio or very eccentric double NS mergers could also eject a similarly large amount of matter (East & Pretorius 2012; Rosswog, Piran & Nakar 2013; Radice et al. 2016a; Dietrich et al. 2017). We investigate whether the kilonova signal associated with the formation of a long-lived remnant in a double NS merger could be distinguished from the kilonova following a BHNS merger with a large mass ejection.
To this aim, we construct the synthetic kilonova signal for a hypothetical BHNS merger ejecting the same amount of material as our fiducial binary NS system, but with the geometry/composition expected for BHNS mergers. More in detail, we assume that 0.05 M⊙ of material is ejected by tidal torques. This material is expected to be very neutron rich and have κ = κred. We assume that the rest of the ejecta originates from the accretion torus formed from the tidal disruption of the NS. Part of the disc outflows, 0.003 M⊙, is in the form of a lanthanide-free neutrino-driven wind, for which we take κ = κblue. An additional 0.127 M⊙ is assumed to be due to the nuclear recombination of the disc. For the latter, we assume similar properties to the viscous disc ejecta from NS mergers: intermediate opacity κ = κpurple and sin 2(θ) angular distribution. The results are shown in Figs 7 (black line) and 8 (black symbols).
We find that, while the kilonova light curves from the two systems share some similarities, they also have important differences that would make them distinguishable. Kilonovae associated with the formation of long-lived remnants peak at a late time in the red bands and are significantly brighter in all bands after the peak times. Furthermore, if the viscous ejecta from the remnant is lanthanide-free, as is assumed to be the case in Figs 7 and 8, then the kilonova peak luminosities in the blue/green bands are significantly larger than those associated with BHNSs. On the other hand, if the viscous ejecta from the remnant are contaminated with lanthanides, then the peak luminosities alone are not sufficient to distinguish long-lived remnants from BHNSs. However, the luminosities after the peak time are still significantly larger in the case of long-lived remnants (Fig. A1) that a determination would still be possible for well-observed kilonovae.
4 CONCLUSIONS
We have studied the outcome of NS mergers by means of numerical relativity simulations focusing on the properties of long-lived or stable remnants. Our calculations employed four microphysical EOSs and an effective treatment of neutrino cooling. We also accounted for heating and compositional changes due to the absorption of neutrinos in one of our simulations. We have compared the properties of long-lived merger remnants to those of rigidly rotating equilibrium configurations.
We have found that the post-merger starts with a short |${\sim }10{-}20\, {\rm ms}$| phase where the evolution is mainly driven by the emission of GWs, as also reported by Bernuzzi et al. (2016), Radice et al. (2016b), and Zappa et al. (2018). Subsequently, the GW luminosity drops substantially. The characteristic time-scale associated with the removal of angular momentum due to GWs exceeds ∼0.5 s, for some binaries by orders of magnitude, and is still growing rapidly at the end of our simulations. This significantly exceeds the time-scale associated with the redistribution of angular momentum operated by the effective turbulent viscosity in the remnant, |$\tau _{\rm visc} \lesssim 0.1\, {\rm s}$| (Hotokezaka et al. 2013; Kiuchi et al. 2018), and it is also likely to exceed the neutrino-cooling time-scale, |$\tau _\nu \sim 2{-}3\, {\rm s}$| (Sekiguchi et al. 2011). Thus, the remnant evolution is mainly driven by the effects of viscosity and neutrino losses. After having reached solid body rotation and over even longer time-scales of many seconds, minutes, or hours, the remnant spins down due to residual GW losses and magnetic torques.
The evolution of the remnants over the viscous time is non-trivial. The reason is that, after the short, GW-driven, post-merger transient, the NS merger remnants are still endowed with too much angular momentum to reach an equilibrium. More precisely, we have shown that there exists no uniformly rotating equilibrium configuration to which the merger remnant can relax under the action of viscosity while conserving baryon mass and angular momentum. Instead, massive and supramassive NSs formed in mergers need to dissipate a significant fraction of their angular momentum within the viscous time-scale. Angular momentum redistribution is likely be accompanied by the emission of massive outflows since GW losses are negligible during this phase of the evolution. These viscous-driven outflows could potentially exceed those typically expected from neutrino-driven winds and from the nuclear recombination of the remnants’ accretion disc. Our results indicate that, for a typical binary, the transition to a uniformly rotating equilibrium could be accompanied by the ejection of up to |${\sim }0.2$| M⊙ of material. The mass ejection is expected to be driven by a combination of effective turbulent viscosity, nuclear recombination, and neutrino heating. However, the details of the ejection process are still not well understood, especially when long-lived remnants are formed. Long-term high-resolution neutrino-radiation GRMHD simulations will be needed to understand how post-merger discs evolve.
Massive and supramassive merger remnants have gravitational masses |${\sim }0.08$| M⊙ larger than those of equilibrium configurations having the same number of baryons. Our results suggest that most of the associated energy is liberated with the emission of neutrinos on a cooling time-scale of few seconds. At the same time Kaplan et al. (2014) showed that, for the temperatures reached in mergers, trapped neutrinos and thermal support yield only minor contributions to the pressure in the core of the remnant. However, hot rigidly rotating equilibrium sequences with increasing angular frequency reach the mass-shedding limit before cold beta-equilibrated sequences. Consequently, the maximum baryonic mass achievable for hot rigidly rotating NSs is |${\sim }0.1$| M⊙ smaller than that of cold rotating NSs. We deduce that the fate of binaries with total masses close to the threshold for the formation of HMNSs depends on a complex interplay between mass ejection and neutrino cooling whose outcome is difficult to predict. For example, remnants with masses below the maximum for cold rigidly rotating NSs could still collapse because of the gravitational mass excess with which they are formed. Conversely, massive remnants could become stable following the ejection of large amounts of material during their viscous evolution. Understanding the long-term evolution of systems with masses close to this threshold is urgent, especially in view of the current efforts to constrain the NS EOSs using the outcome of NS mergers (Margalit & Metzger 2017; Rezzolla et al. 2018; Ruiz et al. 2018). This will be the object of our future work.
Even though our models cannot yet predict the precise path undertaken by the viscous evolution of the remnant, we can nevertheless constrain the spin of the remnant once solid body rotation has been established. This is because, according to our simulations, the end result of the viscous evolution must be close to the mass-shedding limit. This corresponds to spin periods |$P_0 \lesssim 1\, {\rm ms}$|. We have shown that these can be estimated from the final baryonic mass of the remnant using a simple linear fit. The values we found are, however, much smaller than those, around |$10\, {\rm ms}$|, typically inferred from the analysis of SGRBs in the context of the magnetar model (Fan et al. 2013; Gompertz et al. 2013). This tension could be resolved under the assumption that GW losses persist even after the remnant has reached solid body rotation. The spin-down time-scale associated with this persistent emission could be |$\tau _{\rm GW} \sim 100\, {\rm s}$| (Fan et al. 2013; Gao et al. 2016). GW observations of a nearby merger event forming a long-lived remnant might detect this extended signal or severely constrain the magnetar model2 (Fan et al. 2013; Gao et al. 2016).
We have used the model of Perego et al. (2017) to produce synthetic light curves of kilonovae associated with the formation of long-lived NS merger remnants. We have found that the inclusion of viscous-driven ejecta from the merger remnant, in addition to the other outflow components, can boost the peak brightness of the emission by up to one magnitude in all bands. It also significantly broadens the width of the light curves and shifts the peak time in the near-infrared by up to several days. The resulting kilonova is peculiarly bright, blue, and slowly evolving, and would be easily distinguished from kilonovae associated with NS mergers producing BHs or BHNS mergers, despite the fact that the formers can also produce large outflows. Its detection in concomitance with an SGRB or a GW event would constitute smoking gun evidence for the formation of a long-lived remnant.
ACKNOWLEDGEMENTS
It is a pleasure to acknowledge J. Roulet for help with the RNS code, W. Del Pozzo for help with optimizing and improving the kilonova code, and A. Burrows, K. Hotokezaka, and K. Murase for discussions. DR acknowledges support from a Frank and Peggy Taplin Membership at the Institute for Advanced Study and the Max-Planck/Princeton Center (MPPC) for Plasma Physics (NSF PHY-1523261). AP acknowledges support from the INFN initiative ‘High Performance data Network’ funded by CIPE. DR and AP acknowledge support from the Institute for Nuclear Theory (17-2b programme). SB acknowledges support by the EU H2020 under ERC Starting Grant, no. BinGraSp-714626. BZ acknowledges NASA NNX15AK85G for support. Computations were performed on the supercomputers Bridges, Comet, and Stampede (NSF XSEDE allocation TG-PHY160025), on NSF/NCSA Blue Waters (NSF PRAC ACI-1440083), Marconi (PRACE proposal 2016153522), and PizDaint/CSCS (ID 667).
Footnotes
For a fluid with stress energy tensor Tμν this is defined as Tμνnμϕν, where nμ is the normal to the t = const hypersurface and ϕμ is the generator of the rotations in the orbital plane.
REFERENCES
APPENDIX: EFFECT OF REMNANT EJECTA OPACITY
Kilonova light curves (top panel) and dependence of the kilonova peak properties on the ejecta mass (lower panels) for our fiducial binary. Here, we assume the additional ejecta component from the SMNS to be contaminated with lanthanides, with an opacity of |$\kappa _{\rm ej; NS} = 5\, {\rm cm}^2 {\rm g}^{-1}$|. This figure should be contrasted with Figs 7 and 8 that are generated assuming |$\kappa _{\rm ej; NS} = 1\, {\rm cm}^2 {\rm g}^{-1}$|.
