Equation of state of dense matter in the multimessenger era

TL;DR

This paper combines nuclear physics experiments and astrophysical observations to tightly constrain the high-density symmetry energy and the equation of state of dense matter, impacting neutron star properties and gravitational wave signals.

Contribution

It provides new bounds on neutron star radius and tidal deformability by integrating diverse multimessenger data and nuclear constraints.

Findings

Lower bound on neutron star radius: R_{1.4} ≥ 11.1 km

Lower bound on tidal deformability: Λ_{1.4} ≥ 193

Upper bound on radius from GW data: R_{1.4} ≤ 13.3 km

Abstract

While the equation of state (EOS) of symmetric nuclear matter (SNM) at suprasaturation densities has been relatively well constrained from heavy-ion collisions, the EOS of high-density neutron-rich matter is still largely uncertain due to the poorly known high-density behavior of the symmetry energy. Using the constraints on the EOS of SNM at suprasaturation densities from heavy-ion collisions together with the data of finite nuclei and the existence of $2M_\odot$ neutron stars from electromagnetic (EM) observations, we show that the high-density symmetry energy cannot be too soft, which leads to lower bounds on dimensionless tidal deformability of $\Lambda_{1.4} \ge 193$ and radius of $R_{1.4} \ge 11.1$ km for $1.4M_\odot$ neutron star. Furthermore, we find that the recent constraint of $\Lambda_{1.4} \le 580$ from the gravitational wave signal GW170817 detected from the binary neutron star merger by the LIGO and Virgo Collaborations rules out too stiff high-density symmetry energy, leading to an upper limit of $R_{1.4} \le 13.3$ km. All these terrestrial nuclear experiments and astrophysical observations based on strong, EM and gravitational measurements together put stringent constraints on the high-density symmetry energy and the EOS of SNM, pure neutron matter and neutron star matter.

Equations

$$\Lambda_{M}=\frac{2}{3}k_{2}\bigg{(}\frac{c^{2}R}{GM}\bigg{)}^{5},$$

$$E(n,\delta )=E_0(n)+E_{\mathrm{sym}}(n){\delta }^2+\mathcal{O}({\delta }^4),$$

$$v_{i,j}$$

Full text

I Introduction

Dense matter with density comparable to nuclear saturation density $n_{0}$ ($\sim 0.16$ nucleon/fm${}^{3}\approx 2.7\times 10^{14}$ g/cm3) can exist in heavy atomic nuclei, in compact stars, or be produced in heavy-ion collisions. A basic model for understanding such dense matter is the nuclear matter - an ideal static infinite uniform system composed of nucleons (neutrons and protons) with only the strong interaction considered. One fundamental issue in nuclear physics, particle physics and astrophysics is to explore the equation of state (EOS) of nuclear matter Dan02 ; Lat04 ; Oer17 , conventionally defined as energy (or pressure) vs density. Because of the complicated nonperturbative feature of quantum chromodynamics (QCD), it is still a big challenge to determine the nuclear matter EOS from ab initio QCD calculations, especially at suprasaturation densities Bra14 . Therefore, data from terrestrial experiments or astrophysical observations are particularly important to constrain the nuclear matter EOS.

Indeed, the EOS of symmetric nuclear matter (SNM) with equal fraction of neutrons and protons has been relatively well constrained from around $n_{0}$ to about $5n_{0}$ by analyzing the data on giant monopole resonance of heavy nuclei You99 ; Gar18 as well as the kaon production Aic85 ; Fuc06 and collective flow Dan02 in heavy-ion collisions. On the other hand, the EOS of dense neutron-rich matter, especially at suprasaturation densities, remains largely uncertain due to the poorly known high-density behavior of the isospin-dependent part of nuclear matter EOS, characterized by the symmetry energy $E_{\rm sym}(n)$ (see, e.g., Ref. LCK08 ).

