Production and Distribution of ⁴⁴Ti and ⁵⁶Ni in a Three-dimensional Supernova Model Resembling Cassiopeia A

Our simulation of an SN IIb-like case, starting at ∼1400 s after core bounce, is performed with the finite-volume Eulerian multifluid hydrodynamics code Prometheus (Fryxell et al. 1991; Müller et al. 1991a, 1991b). The multidimensional Euler equations are integrated using the dimensional splitting technique of Strang (1968). The code utilizes the piecewise parabolic method (Colella & Woodward 1984) and employs a Riemann solver for real gases (Colella & Glaz 1985). To prevent numerical artifacts created by the odd–even decoupling (Quirk 1994), the AUSM+ Riemann solver of Liou (1996) is applied inside grid cells with strong grid-aligned shocks.

For spatial discretization we utilize the Yin-Yang overlapping grid technique in spherical geometry (Kageyama & Sato 2004) implemented in the Prometheus code to alleviate the restrictive Courant–Friedrich–Lewy (CFL) time-step condition in polar regions (Wongwathanarat et al. 2010a). In our simulation an angular resolution of 2° is used. The radial grid is logarithmically spaced, spanning from $\sim 3.1\times {10}^{5}$ km to 10⁹ km with a relative grid resolution ${\rm{\Delta }}r/r$ of ∼0.5% and ∼1% inside and outside of the progenitor star, respectively. The total number of (static) Eulerian grid zones is $1198\times 47\times 137\times 2$ initially. The innermost radial grid zones are successively discarded as the SN shock propagates to larger radii to relax the CFL condition. The radius of the inner grid boundary, which is treated with an outflow condition, is thus time dependent and is placed at approximately 2% of the minimum SN shock radius. At the end of the simulation the number of radial grid zones has decreased to 466 cells.

Our SN simulation is based on a 15 M_⊙ RSG progenitor, model s15s7b2 of Woosley & Weaver (1995), which was denoted as pre-collapse model W15 in Wongwathanarat et al. (2015). The first seconds of the neutrino-driven explosion of this model were simulated in 3D by Wongwathanarat et al. (2010b, 2013), using the Prometheus-Hotb code with the Yin-Yang grid. The code contains an approximate, gray treatment of neutrino transport (Scheck et al. 2006) connected with a boundary condition for the luminosities of all neutrino species imposed at a contracting Lagrangian radius placed at 1.1 M_⊙ considerably inside the neutrinosphere. By prescribing (time-dependent) neutrino luminosities at this inner grid boundary, it is possible to regulate the strength of the neutrino heating and thus to tune the energy of the neutrino-driven explosion to a desired value.

The particular model considered here is W15-2, which explodes with an energy of about $1.5\times {10}^{51}\,\mathrm{erg}=1.5\,\mathrm{bethe}=1.5$ B (1.13 B after 1.4 s and 1.47 B finally; Wongwathanarat et al. 2010b, 2013, 2015). The evolution until shock breakout from the RSG was tracked as model W15-2-cw by Wongwathanarat et al. (2015).

We stress that the asymmetries of the 3D explosion model developed stochastically mainly by convective overturn in the neutrino-heating layer. Convection, which is the dominant hydrodynamic instability in the postshock region under the conditions of this model, is seeded by initial random perturbations with an amplitude of 0.1% of the radial velocity, imposed as a cell-by-cell random pattern at the beginning of the simulation on a spherically symmetric (1D) post-collapse model about 10 ms after bounce. The large-scale asymmetries that characterize the final distribution of iron-group elements, silicon, and oxygen (see Wongwathanarat et al. 2015) neither were imposed by hand nor grew from preexisting large-scale perturbations, e.g., in the oxygen-burning shell of the progenitor star (for such possible asymmetries, which are not considered here, see Arnett & Meakin 2011; Couch et al. 2015; Müller et al. 2016b).

Since Cas A was diagnosed to be an SN IIb (Krause et al. 2008; Rest et al. 2011), we repeat the 3D long-time simulation in the present work with the W15 progenitor model after having removed the hydrogen envelope down to a rest of ∼0.3 M_⊙; see model W15-IIb in Figure 1, which shows the density profile versus enclosed mass and radius of the modified progenitor model in comparison to the original W15 model. The stellar radius of W15-IIb becomes ${R}_{* }\approx 1.5\times {10}^{7}$ km. Outside of this radius we adopt the treatment by Wongwathanarat et al. (2015) and assume the presence of circumstellar matter with density and temperature profiles following an ${r}^{-2}$ decline. Since the onset and early phases of the explosion are not affected by our modification of the hydrogen envelope, we map initial data from model W15-2-cw at a post-bounce time of ${t}_{\mathrm{pb}}=1431\,{\rm{s}}$ when the SN shock has nearly reached R_*. We then follow the subsequent evolution of the SN ejecta until approximately 1.25 days after core bounce. This new explosion simulation is termed W15-2-cw-IIb.

Zoom In Zoom Out Reset image size

By artificially removing most of the hydrogen, we destroy the self-consistent, hydrostatic structure of the near-surface layers, and the remaining hydrogen envelope becomes unrealistic in its radial profile. This is likely to corrupt the shock propagation through the hydrogen shell and the shock breakout from the stellar surface. However, this manipulation has no effect on the large-scale asymmetries imposed by the explosion mechanism on the innermost SN ejecta such as the ⁵⁶Ni and ⁴⁴Ti yields. Moreover, our modified model is also good enough to permit the demonstration of consequences associated with the absence of the massive hydrogen envelope, in particular for the maximum expansion velocities that can be retained by the mentioned nucleosynthesis products. As shown by Wongwathanarat et al. (2015), the detailed composition-shell structure of the SN progenitor, especially the density profile in the He and H layers, has a strong influence on the propagation of the SN shock, on the development of secondary mixing instabilities and reverse shocks at the composition-shell interfaces, and on the interaction of these instabilities and reverse shocks with the initial ejecta asymmetries imposed by the central explosion mechanism itself. In the case of the lack of a massive hydrogen envelope, for example, the reverse shock forming at the He/H interface after the passage of the SN shock is absent and the corresponding deceleration of fast-moving inner ejecta, especially also of nickel and titanium, does not take place. We will discuss these effects in Section 3.2.

Starting with a neutrino-driven explosion with highly aspherical mass ejection from the beginning, we are unable to quantitatively assess a possible reduction of the strength of Rayleigh–Taylor (RT) instabilities in material of intermediate-mass numbers Z, due to the lack of a strong reverse shock from the He/H interface. Such a reduction is suggested by the works of Ellinger et al. (2012, 2013), where a strong reverse shock had a substantial influence on the RT growth at the C–O/He composition interface. In contrast to our explosion modeling, Ellinger et al. (2012, 2013) launched their explosions in 1D and followed the propagation of a spherical blast wave from the C–O core into the He layer. However, the initial asymmetries created by the neutrino-driven mechanism in collaboration with hydrodynamic instabilities during the very first second of the explosion affect the iron and titanium nucleosynthesis and imprint large-scale asphericities on the innermost ejecta, the outgoing SN shock, and the SN explosion as a whole. They are therefore crucial in seeding the secondary mixing instabilities that grow at the C–O/He and He/H shell interfaces of the progenitor, and the presence of pronounced initial explosion asymmetries causes the radial mixing to be much more efficient and to penetrate deeper (well down into the region of Fe formation) than in the case of explosions initiated in spherical symmetry (see, e.g., Kifonidis et al. 2003, 2006; Wongwathanarat et al. 2015). Since the RT growth factor at the C–O/He-core interface in our model is huge after the passage of the SN shock (see Figure 5 in Wongwathanarat et al. 2015), the development of RT instability at this interface depends primarily on the initial perturbations in the metal-rich ejecta (including the layers containing ⁵⁶Ni and ⁴⁴Ti) and not on the lack or presence of the reverse shock from the He/H interface. As a consequence, we observe the same highly deformed ejecta geometry on the largest scales for model W15-2-cw-IIb and for the RSG explosion W15-2-cw, before in the latter case the reverse shock from the He/H interface decelerates the core of intermediate-mass elements and RT instability at the base of the hydrogen shell triggers the fragmentation of the inner ejecta to smaller-scale structures (compare Figure 7, top row, of Wongwathanarat et al. (2015) with the results reported in Section 3.3 below).

