Combined Constraints on the Equation of State of Dense Neutron-Rich Matter from Terrestrial Experiments and Observations of Neutron Stars

# Combined Constraints on the Equation of State of Dense Neutron-Rich Matter from Terrestrial Nuclear Experiments and Observations of Neutron Stars

Nai-Bo Zhang11affiliation: Department of Physics and Astronomy, Texas AM University-Commerce, Commerce, TX 75429, USA 22affiliation: Shandong Provincial Key Laboratory of Optical Astronomy and Solar-Terrestrial Environment, Institute of Space Sciences, Shandong University, Weihai, 264209, China , Bao-An Li11affiliation: Department of Physics and Astronomy, Texas AM University-Commerce, Commerce, TX 75429, USA **affiliation: Corresponding author: Bao-An.Li@Tamuc.edu , and Jun Xu33affiliation: Shanghai Institute of Applied Physics, Chinese Academy of Sciences, Shanghai 201800, China
###### Abstract

Within the parameter space of equation of state (EOS) of dense neutron-rich matter limited by existing constraints mainly from terrestrial nuclear experiments, we investigate how the neutron star maximum mass M, radius km and tidal deformability of canonical neutron stars all together constrain the EOS of dense neutron-rich nucleonic matter. While the 3-D parameter space of (curvature of nuclear symmetry energy), and (skewness of the symmetry energy and EOS of symmetric nuclear matter, respectively) are narrowed down significantly by the observational constraints, more data are needed to pin down the individual values of , and with quantified uncertainties. The largely controls the maximum mass of neutron stars. While the EOS with is sufficiently stiff to support neutron stars as massive as 2.37 M, to support the hyperthetical ones as massive as 2.74 M (composite mass of GW170817) requires to be larger than its currently known maximum value of about 400 MeV and beyond the causality limit. The upper limit on the tidal deformability of from the recent observation of GW170817 is found to provide upper limits on some EOS parameters consistent with but far less restrictive than the existing constraints of other observables studied.

Dense matter, equation of state, stars: neutron

## 1 Introduction

To summarize what terrestrial experiments and nuclear theories have taught us about the EOS of neutron-rich matter near the saturation density of symmetric nuclear matter (SNM), we first recall that the specific energy in asymmetric nucleonic matter (ANM) can be approximated parabolically in isospin asymmetry as

 E(ρ,δ)≈E0(ρ)+Esym(ρ)δ2 (1)

where is the specific energy in SNM and is the symmetry energy. Around , the and can be Taylor expanded up to as

 E0(ρ) ≈ E0(ρ0)+K02(ρ−ρ03ρ0)2+J06(ρ−ρ03ρ0)3, (2) Esym(ρ) ≈ Esym(ρ0)+L(ρ−ρ03ρ0)+Ksym2(ρ−ρ03ρ0)2+Jsym6(ρ−ρ03ρ0)3 (3)

in terms of several EOS characteristic parameters: the incompressibility and skewness of SNM as well as the slope , curvature and skewness of the symmetry energy. While the above Taylor expansions are not used in our study for neutron stars, they provide asymptotic boundary conditions at for the parameterizations we shall use to describe the EOS of neutron star matter.

Generally speaking, the EOS of neutron-rich nucleonic matter remains very uncertain mostly because of the poorly known nuclear symmetry energy especially at supra-saturation densities (Li et al., 2014). Nevertheless, thanks to the hard work of many people in both astrophysics and nuclear physics over many years, much progress has been made in constraining the EOS of neutron-rich nucleonic matter. In particular, various analyses of terrestrial nuclear experiments and astrophysical observations have constrained the , and in reasonably small ranges around MeV (Shlomo et al., 2006; Piekarewicz, 2010), MeV and MeV. Especially worth noting, the quoted most probable values of and are based on surveys of 53 analyses of different kinds of terrestrial and astrophysical data available up to Oct. 2016 (Li & Han, 2013; Oertel et al., 2017). Moreover, extensive analyses of heavy-ion reactions at intermediate and relativistic energies, especially various forms of nucleon collective flow and kaon production, have provided a reasonably tight constraining band for the EOS of SNM up to about , see, e.g., ref. Danielewicz et al. (2002). However, the coefficients characterizing the high-density behavior of neutron-rich matter, such as the , and are only loosely known to be in the range of MeV, MeV, and MeV mostly based on analyses of terrestrial nuclear experiments and energy density functionals (Tews et al., 2017; Zhang et al., 2017), respectively.

While continuous efforts have been made to constrain the EOS of dense neutron-rich matter using heavy-ion reactions which may involve rare isotopes with large neutron/proton ratios, presently no concensus has been reached yet from analyzing limited data available (Li, 2017). On the other hand, properties of neutron stars and events involving them, such as the mass, radius, moment of inertia, quadrupole deformation, pulsing frequency (Hessels et al., 2006), cooling curve (Yakovlev et al., 2001; Page et al., 2006), frequencies and damping times of various oscillating modes, spin parameter of pulsars, as well as the strain amplitude and phase evolution of gravitational waves from inspiraling neutron star binaries all depend significantly on the EOS of neutron-rich matter, see, e.g., refs. (Lattimer & Prakash, 2016; Oertel et al., 2017; Özel & Freire, 2016; Watts et al., 2016; Grigorian et al., 2016; Newton et al., 2014) for recent reviews. While the observational data of neutron star properties are relatively limited so far, they also provide stringent constraints on the EOS and guide theories of dense neutron-rich matter. In particular, the observed masses around 2M for the two most massive pulsars J1614-2230 (Demorest et al., 2010) and J0348+0432 (Antoniadis et al., 2013) restrict mostly the stiffness of dense SNM EOS, while the radii of neutron stars are known to be more sensitive to the symmetry energy around (Lattimer & Prakash, 2000, 2001). While much progress has been made in measuring the radii of neutron stars, because of the great difficulties involved especially in determining accurately the distance and modeling reliably the spectrum absorptions with different atmosphere models, the reported radii normally suffer from relatively large uncertainties. In fact, some radii extracted from different analyses and observations are still controversial, see, e.g., ref. (Miller & Lamb, 2016) for a recent review. Since we are not in a position to make any judgement on the reliability of any astrophysical observations, in our analysis here we use as an example the radius of M canonical neutron stars () in the range of km inferred from analyzing quiescent low-mass X-ray binaries in ref. (Lattimer & Steiner, 2014). Moreover, we also use the upper limit of the dimensionless tidal deformation from the recent observation of GW170817 by the LIGO+Virgo Collaborations (Abbott et al., 2017; LIGO & Virgo, 2017).