Nuclear data, including those from nuclear structure and heavy-ion collisions, are usually difficult to constrain the high-density $E_{\rm sym}(n)$, although its subsaturation density behavior has been relatively well determined. For example, the nuclear mass can put stringent constraints on $E_{\rm sym}(n)$ around $2/3n_{0}$ (the averaged density of nuclei) Zha13 , while the electric dipole polarizability is mainly sensitive to the $E_{\rm sym}(n)$ around $1/3n_{0}$ since isovector giant dipole resonances are essentially related to the neutrons and protons in nuclear surface Zha15 . Heavy-ion collisions perhaps is the only way in terrestrial labs to produce high density matter but the isospin asymmetry is usually small and thus the current constraints on high-density $E_{\rm sym}(n)$ are strongly model-dependent Xiao09 ; Fen10 ; Rus11 ; Coz13 ; Rus16 ; Coz18 .

In nature, neutron stars (NSs) provide an ideal site to explore dense matter. The discovery of the currently heaviest neutron star PSR J0348+0432 Ant13 with mass $2.01\pm 0.04M_{\odot}$ actually rules out soft NS matter EOSs which are not stiff enough against gravitational collapse. The NS mass-radius ($M$-$R$) relation has been shown to be sensitive to the high-density $E_{\rm sym}(n)$ Lin92 ; Lat01 ; Oze10 ; Ste10 since the averaged density of a NS is about $2.5n_{0}$ and the NS matter is dominated by neutrons with a small fraction ($\sim 10\%$) of protons (and leptons to keep weak equilibrium and charge neutrality). Although the NS mass can be determined precisely, the precise measurement of its radius remains a big challenge Oze16 . A good probe of NS radius is the tidal deformability, i.e., the ratio of the induced quadrupole moment of a neutron star to the perturbing tidal field of its companion, and for a NS with mass $M$, it can be expressed in dimensionless form as Fla08 ; Hin08

[TABLE]

where $k_{2}$ is the tidal Love number and $R$ is the NS radius. The inspiralling binary neutron star (BNS) merger, one important source of gravitational wave (GW) that can be detected by ground-based GW detectors, provides a natural lab to extract information on $\Lambda_{M}$. During the BNS inspiral stage before merger, the tidal effects change the phase evolution of GW waveform compared to that of a binary black hole (BBH) inspiral, and the difference between BNS and BBH inspirals appears from the fifth post-Newtonian order onwards with the leading order contribution proportional to $\Lambda_{M}$ Fla08 ; Dam09 ; Gra18 . The $\Lambda_{M}$ can thus be extracted from the GW signal of BNS inspiral Fla08 ; Hin10 ; Vin11 ; Dam12 .

On 17 August 2017, the first GW signal GW170817 of BNS merger was observed and localized by the LIGO and Virgo observatories Abb17NSMerger , and its electromagnetic (EM) radiation was also detected by many Collaborations (see, e.g., Ref. LAbb17EM ), inaugurating a new era of multimessenger astronomy. Using the GW170817 signal, a large number of studies Mar17 ; Bau17 ; Zho18 ; Fat18 ; De18 ; NBZ18 ; NBZ18b ; Rad18 ; Cou18 ; Ann18 ; Mos18 ; Rez18 ; Rui18 ; Mal18 have been performed to constrain the EOS of NS matter or the properties of NSs. The original analysis of GW170817 suggests an upper limit of $\Lambda_{1.4}\leq 800$ Abb17NSMerger , and a more recent analysis Abb18NSMerger with some plausible assumptions leads to a stronger constraint of $\Lambda_{1.4}=190^{+390}_{-120}$.

In this work, for the first time, by using the same model to analyze simultaneously the data based on strong, EM and gravitational measurements, i.e., the terrestrial data of finite nuclei and heavy-ion collisions, the existence of $2M_{\odot}$ NS from EM observations and the upper limit of $\Lambda_{1.4}\leq 580$ from GW170817, we put stringent constraints on the high-density $E_{\rm sym}(n)$ and the EOS of SNM, pure neutron matter (PNM) and NS matter.