In our simulation of model W15-2-cw-IIb, gravity is taken into account in the hydrodynamics equations. We use the spherical harmonics expansion technique of Müller & Steinmetz (1995) to solve the Poisson equation with source terms for the Newtonian self-gravity of the SN ejecta and the gravitational contribution of a point mass that is placed at the coordinate origin and accounts for the central compact remnant and additional material that starts as ejecta gas but gets absorbed through the inner grid boundary in the course of the evolution by fallback and the successive outward relocation of the inner grid boundary.

The tabulated equation of state of Timmes & Swesty (2000) is employed for describing the stellar plasma as a mix of arbitrarily degenerate and relativistic electrons and positrons, blackbody radiation, and ideal Boltzmann gases of a defined set of fully ionized nuclei, taking into account corrections for Coulomb effects. We consider 11 nuclear species, namely, protons, nine α-nuclei from ⁴He to ⁵⁶Ni (excluding ³²S, ³⁶Ar, ⁴⁸Cr, and ⁵²Fe from the α chain), and a "tracer nucleus" X that (schematically) represents iron-group species synthesized in neutrino-heated ejecta under conditions of neutron excess. The advection of nuclear species with plasma motions is treated with the consistent multifluid advection scheme of Plewa & Müller (1999).

In the discussed simulation of model W15-2-cw-IIb nuclear burning plays no role, because nuclear reactions become unimportant in the expanding SN ejecta during the late stages and on the timescales considered. The simulation is therefore devoted to studying the long-time evolution of the spatial and velocity distribution of the nucleosynthesis products obtained in the earlier phases, i.e., in models W15-2 and W15-2-cw, from which we adopt the initial data for the 3D evolution of model W15-2-cw-IIb. Nuclear fusion processes were taken into account in the computations of models W15-2 and W15-2-cw through a small α-reactions network, connecting the 11 nuclear species mentioned above, details of which were given by Wongwathanarat et al. (2013, 2015). The small network is able to provide rough information about the assembling of heavy nuclei in the SN outflows, but it is unable to give accurate results for the yields of individual nuclear species. The production of ⁴⁴Ti, for example, is massively overestimated by the small network (see Table 3 of Wongwathanarat et al. 2013). For a more accurate assessment of the nucleosythesis of ⁵⁶Ni and ⁴⁴Ti, we therefore perform a post-processing analysis with a large network, employing tracer particle information that we extract from our 3D simulations of models W15-2, W15-2-cw, and W15-2-cw-IIb as described in the next section.

The post-processing for nucleosynthetic yields is performed with a nuclear reaction network of 6300 nuclear species between the proton-drip line and the neutron-drip line, reaching up to the Z = 110 isotopes, as described in Wanajo (2006). The latest reaction library of REACLIB V2.0 (Cyburt et al. 2010) is employed, in which experimentally evaluated values are adopted whenever available. We also test a possible uncertainty of the most influential reaction for ⁴⁴Ti production, ⁴⁴Ti(α, p)⁴⁷V, by dividing the forward and inverse rates in REACLIB V2.0 by a factor of two, according to the new experimental evaluation in Margerin et al. (2014). Using the hydrodynamical trajectories described below, nucleosynthesis calculations are started when the temperature decreases to 10¹⁰ K, assuming initially free protons and neutrons with mass fractions Y_e and $(1-{Y}_{e})$, respectively. Nuclear statistical equilibrium is immediately established at such high temperature, and the detailed initial composition is unimportant. For the tracer particle trajectories with maximum temperatures below 10¹⁰ K, the calculations are started from the beginning of the hydrodynamic simulation (i.e., shortly after core bounce) with the composition of the pre-SN model.

A total number of 131,072 tracer particles is initially (i.e., some 15 ms after core bounce) placed between ∼1000 and ∼6900 km, covering the entire Si layer and the inner shells of the O layer, corresponding to a mass resolution of $3.49\times {10}^{-6}\,{M}_{\odot }$. This setting is suitable to probe the whole volume where ⁵⁶Ni (iron-group material) and ⁴⁴Ti are nucleosynthesized during the explosion (the corresponding peak temperatures and peak densities are shown in Figure 2). It is important to note that we do not co-evolve the tracer particles with the hydrodynamic simulation but use the data files stored from the 3D hydrodynamics run to construct the tracer histories only after the simulation. This is achieved by a fourth-order Runge–Kutta time integrator for the particle trajectories. Linear interpolation of the velocity field is used to obtain velocities at the spatial positions of the particles. The data of our simulation are stored every 100 hydrodynamics time steps. These time intervals between hydrodynamics outputs are further divided into 10 substeps, and the velocity field at the intermediate steps is obtained by linear interpolation in time.

Zoom In Zoom Out Reset image size

Despite the higher-order time integration, the resolution limitations associated with the discrete time sampling by the data outputs (which is much coarser than the small time steps applied in the simulation) can lead to integration errors of the particle evolution. This problem should be kept in mind and will be a matter of discussion in Section 3.2.

The majority of our initial particles, namely, 107,910, are ejected, while 23,162 particles fall through the inner grid boundary and are removed from further analysis. We discriminate between two types of nucleosynthesis-relevant ejecta and the corresponding particles. On the one hand, there are shock-heated ejecta, which are immediately accelerated outward by the expanding SN shock. On the other hand, there are neutrino-processed ejecta, which fall inward below ∼250 km to be affected by electron capture and/or neutrino interactions, which can change the electron fraction, Y_e, and the entropy of the plasma. While the shock-heated matter retains the Y_e from the pre-collapse conditions, the electron fraction of the neutrino-processed material is reset by the neutrino reactions computed with our gray transport approximation. This implies considerable uncertainties not only because our neutrino transport treatment is approximative but also, more generally, because the final neutron-to-proton ratio is very sensitive to many details of the neutrino physics, some of which (such as collective neutrino-flavor oscillations and the lepton-emission self-sustained asymmetry recently found by Tamborra et al. 2014; for a review, see e.g., Mirizzi et al. 2016) are not yet fully understood.

A subset of 75,862 of the expelled particles can be considered as shock-heated ejecta, and another 27,222 particles represent neutrino-processed ejecta. This corresponds to ∼0.265 M_⊙ and $\sim 9.51\times {10}^{-2}\,{M}_{\odot }$ for the shock-heated and neutrino-processed ejecta, respectively. (A small number of 4826 remaining particles did not produce ⁴⁴Ti and ⁵⁶Ni at any significant amounts and therefore are irrelevant for the investigations of this paper.) Because, for the reasons mentioned above, Y_e is uncertain in the neutrino-processed ejecta, we will consider two cases for the nucleosynthetic post-processing of the corresponding particles in order to test the sensitivity of the nucleosynthesis to Y_e: (1) We adopt Y_e as obtained in our 3D simulations with the approximative neutrino treatment at a temperature of $T\simeq 5\times {10}^{9}$ K in outflowing material. The 3D long-time simulation with the corresponding post-processing results for nucleosynthesis will be denoted as W15-2-cw-IIb-Ye_sim. (2) Alternatively, we assume that all of the SN ejecta, including the neutrino-heated ones, inherit their Y_e values directly from the progenitor, i.e., we take the Y_e data from the progenitor profile, where ${Y}_{e}\geqslant 0.4979$ for pre-collapse radii of $r\gtrsim 1160$ km and enclosed masses larger than 1.287 M_⊙, while the more neutron-rich, deeper layers are accreted into the newly formed NS. (The final mass cut is close to 1.35 M_⊙.)

