Two-Dimensional Core-Collapse Supernova Simulations With the Isotropic Diffusion Source Approximation for Neutrino Transport

Kuo-Chuan Pan ({CJK*}UTF8bsmi潘國全), Matthias Liebendörfer, Matthias Hempel and Friedrich-Karl Thielemann Physik Department, Universität Basel, Klingelbergstrasse 82, CH-4056 Basel, Switzerland;

The neutrino mechanism of core-collapse supernova is investigated via non-relativistic, two-dimensional (2D), neutrino radiation-hydrodynamic simulations. For the transport of electron flavor neutrinos, we use the interaction rates defined by Bruenn (1985) and the isotropic diffusion source approximation (IDSA) scheme, which decomposes the transported particles into trapped particle and streaming particle components. Heavy neutrinos are described by a leakage scheme. Unlike the “ray-by-ray” approach in some other multi-dimensional supernova models, we use cylindrical coordinates and solve the trapped particle component in multiple dimensions, improving the proto-neutron star resolution and the neutrino transport in angular and temporal directions. We provide an IDSA verification by performing 1D and 2D simulations with 15 and 20 progenitors from Woosley et al. (2007) and discuss the difference of our IDSA results with those existing in the literature. Additionally, we perform Newtonian 1D and 2D simulations from prebounce core collapse to several hundred milliseconds postbounce with 11, 15, 21, and 27 progenitors from Woosley et al. (2002) with the HS(DD2) equation of state. General relativistic effects are neglected. We obtain robust explosions with diagnostic energies  B for all considered 2D models within approximately milliseconds after bounce and find that explosions are mostly dominated by the neutrino-driven convection, although standing accretion shock instabilities are observed as well. We also find that the level of electron deleptonization during collapse dramatically affect the postbounce evolution, e.g. the ignorance of neutrino-electron scattering during collapse will lead to a stronger explosion.

Subject headings:
hydrodynamics — instabilities —neutrinos — supernovae: general
slugcomment: Accepted for publication in ApJ

1. Introduction

After nearly five decades from the first attempt to obtain a neutrino-driven explosion by Colgate & White (1966), the explosion mechanism of core-collapse supernovae (CCSNe) remains elusive. It is generally believed that neutrino transport and convection are important ingredients to achieve a successful explosion (see recent reviews in Janka 2012; Burrows 2013; Foglizzo et al. 2015). However, modeling CCSNe with Boltzmann transport in three dimensions is numerically expensive and too time consuming with current computing resources, since the neutrino radiation hydrodynamics with Boltzmann transport is essentially a seven-dimensional problem: three spatial coordinates, two angular degrees of freedom, neutrino energy and time. In addition, there are three types of neutrino species and their anti-particles that would require a solution of the Boltzmann equation.

CCSN simulations with Boltzmann transport have been studied in 1D (Mezzacappa & Bruenn, 1993; Liebendörfer et al., 2001, 2004; Sumiyoshi et al., 2005), in 2D (Livne et al., 2004; Ott et al., 2008; Brandt et al., 2011) and recently in 3D with low resolutions and fixed background profiles (Sumiyoshi et al., 2015). However, these 2D and 3D works have ignored the velocity-dependent terms and decoupled these from the energy groups for simplicity, and the spatial resolutions are not sufficient to achieve a correct turbulence cascade. Therefore, approximate methods for the neutrino transport in multi-dimensional simulations are still necessary at this moment.

Simple approximatation schemes include the light-bulb scheme (Murphy & Burrows, 2008; Hanke et al., 2012; Couch, 2013), neutrino leakage (Janka & Mueller, 1996; Rosswog & Liebendörfer, 2003; O’Connor & Ott, 2011; Couch & O’Connor, 2014), and gray transport (Fryer et al., 1999; Scheck et al., 2006). In these schemes, the neutrino transport is relatively cheap and therefore, it is possible to perform high resolution simulations in 3D (with effective angular resolutions ), better describing the turbulence and convection behind the shock and the standing accretion shock instabilities (SASI, Blondin et al. 2003). However the neutrino transport in these schemes is possibly over-simplified because it does not really follow the transport of the neutrino distributions and therefore, cannot describe the neutrino heating self-consistently.

A more sophisticated but still approximated scheme for the neutrino transport is called the Moment scheme: The multi-group flux limited diffusion (MGFLD) scheme (Bruenn et al., 2013; Dolence et al., 2015) takes the zeroth angular moment of neutrino moment expansions (i.e. energy density). The M1 moment scheme additionally evolves the momentum density and considers the higher moment closure with an analytic form (Pons et al., 2000; Kuroda et al., 2012; O’Connor & Ott, 2013; Obergaulinger et al., 2014; O’Connor, 2014; Kuroda et al., 2015) or by a variable Eddington tensor (Burrows et al., 2000; Rampp & Janka, 2002; Thompson et al., 2003; Buras et al., 2006; Müller et al., 2012; Melson et al., 2015). However, as reported by Kuroda et al. (2015), the M1 scheme in 3D general relativity (GR) is still very expensive and difficult to run in long term postbounce simulations with high resolutions (an effective angular resolution ).

Another approximated transport scheme is the isotropic diffusion source approximation (IDSA, Liebendörfer et al. 2009). In the IDSA, the distribution function is separated into an opaque trapped particle component and a transparent streaming particle component. The two components are linked by a source term. Therefore, the transport equation becomes a diffusion problem in the opaque region and a ray tracing problem in the transparent region (Liebendörfer et al., 2009). Multi-dimensional simulations with the IDSA have been performed by Suwa et al. (2011, 2014) in 2D, and Takiwaki et al. (2014, 2012) in 3D. However, most of these multi-dimensional studies with the IDSA or with other transport schemes do not really solve the neutrino transport in multi-dimensions. Instead, they consider the angular distribution to be several 1D problems, i.e. apply the ray-by-ray (RbR) approach. Sumiyoshi et al. (2015) and Dolence et al. (2015) have pointed out that the RbR approach may artificially exaggerate the neutrino flux variations in the angular and temporal components, since the temporal variation of the neutrino fluxes in convective regions is greatly ignored in the RbR approach. Additionally, the RbR approximation overestimates the heating of accretion luminosity on its own ray and underestimates the heating in neighbor rays.

In this Paper, we present two-dimensional CCSN simulations with the IDSA for neutrino transport in cylindrical coordinates, which is in principle easy to extend to 3D and has better resolution and boundary conditions for the description of the proto-neutron star (PNS) than simulations in spherical coordinates. We solve the trapped particle component in multi-dimensions to improve the neutrino transport in angular and temporal directions. Furthermore, we study the neutrino transport effects during collapse by comparing our IDSA solver with a parametrized deleptonization scheme from Liebendörfer (2005). We find that the postbounce explosion dynamics is sensitive to the detailed neutrino interactions before core bounce, such as neutrino-electron scattering (NES) and electron-capture rates. Additional effects such as GR and magnetic fields are also crucial factors on studying CCSN explosion mechanism, in particular for fast-rotating and more massive progenitors, but we leave these parts for future work. In the following section, the numerical methods and the IDSA implementations are described. A verification of our IDSA implementation and a comparison with other neutrino transport schemes is shown in Section 3. In Section 4, we apply our IDSA implementation to multiple progenitors, using a new SN equation of state (EOS). In Section 5, we summarize our results and conclude. A discussion of the different EOS is presented in Appendix A.

2. Numerical Methods and Models

We describe the numerical code and the corresponding setup of our simulations in 2.1. The method and implementation of the IDSA for neutrino transport is demonstrated in 2.2. We present the investigated EOS and supernova progenitors in Sections 2.3 and 2.4.

Neutrino Interaction Reference
1. Bruenn (1985)
2. Bruenn (1985)
3. Bruenn (1985)
4. Bruenn (1985)
5. Mezzacappa & Bruenn (1993); Liebendörfer (2005)

Note. – = free neutrons, = free protons, = free neutrons or protons, = nuclei besides particles, =alpha particles, = electron type neutrinos, = anti-electron type neutrinos, = all type of neutrinos, =electron, and = positron. $\dagger$$\dagger$footnotetext: In this Paper, we include the neutrino-electron scattering by means of using the parametrized deleptonization in Liebendörfer (2005) in the pre-bounce phase, since this reaction is only important during core collapse.

Table 1Neutrino Interactions included in the Simulations

2.1. Numerical Code And Initial Setup

We use FLASH111 version 4 (Fryxell et al., 2000; Dubey et al., 2008; Lee, 2013) to solve the Eulerian equations of hydrodynamics,


where is the gas density, is the velocity vector, is the gas pressure, is the total specific energy, is the gravitational potential, is the electron fraction, is the particle number fraction, and is the particle mean specific energy. The index labels different species of trapped particles (i.e. and ). A detailed description of the neutrino-transport method will be presented in Section 2.2.

FLASH is a parallel, multidimensional hydrodynamics code based on block-structured Adaptive Mesh Refinement (AMR). Our simulation setup is essentially similar to what has been implemented in Couch (2013) and Couch & O’Connor (2014), but replacing their radiation transfer solver by an IDSA solver (Liebendörfer et al., 2009). We use the third-order piecewise parabolic method (PPM, Colella & Woodward 1984) for spatial reconstruction, the HLLC Riemann solver and the Hybrid slope limiter for shock-capture. The “hybrid” Riemann solver is widely used to avoid an odd-even numerical instability when the shock is aligned with the grid (Quirk, 1991). However, we don’t see this instability in our CCSN simulations by comparing simulation results with the hybrid Riemann solver. The HLLC Riemann solver shows a better turbulent cascade based on the implicit large eddy simulations (Radice et al., 2015). Effects from general and special relativity and from magnetic fields are ignored.