II Methods

The nuclear matter EOS, defined as the binding energy per nucleon, can be expressed as the following parabolic approximation form

[TABLE]

where $n=n_{\rm n}+n_{\rm p}$ is the nucleon number density and $\delta=(n_{\rm n}-n_{\rm p})/n$ is the isospin asymmetry with $n_{\rm p}$ and $n_{\rm n}$ denoting the proton and neutron densities, respectively. $E_{0}(n)=E(n,\delta=0)$ is the EOS of SNM, and the symmetry energy is defined by $E_{\mathrm{sym}}(n)=\left.\frac{1}{2!}\frac{\partial^{2}E(n,\delta)}{\partial\delta^{2}}\right|_{\delta=0}.$ At the saturation density $n_{0}$, the $E_{0}(n)$ can be expanded in $\chi=(n-n_{0})/3n_{0}$ as $E_{0}(n)=E_{0}(n_{0})+\frac{1}{2!}K_{0}\chi^{2}+\frac{1}{3!}J_{0}\chi^{3}+{\cal O}(\chi^{4})$, where $K_{0}$ is the incompressibility coefficient and $J_{0}$ is the skewness coefficient. The $E_{\mathrm{sym}}(n)$ can be expanded at a reference density $n_{r}$ in terms of the slope parameter $L(n_{r})$ and the curvature parameter $K_{\mathrm{sym}}(n_{r})$ as $E_{\mathrm{sym}}(n)=E_{\mathrm{sym}}(n_{r})+L(n_{r})\chi_{r}+\frac{1}{2!}K_{\mathrm{sym}}(n_{r})\chi_{r}^{2}+\mathcal{O}(\chi_{r}^{3})$, with $\chi_{r}=(n-n_{r})/(3n_{r})$. Conventionally we have $L\equiv L(n_{0})$ and $K_{\mathrm{sym}}\equiv K_{\mathrm{sym}}(n_{0})$.

In this work, we apply the extended Skyrme-Hartree-Fock (eSHF) model Cha09 ; Zha16 to three systems, i.e., nuclear matter, finite nuclei and neutron stars. As emphasized in Ref. Zha16 , the eSHF model includes additional momentum and density-dependent two-body forces to effectively mimic the momentum dependence of the three-body force and can very successfully describe simultaneously the three systems which involve a wide density region, and thus is especially suitable for our present motivation. The extended Skyrme interaction is expressed as Cha09 ; Zha16

[TABLE]

