A&A 425, 1029-1040 (2004)
DOI: 10.1051/0004-6361:20041108

Nucleosynthesis in multi-dimensional SN Ia explosions

C. Travaglio^1,2 - W. Hillebrandt³ - M. Reinecke⁴ - F.-K. Thielemann⁵

1 - Max-Planck Institut für Astrophysik, Karl-Schwarzschild Strasse 1, 85741 Garching bei München, Germany
2 - Istituto Nazionale di Astrofisica (INAF) - Osservatorio Astronomico di Torino, via Osservatorio 20, 10025 Pino Torinese (Torino), Italy
3 - Max-Planck Institut für Astrophysik, Karl-Schwarzschild Strasse 1, 85741 Garching bei München, Germany
4 - Max-Planck Institut für Astrophysik, Karl-Schwarzschild Strasse 1, 85741 Garching bei München, Germany
5 - Department of Physics and Astronomy, University of Basel, Klingelbergstrasse B2, 4056 Basel, Switzerland

Received 16 April 2004 / Accepted 11 June 2004

Abstract
We present the results of nucleosynthesis calculations based on multi-dimensional (2D and 3D) hydrodynamical simulations of the thermonuclear burning phase in type Ia supernovae (hereafter SN Ia). The detailed nucleosynthetic yields of our explosion models are calculated by post-processing the ejecta, using passively advected tracer particles. The nuclear reaction network employed in computing the explosive nucleosynthesis contains 383 nuclear species, ranging from neutrons, protons, and $\alpha$-particles to ⁹⁸Mo. Our models follow the common assumption that SN Ia are the explosions of white dwarfs that have approached the Chandrasekhar mass ( $M_{\rm ch}\sim 1.39$), and are disrupted by thermonuclear fusion of carbon and oxygen. But in contrast to 1D models which adjust the burning speed to reproduce lightcurves and spectra, the thermonuclear burning model applied in this paper does not contain adjustable parameters. Therefore variations of the explosion energies and nucleosynthesis yields are dependent on changes of the initial conditions only. Here we discuss the nucleosynthetic yields obtained in 2D and 3D models with two different choices of ignition conditions (centrally ignited, in which the spherical initial flame geometry is perturbated with toroidal rings, and bubbles, in which multi-point ignition conditions are simulated), but keeping the initial composition of the white dwarf unchanged. Constraints imposed on the hydrodynamical models from nucleosynthesis as well as from the radial velocity distribution of the elements are discussed in detail. We show that in our simulations unburned C and O varies typically from $\sim$40% to $\sim$50% of the total ejected material. Some of the unburned material remains between the flame plumes and is concentrated in low velocity regions at the end of the simulations. This effect is more pronounced in 2D than in 3D and in models with a small number of (large) ignition spots. The main differences between all our models and standard 1D computations are, besides the higher mass fraction of unburned C and O, the C/O ratio (in our case is typically a factor of 2.5 higher than in 1D computations), and somewhat lower abundances of certain intermediate mass nuclei such as S, Cl, Ar, K, and Ca, and of ⁵⁶Ni. We also demonstrate that the amount of ⁵⁶Ni produced in the explosion is a very sensitive function of density and temperature. Because explosive C and O burning may produce the iron-group elements and their isotopes in rather different proportions one can get different ⁵⁶Ni-fractions (and thus supernova luminosities) without changing the kinetic energy of the explosion. Finally, we show that we need the high resolution multi-point ignition (bubbles) model to burn most of the material in the center (demonstrating that high resolution coupled with a large number of ignition spots is crucial to get rid of unburned material in a pure deflagration SN Ia model).

Key words: hydrodynamics - nuclear reactions, nucleosynthesis, abundances - stars: supernovae: general

1 Introduction

Type Ia supernovae (SN Ia) are known to be stellar explosions with no signs of hydrogen and helium in their spectra, but intermediate mass elements such as Si, S, Ca and Mg near the maximum of their light curves, and many Fe lines at later times. In contrast to massive stars which are the progenitors of type II supernovae (SN II), SN Ia progenitors are thought to be white dwarfs (WDs) in binary systems (see Whelan & Iben 1973; and Hillebrandt & Niemeyer 2000, for a more recent review). In the canonical model the WD, expected to consist mainly of carbon and oxygen, approaches the Chandrasekhar mass ( $M_{\rm ch}$) through a not yet known mechanism, presumably accretion from a companion star, and is then disrupted by a thermonuclear explosion. The declining light is powered by the radioactive decay of ⁵⁶Ni. A strong argument in favor of this scenario is given by the fact that these explosion models fit quite well the observed light curves and spectra (Leibundgut 2001).

Despite the consistency of this general framework with observations the detailed theory of how SN Ia evolve and explode is still subject of considerable efforts. Over the last three decades, one-dimensional spherically symmetric models have been used to study the various channels that may give rise to a successful SN Ia in terms of the predicted spectra, light curves, and nucleosynthesis. Much of this work was centered on the  $M_{\rm ch}$ scenario wherein a C+O white dwarf accretes H or He from a binary companion (Nomoto, Thielemann, & Yokoi 1984) and ignites explosive carbon burning just before it reaches a critical mass of $M_{\rm ch}\sim 1.39$ $M_{\odot }$. The subsequent explosion produces enough ⁵⁶Ni ($\sim$0.6 $M_{\odot }$) and intermediate mass elements to reproduce "normal'' SN Ia lightcurves and spectra, provided that the amount of C+O burned at any given density is suitably chosen. In 1D models this can be achieved by parameterizing the thermonuclear flame speed and, if desired, the density at which a transition to supersonic burning (detonation) occurs (Khokhlov et al. 1999; Niemeyer 1999). Moreover, some mixing of processed matter had to be assumed in order to fit the observed spectra. Alternative scenarios, including sub- $M_{\rm ch}$ explosions and merging white dwarfs (double degenerates), have met with mixed success (see e.g., Arnett & Livne 1994).