Simulations are performed in 1D spherical and 2D cylindrical coordinates. The center of a progenitor star is located at the origin of the simulation box. The simulation box includes the inner  km in radius of a progenitor star in 1D and the full in cylindrical coordinates with  km and  km in 2D. The central  km sphere has the smallest zone size of  km and the AMR mesh is dynamically adjusted based on the magnitude of the second derivatives of gas density, pressure, and entropy. To save computation time, we apply additional AMR decrements based on the distance to the origin. E.g. the first AMR decrement is enforced at  km and the second AMR decrement at  km, the next at  km, and so on. The maximum zone size is  km. We employ the “outflow” boundary condition at the outer boundaries and the “reflect” boundary condition at the inner boundaries. The gravitational potential is solved by the new improved multipole solver in FLASH (Couch et al., 2013) with a maximum Legendre order, .

It should be noted that the “outflow” boundary condition, as reported by Couch (2013), causes a zero-gradient mass accretion at the boundary which will overestimate the inflow and suppress the explosion at late time. We therefore extend our simulation domain to 10,000 km to minimize the effect of boundary conditions.

Progenitor (W02) (W02) (W02) (W02) (W07) (W07)
[g cm]
[g cm]
0.5 0.5 0.5 0.5 0.5 0.5
0.287 0.282 0.279 0.279 0.279 0.273
0.02 0.03 0.017 0.017 0.022 0.017

Note. – Parameters used in the fitting formula from Liebendörfer (2005)

Table 2Parameters for the Parametrized Deleptonization

2.2. Neutrino Transport

In the IDSA (Liebendörfer et al., 2009; Berninger et al., 2013), we decompose the distribution function of transported particle species and the neutrino transport operator into a trapped particle component and a streaming particle component. We assume that these two components evolve separately. With this assumption, we rewrite the transport equation, , where is a collision integral, by two equations,


where is the diffusion source term, which converts trapped particles () into streaming particles () and vice versa.

By using the diffusion limit, the diffusion source term in the trapped transport equation could be expressed by (see Liebendörfer et al. 2009 for a more detailed discussion),




is a non-local diffusion scalar, is the emissivity, is the absorptivity, is the opacity, and is the angle cosine. The distribution functions and can be solved by using equations (8) and (10) in Liebendörfer et al. (2009). Once the distribution functions are known, the net interaction rates of transported particles can be calculated and further updates of the fluid quantities, such as , , , , and , can be derived.

We solve the streaming component in 1D with the original IDSA solver, but solve the trapped component and in multi-dimensions. In each time step, we take the angle-averaged radial bins of fluid and neutrino quantities, such as , , , , , , as 1D inputs for the streaming component. The radial bins contain 600 zones sampled from the center of the proto-neutron star (PNS), defined by the location of maximum density, up to  km. The radial bins are uniformly spaced ( km) up to a radius of  km. Beyond 400 km, the zone spacing logarithmically increases. The streaming component carries only to the location of neutrino spheres and the angular integration of the spectral neutrino flux . The local heating from streaming neutrinos are then determined based on the local neutrino interaction rates based on the multi-dimensional hydro quantities.

It should be noted that the assumption of spherical symmetry in the streaming component may calculate the neutrino flux and heating incorrectly when shocks become highly asymmetric, especially in 2D simulations after the SASI has developed. On the other hand, a RbR approach may artificially enhance the asymmetry and lead to incorrect results as well (Sumiyoshi et al., 2015; Dolence et al., 2015). It is important to have multi-dimensional simulations with both approximations because they are likely to bracket the actual solution: The RbR focusses all heating from the accretion luminosity on its own ray, while the angular integration distributes it over all directions. In reality, a solution in between these extremes is expected.

For the trapped component, we explicitly solve the diffusion source and in multi-dimensions and update locally. Together with the streaming component, the new multi-dimensional interaction rates can be updated. The new neutrino sources are re-evaluated based on multi-dimensional neutrino sources and production rates, and then used for the next streaming step. Since we solve the diffusion part explicitly, a stable solution requires a small neutrino time step . Therefore, in our 1D and 2D simulations, we use a fixed hydro timestep with  sec and do sub-cycling for the neutrino transport. Since the IDSA solver in the streaming component allows larger time step, we adopt two types of sub-cycles. One for the streaming component with  sec, and the other for the trapped component with  sec. Additionally, the neutrino pressure contributes extra momentum and can be expressed by,


We note that Equations (6) and (11) include all thermodynamically important  terms of the Boltzmann transport equation. These  terms have been considered as crucial for CCSN modeling (Mezzacappa & Bruenn, 1993; Liebendörfer et al., 2009; Lentz et al., 2012a).

Our IDSA solver only includes electron flavor neutrinos. We use 20 energy bins which are spaced logarithmically from 3 MeV to 300 MeV. For heavy neutrinos, the energy release is treated by a leakage scheme that is based on the local diffusion and the local production rates (Hannestad & Raffelt, 1998; Rosswog & Liebendörfer, 2003; Perego et al., 2014).

Table 1 summarizes all weak interactions that are included in this work. All weak interactions in the IDSA solver are using the Bruenn description (Bruenn, 1985). Note that the IDSA solver dose not include the NES (interaction 5 in Table 1). NES is important during the collapse phase but provides minor contributions in the postbounce phase (Bruenn, 1989) (see Section 3.2 for a demonstration). Liebendörfer (2005) found a simple relation between the electron deleptonization with the gas density in the collapse phase, where electron fraction and entropy can be parametrized by density and chemical potentials. This parametrized deleptonization (PD) scheme can take interactions effectively into account that are only implemented in the Boltzmann solver that is used to determine the parameters for the PD. We have recalibrated the fitting parameters with different progenitors for the PD scheme by running AGILE-BOLTZTRAN (Liebendörfer et al., 2004), since the original parameters in Liebendörfer (2005) are calibrated for the progenitors G15 (s15s7b2, Woosley & Weaver 1995) and N13 (Nomoto & Hashimoto, 1988). Table 2 summarizes the fitting parameters we use during collapse. In this Paper we include the effect of NES by means of using the PD scheme in the collapse phase. Since the scheme is independent of the neutrino fractions, we update the neutrino fractions and energy through the IDSA solver during collapse. After bounce, we turn off the PD and switch to the IDSA solver. For simulations without NES, we use the IDSA solver from the beginning of core collapse to postbounce.

To verify our IDSA implementation in FLASH, we provide a code comparison of FLASH-IDSA with AGILE-IDSA and with existing results from the literature in Section 3.

2.3. Equation of State

We use the nuclear EOS unit in FLASH which incorporates the finite temperature EOS routines from O’Connor & Ott (2010) and Couch (2013)222 The Lattimer & Swesty EOS (with incompressibility,  MeV; Lattimer & Swesty 1991) and Hempel & Schaffner-Bielich (HS) DD2 EOS are used in this work. The HS(DD2) EOS employs the density-dependent relativistic mean-field interactions of Typel et al. (2010). The description of nuclei in supernova matter is based on (Hempel & Schaffner-Bielich, 2010). This EoS was first applied in core-collapse supernova simulations by Fischer et al. (2014), where further details can be found. LS220 is one of the most common and well-studied EOS in the supernova community. However, it has some deficiencies, for example it is based on the single nucleus approximation for heavy nuclei, and considers only the alpha particle as a degree of freedom of all possible light nuclei. See Hempel et al. (2015) for a comparison of predictions for cluster formation for the HS(DD2) and the LS220 EOS with experiments, where good agreement was found. Furthermore, it was shown by Krüger et al. (2013); Fischer et al. (2014), that the neutron matter EOS of LS220 is in disagreement with constraints from Chiral effective field theory. Furthermore, its low-density symmetry energy deviates from constraints obtained from finite nuclei, see Fig. 9 of Hempel (2014). No multi-dimensional simulations have been performed with HS(DD2) at this moment. We use LS220 for a code verification test in Section 3 and then use DD2 for our main simulations in Section 4. A low-density extension for the EOS is included in the routines from O’Connor & Ott (2010). In Appendex B, we provide a brief discussion of the differences between LS220 and DD2.

2.4. Supernova Progenitors

Figure 1.— Density as function of enclosed mass for the four investigated progenitor models from Woosley et al. (2002) and Woosley & Heger (2007).

In this Paper, we consider four different nonrotating, solar-metallicity progenitors, s11.0, s15.0, s21.0, and s27.0 from Woosley et al. (2002)333 for our multiple progenitor study. We also consider two nonrotating, solar-metallicity progenitors, s15 and s20 from Woosley & Heger (2007) for a comparison study. Figure 1 shows the initial density distribution of these six progenitors. The s11.0 progenitor has the highest core density but the most dilute envelope. The s21.0 and s27.0 and s20 progenitors have a similar density distribution and the most massive envelope among the models, but have a different mass of the iron core and Si/O shell. The locations of regions with a high density gradient correspond to the Si/O interface. For the same progenitor mass, s15 has a denser core and more massive envelope than s15.0. Another common progenitor model used in the literature is the s15s7b2 progenitor from Woosley & Weaver (1995). It has the same progenitor mass as and but the Si/O interface is located much further inside. This difference may make significant changes on the postbounce shock-radius evolutions when the shock reaches the interface due to a different mass-accretion history (Suwa et al., 2014).

To adopt the progenitor models in FLASH, we map the one-dimensional density, temperature, electron fraction, and radial velocity profiles from Woosley et al. (2002) into our 1D/2D grids in FLASH. Other thermodynamic quantities are re-calculated using the EOS solver in FLASH. Neutrino fractions are set to zero at beginning.

3. Idsa Verification

