How tightly is nuclear symmetry energy constrained by unitary Fermi gas?

TL;DR

This paper investigates how the universal equation of state of the unitary Fermi gas constrains the density dependence of nuclear symmetry energy, highlighting the limitations imposed by uncertainties in higher order parameters.

Contribution

It critically assesses the extent to which the unitary Fermi gas EOS constrains the nuclear symmetry energy and emphasizes the importance of precise higher order parameter measurements.

Findings

$E_{UG}( ho)$ provides a lower boundary for $E_{sym}$.

Uncertainties in $K_{sym}$ and skewness parameters significantly affect the $E_{sym}$ and $L$ correlation.

Constraints are limited without precise knowledge of curvature and skewness parameters.

Abstract

We examine critically how tightly the density dependence of nuclear symmetry energy \esym is constrained by the universal equation of state (EOS) of the unitary Fermi gas $E_{\rm{UG}}(\rho)$ considering currently known uncertainties of higher order parameters describing the density dependence of the Equation of State of isospin-asymmetric nuclear matter. We found that $E_{\rm{UG}}(\rho)$ does provide a useful lower boundary for the \esym. However, it does not tightly constrain the correlation between the magnitude $E_{\rm{sym}}(\rho_0)$ and slope $L$ unless the curvature $K_{\rm{sym}}$ of the symmetry energy at saturation density $\rho_0$ is more precisely known. The large uncertainty in the skewness parameters affects the $E_{\rm{sym}}(\rho_0)$ versus $L$ correlation by the same almost as significantly as the uncertainty in $K_{\rm{sym}}$.

Equations

$$E_{\text{sym}}(\mu )\approx E_{\text{PNM}}(\mu )-E_{\text{0}}(\mu ),$$

$$E_{\text{UG}}(\mu )=\frac{3{\hslash }^2{k}_{\text{F}}^2}{10m_{\text{n}}}\xi \equiv E_{\text{UG}}^0{\mu }^{2/3}$$

$$E_{\text{sym}}(\mu )\ge E_{\text{UG}}(\mu )-E_{\text{0}}(\mu )=E_{\text{UG}}^0{\mu }^{2/3}-E_{\text{0}}(\mu ).$$

$$E_{\text{0}}(\mu )=E_0({\rho }_0)+\frac{K_0}{18}(\mu -1{)}^2+\frac{J_0}{162}(\mu -1{)}^3+\mathcal{O}[(\mu -1{)}^4]$$

$$E_{\text{sym}}(\mu )=$$

$$+\frac{J_{\text{sym}}}{162}[(\mu -1{)}^3]+\mathcal{O}[(\mu -1{)}^4]$$

$$E_{\text{sym}}({\rho }_0)\ge$$

$$-\frac{J_n}{162}(\mu -1{)}^3$$

$$L=\frac{2E_{\text{UG}}^0}{{\mu }^{1/3}}-\frac{K_n}{3}(\mu -1)-\frac{J_n}{18}(\mu -1{)}^2.$$

$$E_{\text{sym}}({\rho }_0)\ge$$

$$+\frac{J_n}{81}(\mu -1{)}^3-E_0({\rho }_0).$$

$${\rho }_0=0.157{\text{fm}}^{-3},E_0({\rho }_0)=-15.5\text{MeV},\xi =0.37,$$

$$K_n=370\text{MeV},J_n=-500\text{MeV}$$

Full text

