ABSTRACT
By performing neutrino-radiation hydrodynamic simulations in spherical symmetry (1D) and axial symmetry (2D) with different progenitor models by Woosley & Heger from 12 to 100 M⊙, we find that all 1D runs fail to produce an explosion and several 2D runs succeed. The difference in the shock evolutions for different progenitors can be interpreted by the difference in their mass accretion histories, which are in turn determined by the density structures of progenitors. The mass accretion history has two phases in the majority of the models: the earlier phase, in which the mass accretion rate is high and rapidly decreasing, and the later phase, with a low and almost constant accretion rate. They are separated by the so-called turning point, the origin of which is a change of the accreting layer. We argue that shock revival will most likely occur around the turning point and hence that its location in the plane will be a good measure for the possibility of shock revival: if the turning point lies above the critical curve and the system stays there for a long time, shock revival will obtain. In addition, we develop a phenomenological model to approximately evaluate the trajectories in the plane, which, after calibrating free parameters by a small number of 1D simulations, reproduces the location of the turning point reasonably well by using the initial density structure of progenitor alone. We suggest the application of the phenomenological model to a large collection of progenitors in order to infer without simulations which ones are more likely to explode.
1. INTRODUCTION
A core-collapse supernova is one of the most energetic explosions in the universe. Although the explosion mechanisms are yet to be uncovered, there are a few possible models. Among them, the neutrino heating mechanism (Bethe & Wilson 1985) is the most promising scenario, in which copious neutrinos are emitted in the vicinity of proto-neutron stars (PNSs) and are partially absorbed by postshock material. In this system, neutrinos transfer internal energy from inside to the outside of PNSs and act effectively as a heating source for the postshock layer.
Although this neutrino heating mechanism certainly works, state-of-the-art simulations of neutrino-radiation hydrodynamics cannot produce an explosion in spherical symmetry (Rampp & Janka 2000; Liebendörfer et al. 2001; Thompson et al. 2003; Sumiyoshi et al. 2005). Recently, modern multidimensional simulations became possible, and several exploding simulations have been reported (e.g., in two dimensions (2D), Buras et al. 2006; Marek & Janka 2009; Suwa et al. 2010; Müller et al. 2012a; Bruenn et al. 2013; Pan et al. 2015; in three dimensions [3D], Takiwaki et al. 2012; Lentz et al. 2015; Melson et al. 2015b; Müller 2015). These simulations were limited to a relatively small number of progenitors7 :
- 1.
- 2.
- 3.
- 4.
- 5.
Other progenitor models, 8.1 M⊙ (Müller et al. 2012b) and 9.6 M⊙ (Müller et al. 2013; Melson et al. 2015b), were also investigated. More recently, Bruenn et al. (2013) performed a systematic study using a progenitor series of Woosley & Heger (2007) from 12 to 25 M⊙ and found similar explosions for all progenitors. Dolence et al. (2015) reported, however, that they found that none of them resulted in an explosion.8 More recently, Melson et al. (2015a) performed 2D and 3D simulations with 20 M⊙ of the same series and found an explosion in 2D and failure in 3D with standard neutrino opacities. In this study, we perform two-dimensional simulations for a broader mass range from 12 to 100 M⊙ using the same progenitor series of Woosley & Heger (2007).
The progenitor structure is one of the most important ingredients in the core-collapse supernova explosion mechanism9 because it determines the initial condition and later the accretion rate history. The latter has a strong leverage on the shock wave evolution, since the force balance between the ram pressure of preshocked material and the thermal pressure of postshocked material is the main factor to determine the shock position and the mass accretion rate, , is a good measure of the ram pressure, ρv2. Indeed, it has been demonstrated that the O–Ne–Mg core of an 8.8 M⊙ star can produce an explosion even in spherical symmetry thanks to a rapid decrease in the mass accretion rate onto the shock (Kitaura et al. 2006).
Recently, the progenitor dependence of the supernova dynamics has attracted great attention. Ugliano et al. (2012) performed systematic 1D spherically symmetric simulations for 101 progenitor models from Woosley et al. (2002) with parameterized neutrino luminosities and demonstrated that the explosion energy and characteristics of remnant compact objects (NSs and black holes [BHs]) strongly depend on the initial progenitor structure (see also Ertl et al. 2015). Recently, this study was extended by Nakamura et al. (2014) to 2D self-consistent simulations for the same 101 progenitor models. O'Connor & Ott (2013) performed a similar systematic study based on 32 progenitor models from Woosley & Heger (2007), focusing on the compactness parameter, i.e., the ratio of the mass to radius at a certain mass coordinate, at the bounce. They found that not the zero-age main-sequence (ZAMS) mass but the compactness parameter is a good measure for the neutrino evolution and hence for shock revival in the pre-explosion phase. Couch & Ott (2013) pointed out that, in addition to the density structure, velocity fluctuations in progenitors affect the hydrodynamics of shock revival. They showed in fact that models with velocity fluctuations imposed before collapse can explode more easily than those without them. The conclusion was confirmed later by Müller & Janka (2015), in which more systematic parametric studies were performed. Hence, it is becoming a consensus in society that the initial condition in progenitors is an important ingredient in supernova dynamics.
In this paper, we perform a series of neutrino-radiation hydrodynamic simulations in both spherical symmetry (1D) and axial symmetry (2D) for progenitors with a mass range from 12 to 100 M⊙ in the main-sequence phase, and we pay attention to the postbounce evolutions of mass accretion rates and neutrino luminosities, particularly how they are affected by the progenitor structure. Introducing the turning point on the trajectory in the plane, we argue that its location in the plane will serve as a sufficient condition for shock revival. We also construct a phenomenological model to understand how and when the turning point appears. Applied to a large number of progenitors, the model will be also useful to judge from the progenitor structure alone which progenitors are more likely to produce explosions than others.
The paper begins with the descriptions about the numerical simulations in Section 2. Then, we introduce a new concept, i.e., the turning point on the trajectory in the plane, in Section 3, and we discuss a sufficient condition for shock revival based on the location of the turning point relative to the critical curve. A phenomenological model, which estimates the mass accretion rate and neutrino luminosity from the progenitor structure and gives the location of the turning point qualitatively well, is presented in Section 4. We summarize our results and discuss their implications in Section 5.
2. NUMERICAL SIMULATIONS
2.1. Methods
The numerical methods are basically the same that we used in our previous studies (Suwa et al. 2010, 2011, 2013; Suwa 2014). With the ZEUS-2D code (Stone & Norman 1992) as a base for the hydrodynamics solver, we employ the equation of state of Lattimer & Swesty (1991) with the incompressibility K = 220 MeV, for which the maximum mass of a cold NS is 2.04 M⊙, i.e., more massive than the mass of recently discovered massive NSs (Demorest et al. 2010; Antoniadis et al. 2013). We solve the neutrino transfer equation for νe and by the isotropic diffusion source approximation (IDSA) scheme (Liebendörfer et al. 2009) that splits the neutrino distribution function into two components, which are solved with different numerical techniques. The weak interaction rates for neutrinos are calculated according to Bruenn (1985). The simulations are performed on a grid of 300 logarithmically spaced radial zones extending up to 5000 km, with the smallest grid width being 1 km at the center and 128 equidistant angular zones covering 0 ≤ θ ≤ π for two-dimensional (2D) simulations. For neutrino transport, we use 20 logarithmically spaced energy bins ranging from 3 to 300 MeV.
We conducted 2D simulations in this paper even though fully 3D computations with a spectral neutrino transfer are now becoming possible (Takiwaki et al. 2012; Hanke et al. 2013; Kuroda et al. 2015; Lentz et al. 2015; Müller 2015). The main reason for using 2D simulations is that 3D simulations are computationally too costly and are not suitable for systematic studies. It is true that the nature of turbulence is different between 2D and 3D, but its impact on the critical curve has been disputed by some authors (Murphy & Burrows 2008; Nordhaus et al. 2010; Hanke et al. 2012; Couch 2013b): 1D certainly gives the highest critical luminosity, but the difference between 2D and 3D is subtle. It should be noted that the critical curve seems just shifted vertically, with the shape being almost unchanged as the dimension changes from 1D to 2D to 3D, and we hence believe that the dimensionality will not be a critical factor for the following arguments in this paper. The treatment of neutrino reactions is not as sophisticated as in other simulations (e.g., Müller et al. 2012a), and heavy neutrinos are neglected in this paper. It is stressed, however, that our code is calibrated so that some key quantities, such as luminosities of electron-type neutrinos and anti-neutrinos and shock radii, could be reproduced reasonably well in 1D (see next subsection).
We also performed two additional simulations for model s55 employed in this paper: one with a resolution twice as fine, and the other with the radius of the outer boundary being doubled. In the former computation we confirmed no remarkable change, particularly in neutrino luminosities, whereas in the latter we found that the mass accretion rate was slightly affected at very late times but the shock evolution was essentially intact. We hence conclude that the numerical grid adopted in this paper is adequate at least for the purpose in this paper.
2.2. Code Check
In this subsection, we validate our code by comparison with model N13 in Liebendörfer et al. (2005), in which the authors performed spherically symmetric simulations with weak interactions for νe and alone being implemented as described in Bruenn (1985), and they demonstrated that two different codes (AGILE-BOLTZTRAN and VERTEX) produced consistent numerical solutions. For comparison, we also conduct a spherically symmetric simulation for the same 13 M⊙ model by Nomoto & Hashimoto (1988).
Figure 1 presents the shock radii (left panel) and luminosities of νe and (right panel) as a function of time. Although small differences can be observed, our results are consistent with the other two. The slightly smaller shock radii obtained in our simulation should give more conservative predictions for shock revival.
Figure 1. Comparison of shock radii (left panel) and neutrino luminosities (right panel) of ZEUS (code used in this work; red), AGILE (Liebendörfer et al. 2001; green), and VERTEX (Rampp & Janka 2000; blue). Luminosities of νe are represented by thick lines, and those of are represented by thin lines.
Download figure:
Standard image High-resolution image2.3. Progenitor Structures
We employ progenitors with solar metallicity calculated by Woosley & Heger (2007), who performed stellar evolutionary calculations. They have 12, 15, 20, 30, 40, 50, 55, 80, and 100 M⊙ at ZAMS. Some relevant quantities are presented in Table 1.
Table 1. Properties of Investigated Progenitors
Model | MZAMS | Final Mass | Final Radius | MFe | RFe |
---|---|---|---|---|---|
(M⊙) | (M⊙) | (R⊙) | (M⊙) | (1000 km) | |
s12 | 12 | 10.91 | 638.41 | 1.285 | 1.061 |
s15 | 15 | 12.79 | 831.04 | 1.346 | 1.172 |
s20 | 20 | 15.93 | 1066.68 | 1.540 | 1.591 |
s30 | 30 | 13.89 | 1552.89 | 1.476 | 1.448 |
s40 | 40 | 15.34 | 11.80 | 1.804 | 2.123 |
s50 | 50 | 9.82 | 5.42 | 1.487 | 1.489 |
s55 | 55 | 9.38 | 0.70 | 1.453 | 1.412 |
s80 | 80 | 6.37 | 0.60 | 1.479 | 1.501 |
s100 | 100 | 6.04 | 0.55 | 1.452 | 1.402 |
Download table as: ASCIITypeset image
First, we present the structures of these models. The top two panels in Figure 2 exhibit the density structures as functions of the radius (panel (a)) and enclosed mass (panel (b)). Panel (c) shows the mass–radius relation, in which the free-fall timescales (, where r is the radius, G is the gravitational constant, and M is the enclosed mass) are plotted as dashed lines. One can find that the density structure and ZAMS mass do not correlate with each other in a simple way: the density at the enclosed mass of 2 M⊙ (see panel (b)) becomes the smallest for the model with a ZAMS mass of 12 M⊙ and attains the maximum at 40 M⊙. Models with the ZAMS masses larger than 40 M⊙ have densities in between. This is because strong mass loss during the main-sequence and giant phases yields smaller cores (see also Table 1). From dashed lines in panel (c), one can easily see that the difference of structure leads to different free-fall times of mass elements, which will then result in different mass accretion histories.
Figure 2. Stellar structures for investigated models. The top two panels display the densities as a function of (a) radius and (b) enclosed mass. Panel (c) gives the radii corresponding to the mass and radius relations. The dashed lines show the free-fall times of 0.01, 0.1, and 1 s from bottom to top. Refer to text for details.
Download figure:
Standard image High-resolution imageIn Figure 3, we show the compactness parameter defined by O'Connor & Ott (2011) as
where R(M) denotes the radius for the enclosed mass M. It is evident that the 12 M⊙ model has the smallest ξM at all enclosed masses. ξM increases with the progenitor mass up to 40 M⊙. Interestingly, the models with 50, 55, 80, and 100 M⊙ have smaller ξM than the model with 40 M⊙, which is consistent with the results of O'Connor & Ott (2013), in which they showed that the model with 40 M⊙ gives the maximum values for both ξ1.75 and ξ2.5. In this sense the model s40 is the most compact progenitor, whereas the model s12 is the least compact one. Recently, the temporal evolutions of the ξ parameter during the stellar evolution were studied in detail by Sukhbold & Woosley (2014), and it was shown that ξ is indeed an appropriate quantity to characterize the progenitor structure.
Figure 3. Compactness parameters ξM defined in Equation (1) as a function of mass coordinate M. A larger ξM means a more compact structure: s12 is the least compact progenitor, while s40 is the most compact.
Download figure:
Standard image High-resolution imageThe subsequent subsections present our numerical simulations in 1D and 2D consecutively.
2.4. Spherically Symmetric Simulations
In Figure 4, we show the time evolutions of shock radius (panel (a)) and mass accretion rate at 300 km (panel (b)). It can be seen that all simulations failed to explode owing to insufficient neutrino heating just as expected. We especially pay attention to the evolutions of shock radius and the mass accretion rate, which are intimately related with each other.
Figure 4. Time evolutions of (a) shock radius and (b) mass accretion rate. There are bumps in panel (a), which correspond to the rapid decreases of mass accretion rate (see panel (b)).
Download figure:
Standard image High-resolution imagePanel (b) of Figure 4 shows that the mass accretion rate decreases rapidly until ∼200 ms and then becomes almost constant thereafter for a majority of models, including s20, s30, s50, s55, and s80. This transition in the mass accretion rate originates from the change of accreting layers, i.e., from the silicon layer to the oxygen layer: normally, a large density jump at the boundary of these layers in progenitors can be observed (see Woosley & Heger 2007). As a consequence of this transition, the stalled shock wave expands for a while around the sudden change in the mass accretion rate. We can also recognize some exceptions: model s12 shows a similar change in the mass accretion rate but at a much later time ∼500 ms with a smaller mass accretion rate; in model s100 the mass accretion rate is not settled to a constant value after the change at ∼350 ms but continues to decrease thereafter; in model s15 there is no clear change discernible at all. In the former two cases, the shock expands slightly or stops receding for a while around the times of the changes in the mass accretion rate like the first group. Model s40 is another outlier; although it has a similar pattern in the temporal evolution of the mass accretion rate to those found for the majority of progenitors, the transition occurs at a much later time ∼500 ms with a much high accretion rate. As a result of this high mass accretion and because of a larger ram pressure, the temporary expansion of the stalled shock wave associated with the transition is much less remarkable.
In Figure 5, we show the time evolution of the neutrino luminosity Lν (thick lines) and rms energy (thin lines), which is determined by , where ν is the neutrino energy, μ is the cosine of the angle between radial and neutrino propagation directions, and fν is the distribution function of neutrinos, for νe (solid lines) and (dashed lines). Note that the neutrino heating rate can be written as . In spite of the difference of mass accretion history (see Figure 4(b)), the neutrino luminosities and rms energies are not very different. Among them, model s12 and model s40 exhibit the smallest and the highest values, respectively.
Figure 5. Neutrino luminosities (thick lines) and rms energy (thin lines). Solid lines represent νe, and dashed lines give .
Download figure:
Standard image High-resolution imageIn Figure 6, we show the same evolutions in the plane, where is the mass accretion rate at 300 km and Lν is the total neutrino luminosity (i.e., the sum of the contributions from νe and ). Each model moves from right (high accretion rates) to left (low accretion rates) on these curves. As the mass accretion rate decreases, the neutrino luminosity also diminishes in general, except possibly for the very early phase. One can recognize the steepening of these curves near their left ends for the same majority group of progenitors. This is yet another manifestation of the transition in the mass accretion rate in this plane: the mass accretion rate becomes almost constant after the transition, while the neutrino luminosity continues to decrease. The position of this transition point in the plane is important for shock revival, which will be discussed in detail in the next section. If the point is located more to the top left corner in this plane (i.e., having lower mass accretion rates and higher neutrino luminosities), such a model will be more likely to produce an explosion, particularly in multidimensional simulations, in which the critical curve is supposed to run lower than in 1D. Models s55 and s80 are hence good candidates for exploding models in 2D in this sense.
Figure 6. Model trajectories in the plane. The mass accretion rate is evaluated at 300 km from the center.
Download figure:
Standard image High-resolution imageFigure 7 is the same as Figure 6 but for selected models, i.e., s12, s15, s20, and s55, with some time stamps. The other models not shown in Figure 7 have similar trajectories. Models s20 and s55, members of the majority group, have a clear transition point—which will be referred to as the turning point in the next section—at and , respectively. This is due to rather large jumps in density at the boundary between the silicon and oxygen layers in these models, which are apparent in Figure 2. As is also evident from the time stamps in Figure 7, these models move rapidly from right (high accretion rate) to left (low accretion rate) up to the turning point and then shift downward slowly later. They hence stay near the turning point for a long time, and, as argued later, that is the point where shock revival is most likely to occur. This is why we propose to employ the position of the tuning point as a diagnostic of explosion. It is also noted that the turning point is not always clearly visible and may not exist at all for some models. As mentioned earlier, model s12 does show a transition in the mass accretion rate in Figure 4(b), but the turning point is barely visible near the left end of the trajectory. Model s15 does not show any discernible transition already in Figure 4(b), which is also reflected in Figure 7. It is hence clear that the criterion for shock revival later proposed is a sufficient condition applicable only to those cases with the turning point.
Figure 7. Model trajectories in the plane for selected models with various points representing the times after bounce. The colors of lines are the same as those in Figure 6. Each model evolves from right (high accretion rate) to left (low accretion rate).
Download figure:
Standard image High-resolution image2.5. Axially Symmetric Simulations
In this subsection, we show the results of 2D simulations for the progenitors explored in the previous subsection in 1D. Since all 1D models fail to explode, these 2D simulations will serve as a guide to consider which progenitors are likely to explode.
Figure 8 gives the temporal evolutions of entropy at the north (upper panels) and south poles (lower panels). There are several oscillations in the shock radius for these models, which are the consequences of the standing accretion shock instability (SASI; Blondin et al. 2003; Ohnishi et al. 2006; Blondin & Mezzacappa 2007; Foglizzo et al. 2007; Iwakami et al. 2008; Fernández & Thompson 2009). It is clear that the material in the postshock region is heated up by neutrino irradiation from PNSs (the yellow color represents high entropies). Thanks to this heating, some models (s12, s40, s55, and s80) eventually produce shock re-expansion along the symmetry axis. The other models (s15, s20, s30, s50, and s100) yield no such expansion at least by the end of simulations, even though there is certainly neutrino heating in operation.
Figure 8. Time–space diagrams of specific entropy at poles for two-dimensional simulations. Upper (lower) panels represent the values at the north (south) pole. Models s12, s40, s55, and s80 eventually produce explosions at different times, depending on the initial density structures. The other progenitors, i.e., s15, s20, s30, s50, and s100, failed to produce an explosion at least by the end of simulations.
Download figure:
Standard image High-resolution imageFigure 9 presents time evolutions of the shock radius averaged over a solid angle. It can also be seen in this figure that several progenitors produce a shock revival, which is a necessary condition for a supernova explosion.10 It is important to note that these successful models have the turning points located either more to the left (s12; low accretion rate) or more to the top (s40; high neutrino luminosity) or both (s55, s80; low accretion rate and high luminosity) than unsuccessful models. We use this observation in the next section.
Figure 9. Time evolutions of the angle-averaged shock wave radius. Four of the investigated models, i.e., s12, s40, s55 and s80, clearly show shock expansions.
Download figure:
Standard image High-resolution imageIt seems that the onset of shock revival is delayed from the time of the turning point (see Figure 4(b)). This is because the mass accretion rates in Figure 4(b) are evaluated at 300 km from the center, and it takes some time until it influences the postshock dynamics. Note also that the development of shock oscillations needs some time. It should be mentioned, however, that s40 will probably fail to explode when we take into account general relativity, which is neglected in this paper. O'Connor & Ott (2011) observed in their 1.5D general relativistic simulation that the same progenitor formed a BH at around 550 ms after bounce (similar results were obtained by Sumiyoshi et al. [2006] and Fischer et al. [2009], but with different progenitors). Since this time of BH formation is much earlier than the shock revival time we found in s40, the progenitor leads most likely to a BH formation instead of the very late explosion observed here. Note, however, that s40 is an outlier anyway, having a very large compactness and a very late occurrence of the turning point.
Panel (a) of Figure 10 exhibits the abundance of 28Si (red line) and 16O (green line), as well as the density (blue line) of model s80. One can find that there are two density jumps at 1.66 and 2.17 M⊙ in mass coordinates. Panel (b) of this figure displays as gray lines the trajectories of mass shells at the mass coordinates of 1–1.85 M⊙ with an interval of 0.01 M⊙. Three thin black lines represent the mass coordinates of 1.66, 1.7, and 1.75 M⊙. Note that 1.66 M⊙ corresponds to the interface between Si- and oxygen-burning shells (see also panel (a)). It is interesting to see what happens when this mass shell accretes onto the shock (thick black line). It is evident that several oscillations ensue and the standing shock is finally converted to the expanding shock at ∼400 ms after the bounce. This is a clear demonstration that the transition in the mass accretion rate triggers shock revival.
Figure 10. (a) Initial profiles of density and composition for model s80. The abundances of 28Si (red line) and 16O (green line) and the density (blue line) are given as a function of the mass coordinate. There are two jumps in density, representing the transition of layers. (b) Trajectories of the mass shells with the mass coordinates of 1–1.85 M⊙ with an interval of 0.01 M⊙ are plotted as gray curves for the same model. Thin black lines represent 1.66, 1.7, and 1.75 M⊙ from left to right, respectively. A thick black curve indicates the average shock position. When the mass shell of 1.66 M⊙ runs across the shock, several oscillations ensue in the shock radius. The shock is eventually expanded at ∼400 ms after the bounce.
Download figure:
Standard image High-resolution imageAlthough this is not relevant for the main focus of this paper, we show in Figure 11 for reference the so-called diagnostic energy, which is defined as the integral of the sum of specific internal, kinetic, and gravitational energies over all zones, in which it is positive. Four exploding models (s12, s40, s55, and s80) have indeed nonvanishing diagnostic energies. Some oscillations originate from the shock oscillations. Though the diagnostic energy is gradually increasing, the final value is still much smaller than the typical value of the observed explosion energy, ∼1051 erg. Although even nonexploding models have positive diagnostic energies owing to neutrino heating, it is insufficient to revive the stalled shock wave.
Figure 11. Time evolutions of the diagnostic energy for 2D models. It is defined by the integral of the sum of the specific internal, kinetic, and gravitational energies, over the zones, in which it is positive. The horizontal axis is the postbounce time.
Download figure:
Standard image High-resolution image2.6. Different 15 M⊙ Models
It is a well-known unfortunate fact that stellar evolution calculations by different groups do not agree with one another. In this subsection, we investigate this issue, employing different progenitors with the same typical mass of 15 M⊙ at ZAMS. In addition to model s15 just studied, we use four more models from Nomoto & Hashimoto (1988) (NH88), Woosley & Weaver (1995) (WW95), Woosley et al. (2002) (WHW02), and Limongi & Chieffi (2006) (LC06). The first three of them were also employed in Suwa et al. (2011), in which neutrino oscillation effects on a supernova explosion were investigated. The precollapse density structures are given in Figure 12 (see also Figure 8 of Suwa et al. 2011 for comparison of the density structures at 100 ms after the bounce; in this paper it was argued that the structures are similar among the different models for M < 0.8 M⊙ whereas they are different for M > 0.8 M⊙). It can be observed that even though the initial mass at ZAMS is the same, the density structures prior to collapse become different, depending on both the physics and the numerics implemented in stellar evolutionary calculations. It should be noted in particular that the difference between WW95 and WH07 is substantial for M ≳ 1.1 M⊙ before collapse (compare red and orange lines in Figure 12).
Figure 12. Same as Figure 2, but for progenitors with the ZAMS mass of 15 M⊙. Here we use five models from Nomoto & Hashimoto (1988) (NH88), Woosley & Weaver (1995) (WW95), Woosley et al. (2002) (WHW02), Limongi & Chieffi (2006) (LC06), and Woosley & Heger (2007) (WH07). Owing to the different treatments of physics and numerics for stellar evolutionary calculations, the structures prior to collapse show diversity even if they have the same ZAMS mass. In panel (b), free-fall times are given by dashed lines.
Download figure:
Standard image High-resolution imageFigure 13 presents these models in the plane evaluated for 1D simulations (cf. Figure 6). NH88, WW95, and LC06 have clear turning points, and the former two are located more to the left than the last and are more likely to achieve shock revival. This is a consequence of the density jumps more remarkable for these models as observed in Figure 12. It is noted that all 1D simulations failed to produce an explosion.
Figure 13. Model trajectories in the plane for the 1D simulations of 15 M⊙ progenitors. This is the same as Figure 6, but for different progenitor models. The mass accretion rate is evaluated at 300 km from the center.
Download figure:
Standard image High-resolution imageThe shock evolutions for 2D simulations are given in Figure 14. The two progenitors, NH88 and WW95, indeed succeeded in producing shock revival, whereas the others failed. This is a clear demonstration that not the ZAMS mass but the density structure of the progenitor matters for the dynamics of shock revival. Again, the successful models have turning points that are located more to the left than the unsuccessful models, as seen in Figure 13. This is the same conclusion as in the previous subsection.
Figure 14. Time evolutions of the angle-averaged shock radius for 15 M⊙ progenitors. NH88 and WW95 produce explosion owing to small densities of the envelopes.
Download figure:
Standard image High-resolution image3. TURNING POINT
In this section, we propose a novel idea to diagnose a possibility of shock revival using the trajectory in the or plane (see Figure 15). This plane is often used to discuss the critical curve, which divides this plane into two regions: the region below this line, in which there are steady accretion flows, and the other region above the curve, in which there is no such flow (Burrows & Goshy 1993). The latter is therefore interpreted as the region where shock revival occurs. The question arises, where on the actual trajectory is the critical line crossed from below?
Figure 15. Schematic picture of the critical curve and turning point. If the turning point is located above the critical curve and the luminosity and mass accretion rate stay in the vicinity of the tuning point for a long time, such a model will produce an explosion. The critical curve is expected to be shifted by macrophysics such as dimensionality, and the turning point may be shifted by microphysics, as well as the progenitor structure. The critical curve and turning point are also useful to assess the influence of a particular physics incorporated.
Download figure:
Standard image High-resolution imageIn Figure 15, we present the typical situation we found in the majority of our models in the previous sections as a schematic picture of the trajectory and the critical curve in the plane. The red solid line represents the critical curve, and the black dotted line gives a typical trajectory. As mentioned already in the preceding sections, there is a point on the trajectory at which the slope of the trajectory steepens suddenly as a consequence of the rapid change in the mass accretion rate there. This point is referred to as the turning point in this paper. It is worth noting that the trajectory is shallower than the critical curve before the turning is reached and the order is changed thereafter. Consequently, it is obvious that the trajectory can cross the critical curve if and only if the turning point is located above the critical curve. It should also be clear that shock revival will be fizzled if the system evolves rapidly after the turning point, rolling down the second half of the trajectory and quickly passing the critical point again. Hence, it is important that the system stays for a long time around the turning point.
Since the critical curve is a convex and the monotonically increasing function of the mass accretion rate, the more to the upper left the turning point is located, the more likely shock revival is to obtain. Although the critical curve has been well studied by several groups,11 we emphasize here the importance of the trajectory as well. In principle, multidimensional neutrino-radiation hydrodynamic simulations, or ab initio computations, with detailed neutrino physics and radiative transfer being incorporated, are required to obtain reliable model trajectories. It has been demonstrated, however, that one observed effect of multidimensionality in supernova dynamics is to lower the critical curve (Murphy & Burrows 2008; Nordhaus et al. 2010; Hanke et al. 2012), although the trajectory is also somewhat modified. Hence, it is expected that 1D simulations will be sufficient to find approximate locations of turning points and to infer which models are more likely to explode than others. A 1D model trajectory will also be useful to discuss to what extent particular ingredients included in simulations (e.g., the nuclear equation of state, neutrino interactions, scheme to solve the neutrino transfer) affect the location of the turning point.
In the following, based on the results of our simulations presented so far, we develop a phenomenological model that connects the density structure of the progenitor just prior to the collapse and the model trajectory in the plane.
4. PHENOMENOLOGICAL MODEL
In this section, we construct a phenomenological model. The purpose is twofold: first, it is important to understand qualitatively why and where the turning point appears; second, we aim to expedite the judgment of which progenitors are likely to produce shock revival. As mentioned in the previous sections, the location of the turning point, if any, on the trajectory in the plane may serve as a sufficient condition for shock revival if it is located more to the upper left corner. Although the trajectory evaluated in 1D simulations will be sufficient for this purpose, the procedure may be simplified even further by the employment of the phenomenological model. It is not necessary for the phenomenological models to perfectly reproduce the trajectories obtained numerically. Instead, it is important that the turning points are placed at approximately correct positions and that the relative locations of the turning points for different progenitors are correctly reproduced. The latter point is particularly important, since the numerical results contain systematic errors one way or another.
In this model the mass accretion rate is evaluated as
where tff is the free-fall time, which is defined as a function of the radius by
where α is a parameter introduced to fit to numerical results. Inverting this relation, we regard the radius as a function of tff. Figure 16 shows the mass accretion rates as a function of t, which is identified with tff. The figure should be compared with panel (b) of Figure 4 (note that the vertical scale is different). It can be seen that the model reproduces the characteristic features, namely, that the mass accretion rate is high and rapidly decreasing initially, and when the silicon layer accretes onto the PNS completely, it becomes significantly smaller because of the density drop at the layer boundary and remains almost constant thereafter (see Appendix
Figure 16. Mass accretion rate calculated by the free-falling model.
Download figure:
Standard image High-resolution imageAlthough there have been several approximate functional forms, e.g., (Janka & Mueller 1996), proposed for the total neutrino luminosity as a function of time, we employ the following form based on the diffusion timescale:
where Ldiff = Eint/tdiff is the diffusion luminosity, with being the internal energy stored inside the PNS and tdiff being the diffusion timescale defined shortly, and MPNS and Rν are the mass and radius of the PNS, respectively. Again identifying tff with t results in the following expression:
The diffusion time tdiff can be evaluated as (see Appendix
where σ is the cross section of neutrino-nucleon scattering, which is given as , with σ0 = 1.705 × 10−44 cm2,12 the electron mass me, and the neutrino energy εν; the proton mass is denoted by mp. Moreover, is a characteristic energy of neutrinos inside the PNS.13 The mass of the PNS increases as matter accretes and can be expressed as a function of time by the use of the free-fall time.
The idea underlying Equation (6) is that the material initially located at r falls onto the PNS in its free-fall time and the gravitational energy is converted to the internal energy, which is finally radiated as neutrinos in the diffusion time. We can then evaluate the neutrino luminosity as Equation (6). This phenomenological model is consistent with the finding by Fischer et al. (2009) and Müller & Janka (2014) that the neutrino luminosity seems to be regulated by the smaller of the accretion and diffusion luminosities. Indeed, since and , for , and vice versa.
Figure 17 presents the model trajectories in the plane obtained this way. In this plot, we employ α = 1.5, , and Rν = 50 km in Equations (4) and (7). The turning points, where the slope of trajectory changes rapidly, are clearly produced in most cases. The comparison between the phenomenological model and the numerical results is given in Appendix
should be maximum. Note that, strictly speaking, Equation (8) is valid only for the luminosity of electron-type neutrino, , with the temperature of . Here the Boltzmann constant is denoted by k. We believe, however, that other expressions will not change the following discussions. It can be found that the filled squares coincide with the turning points, whenever they exist. For model s55, for example, this point occurs at s−1, which is consistent with the numerical result (see Figure 6), and Lν is a rapidly increasing function of for smaller mass accretion rates, whereas it increases rather slowly for larger mass accretion rates. As mentioned before, this drastic change comes from the transition of accreting layers, i.e., from silicon to oxygen layers. A comparison between the critical curve of Burrows & Goshy (1993) and turning points presented here is given in Appendix
Figure 17. Phenomenological model for neutrino luminosity as a function of mass accretion rate.
Download figure:
Standard image High-resolution imageGuided by a simple analytical derivation of the critical curve from Janka (2012) (Section 4.3.4 in this review), the critical curve depends not only on but also on NS mass MNS as .14 Note that in this formula MNS denotes the point mass determining the gravitational potential so that it should be replaced with enclosed mass for each mass element. In Figure 18(a), we show the distribution of the ratio between total neutrino luminosity calculated by the phenomenological model (Equation (6)) and the critical curve by Janka (2012) as
The overall factor is chosen to divide exploding and nonexploding models as follows. This quantity, Lν/Lcrit, represents how long model trajectories are distant from the critical curve, which is useful to predict model exploitability. In this figure, exploding models (s12, s40, s55, and s80) are indicated with thick lines, and failing models (s15, s20, s30, s50, and s100) are indicated with thin lines. One can find that three exploding models (s12, s55, and s80) exceed unity and other models stay below unity for the whole regime. The exceptional behavior of s40 originates from omission of neutrino energy dependence. As indicated by more realistic expressions of the critical curve (e.g., Pejcha & Thompson 2012; Müller & Janka 2015), the higher rms energy leads to smaller critical luminosity. The rms energy of s40 is typically higher than in other models, namely, by ∼10% so that its critical luminosity is smaller, by ∼20%. Since the peak value of of s40 is about ∼0.8, the inclusion of rms energy really improves the situation. Note that the rms energy of the other models is very similar (see Figure 5), so that the critical curve of other models does not change considerably. The modeling of rms energy is currently not feasible. However, most models (except a model with a significantly low or high mass accretion rate) imply similar rms energies (see also Liebendörfer et al. 2003), and their dependence on the progenitor structure is lower than luminosities. Thus, the critical curve given by Equation (9) works quite well to predict the explosion/failure from the progenitor model alone. Note also that these peaks in Figure 18(a) coincide with the position of turning points, even though we include the term of enclosed mass in the expression of the critical curve (Equation (9)).
Figure 18. (a) Ratio between model trajectories and the critical curve (Equation (9)) as a function of mass coordinates for each progenitor models; (b) locations of the turning points in the plane. These panels are based on the phenomenological model. In panel (b), additional models (s25, s26, s27, s28, and s29) presented in Appendix
Download figure:
Standard image High-resolution imageAs described in Section 3, a comparison between the turning points and the critical curve gives a useful criterion to judge which progenitors are likely to produce shock revival, which is shown in Figure 18(b). The horizontal axis is replaced from to following Janka (2012). The exploding models are denoted by open circles, and the failed models are denoted by crosses. The comparison works quite well for most of the models, with only the exception of s40, as discussed above. A similar plot is shown in Figure 21 in Appendix
It is true that the trajectories and the critical curve we obtained in this paper will change when the numerics are improved one way or another, but it should be stressed that the following methodology should be applied to any combination of a numerical code and to a collection of progenitors: we first perform a small number of 1D simulations to fix the parameters in the phenomenology and then apply the latter to a large number of progenitors to predict which ones are more likely to explode than others. As an example, we presented the locations of turning points obtained that way for other progenitor models from Woosley & Heger (2007) in Appendix
5. SUMMARY AND DISCUSSIONS
In this paper, we performed neutrino-radiation hydrodynamic simulations in spherical symmetry (1D) and in axial symmetry (2D) for different progenitor models by Woosley & Heger (2007) from 12 to 100 M⊙. We found that all 1D runs failed to produce an explosion and several 2D runs succeeded. The difference in the shock evolutions can be mainly ascribed to different mass accretion histories, which are determined by the density structures of progenitors. For the majority of the models we studied in this paper we found that the mass accretion rate changes its nature suddenly at a certain point: in the earlier phase it is very high and decreasing quickly, whereas it settles to a nearly constant rate in the later phase. In the plane this transition point, which we called the turning point, marks the point where the slope of the trajectory changes rapidly: the trajectory is shallower than the critical curve in the earlier phase and becomes steeper than that thereafter. This means that the trajectory crosses the critical curve from below around the turning point if it ever occurs and that shock revival is most likely to take place there. It is hence obvious that we can employ the position of the turning point as a diagnostic for shock revival. It should be noted, however, that the turning point does not always seem to exist and that the above criterion should be regarded as a sufficient condition. In addition, we developed a phenomenological model to approximately estimate trajectories in the plane. This model utilizes the initial density structure of the progenitor alone and reproduces the locations of turning points reasonably well. Based on these results, we suggest the following usage of the phenomenological model: perform a small number of 1D simulations to fix the free parameters in the model and apply the result to a large number of progenitors to infer which progenitors are more likely to explode. It should be noted that the main effect of the intrinsic multidimensionality of supernova dynamics is to shift the critical curve parallelly downward, and its influences on the trajectory are limited.
The phenomenological model depends on the underlying 1D simulations and hence on the numerics and input microphysics adopted therein. The important thing, however, is that the methodology is applicable to any combination of a numerical code and a collection of progenitor models. Some comments on the simulations by other groups follow. Bruenn et al. (2013, 2014) reported successful explosions for 12, 15, 20, and 25 M⊙ progenitor models of Woosley & Heger (2007), which are identical to those employed in this paper. Dolence et al. (2015) adopted the same progenitors in their simulations and found no explosion. Both of them employed the flux-limited diffusion approximation for neutrino transfer. Dolence et al. (2015) suggested that the reason for this discrepancy might be the use of the ray-by-ray approximation by Bruenn et al. (2013). We obtained an explosion for the 12 M⊙ model, but not for the 15 and 20 M⊙ models, on the other hand, which indicates that our simulations fall somewhere between them. Note that we also employed the ray-by-ray approximation.
Figure 19 compares the trajectories in the plane for these models.15 The solid lines are the trajectories for the models in this paper, the dots are obtained from Figure 7 of Dolence et al. (2015), the dashed lines are drawn based on Figure 1 of Bruenn et al. (2014), and the dot-dashed lines (only for s20) are generated from Figures 3 and 4 of Melson et al. (2015a). Our trajectories are closer to those of Dolence et al. (2015) and Melson et al. (2015a), while the trajectories of Bruenn et al. (2014) are apparently running in the higher-luminosity regime, which may have led to the successful explosions in their simulations. Note that Melson et al. (2015a) also obtained an explosion for s20, which might also be a consequence of higher neutrino luminosity at the turning point. However, the position of the critical curve depends on physical ingredients implemented in each code.
Figure 19. Comparison of model trajectories of this study (solid lines), Dolence et al. (2015) (dots; taken from Figure 7 of their paper), Bruenn et al. (2014) (dashed lines; made from Figure 1 of their paper), and Melson et al. (2015a) (dot-dashed lines, only for s20; made from Figures 3 and 4 of their paper). Note that trajectories of Bruenn et al. (2014) are plotted up to 150 ms after the bounce, while those of this study, Dolence et al. (2015), and Melson et al. (2015a) are plotted up to several hundreds of milliseconds postbounce.
Download figure:
Standard image High-resolution imageRecently, Nakamura et al. (2014) performed 2D simulations for a large number of progenitor models of 10.8–75.0 M⊙ taken from Woosley et al. (2002), which are incidentally different from our progenitor models, and obtained explosions for all the models. Although their code, employing the IDSA and ray-by-ray approximation for neutrino transfer, is quite similar to ours, it is not the same. This is because their hydrodynamics module employs a different shock-capturing scheme, and that may be the source of differences in the outcomes.
Finally, we comment on the assumptions adopted in this study. First, we performed 2D simulations although it is well known that axial symmetry leads to some hydrodynamic features that are qualitatively different from those in three dimensions (3D; Couch 2013b; Hanke et al. 2013; Handy et al. 2014; Takiwaki et al. 2014). The critical curve in 3D is observed to shift from that in 2D (Nordhaus et al. 2010; Hanke et al. 2012), although the magnitude and its dependence on the spatial resolution are still controversial. Incidentally, 3D hydrodynamic simulations with spectral neutrino transfer are currently still computationally too expensive to perform as a systematic study like the one in this paper. Second, the microphysics used in this study is not as elaborate as other numerical studies (Müller et al. 2012a; Bruenn et al. 2013), and our critical curve may be different from theirs. The phenomenological model we proposed in this paper should hence be applied to more sophisticated numerical simulations, which will be carried out as a future project.
We thank A. Heger, M. Limongi, and K. Nomoto for providing progenitor models, M. Tanaka for fruitful discussions, and S. Veith for proofreading. The numerical computations in this study were partly carried out on XT4 and XC30 at CfCA in NAOJ and SR16000 at YITP in Kyoto University. This study was supported in part by the Grant-in-Aid for Scientific Research (Nos. 25103511, 26870823, 23540323, 23340069, 24103006, 26707013, and 24244036), JSPS postdoctoral fellowships for research abroad, MEXT SPIRE, and JICFuS.
APPENDIX A: CONDITIONS FOR CONSTANT
In this section, we give a simple explanation of why we obtain almost constant mass accretion rates at late times for most of the progenitors in this study. We assume the following density structure:
where R is a core radius and ρ0 is the density at r = R. The mass coordinate is given by
where M0 is the mass coordinate at r = R. Then, the mass accretion rate is estimated as
where is the free-fall timescale. Short calculations give
and
Suppose that M(r) ≈ M0, i.e., the central accelerator's mass is dominant and . Then, , and we obtain
which becomes constant if n = 3/2. If the mass at r > R, that is, the mass of accreting matter, is dominant, we get
which leads to
Then we obtain
which is again constant for n = 2. The progenitor models used in this study realize n ≈ 2 for the oxygen layer, so that the latter case is valid for them.
APPENDIX B: DIFFUSION TIMESCALE
Here we obtain a useful expression of the diffusion timescale for neutrinos in a uniform density sphere of radius Rν, which is meant to be a rough approximation to a PNS. The diffusion timescale is given by
where τν is the optical depth of the sphere, which is
Here we used . By combining Equations (19) and (21), we get
APPENDIX C: COMPARISON BETWEEN PHENOMENOLOGICAL AND NUMERICAL MODELS
In this section, we show the comparison between the phenomenological model introduced in Section 4 and the numerical results presented in Section 2.4. Figure 20 presents the model trajectories for models s12, s20, and s80. Solid curves show the model trajectories obtained with the phenomenological models, and dashed curves display the trajectories given by the numerical simulations. It can be found that for s80 these lines agree very well, while for s12 the phenomenological model fails to reproduce the numerical result. As for s20, there is a discrepancy between two lines at high mass accretion rates, whereas the turning point is almost perfectly reproduced. For models s30, s50, and s55, which have clear turning points in the simulations, they are reproduced reasonably well by the phenomenological model. However, for models s15, s40, and s100, we cannot fit well. These results indicate that the phenomenological model is useful for progenitors that have clear turning points, i.e., progenitors with a large density jump between the silicon and oxygen layers.
Figure 20. Comparison between numerical simulation (dashed lines) and the phenomenological model (solid lines) for selected models.
Download figure:
Standard image High-resolution imageAPPENDIX D: TURNING POINTS IN PLANE AND CRITICAL CURVE
In Figure 21 the locations of the turning points obtained by the phenomenological model are plotted for different progenitors as open circles for exploding models and as crosses for nonexploding models in 2D simulations described in Section 2.5. The dashed line is given by , which we fit to roughly divide the exploding models from the nonexploding ones. These critical luminosities are slightly higher than those obtained by Murphy & Burrows (2008) with the light-bulb approximation: . The discrepancy should be ascribed to the difference in the numerical treatments of neutrino transfer. It should be noted, however, that this fit is not perfect and, for instance, s20 lies above the line. We stress that the location of the turning point is indeed a useful measure when one attempts to infer the possibility of explosion from the progenitor structure alone.
Figure 21. Locations of the turning points of the investigated models are represented as open circles for exploding models and crosses for nonexploding models in 2D simulations given in Section 2.5. The dashed line is , which seems to divide the phase space into exploding and nonexploding regions.
Download figure:
Standard image High-resolution imageAPPENDIX E: OTHER PROGENITORS
We show the turning points for all 32 progenitors from Woosley & Heger (2007) in Figure 22. The turning point is defined for each model to be the point in the plane, at which the ratio of takes the maximum value on the trajectory. The values of max() at the turning points are summarized in Table 2, in which the compactness parameters ξ1.5, ξ1.75, and ξ2.5 at the precollapse phase are also given (see Sukhbold & Woosley 2014 for the relations of these quantities with the compactness parameters defined at the bounce). Note that these parameters can be also evaluated only from the density structure of the progenitor and no simulation is required.
Figure 22. Location of turning points for all progenitor models of Woosley & Heger (2007), with the color bar denoting the ZAMS mass. The typical progenitor models are labeled.
Download figure:
Standard image High-resolution imageTable 2. Properties of All Progenitors
Model | MZAMS | Max() | Max() | |||
---|---|---|---|---|---|---|
(M⊙) | ||||||
s12 | 12 | 2.19 | 1.051 | 0.617 | 0.235 | 0.023 |
s13 | 13 | 2.36 | 1.064 | 0.869 | 0.370 | 0.067 |
s14 | 14 | 2.25 | 0.972 | 0.857 | 0.502 | 0.128 |
s15 | 15 | 2.06 | 0.858 | 0.882 | 0.549 | 0.181 |
s16 | 16 | 2.19 | 1.095 | 0.792 | 0.333 | 0.150 |
s17 | 17 | 2.24 | 1.093 | 0.877 | 0.374 | 0.168 |
s18 | 18 | 2.43 | 0.996 | 0.961 | 0.656 | 0.194 |
s19 | 19 | 2.15 | 1.006 | 0.962 | 0.517 | 0.177 |
s20 | 20 | 2.18 | 0.943 | 1.003 | 0.771 | 0.286 |
s21 | 21 | 2.07 | 1.060 | 0.696 | 0.323 | 0.143 |
s22 | 22 | 2.16 | 0.929 | 1.001 | 0.783 | 0.289 |
s23 | 23 | 2.14 | 0.828 | 0.998 | 0.870 | 0.434 |
s24 | 24 | 2.17 | 0.859 | 1.013 | 0.859 | 0.398 |
s25 | 25 | 2.25 | 0.938 | 1.008 | 0.821 | 0.331 |
s26 | 26 | 2.48 | 1.101 | 0.968 | 0.641 | 0.234 |
s27 | 27 | 2.37 | 1.048 | 0.993 | 0.677 | 0.257 |
s28 | 28 | 1.79 | 0.716 | 0.987 | 0.596 | 0.272 |
s29 | 29 | 2.17 | 1.018 | 0.965 | 0.534 | 0.225 |
s30 | 30 | 2.13 | 0.944 | 1.006 | 0.688 | 0.218 |
s31 | 31 | 2.17 | 0.983 | 0.995 | 0.617 | 0.219 |
s32 | 32 | 2.15 | 0.942 | 0.999 | 0.750 | 0.253 |
s33 | 33 | 2.14 | 0.915 | 1.003 | 0.783 | 0.284 |
s35 | 35 | 2.26 | 0.861 | 1.010 | 0.846 | 0.360 |
s40 | 40 | 2.62 | 0.736 | 0.980 | 0.876 | 0.547 |
s45 | 45 | 2.58 | 0.743 | 0.982 | 0.875 | 0.516 |
s50 | 50 | 2.12 | 0.958 | 0.999 | 0.643 | 0.222 |
s55 | 55 | 2.19 | 1.018 | 0.980 | 0.564 | 0.239 |
s60 | 60 | 2.24 | 1.066 | 0.939 | 0.451 | 0.175 |
s70 | 70 | 2.20 | 0.987 | 0.989 | 0.663 | 0.233 |
s80 | 80 | 2.19 | 1.014 | 0.965 | 0.550 | 0.210 |
s100 | 100 | 2.13 | 0.824 | 0.989 | 0.702 | 0.245 |
s120 | 120 | 2.09 | 0.842 | 0.911 | 0.454 | 0.171 |
Download table as: ASCIITypeset image
There is no clear correlation between the maximum values of and the compactness parameters. The critical value of max() that divides exploding from nonexploding models may be set at ∼2.18 because models s55 and s80 (max() = 2.19) explode, while s20 (max() = 2.18) fails. It is hence true that max() is very useful in judging whether a particular model is likely to explode before doing detailed simulations. In addition to max(), we show max() as well in this table. The explosion criterion of this indicator is unity. Therefore, 12 of the 32 models are expected to make explosions with our current numerical code.
To confirm robustness of the criterion, we perform additional simulations of different progenitor models. From Table 2 we pick up from s25 to s29, which include representative models giving maximum and minimum values for max(), i.e., s26 and s28. The resultant shock evolution of 2D simulations is shown in Figure 23, which displays the explosion of three models (s26, s27, and s29) and the failing of two others (s25 and s28). Although the first indicator, max(), does not perfectly fit this result since there is inversion hierarchy between s25 (2.25 and nonexploding) and s29 (2.17 and exploding), the second indicator, max(), solves this contradiction, i.e., s25 is below unity (0.938) and s29 is above unity (1.018).
Figure 23. Shock trajectories of models s25, s26, s27, s28, and s29. Average shock radius is measured in postbounce time.
Download figure:
Standard image High-resolution imageThe value of the second indicator for s29 implies that this model is rather marginal. To check this marginality, we performed additional simulations with different initial perturbations. In addition, we performed additional simulations for s26 as well, because it is expected to explode more robustly, as indicated by max() and max(). Figure 24 presents shock evolutions with different seed perturbations (0.1% in density). One can see that all s26 simulations explode at almost the same time, while half of the s29 simulations explode at very different times. This result implies that models with high values of indicators explode robustly and those with marginal values sometimes fail or explode at different times (see Horiuchi et al. 2014, for the explosion fraction). Note that the diversity of explosion time is narrowed in 3D simulations (Handy et al. 2014; Takiwaki et al. 2014).
Figure 24. Dependence of shock trajectories on initial seed perturbation. Black lines present s26, and gray lines indicate s29.
Download figure:
Standard image High-resolution imageFootnotes
- 7
See also Kotake et al. (2012), which includes the spherically symmetric simulations.
- 8
- 9
There have, of course, been many other important issues raised so far. The most high-profile one for the moment is dimensionality of hydrodynamics (e.g., Ohnishi et al. 2006; Murphy & Burrows 2008; Nordhaus et al. 2010; Burrows et al. 2012; Hanke et al. 2012; Couch 2013b; Takiwaki et al. 2014). See also Suwa et al. (2013) and Couch (2013a) for roles of the the nuclear equation of state in multidimensional hydrodynamic modeling. Influences of various neutrino interactions were also investigated in Suwa et al. (2011) and Müller et al. (2012a).
- 10
Note that a shock revival, or a re-expansion of the stalled shock wave, is not a sufficient condition for a supernova explosion. In fact, shock revival is just a consequence of the dominance of the postshock thermal pressure over the ram pressure in the preshocked region, and the mass accretion to PNSs may continue thereafter, increasing the mass of PNSs. In order to produce a successful explosion, the expanding shock should be strong enough to turn the accretion into an expansion of the envelope. See Suwa et al. (2013) for more details.
- 11
- 12
There are coefficients of O(1), which are neglected for simplicity (see Burrows et al. 2006b, for more details).
- 13
The characteristic value employed here seems rather large compared with the commonly used value ∼10 MeV. This is because the former represents the average energy inside the PNS, where the matter temperature is a few tens of MeV and, as a consequence, the neutrino average energy becomes several tens of MeV. On the other hand, the latter value reflects the matter temperature at the neutrinosphere, O(1) MeV.
- 14
- 15
Note that this is not in the plane.