More recently it has become possible to perform multi-dimensional 2D (Livne 1993; Reinecke et al. 1999; Lisewski et al. 2000) and 3D (Reinecke et al. 2002a,b; Gamezo et al. 2003) simulations of an exploding  $M_{\rm ch}$-white dwarf. The principal difficulty in these models is the fact that the hydrodynamically unstable and turbulent nuclear flame front develops structures on much smaller length scales than can numerically be resolved. However, this problem can be overcome by "large eddy simulations'', i.e. by employing subgrid-scale models for the unresolved scales that provide a guess of the effective turbulent flame speed on the scale of the computational grid (Niemeyer & Hillebrandt 1995a; Röpke et al. 2003). In this flame model, which is well justified in the thin-flame regime and is tested in experiments with premixed turbulent chemical flames, we do not need a detailed prescription of the nuclear reactions. Instead, the fuel consumption rate is roughly propartional to the surface area of the flame front and its normal (turbulent) velocity.

Despite of the need of more detailed studies of such subgrid-scale models, it is important to stress that the multi-dimensional simulations reach a qualitatively different level of predictive power than 1D models. In particular, the amount of material burned at a given density can not longer be fine-tuned but is determined by the fluid motions on the resolved scales and a particular choice of the subgrid model (Reinecke et al. 2002a). Therefore, once the flame model has been fixed numerical simulations of the thermonuclear explosion of a given white dwarf can be done by just choosing the ignition conditions, including the chemical composition of the WD, the only remaining (physical) parameter.

The undeniable influence of SN Ia explosions on, among others, the chemical evolution of galaxies makes the quest for solid theoretical models and nucleosynthetic yields an urgent task. Guided by decades of modeling and nucleosynthesis calculations in spherically symmetric models (the protoype being the W7 model of Nomoto et al. 1984; Iwamoto et al. 1999; Brachwitz et al. 2000), we have begun analyzing the detailed nucleosynthetic yields of our explosion models by post-processing the ejecta. This has been performed adding a "Lagrangian component'' to our Eulerian scheme in the form of tracer particles passively advected with the flow in the course of the Eulerian calculation. Therefore we record their T and $\rho$ history by interpolating the corresponding quantities from the underlying Eulerian grid. A similar method of tracer particles in an Eulerian code to calculate the nucleosynthesis has been adopted in a previous study of multi-dimensional nucleosynthesis in core collapse SNe by Nagataki et al. (1997), and more recently in calculations for very massive stars (Maeda et al. 2002), for core collapse SNe (Travaglio et al. 2004), and for type Ia SNe first preliminary results have been discussed by Niemeyer et al. (2003).

In this paper we present the nuclear yields resulting from several of our multi-dimensional supernova simulations, and we compare them to the standard W7 (Iwamoto et al. 1999; Brachwitz et al. 2000; Thielemann et al. 2003) results. In Sect. 2 we summarize the 2D and 3D SN Ia calculations, discussing different mode of ignition as well as grid resolution of the hydrodynamic code. In Sect. 3 we describe our method to perform nucleosynthesis calculations and the nucleosynthesis network adopted for this work. Finally, in Sect. 4 we present and discuss our nucleosynthesis results. In a first step we have performed resolution studies in 2D consisting of different methods how to distribute the tracer particles, the number of particles used, and grid resolution of the hydrodynamic code. Although 2D simulations cannot be considered to be realistic, as was discussed by Reinecke et al. (2002a), they can serve to guide the more elaborate 3D models. We then discuss the nucleosynthesis resulting from three 3D models, a centrally ignited model and two models with a few and many off-center ignition spots, respectively. It will be shown that the more realistic ignition conditions (central ignition or many ignition spots) also predict nucleosynthesis yields closer to the ones observed in typical SN Ia's.

Concerning nucleosynthesis we will in particular analyze the range and distribution of ⁵⁶Ni masses we are able to produce with our present models, and the sensitivity of the amount and velocity distribution of unburned material (¹²C, ¹⁶O, ²²Ne) to the ignition conditions of the explosions which are still a major uncertainty of SN Ia models.

2 Recent multi-dimensional SN Ia calculations

We have carried out numerical simulations in 2D and 3D, for several different ignition conditions, and for different numerical resolution. Details of these models are given in a series of papers (Reinecke et al. 1999, 2002a,b). A detailed discussion of these models is not the aim of this paper, therefore only a summary of the results essential for our nucleosynthesis calculations will be repeated in this section.

As long as the evolution of the white dwarf before the thermonuclear runaway remains largely unexplored (see recent work by Woosley et al. 2004, and references therein), only very crude constraints can be put on the flame geometry at the onset of the burning. It appears likely that the deflagration sets in at the surface of quietly burning "hot bubbles''. Nevertheless very little is known about the number, size and radial distribution of these hot spots. This is a consequence of the complicated physical processes taking place in the white dwarf's core during the convective smoldering phase prior to ignition lasting for $\sim$1000 years. The long time scales combined with the relatively slow convective motions make numerical simulations of this phase a daunting task which has not been undertaken in its full complexity so far. Theoretical considerations and simplified simulations carried out by Garcia-Senz & Woosley (1995) suggest that fast burning starts on the surface of many small bubbles ($r\le 5$km), within 100 km of the star's center. Central ignition is another possible scenario that has been investigated during the last years using multidimensional calculations (Niemeyer & Hillebrandt 1995b; Khokhlov 1995).

In this work we present four models: c3_2d_512, a 2D model with central ignition, and grid size of 512²; c3_3d_256, a 3D model with central ignition, and grid size of 256³; b5_3d_256, a 3D model with ignition in 5 bubbles, and grid size of 256³; finally, b30_3d_768, a 3D model with ignition in 30 bubbles, and grid size of 768³. This last one is the model with the highest resolution possible to evolve with the computer resources available to us, therefore we will consider it as the "standard'' model for this paper. It achieves a central resolution of 3.33 km, using a grid consisting of 768³ zones. In the simulated octant of our model b30_3d_768, 30 bubbles with a radius of 10 km were distributed randomly. The bubble locations were drawn from a Gaussian probability distribution with a dispersion of $\sigma = 75$ km. Bubbles located more than 2.5$\sigma$away from the center were rejected. In all models we started the simulations with a central ignition density of $2.9\times10$⁹ g/cm³. The simulations have been followed up to 1.5 s. for all the models, except for the b30_3d_768. Due to a very high consumption in computer time the model b30_3d_768 was stopped when no further energy was released.