where the symbols have their conventional meaning Cha09 ; Zha16 . The interaction contains $14$ independent parameters, i.e., the $13$ Skyrme parameters $\alpha$, $t_{0}\sim t_{5}$, $x_{0}\sim x_{5}$, and the spin-orbit coupling constant $W_{0}$. Instead of directly using the $13$ Skyrme parameters, we express them explicitly in terms of the following $13$ macroscopic quantities (pseudo-parameters) Zha16 , i.e., $n_{0}$, $E_{0}(n_{0})$, $K_{0}$, $J_{0}$, $E_{\rm sym}(n_{r})$, $L(n_{r})$, $K_{\rm sym}(n_{r})$, the isoscalar effective mass $m_{s,0}^{\ast}$, the isovector effective mass $m_{v,0}^{\ast}$, the gradient coefficient $G_{S}$, the symmetry-gradient coefficient $G_{V}$, the cross gradient coefficient $G_{SV}$, and the Landau parameter $G_{0}^{\prime}$ of SNM in the spin-isospin channel. The higher-order parameters $J_{0}$ and $K_{\rm sym}$ generally have small influence on the properties of finite nuclei but are critical for the high-density neutron-rich matter EOS and NS properties. In addition, at the subsaturation density $n_{c}=0.11n_{0}/0.16$, the $E_{\rm sym}(n_{c})$ has been precisely constrained to be $E_{\rm sym}(n_{c})=26.65\pm 0.2$ MeV Zha13 by analyzing the binding energy difference of heavy isotope pairs and $L(n_{c})=47.3\pm 7.8$ MeV Zha14 is extracted from the electric dipole polarizability of 208Pb. Therefore, here we fix $J_{0}$ and $K_{\rm sym}$ at various values with $E_{\rm sym}(n_{c})=26.65$ MeV and $L(n_{c})=47.3$ MeV, and the other $10$ parameters are obtained by fitting the data of finite nuclei by minimizing the weighted sum of the squared deviations between the theoretical predictions and the experimental data, i.e., $\chi^{2}(\bm{p})=\sum_{i=1}^{N}\left(\frac{{\cal O}_{i}^{\rm th}({\bm{p}})-{\cal O}_{i}^{\rm exp}}{\triangle{\cal O}_{i}}\right)^{2},$ where the $\bm{p}=(p_{1},...,p_{z})$ define the $z$ dimensional model space, ${\cal O}_{i}^{({\rm th})}$ and ${\cal O}_{i}^{({\rm exp})}$ are the theoretical predictions and the corresponding experimental values of observables, respectively, and $\Delta\mathcal{O}_{i}$ is the adopted error used to balance the relative weights of the various types of observables. We note here that varying $L(n_{c})$ within $L(n_{c})=47.3\pm 7.8$ MeV mainly influences the values of $G_{S}$, $G_{V}$ and $G_{SV}$ which are irrelevant to the nuclear matter EOS, while the parameters $n_{0}$, $E_{0}(n_{0})$, $K_{0}$ characterizing the nuclear matter EOS remain almost unchanged.

In the fitting, we consider the following experimental data of spherical even-even nuclei: (i) the binding energies $E_{\rm{B}}$ of 16O, 40,48Ca, 56,68Ni, 88Sr, 90Zr, 100,116,132Sn, 144Sm and 208Pb Wan17 ; (ii) the charge r.m.s. radii $r_{\rm{c}}$ of 16O, 40,48Ca, 56Ni, 88Sr, 90Zr, 116Sn, 144Sm and 208Pb Ang13 ; Fri95 ; Bla05 ; (iii) the isoscalar giant monopole resonance energies $E_{\rm{GMR}}$ of 90Zr, 116Sn, 144Sm and 208Pb You99 ; (iv) the spin-orbit energy level splittings $\epsilon_{\rm{ls}}^{A}$ for neutron $1p_{1/2}-1p_{3/2}$ and proton $1p_{1/2}-1p_{3/2}$ in 16O, and the proton $2d_{3/2}-2d_{5/2}$, neutron $3p_{1/2}-3p_{3/2}$ and neutron $2f_{5/2}-2f_{7/2}$ in 208Pb Vau72 . To balance the $\chi^{2}$ from each sort of experimental data (see, e.g., Ref. Zha16 ), we assign the errors of $1.0$ MeV and $0.01$ fm to the $E_{\rm{B}}$ and $r_{\rm{c}}$, respectively. In particular, for the $E_{\rm{GMR}}$ we use the experimental error multiplied by $3.5$ to also consider the impact of the experimental error, while for the $\epsilon_{\rm{ls}}^{A}$ a $10\%$ relative error is employed.

