From Supernova to Supernova Remnant: Comparison of Thermonuclear Explosion Models

Our model's basic properties are given in Table 1. All models consider the explosion of a single accreting WD (single degenerate model) of Chandrasekhar mass M = 1.4 M_⊙, initial radius R = 1.96 × 10⁸ cm, and homogeneous composition of ¹²C and ¹⁶O in equal parts by mass. The zero-age main-sequence progenitor is assumed to be of solar metallicity, with electron fraction Y_e = 0.49886.¹¹ The four SN models we have selected cover two explosion mechanisms: there are two delayed detonation models (ddt) from Seitenzahl et al. (2013b) and two deflagration models (def) from Fink et al. (2014).¹² We also compare two numerical setups for the central ignition: two models have many ignition points (N100), leading to rather spherical ejecta, while two models have a few ignition points (N5), leading to more asymmetric ejecta. We note that the N numbers are not to be taken too literally; rather, the ignition spots should be thought of as a convenient parameterization of the initial asymmetry and filling factor of the flame surface. The ignition physics in such highly turbulent flows is not well understood; the most advanced simulations to date find that multiple ignition is unlikely and off-center ignition is likely (Zingale et al. 2009; Nonaka et al. 2012). Our ambition here is not to create a complete mapping between such ignition models and their parameterizations. In this paper, we focus on assessing the impact of existing parameterized ignition models on the eventual morphology of the SNRs.

Table 1.  Summary Data for the Explosion Models

	Ekin (erg)	Mej (M⊙)	${v}_{\max }$ (km s−1)	Composition
N100ddt	1.43 × 1051	1.40	28,700	60% IGE, 30% IME
N5ddt	1.55 × 1051	1.40	27,900	80% IGE, 15% IME
N100def	6.15 × 1050	1.35	14,100	40% IGE, 10% IME
N5def	1.35 × 1050	0.37	13,800	60% IGE,  10% IME

Note. Columns show the explosion's (asymptotic) kinetic energy E_kin, ejected mass M_ej, maximum ejecta velocity ${v}_{\max }$ (100 s after the ignition), and global composition of the ejecta (as an approximate fraction of the unbound mass) where IGE is iron group elements and IME is intermediate group elements (see detailed yields in Seitenzahl et al. (2013b) for the ddt models and in Fink et al. (2014) for the def models).

Download table as:  ASCIITypeset image

In the grid of ddt models from Seitenzahl et al. (2013b), N100ddt produces an SN of normal brightness from the synthesized mass of ⁵⁶Ni ≃ 0.6 M_⊙, and it reproduces well key spectral characteristics of normal SNe Ia (Sim et al. 2013). This is thus the model we first picked, in Paper I, to make simulations of the SNR phase. It is a rather symmetric explosion (Bulla et al. 2016). N5ddt was selected next, because it has fewer ignition kernels and so is more asymmetric. The fewer the ignition kernels, the weaker the initial deflagration, and the more fuel is processed by the detonation, therefore N5ddt is also brighter than N100ddt. Both N100ddt and N5ddt were considered in Williams et al. (2017) for the interpretation of the ejecta distribution in Tycho's SNR. The N3ddt model, which is similar to N5ddt, was discussed in Borkowski et al. (2013) for the SNR G1.9+0.3. From Sim et al. (2013) the best matches for the low-N ddt models are with 1991T-like SNe Ia, although the agreement is not great in detail.

The models in the grid of Fink et al. (2014) are the same as in Seitenzahl et al. (2013b), the only difference being that the DDT is not allowed. Whereas all our ddt models fully unbind the WD, def models are often partial explosions that leave compact remnants behind, that is, a WD with an unusual composition—which may have been observed, see Vennes et al. (2017). A weak explosion such as N5def is thought to be a good model for SN 2002cx-like SNe Ia (Kromer et al. 2013; see also Jordan et al. 2012 for a similar model) and maybe for the larger class of Type Iax SNe (Bulla et al. 2020). N5def is an asymmetric explosion that leaves behind a  ≈1 M_⊙ WD remnant. Compared with N5ddt, the ejected mass is 3.8 times lower while the kinetic energy is 11.5 times lower; the maximum velocity reached is about halved. For completion, we also include the N100def model, although it appears that the strongest def events cannot explain any observed class of SN Ia. Contrary to N5def, N100def is almost all unbound, just a bit less energetic than the ddt version. Compared with N100ddt, the ejected mass is almost the same, while the energy is 2.3 times lower, and the maximum velocity is again about halved.