Typical entropies in the neutrino-heated ejecta at T = 5 GK range from ∼10k_B per nucleon to ∼30k_B per nucleon, while the shock-heated material has entropies between ∼7k_B and ∼12k_B per nucleon during the relevant time of ⁴⁴Ti and ⁵⁶Ni nucleosynthesis.

As we will see in Section 3.1, the nucleosynthetic yields of ⁵⁶Ni and, in particular, of ⁴⁴Ti differ considerably between the two investigated cases. Because of a better match with the Cas A remnant, we will mainly focus on the second case as our preferred and more interesting model, which we term W15-2-cw-IIb.

We have tested the convergence of our post-processing analysis by comparing the yields for our full set of tracer particles with yields estimated on grounds of randomly selected subsets of these particles. We found that typically less than 10% of all particles are sufficient to obtain the produced ⁴⁴Ti and ⁵⁶Ni masses with accuracies of better than a few percent. Of course, good representations of the ⁴⁴Ti and ⁵⁶Ni distributions in 3D space require the use of significantly more than a few thousand tracer particles.

Figure 2 shows the peak temperature and peak density versus "initial" (post-bounce) radius, R_ini, of all particles that contribute to the production of ⁴⁴Ti and ⁵⁶Ni. Shock-heated and neutrino-processed ejecta can be well distinguished. The shock-heated particles reach ${T}_{\mathrm{peak}}\lesssim 9$ GK and ${\rho }_{\mathrm{peak}}\lesssim {10}^{7.6}$ g cm⁻³ and align in rather narrow stripes of temperature and density values as functions of initial radius, because the shock front is only moderately aspherical. The neutrino-processed ejecta form clusters of higher T_peak and ${\rho }_{\mathrm{peak}}$ at post-bounce radii of less than ∼2200 km. While most of these particles are neutrino-heated to temperatures up to about 40 GK in the convective flows within the postshock layer (peak densities ${\rho }_{\mathrm{peak}}\lesssim 2\times {10}^{10}$ g cm⁻³), some of this matter (starting from smaller ${R}_{\mathrm{ini}}\lesssim 1500$ km and reaching to even more extreme peak values of density and temperature) is accreted onto the nascent NS before being blown out again in the neutrino-driven wind.

The electron fraction in the neutrino-processed ejecta, which is set by the competition of ${\nu }_{e}$ and ${\bar{\nu }}_{e}$ absorption and emission processes, sensitively determines the efficiency of producing ⁴⁴Ti and ⁵⁶Ni. For this reason the yields of both nuclei differ strongly between the two considered cases, namely, model W15-2-cw-IIb, where, for a sensitivity test, we adopt the electron fraction, Y_e, from the progenitor conditions, compared to model W15-2-cw-IIb-Ye_sim, where we take Y_e from our 3D simulations with approximative neutrino transport. In the latter case, Y_e is shifted slightly to the neutron-rich side (values at $T\simeq 5\times {10}^{9}$ K lie in the range between 0.456 and 0.498), which strongly suppresses the creation of both ⁴⁴Ti and ⁵⁶Ni as discussed by Magkotsios et al. (2010), where the theoretical framework of the production of these radioactive isotopes was depicted in great detail. As the ⁴⁴Ti yield becomes maximal at ${Y}_{e}\approx 0.5$ with other conditions fixed (Wanajo et al. 2013), our two investigated cases roughly bracket the upper and lower limits of the ⁴⁴Ti production in our model.

Table 1 lists the total yields of models W15-2-cw-IIb and W15-2-cw-IIb-Ye_sim, as well as the individual contributions from shock-heated and neutrino-processed material in both cases, all obtained by the nucleosynthetic post-processing of our tracer particles. In model W15-2-cw-IIb-Ye_sim the shock-heated ejecta component (indicated by the model name extension "shock") dominates the expelled mass of ⁵⁶Ni by far, whereas both ejecta components contribute roughly equally to the making of ⁴⁴Ti. In model W15-2-cw-IIb roughly 95% of the ⁴⁴Ti and about 56% of the ⁵⁶Ni are assembled in the neutrino-processed ejecta (model name extension "νproc"), where Y_e is essentially 0.5. Effectively, this means that model W15-2-cw-IIb produces more than twice as much ⁵⁶Ni as W15-2-cw-IIb-Ye_sim (∼0.096 M_⊙ compared to ∼0.043 M_⊙) and 10 times more ⁴⁴Ti ($\sim 1.6\times {10}^{-4}$ M_⊙ compared to $\sim 1.6\times {10}^{-5}$ M_⊙).

Table 1.  Yields from Nucleosynthetic Post-processing of Tracer Particles with Y_e = 0.5 (Upper Three Lines) and Y_e from the Simulation with Approximate Neutrino Transport (Lower Three Lines), Using Our Standard REACLIB V2.0 Nuclear Rates

Model	M(44Ti) (M⊙)	M(56Ni) (M⊙)
W15-2-cw-IIb	$1.57\times {10}^{-4}$	$9.57\times {10}^{-2}$
W15-2-cw-IIb-shock	$8.66\times {10}^{-6}$	$4.20\times {10}^{-2}$
W15-2-cw-IIb-νproc	$1.49\times {10}^{-4}$	$5.38\times {10}^{-2}$
W15-2-cw-IIb-Yesim	$1.58\times {10}^{-5}$	$4.29\times {10}^{-2}$
W15-2-cw-IIb-shock	$8.66\times {10}^{-6}$	$4.20\times {10}^{-2}$
W15-2-cw-IIb-νproc-Yesim	$7.16\times {10}^{-6}$	$0.10\times {10}^{-2}$

Download table as:  ASCIITypeset image

We note in passing that the calculations with the small α-network that we employ during the hydrodynamic simulations massively overestimate the yield of ⁴⁴Ti, whereas they produce ⁵⁶Ni in the same ballpark as the tracer-particle-based post-processing with the large nuclear network reported above. Table 3 of Wongwathanarat et al. (2013) lists nearly $3\times {10}^{-3}\,{M}_{\odot }$ of ⁴⁴Ti and 0.055–0.139 M_⊙ of ⁵⁶Ni, depending on the undetermined contribution of this nucleus in the neutrino-processed ejecta, for model W15-2.

By far most of the ⁴⁴Ti and ⁵⁶Ni are assembled under α-rich freeze-out conditions, which are experienced by most of the relevant shock-heated and neutrino-heated ejecta, and where the mass fractions of these nuclei reach their highest values. Only subdominant contributions originate from regions of incomplete Si burning (in line with conclusions drawn by Woosley & Hoffman 1991; Magkotsios et al. 2010). This can be concluded from the plots of Figure 3, which show color-coded mass fractions of ⁴⁴Ti (left) and ⁵⁶Ni (right) for all tracer particles of model W15-2-cw-IIb in the peak temperature–peak density plane. The two images should be viewed in relation to Figure 4 of Magkotsios et al. (2010). All of the expelled matter with ${T}_{\mathrm{peak}}\gtrsim 5\times {10}^{9}$ K crosses the α-particle-rich freeze-out region when it expands and cools during ejection, and only a small amount of ⁴⁴Ti and ⁵⁶Ni is produced by incomplete Si burning in matter with ${T}_{\mathrm{peak}}\lt 5\times {10}^{9}$ K. The banana-shaped cluster of particles between peak densities of about 10⁸ and $2\times {10}^{10}$ g cm⁻³ corresponds to neutrino-heated matter in the convective layer between gain radius and SN shock, while the antenna-like narrow extension toward higher densities and temperatures is matter that has transiently been accreted into the outer layers of the proto-NS before it was blown out in the neutrino-driven wind later. The mass ratio of ⁴⁴Ca/⁵⁶Fe as stable daughter nuclei of ⁴⁴Ti and ⁵⁶Ni in our model W15-2-cw-IIb is $1.64\times {10}^{-3}$ and thus compatible with the solar ratio of these elements ($=1.26\times {10}^{-3};$ Lodders 2003), while our less favorable case, model W15-2-cw-IIb-Ye_sim, accounts for a production of only about 1/3 of the solar ⁴⁴Ca/⁵⁶Fe value.