To verify the IDSA implementation in FLASH, we first compare our 1D FLASH simulations with simulations with AGILE-IDSA (Liebendörfer et al., 2009) in Section 3.1. AGILE-IDSA is a 1D spherically symmetric Lagrangian code which is publicly available online444 In Liebendörfer et al. (2009), a nice agreement of AGILE-IDSA with the GR Boltzmann code AGILE-BOLTZTRAN (Liebendörfer et al., 2004) was shown. Since we want to compare our results with the same neutrino transport scheme (IDSA), we run additional simulations in AGILE-IDSA but turn off the GR correction for the gravitational potential. We also turn off the PD during collapse in both codes to test the IDSA solver for neutrino transport before and after bounce.

In Section 3.2, we show 1D and 2D FLASH simulations with the two widely used progenitors, s15 and s20 from Woosley & Heger (2007) and discuss the differences compared to other transport schemes. Furthermore, we run AGILE-BOLTZTRAN simulations to demonstrate the importance of NES during core collapse and discuss the influence of the utilized neutrino opacities and electron capture rates on the shock evolution and neutrino signal.

Figure 2.— 1D radial profiles in mass coordinate during the collapse of the progenitor with LS220 EOS. Different colors represent different times when the central density reaches and at bounce. Solid lines show simulation results from FLASH and dashed lines show simulation results from AGILE-IDSA. Note that there are some slight time offsets between the results from FLASH and AGILE-IDSA because their output files do not exactly match the given central densities.

3.1. Code Comparison with AGILE-IDSA

We perform 1D CCSN simulations of the four investigated progenitors (, , , and ) with the LS220 EOS in both FLASH and AGILE-IDSA. Figure 2 shows the radial profiles of density, electron fraction, entropy and radial velocity of the progenitor during collapse. It is clear to see that the bounce profiles (black lines in Figure 2) are nearly identical in the two codes. Small differences at the center could originate from the fact that AGILE-IDSA is a Lagrangian code with a moving mesh. At the beginning of a simulation, the innermost region of the progenitor star is much less resolved than with the Eulerian grid of FLASH. Therefore, the radial profiles in AGILE-IDSA evolve slightly slower than those in FLASH. Furthermore, the bounce time in FLASH is about  ms later than in AGILE-IDSA, depending on the progenitor star.

Figure 3.— Electron (solid lines) and lepton (dashed lines) fraction profiles at core bounce in model . Red lines show simulation data from Agile-IDSA and blue lines represent data from FLASH. The black line indicates the electron fraction from the fitted PD scheme (Table 2). Both Agile-IDSA and FLASH are 1D simulations with LS220 EOS and without NES.

Figure 3 shows the electron and lepton fractions of the progenitor at core bounce. It is found that the electron fraction is consistent in both codes, but lepton fraction shows a mismatch at  g cm, suggesting a lower electron neutrino distribution at low densities region in FLASH. However, this difference does not alter the postbounce evolution much since the neutrinos are simply free streaming in the low density region. We note that in our multi-dimensional IDSA solver, we have employed a limiter for the diffusion scalar . The limiter enforces the trapped component to have a small value, when , to prevent unphysical oscillation on the , when the trapped neutrino density goes to zero in the free-streaming regime. This limiter leads to small differences of the neutrino luminosity and mean energy, but does not affect the hydrodynamic quantities. Figure 4 shows the particle spectra of the trapped and streaming components at 150 ms postbounce of the progenitor . We show three different regions where the trapped particle component dominates (at  km), trapped and streaming particle components are comparable (at  km), and where the streaming particle component dominates ( km). There are small differences in the particle spectra, but we do not expect exactly identical spectra in the two codes, since the hydrodynamic parts are different.

Figure 4.— Particle spectra at 150 ms postbounce for the trapped particle component (blue lines), streaming particle component (red lines), and their sum (green lines). Solid lines show spectra from FLASH and dashed lines represent spectra from AGILE-IDSA. Left panels show neutrino spectra and right panels show anti-neutrino spectra. Each row indicates spectra at a different radius.

A comparison of radial profiles of the density, electron fraction, entropy, and radial velocity of the progenitor in both codes at different postbounce times is shown in Figure 5. We note that although there are some differences in the radial shock location, the profiles and shock locations are consistent in the mass coordinate. The post bounce shock radius evolution for the four progenitor models are shown in Figure 6. Overall, the evolutions of the shock radius are nearly the same but show slight differences at  ms ( ms) in the model () when the shock reaches the Si/O interface.

Figure 7 shows the time evolution of electron-type neutrino luminosity and mean neutrino energy for the four considered models. The values are sampled at a radius of 500 km in both codes. The electron neutrino luminosities and electron anti-neutrino luminosities are similar for the two codes around bounce and early postbounce. However, at  ms post bounce, the neutrino and anti-neutrino luminosities in FLASH become slightly higher than in AGILE-IDSA for the progenitors , and . The maximum difference is an increase of about in FLASH at  ms post bounce. After  ms post bounce, the neutrino luminosities are back to the same value in both codes but the anti-neutrino luminosity remains slightly higher in FLASH simulations. For the progenitor , the electron neutrino luminosity is the same in both codes but the electron antineutrino luminosity is slightly higher in FLASH. However, we note that the mean neutrino and anti-neutrino energies in FLASH are about lower than in AGILE-IDSA. This difference could originate from the lower neutrino fraction in the low density region in FLASH (see Figure 3), since we measure the neutrino mean energy at  km.

The kernel of our IDSA solver in FLASH is the same as the solver in AGILE-IDSA. The main differences are from the hydrodynamics, the explicit solver for the diffusion scalar in Equation 10, the detailed treatment of the EOS, and potentially, also the gravity solver. The low-density extension for the EOS in FLASH and differences in the internal energy shift may also provide slightly differences. In principle, we should expect identical results in both codes in 1D. As presented above, although some little differences have been observed, most features are consistent and in nice agreement. Liebendörfer et al. (2005) have performed a code comparison of the Boltzmann solver AGILE-BOLTZTRAN with the variable Eddington factor VERTEX code. Both codes have sophisticated physics input in spherical symmetry but different implementations. As pointed out in Liebendörfer et al. (2005), the different grids in the Lagrangian or Eulerian coordinates produce late time differences in the shock evolution when the shock runs through shell interfaces. Since FLASH is also an Eulerian code, similar differences between FLASH and AGILE-IDSA as we have found here can be expected.

Figure 5.— Similar to Figure 2 but for the postbounce evolution.
Figure 6.— Time evolution of the shock radius for the four investigated progenitor models in both FLASH (solid lines) and AGILE-IDSA (dashed lines). Different colors represent different progenitor models. All results are 1D with LS220 EOS and without NES. The differences in the models and after  s are due to a different handling of shell interfaces and a large diffusivity of the implicit adaptive mesh in AGILE-IDSA.
Figure 7.— The electron neutrino (left) and electron antineutrino (right) luminosities and mean energies as a function of time. The quantities are averaged in a  ms time interval. Different colors represent the different progenitor models. Solid lines represent the FLASH simulations and dashed lines indicate the AGILE-IDSA simulations.
Property Models
Progenitor mass () 11 15 21 27
Progenitor compactness () 0.066 0.54 0.68 0.53
Progenitor compactness () 0.004 0.15 0.22 0.23
w PD wo PD w PD wo PD w PD wo PD w PD wo PD
Abbreviation (1D-) DP11 DA11 DP15 DA15 DP21 DA21 DP27 DA27
Time at g cm (ms) 176.3 194.1 267.6 278.2 269.7 284.3 273.0 287.2
Time at g cm (ms) 197.6 224.2 285.7 300.8 286.3 304.0 289.4 306.9
Time at g cm (ms) 201.9 229.3 289.6 305.5 290.2 308.4 293.2 311.4
Time at g cm (ms) 203.0 230.7 290.7 306.8 291.2 309.7 294.3 312.7
Time at bounce (ms) 203.5 231.3 291.2 307.3 291.7 310.2 294.7 313.2
Bounce, central ( g cm) 3.21 3.39 3.14 3.43 3.08 3.35 3.10 3.45
Bounce, central 0.287 0.316 0.282 0.312 0.279 0.310 0.279 0.311
Bounce, shock position () 0.54 0.77 0.56 0.76 0.56 0.75 0.56 0.75
  • Note that the compactness parameter of the s11.0 progenitor is about ten times smaller than others due to a smaller core mass and light envelope.

Table 31D Simulation Results During Core Collapse
Figure 8.— Time evolution of shock radius (left) and electron-type neutrino luminosities (right). Different colors represent different models from FLASH and AGILE-BOLTZTRAN with different neutrino physics. See Section 3.2 for detailed description.
Figure 9.— Time evolution of shock radius (left) and electron-type neutrino signals (right). Different colors represent different progenitor models. In the left panel, thick lines show the evolution of 2D simulations and thin lines represent 1D simulations. In the right panel, thick lines indicate electron-type neutrino luminosities and thin lines show the RMS neutrino mean energies. Solid lines represent electron neutrino signals and dashed lines represent electron antineutrino signals.

3.2. Code Comparisons for the 15 and 20 progenitors

The s15 and s20 progenitors from Woosley & Heger (2007) have been widely studied in the literature, e.g. by Bruenn et al. (2013), Dolence et al. (2015), Suwa et al. (2014), Hanke (2014, PhD Thesis), and Melson et al. (2015). In addition, Lentz et al. (2012a, b) have shown that different neutrino opacities and approximations could lead to different post-bounce evolutions by using spherically symmetric AGILE-BOLTZTRAN simulations with the s15 progenitor. Therefore these two progenitor models are very suitable for comparisons among different SN codes. It is generally agreed that the Vertex-Prometheus (Rampp & Janka, 2002) and the CHIMERA (Bruenn et al., 2013) codes use state-of-the-art physics for neutrino transport. Therefore, it is worth to perform simulations with these progenitor models to have a direct comparison. However, there is still different physics employed in the above mentioned works. For instance, CASTRO simulations (Dolence et al., 2015), ZEUS simulations (Suwa et al., 2014), and our FLASH-IDSA simulations are Newtonian, but AGILE-BOLTZTRAN simulations (Lentz et al. 2012a, b), CHIMERA simulations (Bruenn et al., 2013), and Vertex-Prometheus simulations (Hanke, 2014, and Melson et al. 2015) are GR or post-Newtonian. Furthermore, Dolence et al. (2015) use Shen’s EOS while other groups use the LS220 EOS. Therefore in this Section, we only give a qualitative discussion.

