Constraining dense matter physics using f-mode oscillations in neutron stars

In this undergraduate project, f-mode oscillations in neutron stars are used to constrain the equation of state of dense matter. For the first time, a systematic investigation of the role of nuclear saturation parameters on the mode oscillations is performed. It is found that the uncertainty in the determination of effective nucleon mass plays the most significant role in controlling the f-mode frequencies. Correlations of the frequencies with astrophysical observables relevant for asteroseismology are also investigated. Future detection of f-mode frequencies could then provide a unique way of constraining nuclear empirical parameters and therefore the behaviour of dense matter.


Introduction
With the direct detection of gravitational waves (GW) from binary merger (GW170817) of neutron stars (NSs) [1], a new window of opportunity has opened up to directly probe their interior composition. In conjunction with astrophysical observations using multi-wavelength space-based and ground-based telescopes, gravitational wave detectors now introduce the possibility of multi-messenger astronomy, from which several astrophysical NS observables can be derived [2].
Neutron stars are compact remnants of the evolution of massive stars, associated with core collapse supernovae. These objects have extreme densities up to several times normal nuclear matter density in the core, orders of magnitude beyond matter densities possible to investigate with terrestrial experiments. It is therefore inevitable to resort to theoretical models, in order to extrapolate the known physics at nuclear saturation to higher densities and isospin asymmetries. This extrapolation leads to large uncertainties, reflected by the proliferation of models. Different schemes (ab-initio, phenomenological) have been applied to describe such equations of state (EoS) [3]. The ab initio many-body models use realistic nucleon-nucleon (NN) interaction which is obtained by fitting the NN scattering data. The phenomenological models are based on nuclear density functional theories with effective nucleon-nucleon interactions, which are calibrated by reproducing the nuclear bulk properties around saturation density. Different models show different high density behaviour, resulting in different predictions for NS observables, thus making it difficult to interpret and compare directly with astrophysical data.
NSs are visible throughout the electromagnetic spectrum. A number of global properties of NSs, such as its mass and radius, can be deduced from multi-wavelength astronomical data. Accurate estimations of NS masses using Post-Newtonian effects in relativistic NS binaries indicate maximum masses to be close to 2 M sol [4,5]. Radius estimates obtained from X-ray data are less precise, but the recently launched NICER mission [6] is soon expected to measure radii with a precision of up to 5%. These observables can be obtained from the NS EoS by solving equations for hydrostatic equilibrium. Comparison with astrophysical data then allows us to put important constraints on the EoS models and consequently on the nature of dense matter.
Gravitational waves (GW) are considered one of the most promising tools for constraining dense matter physics, as they can directly probe the interior composition of neutron stars. Quasi-normal modes may be excited in oscillating NSs, producing copious amounts of GWs. These modes, such as fundamental modes (f -modes), pressure modes (p-modes), buoyancy g-modes, rotational r-modes, are classified according to the restoring forces that bring the system back to equilibrium. The most exciting fact is that the mode frequencies and damping timescales contain signatures of the interior composition of neutron stars.
GWs may be emitted from NSs, both isolated or in binary. It was shown that during a merger, the NSs in binary exert strong tidal forces on one another, and the resulting deformation depends on their compactness [1]. This can lead to important constraints on stellar radii and consequently on the dense matter EoS [7]. Further, quasi-normal modes may also be excited during the merger and post-merger phases [8]. Particularly interesting are the f -modes, as they are expected to have frequencies of ∼ 1-2 kHz and produce generous amounts of gravitational radiation through the CFS (Chandrasekhar-Friedman-Schutz) mechanism when unstable [9].
It is therefore of great interest to constrain the NS EoS using studies of f -modes. The goal of NS asteroseismology is to express the f -mode frequency or damping timescale in terms of global NS observables, independent of the underlying EoS. Detection of f -modes would then allow us to invert such relations and obtain constraints for the NS EoS [10,11,12]. However, most of such studies are confined to polytropic or parametrized EoSs due to their simplicity [10,13,14]. There are a few investigations that considered realistic equations of state, but with only a few representative parameter sets [11,12,15]. Although the dispersion of frequencies due to different EoSs is evident, it is impossible to compare between the chosen EoSs as they correspond to very different nuclear matter properties. Further, the EoSs chosen do not span the entire parameter space of empirical quantities allowed by state-of-the-art nuclear experimental data. Fits using such EoSs have been used to derive relations between mode frequencies and global variables (mean density, compactness). It is thus very difficult to extract any direct information about dense matter physics from such correlations.
In this work, which was done as a short exploratory undergraduate project, we investigate the influence of the underlying nuclear matter properties on f -mode frequencies. Within the framework of the Relativistic Mean Field (RMF) model, we systematically explore the parameter space allowed by current nuclear experimental data and perform studies of the sensitivity of the mode frequencies to the uncertainties in each of these parameters. We also investigate the role of these uncertainties on the fit relations relevant for asteroseismology of non-rotating NSs.
In Sec. 2, we describe the microscopic calculation of the EoS as well as the determination of macroscopic NS observables from such an EoS. We also elaborate on the pulsation equations to be solved to obtain f -mode frequencies in non-rotating NSs. We present the results of our investigations in Sec. 3, including the sensitivity study and asteroseismology relations. We summarize our findings in Sec. 4.

