Gravitationalwave signal from binary neutron stars: a systematic analysis of the spectral properties
Abstract
A number of works have shown that important information on the equation of state of matter at nuclear density can be extracted from the gravitational waves emitted by merging neutronstar binaries. We present a comprehensive analysis of the gravitationalwave signal emitted during the inspiral, merger and postmerger of 56 neutronstar binaries. This sample of binaries, arguably the largest studied to date with realistic equations of state, spans across six different nuclearphysics equations of state and ten masses, allowing us to sharpen a number of results recently obtained on the spectral properties of the gravitationalwave signal. Overall we find that: (i) for binaries with masses differing no more than , the frequency at gravitationalwave amplitude’s maximum is related quasiuniversally with the tidal deformability of the two stars; (ii) the spectral properties vary during the postmerger phase, with a transient phase lasting a few millisecond after the merger and followed by a quasistationary phase; (iii) when distinguishing the spectral peaks between these two phases, a number of ambiguities in the identification of the peaks disappear, leaving a simple and robust picture; (iv) using properly identified frequencies, quasiuniversal relations are found between the spectral features and the properties of the neutron stars; (v) for the most salient peaks analytic fitting functions can be obtained in terms of the stellar tidal deformability or compactness. Altogether, these results support the idea that the equation of state of nuclear matter can be constrained tightly when a signal in gravitational waves from binary neutron stars is detected.
pacs:
04.25.D, 04.25.dk, 04.30.Db, 26.60.KpI Introduction
The recent measurement from the advanced interferometric LIGO detectors (Harry et al., 2010) of the first direct gravitationalwave (GW) signal from what has been interpreted as the inspiral, merger and ringdown of a binary system of black holes Abbott et al. (2016) marks, in many respects, the beginning of GW astronomy. As additional advanced detectors such as Virgo (Accadia et al., 2011), and KAGRA (Aso et al., 2013) are going to become operational in the next few years, we are likely to soon witness also signals from the inspiral and postmerger of neutronstar binaries or neutronstar–blackhole binaries. Such systems are not only excellent sources of GWs, but also the most attractive scenario to explain the phenomenology associated with short gammaray bursts (SGRBs). Starting from the first suggestions that merging neutron stars could be behind these phenomena Narayan et al. (1992); Eichler et al. (1989) and supported by the circumstantial evidence coming from numerous astronomical observations (see Berger (2014) for a recent review), numerical simulations have sharpened the contours of this scenario Shibata and Uryū (2000); Baiotti et al. (2008); Anderson et al. (2008a); Liu et al. (2008); Bernuzzi et al. (2012a). What we now know rather reliably is that the merger of a binary neutronstar (BNS) system inevitably leads to the formation of a massive metastable object, which can either collapse promptly to a black hole or survive up to thousands of seconds Ravi and Lasky (2014), emitting gravitational and electromagnetic radiation Zhang and Mészáros (2001); Metzger et al. (2011); Rezzolla and Kumar (2015); Ciolfi and Siegel (2015). Furthermore, if the neutron stars have large magnetic fields and extended magnetospheres, the inspiral can be accompanied by a precursor electromagnetic signal Palenzuela et al. (2013), while the merger can lead to instabilities Siegel et al. (2013); Kiuchi et al. (2015) and to the formation of magnetically confined jet structures once a torus is formed around the black hole Rezzolla et al. (2011); Paschalidis et al. (2015). Hence, the prospects of a multimessenger gravitational and electromagnetic signal are particularly good in the case of merger of binary neutron stars.
In addition to a strong electromagnetic signal, the merger of a neutronstar binary also promises a GW signal that will contain important signatures of the equation of state (EOS) of matter at nuclear densities. These signatures are contained both in the inspiral and in the postmerger signals. The former is reasonably well understood analytically Flanagan and Hinderer (2008); Baiotti et al. (2010); Bernuzzi et al. (2012b, 2015a); Hinderer et al. (2016) and can be tracked accurately with advanced highorder numerical codes Radice et al. (2014a, b, 2015) and over many orbits now Hotokezaka et al. (2016). More importantly, the instantaneous GW frequency at amplitude maximum has been shown to correlate closely with the tidal deformability of the two stars Read et al. (2013); Bernuzzi et al. (2014a); Takami et al. (2015). The latter part of the signal has already been studied in the past Shibata (2005); Oechslin and Janka (2007); Baiotti et al. (2008), but it has become the focus of attention particularly over the last few years Bauswein and Janka (2012); Stergioulas et al. (2011); Bauswein et al. (2012); Hotokezaka et al. (2013); Bauswein et al. (2014); Takami et al. (2014, 2015); Bernuzzi et al. (2015b); Palenzuela et al. (2015); Bauswein and Stergioulas (2015); Dietrich et al. (2015); Foucart et al. (2016); Lehner et al. (2016). This large bulk of work has reached some generally agreedupon conclusions, but also has raised points that are a matter of debate that we hope to clarify here.
Let us therefore start by summarising the aspects of the postmerger GW signal that seem to be robust and confirmed by several groups employing a variety of numerical methods and mathematical approximations. Given a realistic BNS system, namely, with a total mass between and , and a mass difference between the two components that is or less, then the spectrum of the postmerger GW signal will present at least three strong peaks Bauswein and Janka (2012); Stergioulas et al. (2011); Bauswein et al. (2012); Hotokezaka et al. (2013); Bauswein et al. (2014); Takami et al. (2014, 2015). These peaks were dubbed and in Refs. Takami et al. (2014, 2015) and were found to satisfy the following approximate relation: . A simple mechanical toy model was also presented in Takami et al. (2015), that provided an intuitive explanation on the origin of these peaks and why they should be almost equally spaced. In addition to these three peaks, another peak can be identified in the power spectral density (PSD) of the GW signal, although not always. This is given by the coupling between the fundamental mode (which yields the peak) and a quasiradial (fundamental) axisymmetric mode, i.e., with ; this mode was dubbed in Ref. Stergioulas et al. (2011). Finally, Ref. Bauswein and Stergioulas (2015) introduced the concept of the peak frequency and associated it to a “rotating pattern of a deformation of spiral shape”. It is difficult to measure this motion in a numericalrelativity calculation without the possible contamination of spatialgauge effects; however, as we will show in the following, coincides with in the large majority of the cases so that, in the end, the strongest and most robust features of the postmerger signal are confined to four frequencies: , and , respectively.
In addition to selecting the most salient spectral features of the postmerger GW signal, the analyses of Refs. Bauswein and Janka (2012); Stergioulas et al. (2011); Bauswein et al. (2012); Hotokezaka et al. (2013); Bauswein et al. (2014); Takami et al. (2014, 2015); Bernuzzi et al. (2015b) have also tried to associate the values of these frequencies to the stellar properties of the two stars before the merger and hence to their EOS. This was initiated by Refs. Bauswein and Janka (2012); Bauswein et al. (2012), who showed that correlated with the radius of the maximummass nonrotating configuration; the correlation found was rather tight, but restricted to binaries having all the same total mass of Bauswein and Janka (2012); later on, the correlation was explored also when considering a larger sample of masses Bauswein et al. (2012), and it was noted that the correlation was not as tight as previously expected. Indeed, Refs. Takami et al. (2014, 2015) discovered that the correlations of the peaks with the radius of the maximummass nonrotating configuration depend on the total mass of the system and hence are not “universal” in the sense of being only weakly dependent on the EOS (see also Refs. Bauswein et al. (2012); Hotokezaka et al. (2013)). At the same time, Refs. Takami et al. (2014, 2015) showed that the peak is correlated with the total compactness of the stars in a quasiuniversal manner, and highlighted a number of other correlations (24 different ones were presented in Fig. 15 of Takami et al. (2015)); some of these correlations had been already presented in the literature, e.g., in Refs. Read et al. (2013); Bernuzzi et al. (2014a), while most of them were presented there for the first time. More specifically, correlations were found between and frequencies and the physical quantities of the binary system, e.g., the stellar compactness, the average density, or the dimensionless tidal deformability (some of these correlations will be further discussed below).
Some of the results reviewed above are robust, while other are less so. In particular, although the interpretation of the largest peak in the spectrum is normally attributed to the mode of the HMNS, the interpretation of the lowfrequency peak (or peaks) is still subject to debate. More specifically, it is unclear whether such a peak (or peaks) correlates in a “universal” manner with the stellar properties, as shown in Takami et al. (2015) for the peak, or not, as shown in Bauswein and Stergioulas (2015) for the peak.
The purpose of this paper is to try and clarify this debate and, more specifically, to show that for the cases considered here the frequency either coincides with the frequency or falls in a part of the PSD of the signal where no significant power can be found, thus explaining why no universal behaviour was found in Ref. Bauswein and Stergioulas (2015). We reach this conclusion by extending the sample of binaries considered in Refs. Takami et al. (2014, 2015) to include the fully generalrelativistic simulations of the very lowmass binaries with total mass and that were suggested by Bauswein and Stergioulas (2015) to be missing in our sample. In addition, we also consider the largest masses that can be supported for a timescale after the merger sufficiently long to yield an accurate spectrum; depending on the EOS, these masses can be as large as . The large majority of binaries have equal masses, but we consider also four different instances of unequalmass binaries with mass difference of about .
This complete sample of binaries, which counts a total of 56 binaries and doubles the sample presented in Takami et al. (2015), is arguably the largest studied to date with nuclearphysics EOSs and in full general relativity (a sample of comparable size but not in full general relativity was already presented in Ref. Bauswein et al. (2012)). After a systematic analysis of the complete sample, it was possible to sharpen a number of arguments on the spectral properties of the GW signal and that can be summarised as follows:

the GW frequency at amplitude maximum is related to quasiuniversally with the tidal deformability of the two stars; this correlation is strong for equalmass binaries and holds as long as the binaries have masses that do not differ of more than .

the postmerger signal is characterised by a transient phase lasting a few millisecond after the merger, which is then followed by quasistationary phase.

the spectral properties of the GW signal vary during the postmerger phase with a marked difference between the transient and the quasistationary phase.

spectrograms are particularly useful when selecting spectral features in the transient phase because peaks that appear in the short transient, may be subdominant when analysed in terms of the full PSD.

when distinguishing the spectral peaks between the transient and quasistationary phases, a number of ambiguities in the identification of the peaks disappear, leaving a rather simple and robust picture;

the strongest and most robust features of the postmerger signal are confined to four frequencies: , and , where and is the result of a mode coupling.

a number of “universal” relations can be found between the main spectral features and the physical properties of the neutron stars.

for all of the correlations found, simple analytic expressions can be given either in terms of the dimensionless tidal deformability or of the stellar compactness.
When considered as a whole, these results support the idea that the equation of state of nuclear matter can be tightly constrained when a strong postmerger signal in GWs is measured.
The paper is organised as follows. Section II provides a very brief summary of the mathematical and numerical methods used to obtain our results in full general relativity, while Sect. III is dedicated to the illustration of our results. In particular, in Sects. III.1, III.2, and III.3 we concentrate on the waveform properties coming from transient signals, quasistationary signals and from the analysis of the full PSDs, respectively. The analysis of the correlations of the spectral signatures with the stellar properties is instead presented in Sect. III.4, while Sect. IV contains our conclusions and future prospects. Three different appendices offer details on the full set of binaries considered (Appendix A), on the analysis of some modes (Appendix B) or on a twodimensional fit employed in our analysis (Appendix C).
Ii Mathematical and numerical Setup
The mathematical and numerical setup used for the simulations reported here is the same discussed in Takami et al. (2014, 2015) and presented in greater detail in other papers Baiotti et al. (2008, 2009, 2010). For completeness we review here only the basic aspects, referring the interested reader to the papers above for additional information. All of our simulations have been performed in full general relativity using a fourthorder finitedifferencing code McLachlan Brown et al. (2009); Löffler et al. (2012), which solves a conformal traceless formulation of the Einstein equations Nakamura et al. (1987); Shibata and Nakamura (1995); Baumgarte and Shapiro (1999), with a “” slicing condition and a “Gammadriver” shift condition (Alcubierre et al., 2003; Pollney et al., 2007). At the same time, the generalrelativistic hydrodynamics equations are solved using the finitevolume code Whisky (Baiotti et al., 2005), which has been extensively tested in simulations involving the inspiral and merger of BNSs Baiotti et al. (2008, 2009); Rezzolla et al. (2010); Baiotti et al. (2010).
The hydrodynamics equations are solved employing the HartenLaxvan LeerEinfeldt (HLLE) Harten et al. (1983) approximate Riemann solver Harten et al. (1983), which is less accurate but more robust, in conjunction with a Piecewise Parabolic Method (PPM) for the reconstruction of the evolved variables Colella and Woodward (1984). For the time integration of the coupled set of the hydrodynamic and Einstein equations we have used the Method of Lines (MOL) in conjunction with an explicit fourthorder RungeKutta method Rezzolla and Zanotti (2013). In all our simulations we prescribe a CourantFriedrichsLewy (CFL) factor of to compute the size of the timestep.
ii.0.1 Grid structure and extents
We employ an adaptivemesh refinement (AMR) approach that follows closely the one adopted in Baiotti et al. (2008); Kastaun et al. (2013) and where the grid hierarchy is handled by the Carpet meshrefinement driver Schnetter et al. (2004). It implements vertexcentered mesh refinement, also known as the boxinbox method, and allows for regridding during the calculation as well as multiple grid centres. The timestep on each grid is set by the Courant condition and by the spatial grid resolution for that level. Boundary data for finer grids are calculated with spatial prolongation operators employing fifthorder polynomials and with prolongation in time employing secondorder polynomials.
During the inspiral, a grid with the finest refinement and fully covering each star is centred at the position of the maximum restmass density. The grid hierarchy is composed of six refinement levels and a refinement factor for successive levels. The grid resolution varies from (i.e., ) for the finest level, to (i.e., ) for the coarsest level, whose outer boundary is at (i.e., ). To reduce computational costs in the case of equalmass binaries, which represent the large majority of our sample, the grid structure is then replicated employing a symmetric, i.e., a symmetry of degrees around the axis^{1}^{1}1In the case of equalmass binaries we have carried out comparative simulations of the same initial data with and without the symmetry being imposed, finding no appreciable differences in the position (in frequency) of the main peaks of the PSDs. At the same time, the amplitude of such peaks can vary (up to ) and the PSD without symmetry naturally shows more power at higher frequencies as these are not suppressed by the symmetry.. Independently of the mass ratio, the whole grid is set up to be symmetric with respect to the plane both for equal and unequalmass binaries, with a reflection symmetry across the plane, again to reduce computational costs.
The number of grid points across the linear dimension of a star is of the order of , and this is roughly doubled when the merger has taken place and a HMNS has been formed. The boundary conditions are chosen to be “radiative” for the metric in order to prevent GWs (or other numerical perturbations) from scattering back into the grid, and “static” for the hydrodynamical variables.
Finally, we recall that in Ref. Takami et al. (2015) we have carried out a resolution study to assess the influence of the resolution on the spectral properties. By comparing the results with higher and lower resolutions, we have found that the use of a “medium” resolution of is sufficient to provide numerically robust measurements of the peaks, that is, peaks which differ of a few percent only from those obtained with higher resolutions (see Section IV A of Takami et al. (2015)).
ii.0.2 Equations of state
As in our previous work Takami et al. (2014, 2015), we model the stars with five “cold” (i.e., at zero temperature) nuclearphysics EOSs: i.e., APR4 Akmal et al. (1998), ALF2 Alford et al. (2005), SLy Douchin and Haensel (2001), H4 Glendenning and Moszkowski (1991) and GNH3 Glendenning (1985). All of these EOSs satisfy the current observational constraint on the observed maximum mass in neutron stars, i.e., obtained for the pulsar PSR J0348+0432 Antoniadis et al. (2013). In addition, to validate the results also across “hot” EOSs, we consider two binaries described by the LattimerSwesty EOS Lattimer and Swesty (1991), with nuclear compressibility parameter (LS220); these binaries were first studied in Radice et al. (2016), where additional information on their dynamics can be found.
The nuclearphysics EOSs are normally provided in tabular form, but it is more convenient numerically to express them in terms of a number of piecewise polytropes Read et al. (2009); Rezzolla and Zanotti (2013). Four different “pieces” are normally sufficient to reproduce to good precision most of the EOSs, with three of the pieces describing the highdensity core and one the crustal region; we refer to Table I of Ref. Takami et al. (2015) for a list of the properties of the various piecewise polytropes used here.
The cold nuclearphysics EOSs also need to be supplemented by a “hot” contribution that accounts for the considerable increase in the internal energy at the merger. This is normally done through a socalled “hybrid EOS” Rezzolla and Zanotti (2013), in which an idealfluid component that accounts for the shock heating is added to the cold part Janka et al. (1993). In practice the total pressure and specific internal energy are expressed as
(1)  
(2) 
where are given by cold nuclearphysics EOSs (expressed as piecewise polytropes), while the “thermal” part is given by
(3)  
(4) 
where is obtained through the solution of the hydrodynamics equations and is arbitrary, but constrained mathematically to be . After some experimentation carried out in Takami et al. (2015), we have concluded that values do not introduce a significant variance in the spectral properties of the GW signal and we have therefore chosen for consistency with the single polytrope (see discussion in Sect. IV B. in Takami et al. (2015)).
ii.0.3 Initial data
The initial data in our simulations represents quasiequilibrium irrotational BNSs and is computed with the multidomain spectralmethod code LORENE (Gourgoulhon et al., 2001) under the assumption of a conformally flat spacetime metric. All binaries have an initial coordinate separation between the stellar centres of , which yields at least four orbits (or more) before the merger.
The choice of the masses for the binaries is constrained by two considerations. The first one is that, given the substantial computational costs of these simulations we need to consider masses that are realistic, rather than masses that give, for instance, the largest mass difference; in practice, this implies that our masses are around , which indeed we take as our fiducial mass. The second consideration is that in order to model the postmerger reliably we need an HMNS that survives for a sufficiently large amount of time before collapsing to a black hole (e.g., at least ); stated differently, binaries with these EOSs and masses larger than , are not optimal for the postmerger analysis although we use them to have the largest possible sample.
As a result of these considerations, for each cold EOS, we have considered ten equalmass binaries with average (gravitational) mass at infinite separation in the range for the APR4, ALF2, GNH3, H4, and SLy EOSs^{2}^{2}2Note that because the highestmass binaries collapse promptly to a black hole, the only spectral information we can use in these cases is the one relative to the inspiral, i.e., .. Furthermore, as a complement to our set of equalmass binaries, we have also considered four unequalmass binaries with and mass ratio and and mass ratio ; these unequalmass binaries have been chosen for the GNH3 and SLy EOSs as examples of stiff and soft EOSs, respectively. The sample is completed by two equalmass binaries with masses and described by the the hot LS220 EOS. Detailed information on all the models and their properties is collected in Table 1 of Appendix A.
ii.0.4 Gravitationalwave signal
The GW signal is extracted at different surfaces of constant coordinate radius using the NewmanPenrose formalism, so that the GW polarization amplitudes and are related to Weyl curvature scalar by (see Sect. IV of Ref. Baiotti et al. (2008) for details)
(5) 
where the overdot indicates a time derivative and we have introduced the (multipolar) expansion of in spinweighted spherical harmonics Goldberg et al. (1967) of spinweight . In practice, all of our analysis is limited to the dominant mode, i.e., the mode
(6) 
where are the spinweighted spherical harmonics. Following previous work Read et al. (2013); Takami et al. (2015), we align the waveforms at the “time of merger”, which we set to be and define to be correspondent to the time when the GW amplitude
(7) 
reaches its first maximum. As a result, for most binaries we consider GW signals in the time interval . After defining the instantaneous frequency of the GW as , where is the phase of the complex gravitational waveform Read et al. (2013). The time of the merger is also used to define the “frequency at amplitude maximum” (or peak frequency in Ref. Read et al. (2013)) as .
Another quantity used extensively in our analysis is the PSD of the effective amplitude and defined as
(8) 
with
(9) 
and where the and indices refer to the two polarization modes. Using this PSD, we can compute the signaltonoise ratio (SNR) as
(10) 
with being the noise PSD of a given GW detector [e.g., Advanced LIGO url , or the Einstein Telescope (ET) Punturo et al. (2010a, b)].
Iii Results
iii.1 Waveform properties: transient signals
Figure 1 provides a summarising view of some of the waveforms (i.e., of for sources at a distance of ) computed in this paper and that are combined with those of Takami et al. (2015) to offer a more comprehensive impression of the GW signal across different masses and EOSs. The figure is composed of 35 panels referring to the 35 equalmass binaries with nuclearphysics EOSs that we have simulated and that have a postmerger signal of at least ; binaries with shorter postmerger (e.g., ALF2q10M1400) and unequalmass (e.g., GNH3q09M1300) are not reported in the figure but their properties are listed in Table 1. Different rows refer to models with the same mass, while different columns select the five cold EOSs considered and colourcoded for convenience. It is then rather easy to see how small differences across the various EOSs during the inspiral become marked differences after the merger. In particular, it is straightforward to observe how the GW signal increases considerably in frequency after the merger and how lowmass binaries with stiff EOSs (e.g., topleft panel for the GNH3 EOS) show a qualitatively different behaviour from highmass binaries with soft EOSs (e.g., bottomright panel for the APR4 EOS). Also quite apparent is that, independently of the mass considered, the postmerger amplitude depends sensitively on the stiffness of the EOS, with stiff EOSs (e.g., GNH3) yielding systematically larger amplitudes than soft EOSs (e.g., APR4).
What is less evident from Fig. 1 are the features of the transient GW signals emitted a few millisecond after the merger. To this scope, we present in Fig. 2 two representative examples that concentrate on a time window around the merger, i.e., one millisecond in the inspiral (lightgray shaded area) and four milliseconds after the merger. Both panels represent what could be the most realistic reference value for the mass, i.e., , with the left panel referring to a representative soft EOS (i.e., GNH3), while the right panel showing the same but for a representative stiff EOS (i.e., APR4). The top part of each panel reports the gravitational strain for a source at , while the bottom part the corresponding spectrogram, i.e., the evolution of the PSD, where timeseries segments with length and transformed with a Blackman window are overlapped by . Also marked with horizontal line of different type and colour are the various frequencies that have been so far discussed in the literature when describing the spectral properties of the GW signal.
Although we have already mentioned such spectral properties in the Introduction, but we also briefly summarise them below:

the frequencies were first introduced in Ref. Read et al. (2013) and mark the instantaneous GW frequency at the merger, i.e., at GW amplitude maximum. These frequencies were discussed in Ref. Read et al. (2013), where they were first shown to correlate with the tidal deformability of the two stars; similar findings were later reported in Refs. Bernuzzi et al. (2014a); Takami et al. (2015). The values reported here are measured from the data (see also Fig. 7 below).

the frequencies were introduced in Refs. Takami et al. (2014, 2015) and represent the three main peaks of the PSDs measured in those references^{3}^{3}3Other authors, e.g., Bauswein and Janka (2012); Bauswein and Stergioulas (2015), refer to the frequency as to , but we find this convention confusing as there are several “peaks” in the PSD.. The frequencies were found to roughly follow the relation , and a simple mechanical toy model was proposed in Takami et al. (2015) to explain simply this relation.
Figure 3: The same as in Fig. 2, but shown on a much longer timescale, i.e., after the merger. Note that all of the peaks present in the short transient stage (i.e., ) essentially disappear in the quasistationary evolution. The only exception is the peak, which slightly evolves from the frequency. In the case of a stiff EOS (e.g., for the H4 EOS, but not shown here) a trace of the mode is still present, although at very low amplitudes. 
the frequencies correspond to the fundamental mode of the HMNS and hence are equal to twice the rotation frequency of the bar deformation of the HMNS. Their values change slightly in time (by ), and we indicate with the values in the transient phase to distinguish them from the values attained in the subsequent quasistationary evolution of the GW signal.