To evaluate that the PD scheme could effectively represent NES, we perform three Newtonian AGILE-BOLTZTRAN simulations with the s15 progenitor and the LS220 EOS. The first simulation (model AB-NR) includes the NES; the second simulation (model AB-NR-NoNES) ignores the NES; and the third simulation (model AB-NR-Mix) includes NES only in the collapse phase. Figure 8 shows the shock-radius and neutrino luminosity evolutions of these three simulations together with FLASH-IDSA simulations with and without PD. Ignoring NES (models FLASH-IDSA and AB-NoNES in Figure 8) gives a short-period of shock expansion at  ms postbounce. A signature from this can also be seen in the electron anti-neutrino luminosity. This behavior is consistent with the model “N-ReduceOp” in Lentz et al. (2012a), and the shock and electron-type neutrino luminosity evolutions of our FLASH-IDSA simulation without PD is also consistent with our AGILE-BOLTZTRAN simulation (AB-NR-NoNES).

The model AB-NR-Mix is nearly identical to the model AB-NR, demonstrating that NES is important mainly during the collapse phase (Figure 8). Fischer et al. (2012) showed that that NES plays a role in the late time PNS cooling and deleptonization ( s) after the explosion. However, our simulations focus only on the first few hundred milliseconds. The PD scheme in model FLASH-PD greatly improves the post-bounce simulation, making our FLASH-IDSA simulations closer to the AGILE-BOLTZTRAN simulation with full neutrino reactions. Although there is no perfect match, we conclude that the use of the PD scheme in the collapse phase could effectively take into account the NES. It should be noted that heavy neutrinos are treated by a simple leakage scheme which could lead to some difference in our IDSA simulations as well.

In Figure 9, we show the shock-radius evolution and neutrino signatures of our 2D FLASH-IDSA simulations of the s15 and s20 progenitors from Woosley & Heger (2007). To get a fair comparison, we enable the PD scheme during collapse and use the LS220 EOS. Both 2D models explode at  ms postbounce (see Table 4). The explosion time, is defined by the time when the averaged shock radius exceeds  km and never recedes at the end of the simulation. Overall, the shock-radius evolutions are similar to the models B15-WH07 and B20-WH07 in Bruenn et al. (2013) except that the CHIMERA simulations show earlier explosions. It should be noted that we use the old neutrino interactions from Bruenn (1985) and that our simulations are Newtonian. Müller et al. (2012) have shown that GR effects could enlarge the neutrino luminosities and therefore make it easier to explode. The model s20-2007 in Hanke et al. (2014) and Melson et al. (2015) shows a similar explosion time as our model s20, but the shock radius at  ms shrinks to  km due to GR effects. The reasons for the differences between Vertex-Prometheus and CHIMERA are still unclear, but the overall features for the progenitor s20 are still rather similar. However, the model s15-2007 in Hanke et al. (2014) shows a very different result. The shock stalls for  ms and then explodes at around  ms.

On the other hand, 2D CASTRO and ZEUS Newtonian simulations by Dolence et al. (2015) and Suwa et al. (2014) did not obtain explosion with the s15 and s20 progenitors. Our 2D simulations show a fast shock expansion after the prompt convection ( ms, see Figure 9). This is similar to what was observed in Dolence et al. (2015) but somewhat less dramatic. The prompt convection and fast shock expansion coincide with an oscillation of the electron antineutrino luminosity at  ms (see Figure 9 and Figure 6 of Dolence et al. 2015). These could be caused by the reduced opacity or incomplete neutrino interactions as discussed before. Note that Dolence et al. (2015) use the Shen EOS which is considered more difficult to lead to explosions than LS220 (Couch, 2013; Suwa et al., 2013).

Suwa et al. (2014) also use the IDSA (without PD) but with spherical coordinates and the “Ray-by-ray” approach. In principle, we should expect similar results, but the non-explosion of s15 and s20 in Suwa et al. (2014) suggest that the different hydrodynamics code, geometry, resolutions, and multidimensional neutrino transport approximation may also cause significant differences. Suwa et al. (2014) use 300 logarithmically spaced radial zones (from 1km to 5,000 km) and angular resolution. This is roughly three times lower than our simulations. A detailed code comparison is therefore necessary.

Progenitor Abbreviation EOS
(ms) (ms) (ms) (B) () (km)
s15 (W07) 1D-LA15-07 IDSA LS220 237 763 1.85 29.1
s15 (W07) 1D-LP15-07 PD LS220 249 523 1.77 36.1
s20 (W07) 1D-LP20-07 PD LS220 322 678 1.99 30.4
s11.0 (W02) 1D-LA11 IDSA LS220 206 794 1.48 27.5
s15.0 (W02) 1D-LA15 IDSA LS220 273 727 1.84 29.5
s21.0 (W02) 1D-LA21 IDSA LS220 274 726 1.98 30.4
s27.0 (W02) 1D-LA27 IDSA LS220 283 717 1.81 31.2
s11.0 (W02) 1D-DA11 IDSA DD2 231 766 1.47 32.2
s15.0 (W02) 1D-DA15 IDSA DD2 307 693 1.84 34.6
s21.0 (W02) 1D-DA21 IDSA DD2 310 604 1.95 37.1
s27.0 (W02) 1D-DA27 IDSA DD2 313 687 1.81 36.1
s11.0 (W02) 1D-DP11 PD DD2 203 761 1.47 33.6
s15.0 (W02) 1D-DP15 PD DD2 291 709 1.84 36.6
s21.0 (W02) 1D-DP21 PD DD2 292 708 1.98 36.6
s27.0 (W02) 1D-DP27 PD DD2 295 705 1.81 37.1
s15 (W07) 2D-LP15-07 PD LS220 249 312 524 0.298 1.69 37.1
s20 (W07) 2D-LP20-07 PD LS220 324 284 490 0.347 1.86 40.5
s15.0 (W02) 2D-LA15 IDSA LS220 274 209 311 0.464 1.60 46.8
s15.0 (W02) 2D-LA15low IDSA LS220 274 210 369 0.523 1.62 43.5
s11.0 (W02) 2D-DA11 IDSA DD2 232 86 374 0.821 1.31 42.3
s15.0 (W02) 2D-DA15 IDSA DD2 308 186 417 0.506 1.65 43.5
s21.0 (W02) 2D-DA21 IDSA DD2 311 189 411 0.635 1.75 44.8
s27.0 (W02) 2D-DA27 IDSA DD2 314 162 429 0.248 1.66 42.3
s11.0 (W02) 2D-DP11 PD DD2 204 170 376 0.267 1.37 47.4
s15.0 (W02) 2D-DP15 PD DD2 292 275 456 0.255 1.71 45.4
s21.0 (W02) 2D-DP21 PD DD2 292 282 484 0.417 1.82 44.8
s27.0 (W02) 2D-DP27 PD DD2 295 199 484 0.157 1.70 43.5
AAfootnotetext: Neutrino transport scheme during the collapse phase BBfootnotetext: Bounce time CCfootnotetext: Explosion time after bounce DDfootnotetext: Termination time after bounce EEfootnotetext: Diagnostic explosion energy at the end of the simulation FFfootnotetext: PNS mass at the end of the simulation GGfootnotetext: PNS radius (determined at the average radius corresponding to  g cm) at the end of the simulation $\dagger$$\dagger$footnotetext: The model abbreviations are defined by the model’s dimensionality, EOS, neutrino transport scheme, and progenitor mass. See Section 4 for a detailed description. $\dagger\dagger$$\dagger\dagger$footnotetext: The effective angular resolution in this model is a factor of two lower.
Table 4Simulation Parameters and Results

4. Multi-Progenitor Study

We perform 1D and 2D simulations with , , and progenitor models from Woosley et al. (2002). Simulations start from the prebounce core collapse to several hundred milliseconds postbounce with and without the PD in the collapse phase. The former is important to effectively take NES into account. Table 3 shows the core properties of these four progenitors during collapse based on 1D simulations. A summary of all performed simulations is shown in Table 4. The model abbreviations in Tables 3 and 4 are defined by a set of letters and numbers: The first two characters define the dimension of the model (1D or 2D); the first letter after the hyphen denotes the EOS of the model (L for LS220 and D for DD2); the second letter shows the transport scheme during the collapse (A for IDSA and P for PD); and the last two numbers specify the mass of the investigated progenitor model. A “’-07” at the end shows progenitor models from Woosley & Heger (2007), otherwise they are from Woosley et al. (2002). For instance, model 1D-DA15 means a 1D simulation of the progenitor with DD2 EOS and using the IDSA in the collapse phase (i.e. effectively without NES). When we refer to models DA, we consider all models with “DA” in its abbreviations.

4.1. Stellar Collapse and Core Bounce

Simulations are started from the non-rotating, solar-metallicity pre-supernova progenitors from Woosley et al. (2002) without artificial perturbation. Couch & Ott (2013) and Mueller & Janka (2014) show that small perturbations on the Si/O interface during collapse could amplify post-shock turbulence and turn a model that failed to explode towards a successful explosion. In addition, Couch et al. (2015) performed 3D simulations of the final few minutes of iron core growth in a massive star that includes Si shell burning. The results suggest that the non-spherical progenitor structure may have a strong impact on neutrino-driven explosions. In our models, we do not include these non-spherical features, and therefore the 2D simulations during collapse show nice agreement with 1D simulations in all models. However, we will show that spherical variations due to different levels of electron deleptonization or neutrino reactions during collapse may also have a significant impact on neutrino-driven explosions.