For NSs, we consider here the conventional NS model which includes only nucleons, electrons and possible muons ($npe\mu$), and the NS is assumed to contain core, inner crust and outer crust. For the core, the EOS of $\beta$-stable and electrically neutral $npe\mu$ matter is obtained from the eSHF model. For the inner crust in the density region between $n_{\rm out}$ and $n_{\rm t}$, the EOS is constructed by interpolating with $P=a+b{\cal E}^{4/3}$ Car03 where $P$ is pressure and $\cal E$ is energy density. The density $n_{\rm{out}}$ separating the inner and the outer crusts is taken to be $2.46\times 10^{-4}\rm{fm}^{-3}$ while the core-crust transition density $n_{\rm t}$ is evaluated self-consistently by a dynamical approach XuJ09a . For the outer crust, we use the well-known BPS EOS in the density region of $6.93\times 10^{-13}{\rm fm}^{-3}<n<n_{\rm out}$ and the FMT EOS for $4.73\times 10^{-15}{\rm fm}^{-3}<n<6.93\times 10^{-13}{\rm fm}^{-3}$ Bay71 ; Iid97 . It should be noted that all the extended Skyrme interactions used in following NS calculations satisfy the causality condition $dP/d{\cal E}<1$.

III Results and Discussions

Shown in Fig. 1 is pressure vs density for SNM within eSHF in various extended Skyrme parameter sets with $J_{0}$ fixed at ($-300$, $-350$, $-400$, $-450$, $-500$) MeV and $K_{\rm sym}$ in the range of ($-200,60$) MeV. Also included in the figure is the constraint from collective flow data in heavy-ion collisions Dan02 . For all the parameter sets with fixed $J_{0}$ and $K_{\rm sym}$, as expected, the total chi-square $\chi_{\rm tot}^{2}$ falls in the range of $24.45<\chi_{\rm tot}^{2}<36.24$, and the mean $\chi^{2}$ of each sort of experimental data (i.e., $\chi^{2}_{E_{\rm B}}/12$, $\chi^{2}_{r_{\rm c}}/9$, $\chi^{2}_{E_{\rm GMR}}/4$ and $\chi^{2}_{\epsilon_{\rm ls}^{A}}/5$) is approximately equal to $1$. The small variation of $\chi_{\rm tot}^{2}$ suggests that the higher-order parameters $J_{0}$ and $K_{\rm sym}$ indeed have small influence on the properties of finite nuclei. In addition, the pressure of SNM exhibits negligible dependence on the $K_{\rm sym}$, especially for $J_{0}>-500$ MeV. One sees that the pressure of SNM becomes stiffer as the $J_{0}$ increases, and $J_{0}=-300$ MeV predicts a too stiff SNM EOS that violates the flow data. A detailed study indicates an upper limit at $J^{\rm up}_{0}=-342$ MeV.

Using the same parameter sets as used in Fig. 1, we show in Fig. 2 the NS maximum mass $M_{\rm max}$ vs $K_{\rm sym}$. One sees the $M_{\rm max}$ increases sensitively with increasing $J_{0}$ for a fixed $K_{\rm sym}$. The $K_{\rm sym}$ has small influence on the $M_{\rm max}$ when $K_{\rm sym}$ is greater than about $-100$ MeV, but for a fixed $J_{0}$, the $M_{\rm max}$ is drastically reduced with decreasing $K_{\rm sym}$ for $K_{\rm sym}\lesssim-100$ MeV. This behavior is due to the fact that a larger $K_{\rm sym}$ leads to a stiffer high-density symmetry energy, and thus smaller isospin asymmetry in high-density NS matter. The smaller isospin asymmetry in turn suppresses the sensitivity of $M_{\rm max}$ to $K_{\rm sym}$, as the $M_{\rm max}$ is mainly determined by the EOS of high-density NS matter. In particular, for sufficiently large $K_{\rm sym}$ (e.g., $\gtrsim-100$ MeV), the NS matter can become almost isospin-symmetric at high densities, and the NS maximum mass is therefore mainly sensitive to $J_{0}$ which dominates high-density behavior of symmetric nuclear matter. Since the flow data in heavy-ion collisions require $J_{0}\leq-342$ MeV, the parameter sets with $J_{0}=-342$ MeV generally predict largest values of $M_{\rm max}$, and the corresponding results are also included Fig. 2. For $J_{0}=-342$ MeV, it is seen that when $K_{\rm sym}$ is smaller than $-175$ MeV, the predicted $M_{\rm max}$ becomes violating the mass lower limit (i.e., $1.97M_{\odot}$) of the heaviest NS PSR J0348+0432 Ant13 observed so far (its mass $2.01\pm 0.04M_{\odot}$ is shown as a shaded band in Fig. 2), leading to a lower limit at $K^{\rm low}_{\rm sym}=-175$ MeV. For a fixed $J_{0}$, we find that the $\Lambda_{1.4}$ rapidly increases with increasing $K_{\rm sym}$. For a fixed $K_{\rm sym}$, the $\Lambda_{1.4}$ also increases with $J_{0}$ but much weaker than that with $K_{\rm sym}$. Our results indicate that the data of finite nuclei, the flow data in heavy-ion collisions and the existence of $2M_{\odot}$ NS together give the limit of $K_{\rm sym}\geq-175$ MeV, leading a lower limit of $\Lambda^{\rm low}_{1.4}=261$.