Figure 1 shows the energy release for the four models mentioned above. The curves are nearly identical during the first $\sim$0.5 s of the simulation. Owing to the small volume of the bubbles, the initial hydrostatic equilibrium is only slightly disturbed. During the first stages the energy release is therefore lower than in previous simulations. Only after the total flame surface has grown considerably (mostly by deformation of the bubbles), vigorous burning sets in. In the late explosion phase (after about 0.5 s) the total energy differs for the four simulations, and increases moderately with increasing resolution (see Reinecke et al. 2002a,b, for more detailed discussion). We also note that the 3D centrally ignited and the five-bubble models are remarkably similar, even if the centrally ignited has a relatively faster burning between $\sim$0.5 and $\sim$1 s.

	Figure 1: Total energy evolution for the two-dimensional centrally ignited explosion model ( solid line), for the three-dimensional low-resolution centrally ignited ( dashed line) explosion model, for the three-dimensional low-resolution 5 bubbles ( dotted-dotted-dashed line), and for the three dimensional high-resolution 30 bubbles ( dotted-dashed line) explosion model.
Open with DEXTER

It must be noted that the five-bubble model is not identical to the model b5_3d_256 presented by Reinecke et al. (2002b): due to an oversight during the simulation setup the initial positions of the burning bubbles are not the same. As a consequence, the total energy releases of these two simulations are slightly different.

The initial configuration of the front, as well as snapshots of the front evolution at later times are shown in Figs. 2 and 3 for the models c3_3d_256 and b30_3d_768, respectively.

	Figure 2: Snapshots of the front evolution for the centrally ignited model c3_3d_256 at 0 s, 0.2 s, 0.4 s, and 0.6 s.
Open with DEXTER

	Figure 3: Snapshots of the front evolution for the floating-bubble model b30_3d_768 at 0 s, 0.1 s, 0.14 s, and 0.2 s.
Open with DEXTER

In addition to the total energy release, the mass fraction of unburned material in the central region of the remnant appears to be a good criterion for judging the validity of our simulations, because a high amount of C and O in this region would most likely produce a characteristic signature in the late-time spectra which has not yet been observed. In this respect the results of our earlier calculations were not very encouraging since the ashes rose towards the surface in large structures and left nearly pure fuel in the center. Using many initial bubbles, however, seems to alleviate this problem insofar as the statistical isotropy of the initial flame at least delays the development of large-scale turnover motions. As a consequence, C and O is lower than 20% in the central 0.2 $M_{\odot }$ after 0.9 s for model b30_3d_768. For this last model we get $\sim$40% of the total mass stays unburned (we define unburned the material with $T < 1.5 \times 10^9$ K). We will discuss this point in more detail when we will present our nucleosynthesis results.

3 Nucleosynthesis in multi-dimensional SN Ia

The multidimensional SN Ia simulations described in Sect. 2 employed a minimal nuclear reaction network, directly included in the hydrodynamic code, sufficient for a good approximation of the thermonuclear energy release and the predicted gross chemical composition agrees well with the expectations (Reinecke et al. 2002b). It consists of five species ($\alpha$-particles, ¹²C, ¹⁶O, ²⁴Mg, and ⁵⁶Ni) and is intended to model the energy release of the thermonuclear reactions only. No reaction rates are employed: all material behind the flame front is instantaneously converted to a NSE of ⁵⁶Ni and $\alpha$-particles at high densities and to ²⁴Mg at intermediate densities. Below 10⁷ g/cm³, no burning takes place.

Therefore we follow a minimal reaction network directly in the hydrodynamic simulations, and a much more extended network in a post-process step. Here we present the results of the more detailed study of the nuclear abundances in the ejecta obtained by post-processing the output of the four hydrodynamic models discussed above.

Since the multidimensional hydrodynamics scheme used in modeling the explosions is an Eulerian one (i.e. the grid does not move with the fluid), in order to record temperature and density evolution as a function of time (the necessary input for nucleosynthesis calculations) we homogeneously distributed $\sim$10 000 marker particles (in 2D models) and $\sim$20 000 marker particles (in 3D models) and followed their T and $\rho$ evolution. The number of particles in the simulation was chosen in order to reproduce in the best way the resolution of the grid (see discussion below). We then calculated the nucleosynthesis experienced by each marker and computed the total yield as a sum over all the markers, after the decay of unstable isotopes.

3.1 Tracer particles method

In one spatial dimension it is nowadays possible to solve reaction networks consisting of hundreds of species online with the hydrodynamics (see e.g., Rauscher et al. 2002, for explosive nucleosynthesis calculations in core collapse SNe). However, it is more common to use reduced networks in order to obtain the (approximate) energy generation rate for the hydrodynamics and to calculate the detailed chemical composition only afterwards in a post-processing step. This is facilitated by the Lagrangian nature of nearly all 1D codes employed for explosive nucleosynthesis calculations. In Lagrangian schemes, the grid moves with the fluid and therefore it is possible to record the evolution of the temperature and density for different fluid elements (i.e. Lagrangian zones), which is required for the post-processing. In contrast, most grid-based multidimensional hydro schemes are of Eulerian type (i.e. the grid is fixed in space). To obtain the necessary data for the post-processing we added a "Lagrangian component'' to our Eulerian scheme in the form of marker particles that we passively advect with the flow in the course of the Eulerian calculation, recording their T and $\rho$ history by interpolating the corresponding quantities on the underlying Eulerian grid.