Introduction: To understand the nature of neutron-rich nucleonic matter has been a major scientific goal in both nuclear physics and astrophysics. The density dependence of nuclear symmetry energy $E_{\textrm{sym}}(\rho)$ has been a major uncertain part of the equation of state (EOS) of neutron-rich matter especially at high densities, see, e.g., collections in Tesym . Reliable knowledge about the $E_{\textrm{sym}}(\rho)$ has significant ramifications in answering many interesting questions regarding the structure of rare isotopes and neutron stars, dynamics of heavy-ion collisions and supernova explosions as well as the frequency and strain amplitude of gravitational waves from deformed pulsars and/or cosmic collisions involving neutron stars. During the last two decades, significant efforts have been devoted to exploring the $E_{\textrm{sym}}(\rho)$ using both terrestrial laboratory experiments ireview98 ; ibook01 ; ditoro ; LCK08 ; Lynch09 ; Dan09 ; Trau12 ; Tsang12 ; Chuck14 ; Baldo16 ; Li2017 and astrophysical observations Steiner05 ; Lat13 ; Newton14 ; Lida14 ; Pearson ; Far14 ; Fis14 ; Blas16 . Extensive surveys of the extracted constraints on the $E_{\textrm{sym}}(\rho)$ around the saturation density $\rho_{0}$ indicate that the central values of the $E_{\textrm{sym}}(\rho_{0})$ and its slope $L=\left[3\rho(\partial E_{\textrm{sym}}/\partial\rho\right)]_{\rho_{0}}$ scatter around $31.6$ MeV and $58.9$ MeV, respectively Lat13 ; Li-han ; Oer16 . At densities away from $\rho_{0}$, however, the $E_{\textrm{sym}}(\rho)$ remains rather unconstrained especially at supra-saturation densities Li2017 .

Interestingly, recent progresses in another seemingly different field may provide additional information about the density dependence of nuclear symmetry energy. Indeed, theoretical and experimental studies of cold atoms have made impressive progress in recent years, see, e.g., refs. Zwi15 ; Ku12 ; Zu13 ; Endres13 for recent reviews, providing reliable information about the universal EOS ($E_{\textrm{UG}}$) of unitary gas (UG) interacting via pairwise $s$-waves with infinite scattering length but zero effective range. The universal $E_{\textrm{UG}}$ constrains stringently the EOS of pure neutron matter (PNM) at sub-saturation densities, thus provides possibly additional constraints on the nuclear symmetry energy. In fact, it was recently conjectured that the $E_{\textrm{UG}}$ provides the lower boundary of the EOS of PNM ($E_{\textrm{PNM}}$) Kolomeitsev16 . Moreover, using a set of known parameters of symmetric nuclear matter (SNM), and taking zero as an upper bound on the curvature $K_{\textrm{sym}}=\left[9\rho^{2}(\partial^{2}E_{\textrm{sym}}(\rho)/\partial\rho^{2})\right]_{\rho_{0}}$ of $E_{\textrm{sym}}(\rho)$ at $\rho_{0}$, the authors of ref. Kolomeitsev16 obtained a region of $E_{\textrm{sym}}(\rho_{0})$ - $L$ space that is inconsistent with the unitary gas constraints, excluding many $E_{\textrm{sym}}(\rho)$ functionals currently actively used in both nuclear physics and astrophysics.

Our original purposes were to examine several issues not clearly addressed in version-1 of ref. Kolomeitsev16 . We notice that some of these issues are now discussed in more detail in its revised version. Nevertheless, it is still useful to provide our results and opinions on some of these issues. The derivation of the excluded region in $E_{\textrm{sym}}(\rho_{0})$ - $L$ space by ref. Kolomeitsev16 relies on several assumptions NST : the most importance of which is the underlying conjecture that $E_{\textrm{PNM}}(\rho)\geq E_{\textrm{UG}}(\rho)$ where $E_{\textrm{UG}}(\rho)=\xi E_{F}(\rho)$ with $\xi$ being the Bertsch parameter Zwi15 ; Ku12 ; Zu13 ; Endres13 and $E_{F}\propto(\rho/\rho_{0})^{2/3}$ is the energy of a non-interacting non-relativistic degenerate Fermi gas of neutrons. This conjecture is not quite the same as stating that PNM has a greater energy than the UG at all densities, since the UG is experimentally accessible only at low densities. This conjecture merely employs an algebraic expression motivated by the result that apparently the low-density neutron gas has a higher energy than the UG, and by the inference that the repulsive nature of the three nucleon (NNN) interaction coupled with the finite values of the range and scattering length of the two-body (NN) s-wave interaction will ensure that the energy of PNM remains above this algebraic expression at higher densities. A further assumption is that possibly attractive higher order NN interactions (p-wave, d-wave, etc.) are not important relative to the repulsive character of three- and higher-body interactions. Additionally, neutrons are assumed to remain non-relativisitc in the density range considered. Given the important ramifications of the findings in ref. Kolomeitsev16 , we are motivated to critically examine them adopting the same assumptions. We include the third-order terms in density characterized by the skewness coefficients $J_{0}=27\rho_{0}^{3}\partial^{3}E_{0}(\rho)/\partial\rho^{3}|_{\rho=\rho_{0}}$ and $J_{\textrm{sym}}=27\rho_{0}^{3}\partial^{3}E_{\textrm{sym}}(\rho)/\partial\rho^{3}|_{\rho=\rho_{0}}$ in expanding the $E_{\textrm{0}}(\rho)$ and $E_{\textrm{sym}}(\rho)$, respectively. Particularly, we carefully examine the uncertainties of the curvature $K_{\textrm{sym}}$ of the symmetry energy and the skewness coefficients $J_{0}$ and $J_{\textrm{sym}}$, taking into account energy density functionals that are consistent with the PNM EOS derived from microscopic calculations, and examine the effects of those uncertainties on the region of $E_{\textrm{sym}}(\rho_{0})$-$L$ space excluded by the unitary gas constraints.