Zoom In Zoom Out Reset image size

All tracer particles that contribute to the creation of significant amounts of ⁴⁴Ti and ⁵⁶Ni are ejected (i.e., accelerated outward) during the first second of the explosion. In our 3D neutrino-driven explosion model the far dominant part of the ejected ⁴⁴Ti originates from the neutrino-heated, high-entropy postshock layer, whose mass depends on the multidimensional processes that drive nonradial flows in this region. More plasma is channeled through the region of neutrino-matter interactions by such flows than typically experiencing neutrino heating in corresponding 1D explosion models. Moreover, the nucleosynthesis of ⁴⁴Ti is extremely sensitive to the neutron-to-proton ratio and entropy in the ν-processed matter. Spherically symmetric models of neutrino-powered explosions are certainly not able to capture the dynamics of this layer reliably. For example, in our attempts to perform quantitative comparisons of 1D and 3D simulations by varying the parameters that determine the strength of the explosion-driving neutrino source, we found that it is hardly possible to reproduce the shock acceleration and the growth of the explosion energy of our 3D model simultaneously by 1D explosions. The yield of ⁴⁴Ti hinges critically on all characteristic properties (i.e., mass, Y_e, entropy, and expansion velocity) of the ejecta that emerge from the region where the mass cut between compact remnant and SN outflow develops. For this reason 1D models are not qualified to adequately address the problem of ⁴⁴Ti production in SNe. This, in particular, also applies to SN calculations invoking pistons or thermal bombs for initiating the explosions. Both methods are unable to accurately describe the dynamics, thermodynamics, and electron fraction of matter expelled from the vicinity of the mass cut.

It is therefore not astonishing that 3D neutrino-driven explosions can potentially produce much more ⁴⁴Ti than predicted by piston-initiated explosions (as already surmised by Woosley & Hoffman 1991), where this isotope is only made in the shock-heated ejecta, which contribute to our ⁴⁴Ti yield only on the 10% level. Our results also demonstrate that large global asymmetries, e.g., by rapid rotation or even jet-driven (magnetohydrodynamic) explosions (Nagataki et al. 1998), are not needed to explain the ⁴⁴Ti masses observed in Cas A and SN 1987A. The yield of our model W15-3D-cw-IIb is compatible with both observations (see Figure 4),⁷ in particular if we account for uncertainties in the important nuclear reaction rate of the ⁴⁴Ti ${(\alpha ,p)}^{47}{\rm{V}}$ process. Reducing this rate by a factor of two according to the new experimental value of Margerin et al. (2014) boosts the ⁴⁴Ti production by a factor of 1.5 compared to the values listed in Table 1. We did not attempt to regulate the explosion energy of our 3D model for a precise match of the observations. Its value of nearly 1.5 B is close to the one estimated for SN 1987A (around 1.5 B; Utrobin 2005) and only slightly smaller than that of Cas A (roughly 2.3 B; Laming & Hwang 2003; Orlando et al. 2016, and references therein). Despite this energy difference, both SNe seem to have produced ⁴⁴Ti yields in the same ballpark (Figure 4), although the SN 1987A value is not determined very accurately and there is considerable tension between the results obtained from X-ray emission-line measurements (Grebenev et al. 2012; Boggs et al. 2015), late-time spectral analysis (Jerkstrand et al. 2011), and late-time light-curve fitting (Seitenzahl et al. 2014).

Zoom In Zoom Out Reset image size

Neutrino-driven explosions generically eject large amounts of matter with sufficiently high entropies (between ∼10k_B and ∼30k_B per nucleon) to enable an effective production of ⁴⁴Ti, and therefore they can potentially produce considerable yields of this nucleus. The reason is that the explosion energy is carried by matter that must be heated to the point of marginal gravitational binding in the strong field of the newborn NS before outward acceleration of the SN shock and of the postshock plasma can set in (for a detailed discussion, see, e.g., Janka 2001; Scheck et al. 2006). The net energy available to the explosion is then basically provided by the recombination energy of initially free nucleons in the neutrino-heated gas to α-particles and heavy nuclei. This releases at most ∼9 MeV per nucleon, of which under typical conditions 5–7 MeV per nucleon or $\sim (5\mbox{--}7)\times {10}^{18}$ erg g⁻¹ are available to power the SN blast. For the explosion energy one therefore estimates, very roughly,

$$\begin{eqnarray}&&{E}_{\exp }\sim (1.0\,\ldots \,1.4)\times {10}^{51}\left(\displaystyle \frac{{\rm{\Delta }}{M}_{\nu \mbox{-} \mathrm{heat}}}{0.1\,{M}_{\odot }}\right)\mathrm{erg},\end{eqnarray} \tag{ 1 }$$

where ${\rm{\Delta }}{M}_{\nu \mbox{-} \mathrm{heat}}$ is the mass of neutrino-heated material, and the negative gravitational binding energy of the overlying stellar layers exterior to the SN shock and positive energy contributions by explosive nuclear burning of shock-heated ejecta are not taken into account. This expression means that for a net explosion energy of 1 ${\rm{B}}={10}^{51}$ erg, at least ∼0.1 M_⊙ of neutrino-heated matter must be expelled.

Our results for models W15-2-cw-IIb and W15-2-2cw-IIb-Ye_sim, however, demonstrate that the actual mass of ⁴⁴Ti produced in the explosion can vary strongly, depending on the exact value of Y_e in the neutrino-processed ejecta. A reliable determination of the neutron excess requires a more sophisticated treatment of the neutrino transport than applied in our 3D simulations and, moreover, is sensitive to aspects that are uncertain or hard to control with high accuracy. Therefore, since predictions of Y_e with a precision of a few percent are extremely difficult, we do not want to put too much weight on the elemental yields, instead focusing on the distributions of ⁴⁴Ti and ⁵⁶Ni as functions of enclosed mass, radial velocity, and spatial coordinates in the following.

Hydrodynamic instabilities during the neutrino-heating phase determine the overall asphericity of the explosion and the directional asymmetries of the nucleosynthetic ejecta. The initial structures begin to freeze out at the time when the outward acceleration of the SN shock marks the onset of the explosion. The final morphology of the SN blast and the velocity distributions of the chemical elements, however, depend strongly not only on this initial explosion asymmetry but also on the structure of the progenitor star. The latter has an important influence because the outgoing shock accelerates or decelerates when traveling through the different composition shells, depending on whether the density profile is steeper or flatter than ${r}^{-3}$. Deceleration after the shock has passed a composition-shell interface leads to compression of the postshock matter into a dense shell and the formation of a reverse shock that moves backward into the SN ejecta and slows down the expanding flow. Moreover, in regions of crossing pressure and density gradients the growth of RT and Richtmyer–Meshkov instabilities is seeded by the shock asphericity, preexisting perturbations in the progenitor, and, most importantly, the primary explosion asymmetries in the innermost SN ejecta. The secondary instabilities thus enhance the mixing between layers of different chemical composition and cause a fragmentation of the initially large, metal-rich plumes and clumps to smaller, more filamentary structures as discussed in detail for 3D SN simulations of different red and blue supergiant progenitors, including W15-2, by Wongwathanarat et al. (2015).