The first discovery of a neutron star merger using multiple messengers further signifies and sets an excellent example of combining and cross-checking multiple probes of dense neutron-rich matter on earth and in heaven. Given the aforementioned constraints on the EOS parameters of neutron-rich matter mostly based on terrestrial nuclear laboratory experiments, here we study how the astrophysical observations of M, km and all together constrain the high-density EOS in a way consistent naturally with the existing constraints from terrestrial nuclear experiments. For this purpose, fixing the , and at their most probably values mentioned earlier, we explore the intersections of constant surfaces with , km, km, and , respectively, in the 3-dimensional (3-D) parameter space of , and . The 3-D parameter space allowed by all three observational constraints are identified. Moreover, in constructing the EOS of neutron star matter, the crust-core transition density and pressure have to be calculated consistently. While effects of the magnitude and slope of symmetry energy at on the crust-core transition have been extensively studied in the literature, effects of the curvature are less known. We therefore will also explore contours of the and in the versus plane (2-D). The significant role of the is clearly revealed. Furthermore, within the currently known uncertainty ranges of and , by setting in Eqs. (2) and (3) as often done in the literature, we also explore how/what the same three astrophysical constraints may teach us about the EOS of dense neutron-rich matter in the parameter plane. In both the 3-D and 2-D model frameworks, the upper limit of the tidal deformability from the recent observation of GW170817 is found to provide upper limits on some EOS parameters consistent with but less restrictive than the existing constraints on them. Overall, while combing exiting constraints from both terrestrial nuclear experiments and astrophysical observations allows us to limit significantly the EOS parameter space of high-density neutron-rich matter, data of more independent observables are needed to pin down the individual values of , and with quantified uncertainties.

This paper is organized as follows. The construction of the EOS of neutron star matter is presented in Section 2. The Section 3 is devoted to constraining the EOS of dense neutron-rich nucleonic matter with astrophysical observations first in the 3-D space of , and , then in the 2-D plane of with . Finally, a summary is given in Section 4.

## 2 Modeling the equation of state of neutron star matter

The focus of this study is on the EOS of dense neutron-rich matter. We shall thus adopt the NV EOS (Negele & Vautherin, 1973) as that for the inner crust and the BPS EOS (Baym et al., 1971) for the outer crust while focusing on the EOS of the dense core in a large 3-D and 2-D parameter space. In this section, we shall first investigate the crust-core transition density and pressure using a thermal dynamical approach. Our main goal is to examine effects of the symmetry energy, especially its curvature on the transition point. We will then discuss how we construct the EOS for the core. For the purposes of this work, it is sufficient to use the simplest model for non-rotating and charge-neutral neutron stars consisting of only particles at -equilibrium. As we mentioned earlier, essentially all available many-body theories using various interactions have been used to predict both the EOS of SNM and symmetry energy from low to supra-saturation densities. Given their widely different predictions especially at supra-saturation densities, we try to extract parameters charactering the EOS of dense neutron-rich matter from the astrophysical observations using as little as possible predictions of any particular many-body theory and/or interaction. However, we must ensure that the EOS parameters asymptotically become naturally the ones extracted from terrestrial nuclear experiments near .

To avoid misleading the reader, here we discuss in more detail the relationship between the Taylor expansions of Eqs. (1) to (3) used to characterize the EOS near the saturation density in both nuclear theory and in analyzing terrestrial experiments and the parameterizations that we shall use in modeling the EOS of neutron star matter. We also address the question why we can’t just adopt the widely used isospin-independent multi-parameters polytropic EOSs, see, e.g., refs (Topper, 1964; Butterworth, 1976; Read et al., 2009), for the purposes of our work here. Near and , any nuclear energy density functional can be Taylor expanded using Eqs. (1) to (3). The relevant expansion coefficients are determined by the via for the symmetry energy, for the slope parameter of , and for the incompressibility of SNM and the curvature of , as well as and for the skewness of and , respectively. The Taylor expansions become progressively inaccurate for large densities and do not converge when . Therefore, for describing neutron star matter, we parameterize the and in exactly the same form as in Eqs. (1) to (3) but their coefficients are no longer given by the above expressions from any energy density functional . Instead, we treat them as unknown parameters to be extracted from astrophysical observations. More specifically, while the NV+BPS EOSs are used for the crust, starting from the crust-core transition density we parameterize the and as

 E0(ρ)=E0(ρ0)+KNS18(ρρ0−1)2+JNS162(ρρ0−1)3, (4)
 Esym(ρ)=Esym(ρ0)+LNS3(ρρ0−1)+KNSsym18(ρρ0−1)2+JNSsym162(ρρ0−1)3 (5)

using the parameters and . These parameterizations naturally become the Taylor expansions in the limit of . As parameterizations, mathematically they can be used at any density without the convergence issue associated with the Taylor expansions. Any parameterized EOS has to satisfy all the constraints (asymptotic boundary conditions) for the EOS near from terrestrial experiments. The above parameterizations facilitate comparisons with the terrestrial constraints around . When the above parameterizations are used at high densities for modeling the core EOS of neutron stars, the parameters in the above two equations are not necessarily to be the same as the Taylor expansion coefficients in Eqs. (1) to (3). While both mathematically and physically, the and are expected to be very close, the high-order coefficients (parameters) in the Taylor expansions (parameterizations) may be very different. Unfortunately, very little is known about these parameters from experiments/observations and the model predictions for the EOS of dense neutron-rich matter diverge. Since the uncertainty ranges of the Taylor expansion coefficients , and are already so large, it is then meaningful to use their uncertainty ranges as the starting/reference ranges in our search for the high-density EOS parameters and from studying properties of neutron stars.

With the above definitions and explanations, we now introduce/discuss a convention/simplification that we shall use in the following. It is consistent with the existing convention widely used in the literature for the Eq (1). The latter is often referred as the empirical parabolic law/approximation of the EOS of isospin asymmetric nuclear matter (Bombaci & Lombardo, 1991). Namely, it is well known in nuclear physics that the Eq (1) has the dual meanings of being a Taylor expansion in in the limit of on one hand, and on the other hand being a parameterization when it is used in very neutron-rich matter where . Adopting this convention and be consistent with that used for the Eq (1), albeit confusing without the above explanations, it is really not necessary to write both the Taylor expansions of Eqs. (1) to (3) and the parameterizations of Eqs. (4) and (5) as they have the same form although different meanings for the coefficients/parameters involved. In the following discussions, when the Eqs. (1) to (3) are used in describing terrestrial nuclear EOS or predictions of nuclear energy density functional theories, they are Taylor expansions near and . On the other hand, when they are used in describing the EOS of neutron stars, they are parameterizations to be constrained by astrophysical observations. Then, all three equations in Eqs. (1) to (3) have the dual meanings in describing the EOS of isospin asymmetric matter encountered in both terrestrial experiments and neutron stars within broad ranges of and without the convergence problem.

In the sense that because there is no reliable theory for the EOS of neutron-rich matter significantly above the saturation density one has to use parameterizations to describe the EOS of neutron star matter, the spirit of the parameterizations of Eqs. (4) and (5) is essentially the same as many other parameterizations used in the literature, see, e.g., refs. Gandolfi et al. (2009, 2012); Steiner et al. (2016). However, instead of requiring indirectly the involved parameters to reproduce the known properties of symmetric nuclear matter and the symmetry energy near through sometimes complicated equations, we use directly the coefficients of the Eqs. (2) and (3). These coefficients as the asymptotic boundary conditions of the EOS are themselves known near and either through experiments or converged predictions of many reliable models, see, e.g., refs. (Danielewicz et al., 2002; Shlomo et al., 2006; Chen et al., 2009; Steiner et al., 2010; Piekarewicz, 2010; Khan et al., 2012; Dutra et al., 2012, 2014; Li & Han, 2013; Colò et al., 2014; Cai & Chen, 2014; Zhang et al., 2017). In fact, parameterizations similar to the Eqs. (2) and (3), albeit often truncated at the first order in density for and second order for , i.e., using and only, have already been used successfully in studying various properties of neutron stars, see, e.g., refs.(Oyamatsu & Iida, 2007; Sotani et al., 2012). Also, a similar approach of describing approximately the EOS of dense neutron-rich matter was recently proposed in ref. (Margueron et al., 2017a) and successfully used in studying properties of both neutron stars (Margueron et al., 2017b) and finite nuclei (Chatterjee et al., 2017) with Bayesian perspectives. With the above cautions, conventions and justifications, we shall use the parameterized EOS for the core of neutron stars as described in more detail in section 2.3.