Microscopic description
As already mentioned in Sec. 1, we perform this investigation within the framework of the RMF model [16]. Such models have already been applied successfully to describe nuclear matter and nuclei [3]. In order to obtain the EoS of nuclear matter in the RMF Model, we start from the following interaction Lagrangian density: where Ψ N is the Dirac spinor for nucleons N , m is the vacuum nucleon mass. The interaction among the nucleons is mediated by the exchange of the scalar (σ), vector (ω) and the isovector (ρ) mesons. The isoscalar nucleon-nucleon couplings g σ and g ω are determined by fixing them to nuclear saturation properties. The σ meson self-interaction terms b and c ensure the correct description of nuclear matter at saturation density. The isovector and mixed ω-ρ couplings g ρ and Λ ω can be related to empirical quantities such as symmetry energy (J sym ) and its slope (L sym ) [17,18,19,20]. The quartic ω self-coupling ζ is set to zero. We do not consider terms O(3) and above in the expansion with density and asymmetry, as nuclear experimental data to constrain such parameters have large uncertainties.
The uncertainty in the nuclear empirical quantities derived from saturation data is also reflected in the uncertainty in the determination of the RMF model parameters. In order to test our results, we first consider using two commonly used parametrizations [21]: GM1 and GM3, for which the EoSs are well known. Once our numerical scheme is verified, we want to explore the entire RMF parameter space defined by present uncertainties of theoretical and experimental nuclear saturation data. The range of values of empirical RMF parameters explored in this work is summarized in Table 1. For each individual parameter "variation" within the ranges shown in the table, the others are kept at the "fixed" values. Model Given the Lagrangian density Eq. (1), one can solve the equations of motion of the constituent particles as well as those of the mesons. In the mean-field approach, the meson fields are replaced by their mean-field expectation values. One then calculate the EoS (pressure -energy density relationship) using this RMF model. The energy density is given by [20] The pressure P can be derived from the energy density using the Gibbs-Duhem relation [22] where the nucleon chemical potentials are given by

Macroscopic description
Given an EoS, one can obtain the macroscopic NS structure by solving the Tolman Oppenheimer Volkov (TOV) equations of hydrostatic equilibrium [22] dm(r) dr = 4πε(r)r 2 , dp(r) Integrating the TOV equations from the centre of the star to the surface, one can obtain global NS observables, such as mass (M ), radius (R) and compactness (C = M/R). The tidal deformability (Λ) can be obtained by solving a set of differential equations coupled with the TOV equations [23]. These can then be compared to the state-of-the-art limits derived from astrophysical data, in order to impose constraints on the dense matter EoS.  Fig. 1 In Fig. 1, we display all the EoS curves (pressure vs energy density in M eV /f m 3 ) considered in this work, corresponding to the parameter sets in Table 1. The corresponding mass-radius relations are given in Fig. 2. It is evident from Fig. 2 that for all the EoSs considered in this work, the maximum masses lie above 2M sol (horizontal dotted line in the figure).

Solving the mode pulsation equations
The aim of this work is to investigate the influence of the uncertainty in the empirical nuclear parameters on NS f -modes. In general, one must solve coupled fluid and space-time perturbation equations to obtain the mode frequencies [24]. However, the situation simplifies if we consider weak gravitational fields and neglect the metric perturbations. This approach, known as the Cowling approximation, has been widely applied for studying Newtonian as well as relativistic NSs. Although ideally one must employ fully linearized equations in general relativity, it has been shown that the difference on applying the Cowling approximation is less than 20% for f -modes [25]. For this study, we only consider non-rotating NSs (See [15] for investigations of f -modes in rotating NSs).
Within the Cowling approximation (neglecting metric perturbations), the spherically symmetric background spacetime characterized by a line element is given by: ds 2 = −e 2Φ(r) dt 2 + e 2Λ(r) dr 2 + r 2 dθ 2 + r 2 sin(θ)dφ 2 The fluid Lagrangian displacement vector is defined by where W (t, r) and V (t, r) are functions of r and t. The fluid perturbations are decomposed into spherical harmonics Y lm (θ, φ). Using these variables, the perturbation equations for the fluid oscillations can be obtained from the conservation of the energy-momentum tensor δ(∇ ν T µν ) = 0. For a harmonic time dependence, the perturbation functions can be written as W (t, r) = W (r)e iωt and V (t, r) = V (r)e iωt . On simplification, the pulsation equations required to be solved in order to obtain these frequencies are given by [24]: The functions V (r) and W (r), along with the frequency ω, characterize the perturbation vector. We solve coupled Eq. 8 on a fixed background metric from the origin (r = 0), where the solutions behave approximately like [25,26] where C is an arbitrary constant. The other boundary condition that needs to be fulfilled is that the perturbation to the pressure must vanish at the star's surface (r = R). Such condition reads The real part of the complex eigenfrequencies ω give the frequency of the perturbation functions.