the values reported here for , and are measured from the data. In particular, the frequencies are measured from the spectrograms, while the frequencies from the full PSDs; see Sect. III.3). For the frequencies, instead, a first guess is predicted from the analytic expression Eq. (25) of Takami et al. (2015) [i.e., Eq. (18) here]; we then use these guesses and the spectrograms to refine the measured values of the frequencies, which we report in Figs. 2–6 and in Table 2^{4}^{4}4In pratice, we take the analytic guess from Eq. (25) of Takami et al. (2015) to draw a horizontal line in the spectrogram and then correct its vertical position of a few percent till it matches the largest value of the PSD for the longest amount of time.. Such values of the frequencies are also used to obtain an improved estimate of the fitting coefficients (19). Finally, the third frequency is also predicted as .

the frequencies were first introduced in Ref. Stergioulas et al. (2011) and refer to a coupling between the fundamental mode and a quasiradial axisymmetric mode, i.e., with , and that we indicate as . The frequency of the latter mode can be measured, for instance, from the oscillations in the lapse function or the restmass density at the center of the HMNS Bauswein and Stergioulas (2015); Bauswein et al. (2016), so that are effectively measured quantities^{5}^{5}5Note that although the lapse function and the restmass density at the center of the HMNS are both a gauge quantities, they provide a robust and accurate representation of the eigenfrequnecies of oscillating compact objects, as shown in several studies Font et al. (2000, 2002); Baiotti et al. (2005); Radice and Rezzolla (2011)..