The multi-parameters polytropic EOSs, see, e.g., refs (Topper, 1964; Butterworth, 1976; Read et al., 2009) are widely used in modeling the core EOS of neutron stars at supra-saturation densities. Why do not we simply follow this popular approach? This is mainly because the polytropes used so far depend only on the nucleon number (energy) density but explicitly independent of the isospin asymmetry . To extract information about the high-density symmetry energy from neutron stars, we need to know explicitly the underlying isospin composition of the pressure/density. For example, to our best knowledge, all Bayesian analyses done so far using masses and radii of neutron stars have adopted the polytropes. They infer only the total pressures at a few fiducial high densities without any information about the isospin content of the matter at those densities. Indeed, to invert mathematically the TOV equation one only needs to know the total pressure as a function of density regardless of its microphysics composition. However, to understand the entire structure, composition and cooling mechanism of neutron stars, it is necessary to know the isospin dependence of its EOS from the crust to the dense core consistently. Since the pressure has a term proportional to where the isospin-asymmetry profile in neutron stars at -equilibrium is determined uniquely by the as we shall discuss in more detail below, by directly parameterizing the and , separately, and then evaluate the resulting pressure as a function of both and , we shall obtain information about the high-density underlying the pressure-density relation . While directly parameterizing the do not facilitate the extraction of such information about the high-density symmetry energy from astrophysical observations.

### 2.1 Sampling the density dependence of nuclear symmetry energy

First of all, we illustrate the broad variation of by varying independently the coefficients in Eq. (3) within their known uncertainty ranges. It is well known that the isospin asymmetry in neutron stars at a given baryon density is uniquely determined by the in Eq. (1) through the charge neutrality and equilibrium conditions as we shall recall more formally in section 2.3. Generally speaking, because of the term in the EOS, a higher value of will lead to a smaller at equilibrium. As quantitative examples, shown in Figure 1 are the and as functions of reduced density by varying only one coefficient each time while fixing all others: (a) , 50, 60, 70, and 80 MeV, (b) , -300, -200, -100, 0, and 100 MeV, and (c) , 0, 200, and 400 MeV. As their names indicate, the slope , curvature and skewness of symmetry energy play different roles and in order become increasingly more important at higher densities. Obviously, variations of them within their currently known uncertainty ranges allow us to sample very different behaviors of the and the corresponding .

It is worth noting that some combinations of the parameters lead to a decreasing that may even become negative at high densities. As summarized earlier in Szmaglinski et al. (2006) and reviewed very recently in Li et al. (2018), such kind of super-soft at high densities was predicted in a number of theoretical calculations using various interactions. At very high densities, when the short-range repulsive tensor force due to the -meson exchange makes the EOS of SNM increase faster with density than that of pure neutron matter where the tensor force is much weaker, the decreases or even becomes negative at high densities (Pandharipande et al., 1972; Wiringa et al., 1988; Li et al., 2008). To our best knowledge, such a seemingly unusual high-density behavior of the is not excluded by neither any known fundamental physics principle nor experiments/observations so far. Possible consequences of such kind of symmetry energies are discussed in more detail in refs. (Szmaglinski et al., 2006; Li et al., 2018) and references therein. In fact, EOSs with a super-soft can still support massive neutron stars if the SNM parts of the EOSs are sufficiently stiff even without the help from the new light mesons proposed and/or possible modified strong-field gravity for massive objects, see, e.g., refs. (Krivoruchenko et al., 2009; Wen et al., 2009; Lin et al., 2014; Jiang et al., 2015). Interestingly, while not completely settled yet (Li, 2017), there are indeed some circumstantial evidences from intermediate-relativistic energy heavy-ion collisions indicating that the  may become super-soft above about (Xiao et al., 2009). Currently, devoted efforts are being made by the intermediate-relativistic heavy-ion reaction community to pin down the high-density behavior of nuclear symmetry energy, see, e.g., refs. (Li et al., 2014; Russotto et al., 2016; Xu et al., 2016; Tsang et al., 2017; Trautmann and Wolter, 2017).

As we shall discuss in detail next, the around the crust-core transition is mostly controlled by the and parameters when the is fixed at its most probable empirical value of 31.6 MeV. The variation of from 40 to 80 MeV and from -400 to 100 MeV allows us to sample the usual behavior of predicted by various nuclear many-body theories in the sub-saturation density region and within their known empirical constraints at .

### 2.2 Imprints of the density dependence of nuclear symmetry energy on the crust-core transition in neutron stars

Although the crust of neutron stars contributes only a small fraction of the total mass and radius, it plays an important role in various astrophysical phenomena (Baym et al., 1971, 1971b; Pethick & Ravenhall, 1995; Link et al., 1999; Lattimer & Prakash, 2000, 2001; Steiner et al., 2005; Lattimer & Prakash, 2007; Chamel & Haensel, 2008; Sotani et al., 2012; Newton et al., 2013; Pons et al., 2013; Piekarewicz et al., 2014; Horowitz et al., 2015). Critical to many effects of the crust is the transition density between the inner crust and the outer core of neutron stars. Previous studies have found consistently that the transition density is very sensitive to the density dependence of nuclear symmetry energy. In particular, the role of the slope parameter has been extensively studied, see, e.g., refs. (Douchin & Haensel, 2000; Kubis, 2004, 2007; Lattimer & Prakash, 2007; Oyamatsu & Iida, 2007). Often the studies employ predictions of a particular nuclear many-body theory where the values of and are normally correlated. Here we shall first study the individual roles of the and in determining the core-crust transition properties, then contours of the transition density and pressure in the versus plane. Finally, effects of the correlation based on the systematics from analyzing over 500 nuclear energy density functions are examined.

In the present study, we employ the thermodynamical approach (Kubis, 2004, 2007; Lattimer & Prakash, 2007) to estimate the crust-core transition point where the uniform matter becomes unstable against being separated into a mixture of single nucleons and their clusters. The method is known to slightly overestimate the transition density compared to the dynamical approach (Xu et al., 2009; Ducoin et al., 2011; Providência et al., 2014) but sufficiently good for the purposes of this work. Specifically, the transition density is determined by the vanishing effective incompressibility of neutron star matter at equilibrium under the charge neutrality condition (Kubis, 2004, 2007; Lattimer & Prakash, 2007), i.e.,

 Kμ=ρ2d2E0dρ2+2ρdE0dρ+δ2[ρ2d2Esymdρ2+2ρdEsymdρ−2E−1sym(ρdEsymdρ)2]=0. (6)