While Skyrme models consistent with microscopic PNM calculations tend to give $K_{\textrm{sym}}$ in the range -100 to -200 MeV Fattoyev12 ; Brown13 , Relativistic Mean Field (RMF) models consistent with microscopic PNM calculations can give positive values of $K_{\textrm{sym}}$ Fattoyev12 ; Dutra14 , reflecting a difference in the form of these two classes of energy density functionals. Indeed, some reputable non-relativistic and relativistic energy density functionals in the literature, see, e.g., reviews in ref. Lwchen09 ; Dutra12 ; Dutra14 ; Colo14 , predict positive $K_{\textrm{sym}}$ values and meet all existing constraints including the EOS of PNM within their known uncertain ranges. For example, the TM2 RMF interaction has $K_{\textrm{sym}}$ = 50 MeV and passes the PNM test of ref.Dutra14 . To our best knowledge, while the majority of existing models predict negative values for $K_{\textrm{sym}}$, there is no fundamental physics principle excluding positive $K_{\textrm{sym}}$ values. The current situation clearly calls for more studies on the $K_{\textrm{sym}}$ especially its experimental constraints. Hopefully, ongoing experiments at several laboratories Umesh to extract the isospin dependence of nuclear incompressibility and subsequently the $K_{\textrm{sym}}$ from giant resonances of neutron-rich nuclei will help settle the issue in the near future. Interestingly, very recently by analyzing comprehensively the relative elliptical flows of neutrons and protons measured by the ASY-EOS and the FOPI-LAND collaborations at GSI using a Quantum Molecular Dynamics (QMD) model Cozma , the extracted values for the slope and curvature parameters are $L=59\pm 24$ MeV and $K_{\textrm{s}ym}=88\pm 372$ MeV at the $1\sigma$ confidence level. After considering uncertainties of other model parameters including the incompressibility of SNM, neutron-proton effective mass splitting (related to the momentum dependence of the symmetry potential), Pauli blocking and in-medium nucleon-nucleon cross sections, it was concluded that $L=59\pm 24(\textrm{exp})\pm 16(\textrm{th})\pm 10(\textrm{sys})$ MeV and $K_{\textrm{sym}}=0\pm 370(\textrm{exp})\pm 220(\textrm{th})\pm 150(\textrm{sys})$ MeV. To our best knowledge, the latter represents the latest and most accurate constraint on $K_{\textrm{sym}}$.

The lower boundary of nuclear symmetry energy constrained by the universal EOS of unitary Fermi gas: Within the parabolic approximation for the EOS of isospin asymmetric nuclear matter (ANM) in terms of the energy per nucleon $E$, i.e., $E(\rho,\delta)=E_{0}(\rho)+E_{\textrm{sym}}(\rho)\delta^{2}+\mathcal{O}(\delta^{4}),$ the symmetry energy $E_{\textrm{sym}}(\rho)=2^{-1}[\partial^{2}E(\rho,\delta)/\partial\delta^{2}]_{\delta=0}$ can be approximated by

