Physics of CoreCollapse Supernovae in Three Dimensions:
a Sneak Preview
Abstract
Nonspherical mass motions are a generic feature of corecollapse supernovae, and hydrodynamic instabilities play a crucial role for the explosion mechanism. First successful neutrinodriven explosions could be obtained with selfconsistent, firstprinciple simulations in three spatial dimensions (3D). But 3D models tend to be less prone to explosion than corresponding axisymmetric (2D) ones. This has been explained by 3D turbulence leading to energy cascading from large to small spatial scales, inversely to the 2D case, thus disfavoring the growth of buoyant plumes on the largest scales. Unless the inertia to explode simply reflects a lack of sufficient resolution in relevant regions, it suggests that some important aspect may still be missing for robust and sufficiently energetic neutrinopowered explosions. Such deficits could be associated with progenitor properties like rotation, magnetic fields or precollapse perturbations, or with microphysics that could lead to an enhancement of neutrino heating behind the shock. 3D simulations have also revealed new phenomena that are not present in 2D, for example spiral modes of the standing accretion shock instability (SASI) and a stunning dipolar leptonemission selfsustained asymmetry (LESA). Both impose time and directiondependent variations on the detectable neutrino signal. The understanding of these effects and of their consequences is still in its infancy.
Contents
1 Introduction:
NeutrinoDriven Explosions From a
Broader Perspective
A large variety of observational aspects indicates the importance of multidimensional effects in corecollapse supernovae. Width and smoothness of the lightcurve maximum, Doppler shifting and broadening as well as substructures of spectral lines, early emergence of Xray and gammaray emission, and significant levels of polarization in the case of Supernova 1987A and other well observed supernovae testify the presence of largescale asphericity and radial mixing processes [for reviews, see 1, 2, 3]. Also the gaseous remnants of supernovae exhibit global deformation and inhomogeneities as morphological fingerprints of explosion asymmetries. Clumpiness, filaments, and directional variations of the elemental composition are interpreted as consequences and relics of explosion asymmetries whose origin could be linked to the earliest moments of the explosion [e.g., 4, 5, 6, 7, 8, 9]. Measured natal kick velocities of young neutron stars, which cannot be explained by the orbital velocities of disrupted binary systems, provide another empirical hint to considerable asymmetries that are present already in the earliest phase of the explosion [e.g., 10, 11, 12, 13].
Numerical simulations also demonstrate that multidimensionality plays a crucial role to explain why massive stars blow up [for recent reviews on this subject, see 14, 15, 16, 17, 18, 19]. Hydrodynamical instabilities in collapsing stellar cores like convection and the standing accretion shock instability [SASI; 20, 21, 22] have been shown to grow even from small initial perturbations on time scales relevant for the supernova mechanism. They trigger asymmetries and nonradial mass motions on large scales and create turbulent flows on small scales. Also rotation and magnetic fields can have a decisive influence on the development of the explosion, because the initial fields can be strongly amplified in the collapsing stellar core by compression, winding in shear layers, and in particular by the magnetorotational instability in differentially rotating regions [e.g., 23].
Multidimensional effects are not just a byproduct or sideeffect of the explosion, they are essential for the success of the mechanism. The presently known and accepted “standard” physics —i.e., precollapse conditions from stellar evolution calculations, transport opacities for active neutrinos, nuclear reaction rates and the neutronstar equation of state, relativistic plasma dynamics and general relativistic gravity— do not facilitate successful explosions in spherically symmetric (“onedimensional”, 1D) simulations except for stars close to the lowmass end of supernova progenitors [for details, see 14]. This conclusion was unavoidable after elaborate, energy and velocitydependent threeflavor neutrino transport based on the direct solution of the Boltzmann equation [24, 25] as well as iterative solvers of the twomoment equations (i.e. neutrino energy and momentum equations) coupled to a Boltzmann closure [26, 27] had become available and did not bring success to the models.
textbox[b] FLAVOR OSCILLATIONS OF ACTIVE NEUTRINOS Flavor oscillations of the three active neutrino flavors are not selfconsistently included in present supernova models. For solar and atmospheric mixing parameters, flavor changes induced by the matter background according to the MikheyevSmirnovWolfenstein [28, 29] effect take place in the stellar mantle and envelope (at densities 10 g cm) and are suppressed in the dense medium of the supernova core [30, 31]. The situation is less clear for collective flavor transformations caused by neutrinoneutrino interactions. In view of the still incomplete understanding of this extremely complex, highly nonlinear phenomenon, it cannot be excluded that neutrinos propagating through the dominant neutrino background outside of the neutrinosphere might change their flavor identity [for a status report, see 32]. Oscillations of active neutrinos outside of the neutrinosphere are important for the detectable neutrino signal and may have an impact on supernova nucleosynthesis, but it is unlikely that they have a strong influence on the explosion. Since in present supernova models the individual luminosities of muon and tau neutrinos and antineutrinos (, , , , collectively denoted as ) are roughly half of those of electron neutrinos () and antineutrinos () during the postbounce accretion phase, and since the heavylepton neutrino spectra are only moderately harder, even a complete swap of to (or of the antineutrinos) just outside of the neutrinosphere would not enhance the neutrinoenergy deposition behind the stalled shock front. \end@float
Solving the riddle how supernovae explode therefore requires a better understanding of multidimensional phenomena in the collapsing stellar core. Selfconsistent multidimensional numerical simulations are indispensable for that and have enabled considerable progress over the past decade. Supplemented by semianalytic analysis, toy models, parametric studies, and in the case of the SASI even by laboratory experiments [33, 19], such simulations have advanced the field, providing us deeper insights into the nonlinear processes that could aid the revival of the stalled bounce shock and shape the observable asymmetries of supernova explosions.
Of course, astrophysical phenomena in nature involve three spatial dimensions and the ultimate aim of firstprinciple modeling must therefore be 3D simulations. The increasing power of modern supercomputers and the development of highly parallelized numerical codes have made it possible only very recently to perform such 3D simulations with gridbased schemes including energydependent, threeflavor neutrino transport. This progress could thus be achieved more than a decade after the seminal work by C. Fryer and collaborators, who obtained neutrinodriven explosions in fullscale 3D stellar corecollapse models with smoothedparticle hydrodynamics (SPH) including gray, fluxlimited neutrino diffusion [34, 35, 36]. Their results basically confirmed similar, earlier simulations in two dimensions [37, 38] and represent a modern manifestation of the original concepts developed by Colgate & White [39], Arnett [40], and Bethe & Wilson [41]. While these spearheading computations of the first 3D explosions must unquestionably be considered as a breakthrough in supernova modeling, they still granted only a preliminary, diffuse glimpse into the 3D world. The radical simplifications of the neutrino transport and the severe deficiencies of resolution and accuracy in the hydrodynamics treatment did not allow for strong and reliable conclusions on the viability of the neutrinodriven mechanism.
In this review we report the advances of 3D supernova modeling since these early steps and reflect our growing understanding of the physics that determines the explosion mechanism. Because the prompt bounce shock fails to trigger explosions, the energy deposition by the enormously intense neutrino flux radiated from the newly formed neutron star appears as the most plausible explanation for the revival of the stalled supernova shock. Neutrino heating may provide the engine that powers the far majority of “normal” supernovae with energies from less than 10 erg to roughly erg [38, 42, 43, 44, 45].
textbox[b] ALTERNATIVES TO NEUTRINODRIVEN EXPLOSIONS In supernovae with energies higher than 2.5 erg up to the “hypernova” regime, where the explosions can reach several erg for progenitor stars above 20 , rapid rotation and the efficient amplification of magnetic fields are likely to play a crucial role. These explosions are probably driven by magnetohydrodynamic effects around highly magnetized, fastspinning neutron stars (possibly associated with the observed “magnetars”) or around rapidly rotating black holes that accrete infalling stellar gas with high rates from a surrounding torus threaded by strong magnetic fields (see Janka 2012 for a review). Though recent multidimensional models seem to strengthen the case, the paradigm of the neutrinodriven mechanism is far from being convincingly established or even proven. It should not remain unmentioned, that it is still questioned by some because of the remaining problems of selfconsistent models to yield robust explosions and to explain the observed energies of typical supernovae by neutrinoenergy deposition. While these problems may disappear once 3D models become more mature, one should certainly remain openminded as long as there is a lack of solid empirical evidence for the neutrino mechanism. However, the suggested alternatives, e.g. the “jitteringjet mechanism” [46, 47] and “collapseinduced thermonuclear explosions” [48], are based on adhoc assumptions in tension with the presently established understanding of stellar evolution and supernova dynamics. The involvement of speculative ingredients and the missing selfconsistency are neither satisfactory nor convincing and place such suggestions on a level of sophistication far below the current state of the neutrinodriven mechanism. \end@float
2 Advancing From Two to Three Dimensions
In this review we focus on the physics of the neutrinodriven mechanism and provide an update on the developments in 3D modeling after the report given by Janka [14]. Before doing so, this section will briefly summarize what we have learned from 2D simulations and why further progress needs 3D models.
2.1 Status of 2D Modeling
Most of the fullscale modeling of supernovae has been performed, and much of it is still performed, with the assumption of rotational symmetry around a chosen axis, i.e., in two spatial dimensions (2D). Such simulations with stateoftheart neutrino transport lend support to the viability of the neutrinodriven mechanism, although the present results of different groups reveal important and unsatisfactory quantitative and qualitative differences, and no general consensus has been achieved yet.
2.1.1 Results from Different Groups
2D simulations by the Oak Ridge group with “raybyrayplus” (RbR+) fluxlimited diffusion for 12, 15, 20, and 25 progenitors develop neutrinodriven explosions nearly simultaneously [49, 50]. The supernova shock expands essentially identically in all cases, and shock revival (measured by the time when the total energy of some mass in the gain layer becomes positive) occurs after only 200 ms of postbounce accretion. Because of the early onset of the explosions these models possess a big mass in the gain layer and therefore explode fairly strongly with energies between 0.34 and 0.88 B (1 B = 1 bethe = erg), which are still increasing at the end of the simulations. However, these “Series B” models seem to be too optimistic and were (modestly) affected by a numerical deficiency [51].
In contrast, for the same stars the Garching simulations, employing RbR+ twomoment transport with Boltzmann closure, yield later explosions. The onset times differ between the progenitors because the blast wave expands only when the massaccretion rate has dropped sufficiently [52]. The explosions are therefore likely to be less energetic, but the simulations had to be stopped due to small time steps long before the energies could saturate.
O’Connor and Couch [53], using a twomoment (M1) scheme with algebraic closure relation for 2D neutrino transport (but with a subset of neutrino interactions and without energybin coupling) found a dynamical behavior of the four progenitors qualitatively very similar to the results of the Garching group. Quantitative differences in details can probably be traced back to differences in the modeling ingredients, which demand systematic tests.
The Princeton group, applying 2D fluxlimited neutrino diffusion diffusion [54] and M1 transport [55], again for supernova runs of the same progenitors, did not obtain any explosions. However, they used a Newtonian potential, whereas all other groups made use of an approximation of relativistic gravity. The quantitative importance of relativistic effects had already been concluded from 1D and multidimensional models in a variety of works [e.g., 56, 57, 58]. O’Connor and Couch [53] explicitly demonstrated that their models did not explode when Newtonian gravity was employed.
Explosions in 2D for larger sets of progenitors were also reported by Japanese groups [59, 60] for Newtonian simulations with RbR neutrino transport based on the isotropic diffusionsource approximation [IDSA; 61], by the Basel group [62] for Newtonian models with multidimensional IDSA (and a highdensity equation of state different from the ones used in all other works), and by the Garching group for relativistic simulations with RbR+ twomoment transport and Boltzmann closure [15, 57, 63, 64, 65].
2.1.2 Assessment
Although it is assuring that the models in the subset of successful cases exhibit similarities concerning their gross features [with the explosions of Refs. 49, 50, defining outliers in the optimistic direction], much of the agreement could just be incidental. Suwa et al. [59], for example, obtain an explosion for the 12 progenitor but not for the 15 and 20 stars investigated by the other groups mentioned above, and Nakamura et al. [60] report an explosion of an 11.2 model with considerably faster shock expansion and higher energy than found by the Garching team [66, 15, 57]. A detailed comparison of the published results, even for the same progenitor stars, is hampered by the fact that the simulations were not only performed with different transport solvers and approximations but also with different hydrodynamic schemes and computational grids, implying different (numerical) perturbations that seed the growth of instabilities. Moreover, different resolutions, different sets of neutrino reactions with different simplifications for the opacities, different equationofstate descriptions in particular also in the lowdensity regime, and different gravity treatments enhance the difficulties of direct comparisons. Dedicated and coordinated code tests involving the competing groups are highly desirable and currently in preparation.
2.2 Need of 3D Models
While the tension associated with discrepant results of different groups is unsatisfactory, supernova theory and 2D simulations have never been joined in a love marriage but simply a partnership of convenience that was enforced by the limitations of computing resources and the constraints of numerical codes. Twodimensional modeling has played, and may continue to play, an important role for understanding basic effects of nonradial flows and their feedback on the neutrino emission and heating in selfconsistent, multidimensional supernova models. The simplifying assumption of axisymmetry eases the modeling enormously by being computationally much less demanding than 3D simulations. Currently, 3D calculations are at the extreme limit of numerical feasibility because of their tremendous requirements of computing resources, which strongly increase with more elaborate neutrino treatment. Therefore 2D calculations have so far been the only way to perform resolution studies with fullscale supernova models and to explore larger model sets for the dependence of explosions on microphysics and progenitor properties.
However, important questions cannot be finally answered in two dimensions. How robust is the neutrinodriven mechanism? What are the uncertainties associated with remaining approximations and simplifications in the neutrino transport and microphysics? How important are amplitude and mode pattern of perturbations in the precollapse core? When and how does rotation play a role for the explosion, when do magnetic fields become dynamically relevant? Can neutrinodriven explosions explain the observed asymmetries of supernovae and supernova remnants as well as pulsar kicks and spins?
Threedimensional simulations are certainly needed to connect models to real observations. Moreover, they are indispensable to quantitatively assess the dependence of the mechanism on hydrodynamic instabilities for several reasons. First, in 3D new phenomena and new modes of instability can come into play that do not exist in 2D. For example, in addition to axisymmetric () SASI sloshing motions in the 2D case, 3D simulations can develop () spiral SASI modes, too [71, 72, 73, 74, 75, 76]. Second, in agreement with theoretical expectations [70], the energy spectrum of turbulent mass motions was found to differ between 2D and 3D and the turbulent cascading to transport energy in opposite directions [Figure 1; 67, 77, 78, 79, 75]. Third, the imposed axisymmetry constrains SASIdriven or convectively driven largeamplitude shock oscillations to the axial direction, which was claimed to lead to unphysical feedback effects on the neutrino emission, in particular when the RbR approximation is applied [55, 54]. Interestingly, however, the differences in the dynamical evolution for models with M1 and RbR transport reported by Skinner et al. [55] decreased when their simulations were performed with higher resolution. Moreover, these differences are in conflict with the overall consistency between the models with M1 transport by O’Connor & Couch [53] and those by the Garching group [52], which were computed with RbR+ transport.
3 The Explosion Mechanism in Three Dimensions
For the reasons discussed in Sect. 2.2, 3D simulations are the longdesired, natural next step in supernova modeling with gridbased hydrodynamics schemes. First results have been published since 2010 and have already led to interesting insights and even the discovery of new phenomena.
3.1 Parametric and Selfconsistent 3D Modeling
Similar to the spectrum of published 2D simulations reported in Sect. 2.1, 3D calculations are performed with a wide variety of hydrodynamics codes and modeling strategies, using different grids and resolutions, different gravity treatments, different equations of state (ranging from simple idealgas laws to stateoftheart descriptions of the hadronic and leptonic plasma components in all relevant regimes of density, temperature, and composition), and different approximations for the neutrino treatment, neutrino opacities, and included neutrino reactions. All currently applied transport schemes involve approximations, because a rigorous, timedependent solution of the Boltzmann transport equation in sixdimensional phase space (with three spatial coordinates and three momentum components) is not feasible on the available supercomputers but will require exaflop capability.
The neutrino treatments in timedependent 3D models in the literature can be sorted into three basic categories, differing in their degree of sophistication, namely:
All nonleakage transport treatments used in timedependent 3D simulations so far have employed the RbR(+) approximation for the directional variations except the works by Kuroda et al. [58, 99], which, however, are severely limited in resolution and could cover only some 10 to 100 ms after bounce.
Ignoring neutrinos completely or not describing neutrino transport implies a lack of selfconsistency. In particular the omission of feedback effects of the hydrodynamics (e.g. of accretion or rotation) on neutrino emission and heating can be problematic and can seriously limit the conclusions that can be drawn. In the following, such aspects of numerical modeling will not be commented on unless they are relevant for the results and associated discussion.
3.2 The Importance of Nonradial Flows
The impact of multidimensional flows on the neutrino mechanism has been investigated since the publications of the first 2D models [37, 100, 101]. Now with the focus on the recent 3D results, the discussion of the question is still going on, which effects play a role and on which level of importance.
3.2.1 Growth of Convection and Buoyancy
The negative entropy gradient created by neutrino heating can cause the onset of convective overturn in the postshock layer, whose growth rate (inverse growth time scale) is connected to the complex BruntVäisälä frequency: for convective instability. The inward motion of the postshock accretion flow (with negative radial velocity ), however, suppresses the linear growth of buoyant perturbations unless [102]
(0) 
Here, all quantities (, , the gain radius and shock radius ) have to be angleaveraged as discussed, e.g., in Ref. [103], and and are the advection time scale and convective growth time scale in the gain layer, respectively (the angle brackets denote volume averages). The typical accretion velocity of the postshock flow, , scales with the infall velocity of the stellar matter ahead of the shock, which is roughly given by the freefall velocity at the shock: with being the density jump in the shock and being the mass providing the accelerating gravity field (approximated by the neutronstar mass). The threshold defined by , however, can be circumvented when the initial density perturbations (compared to the ambient density ) are in the nonlinear regime already:
(0) 
where is the average value of the gravitational acceleration in the gain layer. In this case a smallscale perturbation is able to rise against the accretion flow. If the whole flow is perturbed, the buoyant motions on small scales affect the situation globally and can allow for the onset of convective overturn also on larger scales [104].
textbox[t] STANDING ACCRETION SHOCK INSTABILITY (SASI) The SASI is a generic, nonradial instability of stagnant accretion shocks that leads to largescale shock deformation with dipolar and quadrupolar () spherical harmonics modes having the biggest growth rates. In the nonlinear regime violent, largeamplitude sloshing and spiral motions of the shock are a consequence [e.g., 22, 71, 105, 72, 104]. The oscillatory growth of the SASI from initial perturbations is caused by an advectiveacoustic cycle in the cavity between stalled shock and accreting protoneutron star [20, 21, 106, 104, 107]. The saturation amplitude of the SASI is determined by dissipative flow effects like parasitic instabilities (e.g., secondary KelvinHelmholtz and RayleighTaylor instability), which extract energy from the largescale coherent flow. A laboratory experiment with a hydraulic jump in a shallow water flow is a physics analogue of the stalled accretion shock that shares basic features with the SASI in supernova cores [apart from neutrino heating and associated buoyancy; 33, 19]. Symmetry breaking of modes to spiral SASI modes requires that the ratio of initial shock radius to neutronstar radius exceeds a critical value of about two [108]. Angular momentum separation and SASI induced explosion asymmetries can be important to explain neutronstar kicks [109, 87] and spins [71, 110]. \end@float
3.2.2 Growth of the SASI
The linear growth rate of SASI shock deformation modes due to an amplifying advectiveacoustic cycle in the accretion flow between stalled accretion shock and neutronstar surface can be coined as [106, 104]:
(0) 
Here, is the cycle efficiency and
(0) 
the duration of the cycle as the sum of sound travel time (second, subdominant term; is the local sound speed) and advection time between shock radius and a cyclecoupling radius , which is located in the flowdeceleration region between neutronstar radius () and gain radius. According to Eq. (3.2.2) the SASI growth rate scales inversely with the cycle period and SASI activity is expected to be strongest during retraction phases of the accretion shock.
3.2.3 Consequences for the Explosion
Dimension is a key to the neutrino mechanism of corecollapse supernova explosions [80]. So far, however, it is not clear which key exactly fits into the keyhole.
3.2.3.1 Improvements Compared to 1D
It is undisputed and supported by all modern simulations that nonradial postshock flows enhance the neutrinoenergy deposition compared to the 1D case. This can be concluded from higher heating rates () and higher heating efficiencies [; is the luminosity of one neutrino species] in the gain layer. Convective buoyancy and SASI motions push the shock to larger radii, thus increasing the mass () in the gain layer, which can be interpreted as a stretching of the effective advection time through the gain layer, with being the progenitorspecific massaccretion rate [111, 66], or as an increase of the dwell time of matter in the gain layer [112, 113].
In a more detailed picture of the flows, accretion downdrafts channel cool matter close to the gain radius, where efficient neutrinoenergy deposition takes place. Simultaneously, rising highentropy plumes carry the neutrinoheated gas outwards, away from the gain radius. The corresponding expansion cooling of this gas diminishes the energy loss by reemission of neutrinos.
3.2.3.2 3D versus 2D: No Consensus
A highly controversial and still not completely settled discussion was instigated by the question how nonradial flows differ in 2D and 3D and what the corresponding consequences for the explosion mechanism could be. In first, parametric (neutrino light bulb) simulations in 3D, Nordhaus et al. [80] found considerably earlier explosions than in 2D, requiring 15–25% lower driving neutrino luminosity than in 2D. Also subsequent, revised work by the Princeton group [114, 78] still showed (slightly) more favorable explosion conditions in 3D. These results could not be reproduced by Hanke et al. [67], who saw little difference between 2D and 3D for their standard resolution but easier explosions with higher angular resolution in 2D and the opposite trend in 3D.
The findings by Hanke et al. [67] were confirmed by selfconsistent 3D simulations of the Garching group with highfidelity neutrino transport [74, 98]. They also received support by other groups: explosions in 3D occur less readily than in 2D [77, 84, 85, 86, 51]. With better numerical resolution, Couch [77] and Abdikamalov et al. [75] obtained a longer delay of the shock revival and higher values of the critical neutrino luminosity for explosions in 3D, in accordance with Hanke et al. [67].
Fernandez [81], in support of Nakamura et al. [76] and Iwakami et al. [115], observed that SASIdominated explosions can be obtained with up to 20% lower driving luminosity in 3D than in 2D because of the ability of spiral modes to store more nonradial kinetic energy than linear sloshing modes. The magnitude of this difference, however, decreases with increasing resolution.
While all of these works were concerned with the question which critical neutrino luminosity is needed to trigger shock revival, Handy et al. (2014) approached the problem from a different angle and fixed the final explosion energy instead of the driving neutrino luminosity in their parametric setups. They found that 3D models require lower neutrino luminosities to produce equally energetic explosions as 2D simulations and concluded that the “convective engine” in their models is 4% more efficient in 3D than in 2D. This result is not in conflict with any of the other ones because the authors explored a different question.
3.2.3.3 Controversial Resolution Effects in 2D
In contrast to the Garching group [116, 67, 52], Couch [77] as well as Skinner et al. [55] observed a trend to later explosions also for better resolved 2D models. Couch & Ott [117] speculated that this discrepancy might be a consequence of massively underresolved turbulence in the Garching models, while their own simulations produce the correct behavior in the inertial range of turbulence. However, it is equally well possible that the disparate resolution dependence of the 2D simulations is simply caused by differences of the numerical schemes.
The code used at Garching [employing spherical polar coordinates or a YinYang grid; 123] retains perfect spherical symmetry for spherical initial data. Therefore the growth of hydrodynamic instabilities must be seeded by artificially imposed perturbations. Adding numerical resolution in these simulations indeed reduces the effects of numerical viscosity and thus can be conducive for the development of hydrodynamic instabilities. The trend to (slightly) earlier 2D explosions with higher angular resolution persists all the way from 3 down to 0.5 bin width [67, 52] with signs of convergence between 1 and 0.5 [116]. This finding can well mean that the flow dynamics that supports the explosion is not fully described by the concept of turbulence. An inversion of the trend was obtained for higher radial resolution by Hanke et al. [67], however, not because of a better representation of turbulence as claimed by Couch & Ott [117], but because of a resolution sensitivity connected to the simple parametrization of approximative neutrino heating and cooling terms. In selfconsistent simulations with neutrino transport the trend of easier 2D explosions with higher angular resolution also holds for models with enhanced radial resolution [52].
In contrast to the Garching code, schemes used by other authors (mostly with cartesian grids) create perturbations for numerical reasons. Runs with higher grid resolution may possess a lower level of such intrinsic noise and correspondingly show less favorable conditions for the growth of hydrodynamic instabilities, delaying the explosion. A verification of this possibility will require detailed comparisons of simulations for well defined testsetups.
3.2.4 Assessment
Some of the conflicting findings reported in the literature may be a consequence of different codes and setups with different levels of numerical noise and different grid resolution in different domains of the simulations. Directly comparing angular resolution of higherorder solvers on polar grids with “effective” angle resolution of cartesian grids [117] is extremely misleading. Some of the seemingly contradictive conclusions may also be traced back to the use of different parametric modeling approaches with different simplifications of complex microphysics and different quantities varied or kept fixed. Finally, also the lack of selfconsistency and the omission of feedback mechanisms that couple hydrodynamics and neutrino emission may be the cause of misleading results. Cautious interpretation is therefore advisable and farreaching conclusions should be avoided.
3.3 Turbulence vs. Convection vs. SASI: an Ongoing Debate
The trend of later explosions in 3D compared to 2D [superimposed by considerable casetocase stochasticity, see 86] is understood as a consequence of the turbulent cascading of energy to small scales in 3D, whereas in 2D energy is pumped into the largest structures [67, 77, 84, 78].
3.3.1 Turbulent Kinetic Energy Spectrum in 2D and 3D
Figure 1 shows energy spectra of nonradial (longitudinal) mass motions based on a decomposition of the azimuthal velocity at a given radius (located in the gain layer) into spherical harmonics :
(0) 
[cf. Ref. [67] for details; for a discussion of differences between spectra for full 3D velocities and velocity components, see Ref. [118]]. As expected from the direction of the energy cascade, in 2D the lowest modes clearly dominate the energy spectrum and the spectral decline towards high spherical harmonics modes (small wavelengths) is steeper (roughly ) in the inertial range of the forward enstrophy cascade. In contrast, in 3D the spectral decay in the inertial range of the forward energy cascade is closer to . While this basic structure of the energy spectra was confirmed by many groups [e.g., 77, 84, 78, 75, 86, 79], vivid twitter started about the exact shape of the powerlaw decline. Some authors interpreted their 3D results by an decay, piecewise power laws, or an exponential slope instead of the cascade expected for Kolmogorov turbulence [68, 69] and surmised that anisotropic turbulence [the radial component of the Reynolds stress is in rough equipartition with the summed tangential components; 119, 120], nonstationarity or flow compressibility might be responsible for nonKolmogorov behavior. However, Radice et al. [118, 82] showed that increasing resolution enables the development of a power law with exponent in an increasing inertial range of wave numbers.
3.3.2 Big Bubbles or Turbulence as Explosion Drivers
The energy cascade in 2D thus feeds the biggest buoyant plumes, which have been recognized as conducive for driving shock expansion [67] because of their better volumetosurface ratio, enabling buoyancy to win over drag effects [77, 103]. 3D is less favorable for the creation of such big bubbles through convective and turbulent processes, for which reason Hanke et al. (2012) emphasized the possible importance of the SASI as a natural mechanism to foster the growth of lowmode asymmetries and to push shock expansion, provided the saturation amplitude of SASI modes (which is limited by parasitic effects like RayleighTaylor and KelvinHelmholtz instabilities) is sufficiently large.
Others continued to favor buoyancydriven convection [114, 119, 83], “penetrative convection” [79], and “turbulent convection” [120, 75] as the main supportive instability in the transition from the stalled shock to outward shock acceleration. Couch & Ott [117], inspired by Murphy et al. [119], proposed that turbulence provides effective pressure that adds to the gas pressure in stabilizing the gain layer, thus favoring the accumulation of energy, mass, and momentum behind the shock and lowering the critical neutrino luminosity for explosion. This hypothesis was motivated by their observation that multidimensional simulations explode with lower neutrinoheating rates than 1D simulations. In this scenario enhanced neutrino heating by the convectively stretched residence time of matter in the gain layer and turbulent Reynolds stresses conspire in fostering shock revival. Couch and Ott [117] conclude a direct correlation between the strength of turbulence and the susceptibility to explosion.
3.3.3 Artifacts of Modeling
The setups investigated by Radice et al. [82], Burrows et al. [114] and Dolence et al. [78] on the one hand and those of Couch & Ott [121, 117] on the other hand are probably not representative for supernovae and their progenitors in general. In the former works the shock stagnation radius is very big (exceeding 300 km) and lingers there in a quasisteady state for many 100 ms. This large, quasistationary shock radius is likely to be an artifact of the parametrized nuclear photodisintegration and neutrino treatment. It is not typical of the shock behavior in more realistic simulations with proper nuclear and neutrino physics, where the shock initially expands to only 150 km, then tends to retreat again, and finally accelerates outwards quickly once it achieves to expand beyond the nucleon dissociation/recombination radius of km ( is the baryon mass).
In the works of Couch & Ott [121, 117], largeamplitude, largescale velocity perturbations in the Si/O layer of the progenitor model were assumed. But even in their simulations without such an imposed velocity field, large nonradial mass motions develop in the gain layer already 20–30 ms after core bounce, in contrast to the simulations by the Garching group, where significant nonradial flows in the gain layer start only 100 ms post bounce. The rapid growth of turbulence points to a sizable amplitude of numerical, probably gridinduced, perturbations in the models of Refs. [121, 117], which could instigate buoyancy in an uncontrolled way and could thus lead to an overestimated importance of turbulent convection. The implications of such purely numerical effects for modeling supernovae must still be clarified.
3.3.4 Relevance of the SASI
SASIdominated phases were observed in selfconsistent 3D supernova simulations [e.g., 74, 94, 95, 96, 98] as well as parametric studies [84, 75, 81] but were weak in others, e.g. Ref. [51] for a fullscale supernova model and Refs. [83, 114, 78, 77, 117] for parametric models. Obviously, unless SASI oscillations were damped by a lack of resolution [cf. 75] or disfavored by the parameters of the modeling setup, the simulated models seem to have passed through different dynamical regimes: the growth of the SASI is favored for small shock radii, in which case the SASI cycle time is short (Eq. 3.2.2) and its growth rate large (Eq. 3.2.2), whereas convective activity is facilitated by bigger shock radii (Eq. 3.2.1) or by greater numerical or imposed perturbations in the infalling stellar core matter (Eq. 3.2.1).
Burrows et al. [114], tentatively seconded by Couch & O’Connor [84], hypothesized that neutrinodriven buoyant convection should almost always dominate the SASI in neutrinopowered supernova explosions. However, Müller et al. [63] and Fernandez et al. [103, 81] demonstrated that neutrinodriven explosions can develop from both the regimes where SASI or convection dominate the dynamics. SASI mass motions, in particular spiral modes, act as a storage of kinetic energy on large scales and create secondary shocks that heat the matter by dissipating kinetic energy. Both effects aid shock expansion.
3.3.5 Assessment
Despite strong opinions in favor of exclusive explanations, a more likely possibility is that all effects of nonradial flows —convective buoyancy, turbulent pressure, and SASI motions— can influence the supernova shock dynamics and could provide support to the postshock layer. The relative importance of different nonradial instabilities can depend on the phase of the postbounce evolution and on the detailed expansion and contraction behavior of the stalled shock, which is sensitive to progenitorspecific properties of the accretion flow. Present results of simulations may also be affected by artifacts connected to the numerical grid and by technical features in the (simplified) modeling setups. Future, wellresolved and fully selfconsistent 3D simulations for larger sets of progenitors and realistic precollapse perturbations in codes with low intrinsic noise level are needed to confirm our expectation that the cores of collapsing stars can evolve through SASIdominated episodes at least transiently.
3.4 Selfconsistent 3D Explosion Models
Successful neutrinodriven explosions were recently obtained in selfconsistent 3D simulations by several groups using different neutrinotransport approximations, all based on RbR or RbR+ treatments of the angle dependence on the employed coordinate grids.
3.4.1 Recent Results
Takiwaki et al. [85, 86], applying an IDSA code for and transport in combination with a leakage scheme for , reported explosions for an 11.2 star, but their models had modest resolution of 300 nonequidistant radial zones and at best 2.8 in the polar and azimuthal directions. Müller [93] also obtained an explosion of this progenitor, using a stationary twostream solution of the RbR transport equation combined with an analytic variable Eddington factor closure [92] and higher resolution (550128256 zones in , , and directions).
The Oak Ridge and Garching groups found explosions with more sophisticated, energydependent, threeflavor transport solvers including more microphysics and intermediate numerical resolution. Oak Ridge, employing the Chimera code with fluxlimited diffusion and 540 logarithmically spaced radial zones and a 180180 zone (nonuniform)  grid, investigated a 15 progenitor [51]. At Garching, explosions were obtained for 9.6 and 20 stars [97, 98] and for a rotating 15 progenitor [Figure 2; 122] with the PrometheusVertex code, using a twomoment transport scheme with Boltzmann closure, 400–600 nonequidistant radial zones (improving the resolution dynamically) and 2 cell size in both angular directions of axisfree YinYang [123] and standard polar grids.
3.4.2 New Messages
The results of all these simulations agree in their basic finding that 3D models are less susceptible to explosion than 2D models. Shock revival, if happening at all, occurs later in 3D than in 2D. This outcome is in line with conclusions drawn from most parametric studies (cf. Sect. 3.2). None of the mentioned models could so far be evolved to the point where the final explosion energy was determined. In earlier simulations of 11.2, 20, and 27 stars the Garching team had not seen successful shock revival until 500 ms after bounce, although the corresponding 2D models had already developed explosions at that time [74, 95, 96, 94]. Diagnostic parameters that signal the proximity to explosion like the maximum shock radius and the ratio of advection to heating time scale suggested marginal failures. Indeed, in the recent 20 simulation [98] a reduction of neutralcurrent neutrinonucleon scattering by effectively 10–20% in the neutrinospheric layers (motivated by possible strangeness contributions to the nucleon spin, affecting the axialvector weak coupling) led to an increase of the neutrino luminosities and of the neutrino heating behind the shock that was sufficient to turn the failed explosion to success. Although for the purpose of demonstration the magnitude of the considered strangeness correction was overestimated compared to the currently best experimental and theoretical limits [124], the result of Melson et al. [98] calls attention to an important sensitivity of 3D supernova models to even only smaller variations of the neutrino opacities.
Interestingly, the 9.6 explosion of Melson et al. [97] developed a slightly higher explosion energy than the corresponding 2D case. This can be explained by differences of the dynamics of convective downflows. In 3D KelvinHelmholtz instability leads to more efficient fragmentation, which decelerates the downdrafts and keeps more matter in the gain layer, thus reducing neutrinoenergy losses below the gain radius. Müller [93] observed a similar, even larger effect in his successful 11.2 model, which exploded earlier and more powerfully in 3D than in 2D [in contrast to the results by 85, 86]. He traced the more favorable 3D situation back to a variety of subtle differences in the accretion and outflow dynamics in 2D and 3D, leading to more efficient driving of neutrinoheated gas outflow in 3D. Whether these interesting 3D effects also apply to more massive progenitors and how they might depend on numerical resolution still needs to be figured out.
3.4.3 Implications
These results suggest that the discussion of 3D effects in comparison to 2D involves at least two separate aspects:

What is the dimensionality dependence of the dynamics prior to the onset of explosion? Do 3D simulations develop shock revival faster or more delayed than in 2D? The majority of current models (with few exceptions) suggests the latter.

What are 3D vs. 2D differences after the onset of the explosion? Can 3D flow dynamics enhance the explosion energy and thus bring 3D models closer to observed supernova energies? Here the results of Refs. [98, 93] suggest interesting advantages of the 3D case. Parametric explosion studies in Ref. [79] seem to yield support.
Overall, one is therefore tempted to conclude that the selfconsistent 3D calculations provide backup for the neutrinodriven mechanism. When the models fail, not much seems to be missing to achieve shock revival. The 2D and 3D simulations of the Oak Ridge group (also those of Series C) are considerably more optimistic than those from Garching, despite similarities in many aspects of the numerical modeling [but also differences in quite a few, cf. 51]. The reason for this discrepancy is unclear and can be clarified only by detailed and direct comparisons.
3.4.4 Caveats
A serious drawback of the current selfconsistent simulations are the constraints of numerical resolution because of the enormous computational demands when detailed microphysics and transport are included. Radice et al. [118, 82], assuming that steadystate turbulence applies and setting up suitable conditions, demonstrated that Kolmogorov scaling can progressively be recovered as the resolution in their toy models is increased. Although they found that the resolution of the published stateoftheart supernova simulations seems to be sufficient to describe gross features of neutrinodriven turbulent convection, the results of Refs. [118, 82] for successive grid refinements exhibit subtle and partly nonmonotonic differences, whose relevance for supernova dynamics is currently unclear. Better resolution in particular diminishes a “bottleneck effect”, in which numerical viscosity hinders the efficient cascading of turbulent energy to small scales and keeps energy on large scales. Abdikamalov et al. [75] expressed concerns that this bottleneck might incorrectly and artificially promote explosions. It will have to be tested whether present selfconsistent supernova simulations are affected by such artifacts as soon as higherresolution calculations can be afforded.
3.5 A Universal Explosion Criterion
Simulations at the onset of runaway shock expansion exhibit a wide range of variation of their diagnostic parameters, e.g. of their neutrino luminosities, average and maximum shock radii, mass accretion rate, total mass in the gain layer and mass fraction of recombined nucleons there, neutrinoheating rate and efficiency, average and maximum entropy, turbulent kinetic energy in the gain layer, etc. The question therefore arises whether there is any individual parameter value or parameter relation that has to be matched when shock revival shall occur? Could such a criterion capture the influence of dimensionality and of multidimensional effects on the susceptibility to explosion?
Burrows & Goshy [125] proposed a critical luminosity condition as a function of the massaccretion rate of the stalled shock, , which defines the threshhold value of the neutrino luminosity that needs to be exceeded to enable shock expansion against the ram pressure of the infalling material. Their semianalytic analysis and followup works (numerical and analytical) by others confirmed the existence of an upper limit of the neutrinoenergy deposition in the gain layer that is compatible with solutions for quasistationary accretion shocks [for a review, see 14].
Although there are interesting proposals of alternative explosion criteria like an antesonic condition [126, 120] and an integral condition for the gain layer [127], the considerations here will focus on a generalization of the critical luminosity relation.
3.5.1 Critical Luminosity in Spherical Symmetry
Janka [14] argued [see also 92] that the radius of the stalled shock in 1D models approximately follows the scaling relation
(0) 
where is defined as the total luminosity of plus , , and denotes the weighted average of the corresponding mean squared energies:
(0) 
Both determine the neutrinoenergy deposition in the gain layer, which is given by
(0) 
where (for ) are defined as squared energies of the neutrinoenergy distributions expressed in terms of the energy moments of the neutrinonumber distributions.
Janka [14] showed that the critical luminosity condition can be deduced from equating the advection time scale,
(0) 
with the heating time scale,
(0) 
Here, is the average massspecific binding energy in the gain layer,
(0) 
where is the total (internal plus gravitational plus kinetic energy) and the mass in the gain layer. In Eq. (3.5.1) the approximations are used that and, roughly, . Employing as well established condition that signals runaway shock expansion [for a detailed discussion, see, e.g., 128] one derives [cf. also 92]:
(0) 
Following Ref. [52] we do not omit and in this relation.
The path to generalize this condition to the multidimensional case was exemplified by Müller & Janka (2015) and Summa et al. (2015). Despite the complex nature of the postshock flow including violent, largeamplitude SASI sloshing, bubble buoyancy, and turbulent convection and flow fragmentation in a highly nonstationary environment, these works found that the overall effects of nonradial mass motions seem to be captured astonishingly well by a simple concept described in the following.
3.5.2 Effects of Turbulence
Guided by Ref. [92] we take multidimensional (“turbulent”) mass motions in the gain layer into account through an isotropic pressure contribution that is coined in terms of the squared Mach number of unordered fluid flows, averaged over the gain region: , ignoring additional complexity, e.g., by turbulent energy transport or centrifugal support [cf. 118, 120]. Including this additional postshock pressure in the shockjump condition, one finds a larger radius of the stalled shock compared to Eq. (3.5.1):
(0) 
with
(0) 
Through a corresponding stretching of the advection time scale (Eq. 3.5.1), the condition leads to a modified critical luminosity [92, 52]:
(0) 
Different from Refs. [92, 52], we empirically choose the total anisotropic kinetic energy, , in the definition of instead of just the lateral component:
(0) 
with being angular averages over spherical shells. This generalization shall make Eq. (3.5.2) applicable to 2D and 3D results with the same proportionality factor. The velocity differences in Eq. (3.5.2) are meant to subtract ordered radial flows like accretion or expansion of the gain layer and coherent angular motions associated with stellar rotation and spiral SASI modes. In contrast to Müller & Janka [92], we do not apply an approximation for the sound speed behind the shock but extract directly from the numerical simulations as massweighted average over the gain layer:
(0) 
For marginal stability of the gain layer to convection in a steadystate situation, the volumeintegrated neutrino heating rate must be balanced by outward “turbulent luminosity”, i.e. [119]. According to Müller & Janka [92], this requirement is reflected by the following relation between and :
(0) 
Here and denote angle averages of the shock and gain radius, respectively. This relation turns out to be well fulfilled when the average shock radius is nearly stationary, but large deviations occur when the shock moves.
3.5.3 Effects of Rotation
The effects of rotation can be included in a similar way. Rotation provides centrifugal support in the infall region ahead of the shock. Instead of the freefall velocity of the matter, , which was used in Eq. (3.5.1), one now gets
(0) 
from solving the equation of motion, , assuming conservation of the specific angular momentum, , during collapse. In Eq. (3.5.3) we have introduced the rotational correction factor
(0) 
Directionaveraging the effect of rotation, one can define as average of the specific angular momentum on spherical shells. Considering in the postshock layer a modest radial increase of the spherically averaged specific angular momentum, [111], centrifugal effects in rotational equilibrium can be absorbed into a correction of the gravitating mass, . The usual equation of hydrostatic equilibrium thus applies with being replaced by . Therefore the density, temperature, and pressure profiles in the relativistic gaspressure dominated gain layer can still be approximated by , , [129]. Making use of this result and of Eq. (3.5.3) for the preshock velocity, a derivation in analogy to the one in Refs. [14, 92] yields for the radius of the stalled shock in the presence of rotation:
(0) 
Accounting for both this change of the shockstagnation radius and for the rotational deceleration of the infall velocity in the postshock advection time scale (Eq. 3.5.1), the condition yields now:
(0) 
Rotation leads to a larger shockstagnation radius (Eq. 3.5.3 with ). It also decreases the critical luminosity because of the factor . In addition, rotational energy in the gain layer provides a positive contribution to , shifting closer to zero, which also decreases the rhs of Eq. (3.5.3).
3.5.4 Universal Critical Luminosity Condition
Including effects of rotation as well as unordered (“turbulent”) mass motions requires the combination of both factors and in the shock stagnation radius (cf. Eqs. 3.5.2 and 3.5.3),
(0) 
and in the critical luminosity criterion (cf. Eqs. 3.5.2 and 3.5.3),
(0) 
where the timedependent quantity subsumes all gainlayer related properties:
(0) 
can be used to correct for variations of the time evolution of gain radius, binding energy, nonradial (turbulent) postshock flows, and rotation, which lead to time and modeldependent variations of the critical luminosity in addition to the basic dependence on and . Doing so we derive a generalized, universal relation for the critical luminosity [cf. 52]:
(0) 
The arbitrary constant is introduced as normalization of the correction factor relative to a chosen reference model, for which is evaluated at the time when .
3.5.4.1 Application to Selfconsistent Supernova Models
Figure 3 shows evolutionary tracks (running from right to left) in the plane spanned by and for the set of 2D models of Ref. [52]. The three panels also include two exploding 2D simulations with rotation as well as exploding and nonexploding 3D models computed at Garching (see Sect. 3.4; the cases with explosions are displayed in Figure 2). Colored symbols on top of the tracks mark the instants when the explosion sets in (i.e. when ) and correspond to the critical luminosities of Eq. (3.5.4). Open symbols in the top panel indicate the locations of the uncorrected values of Eq. (3.5.4).
Obviously, for all 2D models the correction factor successfully compensates for the influence of nonradial mass flows, rotation, and modeltomodel variations of gain radius and specific energy in the gain layer. Consequently, the colored symbols obey a tight correlation as expected from Eq. (3.5.4), and the dashed black line defines a universal critical explosion condition that holds for 2D simulations.
The remaining lowlevel scattering could be linked to approximations that entered the derivation and evaluation of the scaling relation of Eq. (3.5.4), e.g.: simple power laws for the radial structure of the gain layer; assumed scaling of with postshock quantities; use of and as measured at infinity instead of values at the gain radius; omission of gravity contributions from the mass of the gain layer and pressure corrections in computing the preshock velocity (Eq. 3.5.3); simple averaging of rotation effects on spheres; or the assumption of a modelindependent compression ratio in the shock.
The upward bending of the evolutionary tracks at the onset of explosion is caused by a steep drop of in the denominator of , while in the numerator evolves slowly. The decline of (Eq. 3.5.4) occurs because a strong increase of supports the outward acceleration of the shock and, as a consequence, the specific binding energy of the gain layer () plummets.
3.5.4.2 Behavior of 3D Models
The 3D models deserve special discussion because their behavior is less uniform.
Except for model z9.63D, all 3D simulations with successful and failed explosions for nonrotating stars fit the general picture that applies for the 2D cases. The evolution tracks of all nonexploding 3D runs including the ones that marginally fail, remain below the (corrected) critical luminosity curve defined by the 2D simulations (see Figure 3, upper two panels). The nonrotating 3D model s203Ds [98] explodes with about half the value of the average squared Mach number (; Figure 3, bottom panel) of typical 2D models, but its critical point coincides perfectly with the critical luminosity curve.
The lowmass model z9.63D [97] clearly defines an outlier. Its explosion starts () long before the evolution track reaches the critical luminosity curve of Eq. (3.5.4). Here the blastwave acceleration is driven by the outward expansion of buoyant plumes with neutrinoheated matter, and postshock turbulence does not play an important role. Also the specific energy in the gain layer does not increase significantly when runaway shock expansion occurs, because the rising plumes of dilute gas contain only a minor fraction of the mass in the gain layer. In fact, the uncorrected quantity (Eq. 3.5.4) of model z9.63D lies close to the critical curve for all other models (see open star in Figure 3, top panel), and the corrected quantity is farther away from this curve. This unfavorable trend is caused by the normalization with the factor , which is picked from one of the 2D explosion models and reflects the small values of these cases. This observation illustrates that model z9.63D blows up nearly like a 1D explosion. Under such circumstances normalization with from multidimensional models is not appropriate.
The rapidly rotating model m15u63Dartrot also explodes, in contrast to the more slowly spinning case m15u63Drot, which fails to blow up [122]. At the onset of the explosion the value of in model m15u63Dartrot is similar to those of the 2D cases, although its evolutionary track indicates the special, rotationdominated postbounce dynamics of this model by passing the critical luminosity nearly horizontally before bending upwards sharply (Figure 3, middle panel). Shock revival in this case is fostered by a strong spiral SASI mode in agreement with findings in Refs. [76, 115]. It is important to note that preshock rotation plays a negligible role since . The main effect of rotation is its contribution of kinetic energy as part of the total energy in the gain layer. In addition, spiral SASI waves also trigger secondary turbulent mass motions, creating effective turbulent pressure. Both effects are accounted for in our formalism of Sects. 3.5.2 and 3.5.3. Sufficiently rapid rotation can therefore facilitate the explosion and reduces the critical threshold for the neutrino luminosity (Eq. 3.5.4). This can be seen by the open star of model m15u63Dartrot in the top panel of Figure 3. Although our correction procedure of Eq. (3.5.4) moves the critical luminosity of this model closer towards the threshold line, it still lies somewhat below this curve and agrees with the universal relation less well than the other models. The crossing of the evolution track of model m15u63Dartrot therefore happens somewhat later than the equality is fulfilled. However, since both the evolution track and the critical curve are fairly flat, small uncertainties in the vertical location can cause a considerable shift of the crossing point. For this reason it is not clear whether the slight mismatch of m15u63Dartrot is a mere incidence or whether it points to shortcomings or incompleteness of the theoretical framework to describe the effects of rotation in this case.
We conclude that the dynamical evolution of model z9.63D definitely and of m15u63Dartrot possibly is not well captured, or at least not fully represented, by the theoretical concept that underlies the derivation of Eq. (3.5.4). The functional relation given by this equation seems to define a universal critical luminosity condition for the subclass of 2D and 3D models where strong, nonordered (“turbulent”) flows play a major role for the postshock dynamics. Deviations from the relation of Eq. (3.5.4) may signal cases where the approach to explosion is determined by alternative dynamical phenomena like bouyancy or spiral SASI motions that drive the shock expansion. A larger set of 3D models is needed, however, to consolidate this framework.
4 New Phenomena in Three Dimensions
The first selfconsistent 3D simulations have already enabled discoveries of new phenomena, which are associated with the special dynamical behavior of 3D flows in the supernova core and can cause characteristic imprints on the neutrino emission. Predicting such signal features in detail is of crucial importance for the interpretation of the neutrino measurement from a future Galactic supernova.
Based on its particularly low photomultiplier dark noise rates, IceCube is an excellent supernova burst detector by observing a collective rise in all photomultiplier rates on top of the dark noise. Although IceCube will not be able to provide directional information and no eventbyevent energy information, it is unique among all existing SN neutrino observatories in its ability to record the signal with a 2 ms timing resolution. IceCube will thus be sensitive to subtle features in the time structure of the neutrino signal.
4.1 SASI Sloshing and Spiral Modes and NeutrinoEmission Modulations
Largescale, largeamplitude oscillations of the stalled shock and postshock layer associated with SASI motions cause modulations of the massaccretion flow onto the newly formed neutron star. These accretion variations lead to quasiperiodic fluctuations of the neutrino emission, which can well be detected with IceCube for a Galactic supernova [130]. This phenomenon was first observed in 2D simulations with RbR+ neutrino transport [131] as well as multidimensional, multiangle neutrino transport [132, 133].
However, doubts about the physical reality of the phenomenon were expressed because in 2D the imposed axial symmetry constrains the SASI sloshing to the axial direction and the accretion flow occurs mostly close to the equatorial plane, potentially exaggerating the asymmetry. Such a geometrical constraint does not exist in 3D, where in addition symmetrybreaking () spiral modes can occur due to a superposition of phaseshifted bipolar SASI oscillations in different directions [e.g., 71, 72, 73, 74]. Angular momentum in the collapsing stellar core can accelerate the growth of this spiral SASI mode [134].
Tamborra et al. [94, 95] analyzed the 3D simulations of the Garching group for neutrinoemission asymmetries. In order to evaluate the computational models for the observable neutrino signals from different viewing directions, they integrated the neutrino emission over the whole visible hemisphere for all observer positions. In this postprocessing of the numerical data, they applied a method introduced in Ref. [89] to take into account limb darkening and projection effects. This procedure corrects the artificial enhancement of smallscale emission variations associated with the use of RbR+ neutrino transport in the supernova simulations.
During phases of SASI sloshing and spiral motions the emission of and , whose production dominates in the hot accretion flows, exhibits synchronous time modulations with a spatial correlation. On a somewhat lower level such emission variations can be found also for heavylepton neutrinos (). The amplitude of these fluctuations can be up to 20% of the directionaveraged luminosities of and and up to about 5% for . The mean energies of the radiated neutrinos of all species change by about 1 MeV, also in phase with the luminosity variations. The typical frequency of the fluctuations is 50–100 Hz and shows up as a strong peak in the power spectrum of the neutrinodetection rate. This frequency is closely linked to the inverse of the duration of the SASI cycle, (Eq. 3.2.2), which can be expressed in terms of the shock radius and of the neutronstar radius and mass as [65]
(0) 
Since the growth of the SASI is favored during periods of shock retraction (cf. Sects. 3.2.2 and 3.3.4), a measurement of the duration and frequency of SASIinduced time variations in the neutrino signal from a Galactic supernova would yield information about the shock radius and its evolution prior to the onset of the explosion.
The neutrinoemission variations predicted by the 3D models would be well detectable with IceCube for a stellar death in the Milky Way, in particular from observer positions close to the plane of the SASI activity, where the fluctuation amplitude is largest. But even for observers outside of this plane, which is not fixed but can move and differ between different episodes of SASI activity, there is a promising perspective for detection [94, 95].
In evolution phases without SASI shock motions and in models with dominant convective activity in the postshock layer, the neutrino luminosities still exhibit fluctuations caused by the convective variations of the mass accretion by the nascent neutron star. However, these luminosity fluctuations look similar from all observer directions, and they are less regular and have smaller amplitudes of only some percent, associated with a broader power distribution in the Fourier space. IceCube will still be able to detect these signal features, albeit only for a supernova at a distance up to a few kpc [135].
4.2 Leptonnumber Emisson Asymmetry (LESA)
The 3D simulations of the Garching group also revealed a stunning, new, and unexpected neutrinohydrodynamics phenomenon that develops in all models simulated with energydependent, threeflavor neutrino transport handled by the PrometheusVertex code. It was first described by Tamborra et al. [96], who discovered a dipolar asymmetry of the leptonnumber emission in the neutrino data, which they termed LESA: Leptonnumber Emission Selfsustained Asymmetry.
4.2.1 LESA Phenomenology
Within roughly 200 ms after core bounce a large hemispheric asymmetry of the leptonnumber emission builds up from an initial, highorder multipolar pattern (Figure 4). In fact, the leptonnumber ( minus ) flux develops a dominant dipole amplitude that can become stronger than the monopole. It implies that the newly formed neutron star loses its lepton number predominantly in one hemisphere, i.e., the number flux clearly exceeds the number flux on one side whereas the excess is much smaller (or even the emission stronger) on the other side.
This dipole asymmetry persists for hundreds of milliseconds with its amplitude, after passing the maximum, following the gradual decay of the monopole (Figure 5). The dipole direction changes only slowly in models where convective mass motions determine the nonradial flows in the postshock layer, i.e. it drifts on time scales much longer than the typical convective turnover time scale, accounting for an angular velocity of up to 90 over several hundred milliseconds. During phases with violent SASI spiral motions such a longtime drift can be superimposed by a wobbling of the LESA direction by up to some 10 with the high frequency (Eq. 4.1) of the SASI spiraling [96].
The LESA dipole direction and the SASI shockdeformation vector are uncorrelated and not causally connected. The characteristic neutrinoemission properties of SASI and LESA differ fundamentally [95, 96]. SASI asymmetries and timemodulations are synchronized between all neutrino species and can reach 10–20% of the total energy flux. In contrast, the amplitude of the leptonnumber flux dipole can even exceed the monopole, and the and emission maxima peak in opposite hemispheres. Different from these neutrinos, whose individual hemispheric flux asymmetries can be as high as 20%, heavylepton neutrinos () are affected by the LESA effect only on a minor level of order percent. The total energy flux (summed over all neutrino species) therefore exhibits directional variations and a dipole amplitude on the level of several percent with its maximum pointing opposite to the leptonnumber emission dipole (Figure 4, right column).
4.2.2 LESA Origin
The emission asymmetry of the LESA phenomenon originates mostly from the convection layer inside of the protoneutron star. Up to the neutrinosphere about three quarters of the final leptonemission dipole have built up, and another quarter of the total effect is added in the accretion layer that surrounds the neutrinosphere [96].
The emission dipole is mirrored by a pronounced hemispheric asymmetry of the electron distribution in a thick shell just inside of the neutrinosphere: in the hemisphere of the stronger emission the electron fraction (, which is the number of electrons per baryon) is considerably higher than on the opposite side (see Figure 5, where the crosssectional panel shows this shell bounded by the density contours for g cm and roughly g cm). This asymmetry is a consequence of stronger convection in the protoneutron star on the side of the higher and dominant emission. The correspondingly enhanced convective transport of lepton number carries electrons from the leptonrich, deep core (where ) to the more deleptonized layer () enclosed by the neutrinosphere.
Despite the large hemispheric imbalance of the density and pressure distributions remain spherical (in the absence of rotation) as dictated by the monopoledominated gravitational potential. Since the pressure of the nuclear medium is a function of density, temperature, and electron fraction, , the nonspherical distribution of must be compensated by a small asphericity of the temperature (and entropy). Regions with high must be cooler. This explains why the total neutrinoenergy flux is higher in the hemisphere opposite to the leptonemission maximum (Figure 4, right column).
This interior effect is amplified by an exterior feedback cycle. Since are radiated with somewhat harder spectra than , the neutrinoenergy deposition by and absorption in the gain layer is stronger in the hemisphere of relatively higher emission (i.e., on the side opposite to the LESA dipole direction). The stronger heating pushes the accretion shock to a larger radius, thus creating a global dipolar deformation of the shock that persists (on average) as a longtime phenomenon, overlain by shorttime variations associated with convective and SASI activity in the gain layer. Due to the dipolar deformation the accretion shock deflects the incoming massaccretion flow and channels it preferentially to the hemisphere of the smaller shock radius. This effect enhances the mass accretion of the neutron star and therefore the inflow of fresh lepton number on the side of the higher , thus further adding to the asymmetry and amplifying the leptonnumber emission asymmetry. Tamborra et al. [96] speculate that this accretion asymmetry could function as a selfsustaining mechanism that stabilizes the LESA dipole during the accretion phase.
4.2.3 LESA Sensitivity to Model Variations
Figure 5 displays the development of LESA in a set of 3D simulations for different progenitors, varied microphysics, and different rotation rates. The plots show the time evolution of the power in the monopole, dipole, and quadrupole components () of the spherical harmonics decomposition of the leptonnumber emission, defined by
(0) 
where is the local (radial) number flux density of neutrino species . The monopole is the total leptonnumber flux and the dipole is defined in the same way as introduced in Ref. [96]: in coordinates aligned with the dipole direction.
The monopole decays along with the decline of the massaccretion rate. The quadrupole grows faster than the dipole in all models and both appear with the onset of convection in a Ledouxunstable layer inside of the neutron star. In all cases the dipole power begins to increase quasiexponentially at ms after bounce, saturates at a similar magnitude and dominates the quadrupole component and even the monopole at some point, except in the rapidly rotating model m15. Rotation suppresses the growth of the dipole; model m15 with less angular momentum reveals this effect less strongly.
The growth of the dipole takes place during a phase of rapid contraction and deleptonization of the protoneutron star that is a consequence of its fast gain of mass associated with the high massaccretion rate early after bounce. A connection between neutronstar contraction and LESA growth is also suggested by model s20, where due to higher emission because of a reduced neutrinonucleon scattering opacity [98] the neutron star contraction and the dipole growth are faster than in the reference model s20. Inversely, model s20 was computed with a stiffer nuclear equation of state, which slows down the contraction of the nascent neutron star. Indeed, in this simulation the dipole growth is delayed and its climb to maximum strength takes longer.
4.2.4 Explanation of LESA
The LESA phenomenon is not well understood yet and, in particular, the development of a dominant dipole asymmetry still needs to be explained.
Tamborra et al. [96] reasoned that the two opposite hemispheres “communicate” by the accretion asymmetry around the neutron star (cf. Sect.4.2.2). A detailed analysis, however, reveals that also in the convection layer inside of the protoneutron star there is an effective largescale flow between the two hemispheres, which overlies the more local mass motions associated with the smallscale pattern of convective cells [136].
The quasiexponential amplification of the LESA asymmetry might be connected to the fact that convection in the protoneutron star according to the Ledoux criterion,
(0) 
( is the electronic lepton fraction and the entropy per nucleon), depends on entropy and leptonnumber gradients. In this context it is important that changes its sign from negative to positive for values of below a critical limit. This critical threshold is a function of density and temperature and increases for higher and . As the protoneutron star deleptonizes and contracts, the threshold is eventually passed, in which case a negative gradient of becomes stabilizing instead of destabilizing. This effect damps convective activity. As a consequence, electrons are less efficiently dragged upward from the leptonrich central core, while at the same time outward neutrino diffusion (which becomes faster than convective transport near the outer edge of the convection zone) still carries away lepton number and maintains the ongoing deleptonization of the convective shell. With the decreasing lepton fraction the stabilizing influence of the second term in the Ledoux criterion will further increase. This defines an amplifying feedback cycle that could be the underlying reason for the exponential growth of the leptonemission asymmetry, instigated by a sufficiently large seed perturbation.
LESA would thus manifest itself as a neutrinohydrodynamics instability, in which convective and neutrino transport in combination determine the largescale flow dynamics within global, nonspherical modes that encompass a major fraction of the volume of supernova cores.
4.2.5 Is LESA a Numerical Artifact?
All facts reported in the preceding sections seem to provide support for LESA being a physics phenomenon and not a numerical artifact, e.g. as a consequence of the RbR+ transport treatment applied in the 3D simulations performed at Garching. However, a causal connection between LESA and RbR(+) transport appears unlikely, because the RbR(+) approximation tends to produce localized, smallscale extrema, and it is difficult to imagine how this could help assembling a global dipole asymmetry. LESA is certainly not linked to the use of a polar coordinate grid (with its wellknown axis artifacts), because the dipole direction differs from model to model (with its beginning not being correlated with the direction of the polar axis), and because a dipolar leptonemission asymmetry also develops in simulations with an axisfree YinYang grid (e.g. in models z9.6 and m15 of Figures 4 and 5).
Yet, plausibility and consistency is not evidence. Independent confirmation by other groups is needed and a theoretical framework must be developed that assembles the different pieces of the puzzle described above into a consistent picture. Indeed, the presence of a hemispheric asymmetry in the protoneutron star and of a leptonnumber emission dipole with a size of order the monopole magnitude was also observed in 3D simulations with a neutrinoleakage scheme (Evan O’Connor, private communication 2014). Also the 15 3D explosion model of the Oak Ridge group [51], which applied a fluxlimited diffusion solver for the neutrino transport, exhibits a leptonnumber emission dipole. This feature grows at the same time as in the Garching models and reaches peak amplitudes of more than s (Eric Lentz, private communication 2015), which is again in the ballpark of the dipole amplitude in the Garching simulations (see Figure 5). However, in the Oak Ridge model the monopole remains much larger and does not decay below s until 450 ms after bounce, which is a clear difference to the Garching simulations but might be connected to the use of different transport treatments and neutrino opacities by the two groups.
Despite being assuring, the caveat is that all of these 3D simulations were performed with RbR(+) approximations for the neutrino transport. Stronger evidence for the physical nature of LESA requires its confirmation by models with true multidimensional transport treatments.
4.2.6 Consequences of LESA
The LESA phenomenon implies a larger excess of the emission compared to the emission in one hemisphere than in the other. This asymmetry persists until the onset of the explosion and beyond (see models z9.6 and s20, Figures 4 and 5). It leads to different neutrontoproton ratios (and ) in the supernova ejecta expelled in different directions, because in the expanding neutrinoheated matter is set by the competition of and absorptions on neutrons and protons, respectively. The consequences for supernova nucleosynthesis remain to be explored.
It is presently not known how long into the protoneutron star cooling evolution the strong dipole component of the neutrino emission will survive. If a stable dipole radiates the total gravitational binding energy of the newly formed neutron star with an asymmetry of the order of several percent, the corresponding recoil acceleration could explain the observed natal kick velocities of pulsars up to more than 1000 km s. In the 9.6 explosion (the 3D model z9.6) we estimate a (mainly neutrinoinduced) kick of 35 km s. within the simulated time of 450 ms after bounce. Though this is not particularly big, even such a modest kick magnitude would have important astrophysical implications if it defines a lower bound because LESA is a generic phenomenon that plays a role in all newborn neutron stars.
The LESA asymmetry is also of great importance for detailed predictions of the observable neutrino signal from a future Galactic supernova. While a viewingangle dependence of the radiated neutrino luminosities adds complexity to the interpretation of a measured neutrino burst, the picture will become even more complicated if the directiondependent  flux asymmetry leads to selfinduced neutrinoflavor conversions that vary with the observer position [137].
With all these questions being barely scratched, the LESA defines yet another interesting territory that deserves further exploration as simulations penetrate even deeper into the landscape of 3D phenomena happening inside of supernova cores.
Disclosure Statement
The authors are not aware of any affiliations, memberships, funding, or financial holdings that might be perceived as affecting the objectivity of this review.
Acknowledgments
We thank Ewald Müller and Georg Raffelt for illuminating discussions, Evan O’Connor and Eric Lentz for communicating their LESA results to us, Aaron Döring (MPA), Elena Erastova, and Markus Rampp (MPCDF) for visualization support of the 3D simulations, and Georg Stockinger for Figure 5. This research was supported by the Deutsche Forschungsgemeinschaft through grant EXC 153 and by the European Research Council through grant ERCAdG No. 341157COCO2CASA. We acknowledge computing time from the European PRACE initiative on SuperMUC (GCS@LRZ, Germany), Curie (GENCI@CEA, France), and MareNostrum (BSCCNS, Spain). Postprocessing of simulation results was done on Hydra of the Max Planck Computing and Data Facility (MPCDF).
Literature Cited
 [1] Arnett WD, Bahcall JN, Kirshner RP, Woosley SE. Ann. Rev. Astron. Astrophys. 27:629 (1989)
 [2] Hillebrandt W, Hoflich P. Reports on Progress in Physics 52:1421 (1989)
 [3] Maeda K, et al. Science 319:1220 (2008)
 [4] Aschenbach B, Egger R, Trümper J. Nature 373:587 (1995)
 [5] DeLaney T, et al. Astrophys. J. 725:2038 (2010)
 [6] Milisavljevic D, Fesen RA. Astrophys. J. 772:134 (2013)
 [7] Milisavljevic D, Fesen RA. Science 347:526 (2015)
 [8] Grefenstette BW, et al. Nature 506:339 (2014)
 [9] Boggs SE, et al. Science 348:670 (2015)
 [10] Cordes JM, Chernoff DF. Astrophys. J. 505:315 (1998)
 [11] Lai D, Chernoff DF, Cordes JM. Astrophys. J. 549:1111 (2001)
 [12] Arzoumanian Z, Chernoff DF, Cordes JM. Astrophys. J. 568:289 (2002)
 [13] Hobbs G, Lorimer DR, Lyne AG, Kramer M. Mon. Not. R. Astron. Soc. 360:974 (2005)
 [14] Janka HT. Annual Review of Nuclear and Particle Science 62:407 (2012)
 [15] Janka HT, et al. Progress of Theoretical and Experimental Physics 2012:01A309 (2012)
 [16] Kotake K, et al. Progress of Theoretical and Experimental Physics 2012:01A301 (2012)
 [17] Kotake K, et al. Advances in Astronomy 2012:428757 (2012)
 [18] Burrows A. Reviews of Modern Physics 85:245 (2013)
 [19] Foglizzo T, et al. Publ. Astron. Soc. Australia 32:e009 (2015)
 [20] Foglizzo T. Astron. Astrophys. 368:311 (2001)
 [21] Foglizzo T. Astron. Astrophys. 392:353 (2002)
 [22] Blondin JM, Mezzacappa A, DeMarino C. Astrophys. J. 584:971 (2003)
 [23] Akiyama S, Wheeler JC, Meier DL, Lichtenstadt I. Astrophys. J. 584:954 (2003)
 [24] Liebendörfer M, Mezzacappa A, Thielemann FK. Phys. Rev. D 63:104003 (2001)
 [25] Liebendörfer M, et al. Astrophys. J. Suppl. Ser. 150:263 (2004)
 [26] Burrows A, et al. Astrophys. J. 539:865 (2000)
 [27] Rampp M, Janka HT. Astron. Astrophys. 396:361 (2002)
 [28] Mikheyev SP, Smirnov AY. Yadernaya Fizika 42:1441 (1985)
 [29] Wolfenstein L. Phys. Rev. D 17:2369 (1978)
 [30] Wolfenstein L. Phys. Rev. D 20:2634 (1979)
 [31] Hannestad S, Janka HT, Raffelt GG, Sigl G. Phys. Rev. D 62:093021 (2000)
 [32] Mirizzi A, et al. La Rivista del Il Nuovo Cimento 39:1 (2016)
 [33] Foglizzo T, Masset F, Guilet J, Durand G. Physical Review Letters 108:051103 (2012)
 [34] Fryer CL, Warren MS. Astrophys. J. 574:L65 (2002)
 [35] Fryer CL, Warren MS. Astrophys. J. 601:391 (2004)
 [36] Fryer CL, Young PA. Astrophys. J. 659:1438 (2007)
 [37] Herant M, et al. Astrophys. J. 435:339 (1994)
 [38] Fryer CL. Astrophys. J. 522:413 (1999)
 [39] Colgate SA, White RH. Astrophys. J. 143:626 (1966)
 [40] Arnett WD. Canadian Journal of Physics 44:2553 (1966)
 [41] Bethe HA, Wilson JR. Astrophys. J. 295:14 (1985)
 [42] Fryer CL, Kalogera V. Astrophys. J. 554:548 (2001)
 [43] Ugliano M, Janka HT, Marek A, Arcones A. Astrophys. J. 757:69 (2012)
 [44] Ertl T, et al. ArXiv eprints 1503.07522 (2015)
 [45] Sukhbold T, et al. ArXiv eprints 1510.04643 (2015)
 [46] Papish O, Soker N. Mon. Not. R. Astron. Soc. 416:1697 (2011)
 [47] Papish O, Soker N. Mon. Not. R. Astron. Soc. 438:1027 (2014)
 [48] Kushnir D, Katz B. Astrophys. J. 811:97 (2015)
 [49] Bruenn SW, et al. Astrophys. J. 767:L6 (2013)
 [50] Bruenn SW, et al. ArXiv eprints 1409.5779 (2014)
 [51] Lentz EJ, et al. Astrophys. J. Lett. 807:L31 (2015)
 [52] Summa A, et al. ArXiv eprints 1511.07871 (2015)
 [53] O’Connor E, Couch S. ArXiv eprints 1511.07443 (2015)
 [54] Dolence JC, Burrows A, Zhang W. Astrophys. J. 800:10 (2015)
 [55] Skinner MA, Burrows A, Dolence JC. ArXiv eprints 1512.00113 (2015)
 [56] Lentz EJ, et al. Astrophys. J. 747:73 (2012)
 [57] Müller B, Janka HT, Marek A. Astrophys. J. 756:84 (2012)
 [58] Kuroda T, Kotake K, Takiwaki T. Astrophys. J. 755:11 (2012)
 [59] Suwa Y, Yamada S, Takiwaki T, Kotake K. Astrophys. J. 816:43 (2016)
 [60] Nakamura K, Takiwaki T, Kuroda T, Kotake K. Publ. Astron. Soc. Japan 67:107 (2015)
 [61] Liebendörfer M, Whitehouse SC, Fischer T. Astrophys. J. 698:1174 (2009)
 [62] Pan KC, Liebendörfer M, Hempel M, Thielemann FK. Astrophys. J. 817:72 (2016)
 [63] Müller B, Janka HT, Heger A. Astrophys. J. 761:72 (2012)
 [64] Müller B, Janka HT, Marek A. Astrophys. J. 766:43 (2013)
 [65] Müller B, Janka HT. Astrophys. J. 788:82 (2014)
 [66] Marek A, Janka HT. Astrophys. J. 694:664 (2009)
 [67] Hanke F, Marek A, Müller B, Janka HT. Astrophys. J. 755:138 (2012)
 [68] Kolmogorov A. Akademiia Nauk SSSR Doklady 30:301 (1941)
 [69] Landau LD, Lifshit’s EM. Reading, MA: AddisonWesley (1959)
 [70] Kraichnan RH. Physics of Fluids 10:1417 (1967)
 [71] Blondin JM, Mezzacappa A. Nature 445:58 (2007)
 [72] Iwakami W, et al. Astrophys. J. 700:232 (2009)
 [73] Fernández R. Astrophys. J. 725:1563 (2010)
 [74] Hanke F, et al. Astrophys. J. 770:66 (2013)
 [75] Abdikamalov E, et al. Astrophys. J. 808:70 (2015)
 [76] Nakamura K, Kuroda T, Takiwaki T, Kotake K. Astrophys. J. 793:45 (2014)
 [77] Couch SM. Astrophys. J. 775:35 (2013)
 [78] Dolence JC, Burrows A, Murphy JW, Nordhaus J. Astrophys. J. 765:110 (2013)
 [79] Handy T, Plewa T, Odrzywołek A. Astrophys. J. 783:125 (2014)
 [80] Nordhaus J, Burrows A, Almgren A, Bell J. Astrophys. J. 720:694 (2010)
 [81] Fernández R. Mon. Not. R. Astron. Soc. 452:2071 (2015)
 [82] Radice D, et al. ArXiv eprints 1510.05022 (2015)
 [83] Ott CD, et al. Astrophys. J. 768:115 (2013)
 [84] Couch SM, O’Connor EP. Astrophys. J. 785:123 (2014)
 [85] Takiwaki T, Kotake K, Suwa Y. Astrophys. J. 749:98 (2012)
 [86] Takiwaki T, Kotake K, Suwa Y. Astrophys. J. 786:83 (2014)
 [87] Scheck L, Kifonidis K, Janka HT, Müller E. Astron. Astrophys. 457:963 (2006)
 [88] Wongwathanarat A, Janka H, Müller E. Astrophys. J. 725:L106 (2010)
 [89] Müller E, Janka HT, Wongwathanarat A. Astron. Astrophys. 537:A63 (2012)
 [90] Wongwathanarat A, Janka HT, Müller E. Astron. Astrophys. 552:A126 (2013)
 [91] Wongwathanarat A, Müller E, Janka HT. Astron. Astrophys. 577:A48 (2015)
 [92] Müller B, Janka HT. Mon. Not. R. Astron. Soc. 448:2141 (2015)
 [93] Müller B. Mon. Not. R. Astron. Soc. 453:287 (2015)
 [94] Tamborra I, et al. Physical Review Letters 111:121104 (2013)
 [95] Tamborra I, et al. Phys. Rev. D 90:045032 (2014)
 [96] Tamborra I, et al. Astrophys. J. 792:96 (2014)
 [97] Melson T, Janka HT, Marek A. Astrophys. J. Lett. 801:L24 (2015)
 [98] Melson T, et al. Astrophys. J. Lett. 808:L42 (2015)
 [99] Kuroda T, Takiwaki T, Kotake K. ArXiv eprints 1501.06330 (2015)
 [100] Burrows A, Hayes J, Fryxell BA. Astrophys. J. 450:830 (1995)
 [101] Janka HT, Müller E. Astron. Astrophys. 306:167 (1996)
 [102] Foglizzo T, Scheck L, Janka HT. Astrophys. J. 652:1436 (2006)
 [103] Fernández R, Müller B, Foglizzo T, Janka HT. Mon. Not. R. Astron. Soc. 440:2763 (2014)
 [104] Scheck L, Janka HT, Foglizzo T, Kifonidis K. Astron. Astrophys. 477:931 (2008)
 [105] Iwakami W, et al. Astrophys. J. 678:1207 (2008)
 [106] Foglizzo T, Galletti P, Scheck L, Janka HT. Astrophys. J. 654:1006 (2007)
 [107] Guilet J, Foglizzo T. Mon. Not. R. Astron. Soc. 421:546 (2012)
 [108] Kazeroni R, Guilet J, Foglizzo T. Mon. Not. R. Astron. Soc. 456:126 (2016)
 [109] Scheck L, et al. Physical Review Letters 92:011103 (2004)
 [110] Guilet J, Fernández R. Mon. Not. R. Astron. Soc. 441:2782 (2014)
 [111] Buras R, Janka HT, Rampp M, Kifonidis K. Astron. Astrophys. 457:281 (2006)
 [112] Thompson TA, Quataert E, Burrows A. Astrophys. J. 620:861 (2005)
 [113] Murphy JW, Burrows A. Astrophys. J. 688:1159 (2008)
 [114] Burrows A, Dolence JC, Murphy JW. Astrophys. J. 759:5 (2012)
 [115] Iwakami W, Nagakura H, Yamada S. Astrophys. J. 793:5 (2014)
 [116] Scheck L. 2006. Parametric Studies of Hydrodynamik Instabilities in the Supernova Core by Two and Threedimensional Simulations. PhD Thesis, Technische Universität München
 [117] Couch SM, Ott CD. Astrophys. J. 799:5 (2015)
 [118] Radice D, Couch SM, Ott CD. Computational Astrophysics and Cosmology 2:7 (2015)
 [119] Murphy JW, Dolence JC, Burrows A. Astrophys. J. 771:52 (2013)
 [120] Murphy JW, Meakin C. Astrophys. J. 742:74 (2011)
 [121] Couch SM, Ott CD. Astrophys. J. 778:L7 (2013)
 [122] Summa M, Janka HT, Melson T, Marek A. in preparation (2016)
 [123] Wongwathanarat A, Hammer NJ, Müller E. Astron. Astrophys. 514:A48 (2010)
 [124] Hobbs TJ, Alberg M, Miller GA. ArXiv eprints 1601.01729 (2016)
 [125] Burrows A, Goshy J. Astrophys. J. 416:L75 (1993)
 [126] Pejcha O, Thompson TA. Astrophys. J. 746:106 (2012)
 [127] Murphy JW, Dolence JC. ArXiv eprints 1507.08314 (2015)
 [128] Fernández R. Astrophys. J. 749:142 (2012)
 [129] Janka HT. Astron. Astrophys. 368:527 (2001)
 [130] Lund T, et al. Phys. Rev. D 82:063007 (2010)
 [131] Marek A, Janka HT, Müller E. Astron. Astrophys. 496:475 (2009)
 [132] Ott CD, Burrows A, Dessart L, Livne E. Astrophys. J. 685:1069 (2008)
 [133] Brandt TD, Burrows A, Ott CD, Livne E. Astrophys. J. 728:8 (2011)
 [134] Yamasaki T, Foglizzo T. Astrophys. J. 679:607 (2008)
 [135] Lund T, et al. Phys. Rev. D 86:105031 (2012)
 [136] Stockinger G. 2015. Features of LESA – A Study on the LeptonNumber Emission SelfSustained Asymmetry as a new Phenomenon in Supernova Explosions. Master Thesis, Technische Universität München
 [137] Chakraborty S, Raffelt G, Janka HT, Müller B. Phys. Rev. D 92:105002 (2015)
 [138] Arnett WD, Meakin C. Astrophys. J. 733:78 (2011)
 [139] Couch SM, Chatzopoulos E, Arnett WD, Timmes FX. Astrophys. J. 808:L21 (2015)
 [140] Nagakura H, Sumiyoshi K, Yamada S. Astrophys. J. Suppl. 214:16 (2014)
 [141] Sumiyoshi K, Takiwaki T, Matsufuru H, Yamada S. Astrophys. J. Suppl. 216:5 (2015)