In the 3D simulations the star is subdivided into 27³ grid cells equidistant in the integrated mass M(r), azimuthal angle $\varphi$ and cos $\theta$, so that each grid cell contains the same mass. A tracer particle was placed randomly in each of those grid cells, therefore the total number of tracer particles we used is 27 ³ = 19 683. After numerical inversion of the function M(r), the (M, $\varphi$, cos $\theta$) coordinates are mapped onto the Cartesian grid. For the 2D simulation, 100² particles are distributed in rand cos $\theta$ directions, using the same procedure as above. In all cases the simulation covers one octant, therefore to get the total white dwarf mass we mutiply the mass of each marker by 8 and we sum over all the markers. They are distributed in the way to have identical mass, that is therefore calculated as the ratio between  $M_{\rm ch}$ and the total number of tracer particles. The mass of each marker is therefore between $\sim$10 $^{-3}~ M_\odot$ and $\sim$10 $^{-4} ~M_\odot$. The initial distribution of the markers for the b30_3d_768 case is shown in Fig. 4.

	Figure 4: Radial distribution of the tracer particles in the 3D model at the beginning of the simulation.
Open with DEXTER

Finally, we compare the T(t) values for each marker extracted from the hydrodynamical model with the T(t) derived from using the internal energy from the hydro-code (which has contributions from the Boltzmann gas of ionized nuclei, the Planck-spectrum of photons, and relativistic degenerate or non-degenerate electrons positrons). This is done by calculating separately the equation of state for a given maker particle's internal energy, density and composition, and deriving from that the temperature T. As chemical composition for this calculation we use a 16 isotopes network, composed by n, p, ⁴He, ¹²C, ¹⁶O, ²⁰Ne, ²⁴Mg, ²⁸Si, ³²S, ³⁶Ar, ⁴⁰Ca, ⁴⁴Ti, ⁴⁸Cr, ⁵²Fe, ⁵⁶Ni, ⁶⁰Zn. We found that the T(t) obtained directly from the hydrodynamic model is in average lower by 10% up to 20% for markers with high temperatures ( $T \ge 6 \times 10^9$ K) as compared with the temperature derived from the energy. This can be understood from the fact that in the hydro-code the internal energy density is the directly computed variable and, therefore, is more accuratly determined than the temperature. Consequently the more precise T(t) distribution is derived from the internal energy and the actual composition. This is what we used for our nucleosynthesis calculations.

3.2 Nucleosynthesis network

The nuclear reaction network employed in computing our post-process explosive nucleosynthesis calculations contains 383 nuclear species ranging from neutrons, protons, and $\alpha$-particles to molibdenum. A detailed description of the code we used to solve the nuclear network and the reaction rate library utilized is given by Thielemann et al. (1996) and Iwamoto et al. (1999). Weak interaction rates applied in those calculations were taken from Fuller et al. (1985). More recently full large-scale shell model calculations for electron capture and $\beta$-decays became available also for pf-shell nuclei, i.e. the Fe-group (Langanke & Martinez-Pinedo 2000; Martinez-Pinedo et al. 2000). They have already been included in preliminar calculations by Brachwitz et al. (2000) and Thielemann et al. (2003). We also included the new rates in the calculations presented in this paper. As discussed below more in details, the nuclear reaction rates entering the thermonuclear modeling can play an important role. While large portions of the ejecta which experience maximum temperatures in excess of $6\times10^9$ K follow nuclear statistical equilibrium (a chemical equilibrium of all strong and electromagnetic reactions), weak interactions occur on a longer timescale and different choice of Fueller et al. (1985) (as used by Iwamoto et al. 1999) or Langanke & Martinez-Pinedo (2000) (as used for this work; and by Brachwitz et al. 2000; Thielemann et al. 2003), can strongly affect the results.

The initial WD composition we used consists of (mass fraction) 0.475 $M_{\odot }$ of ¹²C, 0.5 $M_{\odot }$ of ¹⁶O, and 0.025 $M_{\odot }$ of ²²Ne (in agreement with the W7 initial composition, Iwamoto et al. 1999). With this initial composition we typically simulate a solar metallicity SN Ia.

When the flame passes through the fuel, ¹²C, ¹⁶O and ²²Ne are converted to ashes with different compositions depending on the intial T and $\rho$. We stop our nucleosynthesis calculations after $\sim$1.5 s. When the temperature in the markers dropped at $\sim$ $1.5 \times 10^9$ K the explosive nucleosynthesis is almost frozen. The distribution of T and $\rho$ as a function of time is shown in Fig. 5 for the model b30_3d_768. The thick lines represent the upper envelope for T and $\rho$, and with the dashed lines we plot T and $\rho$ histories in some markers taken as examples. As one can notice from Fig. 5, the T distribution is not yet below $\sim$ $1.5 \times 10^9$ K for all the tracer particles. At $T \sim 3$- $4\times 10^9$ K (i.e. the upper values shown in Fig. 5) we still expect some explosive C-burning products. Therefore we linearly extrapolate T and $\rho$ until all the particles have $T < 1.5 \times 10^9$ K (that corresponds to $\sim$1 s).

	Figure 5: Temperature ( upper panel) and density ( lower panel) history in the tracer particles for the b30_3d_768 model. The thick lines represent the upper envelope of the distribution and the thin dashed lines represent some of the markers randomly taken as examples.
Open with DEXTER

The combination of T and $\rho$ vs. time in each marker is very important for the nucleosynthesis results (as discussed below). Comparing the model shown in Fig. 5 with the T and $\rho$distribution shown by Iwamoto et al. (1999) for their W7 model, we note differences that can be interesting for the nucleosynthesis calculations. First, the timescale in our models are much smaller ($\sim$1.5 s) with respect to the W7 timescale ($\sim$6 s). Then different combinations of T and $\rho$, in our case rather low T at still high $\rho$, also give us interesting differences in the nucleosynthesis calculations.

4 Discussion and results

In this section we present the results for nucleosynthesis calculations in the model c3_2d_512, c3_3d_256, b5_3d_256, and our "standard'' b30_3d_768. We also compare them with the W7 calculations by Brachwitz et al. (2000) and Thielemann et al. (2003) (note that in figures and tables we better compare with Brachwitz et al. calculations instead of Iwamoto et al. 1999, in order to be consitent with the use of elecetron capture rates). We analyze the consequences of different hydrodynamic resolutions on the nucleosynthesis, we compare their different distribution of burned and unburned material, we discuss the trend of the $Y_{\rm e}$in the markers as a consequence of our nucleosynthesis calculations. Finally, we discuss the velocity distribution of different nuclear species.

4.1 Nucleosynthesis and yields: Comparison between 2D and 3D