the frequencies were first introduced in Ref. Bauswein and Stergioulas (2015) and refer to the contribution to the GW signal coming from a “rotating pattern of a deformation of spiral shape”. Because shapes in gaugedependent quantities such as the restmass density are essentially impossible to measure in numericalrelativity calculations, we cannot measure these frequencies in our calculation. Hence, the values reported are those predicted from the analytic prescription given in Ref. Bauswein and Stergioulas (2015) [cf., Eq. (2) of Bauswein and Stergioulas (2015)]. It was also claimed that the peak can be roughly reproduced in a toy model, but no details were given in Ref. Bauswein and Stergioulas (2015).
The spectrograms in Fig. 2 contain a wealth of information about the transient postmerger phase. To clarify a series of imprecise and sometimes confusing statements recently appeared in the literature, we collect below the main information on the spectral properties of the transient. The points summarised below do not refer only to the models in Figs. 2 and 3, but to all of the models simulated [cf., Fig. 4].

the frequencies (blue dashed lines) are shortlived and evolve into the frequencies as the GW signal reaches its quasistationary phase (i.e., the one after ). This change is of the order of , so that .

the frequencies (light blue solid lines and green dashed lines) are instead shortlived and their amplitude becomes vanishingly small after the transient. This was remarked in Ref. Takami et al. (2015), where a simple toy model was developed to explain these frequencies as a result of the modulation of the (rotating) oscillation of the two stellar cores [cf., Appendix A of Takami et al. (2015)]. Note that although are only predicted analytically, they do coincide with the maximum values of the spectrogram. This provides an important confirmation on the correctness of the interpretation in Ref. Takami et al. (2015).

the frequencies (cyan solid lines) do not have significant power during the transient phase and it is only later that they may produce a contribution and only for rather stiff EOSs (cf., Fig. 4, which we discuss below). Note also that the frequencies are measured from and , respectively, leaving no room for interpretation.

the frequency (red dashed lines) is essentially the same as the frequency for the stiff EOS GNH3, while it is significantly different for the soft EOS APR4. The fact that in some cases but not in all cases, holds true also when considering other EOSs and was found also in Ref. Bauswein and Stergioulas (2015). We will touch on this point further below (see dashed blue and red lines in Figs. 5 and 6 and relative discussion).

because the frequency is predicted from an analytic expression [cf., Eq. (2) of Bauswein and Stergioulas (2015)], its coincidence with the frequency for stiff EOSs suggests that the two frequencies are just the same for stiff EOSs. On the other hand, for soft EOSs the frequencies do select a genuinely different mode, which are however shortlived; as we will show in Fig. 5, the total power stored in these modes is always very small in our data (see also Bauswein and Stergioulas (2015)).
To recap, the analysis of the spectrograms in the transient phase reveals that three modes are clearly visible: . The last two disappear later on, while the first one survives as and with changes of a few percent. The frequencies essentially coincide with the frequencies for stiff EOSs, while marking a different mode for soft EOSs; these latter frequencies are short lived and provide a minimal contribution to the total PSD.
iii.2 Waveform properties: quasistationary signals
We next discuss the spectral properties of the signal when considered over a much longer timescale, which we take to be at least after the merger if the HMNS does not collapse before. This is shown in Fig. 3, again for two representative EOSs and for the fiducial mass of . As anticipated in the previous Section, the only frequency surviving on these timescales is the peak (blue solid lines), which evolves slightly from the peak as the HMNS attains a quasistationary equilibrium. Other peaks, such those associated to essentially vanish after the transient, while the peak retain only small powers.
Figure 4 is rather “dense”, but provides a comprehensive summary in terms of spectrograms of the results discussed in the last two sections. In particular, the left panel reports the spectrograms for the 35 binaries presented in Fig. 1, but concentrating on the transient phase, i.e., for , while the right panel shows the spectrograms for the complete GW signal. In essence, the two panels show that:

for all the EOSs considered here, three frequencies appear in the transient phase: .

the frequencies essentially coincide with the frequencies for stiff EOSs (i.e., GNH3, H4, and ALF2), but differ for soft EOSs (i.e., SLy and APR4).

for soft EOSs, the frequencies are systematically at larger frequencies than the frequencies (see, e.g., the model with for the APR4 EOS), but yield a very contribution to the overall PSD.