This approach has been used extensively in the literature to locate the transition point using various EOSs (see, e.g., refs. Xu et al., 2009; Ducoin et al., 2011; Providência et al., 2014; Routray et al., 2016). Enclosed in the bracket of the above expression for are the first-order and second-order derivatives of the symmetry energy, i.e., quantities directly determining the and . It is thus necessary and interesting to first explore separate roles of the latter on the transition density. Shown in Figure 2 are the transition density as functions of with different values of in the window-a and with different values of in the window-b, respectively. It is clearly seen that the changes more dramatically with the variation of than in their respective uncertainty ranges. This is mainly because the last two terms in the expression for largely cancel out, leaving the curvature of dominate. In addition, the value of is already relatively well constrained in a smaller range than the , making the variation of with look weaker.

Next, we examine the transition density and pressure by varying both the and within their uncertainty ranges continuously. Shown in the two windows of Figure 3 are contours of constant transition densities and pressures in the plane, respectively. For transition densities larger than fm, the required increases monotonically with , while different behaviors are observed for lower values of . The lowest transition density about fm appears around the boundary corner at MeV and MeV. Different from the contours of constant , for a fixed , the required always increases linearly with before reaching the boundary along the line . The latter is used as a limit in exploring properties of neutron stars in the EOS parameter space.

As indicated earlier, in using the Eqs. (2) and (3) we assume the coefficients are independent and intend to constrain them directly from observations using as little as possible predictions of any particular many-body theory. Thus, in determining the crust-core transition point for constructing the EOS of neutron stars, we freely vary the and within their uncertainty ranges specified earlier. Nevertheless, theoretically predicted values of the and are often correlated when model ingredients and/or interactions are varied. Thus, for a consistency check we also study how the predicted correlation between and may limit the transition point. As discussed by many people in the literature (see, e.g., refs. Farine et al., 1978; Ducoin et al., 2011; Providência et al., 2014; Pearson et al., 2014; Danielewicz & Lee, 2009; Chen et al., 2009; Vidaña et al., 2009; Ducoin et al., 2011; Mondal et al., 2017; Tews et al., 2017), based on the systematics of many predictions using various many-body theories and interactions, an approximately universal and linear correlation exists between the and . For example, using totally over 500 energy density functions including 263 Relativistic Mean Field (RMF) models, Hartree-Fok calculations using 240 Skyrme (Dutra et al., 2012, 2014) as well as some realistic and Gogny interactions, Mondal et al. (2017) found the following correlation

 Ksym=(−4.97±0.07)(3Esym(ρ0)−L)+66.8±2.14  MeV. (7)

Using essentially the same sets of energy density functionals but requiring fm, MeV, MeV, and MeV, Tews et al. (2017) rejected some of the energy functionals. They found the following correlation using the remaining 188 Skyrme and 73 RMF interactions

 Ksym=3.501L−305.67±24.26 (±56.59)  MeV, (8)

where the and are error bars including 68.3% and 95.4% of the accepted EOSs, respectively. In Figure 3, the above two correlation functions are shown with the white and red dashed lines, separately. The 68.3% uncertainty range ( MeV) in Eq. (8) is indicated by the red dotted lines while that of Eq. (7) is too small to be shown here. While the parameterizations in Eqs. (7) and (8) are consistent, their mean values have slightly different slopes because of the different selection criteria used. Nevertheless, if the uncertainty range of Eq. (8) is enlarged from 68.3% to 95.4% of the accepted EOSs, the parameterization in Eq. (7) can then be fully covered by that in Eq. (8). Using the above two correlations, the transition density and pressure are then restricted to be about fm and MeV fm, respectively. These values are consistent with the crust-core transition properties often used in the literature. To this end, especially since some of the apparent correlations among EOS parameters from model calculations may be spurious (Margueron et al., 2017a), it is worth noting that while the correlation from the systematics of over 500 energy density functions is very useful for the consistency check, it is still necessary and important to determine the individual values of and from experiments/observations. In our following calculations, we thus use consistently the crust-core transition density and pressure by varying independently the and values within their respective uncertain ranges without using any of the about two correlation functions.

### 2.3 The core EOS of neutron stars

For completeness and the ease of our discussions in the following, we first recall here the formalism for calculating the EOS in the cores of neutron stars. The total energy density of charge neutral matter at -equilibrium can be written as

 ϵ(ρ,δ)=ϵb(ρ,δ)+ϵl(ρ,δ), (9)

where and are the energy density of baryons and leptons, respectively. The can be calculated from

 ϵb(ρ,δ)=ρE(ρ,δ)+ρMN, (10)

where the specific energy of baryons is given in Eq. (1) and is the average rest mass of nucleons. The from the noninteracting Fermi gas model can be expressed as () (Oppenheimer & Volkoff, 1939)

 ϵl(ρ,δ)=ηϕ(t) (11)

with

 η=m4l8π2,  ϕ(t)=t√1+t2(1+2t2)−ln(t+√1+t2), (12)

and

 t=(3π2ρl)1/3ml. (13)

The chemical potential of particle can be calculated from

 μi=∂ϵ(ρ,δ)∂ρi. (14)

The isospin asymmetry and relative particle fractions at different densities in neutron stars are obtained through the -equilibrium condition and the charge neutrality condition The pressure of the system can be calculated numerically from

 P(ρ,δ)=ρ2dϵ(ρ,δ)/ρdρ. (15)

The above expressions allow us to calculate the core EOS which is connected smoothly at the core-crust transition point to the NV EOS (Negele & Vautherin, 1973) for the inner crust followed by the BPS EOS (Baym et al., 1971) for the outer crust.

As mentioned earlier, the EOS parameters , and near are relatively well determined. To investigate how/what high-density EOS parameters are constrained by the three astrophysical observations considered in this work, we construct the EOS of neutron star matter by varying the poorly known , and characterizing the EOS of dense neutron-rich nucleonic matter. In principle, all coefficients used in Eqs. (2) and (3) should be varied simultaneously within a multivariant Bayesian inference (Steiner et al., 2010; Raithel et al., 2016; Margueron et al., 2017b). Such a study is in progress. In the present work, we shall perform traditional analyses first in the 3-D parameter space spanned by , and while fixing all other parameters at their currently known most probable values. Equivalent to re-parameterizing the EOS of SNM and with less parameters as often done in the literature, or expanding the Eqs. (2) and (3) only up to , we shall also explore the EOS in the 2-D parameter space of and by setting while keeping other parameters at their known most probable values. The two cases studied here are similar in spirit to using different numbers of piecewise polytropes or parameters to model the EOS of dense neutron-rich matter. Naturally, the values of the parameters involved may be different in the two cases, but they should all asymptotically approach the same existing constraints on them near .

Certainly, there have been continuous efforts in both astrophysics and nuclear physics to constrain the EOS parameters in both cases. For example, from the pressure of SNM constrained by nucleon collective flow data in relativistic heavy-ion collisions (Danielewicz et al., 2002), a constraint of MeV was obtained by Cai & Chen (2014). Combining it with the mass of neutron star PSR J0348+0432, they further narrowed it down to MeV. By analyzing X-ray bursts, Steiner et al. (2010) extracted a value of or MeV assuming a photospheric radius of or , respectively. All of these constraints on overlap but have different uncertainty ranges. Similarly, the and have not been well constrained either by any experiments/observations so far. Nevertheless, the systematics of over 500 RMF and SHF energy density functionals indicates the following range: MeV (Dutra et al., 2012, 2014), MeV and MeV, respectively (see, e.g., Chen et al., 2009; Dutra et al., 2012, 2014; Colò et al., 2014; Zhang et al., 2017). Thus, we adopt these ranges for the parameters , and to be consistent with both existing experimental and theoretical findings.