The SN explosion simulations were made in 3D using the hydrodynamic code LEAFS (Reinecke et al. 2001). Each simulation was done in two steps. First, a hydro run with a restricted nuclear reaction network, with only the species needed to get the energy budget right: carbon, oxygen, a proxy for intermediate-mass isotopes, nickel, and alpha particles. Then the computation of the nucleosynthesis in post-processing, using one million tracer particles and a full nuclear reaction network with 384 isotopes of elements up to Z = 32. Each model was simulated until a time of 100 s, well after the end of nuclear reactions, and when the ejecta are observed to have settled into the self-similar free expansion phase. The final result has been remapped to Cartesian grids of resolution 200³.

Figure 1 shows the angularly averaged radial profile of the density for the four models. The averaging is done from the explosion center, from which everything expands. The average radial density profiles of the ddt and def models are different. The average density of ddt models is rather uniform in the core and then decreases exponentially, whereas for def models it peaks toward the center (at low velocities) and falls more rapidly at the edge. The obtained profiles are compared with simple analytical profiles commonly used in the SNR community: a power-law profile (e.g., Decourchelle & Ballet 1994) and an exponential profile (Dwarkadas & Chevalier 1998), neither of which is a very good fit. Anyway, such a 1D density profile does not capture the complicated structure of the ejecta, which stems from the dynamics of the deflagration and/or detonation fronts. As shown in Paper I, when reducing the full 3D profile to a single radial profile one does not get the same SNR evolution.

Zoom In Zoom Out Reset image size

The output of the SN simulations described in the previous section was used as input for the SNR simulations described in this section. See Paper I for details. For the def models, the bound remnant has been excised for the simulation of the SNR phase. As explained in Fink et al. (2014), to determine the unbound parts we calculated the asymptotic specific kinetic energy _kin,a = _kin + _grav; the condition to be unbound is _kin,a > 0, then matter will approach velocity ${v}_{{\rm{a}}}=\sqrt{2{\epsilon }_{\mathrm{kin},{\rm{a}}}}$. During the initial self-similar free expansion phase, r ∝ t and ρ ∝ t⁻³. To save on computing time, the ejecta profiles are rescaled from t = 100 s to t = 1 day, and the hydro simulation starts at 1 day. It is observed that the assumption of self-similarity is actually reasonable up to about a year of evolution. We note that we are not concerned here with the SN phase per se, or the modeling of the light curve, although we checked the effect of early energy deposition from radioactive decay on the dynamics (see Section 4.1). The ambient density is assumed to be uniform for simplicity, and to be low: n_ISM = 0.1 cm⁻³ to have overall dynamics similar to Tycho's SNR. Besides this value, the evolution of the SNR depends on the explosion energy ${E}_{\mathrm{SN}}$ and ejecta mass M_ej (Truelove & McKee 1999). From these three quantities one can define three characteristic scales; we are using here the same definitions as in Dwarkadas & Chevalier (1998) and Warren & Blondin (2013):

$\begin{eqnarray}&&{r}_{\mathrm{ch}}={\left(\displaystyle \frac{3{M}_{\mathrm{ej}}}{4\pi {\rho }_{\mathrm{ISM}}}\right)}^{1/3},\end{eqnarray} \tag{ 1 }$

$\begin{eqnarray}&&{u}_{\mathrm{ch}}={\left(\displaystyle \frac{2{E}_{\mathrm{SN}}}{{M}_{\mathrm{ej}}}\right)}^{1/2},\end{eqnarray} \tag{ 2 }$

$\begin{eqnarray}&&{t}_{\mathrm{ch}}=\displaystyle \frac{{r}_{\mathrm{ch}}}{{u}_{\mathrm{ch}}}.\end{eqnarray} \tag{ 3 }$