For the nucleosynthesis calculations the peak temperatures combined with the density distribution achieved during the propagation of the front, are the most important quantities. As shown in Fig. 5, for the b30_3d_768 model the maximum of the Tin the markers covers a large range ( $1.5\times 10^9~{\rm K} < T < 8.4\times 10^9$ K). Note that in this figure each dotted line represents a tracer particle trend (we selected randomly some markers for this plot), and the thick line is the upper envelope of the T and $\rho$ distribution in the total tracer particle sample. At high temperatures all strong interactions and photodisintegrations are so fast that a chemical equilibrium (nuclear statistical equilibrium, NSE) is immediately achieved (in our calculations we assume NSE condition for $T \ge 6.0\times 10^9$ K). The resulting chemical composition is therefore only dependent on $\rho$, T, and the neutron eccess (that is determined by the total amount of electron captures taking place on longer timescales). An example of the behaviour of the chemical abundances in one marker as a function of time is shown in Fig. 6. For this case the marker is originally located at a radius of $\sim$150 km, i.e. in the innermost dense zone; the initial density at the position of the marker is $\sim$ $2.5\times 10^9$ g/cm³ (the central density for this model is $2.9\times 10^9$ g/cm³) and it is heated by the flame front almost immediately (at $\sim$0.1 s). The NSE conditions are achieved very fast and the temperature reaches a quite high peak $T\sim 8.5\times 10^9$ K (Fig. 6, upper left panel). The resulting electron fraction $Y_{\rm e}$ drops rapidly (in $\sim$0.3 s) from the initial value of 0.4989 to 0.468 (Fig. 6, upper right panel). The main O-burning products, ²⁸Si and ³²S are shown in Fig. 6 (left middle panel), together with the abundance of ²⁴Mg resulting from C-burning. Furthermore, the rise of the temperature in excess of $6\times10^9$ K leads to a complete NSE and ²⁸Si exhaustion within 0.1 s. Due to the high density (typical of explosive Si-burning) normal freeze-out occurs. One of the main products is ⁵⁶Fe (Fig. 6, right middle panel). Since T and $\rho$conditions are very high, also neutron-rich nuclei are built up due to electron captures, and ⁵⁶Fe is partly replaced by ⁵⁴Fe (Fig. 6, right middle panel) and ⁵⁸Ni (Fig. 6, left lower panel). In such details we can follow the nucleosynthesis changes along the time in all the tracer particles.

Figure 6: Example of the nucleosynthesis calculation in one tracer particle. T and $\rho$ are plotted in the upper left panel; with $T > 6\times 10^{9}$ K reaches NSE conditions. In the upper right panel the resulting $Y_{\rm e}$ is shown. The other panels give the mass fraction vs. time for 16O, 28Si, 32S, 54,56Fe, 56,57,58Ni, 63Cu, and 64Zn. For 16O, 28Si, 32S the time is plotted only up to 0.2 s since their abundances are zero at later times.

Open with DEXTER

In Fig. 7 we show the distribution of $Y_{\rm e}$ vs. time obtained as a result of the nucleosynthesis calculation in the model b30_3d_768. The two thick lines stand for the upper and lower values of the $Y_{\rm e}$ in the markers. The dashed lines represent the $Y_{\rm e}$ time evolution in some markers randomly taken as examples. The range covers values from $\sim$0.5 (that represent the typical initial composition as well as the composition of markers with rather low T), down to $\sim$0.462, reached by the markers with the highest T (see e.g. the example in Fig. 6).

	Figure 7: The same as Fig. 5 for the history of the $Y_{\rm e}$ in the markers.
Open with DEXTER

	Figure 8: Distribution of the tracer particles vs. radius at $\sim$1.2 s. for b30_3d_768 ( upper panel), c3_3d_256 ( middle panel), and b5_3d_256 ( lower panel) models. The unburned particles ( $T < 1.5 \times 10^9$ K) are plotted with a solid line, and the processed particles with a dashed line.
Open with DEXTER

Table 1 lists the synthesized masses for all the stable isotopes up to ⁶⁸Zn for the models c3_2d_512, c3_3d_256, b5_3d_256, and b30_3d_768. For comparison we also include in Col. 2 the calculations for the W7 model (from Thielemann et al. 2003). For the calculations presented in Cols. 5-7 we include the nucleosynthesis results starting only when the temperature has reached 90% of the peak. In fact, due to finite numerical resolution in the hydrodynamic simulation, the rise of the temperature is not as steep as it would be in reality. Consequently, when markers reach NSE conditions and weak-interactions start to play the most important role, the nucleosynthesis timescales are fast, and even 0.1 s (i.e. typical timescale we have for the rise of T) are crucial for some reactions to give an important contribution. Nevertheless, as one can see from the table comparing Col. 4 (the nucleosynthesis for the c3_3d_256 model has been calculated considering the all rise of the T curve) and Col. 5 (the nucleosynthesis for the same c3_3d_256 model has been calculated only when T has reached 90% of the peak), the differences in the total yields are very small. This is due to the fact that the amount of markers affected by this inaccuracy is a small fraction of the total. In the table we do not include isotopes heavier than ⁶⁸Zn, even if the network we used was extended up to ⁹⁸Mo, since their resulting mass fraction are smaller than $\sim$10 $^{-15} ~M_\odot$. An important thing to notice is the difference in the amount of unburned material (as defined at the beginning of this paper, we consider unburned the material with $T < 1.5 \times 10^9$ K) in the four models. For the 2D model c3_2d_512 $\sim$60% of the total material remains unburned, instead of 40% we obtain for the 3D model c3_3d_256. Both of these models are centrally ignited, therefore the difference in the amount of unprocessed material is mainly a consequence of the difference in the total energy distribution due to multi-dimension effects (see Fig. 1). In Fig. 8 we show the distribution of burned and unburned material for the model b30_3d_768 (upper panel), c3_3d_256 (middle panel), and b5_3d_256 (lower panel), in terms of number of tracer particles. The distribution plotted in Fig. 8 has been taken at 1.1 s for all the three models, that is the final time we reached in the b30_3d_768. For all these three simulations the burned component dominates at a radius of $\sim$ $5.0\times 10^7$ cm, instead the unburned material is more or less uniformly distributed, with a dominant component in the outermost zones and a tail in the central regions. As also discussed by Reinecke et al. (2002b) the distribution of the unprocessed material depends on the initial conditions for the burning. When the model is centrally ignited, as c3_3d_256, the dominant component of unburned material is in the outermost regions. Instead with a floating-bubble model with comparable resolution, like b5_3d_256, unprocessed material can also be concentrated in the center. Nevertheless the amount of unprocessed material in the center for a floating-bubble model depends on the amount of ignition-spots together with the resolution used. In fact (see Fig. 8) for the b30_3d_768 case most of the unburned component in the innermost regions disappeared.