Within the above uncertainty ranges of the EOS parameters, the pressure in neutron stars can be varied within the shaded band shown in the left window of Figure 4. Its upper and lower limit is obtained by using the parameter set of MeV, MeV, MeV and MeV and the set of MeV, MeV, MeV and MeV, respectively. The individual roles of these parameters are examined by varying them independently in the right window. As one expects, the variation of is most effective in modifying the pressure at supra-saturation densities.

## 3 Constraining the EOS of dense neutron-rich matter with observed properties of neutron stars

With the EOSs prepared in the way described above, the mass (M)-radius (R) relationship of neutron stars is obtained by solving the Tolman-Oppenheimer-Volkov (TOV) equations (Tolman, 1934; Oppenheimer & Volkoff, 1939)

 dPdr=−G(m(r)+4πr3P/c2)(ϵ+P/c2)r(r−2Gm(r)/c2), (16)
 dm(r)dr=4πϵr2 (17)

where is the gravitation constant, is the light speed and is the gravitational mass enclosed within a radius . The dimensionless tidal deformability is related to the Love number , neutron star mass M and radius R via

 Λ=23k2⋅(R/M)5. (18)

The is determined by the EOS thorough a differential equation coupled to the TOV equation (Hinderer, 2008; Hinderer et al., 2010). More details about the formalism and code used in this work to calculate the can be found in, e.g., Fattoyev et al. (2013, 2014).

Within the 3-D parameter space of , and and the 2-D parameter space of and under the conditions discussed in the previous section, using the observational data of M, km and , we study how/what high-density EOS parameters are constrained in the following two subsections, separately.

### 3.1 Observational constraints in the J0−Ksym−Jsym EOS parameter space

To be clear, we first emphasize again that in this case the following parameters are fixed at their currently known most probable values based on previous systematic surveys as discussed earlier: MeV, MeV and MeV. We explore constraints on the EOS in the 3-D parameter space with the crust-core transition density consistently determined and the condition that . Technically, in exploring the 3-D parameter space in three loops we change the to reproduce a specific observational data at given values of and which are varied independently. Shown in Figure 5 are the constant surfaces of neutron stars’ maximum mass of M (green), radius of km (magenta) and km (yellow) as well as the upper limit of the dimensionless tidal deformability (orange) of canonical neutron stars, respectively. For clarity, a causality surface limiting M M is not shown.

For ease of the following discussions, it is worth recalling first that, as shown in Eq. (1), the EOS of SNM , the symmetry energy and the isospin asymmetry at equilibrium are the three quantities together determining the total pressure in neutron stars. More specifically, the total pressure is proportional to . While the term is controlled by the parameter with a fixed , the and are determined by the and parameters when the is fixed. Moreover, the symmetry energy contribution to the total pressure is weighted by . When the is softer with smaller or negative and values, the system is more neutron-rich as shown in Figure 1. In particular, for extremely small and values (e.g., MeV and MeV), the becomes negative and the reaches its maximum of 1 at high densities. Then, the necessary contribution to the pressure from the term will require a large value to support massive neutron stars. At the other extreme, however, when both the and are strongly positive (e.g., MeV and MeV at the bottom left corner), the symmetry energy is super-stiff and the is very small as shown in Figure 1. The required to support massive neutron stars is then very small.

It is interesting to note several major features in this rather information-rich 3-D plot summarizing very extensive calculations. Let us first focus on the constant surface of M. From the top right to the bottom left corner, the required first decreases quickly and then stays almost a constant with the increasing values of and from negative to positive. This feature is completely understandable based on the discussions in the previous paragraph. Namely, near the upper right corner, the symmetry energy is super-soft and the resulting is close to 1. The weight of the symmetry energy contribution to the pressure is significant while the value is small and may even be largely negative. To support neutron stars with the maximum mass of M, the value of has to be highly positive to compensate the small or overcome the possibly negative contribution from the symmetry energy. However, near the bottom left corner where the symmetry energy is super-stiff, the resulting at equilibrium becomes so small such that the symmetry energy contribution to the pressure is strongly suppressed. The necessary value of is therefore small and the constant surface of M becomes rather flat at small values. More quantitatively, the required minimum value of is about MeV at MeV and MeV.

For a comparison and see more clearly the role of parameter in determining the masses of neutron stars, a constant surface of M (sky blue) corresponding to the composite mass of GW170817 is also shown. It is seen that the overall shapes of the constant surfaces of M and 2.74 M are rather similar. To support such a hypothetical neutron star would require a value above 400 MeV beyond its current uncertain range and the causality limit. While the fate of GW170817 is not completely determined yet, several analyses using observations of GW170817 combined with theories/simulations and the causality limit under some caveats have placed the upper bounds on neutron star masses in the range of M (Margalit & Metzger, 2017; Shibata et al., 2017; Rezzolla et al., 2018; Ruiz et al., 2018), see, e.g., ref. (Radice et al., 2018) for a very recent review. For instance, Margalit & Metzger (2017) used electromagnetic constraints on the remnant imposed by the kilonova observations after the merger and the gravitational wave information predicted a maximum mass of M with 90% confidence. To support neutron stars with such masses, the just needs to be slightly positive well within its current uncertain range. Thus, once the maximum mass is pinned down, it will put a stringent upper limit on the parameter.

In addition, it is also known that quadrupole deformations of neutron stars depend on the symmetry energy (Krastev et al., 2008a, b; Fattoyev et al., 2013, 2014, 2017; Krastev et al., 2018). Interestingly, Abbott et al. (2017) inferred that the dimensionless tidal deformability of GW170817 has an upper limit of at the 90% confidence level for the low-spin prior. However, it is seen in Figure 5 that the constant surface of (orange) locates far outside the constant surface of km. Thus, limits on the high-density EOS parameters from the constraint alone are presently much looser than the radius constraint extracted from analyzing the X-ray data. Nevertheless, the expected detection of gravitational waves from a large number of neutron star mergers has the potential to improve the situation.

Now, let us move to the EOS constraints from the radii of neutron stars. The radius of canonical neutron stars is known to depend strongly (weakly) on the nuclear symmetry energy (EOS of SNM) (Li & Steiner, 2006). It is thus not surprising that the two constant surfaces of radius at km and km for canonical neutron stars are essentially vertical in Figure 5, indicating a weak dependence on the as one expects. Indeed, they have significant dependences on both the and as indicated by the separation between the two constant-radius surfaces. More quantitatively, the required values of in the constant surfaces of km and km decrease continuously with increasing and . The two surfaces can be approximately described by and , respectively. With a fixed value of , the nuclear pressure becomes stronger with increasing values of and . Thus, the constant surface of km is on the left (having stiffer symmetry energies) of the km surface.

As indicated by the arrows in Figure 5, the space enclosed by the three constant surfaces of M and km are the EOS parameter space allowed by the astrophysical observations of the maximum mass and radii of neutron stars. The space is constrained from the bottom by the maximum mass, its intersections with the two limits on the radii set the boundary lines restricting the and at small values around MeV. At higher values, the EOS parameter space is mainly bounded by the two radius constraints. To show the accepted EOS parameters more clearly, the allowed regions in the versus planes with , , , and MeV, respectively, are shown in Figure 6. The boundaries of these allowed regions in the planes are shown in Figure 7. The regions insides the lines are allowed for a given specified. The shadowed regions are excluded values of and for all values. More specifically, the boundary of the right shadowed region can be divided into two parts based on the slopes. The upper one with a negative slope is obtained by requiring the crust-core transition pressure to always stay positive, i.e, . It can be described approximately by the expression MeV. The lower one with a positive slope is obtained by the intersection line between the surfaces of M and km. It can be fitted by the expression MeV. This boundary sets an upper limit for at about MeV. The left shadowed region is excluded by the intersection line of km and MeV in Figure 5. It can be fitted by the expression MeV.