The hydrodynamic evolution of the SNR is computed up to 500 yr after the explosion (about the age at which we observe Tycho's SNR).

The SNR hydro simulation is performed using a custom version of the hydrodynamics code RAMSES (Teyssier 2002). Since we are interested in the morphology of the remnant, we conduct simulations in 3D. They are performed on a Cartesian grid of resolution 256³. The SN grid is loaded at 1:1, leaving room for the blast wave developing ahead of the ejecta. The grid is co-expanding with the remnant, so that the physical size of the box is changing over time, while the relative resolution of the cells remains constant. Our purpose is to study the morphology of the remnant, independent of its physical size (which for a given SNR at a given age requires knowledge of its distance, which is usually very uncertain).

We focus on the shocked region, which is the part revealed by X-ray observations. To make visible its 3D structure, and to quantify its evolution, we proceed as described in Paper I. At any given time, we extract in the simulation box the surface of the wave fronts: the RS in the ejecta, the CD between the shocked ejecta and the shocked ISM, and the FS in the ISM. We record the radius of each surface, measured from the explosion center, as a function of the direction (θ, ϕ). These data are stored and analyzed using the HEALPix package (Gorski et al. 2005). We interpret these radii (or functions thereof) as functions R(θ, ϕ) on the sphere, which we expand in spherical harmonics ${Y}_{{\ell }}^{m}(\theta ,\phi )$ to obtain their angular power spectrum C_ℓ. The complete set of plots from this analysis, for all models and all wave fronts, is available in an online repository.¹³ We checked that the reconstructed maps are indistinguishable from the original ones by eye, and that the residuals have negligible levels for our analysis (these ancillary plots are present in the repository). The residuals show random variations at the pixel scale, plus some large-scale patterns and/or some hot spots, which may or may not have an obvious counterpart. Across the four models and the 500 years of evolution, the relative residuals have zero mean, standard deviation below 1% for the CD and below 0.1% for the shocks, and local peaks at up to 8% for the CD and 3% for the shocks. The angular spectra presented in Section 3.2 capture quantitatively most of what can be seen from the simulation results.¹⁴ They allow us to see the different evolution of SN-phase modes and SNR-phase modes (see also Paper I) and to make comparisons between the different models (new to this paper).

Furthermore, each SNR simulation is made with two setups: one using the actual SN profile (labeled 3Di) and one using an angularly averaged version of it, so effectively one-dimensional, depending on radius only (labeled 1Di). The results from the 1Di simulations will not all be shown in this paper, but are used as a reference to interpret the results of the 3Di simulations, as well as to assess the possible presence of numerical artefacts (see Paper I for a systematic comparison for the N100ddt model).

Slices. In Figures 3–6 we show slices of the mass density for each model, to show the inner structure of the remnant. On each plot we show a slice in the midplane along three directions (the principal axes of the grid x, y, z), to appreciate the morphological variations, and a slice at three different ages, to appreciate the time evolution: at 1 yr the morphology is still very similar to the initial conditions, at 100 yr the characteristic shell of the SNR phase is visible while the imprint of the SN is still clear, at 500 yr the shell-like structure is undergoing regularization. Movies from 1 yr to 500 yr in steps of 1 yr are also available. In all cases, we can see the competition between two effects shaping the young SNR: the progressive regularization of the SN ejecta and the growth of the RTI. After a few hundred years, all the SN models, despite differences in explosion energy and ejecta distribution, tend to the standard shell structure expected for a young SNR.

For N100ddt, as previously reported, by 500 yr the SNR has regularized in a rather spherical shell, although with some disturbances compared with the case of smooth initial conditions (see the detailed comparison of the 1Di and 3Di setups in Paper I). The ejecta initially have a shape that, although complicated, does not exhibit a preferred direction (this is quantified in Section 3.2); as a result the shell structure at 500 yr looks statistically similar when observed at different angles. The situation is different for N5ddt. This model initially has an elongated structure, which over time leads to a two-sided SNR, with one half of the shell being almost twice as thick as the other half. Looking more carefully at the time evolution, we see that the FS is located at different radii from the very beginning, with the layer of swept-up ISM initially growing uniformly. The more elongated parts of the ejecta are also the less dense, and so are the parts where the RS travels faster inward; this is the main reason for the growing spread between the RS and the FS. The net result at 500 yr is that both the RS and the FS are nearly circular, but the RS is actually centered on the explosion center whereas the FS is off-center. The RT fingers grow in a similar way over the CD, and so are more likely to interact with the FS on the thin side than on the thick side.

This asymmetry between the N100 and N5 models is also visible for the def cases, although it expresses itself in different ways. Because of the way the deflagration fronts propagate, def models tend to have more regular initial contours with a more uniform density distribution over angles. This can be quantified, in a given slice, by comparing the perimeter of the CD with the perimeter of a circle that has a radius equal to the average radius of the CD. At 1 yr this ratio is 5.7 ± 0.5 for N100ddt, 4.8 ± 0.5 for N5ddt, 2.2 ± 0.0 for N100def, and 2.6 ± 0.1 for N5def. So from the start def models have a contour twice as regular as ddt models. As a result def models produce rather symmetric SNRs. But whereas N100def has a fairly isotropic distribution, from the start and so at all times, N5def is clearly biased on one side, leading to an off-center SNR. Compared to N5ddt, the ejecta density is more uniform and so the RS travels in a more symmetric way, leading to a rather balanced shell. But now both shocks, and so the entire SNR shell, are completely off-centered with respect to the explosion center, at all times. Another remarkable property of def models is the presence of an overdensity at the center. This stems from the fact that pure deflagration models fail to fully unbind the WD. A compact remnant remains at the center of the peak. We remind the reader that this remnant is not included in the SNR simulations, only the unbound ejecta are taken into account; we add a word of caution that the cut may not be perfectly modeled, although the presence of low-velocity ejecta is a real feature. For the N5def model, for which only a quarter of the WD mass is ejected, we even see a shock front near the center, which originates in the sloshing motion of the matter that is marginally unbound (Fink et al. 2014). By 500 yr the RS is still far from these structures. We checked that at 1000 yr, the RS has just started its motion inward (in the ISM frame), but has barely reached the central overdense region for the def models. Semianalytical models (assuming spherically symmetric ejecta with a power-law radial density profile), as implemented in the calculator by Leahy & Williams (2017), predict that it will take up to 3000 yr for the RS to reach the center with our range of E_kin and M_ej. For the def models, especially N100def, we also observe dense filaments radiating from the center. Looking back at the time evolution of the N100def SN simulation, we see how these structures are formed as the deflagration develops (see Figures 21 and 22 in the Appendix). These overdense channels are created between large-scale plumes of the deflagration ashes. Being made of unburnt material, they are rich in C and O.

Projections. In Figures 7–10 we show for each model projections of the density squared of the ejecta in the shocked region, which is a proxy for the (broadband) thermal emission from the hot plasma (which falls in the X-ray domain, around keV energies). Again we show for each model the map along three orthogonal directions x, y, z, and at three different times 1 yr, 100 yr, and 500 yr. Movies from 1 yr to 500 yr in steps of 1 yr are also available. Note that, since the emission from the plasma depends on its electron temperature and ionization state, which are time-dependent, a more precise calculation is needed to produce realistic X-ray maps, which we defer to a future study. These maps allow us to assess, at a glance, what kind of morphological differences can be expected on observations of SNRs from the four SN models.

For N100ddt, as in Paper I, maps at 500 yr are rather symmetric, although not as much as simulations made from smooth initial conditions, and so depend somewhat on the direction of observation. The edge of the ejecta generally looks brighter, as expected in projection (limb brightening effect). In comparison, the N5ddt model displays a distinctive two-sided structure, with one half of the SNR significantly brighter than the other (of course, the effect would not be visible if looking right along the dipolar axis). This is a direct consequence of the varying thickness of the shocked region.

The def models, N100def and N5def, have a more angularly symmetric shape in projection (although they may not look exactly the same from every direction). Differences between the N100 and N5 versions are not as obvious as in the ddt case, except for the fact that the N5def SNR is significantly off-center, but that would not be inferred from the map alone without additional information (such as measurements of the radial motion of the ejecta, pointing to the explosion center). A distinctive feature of the def models, present for both N100 and N5 versions, is the appearance of bright lines that cross the SNR surface. We will see them in another form in the next section. Finally, we note that, because these maps in projection show only the matter bounded by the RS and FS, they cannot feature the central overdensity of the def models (otherwise this would be the dominant feature).

Edge of the ejecta. In the left panels of Figures 11–14 we show in a Mollweide projection the surface of the CD (edge of the ejecta). The function plotted is the relative variation in radius from the explosion center, R(θ, ϕ) = (r_CD − 〈r_CD〉)/〈r_CD〉. Regions in tones of red are ahead of the average (R > 0), while regions in tones of blue are lagging (R < 0). Unlike the maps of the previous section, these maps contain information in all directions at once. On the right panels we show the power spectrum resulting from the expansion in spherical harmonics of the function R(θ, ϕ), as a function of angular wavenumber ℓ. The projections allow us to quickly visually assess the shape of the ejecta, while the spectra allow us to quantify the asymmetries. The power C_ℓ plotted is normalized in such a way that each gray bin is the contribution of wavenumber ℓ to the total variance of R(θ, ϕ). As before, the plots are shown at three different times: 1 yr, 100 yr, and 500 yr. Movies from 1 yr to 500 yr in steps of 1 yr are also available.

The CD surface at 1 yr highlights the morphological differences between the different SN models. First we see the difference between N100 models, for which power peaks at intermediate angular scales, and N5 models, for which power peaks at the dipole component ℓ = 1. For N5ddt the dipole takes the form of a protrusion: a relatively small area is ahead of average by 50%, while a larger area is lagging by half this amount. For N5def the dipole is more balanced, with roughly two opposite halves. For N5ddt the dipole strength decreases over time via the mechanism that we have discussed previously; the ejecta become more regular over time, although an asymmetry is still visible at 500 yr. For N5def the dipole remains at all times, as the SNR is effectively off-center as we have seen. Second we see differences between ddt and def models, with the latter having weaker asymmetries (note that the scale for the projection maps has been set to the same value for the ddt and def versions of a given ignition setup).

On top of these large-scale effects, which stem from the chosen SN model, we see the RTI steadily growing over time, at smaller scales. After a couple of hundred years it competes with the initial conditions in shaping the remnant, and at 500 yr it controls most of what is immediately visible by eye.

Another effect, a distinctive feature of the def models, which was not expected, is the progressive appearance over time of contiguous lines, consistently ahead of average. They are particularly clear on the maps for the N100def model (Figure 13), even though they contribute little power. These protruding structures match the lines previously reported on the projected density maps (Figures 9 and 10). They exhibit an irregular pattern that does not match the geometry of our simulation mesh. Looking at the time evolution, we see that some of them grow on top of pre-existing ejecta ridges, while some grow at the bottom of ejecta trenches. They have dynamics completely separate from that of the RTI. After a few hundred years they form a network of seam lines on the ejecta surface, separating the ejecta into patches.

Width of the shocked region. In this paper we are not showing the individual plots for the RS and FS, which show fewer variations than for the CD (these plots are present in the online repository). As shown in Paper I, the RS, being always close to the CD for such young SNRs, shows some echo of the CD fluctuations (feet of the RTI). The FS initially delineates the ejecta surface, and over time detaches from it and becomes more regular, although it may still show some imprint of the turbulent interior (fingers of the RTI). These effects are visible for each of the SN models, to varying extent.

Instead, in Figures 15–18 we are showing the width of the shocked region, which is the radial distance between the RS and the FS, normalized by the average size of the SNR, which is taken to be the average radius of the FS. So explicitly we consider the function R(θ, ϕ) = (r_FS − r_RS)/〈r_FS〉. We now always have R > 0, and we are using tints of red on the map. For each model, we see the growth of the shocked region over time. Initially, it naturally bears the imprint of the shape of the ejecta's edge shown above, from which it emerges. As the shocks travel away, and their surface progressively regularizes through the effect of the high pressure, the shocked region becomes globally more homogeneous in time—except for N5ddt. For all models, small-scale variations visible at 500 yr are mostly caused by a few strongly growing RT fingers that are pushing the FS. At this age, the shocked region comprises about 20%–50% of the radial size of the SNR, with the thickest models being N5ddt (on one side) and N5def (on most of its surface).

In line with what we have seen before on slices, N5ddt has a clearly bipolar morphology, with the shocked region about twice as large on one side as on the other. For the N5 models, comparing the plots for the CD (Figures 12 and 14) and for FS–RS (Figures 16 and 18) reveals that N5ddt and N5def are asymmetric in completely different ways. For N5ddt, the ejecta surface is initially dipolar, but regularizes over time, while the shocked region follows the opposite trend. As explained before, this stems from the fact that the most extended ejecta are also the most diluted, and so the most easily decelerated. This ends up being a two-sided SNR. For N5def, the ejecta are skewed on one side at all times, but the shocked region quickly becomes pretty uniform. This makes a regular-looking but offset SNR.

Finally, for the def models, we note that the seams noticed on the ejecta surface are also present. They are particularly visible around 100 yr. They are mostly caused by the fact that the RS is at larger radii than average in these regions—meaning that it is traveling inward more slowly than average. So the picture that emerges is that these structures are caused by ejecta that, from the beginning, have greater momentum than average: these ejecta move faster, hence the protruding lines on the ejecta surface (Figures 13 and 14), and thus they generate a weaker RS behind them, hence the narrow lines on the width of the shocked region (Figures 17 and 18). Since the initial SN profiles are homologous, differences in momentum at a given radius come from differences in density; indeed the seams show up on the maps of the projected density squared (Figures 9 and 10). For N100def, a closer inspection of the initial mass density at 100 s reveals that a network of overdense ridges is indeed present in the outer ejecta layers, below the ejecta surface. A 3D visualization of this structure is available online.¹⁵ Its origin appears to be connected with the dynamics of the deflagration fronts during the SN phase, previously mentioned and shown in the Appendix. We interpret the network of overdense ejecta as the outer boundaries between burning regions. It is this structure that eventually pierces through and becomes apparent during the SNR phase.

Angular power as a function of time. To finish the presentation of the results, we show the evolution of a global quantity: the summed angular power of the function R(θ, ϕ) = (r − 〈r〉)/〈r〉. At a given time, the power in a given range of angular scales ℓ is the area of the shaded region on the angular spectra between these scales. Its time evolution tells us the contribution of these scales in the shaping of the SNR morphology.

In Figure 19 we plot the total angular power (summed up from ℓ = 1 to ${{\ell }}_{\max }=383$) as a function of time, for the four SN models. On each plot we show the power for the three wave fronts: CD in green (corresponding to the spectra of Figures 11–14), RS in blue, and FS in red. On each plot we show the power for the actual SN model (3Di, solid lines) and for the angularly averaged version of it (1Di, dashed lines). The case of effectively 1D initial conditions only exhibits the effects due purely to the SNR phase, and so, as shown in Paper I, comparing the two setups allows us to isolate the contribution from the SN initial conditions. In the 1Di case not much is expected to happen at the RS and FS, while the curve for the CD shows the continued growth of the RTI. In the 3Di case there is power from the initial conditions. For the FS and RS it only decays in time, but for the CD it starts to decrease as the memory of the initial conditions fades then eventually increases as the RTI makes its presence felt; the long-term behavior at the CD is similar in the 1Di and 3Di cases because it is produced by the same RTI. We therefore see two regimes (except for the N5def model, for reasons explained below).

Zoom In Zoom Out Reset image size

We observed on the power spectra for the 1Di setups that the RTI does not reach below about ℓ = 10 (and peaks around ℓ = 50–60 depending on the models), so we can separate the contributions from the SN phase and the SNR phase. In Figure 20 we show separately the power at small scales (large ℓ, thin lines) and at large scales (small ℓ, thick lines). For the RS and FS this makes little difference since most of the power is at the largest scales. For the CD we see the SN modes that decay in time, versus the SNR modes that grow in time. The time at which the curves cross for the CD (for the 3Di case) is marked by a vertical line. This can serve as an estimate of the timescale at which the SN-phase turbulence and the SNR-phase instabilities have similar impact on the morphology of the SNR. At earlier times the SNR is dominated by the SN initial conditions; at later times the SNR is dominated by processes intrinsic to the SNR phase. The N5def model is a special case. We have seen that it remains off-center, and so exhibits a dipole component at all times. As a result the total power, as well as the power at small scales ℓ < 10, plateaus or even increases in time (bottom right plot of Figures 19 and 20). So we removed the ℓ = 1 contribution (only for the N5def model) in order to be able to see the decay of the other SN modes and apply our timescale-finding procedure. The transition times are 180 yr for N100ddt, 260 yr for N5ddt, 95 yr for N100def, and 280 yr for N5def. We see that the N5 models remain under the SN influence for a longer time than the corresponding N100 models. This is consistent with the N5 ignition setups being more asymmetric than the N100 ignition setups. Also the transition time is twice as large for the ddt as for the def model for the N100 case, but it is similar (and in the opposite direction) for the N5 case. So we cannot readily generalize on the effect of deflagration versus deflagration plus detonation.

Zoom In Zoom Out Reset image size

We first discuss how some effects that were not included in the simulations may affect the interpretation of our results.

Heating from radioactive decay. The most important radioactive species is ⁵⁶Ni, which powers the SN light curve. It decays via beta decay to iron in two steps: ${}_{28}^{56}\mathrm{Ni}\to {}_{27}^{56}\mathrm{Co}\to {}_{26}^{56}\mathrm{Fe}$, with half-lives of respectively 6 days and 77 days. The mass of ⁵⁶Ni synthesized in the models is 0.60 M_⊙ = 43% M_ej in N100ddt, 0.97 M_⊙ = 69% M_ej in N5ddt, 0.35 M_⊙ = 27% M_ej in N100def, and 0.16 M_⊙ = 43% M_ej in N5def. The energy released that may be converted to heat (including photons and electrons/positrons, excluding neutrinos which escape) remains modest compared to the kinetic energy of the ejecta (although somewhat less so for the models newly considered in this paper), as follows. N100ddt: 1.1 × 10⁵⁰ erg = 8% ${E}_{\mathrm{SN}}$, N5ddt: 1.8 × 10⁵⁰ erg = 12% ${E}_{\mathrm{SN}}$, N100def: 6.6 × 10⁴⁹ erg = 11% ${E}_{\mathrm{SN}}$, N5def: 2.9 × 10⁴⁹ erg =22% ${E}_{\mathrm{SN}}$. Another difference between the ddt and def models is the spatial distribution of the ⁵⁶Ni synthesized. For ddt models it is concentrated in the core, with more asymmetries for N5ddt than for N100ddt. For def models it is spread more evenly over radii.

To check the impact of heating from radioactive decay on the SNR dynamics, we re-ran the SNR simulation under the assumption of local deposition of all the energy released, which constitutes an extreme case (see Section 4.2 and Appendix C in Paper I for the details). We had not seen a significant impact for N100ddt; for the other models also we observe only mild differences on the projected maps and power spectra, so that our main results above still hold. The effect of strong heating in the ejecta core is visible on the inner density structures, but does not affect much the dynamics of the wave fronts. Further investigation of this would require diagnostics of the radioactive species still present deep within the unshocked ejecta.

Effects related to ⁵⁶Ni will be strongest in the first few years, and will terminate within a hundred years (so well within the ejecta-dominated phase). Another radioactive element of interest is ⁴⁴Ti; the mass yields and energy deposits are very small, of the order of 10⁻⁸ M_⊙ and 10⁴⁴ erg respectively, but it has a much longer lifetime of about 60 yr. Detecting ⁴⁴Ti in young Type Ia SNRs has been a longstanding objective of hard X-ray/soft γ-ray missions; see Weinberger et al. (2020) and references therein for Tycho and Kepler SNRs.

Magnetic fields and particle acceleration. In Paper I we already discussed the possible effects of the magnetic field and particle acceleration, which we do not expect to change the overall trends observed in this paper, when it comes to comparing different SN models with different ejecta morphology. That being said, we reckon that for a given SNR several physical processes may be at play at the same time at and around the shocks. We have seen that, because of the initial asymmetry of the ejecta, the RS/FS ratio may vary systematically over the SNR, and this is especially marked for the N5ddt model. The location of the shocks also depends on the compressibility of the plasma, which can be substantially affected by the presence of a nonthermal population of particles (e.g., Ferrand et al. 2010; Warren & Blondin 2013). The presence of energetic particles increases the compressibility of the fluid, lowering the effective adiabatic index, which results in a narrower (and denser) shocked region. It remains to be seen how the two effects will work together, and how our findings impact diagnostics of the efficiency of particle acceleration at SNR shock waves.

Impact of the ambient medium. We have observed that the N5 models lead to asymmetric SNRs, either with a dipolar structure for N5ddt or with an offset structure for N5def. These effects in our simulations stem from the 3D distribution of the ejecta density. Similar effects, of systematic variations of the RS/FS ratio across directions and over time, could also be obtained from gradients in the density of the ambient medium (e.g., Orlando et al. 2007). In this study we opted for the reference case of a uniform ISM. What we show is that fluctuations in the ISM are not necessary to obtain SNRs with non-canonical shapes during the ejecta-dominated phase. Older SNRs will most certainly have their shape affected by the ambient medium, which is not uniform on scales larger than a parsec. For recent simulation work on the impact of a turbulent or cloudy ISM on the development of the SNR, relevant to Type Ia SNe, see Wang et al. (2018), Peng et al. (2020), and Villagran et al. (2020).

There are only a few well-known Type Ia SNRs that are young (age  < 1000 yr) and nearby (spatially resolvable). In the Galaxy these are G1.9+0.3 (≃150 yr), the youngest known (Borkowski et al. 2013), Kepler = G4.5+6.8 = SN 1604 (416 yr), with a complicated morphology presumably from the progenitor's mass loss history (Reynolds et al. 2007; Burkey et al. 2013; Sato & Hughes 2017a), Tycho = G120.1+1.4 = SN 1572 (448 yr, Warren et al. 2005; Sato & Hughes 2017b), and SN 1006 = G327.6+14.6 (1014 yr), with a characteristic bilateral morphology from the orientation of the magnetic field (Uchida et al. 2013; Winkler et al. 2014). In the LMC we have SNR 0519-69.0, 450 ± 200 yr estimated from the shock dynamics (Kosenko et al. 2010), 600 ± 200 yr from the light echo (Rest et al. 2005), SNR 0509-67.5 (Warren & Hughes 2004), 400 ± 120 yr from the light echo, and N103B, less than 850 yr from the measured expansion (Williams et al. 2018), consistent with the light echo.

In Paper I we discussed the case of Tycho's SNR with the N100ddt model, which can explain a standard SN Ia. Maps in projection look more realistic in the 3Di setup than the 1Di setup. The SN model leads to power at large scales (small ℓ), which has been reported in observations (Warren et al. 2005), but which had not been reproduced by simulations assuming smooth initial profiles. At age 500 yr, the simulated angular spectrum could be confused with an enhanced RT growth, when really it is the combination of both the decay from the SN modes and the growth from the RTI modes.

N5ddt is the most asymmetric model, from the SN phase into the early SNR phase. The resulting SNR at 500 yr has a clearly two-sided shell. Actually Tycho has a brighter side in X-rays too. Can initial asymmetries of the SN explain ejecta features observed in some young SNRs, like the iron-rich knots at the edge of Tycho (Yamaguchi et al. 2017) and of Kepler (Sato et al. 2020, see in particular their Figure 8)? Such structures are on smaller scales than the dipole of N5ddt. From a purely hydrodynamic perspective, they would require ejecta bullets with a huge momentum to be found at the edge of the SNR (see Tutone et al. 2020 for a study of this aspect, in the context of core-collapse SNRs). In subsequent papers, we will have a closer look at the distribution of elements in space and in time. Polarization studies (Bulla et al. 2016, 2020) indicate that asymmetries are not so large that they would contradict the SN observations, which means that the power to discriminate between models is limited. This makes the striking differences in the SNR phase particularly relevant.

N5def may be a model for sub-energetic SNe events, which are less common than normal SNe Ia. After incompleteness effects are taken into account Foley et al. (2013) estimate that SNe Iax account for about 30% of all SNe Ia. We may or may not have observed a young and nearby remnant belonging to this class. Zhou et al. (2020) just proposed, from a study of the elemental abundances, that SNR Sgr A East = G0.0+0.0 is the remnant of an SN Iax, even though a core-collapse origin had been discussed in the past (Maeda et al. 2002; Park et al. 2005). This appears to be the first such claim. However, the age of this SNR remains uncertain, and it is located in the Galactic center, which is a complicated region. In any case it is interesting that there are specific predictions for the SNR phase for def models. One is the central overdensity, although this part of the ejecta is not readily accessible to observations before it is shocked by the RS (which will take thousands of years). There actually is a class of SNRs with centrally peaked thermal emission in X-rays, the so-called mixed-morphology SNRs (Rho & Petre 1998), also called thermal composites to distinguish them from plerionic composites where the central enhancement is due to a pulsar wind nebula. But these tend to be older SNRs, often interacting with molecular clouds, and other physical mechanisms are invoked such as thermal conduction and evaporation. Besides, at present no SNR in the sample has been unambiguously confirmed to be of Type Ia origin, although there are some candidates, for instance SNR 3C 397 = G41.1-0.3 (Yamaguchi et al 2015; Martínez-Rodríguez et al. 2020). We speculate that the SNR HB 9 = G160.9+02.6 may be one such remnant because of its high Fe/Si ratio (Saito et al. 2020).

Another notable feature of the def models is the presence of a network of seam lines on the ejecta surface, which are overdense and so would be overbright in X-rays. This appears to be an imprint of the deflagration dynamics, and so could possibly serve as an observational signature of this class of explosion models. In the ddt models, such structures are erased by the subsequent detonations. It is currently not clear whether such structures exist in real SNRs.

To make more quantitative comparisons between our simulations and observations, we need to develop the appropriate analysis techniques, which is work in progress. One approach is to use techniques from harmonic analysis, similar to the spherical harmonics expansion used in Section 3.2, but in 2D slab geometry (for the SNR interior) or 1D circular geometry (for the outline of the SNR) as observed in projection. A different approach is to use techniques from topology: Sato et al. (2019) applied the method of "genus statistic" to characterize the clumpy structures in Tycho's SNR. The detailed comparison of observations with our models using such techniques may allow us to probe the explosion mechanisms for SNe Ia.