Equation of state effects in the core collapse of a $20$-$M_\odot$ star

TL;DR

This study explores how uncertainties in dense matter equations of state influence core-collapse supernovae and neutron star properties, highlighting the critical role of nucleon effective mass in supernova dynamics and neutrino signals.

Contribution

We constructed 97 new finite-temperature EOSs consistent with current constraints and systematically analyzed their impact on supernova and neutron star outcomes.

Findings

Nucleon effective mass at high densities is the main uncertainty affecting neutrino signals.

Higher effective mass at nuclear saturation density increases the likelihood of supernova explosion.

Changes in other EOS parameters mainly influence cold neutron star properties, not supernova dynamics.

Abstract

Uncertainties in our knowledge of the properties of dense matter near and above nuclear saturation density are among the main sources of variations in multi-messenger signatures predicted for core-collapse supernovae (CCSNe) and the properties of neutron stars (NSs). We construct 97 new finite-temperature equations of state (EOSs) of dense matter that obey current experimental, observational, and theoretical constraints and discuss how systematic variations in the EOS parameters affect the properties of cold nonrotating NSs and the core collapse of a $20\,M_\odot$ progenitor star. The core collapse of the $20\,M_\odot$ progenitor star is simulated in spherical symmetry using the general-relativistic radiation-hydrodynamics code GR1D where neutrino interactions are computed for each EOS using the NuLib library. We conclude that the effective mass of nucleons at densities above nuclear saturation density is the largest source of uncertainty in the CCSN neutrino signal and dynamics even though it plays a subdominant role in most properties of cold NS matter. Meanwhile, changes in other observables affect the properties of cold NSs, while having little effect in CCSNe. To strengthen our conclusions, we perform six octant three-dimensional CCSN simulations varying the effective mass of nucleons at nuclear saturation density. We conclude that neutrino heating and, thus, the likelihood of explosion is significantly increased for EOSs where the effective mass of nucleons at nuclear saturation density is large.

Tables

Table 1. Table 1: Constraints of nuclear matter properties used in this work grouped in sets defined in Sec.  II.3 . Nuclear matter empirical parameters were compiled in Ref.  Margueron, Hoffmann Casali, and Gulminelli ( 2018a ) , see references therein for details. Meanwhile, nuclear matter pressure at 4 ​ n sat 4 subscript 𝑛 sat 4n_{\rm sat} , P ( 4 ) superscript 𝑃 4 P^{(4)} , for SNM and PNM is from Ref.  Danielewicz, Lacey, and Lynch ( 2002 ) . We use values similar to the ones in Refs.  Margueron, Hoffmann Casali, and Gulminelli ( 2018a ); Danielewicz, Lacey, and Lynch ( 2002 ) , but exclude from our analysis regions of parameter space that fail to reproduce 2 2 2 - M ⊙ subscript 𝑀 direct-product M_{\odot} NSs and, in the case of L sym subscript 𝐿 sym L_{\rm sym} , values that lead to too large radii for NSs Nättilä  et al. ( 2016 ) . We show the averages and one-standard deviations compiled or assumed in this work.

Set	Quantity	Range	This work	Units
$s_M$	$m^{\star }$	$\mathrm{0.75\pm 0.10}$	$\mathrm{0.75\pm 0.10}$	$m_n$
	$\Delta m^{\star }$	$\mathrm{0.10\pm 0.10}$	$\mathrm{0.10\pm 0.10}$	$m_n$
$-$	$n_{\mathrm{sat}}$	$\mathrm{0.155\pm 0.005}$	$0.155$	${\mathrm{fm}}^{-3}$
	${ϵ}_{\mathrm{sat}}$	$-\mathrm{15.8\pm 0.3}$	$-15.8$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$s_S$	${ϵ}_{\mathrm{sym}}$	$\mathrm{32\pm 2}$	$\mathrm{32\pm 2}$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
	$L_{\mathrm{sym}}$	$\mathrm{60\pm 15}$	$\mathrm{45\pm 7.5}$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$s_K$	$K_{\mathrm{sat}}$	$\mathrm{230\pm 20}$	$\mathrm{230\pm 15}$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
	$K_{\mathrm{sym}}$	$-\mathrm{100\pm 100}$	$-\mathrm{100\pm 100}$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$s_P$	$P_{\mathrm{SNM}}^{\left(4\right)}$	$\mathrm{100\pm 50}$	$\mathrm{125\pm 12.5}$	$\mathrm{MeV}{\mathrm{fm}}^{-3}$
	$P_{\mathrm{PNM}}^{\left(4\right)}$	$\mathrm{160\pm 80}$	$\mathrm{200\pm 20}$	$\mathrm{MeV}{\mathrm{fm}}^{-3}$

Table 2. Table 2: Zero-temperature properties of the SLy4 and LS220 EOSs.

Quantity	SLy4	LS220	Units
$m^{\star }$	$0.694$	$1.000$	$m_n$
$\Delta m^{\star }$	$-0.185$	$0.000$	$m_n$
$n_{\mathrm{sat}}$	$0.1595$	$0.1549$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
${ϵ}_{\mathrm{sat}}$	$-15.97$	$-16.00$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
${ϵ}_{\mathrm{sym}}$	$32.00$	$28.61$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$L_{\mathrm{sym}}$	$45.96$	$73.81$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$K_{\mathrm{sat}}$	$229.90$	$219.84$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$K_{\mathrm{sym}}$	$-119.70$	$-24.04$	$\mathrm{MeV}{\mathrm{baryon}}^{-1}$
$P_{\mathrm{SNM}}^{\left(4\right)}$	$127.12$	$107.75$	$\mathrm{MeV}{\mathrm{fm}}^{-3}$
$P_{\mathrm{PNM}}^{\left(4\right)}$	$142.15$	$162.08$	$\mathrm{MeV}{\mathrm{fm}}^{-3}$

Equations

$${ϵ}_B(n,y,T)={ϵ}_{\mathrm{kin}}(n,y,T)+{ϵ}_{\mathrm{pot}}(n,y).$$

$${ϵ}_{\mathrm{kin}}(n,y,T)=\frac{1}{n}\left(\frac{{\hslash }^2{\tau }_n}{2m_n^{\star }}+\frac{{\hslash }^2{\tau }_p}{2m_p^{\star }}\right),$$

$${\tau }_t=\frac{1}{2{\pi }^2}{\left(\frac{2m_t^{\star }T}{{\hslash }^2}\right)}^{\frac{5}{2}}{\mathcal{F}}_{3/2}({\eta }_t),$$

$$\frac{{\hslash }^2}{2m_t^{\star }}=\frac{{\hslash }^2}{2m_t}+{\alpha }_1{n}_t+{\alpha }_2{n}_{-t}.$$

$${\mathcal{F}}_k(\eta )={\int }_0^{\infty }\frac{u^{k}du}{1+\exp (u-\eta )},$$

$${\eta }_t=\frac{{\mu }_t-V_t}{T}.$$

$${ϵ}_{\mathrm{pot}}(n,y)=\sum _{i=0}^N\left[a_i+4b_{i}y(1-y)\right]n^{{\delta }_i}.$$

$${ϵ}_B(n,y)={ϵ}_{\mathrm{is}}(x)+{\delta }^2{ϵ}_{\mathrm{iv}}(x),$$

$${ϵ}_{\mathrm{is}}(x)$$

$${ϵ}_{\mathrm{iv}}(x)$$

$$+\frac{1}{3!}{Q}_{\mathrm{sym}}{x}^3+\dots ,$$

$$s_M$$

$$s_S$$

$$s_K$$

$$s_P$$

$$s_t\approx {\left(\frac{\pi }{\hslash }\right)}^2\frac{m_t^{\star }T}{(3{\pi }^2{n}_t{)}^{2/3}}.$$

$$\underset{T\to 0}{\lim }{\mathcal{F}}_k(\eta )={\int }_0^{\eta }{u}^{k}du+\frac{{\pi }^2}{6}{T}^2{\left(\frac{d(u^k)}{du}\right)}_{\eta }+\dots .$$

$$n_t\simeq$$

$${\tilde{\mu }}_t={\mu }_{tF}\left[1-\frac{{\pi }^2}{12}{\left(\frac{T}{{\mu }_{tF}}\right)}^2\right],$$

$${\tau }_t\simeq {\tau }_{tF}\left[1+\frac{5{\pi }^2}{12}{\left(\frac{T}{{\mu }_{tF}}\right)}^2\right],$$

$$s_B$$

$$\simeq \frac{T}{n}\frac{{\pi }^2}{2}\sum _t\left[\frac{n_t}{{\mu }_{tF}}\right]=\sum _t\frac{n_t}{n}{s}_t,$$

$$P_{th}$$

$$({\alpha }_2+{\alpha }_1)$$

$${\beta }_{\Delta }$$

$$A=\left[\begin{array}{cccccccc}{\mathfrak{a}}_0& {\mathfrak{a}}_0& {\mathfrak{a}}_1& {\mathfrak{a}}_1& {\mathfrak{a}}_2& {\mathfrak{a}}_2& {\mathfrak{a}}_3& {\mathfrak{a}}_3\\ {\mathfrak{a}}_0^{\prime }& {\mathfrak{a}}_0^{\prime }& {\mathfrak{a}}_1^{\prime }& {\mathfrak{a}}_1^{\prime }& {\mathfrak{a}}_2^{\prime }& {\mathfrak{a}}_2^{\prime }& {\mathfrak{a}}_3^{\prime }& {\mathfrak{a}}_3^{\prime }\\ 0& -{\mathfrak{a}}_0& 0& -{\mathfrak{a}}_1& 0& -{\mathfrak{a}}_2& 0& -{\mathfrak{a}}_3\\ 0& -3{\mathfrak{a}}_0^{\prime }& 0& -3{\mathfrak{a}}_1^{\prime }& 0& -3{\mathfrak{a}}_2^{\prime }& 0& -3{\mathfrak{a}}_3^{\prime }\\ {\mathfrak{a}}_0^{\prime \prime }& {\mathfrak{a}}_0^{\prime \prime }& {\mathfrak{a}}_1^{\prime \prime }& {\mathfrak{a}}_1^{\prime \prime }& {\mathfrak{a}}_2^{\prime \prime }& {\mathfrak{a}}_2^{\prime \prime }& {\mathfrak{a}}_3^{\prime \prime }& {\mathfrak{a}}_3^{\prime \prime }\\ 0& -{\mathfrak{a}}_0^{\prime \prime }& 0& -{\mathfrak{a}}_1^{\prime \prime }& 0& -{\mathfrak{a}}_2^{\prime \prime }& 0& -{\mathfrak{a}}_3^{\prime \prime }\\ {\mathfrak{b}}_0& {\mathfrak{b}}_0& {\mathfrak{b}}_1& {\mathfrak{b}}_1& {\mathfrak{b}}_2& {\mathfrak{b}}_2& {\mathfrak{b}}_3& {\mathfrak{b}}_3\\ 0& -{\mathfrak{b}}_0& 0& -{\mathfrak{b}}_1& 0& -{\mathfrak{b}}_2& 0& -{\mathfrak{b}}_3\end{array}\right],$$