Figures 10-13 show the radial density, electron fraction, entropy, and radial velocity evolutions of 2D models with HS(DD2) EOS based on the IDSA or the PD at different time during the collapse phase. The 1D data are not shown because they are barely distinguishable from the 2D data before bounce. In the bounce profiles, the models DP show a slightly lower central density but a higher density distribution outside of the iron core, due to an earlier bounce time.

The core evolutions of different progenitors behave qualitatively similar during collapse, but show quantitive differences for the central electron fraction, central density, and bounce shock locations. Table 3 summarizes important quantities for the four different progenitors during the collapse phase. Due to a higher initial central density of the progenitor , the central density of reaches  g cm  ms earlier than other progenitors. Once the core densities in different progenitors reach the same , the more massive progenitors collapse faster than other progenitors.

In addition, since the electrons are highly degenerate and can only gain energy, this means that neutrinos lose their energy through NES and therefore escape more easily, accelerating the collapse process (Bethe, 1990). Therefore models DP (with effective NES) collapse faster than models DA (without NES) and have a lower (more captures), (escape) and at core bounce.

Core bounce is defined here by the first time when the maximum density in the core exceeds  g cm and the maximum peak entropy is above  k baryon. At bounce, the bounce shock emerges at (defined by the enclosed mass within the shock front at bounce) in models DP and at in models DA. Since the core mass is proportional to (Yahil & Lattimer, 1982), models DA have higher core mass than models DP at bounce. The highest infall velocity at bounce is about  km s in all progenitor models. It should be noted that the Si/O interface, which corresponds to the high entropy gradient region at about  km in Figures 10-13, is at different radii for models DA and DP due to different bounce times, since the bounce time directly sets the location of the Si/O interface. These difference will strongly influence the shock evolution after bounce.

Figure 10.— Angle averaged radial profiles of progenitor at different times during core collapse. All models are using the DD2 EOS. Different colors indicate different averaged profiles at certain central densities and bounce. The solid lines show simulations with the PD, and the dashed line show simulations with the IDSA.
Figure 11.— Similar to Figure 10 but for progenitor .
Figure 12.— Similar to Figure 10 but for progenitor .
Figure 13.— Similar to Figure 10 but for progenitor .
Figure 14.— Average shock radius vs. postbounce time for different progenitors and neutrino transport approximations. All models are using the DD2 EOS. Different colors indicate different progenitor models. The left panel represents simulations without the PD, and the right panel represents simulation with an effective inclusion of NES by the PD. The solid lines show 2D simulations, and the dashed lines show 1D simulations.

4.2. Shock Expansion and Instabilities

Figure 14 shows the angle averaged shock radius as a function of time in our 1D and 2D simulations with the DD2 EOS. The shock radius of a sample ray is defined here by the highest radius where the entropy is above   baryon (  baryon for progenitor models and ) after the postbounce time was reached  ms. In 2D simulations, we sample 180 radial rays. Before  ms postbounce, the shock front is defined by the minimum radial infall velocity.

1D and 2D models behave very similarly in the first few milliseconds until the bounce shock passes through the neutrino sphere at  ms. At that time, a prompt entropy- and electron-driven convection occurs and makes the 2D simulations different from 1D. Later on, this prompt convection causes a fast shock expansion, and changes the shock radius from  km to  km at  ms (see Figure 14). While the prompt convection is caused by the negative entropy and gradient, it should be noted that this early fast shock expansion during  ms may be amplified by our incomplete neutrino interactions and the old neutrino opacities, and by the grid-effect in cylindrical coordinates in consideration with the multidimensional treatment of neutrino diffusion, since a Rayleigh-Taylor bubble along the diagonal axis is observed at that time. A similar fast shock expansion is observed in Dolence et al. (2015) as well, but is not obvious in any other 2D simulation with spherical coordinates. However, unlike the results in Dolence et al. (2015), the shock radius in our simulations shrinks back to  km after this prompt convection at  ms.

Figure 15.— The color maps show the time evolution of the angle-averaged Brunt-Väisälä (BV) frequency in units of ms (left) and the anisotropic velocity in units of the speed of light c (right) in our models 2D-DP. The red lines show the time evolution of the maximum shock radius and the black (left) and white (right) lines show the time evolution of the averaged shock radius.

The prompt and late convection can be understood from the local stability analysis via the Ledoux criterion (Ledoux, 1947),


where we assume for simplicity. The Brunt-Väisälä (BV) frequency, , describes the linear growth frequency for convection. Follow the definition of (Buras et al., 2006; Ott et al., 2013), one obtains


where is the local gravitational potential and the approximation was used. Once convection is active, another useful quantity to describe the strength of convention is the anisotropic velocity. We follow the definition from Takiwaki et al. (2012) and define the anisotropic velocity as


In Figure 15, we show the evolution of the angle-averaged BV frequency and anisotropic velocity of our models DP. The BV frequencies are positive and high when the shock break through the neutrinosphere and the prompt convection happens at  ms. Starting from the core regions, the anisotropic velocities become very strong after the prompt convection has developed, and therefore drive the fast-shock expansion at  ms. However, the convection stops and the anisotropic velocity returns to small values at  ms. Models DA behave similar but on a different time scale due to different shock evolutions. We note that our estimate of may be incorrect at small radii, where neutrinos are trapped, due to the ignorance of the contribution from neutrinos in , but the overall features should be the same.

Figure 16.— Angle averaged radial profiles of progenitor at different postbounce times. All models are evolved in 2D and are using the DD2 EOS. Different colors indicate different postbounce times. The solid lines show simulations with the PD, and the dashed lines show simulations without the PD.
Figure 17.— Similar to Figure 16 but for progenitor .
Figure 18.— Similar to Figure 16 but for progenitor .
Figure 19.— Similar to Figure 16 but for progenitor .

For a given progenitor model, the bounce shocks have a similar shock strength in both models DA and DP. The difference of shock radius in models DA and DP during the prompt convection is mainly due to the difference of the bounce time, since different bounce time leads to a different shock strength when the bounce shock reaches the location of the negative entropy gradient. Since convection is suppressed in 1D, we don’t observe this prompt convection in 1D models. However, the prompt convection in 2D models leads to a larger shock radius than in 1D models at  ms.

After  ms, the shock stalls at  km for another  ms (except for model 2D-DA11). It is no surprise that all 1D models fail to explode due to the lack of convection. However, at first it seems surprising that all our 2D models explode within milliseconds postbounce regardless of the progenitor mass and the inclusion of NES during collapse. But on second thought, one has to consider that our simulations are Newtonian (comparing with Bruenn et al. 2001; Liebendörfer et al. 2001; Müller et al. 2012) and do not include the most recent electron capture rates (Langanke et al., 2003; Hix et al., 2003). In general, Newtonian simulations disfavor explosions compared to general relativistic simulations (Müller et al., 2012), and also the old neutrino interactions and opacities could make our simulations too optimistic with respect to explosions. Furthermore, we are using the DD2 EOS, which has more realistic nuclear matter properties than the LS220 EOS that was used in the above mentioned simulations. In Appendix B, we show that DD2 leads in our simulations to more favorable conditions for explosions than LS220. Table 4 summarizes the essential parameters of our 1D and 2D simulations at the end of the simulations. Model 2D-DA11 is the fastest explosion model that explodes already at  ms by neutrino-driven convection. For the other DA models, the explosion time is  ms. In addition, because of the effective NES, models DP evolve more slowly and explode at a late  ms ( ms for model 2D-DP11, see Table 4 for detailed information).

In model 2D-DP11, the shock radius first expands to  km at  ms postbounce when it reaches the Si/O interface, but temporarily drops back to  km at  ms postbounce before it explodes. Progenitor and explode later than progenitor and but at a similar time in both models 2D-DA and 2D-DP. In Figures 16-19, we compare the average radial profiles of density, entropy, electron fraction, and radial velocity at 1, 3, 50, 100, and 200 ms postbounce for models DA and DP. We terminate the simulations at  ms postbounce in 1D, and  ms postbounce in 2D.

The 2D entropy distribution of models 2D-DP are shown in Figure 20 at 150, 250, and 300 ms postbounce. At  ms, the convection is getting stronger (see Figure 15) but the shock radius distribution is still spherically symmetric. Later on, at  ms the SASI starts in models 2D-DP15, 21, and 27 but the models of progenitor show mainly convection without SASI. The SASI activities can be seen in Figure 21, where we show the normalized coefficients of the decomposition of the shock radius into Legendre polynomials (Müller et al., 2012). can be calculated by:


and corresponds to the averaged shock radius. For the progenitor , there is no obvious evidence of SASI activities in both models 2D-DA11 and 2D-DP11. The amplitude of the normalized coefficients for , 2 and 3 modes are small and within the same order of magnitude. For progenitors , , and , the SASI activities could be seen and start to grow at  ms postbounce. After  ms postbounce, the dipole () and quadrupole () modes grow to in the progenitors and , and in progenitor . The amplitudes do not show significant differences between models DA and DP, but the starting time of the high growth rates of the amplitudes corresponds to the time of fast shock expansion. SASI activities could also be seen in Figure 22 for the entropy distribution along the north and south poles.

To verify the code convergence, we have also performed a low-resolution run by reducing the angular resolution by a factor of 2 (model 2D-LA15low). We find that the explosion time is 1ms delayed and the shock expansion evolves slightly slower than for the standard resolution run (model 2D-LA15). Overall, we find no significant differences between the low-resolution run and the standard run.