We also tested consequences of burning for longer times. As mentioned in Sect. 3.2, at $\sim$1.2 s the upper limit of the temperature is $\sim$ $4\times 10^9$ K (see also Fig. 5), and explosive C- and a bit of Ne-burning can still occur. Therefore, to calculate the nucleosynthesis we extrapolated temperature and density for $\sim$1 s further (when the upper values for T are not higher that $1.5 \times 10^9$ K). If we compare the nucleosynthetic yields calculated at the end of the hydro simulation with the nucleosynthetic yields calculated using the extrapolated values of T and $\rho$, we obtain that only few % of ¹²C burns. The consequences can be relevant only for few isotopes, the main products of explosive C-burning, i.e. ^20,21Ne, ²³Na, ^25,26Mg, ²⁷Al. Therefore only for those isotopes we expect important changes if we could follow the hydrodynamic simulation for longer timescales.

It is interesting to notice the difference in the ¹²C and ¹⁶O abundance of the W7 model and our multi-dimensional models. The ¹²C of the W7 is about a factor of $\sim$10 lower than in our cases, instead ¹⁶O of W7 is only a factor of $\sim$3 lower. While ¹⁶O and ¹²C we obtain are built up by unprocessed particles, in the W7 C-burning is more efficient and burns a significant amount of ¹²C at low T and $\rho$, with a resulting different C/O ratio. Also, the ⁵⁶Fe mass (mainly deriving from the decay of the long-lived ⁵⁶Ni) obtained in the W7 model (0.696 $M_{\odot }$, Thielemann et al. 2003) is by far higher than the ⁵⁶Fe mass resulting from the multi-dimensional SN Ia models. The highest value we can reach in our models is obtained with the highest resolution floating-bubble model b30_3d_768 (0.44 $M_{\odot }$). We note that the initial conditions (in this case 30 ignition spots) are crucial for a more precise study of the nucleosynthesis, in particular of the innermost regions. Possibly, a model with even higher number of ignition bubbles would give us a still higher ⁵⁶Fe mass.

We note that neutron-rich isotopes like ⁴⁸Ca, ⁵⁰Ti, ⁵⁴Cr, ⁵⁴Fe, ⁵⁸Ni are strongly underproduced with respect to the W7 model presented by Iwamoto et al. (1999), instead are in good agreement with the Branchwitz et al. (2000) and Thielemann et al. (2003) predictions. This is due to the differences in the electron capture rates adopted (as just discussed above). As a consequence, electron capture rates of nuclei affects directly the electron fraction $Y_{\rm e}$(Thielemann et al. 1986). As shown in Fig. 7, the lowest $Y_{\rm e}$ reached in the markers is 0.462, instead the lowest value in the W7 model by Iwamoto et al. (1999) is $\sim$0.446. For the model by Iwamoto et al. (1999) low $Y_{\rm e}$ values (<0.46) are reached in the innermost zones, i.e. $M < 0.03 ~M_\odot$, where the highest temperatures are reached ($\sim$ $9.0 \times 10^9$ K). The highest T reached in our models, as shown in Fig. 5, are $\sim$ $8.4\times 10^9$ K with consequently higher $Y_{\rm e}$.

Finally, in Table 2 we report the synthesized masses of the main radioactive species from ²²Na up to ⁶³Ni. Bigger differences between our models and W7 are the abundances of ⁴⁸Ca and ⁶⁶Zn (in this case we use for comparison the W7 data from Iwamoto et al. 1999; since the Thielemann et al. (2003) data are not yet available for the long-lived isotopes). As discussed in Thielemann et al. (2003), these isotopes are very sensitive to small variations in the central density of the model. We used $2.9\times 10^9$ g/cm³ instead of the $2.0\times 10^9$ g/cm³ of the model B2C20 presented by Thielemann et al. (2003) and used in the current paper as a comparison to our model.

	Figure 9: Nucleosynthetic yields (in mass fraction normalized to the solar value and to the corresponding solar ratio) obtained using 19 683 tracer particles in the 3D model b30_3d_768 compared to the W7 yields given by Thielemann et al. (2003).
Open with DEXTER

In Fig. 9 we show the yields obtained for our "standard'' model b30_3d_768 compared with the W7 yields (Brachwitz et al. 2000; Thielemann et al. 2003), scaled to their relative solar abundances and to the ⁵⁶Fe abundance. As one can see from the figure, a part from the difference in the relative ⁵⁶Fe mass (0.44 $M_{\odot }$ for the b30_3d_768 and 0.669 $M_{\odot }$ for W7), as well as in the unburned material (i.e. ¹²C, ¹⁶O and ²²Ne), the trend for the production of different isotopes is quite similar.

In Fig. 10 we plotted the yields of the c3_2d_512 and c3_3d_256 (upper panel), and b5_3d_256 (lower panel), normalized to the "standard'' model b30_3d_768. From Fig. 11 one can notice that with a similar amount of burned material (the difference between the two models in the total burned material it is not more than 10%), the b5_3d_256 produces much more $\alpha$-elements. Infact the higher efficiency of ¹²C burning is clearly shown by a higher production of ²³Na and ⁴⁰Ca. Also the b5_3d_256 model has a lower ⁵⁶Fe production (0.34 $M_{\odot }$, instead of 0.44 $M_{\odot }$ of the b30_3d_768), and Fe-group elements.

	Figure 10: Nucleosynthetic yield ratio, comparing the model c3_2d_512 ( dotted line), and c3_3d_256 ( solid line), with our standard model b30_3d_768 ( upper panel). Nucleosynthetic yield ratio, comparing the model b5_3d_256 with our standard model b30_3d_768 ( lower panel).