$${\mathfrak{a}}_i$$

$${\mathfrak{a}}_i^{\prime }$$

$${\mathfrak{a}}_i^{\prime \prime }$$

$${\mathfrak{b}}_i$$

$$B=\left[\begin{array}{c}{ϵ}_{\mathrm{sat}}-{ϵ}_{\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ n_{\mathrm{sat}}^{-1}{P}_{\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ {ϵ}_{\mathrm{sym}}-{ϵ}_{\mathrm{sym},\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ L_{\mathrm{sym}}-L_{\mathrm{sym},\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ K_{\mathrm{sat}}-K_{\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ K_{\mathrm{sym}}-K_{\mathrm{sym},\mathrm{kin}}(n_{\mathrm{sat}},0.5)\\ P_{\mathrm{SNM}}^{(4)}-P_{\mathrm{kin}}(4n_{\mathrm{sat}},0.5)\\ P_{\mathrm{PNM}}^{(4)}-P_{\mathrm{kin}}(4n_{\mathrm{sat}},0)\end{array}\right].$$

$$P_{\mathrm{kin}}(n,y)$$

$$K_{\mathrm{kin}}(n,y)$$

$$K_{\mathrm{sym},\mathrm{kin}}(n,y)$$

$${ϵ}_{\mathrm{sym},\mathrm{kin}}(n,y)$$

$$L_{\mathrm{sym},\mathrm{kin}}(n,y)$$

Full text

I Introduction

Stars with masses above roughly eight times the mass of the Sun ($M_{\odot}$), end their lives in a core collapse event, in many cases leading to a core-collapse supernova (CCSN) explosion. Core collapse sets in once electron degeneracy pressure in the nickel-iron core of a massive star can no longer support it against gravity Woosley, Heger, and Weaver (2002).

Core collapse proceeds until the inner core reaches nuclear saturation density, $\rho_{\rm sat}\gtrsim 2.7\times 10^{14}\,\mathrm{g\,cm}^{-3}$, at a temperature of $10-20\,\mathrm{MeV}$. At this point, the residual nuclear force prevents the inner core from contracting any further and it rebounds into the still infalling outer core, creating a shock wave. As the shock wave propagates through the outer core it eventually stalls because of energy losses resulting from dissociation of heavy nuclei and to lesser extent due to neutrino losses from behind the shock.

A few mechanisms that revive the shock and lead to successful CCSNe have been suggested, see discussion in Ref. Burrows (2013) and references therein. Simulations have shown that it is likely that a multitude of macroscopic (e.g., progenitor structure, large-scale convection, magnetohydrodynamic forcing) and microscopic properties and processes (e.g., neutrino heating) couple non-linearly to drive an explosion. Still, it is believed that the main contributor to shock revival is the neutrino heating mechanism Wilson (1985); Bethe and Wilson (1985), whereby $\sim 10\%$ of the outgoing electron-flavor neutrino luminosity is deposited behind the shock. This provides the shock with thermal support, drives turbulence, and aids in shock runaway Janka, H.-Th. (2001); Couch and Ott (2015); Burrows (2013).

One of the fundamental ingredients to understand the dynamics of core collapse events is the equation of state (EOS) of dense matter. The density at which the collapse halts, how many protons are converted into neutrons during the collapse, the spectra of neutrinos, how much energy is deposited behind the shock and its expansion rate, the ejecta mass and its composition, the proto neutron star (PNS) mass, its radius, cooling rate, and whether it later collapses into a black-hole (BH) as well as the gravitational wave (GW) signal, are all dependent on the EOS. In a CCSN, and also in NS mergers, matter exists in a wide range of temperatures, $0\lesssim{T}\lesssim\mathcal{O}(100\,\mathrm{MeV})$, densities, ${\rho}\lesssim 10^{15}\,\mathrm{g\,cm}^{-3}$, and proton fractions, $0.0\lesssim{y}\lesssim 0.5$. Some of these conditions are so extreme they are not readily available to laboratory experiments and, thus, such regions of parameter space can only be probed indirectly from observations in consent with computational and theoretical models.

Recently Ref. Margueron, Hoffmann Casali, and Gulminelli (2018a) introduced the concept of meta-modeling for the nuclear EOS (see also Ref. Steiner et al. (2005)). In their model, the EOS is parametrized in terms of empirical parameters, i.e., nuclear matter binding energy, saturation density, incompressibility, symmetry energy, and so on. The average values of the empirical parameters and their uncertainties are estimated based on experimental and theoretical nuclear physics constraints. In follow-up studies meta-modeling was used to study the effects of uncertainties in the empirical parameters on NS properties Margueron, Hoffmann Casali, and Gulminelli (2018b), finite size effects in the description of nuclear masses and radii of ground state nuclei Chatterjee et al. (2017), and to compute correlations between empirical parameters from known constraints Margueron and Gulminelli (2019).

We follow the meta-modeling approach Margueron, Hoffmann Casali, and Gulminelli (2018a, b) and analyze how uncertainties in properties of nuclear matter affect cold NS properties and the core collapse of a $20$-$M_{\odot}$ progenitor star. We use the meta-modeling formalism to construct a family of finite temperature EOSs of dense matter. The EOSs are built using the recently developed open-source SROEOS code Schneider, Roberts, and Ott (2017), which is itself based on the Lattimer and Swesty liquid-drop model of nuclei Lattimer and Swesty (1991), with a few improvements. The main improvements relevant to this work are the possibility to compute EOSs where (1) the effective mass of nucleons is different from their vacuum values and (2) for any desired value of the incompressibility of nuclear matter $K_{\rm sat}$, instead of the canonical values of $180\,\mathrm{MeV}$, $220\,\mathrm{MeV}$, and $375\,\mathrm{MeV}$ to which the code of Lattimer and Swesty (1991) is essentially limited. The SROEOS model has also been extended to transition to a description of many nuclear species in nuclear statistical equilibrium (NSE) at low densities.

The main goal of this study is to separately determine how each empirical parameter of the EOS may affect a core collapse event and the resulting PNS. This is only possible using many EOSs obtained within a single formalism. Previous studies have studied the effect of the EOS on CCSNe and their observables, e. g., Refs. O’Connor and Ott (2011); Hempel et al. (2012); Steiner, Hempel, and Fischer (2013); Couch (2013); Suwa et al. (2013); Char, Banik, and Bandyopadhyay (2015); Fischer et al. (2014); Furusawa et al. (2017); Richers et al. (2017); Nagakura et al. (2018); Morozova et al. (2018). The main drawback of these studies is that often the EOSs being compared were obtained with distinct approaches, used different prescriptions to describe low density matter, and, with the exception of Ref. Richers et al. (2017) which analyzed changes resulting from using 18 different EOSs in their simulations, the number of EOSs investigated was rather small. Thus, in many cases, it was challenging to disentangle how a parameter of the EOS contributed to a given observable.

In this paper, we focus on EOS effects on the neutrino heating mechanism and delay the study of GW signals to future work. We simulate the core collapse of a single non-rotating $20$-$M_{\odot}$ progenitor star taken from Woosley and Heger (2007) using 97 distinct EOSs that each vary in at most two different empirical parameters from a baseline EOS. The SROEOS code is ideal for this type of sensitivity study as it allows one to compute many EOSs within the same framework using arbitrary Skyrme-type parametrizations of the nuclear forces. Furthermore, to limit our assessment only to the effects of the high-density part of the EOS, we use the same nuclear surface parametrization for all EOSs and the same NSE EOS at low densities for all simulations. The CCSN simulations are performed using the open-source general-relativistic multi-group radiation-hydrodynamics code GR1D O’Connor and Ott (2010, 2011) and the NuLib neutrino transport libraries O’Connor and Ott (2013); O’Connor (2015). Since the GR1D code is limited to spherical symmetry, we also perform six three-dimensional (3D) simulations, limited to an octant of the 3D cube to keep computational demands manageable. For this, we employ the open-source 3D general-relativistic radiation-hydrodynamics code Zelmani Roberts et al. (2016); Ott et al. (2018), which is based on the Einstein Toolkit Löffler et al. (2012); Mösta et al. (2014). We perform the octant 3D runs for five variations of the SLy4 EOS Chabanat et al. (1998) and the LS220111LS220 is the Lattimer & Swesty EOS with incompressibility $K_{\rm sat}=220\,\mathrm{MeV}$ EOS. In this work the LS220 EOS was recomputed using the SROEOS code Schneider, Roberts, and Ott (2017). EOS Lattimer and Swesty (1991).

This paper is structured as follows. In Section II, we discuss a variant of the meta EOS model of Ref. Margueron, Hoffmann Casali, and Gulminelli (2018a) that suits our needs. We proceed to discuss how each of the empirical parameters affects the properties of cold beta-equilibrated NSs in Sec. III and spherically-symmetric core collapse in Sec. IV. In Sec. V, we discuss 3D runs with octant symmetry. We conclude in Sec. VI.

II Meta EOS

Motivated by Ref. Margueron, Hoffmann Casali, and Gulminelli (2018a), we use a metamodeling formalism to compute Skyrme parameters for the nucleonic EOS in terms of empirical nuclear parameters. In this work, matter is assumed to be made solely of nucleons, electrons, positrons, and photons222It is expected that the EOS softens at very high temperatures and densities due to the appearance of heavy leptons, hyperons, condensates, and quark-gluon plasmas Sagert et al. (2009). They are not explicitly included here since we take a parameterized approach to the high-density EOS. . Electrons, positrons, and photons are treated as uniform free gases and charge neutrality is assumed. Therefore, their contributions to the EOS decouple from the nucleon contributions. Our treatment of these components of the EOS is discussed in detail in Appendix A of Ref. Schneider, Roberts, and Ott (2017).

II.1 Skyrme Model

The bulk nuclear contribution to the EOS is computed assuming non-relativistic effective Skyrme-type nucleon-nucleon interactions. In this approach the energy per baryon $\epsilon_{B}$ of nucleonic matter with number density $n$ and proton fraction $y$ can be separated into its kinetic and potential energy density contributions i.e.,

[TABLE]

The kinetic energy density term is

[TABLE]

where

[TABLE]

and the density dependent effective nucleon masses $m_{t}^{\star}$ are given by

[TABLE]

Above, $n_{t}$ and $m_{t}$ are, respectively, the density and vacuum mass of a nucleon with isospin $t$, where $t=n$ for neutrons and $t=p$ for protons, and, if $t=n$ then $-t=p$ and vice-versa. The neutron and proton densities are related to the proton fraction $y$ and the nucleon density $n$ by $n_{n}=(1-y)n$ and $n_{p}=yn$, respectively. The quantities $\alpha_{1}$ and $\alpha_{2}$ are parameters of the model and establish a simple dependence of the nucleon effective masses on the density and proton fraction of the system. We stress that the Skyrme model treatment of effective masses is rudimentary, other models allow for much more complex dependencies of $m^{\star}$ Constantinou et al. (2015). Nevertheless, we use this model as a guide to teach us how each piece of the EOS affects neutron star (NS) properties and the dynamics of CCSNe.

The Fermi integral in Eq. (3) is defined as

[TABLE]

and is a function of the degeneracy parameter

[TABLE]

Here, $\mu_{t}$ is the nucleon chemical potential and $V_{t}$ is the single-particle potential (see Schneider, Roberts, and Ott (2017) for more details).

The temperature-independent potential energy density term in Eq. (1) has the form

[TABLE]

where $a_{i}$, $b_{i}$, and $\delta_{i}$ are constant parameters of the Skyrme model. The $i=0$ term is chosen to represent two-body nucleon interactions. Therefore, we fix $\delta_{0}=1$ for all models. Meanwhile, the $i>0$ terms approximate effects of many-body interactions Lattimer and Swesty (1991). The summation in most Skyrme parameterizations ends at $N=1$, while only a small number of studies in the literature consider $N>1$ Dutra et al. (2012). To allow for more flexibility in our empirically fitted models, we choose to fix $N=3$ and $\delta_{0}=1$, $\delta_{1}=4/3$, $\delta_{2}=2$, and $\delta_{3}=7/3$ (the last three terms amount to an expansion in terms of the Fermi momenta of the nucleons $k_{t}\propto n_{t}^{2/3}$ Lim and Holt (2017); Tews et al. (2017)). Therefore, the EOS model contains ten free parameters $\{a_{0},b_{0},a_{1},b_{1},a_{2},b_{2},a_{3},b_{3},\alpha_{1},\alpha_{2}\}$ that we fit using a set of empirical properties of nuclear matter.

II.2 Empirical Parameters

Now, we would like to define a set of empirical properties with which to constrain our Skyrme EOS parameters. First, we consider measurable properties of nearly symmetric nuclear matter near nuclear saturation density. In these conditions, the zero-temperature nuclear EOS can be expanded about nuclear saturation density, $n=n_{\rm sat}\simeq 0.155\,\mathrm{fm}^{-3}$, for symmetric matter ($y=1/2$) in a Taylor series, giving rise to a set of expansion parameters that can be empirically constrained. This expansion is written as

[TABLE]

where $x=(n-n_{\rm sat})/(3n_{\rm sat})$ and $\delta=1-2y$ is the isospin asymmetry. Here, the isoscalar (is) and isovector (iv) expansion terms are Margueron, Hoffmann Casali, and Gulminelli (2018a); Piekarewicz and Centelles (2009)

[TABLE]

shown here explicitly up to third order. The empirical parameter $\epsilon_{\rm sat}$ is the energy per baryon at nuclear saturation density $n_{\rm sat}$, $K_{\rm sat}$ is the isoscalar incompressibility modulus, and $Q_{\rm sat}$ the isoscalar skewness. Similarly, $\epsilon_{\rm sym}$ is the symmetry energy, $L_{\rm sym}$ is related to the slope of symmetry energy in the direction of increasing density, $K_{\rm sym}$ is the isovector incompressibility modulus, and $Q_{\rm sym}$ is the isovector skewness. By definition of the saturation density $n_{\rm sat}$, the linear term in $x$ of $\epsilon_{\rm is}$ vanishes. In principle, all of these expansion parameters can be determined experimentally, with varying degrees of difficulty. Nevertheless, the lower-order parameters are substantially easier to constrain. Therefore, we only include the well constrained saturation density empirical parameters $\{n_{\text{sat}},\epsilon_{\text{sat}},K_{\text{sat}},\epsilon_{\text{sym}},L_{\text{sym}},K_{\text{sym}}\}$ in our Skyrme model fits described below.

Although this expansion is useful near saturation density, it cannot accurately describe the behavior of the nuclear EOS at densities larger than a few times saturation density since $x$ is no longer small and the expansion breaks down. Densities this large are reached in CCSNe and in the cores of NSs. Therefore, we also require empirical constraints at higher density. Most experiments probe densities near saturation density, but there are some results available for higher densities. Using measurements of flow in heavy ion collisions and theoretical transport models, Danielewicz, Lacey, and Lynch (2002) constrained the baryonic pressure $P_{B}=n^{2}\partial\epsilon_{B}/\partial n$ of symmetric nuclear matter (SNM) and pure neutron matter (PNM), albeit in a model dependent way, at four times nuclear saturation density, $P_{\text{SNM}}^{(4)}=P_{B}(n=4n_{\text{sat}},y=1/2)$ and $P_{\text{PNM}}^{(4)}=P_{B}(n=4n_{\text{sat}},y=0)$. Constraints on these pressures have recently been made sharper by combining the results of these flow experiments with constraints on the tidal deformability of NSs inferred from GW170817 Tsang et al. (2018).

Finally, although they do not enter into the expansion above, the nucleon effective masses at saturation density can also be considered a quasi-empirical parameter Margueron, Hoffmann Casali, and Gulminelli (2018a). However, there is considerable complexity involved in extracting this property of the single quasi-particle energies. Nevertheless, the nucleon effective masses are particularly important for determining the temperature dependence of the nuclear EOS [see Eq. (2) above]. Therefore, we include the nucleon effective mass at saturation density in SNM, $m^{\star}\equiv m_{n}^{\star}(n=n_{\rm sat},y=1/2)$, and the neutron-proton effective mass splitting in PNM, $\Delta m^{\star}\equiv m_{n}^{\star}(n=n_{\rm sat},y=0)-m_{p}^{\star}(n=n_{\rm sat},y=0)$, in our list of empirical parameters.

In total, this gives ten empirical parameters that we consider in this work, $\{n_{\text{sat}},\epsilon_{\text{sat}},K_{\text{sat}},\epsilon_{\text{sym}},L_{\text{sym}},K_{\text{sym}},m^{\star},\Delta m^{\star},P_{\text{PNM}}^{(4)},P_{\text{SNM}}^{(4)}\}$. Due to their small uncertainties, we fix the values of the nuclear saturation number density $n_{\rm sat}=0.155\,\mathrm{fm}^{-3}$ (mass density $\rho_{\rm sat}=2.7\times 10^{14}\,\mathrm{g\,cm}^{-3}$) and of the energy at nuclear saturation density $\epsilon_{\rm sat}=-15.8\,\mathrm{MeV}$. Other saturation density quantities are allowed to vary within their experimental or theoretical uncertainties (as compiled in Margueron, Hoffmann Casali, and Gulminelli (2018a)) as long as they are able to produce $2$-$M_{\odot}$ NSs Demorest et al. (2010); J. Antoniadis et al. (2013); Fonseca et al. (2016). The exception to this choice is the slope of the symmetry energy $L_{\rm sym}$. Instead of using the average values of Ref. Margueron, Hoffmann Casali, and Gulminelli (2018a), $L_{\rm sym}=60\pm 15\,\mathrm{MeV}$, we set $L_{\rm sym}=45\pm 7.5\,\mathrm{MeV}$. Although this choice only probes the lower half of possible values compiled in Ref. Margueron, Hoffmann Casali, and Gulminelli (2018a), we choose these limits so that the mass-radius relationships of NSs in this work are centered near the center of the constraints computed from observations of x-ray bursts Nättilä et al. (2016). These limits also agree with combined theoretical calculations of pure neutron matter and astrophysical observations Lattimer and Lim (2013); Tews et al. (2017); Oertel et al. (2017). Even though $L_{\rm sym}$ is correlated with radii of low mass NSs Alam et al. (2016), for the systems we study, our limited choice for $L_{\rm sym}$ has little effect on PNS properties in the first second after core collapse. Finally, we ignore existing correlations between the different empirical nuclear matter parameters Alam et al. (2016); Tews et al. (2017); Margueron and Gulminelli (2019). Note, however, that the allowed ranges for empirical parameters contain EOSs that do not fulfill expected correlation between $\epsilon_{\rm sym}$ and $L_{\rm sym}$ determined on the basis of unitary gas considerations Tews et al. (2017). We justify our choice with our primary interest in how different parameters of the EOS affect CCSNe. Our focus is less on particularly intricate details of the EOS. In Tab. 1 we summarize the constraints used in this work.

II.3 Empirically Constrained Skyrme EOS Models

For a given set of Skyrme parameters, the empirical parameters described in the last section can be calculated from the Skyrme energy density [Eq. (1)], its derivatives, and the Skyrme expression for the effective masses [Eq. (4)]. Conversely, for a given choice of the ten empirical parameters given above, the ten Skyrme parameters are fixed. Our method for finding the Skyrme parameters from the empirical parameters is given in Appendix B. We stress that the fitted Skyrme parameterization only matches the saturation density expansion [Eq. (8)] at saturation density since the Skyrme model has a different functional form from the polynomial expansion.

To investigate the impact of EOS uncertainties on cold NSs and core collapse, we build a set of 97 Skyrme EOSs by picking 97 sets of the empirical parameters in the ranges given in Tab. 1. We initially set the quantities used to obtain the Skyrme parametrization to their average values. Then, two-sigma variations in the nuclear properties are implemented for four sets of nuclear properties with two quantities each. The sets are

[TABLE]

Thus, for set $s_{M}$ the values of $m^{\star}$ and $\Delta m^{\star}$ can be their average values ($m^{\star}=0.75$ and $\Delta m^{\star}=0.10$), or their average values plus or minus one standard deviation ($m^{\star}=0.75\pm 0.10$ and $\Delta m^{\star}=0.10\pm 0.10$) or two standard deviations ($m^{\star}=0.75\pm 0.20$ and $\Delta m^{\star}=0.10\pm 0.20$). Similar variations are implemented for all other sets, leading to a total of 97 different parametrization for the EOS333There are 25 EOSs in each set $s$. However, the baseline EOS with the average values of the observables is the same for all 4 sets.. For each of the parametrizations we build an EOS table using the open-source SROEOS code we have recently developed Schneider, Roberts, and Ott (2017).

II.4 Non-uniform and low density matter

To limit our focus to the effects of the empirical parameters on CCSNe, we set the same parametrization of the nuclear surface for all EOSs. This is different from what we presented in Ref. Schneider, Roberts, and Ott (2017), where the parametrization of the surface properties was computed self-consistently based on the Skyrme parameters. We defer to future work a detailed study of nuclear surface effects on CCSNe. Here, the surface parameters are chosen to be $\sigma_{s}=1.15\,\mathrm{MeV\,fm}^{-2}$, $q=16$, $\lambda=3.0$, and $p=1.5$, see Eqs. (19) and (20) in Ref. Schneider, Roberts, and Ott (2017). The surface parametrization chosen here leads to a surface symmetry energy $S_{S}=57.8\,\mathrm{MeV}$, in agreement with the value $S_{S}=58.9\pm 1.1\,\mathrm{MeV}$ of Ref. Jiang et al. (2012), and a surface level density $A_{S}=0.13\,\mathrm{MeV\,fm}^{-1}$.

Once empirical and surface parametrizations are set, we use the SROEOS code to obtain the EOS table. The EOSs in the Skyrme model are obtained in the single nucleus approximation (SNA) Lattimer and Swesty (1991); Schneider, Roberts, and Ott (2017) although extensions to accommodate multiple nuclear species have recently been proposed Grams et al. (2018); Furusawa and Mishustin (2018). We take the same approach discussed in our previous work and match our Skyrme-type EOSs to an EOS of 3,335 nuclei in nuclear statistical equilibrium (NSE) Schneider, Roberts, and Ott (2017). A unified method to connect SNA and NSE EOSs is the subject of Refs. Gulminelli and Raduta (2015); Raduta and Gulminelli (2019). Here we follow the simple prescription to transition between SNA to NSE EOSs using a density dependent function as discussed in Sec. VII A of Ref. Schneider, Roberts, and Ott (2017). Here, we set the transition parameters $n_{\rm tr}=10^{-3}\,\mathrm{fm}^{-3}$ and $n_{\delta}=0.33$, see Eqs. (57) and (58) of Ref. Schneider, Roberts, and Ott (2017). Note that the parameter $n_{\rm tr}$ is different from $n_{\rm tr}=10^{-4}\,\mathrm{fm}^{-3}$ used in Ref. Schneider, Roberts, and Ott (2017). The reason for this change is that the time to bounce in core collapse is insensitive to $n_{\rm tr}$ in the range $10^{-2}\,\mathrm{fm}^{-3}\gtrsim n_{\rm tr}\gtrsim 10^{-3}\,\mathrm{fm}^{-3}$, while it is a function of $n_{\rm tr}$ for $n_{\rm tr}<10^{-3}\,\mathrm{fm}^{-3}$. We note that setting $n_{\delta}\lesssim 0.5$ has little effect on CCSN simulations. However, larger values may have an effect since the SNA (NSE) EOS will have significant contributions at low (high) densities.

III Cold neutron stars

We study how variations in the empirical parameters of the EOS and of the pressure of nuclear matter at high densities affects the zero-temperature EOS and properties of cold nonrotating beta-equilibrated NSs using the suite of EOSs discussed in Section II. We consider each set of empirical parameter variations [see Eqs. (11)] separately.

III.1 Effective Mass

The Tolman-Oppenheimer-Volkoff (TOV) equations of NS structure only depend on the relationship between the pressure and energy density for the cold, beta-equilibrated EOS, $P_{\beta-\text{equil},T=0}(\epsilon_{B}n)$, where $n$ is the baryonic number density and $\epsilon_{B}$ the energy per baryon defined in Eq. (1). Since $\epsilon_{B}$ and its first few derivatives are fixed at saturation density by the empirical expansion parameters, varying only the effective masses, set $s_{M}$ defined in Eq. (11a), has a limited impact on $P_{\beta-\text{equil},T=0}(\epsilon_{B}n)$ and one expects small variations in the nonrotating NS mass radius relation444 Due to our choice of fixing the empirical parameters of order 2 and lower in Eq. (8) as well as the baryonic pressures for SNM and PNM at $4n_{\rm sat}$, the zero-temperature baryonic pressure, $P_{B}=n^{2}\partial\epsilon_{B}/\partial n$, is almost independent of $m^{\star}$ and $\Delta m^{\star}$. Small variations in the cold EOS for distinct choices of $m^{\star}$ and $\Delta m^{\star}$ result from how the Skyrme parameters, and, thus, the empirical parameters of order 3 and higher in Eq. (8), adjust to reproduce the fixed empirical parameters and the pressure at $4n_{\rm sat}$. Our method contrasts with the one in Ref. Yasin et al. (2018), where a large effect in the EOS and mass-radius relations of cold beta equilibrated NSs due to variations of the effective mass is observed. .

The limited impact of the effective masses on the zero-temperature EOS is visible in the first column of Fig. 1, where we plot the zero-temperature pressures of SNM (top) and PNM (bottom) as a function of density. No perceptible differences are seen for the EOS of SNM as the effective masses are changed. Meanwhile, only minor changes in the EOS of PNM occur for the different effective masses. As in the SNM case, the EOS of PNM is, by construction, within the bounds determined from flow experiments Danielewicz, Lacey, and Lynch (2002), since we fix the pressure of PNM at four times saturation density. There are two bands shown for the pressure of PNM where the lower (higher) pressure band represents the pressure of PNM considering the softest (stiffest) density dependence of the PNM EOS proposed in Ref. Prakash, Ainsworth, and Lattimer (1988). Our results cross the two different bands and, at the highest densities, coincide with the upper limit of the range obtained in Ref. Danielewicz, Lacey, and Lynch (2002). The explored range agrees with results from Ref. Tsang et al. (2018), which compares results from flow experiments Danielewicz, Lacey, and Lynch (2002) with the tidal deformability computed for the NS merger event GW170817 B. P. Abbott et al. (2017). We add to our comparisons the pressure of PNM obtained from chiral effective field theory (EFT) Tews et al. (2018). For densities up to $n\simeq 1.5n_{\rm sat}$, the values from the Skyrme EOSs are within the constraints of chiral EFT, although they are slightly above the limits for higher densities.

In the second column of Fig. 1, we plot the mass-radius relations of cold beta-equilibrated NSs obtained solving the TOV equations (top) and the baryonic mass of the cold NS with proton fraction $y$ above a critical value set to $y_{\rm crit}=0.11$ (bottom) as is the condition necessary for direct Urca processes to take place inside a NS Brown et al. (2018). Because we limit our analysis to EOSs that predict a large pressure at high densities, see Tab. 1, all EOSs satisfy the observational constraints for the mass of PSR J0348+0432, $2.01\pm 0.04\,M_{\odot}$ J. Antoniadis et al. (2013). A similarly large NS mass, $M=1.93\pm 0.02\,M_{\odot}$, has been observed for PSR J1614-2230 Fonseca et al. (2016). Furthermore, our choices of the other empirical parameters guarantee that the mass-radius relations are within the $1\sigma$ range of “model A” of Ref. Nättilä et al. (2016) obtained from observations of x-ray bursts. The EOSs also obey the constraints for the radius of a $1.4$-$M_{\odot}$ NS, $12.00\,\mathrm{km}<R_{1.4}<13.45\,\mathrm{km}$, computed from the data for the NS merger observation GW170817 Most et al. (2018). This constraint is more stringent than obtained by others for the same event, e. g., Ref. De et al. (2018) constrain radii of NSs to be in the $8.9\,\mathrm{km}<\bar{R}<13.2\,\mathrm{km}$ range while results from the LIGO and Virgo Collaborations suggest $R=11.9\pm 1.4\,\mathrm{km}$ B. P. Abbott et al. (2018). The constraint of Ref. De et al. (2018) was computed assuming hadronic EOSs for high density matter and from inference of the dimensionless tidal deformability deduced from the GW170817 event that suggests $\tilde{\Lambda}<800$ B. P. Abbott et al. (2017). Meanwhile, the LIGO/Virgo results require that both bodies that generated the GW170817 event are NSs described by the same EOS with spins within the range observed in Galactic binary NSs and are able to produce $1.97$-$M_{\odot}$ NSs. We notice only minor differences in the mass-radius relations as a function of the effective masses, mostly in the mass range $0.5\,M_{\odot}\lesssim M\lesssim 1.5\,M_{\odot}$.

Recently, it has been shown that the cooling rate of the NS in the transient system MXB 1659-29 while in quiescence is consistent with direct Urca reactions occurring in a small fraction of the core, $\approx 0.03\,M_{\odot}$ Brown et al. (2018). Assuming hadronic matter, this is only possible if nucleons in the core are unpaired and the proton fraction exceeds a critical value $y_{\rm crit}$ in the range $0.11-0.15$ Lattimer et al. (1991). Here we set $y_{\rm crit}=0.11$ and compute for each NS the total baryonic mass in the core which exceeds $y_{\rm crit}$, $M_{\rm baryon}(y>y_{\rm crit})$. We define $M_{\rm baryon}(y>y_{\rm crit})$ as the integrated baryonic mass in regions of the star where $y\geq y_{\rm crit}$ excluding the crust, i.e., the outer $\simeq 1\,\mathrm{km}$ of the star, as densities there are to low to induce direct Urca reactions. If the values chosen for the empirical parameters hold, the EOS described by those parameters implies that the NS in the MXB 1659-29 system has a mass in the range $1.6$ to $1.8\,M_{\odot}$ as lower mass values would imply that the proton fraction in the core never reaches the critical value $y_{\rm crit}$ to start the direct Urca process. Meanwhile, NSs with larger masses would cool at a much faster rate through direct Urca processes. Thus, under the assumption that matter in the core of a NS is made of unpaired nucleons, combined measurements of NS masses and cooling rates may be used to improve constraints on the EOS of dense matter.

Finally, in the last column of Fig. 1, we compare the interior properties of a canonical $1.4$-$M_{\odot}$ NS for the different EOSs. Although there are no clear visible changes for the density as a function of NS radius, we notice that there are, as in the case of the gravitational mass with proton fraction above $y_{\rm crit}$, small changes in the proton fraction in the core region as a function of the nucleon effective masses. These variations in proton fraction in the inner core of a $1.4$-$M_{\odot}$ NS are inversely (directly) correlated with $m^{\star}$ ($\Delta m^{\star}$). However, these changes are small, and the nucleon effective masses affect the central proton fraction $y_{1.4}$ at the center of a $1.4$-$M_{\odot}$ NS by at most 0.02. Nevertheless, a clear trend is observed here: EOSs that predict smaller radii for the same mass NS also predict a larger isospin asymmetry in their cores.

III.2 Symmetry energy and its slope

We now discuss the variation set $s_{S}$ [Eq. (11b)], where the symmetry energy $\epsilon_{\rm sym}$ and its logarithmic derivative with respect to density $L_{\rm sym}$ at saturation density are varied.

In Fig. 2, we plot the pressure as a function of density and properties of cold beta-equilibrated NSs considering changes in these quantities according to Tab. 1. Because we are only modifying parameters of the symmetry energy, the pressure of SNM remains unchanged across EOSs, see the top left panel in Fig. 2. Meanwhile, there are some variations in the pressure of PNM, as depicted by the bottom left panel of Fig. 2. The differences between the EOSs are largest below $\simeq 2n_{\rm sat}$ since the higher density behavior of the symmetry energy is strongly constrained by the fixed values of $P^{(4)}_{\rm SNM}$, $P^{(4)}_{\rm PNM}$, and $K_{\rm sym}$ for all EOSs in the variation set $s_{S}$ [Eq. (11b)]. Therefore, all EOSs obey the flow constraints from Danielewicz et al. Danielewicz, Lacey, and Lynch (2002) across a wide range of densities. In comparison, some of the $s_{S}$ EOSs become slightly inconsistent with the sub-saturation density chiral effective field theory constraints Tews et al. (2018) at low density.

The mass-radius curve of cold beta-equilibrated NSs, the top center plot in Fig. 2, is most impacted by symmetry energy variations at lower NS mass. For NSs with mass $M\lesssim 1.5\,M_{\odot}$, larger symmetry energies at saturation density $\epsilon_{\rm sym}$ and symmetry energy slopes $L_{\rm sym}$ result in larger NS radii. This is consistent with the results of Refs. Lattimer and Prakash (2001); Steiner et al. (2005), which highlight the impact of the density dependence of the symmetry energy on the NS radius. However, there are only minor changes in the mass-radius relationship in the region $M\gtrsim 2\,M_{\odot}$, as these NSs reach quite high densities in their cores where the pressure is fixed by $P^{(4)}_{\rm SNM}$ and $P^{(4)}_{\rm PNM}$. Nevertheless, massive NSs with approximately the same radius have very different inner compositions. See the bottom center panel of Fig. 2. For the variations considered here, we observe an inverse relationship between the NS radius and the amount of matter with proton fraction larger than the critical value $y_{\rm crit}=0.11$, i.e., the isospin asymmetry. This is also clearly seen in the composition of the $1.4$-$M_{\odot}$ NS, see bottom right plot in Fig. 2. At densities near or above $n_{\rm sat}$, the density profile of $1.4$-$M_{\odot}$ NSs is similar for all EOS parametrizations that differ only in $\epsilon_{\rm sym}$ and $L_{\rm sym}$, top right of Fig. 2. However, these NS radii may differ by up to $800\,\mathrm{m}$ due to different density profiles at densities lower than $n_{\rm sat}$.

III.3 Incompressibility

We now consider set $s_{K}$ [Eq. (11c)], where we analyze variations in the isoscalar incompressibility $K_{\rm sat}$, which is well constrained, and the isovector incompressibility $K_{\rm sym}$, which is poorly known (see Section II).

In Fig. 3, we plot the pressure of SNM (top left) and of PNM (bottom left). Small differences are evident in SNM for different $K_{\rm sat}$, while the differences in PNM are substantial due to the large range of values allowed for $K_{\rm sym}$. Since we keep the pressure of SNM and PNM at $n=4n_{\rm sat}$ fixed for all EOSs, the curves for the pressures cross at this value and at $n=n_{\rm sat}$. This limits the effect of both $K_{\rm sat}$ and $K_{\rm sym}$ at high density.

Variations in the incompressibilities cause drastic differences in the mass-radius relationships and compositions of cold NSs (see the center upper and center lower panels of Fig. 3, respectively). There is an inverse correlation between the radius of a NS predicted by a given EOS and its isospin asymmetry, which is similar to what we see for variation sets $s_{M}$ and $s_{S}$, Secs. III.1 and III.2, respectively. This is particularly obvious in the rightmost panels of Fig. 3, which show the internal properties of $1.4$-$M_{\odot}$ NSs.

We also observe different qualitative behaviors in the core composition that relate to the isovector incompressibility $K_{\rm sym}$. While for $K_{\rm sym}\lesssim-200\,\mathrm{MeV}$ the proton fraction in the NS core is almost constant, for $K_{\rm sym}\gtrsim-200\,\mathrm{MeV}$ the core asymmetry decreases with $K_{\rm sat}$. Similar properties are found across NSs with the same mass but different EOSs except for the most massive ones, $M\gtrsim 2\,M_{\odot}$.

III.4 Pressure at high-density

Finally, we consider the variation set $s_{P}$ [Eq. (11d)], where the pressures of SNM and PNM are varied at four times nuclear saturation density, while leaving all other empirical parameters constant. These variations begin to have an impact at densities $\gtrsim 2n_{\rm sat}$, since the saturation density properties of the EOSs are unaltered. This is clearly visible in the top and bottom left panels of Fig. 4.

Changes in the pressure at high densities translate directly into variations in the mass-radius relationship of high-mass NSs, which probe these high densities in their cores (see the center panels of Fig. 4). Since the pressure in the NS core is somewhere between the SNM and PNM pressures, increasing either one stiffens the EOS and gives rise to a larger radius for a fixed NS mass. Additionally, increasing the pressure of either SNM or PNM increases the predicted maximum NS mass.

Varying these pressures also impacts the predicted lepton richness of NSs. In the lower left panel of Fig. 4, it can be seen that the pressure of PNM is anti-correlated with the isospin asymmetry in the NS core while the pressure of SNM is correlated with the isospin asymmetry.

IV Spherically-symmetric core collapse

We now focus on how variations in the empirical parameters of the EOS and of the pressure of nuclear matter at high densities affects the core collapse of a massive star and its CCSN evolution. We will mainly investigate the impact of the EOS on neutrino emission during the postbounce phase. The details of neutrino emission from high-density matter in a CCSN is interesting both because these neutrinos can be directly detected from a galactic CCSN (e.g., Scholberg (2012)) and because these neutrinos can be re-absorbed in the lower density matter behind the CCSN shock and play a role in powering the explosion Wilson (1985). Uncertainties in the nuclear EOS translate into uncertainties in predictions of CCSN neutrino fluences, which in turn introduce uncertainty in the detectability of the neutrino emission and into the CCSN mechanism itself. Both the explosion mechanism and detectability are sensitive to changes in the neutrino energy spectra, which we will characterize by the root-mean-square (RMS) neutrino energy, $\sqrt{\langle\epsilon_{\nu}^{2}\rangle}$, and in the neutrino luminosities, $L_{\nu}$. Larger luminosities and RMS energies of electron neutrinos and antineutrinos result in higher predicted neutrino detection rates and more favorable conditions for explosion due to the quadratic energy dependence of neutrino interaction cross-sections.

Specifically, we study the collapse and bounce of a $20$-$M_{\odot}$ progenitor star (s20WH07 of Woosley and Heger (2007)) in spherical symmetry using the radiation-hydrodynamics code GR1D O’Connor and Ott (2011); O’Connor (2015). We study this progenitor star since it (1) has been studied by many other groups Bruenn et al. (2013); Char, Banik, and Bandyopadhyay (2015); Dolence, Burrows, and Zhang (2015); Melson, Janka, and Marek (2015); Bruenn et al. (2016); Pan et al. (2016); Suwa et al. (2016); Summa et al. (2016); O’Connor et al. (2017); Bandyopadhyay (2017); Ott et al. (2018); Just et al. (2018); Glas et al. (2019); O’Connor and Couch (2018); O’Connor et al. (2018), so comparisons can be readily made, (2) produces a massive PNS, (3) does not collapse into a black-hole within the first second after bounce, and (4) often exhibits the onset of an explosion in multi-dimensional simulations Bruenn et al. (2013); Melson, Janka, and Marek (2015); Summa et al. (2016); Ott et al. (2018); O’Connor and Couch (2018) soon after the density discontinuity from the Si/Si-O shell boundary crosses the shock radius. Furthermore, (5) the PNS central number density during the first second after bounce is in the range $2-3n_{\rm sat}$. Since we constrain our EOS with empirical properties at saturation density and at four-times saturation density, this maximum density does not go beyond the range of densities over which the EOSs have been fit.

For each EOS table described in Section II, a consistent set of neutrino opacities is generated using the NuLib library O’Connor (2015). We then run a core collapse simulation until $800$-$1000\,\mathrm{ms}$ after bounce. In the simulation, we consider electron neutrinos and electron antineutrinos separately and group the heavy flavored neutrinos and anti-neutrinos into a single composite species. For each species, we follow 24 logarithmically spaced neutrino energy groups running from 1 $\,\mathrm{MeV}$ to \simeq 269$$\,\mathrm{MeV}. The computational grid is set to have $1\,500$ grid cells, constant cell size of $100\,\mathrm{m}$ out to a radius of $20\,\mathrm{km}$, and then geometrically increasing cell size to an outer radius of $20\,000\,\mathrm{km}$. We map stellar mass rest-mass density $\rho$, proton fraction $y$, and pressure $P$ from the progenitor star to GR1D as described in Schneider, Roberts, and Ott (2017).

IV.1 Effective Mass

First, we consider the impact of variation set $s_{M}$ [Eq. (11a)] on core collapse, where the effective mass $m^{\star}$ and the effective mass splitting $\Delta m^{\star}$ are varied. Since the temperature enters only through the factor $m^{\star}T$ in the Skyrme model we use (see Eq. 1), one expects the finite temperature behavior of the EOS to be substantially impacted by changes in the effective mass. As shown in Sec. III, varying the effective mass in our EOS fitting procedure has a negligible impact on the zero-temperature EOS and therefore a negligible impact on cold-NS structure. On the other hand, in CCSNe, temperatures of tens of MeV can be reached and the finite-temperature properties of the EOS may have a substantial impact.

The high-density EOS impacts the neutrino emission by changing the structure and thermodynamic state of the region from which most neutrinos are emitted, the neutrinospheres. The position of the neutrinosphere depends on both the energy and neutrino species (flavor, neutrino/antineutrino). Here, we consider the properties of a neutrino-energy averaged neutrinosphere, which qualitatively captures the state of the material from which the bulk of the neutrinos are emitted. The neutrinosphere is defined as the location where the opacity is equal to $\tau_{\nu}=2/3$ O’Connor and Ott (2013). Generally, before explosion, the neutrinospheres move to smaller radius, higher density, and higher temperature over time. The electron neutrinospheres and antineutrinospheres also stay nearly in neutrino-free beta-equilibrium since they can efficiently lose lepton number by definition.

In Fig. 5, the influence of varying the effective mass on the neutrinosphere properties is shown. Increasing the SNM effective mass at saturation density, $m^{\star}$, increases the temperature of the neutrinosphere for all flavors and at all times. On the other hand, increasing $m^{\star}$ decreases the neutrinosphere radii for all flavors. For electron neutrinos and antineutrinos, increasing $m^{\star}$ causes a decrease in the neutrinosphere density. Higher temperatures result in larger values of the beta-equilibrium $y_{e}$. For heavy-lepton neutrinos, increasing $m^{\star}$ increases the neutrinosphere density slightly. The impact of variations in $\Delta m^{\star}$ on the neutrinosphere properties is relatively small.

It is somewhat counter intuitive that the neutrinosphere temperatures increase with the effective mass, since the effective mass enters the EOS in the combination $m^{\star}T$. Nevertheless, it is easy to understand this behavior. First, note that deviations in the nucleon effective masses from their vacuum mass depends linearly on the density. Since the density of the neutrinospheres is less than a hundredth nuclear saturation density this means that $m^{\star}$ at the neutrinosphere is essentially the bare nucleon mass. Therefore, the impact of varying the effective mass on the neutrinosphere properties must be indirect. For small temperatures where the Sommerfeld expansion is valid, the entropy in nucleon species $t$ is given by, see Appendix A,

[TABLE]

In the same approximation, the temperature dependent contribution from species $t$ to the nucleon pressure is given by $P_{\text{th},t}=Tn_{t}s_{t}/3$. Therefore, in so far as the density and entropy throughout the outer layers of the PNS are not impacted by changes in the effective mass, the pressure of material below the neutrinosphere goes down with increasing effective mass since $T\propto(m^{\star})^{-1}$. This suggests that increasing the effective mass results in more compact outer layers of the PNS. This is consistent with what our simulations show. As can be seen in Fig. 5, increasing the effective mass results in a smaller radius neutrinosphere which, in turn, results in a larger virial temperature for the neutrinosphere.

Variations in the neutrinosphere properties are directly imprinted in the CCSN neutrino emission itself. In Fig. 6 we plot the RMS energy (top) and luminosity (bottom) of the three neutrino species considered, i.e., $\nu_{e}$, $\bar{\nu}_{e}$, and $\nu_{x}=\nu_{\mu/\tau}=\bar{\nu}_{\mu/\tau}$. Soon after core bounce, $t-t_{\rm bounce}\lesssim 200\,\mathrm{ms}$, all EOSs predict RMS energy and luminosity of neutrinos emitted that differ only by $\lesssim 5\%$ in the most extreme cases. However, after the first $\simeq 200\,\mathrm{ms}$, neutrino energies and luminosities start to diverge. The average RMS energy of all neutrino flavors and the luminosity for $\nu_{e}$ and $\bar{\nu}_{e}$ neutrinos is higher the larger the nucleon effective mass $m^{\star}$ at saturation density is. Meanwhile, there is barely any change in the neutrino luminosity for the heavy-lepton neutrinos $\nu_{x}$ as the effective mass changes. Moreover, differences in neutrino properties are only affected at the $\simeq 1\%$ level by the change of the nucleon effective mass difference $\Delta m^{\star}$.

In Fig. 6 we see that the largest variation in the RMS energies occurs for the heavy-lepton neutrinos $\nu_{x}$ after $\simeq 400\,\mathrm{ms}$ after core bounce, although the heavy lepton neutrino luminosities are barely affected. Nevertheless, supernova electron neutrinos and antineutrinos have a larger impact on the supernova explosion mechanism with the latter being easier to detect Mirizzi et al. (2016). We observe that an increase in the effective mass $m^{\star}$ also leads to an increase in the RMS electron neutrino and antineutrino energies by about $2$ to $3\,\mathrm{MeV}$ soon after core bounce, $t-t_{\rm bounce}\gtrsim 200\,\mathrm{ms}$, while luminosities increase by up to $30\%$. An interesting question is whether different neutrino interactions, e.g., inelastic neutrino-nucleus scattering, will result in the same qualitative and quantitative differences. We postpone an investigation of this question to future work.

An increasing effective mass increases the luminosity and average energy of electron neutrinos and antineutrinos and thereby increases the rate of neutrino heating behind the SN shock. Therefore, it might be expected that a higher effective mass makes conditions more favorable for shock runaway. Nevertheless, we find that larger effective masses result in smaller shock radii in spherically-symmetric runs. In Fig. 7 we observe that the shock radius $R_{\rm shock}$ follows the PNS radius $R_{12}$. In these spherically-symmetric simulations, the impact of the reduced PNS radius on the shock overwhelms the increased neutrino heating rate when the effective mass is increased. Nevertheless, in multi-dimensional simulations, the larger neutrino luminosities and average energies may instead lead to shock radii that expand faster for larger nucleon effective masses $m^{\star}$. This is discussed in Sec. V.

Besides neutrinos emitted during core collapse, we also discuss the hot PNS evolution during the first second after collapse, see Fig. 7. In Ref. Summa et al. (2016), it is argued in the context of 2D simulations, that the LS220 EOS leads to fast contracting PNSs because this EOS generates compact cold beta equilibrated NSs. In our simulations we see that the collapse of a massive star simulated using EOSs that differ only in their effective masses predict very similar mass-radius relations for cold NSs, see Fig. 1. Although all these EOSs produce very similar cold beta equilibrated NSs, they predict distinct behaviors for the PNSs formed in core collapse. In Fig. 7, we plot the core temperature $T_{c}$ and density $\rho_{c}$ as well as shock radius $R_{\rm shock}$ and core radius $R_{12}$, the latter defined as the radius where mass density is $\rho=10^{12}\,\mathrm{g\,cm}^{-3}$. There is a clear correlation between the effective mass $m^{\star}$ and the core density after bounce as well as how fast the PNS radius and shock contract after reaching their maximum values. The core temperature, on the other hand, is higher (lower) the lower (higher) $m^{\star}$ is.

Density, temperature, and proton fraction profiles of the PNS at $500\,\mathrm{ms}$ after bounce are plotted in Fig. 8. It is clear that EOSs with higher $m^{\star}$ produce less thermal support since their temperatures are lower in most of the PNS interior and hot mantle, although the temperature is higher in the region where it peaks. Thus, we deduce the reason the LS220 EOS leads to faster contraction when compared to other EOSs is better explained by its assumptions about its effective mass, set by $m^{\star}=m_{n}$, rather than by the mass-radius relation it predicts for cold beta equilibrated NSs, which is barely affected by the effective mass.

Fig. 4 of Ref. Steiner, Hempel, and Fischer (2013) shows that the PNS radius that follows from the core collapse of a $11.2$-$M_{\odot}$ progenitor star simulated with the LS220 EOS contracts significantly faster than the radius of PNSs simulated with other EOSs that have $m^{\star}/m_{n}\simeq 0.61-0.76$. However, the EOSs in that work use diverse prescriptions to compute the EOSs at low and high densities, which makes a direct comparison between our results and their results non-trivial. In this work, by unifying the formalism used for all EOSs, we are able to draw stronger conclusions about the effect of each parameter of the EOS on the core collapse, and specifically on the role of the effective mass.

IV.2 Symmetry energy and its slope

We perform core collapse simulations using variation set $s_{S}$, where the symmetry energy and its slope are varied. We observe that for the range of variations considered for $\epsilon_{\rm sym}$ and $L_{\rm sym}$, the changes in the neutrino spectra and the PNS properties are rather small. They are of comparable in magnitude to the changes seen from varying the nucleon effective mass splitting, $\Delta m^{\star}$. Thus, for the purpose of simulations of CCSNe, these two quantities are rather well constrained and we expect that even substantial variations around the current best estimates for these two observables will not affect simulation results significantly.

It may be the case, however, that if we were to simulate these CCSNe for longer timescales, including into the cooling phase, that larger differences between EOSs could become apparent. We defer this, as well as CCSN simulations of different progenitors, to future work.

IV.3 Incompressibility

We now discuss effects in CCSN simulations due to changes in the incompressibility parameters $K_{\rm sat}$ and $K_{\rm sym}$ of the EOS. As in the case of variations in the symmetry energy $\epsilon_{\rm sym}$ and its slope $L_{\rm sym}$, the relative changes in the neutrino spectra are rather small and at most twice those observed for changes in the nucleon effective mass splitting, $\Delta m^{\star}$. Nevertheless, it is likely that this is the case only for the short times we evolved the collapsing star, $t\lesssim 1.0\,\mathrm{s}$. For longer evolutions or more massive progenitors, larger differences between the EOSs are likely. This may be inferred from Figs. 9 and 10. The former shows the evolution of the core density, core temperature, shock radius, and PNS radius. The latter shows the PNS density, temperature, and proton fraction as a function of radius. The central density and temperature of the PNS at $\simeq 1\,\mathrm{s}$ after bounce differ by $\simeq 20\%$ between the most extreme cases. Meanwhile, changes in the shock radius and PNS radius, are affected only in the $\simeq 5\%$ range. At that time, the maximum PNS mantle temperature is correlated with both $K_{\rm sat}$ and $K_{\rm sym}$. On the other hand, the PNS and shock radius are anti-correlated with these quantities. We observe that despite the much larger error bar in $K_{\rm sym}$ when compared to $K_{\rm sat}$, both lead to uncertainties in PNS structure of similar magnitudes. We expect these differences to be amplified in multi-dimensional simulations due to the interplay between neutrino heating and hydrodynamic instabilities that can lead to shock revival Müller (2016). Hence, it is important for realistic simulations that these two parameters are constrained further in the future.

IV.4 Pressure at high-density

We also study the differences in the neutrinos spectra and in the PNS evolution during the first second of collapse for the $20$-$M_{\odot}$ progenitor star due to changes in the pressure of SNM and PNM at $n=4n_{\rm sat}$, set $s_{P}$ in Eq. (11d). As expected, by the end of our runs neither the emitted neutrinos nor the PNS properties were significantly altered by changes in the pressure at high densities. Except for changes of $\simeq 5\%$ with respect to the baseline EOS for the density and temperature in the core near the end of the runs, none of the other quantities studied (neutrino luminosity and RMS energy, and shock and PNS radii) differed by more than 1% during the run. This is due to the maximum density in the PNS still being below $2.5n_{\rm sat}$ at $t-t_{\rm{bounce}}\simeq 1\,\mathrm{s}$ and, thus, the EOSs used in all runs did not reach regions were the differences become large. Lower pressures at high densities cause densities (temperatures) in the core to increase faster (slower). As in the cases of changes in the symmetry energy $\epsilon_{\rm sym}$ and its slope $L_{\rm sym}$, we expect that longer evolutions will show differences for the different EOSs, as the densities reached throughout the PNS will be higher. Furthermore, we expect the pressure at high densities to play a significant role in setting the time of collapse of the PNS to a BH. Such a study is currently underway Schneider et al. .

V Three-dimensional CCSN Simulations

In order to further investigate the insights gained from performing spherically-symmetric (1D) core collapse simulations discussed in Sec. IV for different EOSs, we perform six three-dimensional (3D) octant runs, i.e., limited to one octant of the 3D cube, for the same non-rotating $20$-$M_{\odot}$ presupernova model s20WH07 Woosley and Heger (2007). In our 1D simulations, we find that increasing the nucleon effective mass makes the PNS atmosphere more compact and increases the neutrino energies and luminosities. In the spherically-symmetric simulations, the impact of a reduced PNS radius overwhelmed the impact of increased neutrino heating. Hence, larger effective masses result in smaller maximum shock radii. Nevertheless, spherical symmetry inhibits hydrodynamic instabilities that may be present behind the shock and these conclusions may not hold in more realistic three-dimensional simulations. Five of the 3D runs are performed using variants of the finite temperature SLy4 EOS Chabanat et al. (1998); Schneider, Roberts, and Ott (2017). Additionally, we perform one run with the often used Lattimer & Swesty EOS with $K_{\rm sat}=220\,\mathrm{MeV}$, LS220. The SLy4 and LS220 EOS properties at $T=0$ are listed in Tab. 2.

The variants of the SLy4 EOS are computed using the methods described in Sec. II and Appendix B by keeping all empirical quantities except the effective mass for SNM at saturation density $m^{\star}$ constant. The values used for the effective mass are $m^{\star}/m_{n}=0.6$, $0.7$, $0.8$, $0.9$, and $1.0$. In the discussion that follows we differentiate between the different SLy4 EOSs by adding a subscript that corresponds to the effective mass used, SLy4${}_{m^{\star}/m_{n}}$. As in Sec. IV, the SLy4 EOSs as well as the LS220 EOS are connected to a low-density EOS of 3 335 nuclei in NSE using the prescription outlined in Sec. II.4.

Following Sec. IV, we simulate the collapse of the progenitor star using the GR1D code O’Connor and Ott (2013); O’Connor (2015). In this phase, the neutrino reactions are considered in the exact same manner as discussed in the previous section. Following Ref. Roberts et al. (2016), we map the spherically-symmetric collapsing progenitor 20$\,\mathrm{ms}$ after core bounce to a high-resolution octant 3D geometry with reflecting boundary conditions on the $xy$, $yz$, and $zx$ planes. The remainder of the simulation is performed using the general-relativistic radiation-hydrodynamics code Zelmani Roberts et al. (2016), which is itself based on the Einstein Toolkit Löffler et al. (2012); Mösta et al. (2014). At this point we modify the neutrino transport and consider only 16 energy groups. As in Ref. Ott et al. (2018), we employ the subset of neutrino opacities from Ref. Bruenn (1985), but now leave out velocity dependence and inelastic neutrino-electron scattering.

In Fig. 11, we plot the neutrino RMS energies, $\sqrt{\langle\epsilon_{\nu}^{2}\rangle}$, and luminosities, $L_{\nu}$, after core bounce for the three considered neutrino species. As in the spherically-symmetric case, both neutrino energies and luminosities, for the Skyrme-type EOSs, increase as the effective mass is increased. In the range of effective masses studied, differences in neutrino RMS energies are approximately $1.5\,\mathrm{MeV}$ for all neutrino species. We observe that neutrino energies and luminosities, especially for the heavy-lepton neutrinos $\nu_{x}$, computed for the LS220 EOS are higher than for the SLy41.0 EOS, even though both have the same effective mass for SNM at saturation density, $m^{\star}=m_{n}$. The reason for this is that most of the empirical parameters that differ between the two EOSs, see Tab. 2, shift neutrino luminosities and energies to higher values for the LS220 EOS with respect to the SLy41.0 EOS. The exception is $K_{\rm sym}$, which slightly decreases the neutrino output for the LS220 when compared to the SLy41.0 EOS. The pressure at high densities, represented by $P^{(4)}_{\rm SNM}$ and $P^{(4)}_{\rm PNM}$, meanwhile, does not have a significant effect for this progenitor within the first second of core bounce.

Although variations of the effective mass have a similar impact on the RMS neutrino energies in the 3D simulations as they had in the 1D simulations, the resulting shock radius evolutions differ substantially. The s20WH07 progenitor has a steep density and specific entropy discontinuity at the Si/Si-O shell interface. In the full 3D simulations of Ref. Ott et al. (2018) for the same progenitor star but using the SFHo EOS Steiner, Hempel, and Fischer (2013), the abrupt decrease in the ram pressure at the shock as the discontinuity is accreted results in shock runaway. In Fig. 12, we plot the shock radius and accretion rate for our six octant 3D simulations. The accretion rates for all octant runs agree within 1% or less, while the shock radius after the shock crosses the Si/Si-O is very sensitive to the EOS, with only the LS220 EOS predicting shock runaway.

In this paper, we choose not to carry out a direct comparison between our results and that of Ref. Ott et al. (2018). We do so for a number of reasons. First, full 3D runs appear to more readily lead to shock runaway than octant runs Roberts et al. (2016). Second, when setting the initial conditions of the run we choose to preserve density $\rho$, proton fraction $y$, and pressure $P$, while in Ref. Ott et al. (2018) chose density $\rho$, proton fraction $y$, and temperature $T$. This leads to different times of core bounce and a different accretion history. Finally, the SFHo EOS, including its low-density part, is generated using a relativistic mean-field approach and not a Skyrme model. Fig. 15 of Ref. Schneider, Roberts, and Ott (2017) shows how changes in the low density EOS affect the postbounce accretion rate. Understanding how the difference in the low density EOS as well as in the initial conditions lead to differences in the PNS profile and CCSNe evolution is beyond the scope of the present work.

With respect to the shock radius evolutions resulting from the different EOSs, we note that for the octant runs EOSs with higher effective masses for SNM at saturation density $m^{\star}$ generally lead to larger shock radius after bounce. In the LS220 run, the shock runs away approximately $350\,\mathrm{ms}$ after core bounce reaching, on average, $500\,\mathrm{km}$ by the end of the run. In the SLy41.0 run, on the other hand, the average shock radius grows up to $220\,\mathrm{km}$ at $320\,\mathrm{ms}$ after core bounce, only slightly lower than what is predicted for the LS220 EOS, but then recedes. Although this is opposite to the pattern seen for the shock radii in 1D runs, see Fig. 7, this is expected in 3D simulations due to the higher neutrino luminosities and RMS energies for EOSs that have higher $m^{\star}$. Compare Fig. 6 for 1D runs and Fig. 11 for the 3D octant runs. An exception is the SLy40.6 EOS, whose 3D simulation predicts shock radius behavior similar to the SLy40.8 run and higher radii than what we observe in the SLy40.7 run, despite its lower neutrino luminosities and average energies. This is likely due counteracting effects of lower neutrino production, but larger initial mass in the gain region for EOSs with lower effective masses, see Fig. 13.

In Fig. 13, we present diagnostics that help us understand variations in the results for the different EOSs. First, higher neutrino energies and luminosities lead to higher integrated neutrino heating, heating minus cooling $\dot{Q}$, and higher heating efficiency, $\eta=\dot{Q}(L_{\nu_{e}}+L_{\bar{\nu}_{e}})^{-1}$, in the gain layer. Ref. Ott et al. (2018) showed that for the first $80-100\,\mathrm{ms}$ after bounce, the heating efficiency $\eta$ is almost independent of the progenitor. Here we observe that $\eta$ is also almost completely EOS independent early after bounce. However, it is clearly correlated with the effective mass $m^{\star}$ at later postbounce times. At the time when the Si/Si-O interface reaches the shock, $\eta$ is $\simeq 50\%$ higher for EOSs with $m^{\star}=m_{n}$ compared to the ones with $m^{\star}\geq 0.6m_{n}$.

Next, from Fig. 13, we see that the mass in the gain layer $M_{\rm gain}$ is mostly EOS independent until the Si/Si-O shell crosses the shock radius. After this occurs, EOSs that predict higher PNS compactness, $(GM_{\rm PNS})/(R_{\rm PNS}c^{2})$, also predict larger mass in the gain layer, another indicator of favorable conditions for shock runaway. The ratio between the timescales $\tau_{\rm adv}\simeq M_{\rm gain}\dot{M}^{-1}$ for material to advect through the gain layer and $\tau_{\rm heat}\simeq|{E_{\rm gain}}|{\dot{Q}^{-1}}$ for neutrino heating is another such indicator Janka, H.-Th. (2001); Thompson, Quataert, and Burrows (2005). Following implementation details of Ref. Müller, Janka, and Marek (2012), we find that two of the EOSs, LS220 and SLy41.0, cross the $\tau_{\rm adv}/\tau_{\rm heat}\gtrsim 1$ threshold set as a condition that favors shock runaway, while SLy40.9 comes very close to it. While the LS220 EOS results in shock runaway, none of the simulations using variants of the SLy4 EOS lead to shock runaway within $400\,\mathrm{ms}$ of core bounce. Not even the SLy41.0 EOS, despite reaching a ratio between advection and heating timescales $\tau_{\rm adv}/\tau_{\rm heat}\gtrsim 1.5$. As discussed in Ref. Ott et al. (2018), $\tau_{\rm adv}/\tau_{\rm heat}$ serves more as a diagnostic of shock runaway than a condition for explosion. Even at times where $\tau_{\rm adv}/\tau_{\rm heat}\gtrsim 1$ for the simulations employing the SLy41.0 EOS, the mass in the gain layer continues to decrease and the shock stabilizes at $\langle{R_{\rm shock}}\rangle\simeq 200\,\mathrm{km}$ before receding. In the simulation using the LS220 EOS, the mass in the gain layer stabilizes and then grows once explosion sets in.

Finally, we also plot in Fig. 13 the average radial and angular turbulent energies as defined in Ref. Müller et al. (2017). As argued in Refs. Murphy, Dolence, and Burrows (2013); Couch and Ott (2015); Radice et al. (2016); Roberts et al. (2016) we find that the total turbulent energy is anisotropic on large scales, i.e., $\langle\epsilon_{\rm{turb}}\rangle_{r}\simeq\langle\epsilon_{\rm{turb}}\rangle_{\theta+\phi}$. Furthermore, $\langle\epsilon_{\rm{turb}}\rangle$ is mostly EOS independent until $\simeq 260\,\mathrm{ms}$ after core bounce, when shock behavior becomes very sensitive to the EOS. For simulations using the SLy4 EOSs with $m^{\star}\leq 0.9m_{n}$ shock radius recedes quickly in the late stages of the run while $\langle\epsilon_{\rm turb}\rangle$ increases. On the other hand, using the LS220 EOS leads to saturation of $\langle\epsilon_{\rm turb}\rangle$ as its shock runs away, Fig. 12. The SLy41.0 EOS predicts a behavior for the turbulent energy density that is a mix of the predictions by the simulations using the LS220 and the other SLy4 EOSs: a momentary stabilization of $\langle\epsilon_{\rm turb}\rangle$ is achieved while $\langle{R_{\rm shock}}\rangle\simeq 200\,\mathrm{km}$ followed by a fast rise as the shock radius recedes.

VI Conclusions

We carried out a detailed study of the impact of variations of different experimentally accessible parameters of the nuclear matter EOS on the properties of cold beta-equilibrated neutron stars (NSs) and on the core collapse and postbounce evolution of a massive star.

Using the SROEOS code Schneider, Roberts, and Ott (2017), we constructed 97 finite-temperature EOSs in which we systematically varied the empirical parameters of the EOS based on the experimental and theoretical constraints compiled in Refs. Margueron, Hoffmann Casali, and Gulminelli (2018a, b); Danielewicz, Lacey, and Lynch (2002). We then used these EOSs to compute the properties of cold beta-equilibrated NSs and to simulate the core collapse of the $20$-$M_{\odot}$ presupernova stellar model of Ref. Woosley and Heger (2007). We carried out core-collapse supernova (CCSN) simulations using the spherically-symmetric general-relativistic radiation-hydrodynamics code GR1D O’Connor and Ott (2011, 2010); O’Connor (2015). We carried out the simulations to approximately and investigated the neutrino signals and protoneutron star (PNS) evolution for each EOS.

Although the uncertainty in the effective nucleon mass at saturation density has a negligible impact on the properties of cold NSs in our EOS model, we find variations in the effective mass have a substantial impact on the postbounce evolution of our CCSN models. The effective nucleon mass mainly regulates the temperature dependence of the Skyrme-type EOSs we consider, so it impacts the structure of the shock heated material in the PNS. Specifically, we found that the effective mass of nucleons in SNM at saturation density, $m^{\star}$, impacts the interior structure of the PNS, the PNS radius, the CCSN neutrino emission, and the evolution of the CCSN shock. Increasing the effective mass increases the average neutrino energies for all neutrino types and their total luminosity. This is because increasing the effective mass $m^{\star}$ leads to more compact PNSs with hotter neutrinospheres, although the larger effective masses result in lower PNS core temperatures. Recently, similar conclusions regarding the impact of the effective mass were reported from spherical-symmetric simulations of a $15$-$M_{\odot}$ progenitor star Yasin et al. (2018).

Variations in other parameters of the EOS, such as changes in the neutron-proton effective mass splitting in PNM, have a small impact on CCSN evolution. Moreover, changes in the isoscalar part of the incompressibility, $K_{\rm sat}$, affects temperature and density in the core of PNSs, but has limited impact on the neutrino signal, and the outer regions of the PNS. Although it is more weakly experimentally constrained, varying the isospin incompressibility, $K_{\rm sym}$, leads to variations in neutrinos signal and PNS evolution of the same order of magnitude as the isoscalar incompressibility, $K_{\rm sat}$. Furthermore, for the purpose of CCSN evolution, symmetry energy terms and the pressure at high densities, $n\gtrsim 4n_{\rm sat}$, have even smaller impact on the outcome of the core collapse than changes in the incompressibility. Based on the spherically-symmetric simulation results, we conclude that most of the uncertainty introduced into simulations of core collapse evolution and its neutrino signal by uncertainties in the EOS is due to the temperature dependence of the EOS and, to a lesser degree, due to the nuclear incompressibility.

To confirm these spherically-symmetric results, we performed six octant 3D simulations using the LS220 EOS and five variants of the SLy4 EOS where the effective mass of nucleons for SNM at saturation density was varied in the $m^{\star}=0.6-1.0\,m_{n}$ range. The runs were performed using the same set-up as the spherically-symmetric runs up to $20\,\mathrm{ms}$ after bounce and using the Zelmani code Roberts et al. (2016) leaving out the velocity dependence and inelastic neutrino-electron scattering in the neutrino transport.

Among the octant runs, lower $m^{\star}$ causes lower neutrino average energies and luminosities, as was the case in the spherically-symmetric runs. The lower neutrino energies result in less neutrino heating of the gain layer which subsequently leads to lower shock radii and failed explosions. Only the simulation using the LS220 EOS ($m^{\star}=m_{n}$) shows shock runaway at $\sim 350\,\mathrm{ms}$ after core bounce. For the SLy4 EOS variants there is a strong correlation between the shock radii and the value of $m^{\star}$. For runs employing the SLy4 EOS variant with $m^{\star}/m_{n}=1.0$, SLy41.0, and $0.9$, SLy40.9, the average shock radius reaches $\simeq 220\,\mathrm{km}$ and $180\,\mathrm{km}$, respectively, before starting to recede. For the other SLy4 EOS variants, the maximum average shock radius is limited to $160\,\mathrm{km}$. Analysis of our simulations shows that the run using the SLy41.0 EOS reached conditions very close to those that induce shock runaway. Specifically, the ratio between the advection and heating time scales is well above the limit $\tau_{\rm adv}/\tau_{\rm heat}\gtrsim 1$, usually indicative of impending shock runaway. It is likely that the small differences in nuclear saturation density properties between SLy41.0 and LS220, which play only a secondary role in our spherically-symmetric runs, determine that the shock runs away in the latter simulation while it does not in the former. We expect full 3D simulations to more easility lead to shock runaway than the octant simulations considered here Roberts et al. (2016). Thus, it is likely that for such conditions, the SLy41.0, and maybe even some of the other SLy4 EOS variants with lower $m^{\star}$, may experience shock runaway in full 3D.

Our octant runs may be compared to the full 3D run of Ott et al. for the same progenitor Ott et al. (2018). That run used the SFHo EOS Steiner, Hempel, and Fischer (2013), which has $m^{\star}=0.76m_{n}$. Nevertheless, despite the relatively low value of $m^{\star}$, that simulation saw shock runaway. It is likely that full 3D, differences in the high and low-density EOS, and differences in the setup of the initial conditions all played a role in the outcome of that simulation. This highlights the difficulty of comparing the role of the EOS between simulations that differ in many ways.

Understanding the effects each element of the EOS has on the outcome of a core collapse event is a long standing problem in nuclear and computational astrophysics. Using the SROEOS code Schneider, Roberts, and Ott (2017) we have, for the first time, determined in a consistent manner the pieces of the EOS that most significantly affect core collapse dynamics and PNS evolution. We demonstrated that uncertainties in the temperature dependence of the EOS affect neutrino energies and luminosities and play an important role in determining whether shock runaway takes place. We stress the need to extend our study to understand the EOS effects with different progenitors, full 3D simulations, and using other CCSNe simulation codes O’Connor et al. (2018) to confirm our findings.

Acknowledgements.

We acknowledge helpful discussions with H. Nagakura, I. Tews, C. Constantinou, M. Prakash, C. J. Horowitz, S. Couch, and MK.L. Warren. This research was funded by the National Science Foundation under award No. AST-1333520, CAREER PHY-1151197, PHY-1404569, OAC-1550514, and by the Sherman Fairchild Foundation. This research used resources of the National Energy Research Scientific Computing Center (NERSC), a U.S. Department of Energy Office of Science User Facility operated under Contract No. DE-AC02-05CH11231. The authors acknowledge the Texas Advanced Computing Center (TACC) at The University of Texas at Austin for providing HPC resources that have contributed to the research results reported within this paper.

Appendix A Sommerfeld Expansion

To compute the Sommerfeld expansion we make use of

[TABLE]

Some algebra leads to

[TABLE]

where we defined $\kappa_{t}=(1/2\pi^{2})(2m_{t}^{\star}/\hbar^{2})^{3/2}$, $\mu_{tF}=(\hbar^{2}/2m_{t}^{\star})(3\pi^{2}n_{t})^{2/3}$, and $\tilde{\mu}_{t}=T\eta_{t}$. We may invert Eq. (14) to obtain

[TABLE]

where $\mu_{tF}=\hbar^{2}k_{tF}^{2}/2m_{t}^{\star}$ is the Fermi chemical potential with $k_{tF}=(3\pi^{2}n_{t})^{1/3}$ the Fermi momentum.

A similar procedure implies that the kinetic energy density is

[TABLE]

with $\tau_{tF}=\frac{3}{5}k_{tF}^{2}n_{t}$. Thus, the low temperature limit of the specific entropy

[TABLE]

where $s_{t}$ is given in Eq. (12). The thermal contribution to the pressure is readily obtained from

[TABLE]

which reduces to $P_{th}\simeq\tfrac{1}{3}Tns_{B}$. Expressions containing higher order terms can be found in Ref. Prakash et al. (1997).

Appendix B Linear Equations

We present the linear equations discussed in Sec. II used to obtain the Skyrme parametrization given the set of EOS properties in Tab. 1.

The $\alpha_{1}$ and $\alpha_{2}$ parameters are computed from the properties of the effective masses $m^{\star}(n,y)$ at two distinct points in the $n$, $y$ phase space. We set the neutron effective mass value at $m_{n}^{\star}(n_{\rm sat},1/2)$ and $\Delta m^{\star}(n_{\rm sat},0)=m_{n}^{\star}(n_{\rm sat},0)-m_{p}^{\star}(n_{\rm sat},0)$ and compute the $\alpha$ parameters from the coupled equations:

[TABLE]

Here $\beta_{t}=\hbar^{2}/2m_{t}$ and $\beta_{\Delta}=\hbar^{2}/2(m_{n}^{\star}-m_{p}^{\star})$. Eq. (20) reduces to $(\alpha_{2}-\alpha_{1})n_{\rm sat}=(\beta_{n}-\beta_{p})$ when $\Delta m^{\star}(n_{\rm sat},0)=0$. We decided to compute the parameters $\alpha$ that set the effective mass of nucleons separately from the other Skyrme parameters to avoid negative effective masses at high densities and/or large isospin asymmetries.

The parameters $a_{i}$ and $b_{i}$ in Eq. (7) are computed by solving the system of linear equations $Ax=B$ where

[TABLE]

where we defined

[TABLE]

$x=(a_{0},b_{0},a_{1},b_{1},a_{2},b_{2},a_{3},b_{3})^{T}$, and

[TABLE]

In Eq. (26) $n_{\rm sat}$, $\epsilon_{\rm sat}$, $\epsilon_{\rm sym}$, $K_{\rm sat}$, $K_{\rm sym}$, $L_{\rm sym}$, $P^{(4)}_{{\rm SNM}}$, and $P^{(4)}_{{\rm PNM}}$, are, respectively, the nuclear saturation density, energy at saturation, symmetry energy at nuclear saturation density, isoscalar incompressibility, isovector incompressibility, the slope of the symmetry energy, and the pressures of SNM and PNM at $4n_{\rm sat}$. Furthermore, $\epsilon_{\rm kin}(n,y)$ is the kinetic energy term of the specific energy and was defined in Eq. (2) while

[TABLE]