Testing the numerical scheme
In order to test the numerical scheme, we first reproduce well known results for f -mode frequencies for the GM1 parameter set (see e.g., [25,26]) and also for the GM3 set. In Figure 3 we display the fundamental f -mode frequencies for these two reference parameter sets. As expected, the f -modes have frequencies within 1-3 kHz compatible with previous results reported in the literature.

Calculation of l = 2 f -modes
In Sec. 2, we discussed the uncertainties associated with the nuclear saturation parameters, which in turn result in the uncertainty in the EoS. Having tested the numerical scheme of the f -mode frequencies in Sec. 3.1, we now extend it to investigate the entire parameter space. We vary each of the nuclear saturation parameters individually within their known uncertainties listed in Table 1 and study the sensitivity of the f -mode frequencies to each of the variations. It is obvious from the panels in Fig. 4 that the influence of varying the isoscalar parameters, energy per particle at saturation E sat and compressibility K sat , is negligible. Similarly, the variation of the isovector parameters, symmetry energy J sym and its slope L sym , do not vary the f -mode frequencies significantly, as can be seen from the same figure. In Fig. 5, we find a small non-zero effect of the variation in the saturation density ρ sat . However, the most important influence on the frequencies   Fig. 6, we can clearly differentiate between the different EoSs with varying effective mass in the f -modes as a function of the stellar mass. The frequencies as a function of NS masses are seen to vary between 2 to 2.6 kHz, and the variation is monotonic with increasing m * /m. This could have interesting consequences for extracting dense matter physics from the detection of f -mode frequencies for known stellar masses.

f -modes and tidal deformability
Among the various NS astrophysical observables that can help to constrain the nuclear EoS, one of the most promising quantities that has recently emerged is the tidal deformability. With the discovery of gravitational waves from mergers of NSs, it was seen that during the inspiral phase the NSs exert strong gravitational forces on each other, and the deformation produced depends on their EoS [1]. The relation between the dimensionless tidal deformability and the NS compactness C = M/R is given by Λ = 2k 2 3C 5 , where k 2 is the second tidal Love number. In order to estimate the influence of the m * /m, the dominant empirical parameter that affects the f -mode frequencies, we first demonstrate its effect on the EoS and tidal deformability. In Fig. 7, we plot the dimensionless tidal deformability Λ as a function of NS mass (in M sol ). We see that the variation of effective nucleon masses m * /m influences Λ as a function of the NS mass, and the range of values is consistent with recent observations [27,28]. Now, we investigate the influence of variation in m * /m on the f -mode frequencies as a function of Λ in Fig. 8. We find that the variation in effective nucleon masses causes a change in f -mode frequencies in the range (2-2.6 kHz), for a corresponding change in dimensionless tidal deformability in the range ∼ 1-10. From the curves, one can obtain the frequencies corresponding to the lower limit of tidal deformability (∼ 160, constraints from terrestrial nuclear experiments), and the upper limit of tidal deformability of NSs extracted from the GW170817 event by LIGO and VIRGO Collaborations (∼ 580) which can provide an interesting constraint for the effective nucleon mass and hence the nuclear EoS [13]. We already mentioned in the previous section that the tidal deformability is related to the stellar compactness. In order to use GW observations to estimate the NS mass and radius and to differentiate between different families of EoS, empirical relations between the frequency of f -modes and the compactness of the star may be useful [10,11,29]. We study in Fig. 9 the effect of variation of m * /m on the f -mode frequencies as a function of C = M/R. Alternatively, one may also obtain the compactness from observations of the gravitational redshift Z from spectral lines, as they are related by Z = (1 − 2C) −1/2 − 1. The corresponding variation of f -mode frequencies as a function of Z are also shown in Fig. 10.