Figure 20.— The entropy distributions of our 2D models with the PD and the DD2 EOS (See Table 4). Each frame shows a section of the domain spanning 800 km. The color scale indicates the entropy in baryon.
Figure 21.— First, second, and third coefficients , , and of the spherical decomposition of the shock radius into Legendre polynomials, normalized to the average shock radius () for different progenitors.
Figure 22.— Entropy along the north and the south poles as a function of time for different progenitor masses and for models DA (left) and DP (right) in Table 4.
Figure 23.— The diagnostic explosion energies defined in Section 4 as a function of time. Different colors represent the different progenitor models. Solid lines represent the models 2D-DP and dashed lines indicate the models 2D-DA in Table 4. The magenta and yellow lines use LS220 EOS, while other lines use DD2 EOS.
Figure 24.— Similar to Figure 7 but for the models DP (solid lines) and the models DA (dashed lines) in Table 4.
Figure 25.— The net heating (red) and cooling (blue) distributions of our models 2D-DP (See Table 4). The solid black lines show the radius of the PNS.
Figure 26.— The integral quantities of net heating (top), mass in the gain region (middle), and heating efficiency (bottom) as functions of time. See Section 4 for the definition of these integral quantities. Different colors represent different progenitor models. Solid lines represent the models 2D-DP and dashed lines indicate the models 2D-DA in Table 4.

4.3. Neutrino Heating and Explosion Energy

The energy released in neutrinos is about  erg in a typical CCSN (Colgate & White, 1966). To power the observed kinetic energy of a CCSN ( erg ), the baryonic matter has to absorb of the neutrino energy (Bethe & Pizzochero, 1990). However, the explosion energies in most published 2D models are still lower than the standard (except for those of Bruenn et al. 2013, 2014). The real explosion energy is not straightforward to calculate and the definition may differ from group to group. Figure 23 shows the “diagnostic energy” of our models 2D-DA and 2D-DP. The diagnostic energy, , is defined by


where the volume integration is performed over the region where the total specific energy, ,


is positive, represents the specific internal energy (thermal energy plus binding energy), and and are the specific kinetic and gravitational energies, respectively. is the reference energy value that is defined by the minimum specific internal energy at the beginning of the simulation, which leads to a negligible diagnostic explosion energy at beginning. We have checked that alternative values for the reference energy, that are suggested in the literature, do not have a significant effect on our final results. Depending on progenitor and simulation domain, is about  erg g. We find that the diagnostic energy is insensitive to if we make it a few times larger or smaller. Note that this expression does not take into account properly the different binding energy contributions from different nuclei to thermal energy. Therefore the so-called recombination energy is not included. corresponds only to the thermal energy that actually should be used in , if the composition was dominated by , or if the composition has a similar average binding energy like . The diagnostic explosion energies (Figure 23) of our models 2D-DP are between and at  ms postbounce and are still increasing at the end of our simulations. The models 2D-DA, which ignored the NES in the collapse phase, have a much stronger explosion energy, at  ms postbounce. This result is consistent with the 1D results in Bruenn (1989). It should be noted that the only difference between models 2D-DA and 2D-DP is the deleptonization method during collapse. The postbounce physics employed are identical in both models 2D-DA and 2D-DP, suggesting that the initial conditions at bounce may dramatically affect the postbounce evolutions. Furthermore, note that it is difficult to compare the results of models DP with other groups because there are differences in the input physics and methods. Additionally, we show the diagnostic explosion energies of models 2D-LA15 and 2D-LA15low for a comparison with LS220 EOS and with low resolution. We find that the model with LS220 EOS leads to a slight lower diagnostic explosion energy due to a delay of the explosion. The low resolution model also gives a similar results, suggesting that the fast explosion in our models is not due to insufficient resolution (Figure 23).

Ugliano et al. (2012); Nakamura et al. (2014), and Perego et al. (2015) suggest that there is a correlation of the compactness parameter with the explosion energy. We do see this trend in our simulations except for the progenitor , if we define the compactness parameter at an enclosed mass of , i.e. (O’Connor & Ott, 2011). In our simulations, the progenitor , which has the lowest compactness parameter, has the second highest diagnostic explosion energy (see Tables 3 and 4). However this trend in our simulations will disappear if we use . It should be noted that we only have four progenitor models in our simulations, and the correlation found in Ugliano et al. (2012) and Perego et al. (2015) show a huge scatter between different models, indicating that the compactness parameter oversimplifies the complexity of the explosion mechanism.

Figure 24 shows the neutrino luminosity and mean energy of electron type neutrinos and anti-neutrinos as functions of time. The neutrino mean energy, is defined by , where is the distribution function. The neutrino luminosities in models DA after  ms grow faster than in models DP but the peak luminosities after  ms are similar in both DA and DP models. The neutrino mean energies show a similar trend but the differences are less pronounced. We note that the neutrino luminosities and mean energies in our simulations are of a similar order of magnitude as the typical results in the literature (Janka, 2012). Therefore, the high explosion energies in our simulations require a larger mass in the gain region or a higher heating efficiency than other investigators.

Figure 25 presents the 2D net heating and cooling rates in models 2D-DP. The red color shows the heating region and the blue color indicates the cooling region. The PNS radius (defined by the average radius of density  g cm) is plotted as well in Figure 25. The gain region is defined as the post-shock region with positive net heating. In Figure 26, we show the integrated quantities of net heating (), total mass in the gain region, and the heating efficiency as functions of time. The heating efficiency is defined by , since the heating originates mainly from the electron type neutrinos.

The early prompt convection shows a very high heating efficiency peak at  ms, which could be associated with a grid-effect in cylindrical coordinates. However, this prompt heating does not affect the explosion energy directly since matter is still bound at that time. The shock expansion and the early convection at  ms enlarge the gain radius and the mass in the gain region. After  ms, the mass in the gain region continues to increase in time. We noted that the masses in the gain region () in our simulations are higher than the values that are reported by other groups, explaining the high diagnostic explosion energy in our simulations.

The models 2D-DA11 and 2D-DP11 have the lowest net heating and heating efficiency among the models (except during the prompt convection period), but have the second highest mass in the gain region and diagnostic energy. This is probably due to the fast shock expansion (see Figure 14) and low mass-accretion in the progenitor .

Figure 27.— The time evolution of the PNS mass (left) and radius (right). The PNS mass and radius are determined at the average radius corresponding to  g cm. Different colors represent different progenitor stars. Solid lines show the models 2D-DP and dashed lines show the models 2D-DA in Table 4.
Figure 28.— The time evolution of central density (, upper left), central electron fraction () and lepton fraction (, upper right), central entropy (lower left), and central temperature (lower right). Different colors represent different progenitor stars. Solid lines show the models 2D-DP and dashed lines show the models 2D-DA in Table 4. In the upper right panel, the thin lines represent the central lepton fraction and the thick lines represent the central electron fraction.

4.4. Proto Neutron Stars

During its evolution, the PNS shrinks due to the neutrino cooling behind the gain radius. In Figure 27, we show the mass and radius evolutions of the PNSs after bounce. An initial oscillation of the PNS is taking place at  ms, when the shock breaks through the neutrino sphere. Later on, mass continuously accretes onto the PNS for several hundred milliseconds so that it reaches at the end of the simulations. At that time, the PNS radii have shrunken to  km. The mass and radii of the PNSs at the end of the simulations are summarized in Table 4.

The convection and SASI activities help the PNS to accrete mass during the setup of an explosion. Accretion funnels are created along the high entropy convective plumes. This can be seen in Figure 20. Therefore, the progenitor which has the less SASI activities and accretion funnels, shows the lowest PNS mass increment and the growth rate of the PNS mass is the lowest as well. Furthermore, the progenitor has the lowest density distribution outside the iron core, which gives the lowest mass accretion rate. After  ms, the mass increment of the PNS becomes small but the mass of the PNS is still increasing at the end of the simulations. We find that although the growth rate of the PNS mass in models 2D-DA and 2D-DP are different, the final PNS mass in models 2D-DP are only slightly larger (), and the radii of PNSs in models 2D-DA are slightly lower than that in models 2D-DP, giving a more compact PNS core in models 2D-DA. This can also be related to the earlier explosion times of models DA.

Figure 28 shows the central density, central electron, lepton fraction, central entropy and central temperature as functions of time. Models 2D-DA show higher central densities and higher central (higher electron pressures), but lower central entropies and central temperatures than the models 2D-DP, giving a more compact PNS core. However, the density in the outer layer of the PNS is higher and more convective in models 2D-DP than that in models 2D-DA. Therefore, models 2D-DP have a higher PNS mass and radius (see Figure 27).

5. Conclusions

We have performed self-consistent CCSN simulations for the four nonrotating, solar-metallicity progenitors, s11.0, s15.0, s21.0, and s27.0 from Woosley et al. (2002) and s15 and s20 from Woosley & Heger (2007) using the AMR code FLASH with an IDSA for neutrino transport. A very good agreement of FLASH-IDSA with AGILE-IDSA is shown in Section 3.1, though some small differences still exist. In addition, we have provided a comparison of our multi-dimensional IDSA results for the s15 and s20 progenitors from Woosley & Heger (2007) with results in the existing literature. However, a detailed comparison with an exchange of data and collaboration among groups would be necessary to understand remaining different results with regard to the underlying physics.

The new SN EOS HS (DD2) is employed in this study and a comparison with the standard LS220 EOS is discussed in Appendix A. It is found that the DD2 EOS shows a faster shock expansion, and a higher mass in the gain region, while the neutrino luminosities are similar, causing the DD2 EOS to lead to slightly earlier and stronger explosions than the LS220.