Open with DEXTER

4.2 Radial velocity distribution

In Fig. 11 we show the mass fractions of selected isotopes as a function of the radial velocity, taken at the end of our simulation (i.e. $\sim$1.2 s) for the model b30_3d_768 and c3_3d_256. With thick lines we plot the unprocessed material in the form of ¹²C, ¹⁶O, and ²²Ne. We note that at the lowest velocities ($\sim$1000 km s^-1) the dominant component is represented by the unburned material (i.e. ¹²C and ¹⁶O) for the floating-bubble high-resolution model, and by ⁵⁶Fe for the centrally ignited model. On the opposite site, at the highest velocities (>10 000 km s^-1) the unburned material dominates in the centrally ignited model, instead is comparable to the ⁵⁶Fe component in the floating-bubble high resolution model. Maximum velocities reached are about 12 000 km s^-1 in both cases. As mentioned in the previous section, our models at 1.2 s are not yet in homologous expansion, i.e. pressure and gravity still play a role, changing the velocity distribution and possibly also the density, therefore the velocity distribution of the elements has to be taken with care. However the distribution in velocitites might leave observable features in the spectra that could be used for diagnostic purposes. Finally, by projecting our 3D results on previously obtained 1D ones, the observed evolution of the Si, S, and Ca lines might give us the possibility to understand whether normal SN Ia are well mixed deflagrations.

5 Summary and conclusions

In this paper we presented the results of nucleosynthesis calculations obtained coupling a tracer particle method to two-dimensional and three-dimensional Eulerian hydrodynamic calculations of SN Ia explosion. The multidimensional SN Ia simulations described in this work employed a minimal nuclear reaction network, sufficient for a good approximation of the thermonuclear energy release. Although the predicted chemical composition agrees well with the expectations, we presented here the results of a very detailed study of the nuclear abundances in the ejecta obtained by post-processing the output of the hydrodynamic models. Since the multidimensional hydrodynamical scheme applied is of Eulerian type (i.e. the grid does not move with the fluid), we added a Lagrangian component to the calculations in the form of tracer particles. In order to record temperature and density evolution as a function of time (necessary input for the nucleosynthesis calculations) we homogeneously distributed $\sim$20 000 marker particles and followed their T and $\rho$ evolution. We then calculated the nucleosynthesis experienced by each marker and computed the total yield as a sum over all the markers including the decays of unstable isotopes.