The above analyses mean the $E_{\rm sym}(n)$ cannot be too soft. When the $K_{\rm sym}$ increases, the $E_{\rm sym}(n)$ becomes stiffer and the $\Lambda_{1.4}$ increases accordingly. The most recent limit of $\Lambda_{1.4}\leq 580$ Abb18NSMerger thus can put an upper limit for $K_{\rm sym}$ for each $J_{0}$ as indicated Fig. 2. The limit of $\Lambda_{1.4}\leq 580$ together with the data of finite nuclei, the flow data in heavy-ion collisions and the existence of $2M_{\odot}$ NS thus give an allowed region for the higher-order parameters $J_{0}$ and $K_{\rm sym}$, as shown by green region in Fig. 2, which leads to $-464~{}{\rm MeV}\leq J_{0}\leq-342~{}{\rm MeV}$ and $-175~{}{\rm MeV}\leq K_{\rm sym}\leq-36~{}{\rm MeV}$. Moreover, the largest NS mass is determined to be $2.28M_{\odot}$ at $(J_{0},K_{\rm sym})=(-342,-62)$ MeV as indicated in Fig. 2.

Shown in Fig. 3(a) is $\Lambda_{1.4}$ vs $K_{\rm sym}$ within eSHF using the same parameter sets as used in Fig. 2. The corresponding results for $\Lambda_{1.4}$ vs $R_{1.4}$ are shown in Fig. 3(b) and the inset in Fig. 3(b) displays the results for $k_{2,1.4}$ vs $R_{1.4}$. The allowed region for $J_{0}$ and $K_{\rm sym}$ is also included in Fig. 3(a). As already mentioned, one indeed sees that $\Lambda_{1.4}$ is sensitive to $K_{\rm sym}$ but is less affected by the $J_{0}$. From Fig. 3(b) and the inset, one sees both $\Lambda_{1.4}$ and $k_{2,1.4}$ exhibit very strong correlation with $R_{1.4}$, and they can be nicely fitted by the formulae $\Lambda_{1.4}=a_{\Lambda}R_{1.4}^{\alpha_{\Lambda}}$ and $k_{2,1.4}=a_{k_{2}}R_{1.4}^{\alpha_{k_{2}}}$, respectively, with $a_{\Lambda}=(1.41\pm 0.14)\times 10^{-6}$, $\alpha_{\Lambda}=7.71\pm 0.04$, $a_{k_{2}}=(8.25\pm 0.58)\times 10^{-5}$ and $\alpha_{k_{2}}=2.69\pm 0.03$. The correlation coefficient is $r_{\Lambda}=0.999$ for $\Lambda_{1.4}=a_{\Lambda}R_{1.4}^{\alpha_{\Lambda}}$ and $r_{k_{2}}=0.996$ for $k_{2,1.4}=a_{k_{2}}R_{1.4}^{\alpha_{k_{2}}}$. These relations together with $261\leq\Lambda_{1.4}\leq 580$ lead to the stringent constraints of $R_{1.4}\in[11.8,13.1]$ km and $k_{2,1.4}\in[0.064,0.085]$.