The size of such effects, in particular the efficiency of the metal penetration into the hydrogen envelope and the corresponding fragmentation of the mixed layer, depends on the shock velocity after the shock has crossed the C–O core, the shock deceleration in the He layer, and the density gradient at the He/H interface. In RSGs the shock accelerates strongly at a steep density decline between the He and H shells and experiences a similarly dramatic subsequent deceleration when traveling through the extended hydrogen envelope. This leads to large growth rates of the RT instability at the shell interface and therefore efficient metal and helium mixing into the hydrogen envelope. In blue supergiants the density gradient at the He/H boundary is less pronounced, the shock acceleration at this location is therefore far more moderate, and the shock deceleration in the more compact H envelope is less extreme. Because the RT growth rates at the He/H transition are smaller, the intrusion of metals into the hydrogen envelope is impeded, unless the fastest metal plumes stay in close touch with the outgoing shock and thus can enter the hydrogen envelope before they experience significant deceleration by interaction with the reverse shock and with a dense "helium wall" that builds up from decelerated helium. If such a situation holds, extended metal-enriched fingers can reach with high velocities deep into the surrounding hydrogen also in the case of blue supergiants (see Wongwathanarat et al. 2015).

These results demonstrate that the presence of the hydrogen envelope in SNe IIP has important consequences for the final velocities and the mixing and fragmentation of the distribution of metal-containing ejecta. In Figure 5 we present the mass distributions as functions of the radial velocity for a set of nucleosynthesis components for model W15-2-cw at the time of shock breakout (bottom panel; data from Figure 12 of Wongwathanarat et al. 2015) compared to results of our present simulation of model W15-2-cw-IIb, i.e., after having removed the hydrogen envelope except for a relic mass of ∼0.3 M_⊙.

Zoom In Zoom Out Reset image size

In model W15-2-cw the mass distributions of ⁵⁶Ni and ⁴⁴Ti reach up to about 5000 km s⁻¹ for a small faction of this radioactive material, which indicates the extent of mixing in velocity space. The overall shapes of the distributions of all elements (iron-group, intermediate-mass, and light elements) resemble each other, with the bulk ejecta velocities ranging from ∼1000 to ∼3000 km s⁻¹, and the maxima and the positive slopes below the maxima exhibiting the order expected from homologous expansion ($v\propto r$) of spherically symmetric explosions, where lighter elements reside at larger radii with higher velocities. The tails of the fastest ejecta above the distribution maxima stretch out to roughly 7000 km s⁻¹ and contain, besides hydrogen and helium, also intermediate and light elements (silicon to carbon) corresponding to the chemical composition of the progenitor envelope.

For a given explosion energy E_exp, the average velocity of the ejecta (mass M_ej) is approximately given by ${\bar{v}}_{\mathrm{ej}}\sim \sqrt{2{E}_{\exp }/{M}_{\mathrm{ej}}}$, and also the maximum ejecta velocity (ignoring a small amount of mass accelerated off the stellar surface at shock breakout) follows roughly the scaling with the ratio $\sqrt{{E}_{\exp }/{M}_{\mathrm{ej}}}$. Our RSG model W15-2-cw (${E}_{\exp }\approx 1.5\,{\rm{B}}$, ${M}_{\mathrm{ej}}\approx 13.6\,{M}_{\odot }$) yields ${E}_{\exp }/{M}_{\mathrm{ej}}\approx 0.11\,{\rm{B}}/{M}_{\odot }$, corresponding to ${\bar{v}}_{\mathrm{ej}}\sim 3300$ km s⁻¹. Without the decelerating influence of the massive hydrogen envelope, model W15-2-cw-IIb possesses a higher value of ${E}_{\exp }/{M}_{\mathrm{ej}}\approx 1.5\,{\rm{B}}/3.3\,{M}_{\odot }\approx 0.45$ B/M_⊙, for which reason one expects about twice the ejecta velocities of the RSG model.

Indeed, the main ejecta mass (helium for W15-2-cw-IIb) exhibits velocities around 6000–8000 km s⁻¹, and small amounts of material expand with more than 10,000 km s⁻¹ (Figure 5, top panel). The mass distributions of different elements are much more dissimilar in model W15-2-cw-IIb than in model W15-2-cw. This suggests considerably less mixing, fully in line with the absence of a reverse shock and RT-unstable conditions at the He/H interface. Hydrogen and helium naturally show the tendency of having the highest velocities for the bulk of their mass, followed by carbon and silicon, whereas the dominant masses of oxygen, titanium, and nickel are distributed in a broad maximum below ∼7000 km s⁻¹. Only minor fractions of the radioactive elements (of order ≲1% of the total ⁴⁴Ti and ⁵⁶Ni mass) are contained in high-velocity tails that extend from 7000 km s⁻¹ up to ∼9000 km s⁻¹.

The highest velocities of significant amounts of ⁴⁴Ti (around 7000 km s⁻¹) in our simulations are in good agreement with the upper bound on the fastest such material determined by Cas A observations (5350 ± 1610 km s⁻¹, Grefenstette et al. 2014; ∼6300 ± 1250 km s⁻¹, Grefenstette et al. 2017).⁸ Besides the question how much dilute, high-velocity ⁴⁴Ti might be present above this velocity range, but hard to detect observationally, our model results need to be considered with some caution, too. First, we remind the reader once again of the fact that we did not attempt to fine-tune our simulations for optimal comparison with the observations. Detailed analysis of Cas A suggests an explosion energy of ∼2.3 B for an ejecta mass of ∼4 M_⊙ (Vink 2004; Orlando et al. 2016), which implies a value of ${E}_{\exp }/{M}_{\mathrm{ej}}\approx 0.575$ B/M_⊙. Therefore, an optimal Cas A model might possess bulk and maximum velocities of ⁴⁴Ti that are even 10%–15% higher than in our model W15-2-cw-IIb. However, our simulations currently track the evolution of the SN explosion only for the first day. Although the ejecta are already close to homology at this time, our models ignore important future effects such as the heating associated with radioactive decays and the deceleration of the fastest ⁴⁴Ti by the reverse shock that originates from the interaction of the ejecta with the circumstellar environment decades later. The long-time evolution from the first day to the remnant stage as observed at the present epoch will be the subject of future studies, for which reason we postpone a more detailed quantitative comparison of the velocity information with Cas A observations to a later stage.

The middle panel of Figure 5 displays the mass distributions in velocity space for ⁴⁴Ti and ⁵⁶Ni as obtained by our post-processing of tracer particles with a large nucleosynthesis network (as described in Section 2.4), in contrast to the top and bottom panels of this figure, where results from hydrodynamic simulations with the small nucleosynthesis network are presented. While the distributions for model W15-2-cw-IIb in the top and middle panels look roughly similar for the bulk of the radioactive material, they nevertheless exhibit an important difference: the tracer particle distributions possess a considerably more populated high-velocity tail, which also reaches up to higher velocities. While the hydrodynamic simulation of model W15-2-cw-IIb with the small network yields ≲1% of the ejected ⁴⁴Ti and ⁵⁶Ni with velocities of more than 7000 km s⁻¹, the tracer particle analysis predicts such a small fraction of the radioactive nuclei to be faster than ∼8500 km s⁻¹, a tiny amount to be ejected with velocities up to 9700 km s⁻¹ instead of ∼9000 km s⁻¹ for the hydrodynamic run, and some 5% (8%) of the ⁵⁶Ni (⁴⁴Ti) to be expelled with over 7000 km s⁻¹.