Asteroseismology relations for f -modes
In the previous Sec. 3.2.2, we studied the correlations between the f -mode frequency and stellar compactness M/R. Now we investigate some universal relations relating the scaled f -mode frequency versus compactness of NSs [11]. Correlations of the f -mode frequency and damping time (or appropriately scaled by stellar mass and radius) with compactness and tidal deformability have been recently investigated [13,14]. In Fig. 11 we plot the frequency (scaled by the NS mass in km) and Fig. 12 the frequency ω (scaled by the NS radius in km) as a function of compactness. We find that ωM as a function of C retains a universal behaviour while ωR shows slight deviation from universal behaviour with variation in effective nucleon mass.
If we want to use GW data from NSs to estimate global properties like radius or mass, we would require empirical relations as much independent of the underlying EoS as possible. It was already pointed out in previous studies [12,15,29] that f -mode frequencies show a linear behaviour as a function of square-root of the mean density (M/R 3 ) 1/2 . However such studies considered polytropic EoSs [10] or realistic EoSs with only a few representative parameter sets [11,15]. Many of the chosen EoSs are not consistent with recent astrophysical data, such as the maximum mass constraint of 2M sol . Linear fits derived from such studies would then introduce an error in obtaining relations useful for asteroseismology.
In Sec. 3.2, we already determined the key nuclear parameters that control the f -mode amplitudes. In Figs. 13 and 14, we now plot the l = 2 f -mode frequencies as a function of (M/R 3 ) 1/2 for varying   [15]. By inverting such relations, the mass and the radius of each NS can then be calculated upto the accuracy of the average density and the stellar compactness. Recent calculations [30,31] have investigated the instability window relevant to f -modes and concluded that higher order modes (e.g. l = 3, 4 modes) could be more dominant than the quadrupole l = 2 mode. In analogy with l = 2 case, we perform a study of the variation of the f -mode frequencies as a function of (M/R 3 ) 1/2 for the variation in ρ sat for higher orders in l = 3, 4. These are given in Figs. 15 and 17. The variation of f -mode frequencies with effective mass for higher orders are similarly plotted in Figs. 16 and 18 respectively. We also obtained the linear fits to the curves for the l = 3 case: (i) For variation in saturation density y = 0.7628 + 2.9039x

Higher order f -modes
(ii) For variation in effective mass: where y is the l = 3 f -mode frequency (in kHz) as a function of x = (M /R 3 ) 1/2 defined previously.

Discussions
In this work, we studied the role of the underlying dense matter physics on f -modes. Within the framework of the RMF model, we performed a systematic investigation of the influence of the uncertainties in empirical nuclear EoS parameters, consistent with recent nuclear experimental data, on the f -mode frequencies. The equations of state obtained from these parameter sets were all compatible with the 2M sol maximum mass limit. We found the nucleon in-medium mass to be the most dominant parameter, while the saturation density has a small non-zero effect. Variation in effective mass between 0.55-0.75 resulted in a corresponding variation of f -mode frequencies between 2-2.6 kHz as a function of NS masses.
We further investigated the effect of uncertainty in effective nucleon mass on the relation between f -mode frequencies and tidal deformability, compactness and gravitational redshift. These quantities can be derived and will be further constrained with the rapid improvement in multi-messenger astrophysical observations. We investigated whether the uncertainty in effective mass affects the universality in asteroseismology relations. We also obtained linear fits to correlations between f -mode frequencies and square root of the average density, for quadrupole (l = 2) and higher order (l = 3, 4) modes.
In order to extract information about f -modes in NS merger remnants, one must include rotation effects. We have not considered rotating NSs in this investigation, as it is beyond the scope of this undergraduate project. One must also consider the influence of exotic constituents of matter (hyperons, kaons, deconfined quark matter) on f -modes. Several works have attempted such an investigation, but again for a few selected EoSs [25,26]. We are also currently performing such an investigation and the results will be reported in a forthcoming publication.
With the recent detection of gravitational waves from binary compact objects, particularly mergers of binary neutron stars, prospects of constraining nuclear physics using gravitational waves have become attractive. It has been speculated that f -modes are among the most significant sources of GWs due to the Chandrasekhar-Friedman-Schutz (CFS) mechanism, for both isolated NSs or in the post-merger scenario. Recent studies suggest that GWs produced by unstable l = m = 2 and the l = m = 4 f -modes could be detectable by the future Einstein Telescope for sources in the Virgo cluster or for l = m = 3 modes even by the LIGO/VIRGO [31]. Then the possibility of distinguishing between NS EoSs using information about f -modes and NS global observables could prove be very interesting.