[TABLE]

where $\mu=\rho/\rho_{0}$ is the reduced density and $\delta=(\rho_{\textrm{n}}-\rho_{\textrm{p}})/\rho$ is the isospin asymmetry of ANM. Using the conjecture $E_{\textrm{PNM}}(\rho)\geq E_{\textrm{UG}}(\rho)$ Kolomeitsev16 and the EOS of unitary gas

[TABLE]

where $k_{\textrm{F}}$ is the neutron Fermi momentum, the lower boundary of symmetry energy can be obtained from

[TABLE]

It is necessary to caution that the above lower boundary of $E_{\textrm{sym}}(\mu)$ is estimated based on the conjecture $E_{\textrm{PNM}}(\rho)\geq E_{\textrm{UG}}(\rho)$ and the assumption that $\xi$ is a constant in the density range we study. As emphasized in ref. Kolomeitsev16 , the conjecture is empirical in nature. While there are strong supports for the conjecture by comparing the EOSs of PNM calculated from various microscopic many-body theories with the $E_{\textrm{UG}}(\mu)$ using a constant $\xi\approx 0.37$ (see Fig.1 of ref. Kolomeitsev16 and Fig.2 of ref. CaiLi15 ) up to about $\rho_{0}$, the rigorous condition for unitarity is expected to be reached in PNM only at very low densities. Although one can not prove the validity of the conjecture at high densities, strong physical arguments were made to justify and use it up to about $1.5\rho_{0}$ in ref. Kolomeitsev16 . Thus, the results of our study should be understood with the caveat that they are obtained under the above reasonable but not rigorously proven conjecture and assumptions. Nevertheless, they are useful for comparing with the results of ref. Kolomeitsev16 obtained using the same assumptions.

The EOS of SNM around $\rho_{0}$ can be expanded to the third order in density as

[TABLE]

in terms of the incompressibility $K_{0}$ and skewness $J_{0}$. At the saturation point of SNM, we adopt $E_{0}(\rho_{0})=-15.9$ MeV and $\rho_{0}=0.164$ fm*-3* Brown13 . The lower boundary of $E_{\textrm{sym}}(\rho)$ thus depends on the values of $K_{0}$, $J_{0}$ and $\xi$.

The incompressibility $K_{0}$ of SNM has been extensively investigated Colo14 ; Stone14 , and the most widely used values are $K_{0}=240\pm 20$ MeV Shlomo06 ; Piekarewicz10 or $230\pm 40$ MeV Khan12 . However, the skewness coefficient $J_{0}$ is still poorly known Meixner13 ; Chen11 ; Farine97 ; Klahn06 ; Mas16 ; Cai14 . In Fig. 1(a), we show $K_{0}$ and $J_{0}$ from 274 parameterizations of the Skyrme and RMF models that pass the PNM tests of Dutra et al Dutra12 ; Dutra14 and satisfy $K_{0}=230\pm 40$ MeV. The spread in values for $J_{0}$ is very large, covering the range $\approx-800<J_{0}<$400 MeV.

As reviewed recently in refs. Endres13 ; Kolomeitsev16 , currently the best estimate for the Bertsch parameter $\xi$ from lattice Monte Carlo studies is $\xi=0.372(5)$ consistent with the most accurate experimental value of $\xi=0.376(4)$. While its values from various other models and experiments have scattered between 0.279 and 0.449(9) within the last decade, it appears that it now has converged to $\xi=0.37\pm 0.005$ which we adopt in this work. Shown in Fig. 2 with the red dashed lines are the variation of $E_{\textrm{sym}}(\rho)$ with $\xi=0.37\pm 0.005$, $J_{0}=0$  and $K=230$ MeV. Effects of varying the $\xi$ value are very small within the range considered.