Tests showed that this discrepancy is a consequence of numerical inaccuracies in our calculation of the tracer particle histories. As mentioned in Section 2.4, we do not evolve the tracer particles online with the hydrodynamic model but reconstruct the particle velocities and positions by higher-order time integration based on a finite number of output files from the hydrodynamic run. Since these output files are dense in the early phase (the first seconds of the explosion), we do not consider this approach as problematic during the time interval when the nucleosynthesis takes place. Indeed, inspecting the mass distributions versus radial velocity at 10 s, the results of the calculations with small network (during the hydro run) and large network (based on particle post-processing) exhibit much closer similarities with respect to their shape and spread in velocity space, apart from effects that are expected as consequences of nucleosynthesis differences between the approximative α-network treatment and the more detailed calculation (e.g., a shift of the ⁴⁴Ti production in the latter case to less strongly heated, higher-velocity material instead of a distinctly synchronous production with ⁵⁶Ni by the small α-network). Integrating the particle histories to late times, however, numerical errors accumulate, in particular since data files for the later stages are stored only with larger time intervals.

In spite of this obvious caveat, we will further analyze the nucleosynthesis-related aspects of our results on the basis of the particle information, because only the big-network calculation provides reliable quantitative information on the production of ⁴⁴Ti and ⁵⁶Ni for their mutual comparison (though the treatment of the latter is less problematic with the small network; see Section 3.1). Yet we should keep in mind that this implies an overestimation of the abundances of these nuclei in the highest-velocity tails of their late-time distributions. However, since the problem concerns only a small amount of matter (some percent of the total ejected ⁴⁴Ti and ⁵⁶Ni masses), it is not of great relevance for our discussion of the spatial asymmetries in Section 3.3.

Figure 6 shows the initial and final distributions of the bulk of ⁴⁴Ti and ⁵⁶Ni in mass space. The mass coordinate of a tracer particle is defined by the stellar mass contained in a sphere whose radius is given by the particle's radial position. All of the nickel and the far dominant part of the titanium originate from a mass range between 1.35 and 1.6 M_⊙ in the progenitor star, which encompasses neutrino-processed and shock-heated ejecta (left panel of Figure 6). A small fraction (of order 1%) of the ⁴⁴Ti is made in matter associated with an initial mass interval from 1.6 to 1.8 M_⊙, where the outgoing shock raises the peak temperatures to $(4\mbox{--}5)\times {10}^{9}$ K and allows some ⁴⁴Ti to be produced during incomplete silicon burning. The dip between the two maxima is connected to a depletion region called the "chasm" by Magkotsios et al. (2010).

Zoom In Zoom Out Reset image size

At nearly a day the distributions of ⁴⁴Ti and ⁵⁶Ni are smeared over the whole range of ejecta masses, indicating efficient radial transport of the innermost SN ejecta into the carbon-oxygen and helium layers of the exploding star (right panel of Figure 6; corresponding cumulative distributions versus radius and enclosed mass are displayed in Figure 7). While the main formation of nickel and titanium takes place in the same mass region, but with considerable variation of the relative production (left panel), the close similarity of the normalized mass distributions of both radioactive nuclei after 1 day does not imply that both species are well mixed to equal composition in all regions. As we will see from our analysis of the spatial distributions in Section 3.3, there are regions that survive with ⁴⁴Ti being enhanced relative to ⁵⁶Ni and vice versa, compared to the average production ratio of the two species.

Zoom In Zoom Out Reset image size

Finally, we again allude to the fact that the tracer particles with the highest velocities are not very reliable because of integration errors in the particle histories. For the mass distributions in the right panel of Figure 6, this concerns a subdominant mass of ⁴⁴Ti and ⁵⁶Ni accumulating in a local maximum between ∼4.2 and ∼4.7 M_⊙. This feature is absent in our results from the small α-network calculation. The peak contains a small fraction (a few percent) of the material that originates from the initial mass interval of 1.35–1.55 M_⊙. In the case of nickel, both shock-heated and neutrino-heated matter contributes, whereas in the case of titanium the peak is supplied only by tracer particles carrying yields from the far dominant neutrino-processed ejecta component.

Spatial variations of the abundance mix of ⁴⁴Ti and ⁵⁶Ni can be concluded from Figure 8, where we show column densities of both nuclei as observed from a direction perpendicular to the kick direction of the newly formed NS, which is indicated by a white arrow. The orientation of our image is chosen such that direct comparison with Figures 2 and 3 and extended data Figure 1 of Grefenstette et al. (2014) (see Figure 9) is possible and exhibits the closest resemblance that we can achieve. Although nickel and titanium can be found in coexistence in the whole mass-filled volume, the color variations reflect considerable local differences of the abundance ratios of both nuclei, with more intense blue colors signaling relatively higher concentrations of ⁴⁴Ti, and more intense green picturing higher concentrations of ⁵⁶Ni.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In the right panel of Figure 8, we show the nickel ejecta only at large radii but do not display the 50% of the expelled ⁵⁶Ni mass that is contained in an inner, spherical volume. This is motivated by investigations of the Cas A SNR by Milisavljevic & Fesen (2015) and Orlando et al. (2016), who concluded that about half of the iron produced by the explosion could reside in the central region that has not yet been heated to X-ray emission temperatures by the reverse shock from the ejecta-environment interaction. After ⁵⁶Ni has decayed to stable ⁵⁶Fe, ionization by the reverse shock is necessary to enable Fe X-ray emission from atomic transitions. In contrast, ⁴⁴Ti continues to be radioactive well into the SNR stage, even at the age of Cas A, and therefore radiates X-rays by nuclear transitions. We note that somewhat less than half of the total ⁴⁴Ti mass (namely, roughly 43%) is located within the 50% nickel sphere (Figure 7), consistent with the fact that the velocity distribution of ⁴⁴Ti is slightly shifted to higher velocities compared to the ⁵⁶Ni distribution (see middle panel of Figure 5).

The overall morphological similarity between the NuSTAR image of Cas A (Figure 2 of Grefenstette et al. 2014; ⁴⁴Ti distribution in Figure 9) and our model results (Figure 8) is astounding. We emphasize once more that this similarity is not by construction but by chance, resulting as a consequence of stochastic asymmetries that develop at the origin of the explosion. For this reason the focus here should be on basic properties, not structural details. It is obvious that the NS kick vector points away from the direction in which the highest concentrations of iron-group material and ⁴⁴Ti are expelled (see also Grefenstette et al. 2017). This is fully in line with the gravitational tugboat mechanism for NS acceleration in asymmetric SN explosions as discussed by Wongwathanarat et al. (2010b, 2013) and Scheck et al. (2006) (see also Nordhaus et al. 2010, 2012).

This kick mechanism accelerates the NS over timescales of seconds mainly by gravitational forces between the anisotropic and clumpy innermost SN ejecta rather than by direct hydrodynamic interaction. In the context of the neutrino-driven explosion mechanism, which was investigated by 3D SN simulations by Wongwathanarat et al. (2010b, 2013), hydrodynamical instabilities such as convective overturn and the SASI naturally create large-scale asymmetries behind the SN shock in the prelude of the explosion. Due to the aid of these instabilities, neutrino-energy deposition revives the stalled SN shock and initiates its outward expansion anisotropically. In the direction of the more powerful shock acceleration, driven by the largest high-entropy plumes of neutrino-heated matter, more neutrino-processed plasma is expelled and, in addition, stronger shock heating of the stellar shells swept up by the expanding blast wave allows for more efficient explosive nucleosynthesis. This leads to enhanced production of elements from roughly ²⁸Si to the iron group (and beyond) in the direction of the stronger explosion. The corresponding hemispheric differences in the total yields of these elements can be as large as factors of 3–4 for highly asymmetric explosions, and even greater differences might well be possible for more extreme cases than those obtained in the set of 3D simulations of Wongwathanarat et al. (2013), or for nuclear species not included in their nucleosynthesis treatment with a small α-network.