We have presented two sets of simulations which compared the neutrino transport with the IDSA (but ignoring NES; abbreviated DA) with the parametrized deleptonization scheme (which includes NES effectively; abbreviated DP) during collapse. The results show clearly that the treatment of neutrino weak interactions and the level of deleptonziation during the collapse phase have a significant impact on the neutrino-driven explosions. All our 2D models explode within  ms. Models without NES explode much easier, stronger and faster than models with NES, The diagnostic explosion energy in our simulations are around in models 2D-DP (see Table 4) and around  B in models 2D-DA. Our explosion energies are likely to decrease when we include general relativistic effects and better electron-capture rates in future models.

In our simulations, we do not use the typical RbR approach for the neutrino transport. Instead, we solve the diffusion part in multi-dimensions, improving the neutrino transport in angular and temporal directions. With this approach, the prompt convection after the core bounce causes a fast shock expansion ( ms), enlarging the gain region at late time, and therefore helps to increase the explosion energy at late time. However, this prompt convection may be amplified by grid-effects in cylindrical coordinates. To address this issue, a full 3D simulation will be necessary or a detailed comparison with simulations in spherical coordinates. Furthermore, we have found SASI activities at late time in progenitors , , and , but in most cases, the explosions are mainly due to neutrino-driven convection. Future work will relax the constraints imposed by axisymmetry. Models in 3D will permit to study a more realistic turbulence cascade and rotation.

We thank the anonymous referee for his/her valuable comments and suggestions. This work was supported by the European Research Council (ERC; FP7) under ERC Advanced Grant Agreement N 321263 - FISH, by the PASC High Performance Computing Grant DIAPHANE, and by the Swiss National Science Foundation (SNF). The Basel group is member of the COST Action New Compstar. Part of this work was inspired during the MICRA2013 meeting in ECT*. We thank Albino Perego for the development of the leakage scheme for heavy neutrinos. We thank Takami Kuroda, Rubén Cabezón, Kei Kotake, Tomoya Takiwaki, Ko Nakamura, and Yudai Suwa for useful discussions. K.C.P. acknowledges Sean Couch for releasing his supernova setup in the public version of FLASH. FLASH was developed largely by the DOE-supported ASC/Alliances Center for Astrophysical Thermonuclear Flashes at the University of Chicago. The simulations have been carried out at the CSCS Monte-Rosa under grant No. 412. Analysis and visualization of simulation data were completed using the analysis toolkit yt (Turk et al., 2011).

Appendix a
A Comparsion Between Dd2 and Ls220

Figure 29.— Similar to Figure 10 but for models 2D-DA15, and 2D-LA15 in Table 4.
Figure 30.— Similar to Figure 24 but for models 1D-DA15, 1D-LA15, 2D-DA15, and 2D-LA15 in Table 4.

Several simulations in this paper use the new DD2 EOS, while most multi-dimensional simulations found in the literatures use LS220. We therefore briefly describe the qualitative differences between the simulations with the DD2 and the LS220 EOS in this section. We perform 1D and 2D simulations of the progenitor with DD2 and LS220 EOS (models 1D-DA15, 2D-DA15, 1D-LA15 and 2D-LA15 in Table 4). Hereafter, we use DA15 to refer to both models 1D-DA15 and 2D-DA15, and similarly we use LA15 to refer to the models 1D-LA15 and 2D-LA15. Their main features are summarized in Table 4). The PD scheme is turned off in order to study the impact of the EOS on the neutrino transport.

Figure 31.— Similar to Figure 28 but for models 1D-DA15, 1D-LA15, 2D-DA15, and 2D-LA15 in Table 4.

The angle-averaged radial profiles of density, electron fraction, entropy, and radial velocity are shown in Figure 29 for the prebounce phase. As discussed in Section 4, the radial profiles in the collapse phase are nearly identical in both 1D and 2D. In addition, the density and radial velocity evolutions in models DA15 and LA15 are also similar during the collapse but models DA15 collapse slower than models LA15. Therefore, models DA15 reach core bounce about  ms later than models LA15.

Figure 32.— Similar to Figure 27 but for models 2D-DA15, and 2D-LA15 in Table 4.

While there are no significant differences of the neutrino luminosity and mean energy between models DA15 and LA15 (Figure 30), models DA15 show a lower central density and higher PNS radius than LA15. The time evolutions of central density and electron fraction of models DA15 and LA15 are shown in Figure 31. Furthermore, models DA15 have a higher central electron and lepton fraction than models DL15. 2D models also show a lower PNS mass and a higher PNS radius than 1D (Figure 32), but the PNS masses are similar in both DD2 and LS220 (see Table 4).

Figure 33.— Similar to Figure 26 but for models 2D-DA15, and 2D-LA15 in Table 4.
Figure 34.— Similar to Figure 14 but for models 1D-DA15, 1D-LA15, 2D-DA15, and 2D-LA15 in Table 4.
Figure 35.— Similar to Figure 16 but for models 1D-DA15, 1D-LA15, 2D-DA15, and 2D-LA15 in Table 4.

Since the bounce time in models DA15 take longer than models LA15, the Si/O interface reaches a smaller radius in models DA15 ( km) than in LA15 ( km, see Figure 29). The bounce shock hits the Si/O interface earlier in models DA15 than in LA15, generating a prompt convection and fast shock expansion at  ms postbounce. The earlier interaction between the bounce shock with the Si/O interface in model 2D-DA15 makes the shock radius in model 2D-DA15 larger than in simulation 2D-LA15. Furthermore, the larger shock radius in model 2D-DA15 produces a higher mass in the gain region (see Figure 33), suggesting that the models with the DD2 EOS explodes more easily than the models with the LS220 EOS. The averaged shock radius evolution of models DA15 and LA15 can be seen in Figure 34. The explosion time () in 2D-DA15 is about  ms earlier than in 2D-LA15. In Figure 35, we compare the radial density, electron fraction, entropy, and velocity profiles in 1D and 2D for DD2 and LS220 at 1, 3, 50, 100, and 200 ms postbounce. Before  ms postbounce, models LS15 show a faster shock expansion than models DA15 but the shock radius in DA15 become larger after the Si/O interface has reached the shock at  ms.

Figure 36.— Similar to Figure 21 but for models 2D-LA15 (Left) and 2D-LD15 (Right) in Table 4.
Figure 37.— Similar to Figure 22 but for models 2D-LA15 (Left) and 2D-DA15 (Right) in Table 4.