The nuclear reaction network employed in computing the explosive nucleosynthesis contains 383 nuclear species ranging from neutrons, protons, and $\alpha$-particles to molibdenum. For this work, the initial mixture we used consists of (mass fraction) 0.475 $M_{\odot }$ of ¹²C, 0.5 $M_{\odot }$ of ¹⁶O, and 0.025 $M_{\odot }$ of ²²Ne. When the flame passes through the fuel, C, O and Ne are converted to heavier elements, with different compositions depending on the T and $\rho$ history. Nuclear statistical equilibrium conditions are assumed in the marker particles with $T > 6\times 10^{9}$ K. At such temperatures ( $T > 6\times 10^{9}$ K) a mixture of ⁵⁶Ni and $\alpha$-particles in NSE is synthesized. Below that temperature burning only produces intermediate mass elements. Once the temperature drops $T < 1.5 \times 10^9$ K, no burning takes place during the short timescale ($\simeq$1.5 s) of the explosion ("unburned'' material).

The current research focused on the sensitivity of the explosion on the ignition conditions and on the detailed nucleosynthetic yields that they predict. We could demonstrate that multi-dimensional explosion models allow us to predict their nucleosynthesis yields with some confidence. It was shown that only 3D models are potentially able to produce enough ⁵⁶Ni to explain the light curves of "normal'' type Ia supernovae, and that also the ignition conditions (central ignition vs. several off-center ignition spot) affect the nucleosynthesis yields. Since the number of ignition spots we can put into the numerical models depend on the spatial resolution and since the explosion energy as well as the Ni-mass increase with increasing number of spots, we expect that our best resolved b30_3d_768 model is closest to what a "typical'' pure-deflagration supernova might eject. The general nucleosynthesis outcome of SNeIa is dominated by Fe-group elements, involving also sizable fractions of Si-Ca and minor amounts of unburned (C and O) or pure C-burning products (e.g. Na, Ne, Mg). Despite of the fact that differences with respect to the standard W7 nucleosynthesis (Iwamoto et al. 1999; Brachwitz et al. 2000; Thielemann et al. 2003) are found, in particular in the ⁵⁶Ni mass produced, as well as in the final amount of unburned material, in general the nuclear yields are consistent with expectations. We can therefore say with some confidence and without parameterization, that the Chandrasekhar mass scenario with a pure turbulent deflagration is a viable candidate for SN Ia explosions. We also note that the significant amount of unburned material ejected by our SN Ia models may have an interesting impact on the role of SN Ia in the context of Galactic chemical evolution of C (giving a contribution of the order of $\sim$20% to the total C at the solar composition). In the case of O still the main sources are massive stars.

Comparing the nucleosynthesis presented in this paper to observed SN Ia spectra, the reader should keep in mind that our models do not reach the homologous expansion phase. We are currently working to modify the combustion hydrocode, using a moving grid that will allow us to follow the evolution much longer. The results will be published elsewhere. We are also performing a detailed parameter study of the variation of the central density and of the initial carbon/oxygen ratio of the SN Ia progenitor (Röpke et al., in preparation). Finally, recent calculations by Timmes et al. (2003) indicate large variations of the ⁵⁶Ni mass as a function of metrallicity (measured by the original ²²Ne content). An investigation of the metallicity effect on the nucleosynthesis and yields is also in progress (Travaglio et al., in preparation).

	Figure 11: Mass fractions of selected isotopes as a function of the radial velocity of the markers (taken at $\sim$1.2 s) for the b30_3d_768 model ( upper panel) and for the c3_3d_256 ( lower panel). The width of each velocity bin is 300 km s-1. For each isotope we sum its abundance over all markers in a certain velocity bin.
Open with DEXTER

Acknowledgements

C.T. thanks the Alexander von Humboldt Foundation, the Federal Ministry of Education and Research, and the Programme for Investment in the Future (ZIP) of the German Governmen, and the Max-Planck Institute für Astrophysik (Garching bei München), for their financial support. W.H. and M.R. acknowledge support by the DFG under grant SFB-375/C5 and by the European Commission under contract HPRN-CT-2002-00303. F.-K.T. is supported by the Swiss National Science Foundation.

References

• Arnett, W. D., & Livne, E. 1994, ApJ, 427, 31 [NASA ADS] (In the text)
• Brachwitz, F., Dean, D. J., Hix, W. R., et al. 2000, ApJ, 536, 934 [NASA ADS] [CrossRef] (In the text)
• Fuller, G. M., Fowler, W. A., & Newman, M. J. 1985, ApJ, 293, 1 [NASA ADS] (In the text)
• Gamezo, V. N., Khokhlov, A. M., Oran, E. S., Chtchelkanova, A. Y., & Rosenberg, R. O. 2003, Science, 299, p.77 (In the text)
• Garcia-Senz, D., & Woosley, S. E. 1995, ApJ, 454, 895 [NASA ADS] [CrossRef] (In the text)
• Hillebrandt, W., & Niemeyer, J. C. 2000, ARA&A, 38, 191 [NASA ADS] (In the text)
• Iwamoto, K., Brachwitz, F., Nomoto, K. I., et al. 1999, ApJS, 125, 439 [NASA ADS] [CrossRef] (In the text)
• Khokhlov, A. M. 1995, ApJ, 449, 695 [NASA ADS] [CrossRef] (In the text)
• Khokhlov, A. M., Höflich, P. A., Oran, E. S., et al. 1999, ApJ, 524, L107 [NASA ADS] [CrossRef] (In the text)
• Langanke, K., & Martinez-Pinedo, G. 2000, Nucl. Phys. A, 673, 481 [NASA ADS] [CrossRef] (In the text)
• Leibundgut, B. 2001, ARA&A, 39, 67 [NASA ADS] (In the text)
• Lisewski, A. M., Hillebrandt, W., Woosley, S. E., Niemeyer, J. C., & Kerstein, A. R. 2000, ApJ, 537, L405 [NASA ADS] [CrossRef] (In the text)
• Livne, E. 1993, ApJ, 406, L17 [NASA ADS] [CrossRef] (In the text)
• Maeda, K., Nakamura, T., Nomoto, K., et al. 2002, ApJ, 565, 405 [NASA ADS] [CrossRef] (In the text)
• Martinez-Pinedo, G., Langanke, K., & Dean, D. J. 2000, ApJS, 126, 493 [NASA ADS] [CrossRef] (In the text)
• Nagataki, S., Hashimoto, M. A., Sato, K., & Yamada, S. 1997, ApJ, 486, 1026 [NASA ADS] [CrossRef] (In the text)
• Niemeyer, J. C. 1999, ApJ, 523, L57 [NASA ADS] [CrossRef] (In the text)
• Niemeyer, J. C., & Hillebrandt, W. 1995a, ApJ, 452, 769 [NASA ADS] [CrossRef] (In the text)
• Niemeyer, J. C., & Hillebrandt, W. 1995b, ApJ, 452, 779 [NASA ADS] [CrossRef]
• Niemeyer, J. C., Reinecke, M., Travaglio, C., & Hillebrandt, W. 2003, From Twilight to Highlight: The Physics of Supernovae, Proc. ESO/MPA/MPE Workshop, Garching, Germany, 151 (In the text)
• Nomoto, K., Thielemann, F. K., & Yokoi, K. 1984, ApJ, 286, 644 [NASA ADS] [CrossRef] (In the text)
• Rauscher, T., Heger, A., Hoffman, R. D., & Woosley, S. E. 2002, ApJ, 576, 323 [NASA ADS] [CrossRef] (In the text)
• Reinecke, M., Hillebrandt, W., & Niemeyer, J. C. 1999, A&A, 347, 739 [NASA ADS] (In the text)
• Reinecke, M., Hillebrandt, W., & Niemeyer, J. C. 2002a, A&A, 386, 936 [EDP Sciences] [NASA ADS] [CrossRef] (In the text)
• Reinecke, M., Hillebrandt, W., & Niemeyer, J. C. 2002b, A&A, 391, 1167 [EDP Sciences] [NASA ADS] [CrossRef]
• Röpke, F. K., Niemeyer, J. C., & Hillebrandt, W. 2003, ApJ, 588, 952 [NASA ADS] [CrossRef] (In the text)
• Thielemann, F.-K., Nomoto, K., & Yokoi, K. 1986, A&A, 158, 17 [NASA ADS] (In the text)
• Thielemann, F.-K., Nomoto, K., & Hashimoto, M. 1996, ApJ, 460, 408 [NASA ADS] [CrossRef] (In the text)
• Thielemann, F.-K., Argast, D., Brachwitz, F., et al. 2003, From Twilight to Highlight: The Physics of Supernovae, Proc. ESO/MPA/MPE Workshop, Garching, Germany, 331 (In the text)
• Timmes, F. X., Brown, E. F., & Truran, J. W. 2003, ApJ, 590, L83 [NASA ADS] [CrossRef] (In the text)
• Travaglio, C., Kifonindis, K., & Müller, E. 2004, New Astron. Rev., 48, 25 [NASA ADS] [CrossRef] (In the text)
• Whelan, J., & Iben, I. 1973, ApJ, 186, 1007 [NASA ADS] [CrossRef]
• Woosley, S. E., Wunsch, S., & Kuhlen, M. 2004, ApJ, 607, 921 [NASA ADS] [CrossRef] (In the text)

 

  
Online Material

Table 1: Synthesized Mass ($M_{\odot }$) in SN Ia models.

Table 2: Synthesized Mass ($M_{\odot }$) for radioactive species in SN Ia models.