According to the allowed parameter space for $J_{0}$ and $K_{\rm sym}$ as shown in Fig. 2, we can determine the EOS of dense matter. The obtained results for $E_{\rm sym}(n)$ are shown in Fig. 4(a), and the pressure vs density for SNM, PNM and NS matter is exhibited in Figs. 4(b), (c) and (d), respectively. Also included in Fig. 4(a) are the constraints at subsaturation densities from midperipheral heavy-ion collisions of Sn isotopes Tsa09 , the isobaric analog states (IAS) and combining the neutron skin data (IAS + NSkin) Dan14 , and the electric dipole polarizability ($\alpha_{D}$) in 208Pb Zha15 . In addition, the constraints on pressure for SNM from flow data in heavy-ion collisions Dan02 and that for NS matter from GW170817 Abb18NSMerger are also included in Figs. 4(b) and (d), respectively. Furthermore, we include the corresponding results with $L(n_{c})=55.1$ MeV and $39.5$ MeV to display the uncertainty due to the $L(n_{c})$ (i.e., $L(n_{c})=47.3\pm 7.8$ MeV Zha14 ). One sees that our results are consistent with the existing constraints but with much higher precision due to the simultaneous consideration of the data of finite nuclei, the flow data in heavy-ion collisions, the observed heaviest NS and the GW170817 signal. We would like to point out that the high density $E_{\rm sym}(n)$ still has large uncertainty and it could be negative at high densities, which will cause isospin instability and thus the presence of PNM core in NSs.

Furthermore, our results indicate that the $L(n_{c})=39.5$ MeV gives smaller lower limits of $K^{\rm low}_{\rm sym}=-203$ MeV, $\Lambda^{\rm low}_{1.4}=193$, $J^{\rm low}_{0}=-475$ MeV and a heavier $M_{\rm max}=2.30M_{\odot}$, while $L(n_{c})=55.1$ MeV gives lager lower limits of $K^{\rm low}_{\rm sym}=-138$ MeV, $\Lambda^{\rm low}_{1.4}=380$, $J^{\rm low}_{0}=-455$ MeV and a lighter $M_{\rm max}=2.26M_{\odot}$. In addition, the $E_{\rm sym}(2n_{0})$ is found to be [46.9, 57.6] MeV, [39.4, 54.5] MeV and [33.0, 51.3] MeV for $L(n_{c})=39.5$ MeV, $47.3$ MeV and $51.1$ MeV, respectively. Using $L(n_{c})=47.3\pm 7.8$ MeV, therefore, we obtain $J_{0}\in[-464_{-11}^{+9},-342]$ MeV, $K_{\rm sym}\in[-175_{-28}^{+37},-36\mp 2]$ MeV, $E_{\rm sym}(2n_{0})\in[39.4_{+7.5}^{-6.4},54.5_{+3.1}^{-3.2}]$ MeV, $\Lambda_{1.4}\in[261_{-68}^{+119},580]$, $R_{1.4}\in[11.8_{-0.7}^{+0.8},13.1\pm 0.2]$ km, and $M_{\rm max}=2.28\mp 0.02~{}M_{\odot}$.