In the direction where the explosion and the nucleosynthesis are weaker, the shock and the postshock material accelerate outward more slowly. On this side the nascent NS can therefore accrete for a longer period of time before the mass infall is quenched by the acceleration of the SN explosion. The momentum transfer both by hydrodynamic accretion flows and by the gravitational attraction from more inert, typically denser and more massive innermost SN ejecta pulls the NS to the side of the weaker blast wave. Thus, the NS receives a kick opposite to the direction of the stronger shock expansion, consistent with momentum conservation during the explosion (for detailed discussions, see Scheck et al. 2006; Wongwathanarat et al. 2013).

This means that the NS experiences a recoil acceleration that points away from the hemisphere where the SN ejects more heavy elements with atomic numbers $Z\sim 14$ and higher. The kick velocity of the NS depends on the stochastic explosion asymmetry, the explosion energy, and the density around the newly formed NS (see Janka 2017). The environmental density of the NS, in turn, depends on the compactness of the progenitor core and determines the amount of matter that is neutrino-processed and shock-heated in the region of explosive nucleosynthesis.

From the set of 3D simulations by Wongwathanarat et al. (2013), model W15-2, which we consider in the present paper, develops a fairly large explosion asymmetry, and the NS receives a kick of 575 km s⁻¹ until ∼3 s after the onset of the SN blast (with acceleration continuing on a low level for an even longer time). Information on hemispheric differences of the yields of nuclei included in the α-network can be found in Table 3 of Wongwathanarat et al. (2013), and visualizations of the 3D distribution of the ejected nickel (compared to other models) can be found in Figures 14 and 15 of that paper. Even for a case with rather high NS kick like model W15-2, radionuclei as nucleosynthesis products in the innermost SN ejecta, where the explosion asymmetry is most extreme, are not just expelled in one hemisphere, but some of this material can be ejected also on the side of the kicked NS. The exact geometry, however, strongly varies from case to case, and lower NS kicks go hand in hand with more isotropic ejection of radioactive material (compare the cases displayed in Figure 14 and 15 of Wongwathanarat et al. 2013).

In Figure 10 we provide a closer comparison of the 3D distributions of ⁴⁴Ti and ⁵⁶Ni including volumetric information. The images of this figure show, for both nuclear species, isosurfaces corresponding to different values of the mass fraction. These isosurfaces were determined such that they enclose 50%, 75%, 90%, and 97% of the total mass ejected of the considered nucleus.

Zoom In Zoom Out Reset image size

The plots confirm our previous conclusions drawn on the basis of the mass distributions in velocity space (Section 3.2): since ⁴⁴Ti and ⁵⁶Ni are nucleosynthesized in close spatial proximity in regions of (incomplete) silicon burning and the α-particle-rich freeze-out (Section 3.1), they are expelled in close connection and there is no process at work that could decouple or decompose them in the ejecta. Their distributions therefore closely resemble each other, and the two nuclei, overall, trace the same 3D geometry.

The bulk of the nickel and titanium is concentrated in relatively small, highly enriched clumps and knots that contain half of the ejecta masses of these nuclei (see left panels of Figure 10). The majority of these clumps expand away from the center of the explosion in the hemisphere opposite to the direction of the NS motion. A close inspection of the two left panels in comparison reveals that the lumps containing 50% of the ejected ⁴⁴Ti are considerably more extended than those of ⁵⁶Ni, meaning that titanium is clearly more diluted. Moreover, one can find regions of high ⁵⁶Ni concentration and little ⁴⁴Ti and vice versa. Note that all of the clumps with radial velocities of ≲4000 km s⁻¹ (blue and light-blue colors) in the left two columns of Figure 10 are located within the volume of the inner sphere that is assumed not to be reverse-shock heated in the right panel of Figure 8. The radius of this sphere is about half the distance from the center to the outermost tips of the largest Ni and Ti fingers. The maximum velocities of titanium in this inner, unshocked sphere are in the ballpark of the fastest material seen by NuSTAR (Grefenstette et al. 2014). Our model suggests that a considerable amount of ⁴⁴Ti should exist in Cas A with higher velocities outside of the reverse-shock radius, which is in line with the 3D data published recently by Grefenstette et al. (2017).

The isosurfaces enclosing 90% and 97% of the ejecta masses for the two nuclei (right two columns in Figure 10) exhibit differences only in details. One example is a larger nickel feature at the 3:30 o'clock position (i.e., right of the root point of the kick vector) that contains little titanium. This feature is even better visible in the second column of Figure 10. Another example are the tips of the most extended, longest finger-like structures, which are more inflated and show higher velocities for titanium. Indications of this trend can already be observed in the images of the left and second columns.

The series of plots in both rows also show significant amounts of low-velocity radionuclei distributed more homogeneously and in a more dilute, volume-filling manner in the central part of the 3D structure formed by the ⁴⁴Ti and ⁵⁶Ni ejecta. Again, most of this material expands away from the center in the hemisphere pointing opposite to the NS kick. A major exception is the big feature between the six o'clock direction and the 9:30 o'clock direction, which mostly shares the same hemisphere with the NS kick velocity vector. One might diagnose a more knotty ⁴⁴Ti NuSTAR map instead of the more homogeneous distribution of titanium with radial velocities of ≲4000 km s⁻¹. Note, however, that our simulations do not account for the long-time effects of radioactive decay heating of nickel and titanium, which might have an influence on porosity and clumping in the titanium distribution 340 yr later.

The "three-finger" geometry of the radioactive ejecta displayed in Figures 10 and 11 reminds one of the three extended iron-rich structures (called "shrapnels" by Orlando et al. 2016) that are seen outside of the reverse-shock radius in Cas A (e.g., Hwang et al. 2004). In order to explore this morphological resemblence of our 3D model W15-2-cw-IIb with Cas A in more detail, we compare the properties of the innermost SN ejecta with the 3D visualization of the observational data of Cas A available at the Smithsonian X3D website http://3d.si.edu/explorer?modelid=45 (see Figure 12) and at http://www.cfa.harvard.edu/~dmilisav/-casa-webapp/model.html. Our main focus in this section is on the basic aspects of the ⁵⁶Ni (or, after decay, Fe) distribution, which carries direct information on the geometry of the explosion in its earliest phase.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In Figure 11 we present images of the iron and ⁴⁴Ti distributions from four different observational perspectives corresponding to the "front view" (left panels), "top view" (second column), "left view" (third column), and "back view" (right panels) defined by the Smithsonian visualization (Figure 12). For optimal resemblance we chose the front view similar to the observer direction of Figures 8 and 10, but with a slightly different inclination of the NS kick vector relative to the line of sight. Instead of being perpendicular to the line of sight as in the previous figures, the kick vector points toward the observer and out of the plane of the page with an angle of about 37°. We note in passing that the corresponding velocity component of the NS motion perpendicular to the line of sight is around 460 km s⁻¹. This is in the ballpark of the transverse velocity of the NS in Cas A, which was estimated to be about 350 km s⁻¹ (Fesen et al. 2006b) and thus typical of young pulsars (e.g., Lyne & Lorimer 1994; Arzoumanian et al. 2002).