Overall, within the framework of our analyses, using the maximum mass of neutron stars as well as the upper and lower limits of the radii of canonical neutron stars, the three EOS parameters are only limited in a space shown in Figure 5 not completely closed with its boundaries partially given in Figure 7. Obviously, data of more independent observables from either or/both astrophysics and nuclear physics are needed to determine the individual values of , and .

### 3.2 Observational constraints in the L−Ksym EOS parameter plane

Within the currently known uncertainty ranges of and , one can parameterize the EOS with less parameters by setting in Eqs. (2) and (3) as often done in the literature. Then, the two most poorly known parameters are the and . The latter is known to be around MeV as mentioned earlier. In this 2-D model framework, here we explore how/what the same three astrophysical constraints may teach us about the EOS of dense neutron-rich matter.

Shown in Figure 8 are contours of constant maximum masses of neutron stars in the plane. The red, magenta, blue, and cerulean shadows represent regions where M, , km, and , respectively. In the white excluded region, the crust-core transition pressure . It is seen that the maximum mass of neutron stars increases with increasing as one expects. However, it is rather insensitive to the variation of in the region considered. Again and obviously, the tidal polarizability is much less restrictive than the radius constraint of km. It is also seen that the km constraint is covered by the constraint M. Consequently, the area bounded by the curves of M and km is the region allowed by the astrophysical observations considered. In this region, the maximum value of is about 52 MeV consistent with that found in the 3-D analyses in the previous subsection. Also in this region, the upper limit of the maximum mass is about 2.37 M reached at MeV. Probably incidentally, the latter is also currently the most probable value of based on the surveys of 53 analyses of existing data (Li & Han, 2013; Oertel et al., 2017). Interestingly, the observed upper limit of the maximum mass is in good agreement with the findings by Fryer et al. (2015), Lawrence et al. (2015) and Alsing et al. (2017). They found that the upper limit of the maximum mass of neutron stars is between and M. However, it is necessary to caution here that we have taken in the 2-D study. As we have discussed in the previous subsection, the value of affects significantly the maximum mass of neutron stars. Thus, without more precise knowledge about the parameter, the absolute maximum mass of neutron stars can not be pinned down. Based on our 2-D model analyses here, neutron stars more massive than 2.37 M would require a positive value of .

## 4 Summary and outlook

In summary, within both the 2-D and 3-D EOS parameter spaces limited by the existing constraints from terrestrial nuclear experiments, we studied how the astrophysical observations of M, km and all together constrain the EOS parameters of dense neutron-rich nucleonic matter. We also investigated effects of the curvature of nuclear symmetry energy on the crust-core transition in neutron stars. The consistently calculated transition density in the plane is used in constructing the EOS of neutron star matter from the surface to the core. The is found to affect significantly the crust-core transition density and pressure. Fixing the , and at their most probably values determined mainly by terrestrial nuclear experiments, we explored the intersections of constant surfaces with , km, km, and , respectively, in the 3-D parameter space of , and . The 3-D parameter space narrowed down significantly by the observational constraints is clearly identified. This helps guide nuclear theories for dense neutron-rich matter and related studies in terrestrial experiments. However, to pin down the individual values of , and , data of additional independent observables from either astrophysical observations and/or laboratory experiments are needed. In particular, the skewness parameter of SNM largely controls the maximum mass of neutron stars. The 2-D EOS with is found sufficiently stiff to support neutron stars as massive as 2.37 M, while to support a hypothetical neutron star as massive as 2.74 M (composite mass of GW170817) would require to be larger than about 400 MeV beyond its known limit. In both the 2-D and 3-D model frameworks considered in this work, the upper limit of the tidal deformability from the recent observation of GW170817 is found to provide upper limits on some EOS parameters consistent with but less restrictive than the existing constraints. In particular, its constraints on the symmetry energy parameters are far less restrictive than the observation of km from analyzing the X-ray data.

While the analyses and results presented here are useful in their own rights, some aspects of our present work should be improved in future studies. In particular, more data and/or a better approach are necessary to determine precisely individual values of the high-density EOS parameters with quantified uncertainties. In the era of gravitational wave astronomy accompanied by the planned new experiments using advanced radioactive beam facilities around the world, we are very hopeful that more data will come soon. Moreover, even with the limited data available, more quantitative information about the high-density EOS parameters may be obtained by using a more robust statistical approach. For example, we fixed the parameters describing the EOS near the saturation density at their most probable values mostly from terrestrial experiments. By doing so we eliminated any possible correlation between these (fixed) parameters and the high-density parameters , and that we want to determine from astrophysical observations. A Bayesian analysis where the prior distribution of parameters encodes what is already known and the new observational data refines such prior knowledge is particularly well suited for inferring the posterior probability density distributions of all the EOS parameters. Such a study using the masses and radii of 14 neutron stars extracted from Chandra X-Ray observations, the maximum masses of neutron stars from GW170817 and earlier observations of neutron stars as well as all the information from terrestrial experiments about the low-order EOS parameters will be carried out and reported elsewhere.

We thank F. J. Fattoyev for providing us the code to calculate the dimensionless tidal deformability. We would also like to thank Alessandra Corsi, Benjamin Owen and Renxin Xu for very helpful discussions as well as Ang Li and Bing Zhang for useful communications. NBZ is supported in part by the China Scholarship Council. BAL acknowledges the U.S. Department of Energy, Office of Science, under Award Number DE-SC0013702, the CUSTIPEN (China-U.S. Theory Institute for Physics with Exotic Nuclei) under the US Department of Energy Grant No. DE-SC0009971 and the National Natural Science Foundation of China under Grant No. 11320101004. JX is supported in part by the Major State Basic Research Development Program (973 Program) of China under Contract Nos. 2015CB856904 and 2014CB845401, the National Natural Science Foundation of China under Grant Nos. 11475243 and 11421505, the “100-talent plan” of Shanghai Institute of Applied Physics under Grant Nos. Y290061011 and Y526011011 from the Chinese Academy of Sciences, and the Shanghai Key Laboratory of Particle Physics and Cosmology under Grant No. 15DZ2272100.

## References