Secondly, effects of the skewness parameter are shown by varying the value of $J_{0}$ between $-800$ MeV and $400$ MeV, the range covered by the models plotted in Fig. 1(a). Although the range is very large, it translates to a range of uncertainty for $E_{\textrm{sym}}(\rho)$ that is equivalent to the range of uncertainty in $K_{0}$ (80 MeV). This is easy to understand as the expansion of SNM’s EOS converges quickly around the normal density by design (the $J_{0}$ contribution of $J_{0}/162$ is a factor of 9 less than the $K_{0}/18$ term).

Considering the uncertainties of all relevant parameters involved, the most conservative lower boundary of $E_{\textrm{sym}}(\rho)$ shown as the shadowed region in Fig. 2 is obtained by using $\xi=0.37$, $\rho_{0}=0.157$ fm*-3*, $E_{0}(\rho_{0})=-15.5$ MeV, and $K_{0}$=270 MeV; for $\mu\leq 1$, $J_{0}=-800$ MeV and for $\mu>1$, $J_{0}=400$ MeV. Overall, our observations and results are consistent with the findings in ref. Kolomeitsev16 .

Constraining the $E_{\textrm{sym}}(\rho_{0})$ versus $L$ boundary: The symmetry energy $E_{\textrm{sym}}(\mu)$ can be expanded around $\rho_{0}$ to third order in density as

[TABLE]

in terms of its magnitude $E_{\textrm{sym}}(\rho_{0})$, slope $L$, curvature $K_{\textrm{sym}}$ and skewness $J_{\textrm{sym}}$ at $\rho_{0}$. Inserting the above equation into Eq. (3), the lower boundary of $E_{\textrm{sym}}(\rho_{0})$ can be expressed as

[TABLE]

where $K_{n}=K_{\textrm{sym}}+K_{0}$ and $J_{n}=J_{\textrm{sym}}+J_{0}$. Taking the derivative of the above equation with respect to density on both sides, one can readily get an expression for the lower boundary of $L$

[TABLE]

Then, putting the above expression back to Eq. (How tightly is the nuclear symmetry energy constrained by a unitary Fermi gas?) the latter can be rewritten as

[TABLE]

These two equations reveal the correlation between the $E_{\textrm{sym}}(\rho_{0})$ and $L$ along their lower boundaries through the arbitrary density $\mu$. Setting $J_{n}$=0, the Eqs. (7) and (How tightly is the nuclear symmetry energy constrained by a unitary Fermi gas?) reduce exactly to the parametric equations of $E_{\textrm{sym}}(\rho_{0})$ and $L$ derived slightly differently in version-1 of ref. Kolomeitsev16 . We note that the quantities that determine the boundary of allowed values of $E_{\textrm{s}ym}(\rho_{0})$ and $L$ are the total curvature parameter $K_{n}$ and total skewness parameter $J_{n}$ also emphasized in version-2 of ref. Kolomeitsev16 .

While having noted that $K_{\textrm{sym}}$ is experimentally and theoretically poorly known, the $E_{\textrm{sym}}(\rho_{0})$ versus $L$ correlation along their boundaries was obtained in ref. Kolomeitsev16 by setting $K_{\textrm{sym}}=0$ based on the prediction of a chiral effective field theory. It was found that the resulting correlation excludes many of the currently actively used models for $E_{\textrm{sym}}(\rho)$ . We reexamine this correlation by varying the $\xi$, $J_{n}$ and $K_{\textrm{sym}}$ within their known uncertain ranges. Again, the value of $\xi$ is now well settled around $0.37\pm 0.005$. Taking $K_{\textrm{sym}}=0$, $J_{n}=0$ and $K_{0}=230$ MeV, the two red dashed lines obtained with $\xi=0.37\pm 0.005$ in both (a) and (b) of Fig. 3 show the resulting lower boundaries of the $E_{\textrm{sym}}(\rho_{0})$ versus $L$ correlation.