after the transient, only the frequencies survive as an adjustment of the frequencies produced during the transient.

the frequency can be easily measured from the oscillations of the lapse function but the associated power in the spectrograms is always extremely small and appreciable only after the transient and for a limited period of time (see, e.g., model with for stiff EOSs GNH3, H4 and ALF2).
iii.3 Waveform properties: analysis of the full PSDs
The discussion has so far been focused on the analysis of the spectral properties of the GW signal as deduced when looking at the spectrograms. We next investigate how the spectral properties appear when analysing the full PSDs of the GW signal.
We start this discussion by considering two representative examples in Fig. 5, which reports the PSDs of a series of lowmass (i.e., ) and mediummass binaries (i.e., ); solid lines refer to the postmerger signal only, while the dotted lines report also the power during the short inspiral. The frequencies relative to the peaks , , , and are shown with vertical dashed lines of different colours. Additionally, the frequency of Ref. Bauswein and Stergioulas (2015), which was reported for and [cf., Eq. (2) of Bauswein and Stergioulas (2015)], are also shown as a reference. In essence, the PSDs in Fig. 5 reveal that:

for all the EOSs considered, the peak corresponding to the frequency is rather easy to recognise and is reasonably well reproduced by an analytic expression that we will discuss in Sect. III.4.

the frequencies as measured from the spectrograms do not correspond to any visible peak in the total PSDs; this is to be expected given that these frequencies are only short lived and their contribution to the total PSD is much smaller than that of the frequencies.

smaller but still clearly visible are the contributions of the and frequencies (the latter are not reported in Fig. 5 for clarity). This behaviour too is not surprising and is due to the short duration of these modes. Note that the frequencies in Fig. 5, which we recall are predicted analytically, also mark the presence of a local maximum in the PSD.

for stiff EOSs, e.g., GNH4 and H4, the peaks corresponding to the and frequencies are very similar, but they become distinct for soft EOSs, e.g., ALF2, SLy and APR4.

when not being comparable to , the frequencies do not seem to mark any local maximum in the PSDs, see, for example the BNS with and the SLy EOS, or the BNS with and the EOSs SLy and APR4.