Figure 36 shows the normalized decomposition into Legendre polynomials (Equation 15) of the shock radius of models 2D-DA15 and 2D-LA15. The corresponding entropy distribution at the north and south pole is shown in Figure 37. The SASI activities start at  s in both 2D-DA15 and 2D-LA15 after the Si/O interface reach the shock and the amplitude of the normalized coefficients do not show significant difference. We remark that a more through investigation of the EOS effects will be the subject of a future study.


  • Abdikamalov et al. (2014) Abdikamalov, E., Ott, C. D., Radice, D., et al. 2014, arXiv:1409.7078
  • Bethe (1990) Bethe, H. A. 1990, Reviews of Modern Physics, 62, 801
  • Bethe & Pizzochero (1990) Bethe, H. A., & Pizzochero, P. 1990, ApJ, 350, L33
  • Berninger et al. (2013) Berninger, H., Frenod, E., Gander, M., Liebendorfer, M., & Michaud, J. 2013, SIAM Journal on Mathematical Analysis, 45, 3229
  • Blondin et al. (2003) Blondin, J. M., Mezzacappa, A., & DeMarino, C. 2003, ApJ, 584, 971
  • Brandt et al. (2011) Brandt, T. D., Burrows, A., Ott, C. D., & Livne, E. 2011, ApJ, 728, 8
  • Bruenn (1985) Bruenn, S. W. 1985, ApJS, 58, 771
  • Bruenn (1989) Bruenn, S. W. 1989, ApJ, 340, 955
  • Bruenn (1989) Bruenn, S. W. 1989, ApJ, 341, 385
  • Bruenn et al. (2001) Bruenn, S. W., De Nisco, K. R., & Mezzacappa, A. 2001, ApJ, 560, 326
  • Bruenn et al. (2013) Bruenn, S. W., Mezzacappa, A., Hix, W. R., et al. 2013, ApJ, 767, LL6
  • Bruenn et al. (2014) Bruenn, S. W., Lentz, E. J., Hix, W. R., et al. 2014, arXiv:1409.5779
  • Buras et al. (2006) Buras, R., Rampp, M., Janka, H.-T., & Kifonidis, K. 2006, A&A, 447, 1049
  • Burrows (2013) Burrows, A. 2013, Reviews of Modern Physics, 85, 245
  • Burrows et al. (2000) Burrows, A., Young, T., Pinto, P., Eastman, R., & Thompson, T. A. 2000, ApJ, 539, 865
  • Colella & Woodward (1984) Colella, P., & Woodward, P. R. 1984, Journal of Computational Physics, 54, 174
  • Colgate & White (1966) Colgate, S. A., & White, R. H. 1966, ApJ, 143, 626
  • Couch (2013) Couch, S. M. 2013, ApJ, 765, 29
  • Couch et al. (2013) Couch, S. M., Graziani, C., & Flocke, N. 2013, ApJ, 778, 181
  • Couch & Ott (2013) Couch, S. M., & Ott, C. D. 2013, ApJ, 778, LL7
  • Couch & O’Connor (2014) Couch, S. M., & O’Connor, E. P. 2014, ApJ, 785, 123
  • Couch et al. (2015) Couch, S. M., Chatzopoulos, E., Arnett, W. D., & Timmes, F. X. 2015, arXiv:1503.02199
  • Dolence et al. (2015) Dolence, J. C., Burrows, A., & Zhang, W. 2015, ApJ, 800, 10
  • Dubey et al. (2008) Dubey, A., Reid, L. B., & Fisher, R. 2008, Physica Scripta Volume T, 132, 014046
  • Fischer et al. (2012) Fischer, T., Martínez-Pinedo, G., Hempel, M., & Liebendörfer, M. 2012, Phys. Rev. D, 85, 083003
  • Fischer et al. (2014) Fischer, T., Hempel, M., Sagert, I., Suwa, Y., & Schaffner-Bielich, J. 2014, European Physical Journal A, 50, 46
  • Foglizzo et al. (2015) Foglizzo, T., Kazeroni, R., Guilet, J., et al. 2015, arXiv:1501.01334
  • Fryer et al. (1999) Fryer, C., Benz, W., Herant, M., & Colgate, S. A. 1999, ApJ, 516, 892
  • Fryxell et al. (2000) Fryxell, B., Olson, K., Ricker, P., et al. 2000, ApJS, 131, 273
  • Hanke et al. (2012) Hanke, F., Marek, A., Müller, B., & Janka, H.-T. 2012, ApJ, 755, 138
  • Hannestad & Raffelt (1998) Hannestad, S., & Raffelt, G. 1998, ApJ, 507, 339
  • Hempel & Schaffner-Bielich (2010) Hempel, M., & Schaffner-Bielich, J. 2010, Nuclear Physics A, 837, 210
  • Hempel et al. (2012) Hempel, M., Fischer, T., Schaffner-Bielich, J., & Liebendörfer, M. 2012, ApJ, 748, 70
  • Hempel (2014) Hempel, M. 2014, arXiv:1410.6337
  • Hempel et al. (2015) Hempel, M., Hagel, K., Natowitz, J., Röpke, G., & Typel, S. 2015, arXiv:1503.00518
  • Hix et al. (2003) Hix, W. R., Messer, O. E., Mezzacappa, A., et al. 2003, Physical Review Letters, 91, 201102
  • Janka (2012) Janka, H.-T. 2012, Annual Review of Nuclear and Particle Science, 62, 407
  • Janka & Mueller (1996) Janka, H.-T., & Mueller, E. 1996, A&A, 306, 167
  • Krüger et al. (2013) Krüger, T., Tews, I., Hebeler, K., & Schwenk, A. 2013, Phys. Rev. C, 88, 025802 & Liebendörfer, M. 2011, ApJS, 195, 20
  • Kuroda et al. (2012) Kuroda, T., Kotake, K., & Takiwaki, T. 2012, ApJ, 755, 11
  • Kuroda et al. (2015) Kuroda, T., Takiwaki, T., & Kotake, K. 2015, arXiv:1501.06330
  • Langanke et al. (2003) Langanke, K., Martínez-Pinedo, G., Sampaio, J. M., et al. 2003, Physical Review Letters, 90, 241102
  • Lattimer & Swesty (1991) Lattimer, J. M., & Douglas Swesty, F. 1991, Nuclear Physics A, 535, 331
  • Ledoux (1947) Ledoux, P. 1947, ApJ, 105, 305
  • Lee (2013) Lee, D. 2013, Journal of Computational Physics, 243, 269
  • Lentz et al. (2012a) Lentz, E. J., Mezzacappa, A., Messer, O. E. B., et al. 2012, ApJ, 747, 73
  • Lentz et al. (2012b) Lentz, E. J., Mezzacappa, A., Messer, O. E. B., Hix, W. R., & Bruenn, S. W. 2012, ApJ, 760, 94
  • Liebendörfer et al. (2001) Liebendörfer, M., Mezzacappa, A., & Thielemann, F.-K. 2001, Phys. Rev. D, 63, 104003
  • Liebendörfer et al. (2001) Liebendörfer, M., Mezzacappa, A., Thielemann, F.-K., et al. 2001, Phys. Rev. D, 63, 103004
  • Liebendörfer et al. (2004) Liebendörfer, M., Messer, O. E. B., Mezzacappa, A., et al. 2004, ApJS, 150, 263
  • Liebendörfer (2005) Liebendörfer, M. 2005, ApJ, 633, 1042
  • Liebendörfer et al. (2005) Liebendörfer, M., Rampp, M., Janka, H.-T., & Mezzacappa, A. 2005, ApJ, 620, 840
  • Liebendörfer et al. (2009) Liebendörfer, M., Whitehouse, S. C., & Fischer, T. 2009, ApJ, 698, 1174
  • Livne et al. (2004) Livne, E., Burrows, A., Walder, R., Lichtenstadt, I., & Thompson, T. A. 2004, ApJ, 609, 277
  • Melson et al. (2015) Melson, T., Janka, H.-T., & Marek, A. 2015, ApJ, 801, L24
  • Melson et al. (2015) Melson, T., Janka, H.-T., Bollig, R., et al. 2015, ApJ, 808, L42
  • Mezzacappa & Bruenn (1993) Mezzacappa, A., & Bruenn, S. W. 1993, ApJ, 405, 669
  • Mezzacappa & Bruenn (1993) Mezzacappa, A., & Bruenn, S. W. 1993, ApJ, 410, 740
  • Mueller & Janka (2014) Mueller, B., & Janka, H.-T. 2014, arXiv:1409.4783
  • Müller et al. (2012) Müller, B., Janka, H.-T., & Marek, A. 2012, ApJ, 756, 84
  • Müller et al. (2012) Müller, B., Janka, H.-T., & Heger, A. 2012, ApJ, 761, 72
  • Murphy & Burrows (2008) Murphy, J. W., & Burrows, A. 2008, ApJ, 688, 1159
  • Nakamura et al. (2014) Nakamura, K., Takiwaki, T., Kuroda, T., & Kotake, K. 2014, arXiv:1406.2415
  • Nomoto & Hashimoto (1988) Nomoto, K., & Hashimoto, M. 1988, Phys. Rep., 163, 13
  • Obergaulinger et al. (2014) Obergaulinger, M., Janka, H.-T., & Aloy, M. A. 2014, MNRAS, 445, 3169
  • O’Connor & Ott (2010) O’Connor, E., & Ott, C. D. 2010, Classical and Quantum Gravity, 27, 114103
  • O’Connor & Ott (2011) O’Connor, E., & Ott, C. D. 2011, ApJ, 730, 70
  • O’Connor & Ott (2013) O’Connor, E., & Ott, C. D. 2013, ApJ, 762, 126
  • O’Connor (2014) O’Connor, E. 2014, arXiv:1411.7058
  • Ott et al. (2008) Ott, C. D., Burrows, A., Dessart, L., & Livne, E. 2008, ApJ, 685, 1069
  • Ott et al. (2013) Ott, C. D., Abdikamalov, E., Mösta, P., et al. 2013, ApJ, 768, 115
  • Perego et al. (2014) Perego, A., Rosswog, S., Cabezón, R. M., et al. 2014, MNRAS, 443, 3134
  • Perego et al. (2015) Perego, A., Hempel, M., Fröhlich, C., et al. 2015, arXiv:1501.02845
  • Pons et al. (2000) Pons, J. A., Ibáñez, J. M., & Miralles, J. A. 2000, MNRAS, 317, 550
  • Quirk (1991) Quirk, J. J. 1991, Ph.D. Thesis,
  • Radice et al. (2015) Radice, D., Couch, S. M., & Ott, C. D. 2015, arXiv:1501.03169
  • Rampp & Janka (2002) Rampp, M., & Janka, H.-T. 2002, A&A, 396, 361
  • Rosswog & Liebendörfer (2003) Rosswog, S., & Liebendörfer, M. 2003, MNRAS, 342, 673
  • Scheck et al. (2006) Scheck, L., Kifonidis, K., Janka, H.-T., Müller, E. 2006, A&A, 457, 963
  • Shen et al. (1998) Shen, H., Toki, H., Oyamatsu, K., & Sumiyoshi, K. 1998, Nuclear Physics A, 637, 435
  • Sumiyoshi et al. (2015) Sumiyoshi, K., Takiwaki, T., Matsufuru, H., & Yamada, S. 2015, ApJS, 216, 5
  • Sumiyoshi et al. (2005) Sumiyoshi, K., Yamada, S., Suzuki, H., et al. 2005, ApJ, 629, 922
  • Suwa et al. (2011) Suwa, Y., Kotake, K., Takiwaki, T., Liebendörfer, M., & Sato, K. 2011, ApJ, 738, 165
  • Suwa et al. (2013) Suwa, Y., Takiwaki, T., Kotake, K., et al. 2013, ApJ, 764, 99
  • Suwa et al. (2014) Suwa, Y., Yamada, S., Takiwaki, T., & Kotake, K. 2014, arXiv:1406.6414
  • Takiwaki et al. (2014) Takiwaki, T., Kotake, K., & Suwa, Y. 2014, ApJ, 786, 83
  • Takiwaki et al. (2012) Takiwaki, T., Kotake, K., & Suwa, Y. 2012, ApJ, 749, 98
  • Thompson et al. (2003) Thompson, T. A., Burrows, A., & Pinto, P. A. 2003, ApJ, 592, 434
  • Turk et al. (2011) Turk, M. J., Smith, B. D., Oishi, J. S., et al. 2011, ApJS, 192, 9
  • Typel et al. (2010) Typel, S., Röpke, G., Klähn, T., Blaschke, D., & Wolter, H. H. 2010, Phys. Rev. C, 81, 015803
  • Ugliano et al. (2012) Ugliano, M., Janka, H.-T., Marek, A., & Arcones, A. 2012, ApJ, 757, 69
  • Woosley & Weaver (1995) Woosley, S. E., & Weaver, T. A. 1995, ApJS, 101, 181
  • Woosley et al. (2002) Woosley, S. E., Heger, A., & Weaver, T. A. 2002, Reviews of Modern Physics, 74, 1015
  • Woosley & Heger (2007) Woosley, S. E., & Heger, A. 2007, Phys. Rep., 442, 269
  • Yahil & Lattimer (1982) Yahil, A., & Lattimer, J. M. 1982, NATO Advanced Science Institutes (ASI) Series C, 90, 53
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

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