The skewness coefficients $J_{0}$ and $J_{\textrm{sym}}$ in $J_{n}$ are both poorly known. We show in Fig 1(b) values of $J_{n}$ against $K_{n}$ for the 275 Skyrme and RMF models. $J_{n}$ varies approximately in the range $-500$ MeV $\leq J_{n}\leq 1000$MeV. To our best knowledge, there is no experimental constraint available on this quantity. $K_{n}$ varies approximately in the range $-150$ MeV $\leq K_{n}\leq 370$ MeV. Most of this comes from the big uncertainties in determining the value of $K_{\textrm{sym}}$, which are discussed in detail in ref. Li2017 . This is partially because the $K_{\textrm{sym}}$ depends on not only $L$ but also its derivative $(\mathrm{d}L/\mathrm{d}\rho)_{\rho_{0}}$ by definition. Microscopically, it depends on not only the nucleon isoscalar effective mass $m^{*}_{0}$ and neutron-proton effective mass slitting $m_{n}^{*}-m_{p}^{*}$ but also their momentum and density dependences that are all essentially completely unknown Li2017 . The latest calculations within many Skyrme Hartree-Fock and/or relativistic mean-field models indicate that $-400\leq K_{\textrm{sym}}\leq 100$ MeV Lwchen09 ; Dutra12 ; Dutra14 ; Colo14 .

The results using $K_{\textrm{sym}}=0$, $K_{0}=230$ MeV (so $K_{n}=$ 230 MeV) and $\xi=0.37$ for different values for $J_{n}$ are shown in Fig. 3 (a). It is seen that $L$ becomes larger as $J_{0}$ decreases, and that the upper boundary of the allowed region will correspond to the lower limit of $J_{n}$. At $E_{\textrm{sym}}(\rho_{0})=40$ MeV, for example, the increase is only about $20\%$.

By setting $K_{\textrm{sym}}=-200,-100$, [math], and $100$ MeV with $K_{0}=230$ MeV (corresponding to $K_{n}=30,130,230$ and $330$ MeV) we can illustrate effects of the $K_{\textrm{sym}}$ on the lower boundary of $E_{\textrm{sym}}(\rho_{0})$ versus $L$ correlation in Fig. 3 (b). It is seen that the $K_{\textrm{sym}}$ affects the results significantly; as $K_{\textrm{s}ym}$ (and hence $K_{n}$) increases, $L$ increases, and the upper boundary of the allowed region will correspond to the largest value of $K_{\textrm{s}ym}$ (and hence $K_{n}$). The overall uncertainty in $K_{n}$ leads to a variation of the upper boundary of the allowed region that is about twice that caused by the uncertainty in $J_{n}$. We note that many of the models allowed by $K_{\textrm{sym}}=0$ MeV would be excluded by using $K_{\textrm{sym}}=-200$ MeV.

Adopting the following values for the five parameters after taking into account the full range presented by the 275 Skyrme and RMF models, i.e.,

[TABLE]

a lower boundary excluding only the shadowed region in Fig. 3 is obtained. It is seen that only the TMA and NL$\rho\delta$, NL3 and LS220 may be excluded, while the STOS, TM1, NL$\rho$, LS220, and KVR, which have been surely excluded previously in ref. Kolomeitsev16 may be allowed. This is mainly a result in extending the upper bound in the uncertainty region of $K_{\textrm{s}ym}$, since the additional uncertainty in $J_{n}$ moves the excluded region to the right, (because $J_{n}\gtrsim 0$). As pointed out in the final analyses of ref. Kolomeitsev16 , we agree that strong empirical correlations among $J_{\textrm{sym}},K_{\textrm{sym}}$ and $L$ exist. These correlations can be used to refine the $E_{\textrm{sym}}(\rho_{0})$-L constraint shown in Fig. 3.

It is well known that the detailed density dependence of nuclear symmetry energy $E_{\textrm{sym}}(\rho)$ contains many interesting and some unknown physics. It is probably not surprising that the estimation of the correlation between its zeroth-order and first-order density expansion coefficients $E_{\textrm{sym}}(\rho_{0})$ and $L$ depends on what we assume about the immediate next high-order term characterized by the curvature $K_{\textrm{sym}}$ at $\rho_{0}$.