the behaviour of the frequencies is far less clear. In those cases where it is comparable with the frequencies (e.g., for the SLy EOS), these frequencies can be associated to the same power excess attributed to the frequencies. In other cases, however, they are either associated to peaks with very limited power^{6}^{6}6This is the case, for instance, for the binaries with with EOS GNH3, H4 and ALF2. or are associated to peaks^{7}^{7}7This is the case, for instance, for the binaries with with EOS GNH3, H4 and ALF2, or for the binaries with with EOS APR4.. This is not surprising since the peaks result from a mode coupling and are therefore expected to be less energetic.
The properties of the PSDs listed above and illustrated in Fig. 5 are not limited to the cases of the binary masses reported in that figure. This conclusion can be reached after inspecting Fig. 6, which is the same as Fig. 5, but reports also all the other masses considered. For compactness we do not report here the PSDs relative to the unequalmass binaries considered, and whose PSDs show a very similar behaviour to the ones discussed so far.
iii.4 Waveform properties: correlations with the stellar properties
In what follows we make use of the results discussed so far to correlate the spectral properties of the GW signal with the properties of the progenitor stellar models. As discussed extensively in Refs. Bauswein and Janka (2012); Bauswein et al. (2012); Read et al. (2013); Bernuzzi et al. (2014a); Takami et al. (2014, 2015); Bernuzzi et al. (2015b), some of these correlations appear to be “universal”, i.e., only slightly dependent on the EOS, and can therefore be used to constrain the physical properties of the progenitor stars and hence the EOS.
iii.4.1 Inspiral and merger
We start by considering the inspiral part of the signal. Figure 7 reports the massweighted frequencies at amplitude maximum as a function of the tidal polarizability parameter for a generic unequalmass binary. We recall that the latter is defined as (see, e.g., Bernuzzi et al. (2014a))
(11) 
where and refer to the primary and secondary stars in the binary
(12) 
are the dimensionless tidal Love numbers, and are the compactnesses. In the case of equalmass binaries, , and expression (11) reduces to
(13) 
where the quantity
(14) 
is another commonly employed way of expressing the tidal Love number for equalmass binaries Read et al. (2013), while is its dimensionless counterpart and was employed in Takami et al. (2015).
The black solid line in the left panel of Fig. 7 shows the fit as given by Eq. (24) of Takami et al. (2015), where five realistic EOSs (i.e., APR4, ALF2, SLy, H4, GNH3) and an idealfluid EOS with were employed across a range of five different masses for each EOS. Such a fit was expressed as [cf., Eq. (24) of Takami et al. (2015)]
(15) 
where^{8}^{8}8Note that in Takami et al. (2015) the fit was actually done using the quantity , so that the coefficient reported there is different by a factor .
(16) 
Making use of the larger sample of binaries and excluding from the fit the idealfluid EOS, we can further refine the fit and obtain new and slightly modified coefficients
(17) 
The new fit is indicated with a black dashed in the left panel of Fig. 7, which also reports in an inset the same data but when represented via the dimensionless tidal deformability to highlight the essentially linear dependence in terms of this variable.
In the right panel of Fig. 7, we instead show the same correlation as in the left panel (and the corresponding fitting expressions) but when using the data taken from Refs. Read et al. (2013); Bernuzzi et al. (2014a). Note that in these cases, the frequencies reported are measured at the peak of amplitude of one of the polarization modes of the strain, i.e., at the maximum amplitude of rather than of ; we believe this slight differences is responsible for the larger variance in the correlation (see discussion in Sect. V C. in Takami et al. (2015)).
Overall, the data reported in Fig. 7 confirms what was first pointed out in Ref. Read et al. (2013), namely, that a rather tight “universal” correlation exists between the frequency at peak amplitude and the tidal deformability; for equalmass binaries, the largest difference between the values measured for and those predicted by the fit are , but the average deviation is much smaller and only.
The correlation becomes weaker when considering unequalmass binaries and this is clearly shown by the data marked with empty squares at and at . Both points refer respectively to the APR4 and GNH3 binaries with the smallest mass ratio of and this can be interpreted as the “breaking” of the universality for small (and possibly unrealistic) mass ratios. Given that the dynamics and GW signal in these cases is rather different from the corresponding equalmass binaries (the merger necessarily happens at lower frequencies as the tidal interaction is amplified and the lowermass star disrupted), this is perfectly plausible; however additional simulations will be needed to confirm this conjecture.
iii.4.2 Postmerger
Next, we consider the correlations in the postmerger part of the signal and concentrate initially on the lowfrequency peaks. The left panel of Fig. 8 reports the values of such frequencies as a function of the stellar compactness. The solid black line represents the “universal” relation first reported in Ref. Takami et al. (2015) and expressed as a cubic polynomial [cf., Eq. (25) of Ref. Takami et al. (2015)]
(18) 
where the fitting coefficients given by Takami et al. (2015) have been further refined after using the larger sample of data considered here and are given by
(19) 
Also indicated with a gray shaded band is the overall uncertainty in the fit related either to the truncation error or to the determination of the frequencies Takami et al. (2015).
Note that the data relative to the new binaries is in full agreement with our previous results in Takami et al. (2015) and that even the very lowmass binaries, i.e., with , support the universal relation. Hence, we believe that the mismatch found in Bauswein and Stergioulas (2015) is due to the incorrect association of the frequency with the frequencies.
Also reported in the figure are the putative frequencies as taken from Refs. Dietrich et al. (2015); Foucart et al. (2016), respectively for binaries with and and , evolved with EOSs MS1, H4, ALF2, and SLy Dietrich et al. (2015), as well as for binaries with and evolved with EOSs LS220, DD2, and SFHo Foucart et al. (2016). In particular, in the case of the SFHo EOS, we believe that the value reported in Table II of Foucart et al. (2016) (i.e., ) actually refers to the frequency and that the correct value for the frequency is the lowfrequency one which is clearly marked in the left panel of Fig. 6 in Foucart et al. (2016) (i.e., ). We have already remarked that in the case of soft EOSs and lowmass binaries it is very easy to confuse the frequency with the one; we believe this is one of those instances. Indeed, when this correction is made, all of the data of Refs. Dietrich et al. (2015) and Foucart et al. (2016) appears in good agreement with the expected “universal” behaviour and its variance.
The frequencies also show a very good tight correlation with the the tidal deformability. This was already reported in Ref. Takami et al. (2015) (cf., Fig. 15 of Takami et al. (2015)) and we further remark on this point in the right panel of Fig. 8, where we employ the same data of the left panel and fit it with a thirdorder polynomial expansion of the dimensionless tidal deformability
(20) 
with the coefficients being given by
(21) 
When considering only equalmass binaries, the largest deviation from the fit (20) is , while the average is only .
We next consider the correlations of the and frequencies with the stellar properties. We recall that the and frequencies correspond to the same fundamental mode of oscillation and the different denomination refers to whether the frequencies are measured during the transient phase () or in the quasiperiodic one (). These frequencies are reported as a function of the dimensionless tidal deformability in the left panel Fig. 9, where we report the data taken from all of our simulations, including those that produce a black hole before 25 ms; a similar figure was presented in Ref. Takami et al. (2015) (cf., Fig. 15 of Takami et al. (2015)), but here we also provide a fitting function in terms of the dimensionless tidal deformability, both for the and the frequencies (see inset in the left panel of Fig. 9) as
(22)  
(23) 
Relative to the equalmass binaries, the largest (average) deviations from the fit are for the frequencies and for the frequencies.
In addition, we present in the right panel of Fig. 9 the correlation between the frequencies and the tidal deformability. We recall that the frequencies are simply the difference between the frequencies, which are easy to measure from the PSDs and the frequencies that can be estimated from the analysis of the data of the simulations (e.g., via the oscillation frequencies of the central restmass density or of the lapse function). In this sense, the frequencies are straightforward to compute; yet, as discussed above, these do not always correspond to a clearly visible peak in the total PSDs (cf., Fig. 6). Using all of this data we obtain a linear fit in terms of shown as a black solid line
(24) 
A similar correlation is exhibited by the compactness (see Bauswein and Stergioulas (2015) and also the inset in the right panel of Fig. 9) and this is not particularly surprising given the fundamental nature of the frequencies. Because the correlation is weaker, the maximum and average deviations from the fit are of and , respectively.
We conclude this Section by reporting the fitting expressions for the frequencies as given by Eq. (2) in Bauswein and Stergioulas (2015) for equalmass binaries with ,
(25) 
No fitting expression was provided in Bauswein and Stergioulas (2015) for the binaries with , but it was not difficult to reconstruct the behaviour of through a quadratic twodimensional fit in terms of the compactness and average gravitational mass of the binary, i.e.,
(26) 
and where
(27) 
Expression (III.4.2) is the one used when indicating the expected frequencies in all the figures of this paper; a comparison of the twodimensional fit (III.4.2) with the various fits suggested by Ref. Bauswein and Stergioulas (2015) will be presented in Appendix C.
Iv Conclusions
The merger of a neutronstar binary is expected to be accompanied not only by a highly energetic electromagnetic signal, but also by a strong signal in GWs that will contain important signatures of the EOS of matter at nuclear densities. These signatures are contained both in the inspiral and in the postmerger signals and, following previous work of a number of authors Baiotti et al. (2010); Bernuzzi et al. (2012b); Radice et al. (2014a, b, 2015); Hotokezaka et al. (2016); Read et al. (2013); Bauswein and Janka (2012); Stergioulas et al. (2011); Bauswein et al. (2012); Hotokezaka et al. (2013); Bauswein et al. (2014); Takami et al. (2014, 2015); Bernuzzi et al. (2015b); Palenzuela et al. (2015); Bauswein and Stergioulas (2015); Dietrich et al. (2015); Foucart et al. (2016) we have here analysed the GW signal from merging neutron star binaries.
A particularly important goal of this work has been the analysis of the postmerger signal, which presents a number of quite strong spectral features that, just like spectral lines from the atmospheres of normal stars, can be used to extract information on the physical properties of the progenitor neutron stars. More specifically, we have tried to clarify whether or not some of these spectral features can be correlated with the properties of the stars before the merger in a way that is quasiuniversal, that is, only weakly dependent on the EOS.
To this scope we have simulated a large number of neutronstar binaries in full general relativity and relativistic hydrodynamics. The binaries, which are mostly with equal masses, have total masses as low as and as large as , spanning six different nuclearphysics EOSs. Because of the close similarities in the EOSs, the new set of binaries has been analysed together with the one already presented in Refs. Takami et al. (2014, 2015), thus providing a total sample of 56 binaries, arguably the largest studied to date in general relativity and with nuclearphysics EOSs.
A systematic analysis of the complete sample has allowed us to obtain a rather robust picture of the spectral properties of the GW signal and hopefully clarify a number of aspects that have been debated recently in the literature. In essence, our most important findings can be summarised as follows:

the instantaneous GW frequency when the amplitude reaches its first maximum is related quasiuniversally with the tidal deformability of the two stars.

this correlation is observed for binaries with masses that do not differ of more than .

the spectral properties vary during the postmerger phase with a marked difference between a transient phase lasting a few millisecond after the merger and a following quasistationary phase.

the most robust features of the postmerger signal pertain four frequencies: , and , where and is the result of a mode coupling.

when distinguishing the spectral peaks between these two phases, a number of ambiguities in the identification of the peaks disappear, leaving a rather simple and robust picture.

“universal” relations are found between the spectral features and the physical properties of the neutron stars.

for all of the correlations between the spectral features and the stellar properties, simple analytic expressions can be found either in terms of the dimensionless tidal deformability or of the stellar compactness.
When considered as a whole and in the light of recent direct detection of GWs Abbott et al. (2016), these results open the exciting and realistic prospects of constraining the EOS of nuclear matter via GW observations of merging BNSs. At the same time, the robustness of these results needs to be validated when different physical conditions are assumed for the merging neutron stars. These include accounting a nonzero amount of initial stellar spin Kastaun et al. (2013); Tsatsin and Marronetti (2013); Bernuzzi et al. (2014b); Tsokaros et al. (2015); Kastaun and Galeazzi (2015), evaluating the modifications introduced by ideal and resistiveMHD effects Giacomazzo et al. (2011); Dionysopoulou et al. (2015), and assessing the impact that the shear instability Ou and Tohline (2006); Corvino et al. (2010); Anderson et al. (2008b); East et al. (2016) may have on the postmerger spectrum when the instability is able to fully develop and persist on secular timescales. Furthermore, while some examples of unequalmass binaries or of genuinely threedimensional EOSs have been considered here, a more systematic analysis needs to be performed; this will be part of our future work.
Acknowledgements.
We thank L. Baiotti, L. Bovard, M. Hanauske for useful comments and discussions and D. Radice for help with the LS220 binaries. Partial support comes from JSPS KAKENHI Grant Number 15H06813, from “NewCompStar”, COST Action MP1304, from the LOEWEProgram in HIC for FAIR, and the European Union’s Horizon 2020 Research and Innovation Programme under grant agreement No. 671698 (call FETHPC12014, project ExaHyPE). The simulations were performed on SuperMUC at LRZMunich and on LOEWE at CSCFrankfurt.Appendix A Summary of stellar properties
We report in Table 1 a summary of the main physical properties of the binaries simulated in this work; some of the data has already been provided in Ref. Takami et al. (2015), but we report it also here for completeness. The various columns denote the gravitational mass ratio at infinite separation, the average gravitational mass at infinite separation, the average radius at infinite separation, the ADM mass of the system at initial separation, the baryon mass , the compactness , the orbital frequency at the initial separation, the total angular momentum at the initial separation, the dimensionless moment of inertia at infinite separation (which is tightly correlated with the compactness Breu and Rezzolla (2016)), the dimensionless tidal Love number at infinite separation, the dimensionless tidal deformability , and the contact frequency . All quantities with a bar are defined as averages, i.e., .
model  EOS  

GNH3q10M1200  GNH3  
GNH3q10M1225  GNH3  
GNH3q10M1250  GNH3  
GNH3q10M1275  GNH3  
GNH3q10M1300  GNH3  
GNH3q10M1325  GNH3  
GNH3q10M1350  GNH3  
GNH3q10M1375  GNH3  
GNH3q10M1400  GNH3  
GNH3q10M1500  GNH3  
GNH3q09M1300  GNH3  
GNH3q08M1275  GNH3  
H4q10M1200  H4  
H4q10M1225  H4  
H4q10M1250  H4  
H4q10M1275  H4  
H4q10M1300  H4  
H4q10M1325  H4  
H4q10M1350  H4  
H4q10M1375  H4  
H4q10M1400  H4  
H4q10M1500  H4  
ALF2q10M1200  ALF2  
ALF2q10M1225  ALF2  
ALF2q10M1250  ALF2  
ALF2q10M1275  ALF2  
ALF2q10M1300  ALF2  
ALF2q10M1325  ALF2  
ALF2q10M1350  ALF2  
ALF2q10M1375  ALF2  
ALF2q10M1400  ALF2  
ALF2q10M1500  ALF2  
SLyq10M1200  SLy  
SLyq10M1225  SLy  
SLyq10M1250  SLy  
SLyq10M1275  SLy  
SLyq10M1300  SLy  
SLyq10M1325  SLy  
SLyq10M1350  SLy  
SLyq10M1375  SLy  
SLyq10M1400  SLy  
SLyq10M1500  SLy  
SLyq09M1300  SLy  
SLyq08M1275  SLy  
APR4q10M1200  APR4  
APR4q10M1225  APR4  
APR4q10M1250  APR4  
APR4q10M1275  APR4  
APR4q10M1300  APR4  
APR4q10M1325  APR4  
APR4q10M1350  APR4  
APR4q10M1375  APR4  
APR4q10M1400  APR4  
APR4q10M1500 