The large angle of $53^\circ$ between the NS kick and the observer direction for the front view of our model implies that the hemispheric asymmetry of the ejecta is mapped well from this perspective, where most of the ⁴⁴Ti and ⁵⁶Ni mass is seen in the northern part, while the NS velocity vector points southward. The distribution of the radionuclei along the line of sight is far less one-sided. This can be concluded from the "top view" panels in the second column and the "left view" panels in the third column of Figure 11, where one of the lower images permits the view on the NS kick vector. There are roughly similar amounts of matter moving toward the observers (who are located in the $(-y)$-direction) and away from them. This seems to be compatible with the observed redshift of the ⁴⁴Ti line in Cas A. The corresponding bulk line-of-sight Doppler velocity was determined to be 1100–3000 km s⁻¹ (Grefenstette et al. 2014). Boggs et al. (2015) estimated a "look-back" redshift velocity (associated with the finite lifetime of ⁴⁴Ti and a longer light-travel time from the far side of Cas A) of ∼1400 km s⁻¹ for the age and expansion speed of the SNR. They reasoned that the look-back effect is consistent with the measured redshift of the ⁴⁴Ti line and that the spatially integrated ⁴⁴Ti spectrum alone does not support any statistically significant asymmetry. Rather than having a large intrinsic mass distribution asymmetry along the line of sight, Cas A reveals major asymmetry by the spatial brightness distribution of the ⁴⁴Ti emission (Boggs et al. 2015). This look-back argument has been revised by the recent 3D analysis of Grefenstette et al. (2017). Nevertheless, the overall conclusions of these authors on the ejecta asymmetry in Cas A are in line with those of Boggs et al. (2015) (see below).

Besides the similarities between our model and Cas A discussed already, the different views in Figures 11 and 12 reveal yet another resemblance. The three extended, iron- and titanium-containing fingers lie, by chance, in the same plane. For this reason, the ejecta fill a much narrower volume when looked at from the top and left views (panels in the two columns in the middle of Figure 11). The configuration is reminiscent of Cas A, where reverse-shock-heated metals (Si, Ar, Ne, S, O) are arranged in a thick, ring-like belt girding the reverse-shock sphere, and the unshocked interior ejecta (seen in infrared via Si ii; DeLaney et al. 2010) are apparently assembled in a "tilted thick disk" (Grefenstette et al. 2017). This outer belt also encompasses the iron-rich regions and the fiducial, Si- and Mg-rich "jets" (see the 3D visualizations of the Cas A observations on the websites mentioned above). Since this peculiar structure of the Cas A ejecta is not perpendicular to the "jets" but, quite the contrary, its central plane cuts the northeast and southwest "jets" of Cas A, rapid rotation seems highly unlikely as an explanation of the Cas A morphology. Adding to this argument, the NS kick in Cas A is perpendicular to the NE–SW direction, which contradicts expectations for a jet-driven explosion (Hwang et al. 2004).

Our long-time simulation demonstrates that the origin of the three iron fingers in our (nonrotating) 3D model W15-2-cw-IIb goes back to three dominant high-entropy plumes arising from postshock convective instability at the time when the explosion was launched (see left panel in the top row of Figure 7 of Wongwathanarat et al. 2015). The morphological resemblance to the iron "shrapnels" of Cas A supports the interpretation that this remnant carries clear fingerprints of the initial explosion asymmetries, and neither a dynamically relevant amount of rotation nor a jet-driven explosion needs to be invoked to explain the overall structure of the iron distribution in the Cas A remnant.

After the submission of our paper, a 3D analysis of the ⁴⁴Ti distribution in Cas A was posted on arXiv by Grefenstette et al. (2017). In basic aspects their findings and conclusions are, again, in remarkably good agreement with our 3D results for model W15-2-cw-IIb as discussed above:

• (i)  
  The estimated initial mass of ⁴⁴Ti, $(1.54\pm 0.21)\,\times {10}^{-4}\,{M}_{\odot }$, with an Fe/Ti ratio around 500 is in the ballpark of our model yields, although our nucleosynthesis calculations suggest that the iron content of Cas A could be somewhat higher than inferred by Grefenstette et al. (2017).
• (ii)  
  In excellent correspondence to our simulated remnant geometry visualized in Figures 10 and 11, the average momentum of the ⁴⁴Ti in Cas A (estimated on grounds of the flux-weighted average of the ejecta velocities) is oriented in the plane of the sky almost precisely opposite to the direction of motion of the central compact object, and the NS is determined to have a significant line-of-sight velocity component toward the observer, even with an angle similar to the one shown in Figure 11.
• (iii)  
  ⁴⁴Ti knots near and exterior to the reverse shock are associated with emission from shock-heated iron (which should mostly be the product of radioactive ⁵⁶Ni synthesized along with ⁴⁴Ti); while the existence of titanium-rich regions without iron emission exterior to the reverse shock is not finally assured, the presence of shock-heated iron without observed titanium requires that the ⁴⁴Ti yield is suppressed relative to the iron yield by at least a factor of two. This is easily compatible with predictions from our model, for which a detailed analysis reveals that the color variations visible in Figure 8 correspond to large-scale variations of the titanium-to-iron ratio up to factors of 2–4 (and small-scale variations with even larger amplitudes).
• (iv)  
  Significant amounts of ⁴⁴Ti and ⁵⁶Ni are diagnosed to lie within the reverse shock of Cas A, as suggested by our simulations. While our preferred model implies that ∼43% of the titanium could be in the unshocked inner sphere (Figure 7), Grefenstette et al. (2017) found roughly 40% to be located clearly interior to the reverse shock. But current observational constraints (in particular of the part of the iron that is still hidden, i.e., unshocked, at the age of the Cas A remnant) are very uncertain, and better information from explosion models requires late-time simulations that follow the effects of radioactive decay heating by ⁵⁶Ni and ⁴⁴Ti, as well as the evolution of the reverse shock from the ejecta interaction with the circumstellar medium.
• (v)  
  Grefenstette et al. (2017) did not see any evidence for ⁴⁴Ti ejecta associated with the blueshifted half of the thick disk of interior, unshocked ejecta, though redshifted titanium regions were found to lie in the redshifted half of the disk. However, the identity and location of the blueshifted knot 20b are disputable, as discussed by Grefenstette et al. (2017) and mentioned in a footnote in Section 3.2. Maybe some part of this matter belongs to the blueshifted half of the thick disk. Different from Boggs et al. (2015), Grefenstette et al. (2017) did not diagnose any significant "look-back" effect (associated with the difference in light-travel time between blueshifted and redshifted ejecta) affecting the results but concluded a moderate asymmetry of the ⁴⁴Ti distribution with some excess of redshifted material (mean line-of-sight velocity of 920 ± 510 km s⁻¹). This seems to be basically compatible with the asymmetric distribution of titanium in our model, as reflected by the line-of-sight component of the NS velocity toward the observer in Figure 11.

In future work we plan to perform a closer comparison of our model predictions and the Cas A data with respect to a wider set of nucleosynthesis products and their spatially inhomogeneous distribution, e.g., in the fiducial Si/Mg-rich, high-velocity, wide-angle "jets" (e.g., Fesen et al. 2006b), in the form of large cavities and "bubbles" in the unshocked inner sphere (Milisavljevic & Fesen 2015), and as crown-like and ring-like features (containing Ar) that also surround the iron-rich "shrapnels" in the Cas A SNR's optical main shell (e.g., Fesen et al. 2001; DeLaney et al. 2010; Milisavljevic & Fesen 2013). For this goal we will continue the present 3D explosion simulations toward the remnant stage. Such calculations will have to include the effects of radioactive decay heating, which may inflate the regions containing the radionuclei against their surroundings. They will also have to take into account the reverse shock from the ejecta-environment interaction, which will have an impact on the radial velocity profile and is likely to modify the detailed shape of the shock-compressed ejecta regions.