• Abbott et al. (2017) Abbott, B. P., et al. 2017, Phys. Rev. Lett., 119, 161101
• Alsing et al. (2017) Alsing, J., Silva, H. O., & Berti, E. 2017, arXiv:1709.07889
• Antoniadis et al. (2013) Antoniadis, J., et al. 2013, Science, 340, 448
• Baym et al. (1971) Baym, G., Pethick, C. J., & Sutherland, P. 1971, APJ, 170, 299
• Baym et al. (1971b) Baym, G., Bethe, H. A., & Pethick, C. J. 1971b, Nucl. Phys. A, 175, 225
• Bombaci & Lombardo (1991) Bombaci, I., & Lombardo, U. 1991, Phys. Rev. C44, 1892
• Butterworth (1976) Butterworth, E. M. 1976, APJ, 204, 561
• Cai & Chen (2014) Cai, B. J., & Chen, L. W. 2017, Nucl. Sci. Tech., 28, 185
• Chamel & Haensel (2008) Chamel, N., & Haensel, P. 2008, Living Rev. Relat., 11, 10
• Chatterjee et al. (2017) Chatterjee, D., Gulminelli, F., Raduta, Ad. R. & Margueron, J. 2017, Phys. Rev. C, 96, 065805.
• Chen et al. (2009) Chen, L. W., Cai, B. J., Ko, C. M., Li, B. A., Shen, C., & Xu, J. 2009, Phys. Rev. C, 80, 014322
• Colò et al. (2014) Colò, G., Garg, U., & Sagawa, H. 2014, Eur. Phys. J. A, 50, 26
• Danielewicz et al. (2002) Danielewicz, P., Lacey, R., & Lynch, W. G. 2002, Science, 298, 1592
• Danielewicz & Lee (2009) Danielewicz, P., & Lee, J. 2009, Nucl. Phys. A, 818, 36
• Demorest et al. (2010) Demorest, P. B., Pennucci, T., Ransom, S. M., Roberts, M. S. E., & Hessels, J. W. T. 2010, Nature, 467, 1081
• Douchin & Haensel (2000) Douchin, F., & Haensel, P. 2000, Phys. Lett. B, 485, 107
• Ducoin et al. (2011) Ducoin, C., Margueron, J., Providência, C., & Vidaña, I. 2011, Phys. Rev. C, 83, 045810
• Dutra et al. (2012) Dutra, M., Loureno, O., Martins, J. S. S., Delfino, A., Stone, J. R., & Stevenson, P. D. 2012, Phys. Rev. C, 85, 035201
• Dutra et al. (2014) Dutra, M., Loureno, O., Avancini, S. S., Carlson, B. V., Delfino, A., Menezes, D. P., Providncia, C., Typel, S., & Stone, J. R. 2014, Phys. Rev. C, 90, 055203
• Fattoyev et al. (2013) Fattoyev, F. J., Carvajal, J., Newton, W. G., & Li, B. A. 2013, Phys. Rev. C, 87, 015806
• Fattoyev et al. (2014) Fattoyev, F. J., Newton, W. G., & Li, B. A. 2014, Eur. Phys. J. A, 50, 45
• Fattoyev et al. (2017) Fattoyev, F. J., Piekarewicz, J. & Horowitz, C.J., 2017, arXiv:1711.06615.
• Farine et al. (1978) Farine, M., Pearson, J. M., & Rouben, B. 1978, Nucl. Phys. A, 304, 317
• Fryer et al. (2015) Fryer, C. L., Belczynski, K., Ramirez-Ruiz, E., Rosswog, S., Shen, G., & Steiner, A. W. 2015, APJ, 812, 24
• Gandolfi et al. (2009) Gandolfi, S., Illarionov, A. Yu., Schmidt, K. E., Pederiva, F., & Fantoni, S. 2009, Phys. Rev. C, 79, 054005
• Gandolfi et al. (2012) Gandolfi, S., Carlson, J., & Reddy, S. 2012, Phys. Rev. C, 85, 032801
• Grigorian et al. (2016) Grigorian, H., Voskresensky, D. N. & Blaschke, D. 2016, Eur. Phys. J. A, 52, 67
• Hessels et al. (2006) Hessels, J. W. T., et al. 2006, Science, 311, 1901
• Hinderer (2008) Hinderer, T. 2008, APJ., 677, 1216
• Hinderer et al. (2010) Hinderer, T., Lackey, B. D., Lang, R. N., & Read, J. S. 2010, Phys. Rev. D, 81, 123016
• Horowitz et al. (2015) Horowitz, C. J., et al. 2015, Phys. Rev. Lett., 114, 031102
• Jiang et al. (2015) Jiang, W.Z., Li, B.A. & Fattoyev, F.J. 2015, Eur. Phys. J. A, 51, 119
• Khan et al. (2012) Khan, E., Margueron, J., & Vidaña, I. 2012, Phys. Rev. Lett., 109, 092501
• Krastev et al. (2008a) Krastev, P. G., Li, B.A. and Worley, A. 2008, Phys. Lett. B, 668, 1
• Krastev et al. (2008b) Krastev, P. G., Li, B.A. and Worley, A. 2008, APJ, 676, 1170
• Krastev et al. (2018) Krastev, P. G., Li, B.A., 2018, arXiv:1801.04620
• Kubis (2004) Kubis, S. 2004, Phys. Rev. C, 70, 065804
• Kubis (2007) Kubis, S. 2007, Phys. Rev. C, 76, 025801
• Krivoruchenko et al. (2009) Krivoruchenko, M.I., Simkovic, F. & Faessler, A. 2009, Phys. Rev. D, 79, 125023
• Lattimer & Prakash (2000) Lattimer, J. M., & Prakash, M. 2000, Phys. Rep., 333, 121
• Lattimer & Prakash (2001) Lattimer, J. M., & Prakash, M. 2001, APJ, 550, 426
• Lattimer & Prakash (2007) Lattimer, J. M., & Prakash, M. 2007, Phys. Rep., 442, 109
• Lattimer & Steiner (2014) Lattimer, J. M., & Steiner, A. W. 2014, Eur. Phys. J. A, 50, 40
• Lattimer & Prakash (2016) Lattimer, J.M. & Prakash, M. 2016, Phys. Rep, 621, 127
• Lawrence et al. (2015) Lawrence, S., Tervala, J. G., Bedaque, P. F., & Miller, M. C. 2015, APJ, 808, 186
• Li & Steiner (2006) Li, B.A. & Steiner, A.W. 2006, Phys. Lett. B, 642, 436
• Li et al. (2008) Li, B.A., Chen, L.W. & Ko, C.M. 2008, Phys. Rep. 464, 113
• Li & Han (2013) Li, B. A., & Han, X. 2013, Phys. Lett. B, 727, 276.
• Li et al. (2014) Li, B. A., Ramos, À., Verde G., & Vidana, I. (Eds.) 2014, Eur. Phys. J. A, 50, 9
• Li (2017) Li, B. A. 2017, Nuclear Physics News, 27, 7
• Li et al. (2018) Li, B. A., Cai, B. J., Chen, L. W., & Xu, J. 2018, Progress in Particle and Nuclear Physics, 99, 29
• LIGO & Virgo (2017) LIGO Scientific Collaboration, Virgo Collaboration, et al. 2017, ApJL, 848, L12
• Link et al. (1999) Link, B., Epstein, R., & Lattimer, J. M. 1999, Phys. Rev. Lett., 83, 3362
• Lin et al. (2014) Lin, W.K., Li, B.A., Chen, L.W., Wen, D.H. & Xu, J. 2014, J. of Phys. G, 41, 075203
• Margalit & Metzger (2017) Margalit, B., & Metzger, B. D. 2017,APJ, 850, L19
• Margueron et al. (2017a) Margueron, J., Hoffmann Casali, R., Gulminelli, F., 2018, Phys. Rev. C 97, 025805
• Margueron et al. (2017b) Margueron, J., Hoffmann Casali, R., Gulminelli, F., 2018, Phys. Rev. C 97, 025806
• Mondal et al. (2017) Mondal, C., Agrawal, B. K., De, J. N., Samaddar, S. K., Centelles, M., & Viñas, X. 2017, Phys. Rev. C, 96, 021302(R)
• Miller & Lamb (2016) Miller, M.C. & Lamb, F. K. 2016, European Physical Journal A 52, 63.
• Negele & Vautherin (1973) Negele, J. W., & Vautherin, D. 1973, Nucl. Phys. A, 207, 298
• Newton et al. (2013) Newton, W. G., Gearheart, M., & Li, B. A. 2013, APJ Suppl. Ser., 204, 9
• Newton et al. (2014) Newton, W.G., et al., 2014, Euro Phys. J. A50, 41.
• Oertel et al. (2017) Oertel, M., Hempel, M., Klähn, T., & Typel, S. 2017, Rev. Mod. Phys., 89, 015007
• Oppenheimer & Volkoff (1939) Oppenheimer, J., & Volkoff, G. 1939, Phys. Rev., 55, 374
• Oyamatsu & Iida (2007) Oyamatsu, K., & Iida, K. 2007, Phys. Rev. C, 75, 015801
• Özel & Freire (2016) Özel, F. & Freire, P. 2016, Annual Reviews of Astronomy and Astrophysics, 54, 401
• Page et al. (2006) Page, D., Geppert, U., & Weber, F. 2006, Nucl. Phys. A, 777, 497
• Pandharipande et al. (1972) Pandharipande, V. R., et al. 1972, Phys. Lett. B, 39, 608
• Pearson et al. (2014) Pearson, J. M., Chamel, N., Fantina, A. F., Goriely, S. 2014, Eur. Phys. J. A, 50, 43
• Pethick & Ravenhall (1995) Pethick, C. J., & Ravenhall, D. G. 1995, Ann. Rev. Nucl. Part. Sci., 45, 429
• Piekarewicz (2010) Piekarewicz, J. 2010, J. Phys. G, 37, 064038
• Piekarewicz et al. (2014) Piekarewicz, J., Fattoyev, F. J., & Horowitz, C. J. 2014, Phys. Rev. C, 90, 015803
• Pons et al. (2013) Pons, J. A., Viganó, D., & Rea, N. 2013, Nat. Phys., 9, 431
• Providência et al. (2014) Providência et al., C. 2014, Eur. Phys. J. A, 50, 44
• Radice et al. (2018) Radice, D., Perego, A., Zappa, F., & S. Bernuzzi 2018, APJ, 852, L29
• Raithel et al. (2016) Raithel, C. A., Özel, F., & Psaltis, D. 2016, APJ, 831, 44; ibid, arXiv:1704.00737
• Read et al. (2009) Read, J. S., Lackey, B. D., Owen, B. J., & Friedman, J. L. 2009, Phys. Rev. D, 79, 124032
• Rezzolla et al. (2018) Rezzolla, L., Most, E. R., & Weih, L. R. 2018, APJ, 852, L25
• Ruiz et al. (2018) Ruiz, M., Shapiro, S. L., & Tsokaros, A. 2018, PRD, 97, 021501(R)
• Routray et al. (2016) Routray, T. R., Viñas, X., Basu, D. N., Pattnaik, S. P., Centelles, M., Robledo, L. B., & Behera, B. 2016, J. Phys. G, 43, 105101
• Russotto et al. (2016) P. Russotto et al. (ASY-EOS Collaboration) 2016, Phys. Rev. C, 94, 034608
• Shibata et al. (2017) Shibata, M., Fujibayashi, S., Hotokezaka, K., et al. 2017, Phys. Rev. D, 96, 123012
• Shlomo et al. (2006) Shlomo, S., Kolomietz, V. M., & Colò G. 2006, Eur. Phys. J. A, 30, 23
• Sotani et al. (2012) Sotani, H., Nakazato, K., Iida, K., & Oyamatsu, K. 2012, Phys. Rev. Lett., 108, 201101
• Steiner et al. (2005) Steiner, A. W., Prakash, M., Lattimer, J. M., & Ellis, P. J. 2005, Phys. Rep., 410, 325
• Steiner et al. (2010) Steiner, A. W., Lattimer, J. M., & Brown, E. F. 2010, APJ, 722, 33
• Steiner et al. (2016) Steiner, A. W., Lattimer, J. M., & Brown, E. F. 2016, Eur. Phys. J. A, 52,18
• Szmaglinski et al. (2006) Szmaglinski, A., Wojcik, W. & Kutschera, M. 2006, Acta Phys. Polon. B, 37, 277
• National Academies (2011) The National Academies, New Worlds, New Horizons in Astronomy and Astrophysics, 2011,
https://www.nap.edu/catalog/12951/new-worlds-new-horizons-in-astronomy-and-astrophysics
• National Academies (2012) The National Academies, Nuclear Physics: Exploring the Heart of Matter, 2012,
• U.S. LRP (2015) The 2015 U.S. Long Range Plan for Nuclear Science, Reaching for the Horizon,
• NuPECC LRP (2017) The Nuclear Physics European Collaboration Committee (NuPECC) Long Range Plan 2017, Perspectives in Nuclear Physics,
• Tews et al. (2017) Tews, I., Lattimer, J. M., Ohnishi, A., & Kolomeitsev, E. E. 2017, APJ, 848, 105
• Tsang et al. (2017) Tsang, M.B. et al. (SRIT Collaboration) 2017, Phys. Rev. C, 95, 044614
• Tolman (1934) Tolman, R. C. 1934, Proc. Natl. Acad. Sci. U.S.A., 20, 3
• Topper (1964) Topper, R. F. 1964, APJ, 140, 434
• Trautmann and Wolter (2017) Trautmann, W. & Wolter, H.H. 2017, arXiv:1712.03093
• Vidaña et al. (2009) Vidaña, I., Providência, C., Polls, A., & Rios, A. 2009, Phys. Rev. C, 80, 045806
• Watts et al. (2016) Watts, A.L. et al, 2016, Rev. Mod. Phys., 88, 021001
• Wen et al. (2009) Wen, D.H., Li, B.A., & Chen, L.W. 2009, Phys. Rev. Lett., 103, 211102
• Wiringa et al. (1988) Wiringa, R. B., Fiks, V., & Fabrocini, A. 1988, Phys. Rev. C, 38, 1010
• Xu et al. (2009) Xu, J, Chen, L. W., Li, B. A., & Ma, H. R. 2009, APJ, 697, 1549
• Xu et al. (2016) Xu, J. et al. (Transport Model Comparison Project) 2016, Phys. Rev. C, 93, 044609
• Xiao et al. (2009) Xiao, Z.G., Li, B.A., Chen, L. W., Yong, G. C. & Zhang, M. 2009, Phys. Rev. Lett., 102, 062502
• Yakovlev et al. (2001) Yakovlev, D. G., Kaminker, A. D., Gnedin, O. Y., & Haensel, P. 2001, Phys. Rep., 354, 1
• Zhang et al. (2017) Zhang, N. B., Cai, B. J., Li, B. A., Newton, W. G., & Xu, J. 2017, Nucl. Sci. Tech., 28, 181.
You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters

169779

How to quickly get a good answer:
• Keep your question short and to the point
• Check for grammar or spelling errors.
• Phrase it like a question
Test
Test description