For the higher-order parameters $J_{0}$ and $K_{\rm sym}$, the $J_{0}\in[-464_{-11}^{+9},-342]$ MeV gives the strongest constraint compared to the existing ones Cai14 , and the $K_{\rm sym}\in[-175_{-28}^{+37},-36\mp 2]$ MeV is also consistent with those extracted from the symmetry energy systematics with some correlations ChenLW11 ; ChenLW15 ; Mon17 ; Mal18 or from heavy-ion collisions Coz18 . Moreover, the $E_{\rm sym}(2n_{0})\in[39.4_{+7.5}^{-6.4},54.5_{+3.1}^{-3.2}]$ MeV is in good agreement with those extracted from the symmetry energy systematics ChenLW15 , heavy-ion collisions Rus16 and the recent analyses on the NS observation and the GW170817 signal NBZ18 ; NBZ18b . As for the NS properties, the $\Lambda^{\rm low}_{1.4}=261_{-68}^{+119}$ agrees with the constraints from analyzing the GW170817 signal Ann18 ; Mos18 or its EM signals Rad18 ; Cou18 , the $R_{1.4}\in[11.8_{-0.7}^{+0.8},13.1\pm 0.2]$ km agrees with those from analyzing the GW170817 Ann18 ; Mos18 ; Mal18 , and the $M_{\rm max}=2.28\mp 0.02~{}M_{\odot}$ is consistent with the results from analyzing GW170817 Rez18 ; Rui18 .

Although the polytropic form of $P=a+b{\cal E}^{4/3}$ has been extensively used to approximate the EOS of the NS inner crust Lat01 ; Car03 ; XuJ09a ; Zha16 , where the nuclear pasta may exist, the sensitivity of our results to the choice of polytropic index needs to be studied. As an example, based on the EOS with $J_{0}=-342$ MeV and $K_{\rm sym}=-175$ MeV, which gives $M_{\rm max}=1.97M_{\odot}$, we obtain $\Lambda_{1.4}=261.76~{}(260.86,260.50)$, $R_{1.4}=12.52~{}(11.82,11.55)$ km and $k_{2,1.4}=0.048~{}(0.064,0.071)$ with the polytropic index $1~{}(4/3,5/3)$. These results indicate that the increase of the polytropic index leads to considerable increase of the $k_{2,1.4}$, small reduction ($<10\%$) of $R_{1.4}$, and negligible decrease of $\Lambda_{1.4}$, suggesting the choice of the polytropic index has negligible effects on $\Lambda_{1.4}$.

IV Summary

Using the eSHF model to simultaneously analyze the data from terrestrial nuclear experiments and astrophysical observations based on strong, EM and gravitational measurements, we have put stringent constraints on the high-density $E_{\rm sym}(n)$ and the pressure of SNM, PNM and NS matter. We have found that the nuclear data and the existence of $2M_{\odot}$ NS rule out too soft high-density $E_{\rm sym}(n)$, leading to lower limits of $\Lambda_{1.4}\geq 193$ and $R_{1.4}\geq 11.1$ km. Further combining the upper limit of $\Lambda_{1.4}\leq 580$ from GW170817 excludes too stiff high-density $E_{\rm sym}(n)$, leading to an upper limit of $R_{1.4}\leq 13.3$ km. Using $L(n_{c})=47.3\pm 7.8$ MeV, we have obtained $J_{0}\in[-464_{-11}^{+9},-342]$ MeV, $K_{\rm sym}\in[-175_{-28}^{+37},-36\mp 2]$ MeV, $E_{\rm sym}(2n_{0})\in[39.4_{+7.5}^{-6.4},54.5_{+3.1}^{-3.2}]$ MeV, and $M_{\rm max}=2.28\mp 0.02~{}M_{\odot}$. In future, more precise limit on $L(n_{c})$, possible discovery of heavier NS and tighter bound on $\Lambda_{1.4}$ from BNS merger would put stronger constraints on the high-density $E_{\rm sym}(n)$ and thus the EOS of dense neutron-rich matter.

Acknowledgments

The authors thank Tanja Hinderer, Ang Li, Bao-An Li and Jorge Piekarewicz for helpful discussions. This work was supported in part by the National Natural Science Foundation of China under Grant No. 11625521, the Major State Basic Research Development Program (973 Program) in China under Contract No. 2015CB856904, the Program for Professor of Special Appointment (Eastern Scholar) at Shanghai Institutions of Higher Learning, Key Laboratory for Particle Physics, Astrophysics and Cosmology, Ministry of Education, China, and the Science and Technology Commission of Shanghai Municipality (11DZ2260700).