Concluding remarks: The universal EOS $E_{\textrm{UG}}$ of the unitary Fermi gas was conjectured in ref. Kolomeitsev16 to provide the lower boundary of the EOS $E_{\textrm{PNM}}$ of PNM, and thus a constraint on nuclear symmetry energy. Although unproven, the conjecture has strong empirical supports, and its implications are important enough that they should be examined rigorously. We found that the $E_{\textrm{UG}}$ does provide a useful lower boundary of nuclear symmetry energy. Moreover, this boundary is essentially not affected by the known uncertainty in the skewness coefficient $J_{0}$ of SNM. However, it does not tightly constrain the correlation between the magnitude $E_{0}(\rho_{0})$ and slope $L$ unless the curvature $K_{\textrm{sym}}$ of the $E_{\textrm{sym}}(\rho)$ at $\rho_{0}$ and the skewness parameters $J_{0}$ and $J_{\textrm{s}ym}$ are better known. Of these, $K_{\textrm{sym}}$ is the more important quantity, its uncertainty affecting the lower boundary of $E_{\textrm{sym}}(\rho)$ by up to twice as much as the uncertainty in $J_{n}$. Most of the previously excluded $E_{\textrm{sym}}(\rho)$ functionals by the universal EOS of unitary Fermi gas assuming $K_{\textrm{sym}}=0$ may not be excluded considering the currently known big uncertainties of the $K_{\textrm{sym}}(\rho_{0})$.

For many purposes in both nuclear physics and astrophysics, it is necessary to map out precisely both the $E_{\textrm{sym}}(\rho)$ and the EOS $E_{0}(\rho)$ of SNM in a broad density range. Near the saturation density $\rho_{0}$, this requires accurate knowledge of the $K_{\textrm{sym}}$ and $J_{0}$ besides the $E_{0}(\rho_{0})$, $L$ and $K_{0}$. To this end, it is interesting to mention briefly quantities that are sensitive to the higher-order EOS parameters and current efforts to determine them. For example, the skewness coefficient $J_{0}$ characterizes the high-density behavior of $E_{0}(\rho)$, and it has been found to affect significantly the maximum mass of neutron stars Cai14 . Moreover, at the crust-core transition point where the incompressibility of neutron star matter at $\beta$-equilibrium vanishes, the value of $J_{0}$ influences significantly the exact location of the transition point NBZ . Thus, astrophysical observations of neutron stars can potentially constrain the $J_{0}$ albeit probably not before other EOS parameters are well determined. On the other hand, in terrestrial laboratory experiments, there have been continued efforts to determine the $K_{\textrm{sym}}$ Colo14 . One outstanding example is the measurement of the isospin dependence of nuclear incompressibility $K(\delta)\approx K_{0}+K_{\tau}\delta^{2}+\mathcal{O}(\delta^{4})$ where $K_{\tau}=K_{\textrm{sym}}-6L-J_{0}L/K_{0}$ using giant resonances of neutron-rich nuclei TLi ; Jorge-TLi . While the current estimate of $K_{\tau}\approx-550\pm 100$ MeV Colo14 from analyzing many different kinds of terrestrial experiments is still too rough to constrain tightly the individual values of $J_{0}$ and $K_{\textrm{sym}}$, new experiments with more neutron-rich beams have the promise of improving significantly the accuracy of the measured $K_{\tau}$ Umesh . Thus, we are hopeful that not only the zeroth and first-order parameters $K_{0}$, $E_{\textrm{sym}}(\rho_{0})$ and $L$ but also high-order coefficients $J_{0}$ and $K_{\textrm{sym}}$ can be pinned down in the near future by combining new analyses of upcoming astrophysical observations and terrestrial experiments.

Acknowledgements: We would like to thank Umesh Garg for helpful 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, the National Natural Science Foundation of China under Grant No. 11320101004 and the Texas Advanced Computing Center. 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, the Shanghai Key Laboratory of Particle Physics and Cosmology under Grant No. 15DZ2272100, and the Shanghai Pujiang Program under Grant No. 13PJ1410600.