EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 642 (October 2020) A&A, 642 (2020) A135 Full HTML
Free Access
Issue
A&A
Volume 642, October 2020
Article Number A135
Number of page(s) 15
Section Astrophysical processes
DOI https://doi.org/10.1051/0004-6361/202038116
Published online 13 October 2020

© ESO 2020

1. Introduction

Molecules in supernovae (SNe) must form in hostile conditions. The initial supernova explosion releases a large amount of energy (∼1051 erg for a typical type II-P SN), most of which is converted into kinetic energy that ejects stellar gas at velocities up to a few percent of the speed of light. The decay of radioactive species produces a population of high-energy Compton electrons that ionize and dissociate molecules. A strong radiation field in the hot ejecta provides another destruction mechanism.

Despite these destruction mechanisms, there have been molecules detected as early as ∼100 days after the initial explosion. The first molecule identified in any supernova is carbon monoxide (CO) in the nearby SN 1987A at 136 days after the explosion (Catchpole et al. 1988; Spyromilio et al. 1988), followed by silicon monoxide (SiO) from day 192 (Aitken et al. 1988; Meikle et al. 1989; Roche et al. 1991). The possible detection of molecules other than CO and SiO at these epochs were also reported, such as CO+ (Meikle et al. 1989; Spyromilio et al. 1988) and CS (Meikle et al. 1989). However, these detections typically relied on a single spectral feature and, therefore, they are more uncertain. Contemporary observations of SN 1987A have further identified H2 (Fransson et al. 2016), HCO+, and SO (Matsuura et al. 2017), in addition to a continued growth of the CO mass over decades (Matsuura et al. 2017).

Subsequent observations have confirmed that molecular formation similar to that observed in SN 1987A is likely to be commonplace in SNe. Both CO and SiO have been identified in several more core-collapse SNe of Type IIP/IIL (Spyromilio & Leibundgut 1996; Kotak et al. 2005, 2006, 2009; Pozzo et al. 2006; Yuan et al. 2016; Banerjee et al. 2018; Tinyanont et al. 2019), Type IIb/Ib (Ergon et al. 2015; Drout et al. 2016), and Type Ic (Gerardy et al. 2002; Hunter et al. 2009). Despite the harsh environment, there seems to be substantial evidence to support a relatively rich molecular chemistry in SNe. The first models for the creation and destruction of CO and SiO predicted a condensation fraction of the atomic species into molecules on the order of order 10−3 (Petuchowski et al. 1989; Lepp et al. 1990), which is in rough agreement with observations.

The study of molecules in supernovae is important for several reasons. From molecular rovibrational emission in the infrared, we can infer physical conditions and supernova parameters, such as temperature and density (e.g., Liu & Dalgarno 1995). Recently, the Atacama Large Millimeter/submillimeter Array (ALMA) has allowed for this kind of analysis to extend into the radio regime and for the nearby SN 1987A spatially resolved maps of the molecules to allow new kinds of constraints to be set (Abellán et al. 2017; Cigan et al. 2019). Forthcoming observations with the James Webb Space Telescope (JWST) of SNe in the infrared (IR), where the molecular contribution is significant, will also provide a crucial new observational viewpoint.

Beyond their observational signatures, molecules are important for the physical conditions in SNe, which, in turn, govern the atomic emission. The emission from CO in the IR may dominate the cooling of the oxygen and carbon-rich gas after around a hundred days, which then affects the ionization state (see e.g., Liu & Dalgarno 1995). Thus, including the effects of molecules on the temperature and ionization is paramount for understanding and modeling the thermal evolution of SNe, which is, in turn, vital for modeling the optical and near infrared (NIR) emission by atoms and ions.

The status of supernova optical/NIR spectral modeling at nebular phases is now at a level where attempts to determine detailed nucleosynthesis yields and supernova progenitor masses are being made (e.g., Dessart et al. 2013; Jerkstrand et al. 2014). This has led to a renewed impetus in terms of the red supergiant problem (Smartt et al. 2009; Smartt 2015) and in gaining an understanding of the landscape of successful and failed explosions as the stellar cores collapse (O’Connor & Ott 2011; Ertl et al. 2016). In the majority of such efforts, either the molecules are ignored or treated in a parameterized fashion (e.g., Jerkstrand et al. 2012, 2014).

Research that has focused on the molecular formation chemistry (e.g., Petuchowski et al. 1989; Lepp et al. 1990; Liu et al. 1992; Liu & Dalgarno 1995, 1996; Gearhart et al. 1999) has not considered, on the other hand, the impact of the molecules on the atomic and ionic emission or the spectral synthesis. Several works also concentrate on molecule formation as a pathway to dust condensation (e.g., Clayton et al. 1999, 2001; Cherchneff & Dwek 2009, 2010; Biscaro & Cherchneff 2014, 2016; Sluder et al. 2018) given the fact that molecules can be predecessors to dust. The dust production that is yielded from different types of SNe at different metallicities is still under debate, partly because it is unknown to what extent dust survives shocks produced in later circumstellar interaction stages. It is unclear how the SN dust yields compare to other major dust production sites, such as the winds of asymptotic giant branch stars and dust growth in the ISM, especially at early times (see e.g., Zhukovska et al. 2008; Valiante et al. 2009; Dwek & Cherchneff 2011; Triani et al. 2020). The main goal of the aforementioned works investigating theories of dust production in SNe is typically to investigate dust formation and provide dust yields for different types of SNe, without necessarily addressing the influence of either molecules or dust on spectra.

The work we initiate here aims to join these two strands (so far) of parallel research and to produce, in a self-consistent way, predictions for the near-UV to far-infrared (UVOIR) spectra of supernovae that take molecules into consideration. In this first paper, we describe the methodology used for implementing molecular formation and cooling in the SN spectral synthesis code SUMO, we produce test calculations, which provide estimates of molecular masses, and perform sensitivity tests. While several different molecular species are present in this first test model, we focus our analysis specifically on CO, as it is known from observations to be abundant in SNe and acts an important coolant, possibly cooling the carbon and oxygen-rich areas of the ejecta by several thousand degrees. Before embarking on realistic multi-zone models in later studies, here we study the most important reactions, any sensitivity to uncertain chemical reaction rates, and the kind of ejecta properties that could be most promisingly diagnosed by observed molecular emission.

This paper is organized as follows. Section 2 presents a description of the spectral modeling methods and the implementation of molecular formation. A first test model of the CO production in SN 1987A is presented in Sect. 3. The sensitivity of this model to input parameters is investigated in Sect. 4 and the possibility of using inferred molecular masses to constrain the density and deposition energy is explored in Sect. 5. In Sect. 6, we provide our summary and conclusions.

2. Modeling methods

2.1. Supernova model

The SUMO code (Jerkstrand et al. 2011, 2012) takes a hydrodynamic model for supernova ejecta (density and composition as a function of velocity) as input and calculates the physical conditions and emergent spectra by solving the temperature equation, the statistical equilibrium equations for both ion abundances and level populations, and the radiation field. It is tailored for the post-peak phases and gives snapshot steady-state solutions at any specified epoch. It currently includes 22 elements between hydrogen and nickel, with about 8000 levels, and 300 000 lines.

2.2. Chemical model

For modeling molecule formation, we incorporate a chemical kinetic description into the SUMO code. For the chemical network we use Ns species and Nr reactions, with species referring to all included neutral atoms, atomic ions, neutral molecules, molecular ions, and free thermal electrons. As chemical reactions occur, number fractions either decrease (if reactions destroy a species) or increase (if reactions create a species). The rate of change of the number density [X] (cm−3) of a species X is solved for in steady state as

d ( [ X ] ) d t = ∑ C i − ∑ D i = 0 , $$ \begin{aligned} \frac{\mathrm{d}(\mathrm{[X]})}{\mathrm{d}t} = \sum \mathcal{C} _i - \sum \mathcal{D} _i = 0 ,\end{aligned} $$(1)

where 𝒞i and 𝒟i are the rates of processes that create and destroy the species, respectively.

The set of Eq. (1), one for each included species, form a system of coupled, non-linear, algebraic equations, which is solved with a Newton-Raphson scheme. At each epoch, this is done several times in an outer loop where temperature, radiation field, and ionization balance are iterated over. The assumption of a steady state is expected to be a good approximation for the range of epochs investigated in this work (150–600 d). This is discussed in Cherchneff & Dwek (2009) and examined in detail in Sect. 4 in Gearhart et al. (1999), where relevant reaction timescales are found to be shorter than the dynamic and radioactive timescales. The validation of the steady-state approximation is also revealed from the results of our calculations in Sect. 3.1 It should, however, be noted that steady-state does not necessarily hold at later epochs, as noted in Cherchneff & Dwek (2009).

Reactions in a chemical network can be categorized according to how many reactants are involved; unimolecular for one reactant, bimolecular for two reactants, and termolecular for three reactants, and so on. In the following examples, the quantities Q1 and Q2 refer to the reacting species, and P1 and P2 to the products of a reaction in the chemical network, which could be either molecules, atoms or ions depending on the specific reaction. For a unimolecular reaction (u),

Q 1 → P 1 + P 2 , $$ \begin{aligned} \mathrm{Q}_1 \rightarrow \mathrm{P}_1 + \mathrm{P}_2 ,\end{aligned} $$(2)

where the reactant Q1 turns into the products P1 and P2, for example, as a result of thermal decomposition. The rates of change of each involved species due to this reaction are:

R u = − ( d ( [ Q 1 ] ) d t ) u = ( d ( [ P 1 ] ) d t ) u = ( d ( [ P 2 ] ) d t ) u = k u [ Q 1 ] , $$ \begin{aligned} R_u = -\left( \frac{\mathrm{d}([\mathrm{Q}_1])}{\mathrm{d}t} \right)_u =\left(\frac{\mathrm{d}(\mathrm{[P}_1])}{\mathrm{d}t} \right)_u= \left(\frac{\mathrm{d}(\mathrm{[P}_2])}{\mathrm{d}t} \right)_u= k_u \mathrm{[Q}_1] ,\end{aligned} $$(3)

with k being the temperature-dependent rate coefficient, which is specific for each reaction and denoted by a suitable subscript.

Similarly for a bimolecular reaction (b), for example, a collision between two species, Q1 and Q2, creating products, P1 and P2, as

Q 1 + Q 2 → P 1 + P 2 , $$ \begin{aligned} \mathrm{Q}_1 + \mathrm{Q}_2 \rightarrow \mathrm{P}_1 + \mathrm{P}_2 ,\end{aligned} $$(4)

the rates of change are

R b = − ( d ( [ Q 1 ] ) d t ) b = − ( d ( [ Q 2 ] ) d t ) b = ( d ( [ P 1 ] ) d t ) b = ( d ( [ P 2 ] ) d t ) b = k b [ Q 1 ] [ Q 2 ] . $$ \begin{aligned}&R_b = -\left( \frac{\mathrm{d}([\mathrm{Q}_1])}{\mathrm{d}t} \right)_b =-\left( \frac{\mathrm{d}([\mathrm{Q}_2])}{\mathrm{d}t} \right)_b \nonumber \\&\quad =\left(\frac{\mathrm{d}(\mathrm{[P}_1])}{\mathrm{d}t} \right)_b = \left(\frac{\mathrm{d}(\mathrm{[P}_2])}{\mathrm{d}t} \right)_b = k_b \mathrm{[Q}_1] \mathrm{[Q}_2]. \end{aligned} $$(5)

Reactions with more reactants, as with termolecular reactions such as a three-body association, are disregarded here as such reactions only become significant at epochs with higher densities than those examined in this work.

The units of the rate coefficients depend on the number of reactants: s−1 and cm3 s−1 for unimolecular and bimolecular reactions, respectively. The reaction rates at a temperature, T, are typically expressed in an Arrhenius-type form (McElroy et al. 2013) as:

k ( T ) = α × ( T 300 K ) β × exp ( − γ / T ) , $$ \begin{aligned} k(T) = \alpha \times \left( \frac{T}{300\,\text{ K}} \right)^{\beta } \times \exp ({-\gamma / T}) ,\end{aligned} $$(6)

where α, β, and γ are parameters that are specific to each reaction.

The various reaction types included and the form of the different rate equations are discussed below and summarized in Table 1. The primary source of rate coefficients used here is the UMIST Database for Astrochemistry (McElroy et al. 20131). This is supplemented by reaction rate coefficients from Sluder et al. (2018) and Cherchneff & Dwek (2009), where reaction networks in supernova remnants were investigated in a way that is similar to the present work. See Appendix B for a complete reference list.

Table 1.

Reaction types included which pertain to molecules.

2.3. Included reaction types

In the following section, we describe the reaction types used in this work. Throughout this section, A, B, C, and D are used to represent (neutral) atoms, while AB and BC represent molecules made up by either A and B or B and C. While the examples here only include diatomic molecules, the notation scheme can be straightforwardly generalized to triatomic or larger molecules.

2.3.1. Thermal collision reactions

Collisions between molecules and atoms, or between atoms and atoms, can change both the composition of the gas by forming or destroying molecules and change the ionization state by charge transfer.

Molecules may form directly by radiative association, where two species, A and B, combine to an initially energized complex. As the gas densities generally are too low for collisional stabilization, where a collision with a third species carries away energy before the energized complex can dissociate, the main channel is stabilization by radiation (Smith 2011). Then A and B form the new species, AB, and emit a photon as

A + B → AB + h ν $$ \begin{aligned} \mathrm{A} + \mathrm{B} \rightarrow \mathrm{AB} + \mathrm{h}\nu \end{aligned} $$(7)

for two neutral reactants A and B, or

A + + B → AB + + h ν , $$ \begin{aligned} \mathrm{A}^+ + \mathrm{B} \rightarrow \mathrm{AB}^+ + \mathrm{h}\nu ,\end{aligned} $$(8)

in the case of an ion and a neutral.

Species can also rearrange through ion-neutral or neutral-neutral interactions. During ion-neutral interactions, the neutral species is polarized by the electric field of the ion, which induces an electric dipole moment and causes an attractive force between the ion and the neutral. The outcome of such interaction depends on what is energetically favorable for the involved species. One such possible reaction is an ion-neutral exchange of an atomic or molecular ion colliding with a neutral species exchanging one of the components in the molecule as

AB + + C → A + BC + . $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{C} \rightarrow \mathrm{A} + \mathrm{BC}^+ .\end{aligned} $$(9)

A second possibility is ion-neutral charge exchange, which involves ions colliding with a neutral species transferring charge as

A + + B → A + B + . $$ \begin{aligned} \mathrm{A}^+ + \mathrm{B} \rightarrow \mathrm{A} + \mathrm{B}^+.\end{aligned} $$(10)

Similarly, there can be charge exchange between a molecule and an atom as

AB + + C → AB + C + $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{C} \rightarrow \mathrm{AB} + \mathrm{C}^+ \end{aligned} $$(11)

or between two molecules. If the neutral species does not possess a dipole moment, the reaction rates of ion-neutral interactions can be assumed to be temperature independent, however, if the neutral has a dipole moment, the average attraction increases for lower temperatures (Larsson et al. 2012). Consequently, these reactions are very important for low-temperature environments, such as ISM and late-stage supernovae.

Similarly, in a neutral-neutral exchange, two neutral species collide and react as

AB + C → A + BC . $$ \begin{aligned} \mathrm{AB} + \mathrm{C} \rightarrow \mathrm{A} + \mathrm{BC} .\end{aligned} $$(12)

Such interactions are typically weakly attractive at larger distances due to van der Waal forces and repulsive at small distances. As a result, there is a higher energy barrier to the reactions, with reaction rates increasing with temperature (Smith 2011). The temperatures in supernovae are, however, sufficiently high and these types of reactions are included in the network.

Due to low density and comparably low temperatures, we disregard three-body reactions and thermal fragmentation reactions. The CO thermal fragmentation rate (using rate coefficients from Appleton et al. 1970) at the temperature and densities here is always at least five orders of magnitude smaller than other important destruction processes. We allow for up to two reactants and three products in the reactions. All collision reactions are then treated as bimolecular reactions, described in Eqs. (4) and (5), and the rates used for the included collisions can be seen in Table B.1.

2.3.2. Recombination

Recombination with electrons can recombine an atomic ion as well as recombine or dissociate a molecular ion (or both). Radiative recombination between an atomic ion and a thermal electron results in neutralization as

A + + e − → A + h ν . $$ \begin{aligned} \mathrm{A}^+ + \mathrm{e}^- \rightarrow \mathrm{A} + \mathrm{h}\nu .\end{aligned} $$(13)

For the case of molecules more possible outcomes of recombination are feasible, such as

AB + + e − → AB + h ν . $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{e}^-&\rightarrow \mathrm{AB} + \mathrm{h}\nu .\end{aligned} $$(14)

AB + + e − → A + B . $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{e}^-&\rightarrow \mathrm{A} + \mathrm{B}.\end{aligned} $$(15)

AB + + e − → A + + B + e − . $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{e}^-&\rightarrow \mathrm{A}^+ + \mathrm{B} + \mathrm{e}^- .\end{aligned} $$(16)

AB + + e − → A + B + + e − . $$ \begin{aligned} \mathrm{AB}^+ + \mathrm{e}^-&\rightarrow \mathrm{A} + \mathrm{B}^+ + \mathrm{e}^- . \end{aligned} $$(17)

Here, the reaction in Eq. (14) is analogous to the atomic case in Eq. (13). For molecules, the dissociative recombination pathway in Eq. (15) will typically dominate by a large factor (Larsson et al. 2012). The rates for the molecular dissociative recombination reactions (Eq. (15)) used in this work are listed in Table B.2. We also assume that radiative recombination reactions of the type shown in Eq. (14) occurs with a small fraction (0.1%) of the corresponding dissociative recombination reaction rate and that dissociative recombination reactions resulting in ions (Eqs. (16) and (17)) are negligible.

2.3.3. Destruction by Compton electrons

Supernovae provide one quite particular component to the physical environment; a population of high-energy Compton electrons (which we label e C − $ \rm{e}_\mathrm{{C}}^{-} $ here) that give a significant and often dominant molecular destruction rate. These Compton electrons originate from γ-photons produced in radioactive decay that repeatedly Compton scatter on bound and free thermal electrons, giving them high energies. These high-energy electrons can ionize atomic species, and ionize or dissociate molecules. For a diatomic molecule the following destruction processes are considered:

AB + e C − → AB + + e − + e C − . $$ \begin{aligned} \mathrm{AB} + \mathrm{e}^-_{\rm C}&\rightarrow \mathrm{AB}^+ + \mathrm{e}^- + \mathrm{e}^-_{\rm C}. \end{aligned} $$(18)

AB + e C − → A + B + e C − . $$ \begin{aligned} \mathrm{AB} + \mathrm{e}^-_{\rm C}&\rightarrow \mathrm{A} + \mathrm{B} + \mathrm{e}^-_{\rm C}. \end{aligned} $$(19)

AB + e C − → A + + B + e − + e C − . $$ \begin{aligned} \mathrm{AB} + \mathrm{e}^-_{\rm C}&\rightarrow \mathrm{A}^+ + \mathrm{B} + \mathrm{e}^- + \mathrm{e}^-_{\rm C}. \end{aligned} $$(20)

AB + e C − → A + B + + e − + e C − . $$ \begin{aligned} \mathrm{AB} + \mathrm{e}^-_{\rm C}&\rightarrow \mathrm{A} + \mathrm{B}^+ + \mathrm{e}^- + \mathrm{e}^-_{\rm C}. \end{aligned} $$(21)

The rate coefficients for these reactions can be calculated by assuming energy from Compton electron collisions is deposited in the zone at a rate of L [erg s−1]. For Ntot total number of molecules or atoms, the energy deposition rate per particle is then L/Ntot. The rate coefficient for a specific Compton destruction reaction can then be defined as

k C = L N tot W i , $$ \begin{aligned} k_{C} = \frac{L}{N_{\rm tot} W_i} ,\end{aligned} $$(22)

where Wi is the mean energy per ion pair for a given species. It generally depends on the composition and ionization state. The rate of change per volume for a reaction involving species AB, for example AB+ eC − → AB + + e − + eC − $ {\rm AB} + {\rm e}^-_{\rm C} \rightarrow {\rm AB}^+ + {\rm e}^- + {\rm e}^-_{\rm C} $, then RC = kC[AB], effectively treating it as a unimolecular reaction, as described in Eqs. (2) and (3).

The Compton electron destruction rates for CO were investigated in Liu & Dalgarno (1995), presenting a relationship between the electron fraction and Wi, for a fixed composition (49.5% O I, 49.5% C I, 1% CO I). Wi settles to a constant value at low values of electron fraction, xe ≲ 0.01. A constant Wi is often assumed (e.g., Cherchneff & Dwek 2009; Sluder et al. 2018), but instead, we fit a cubic polynomial to the data presented in Fig. 2 in Liu & Dalgarno (1995), therefore allowing Wi to vary with the fractional ionization. For details about the cubic polynomial fits, see Appendix A.

We note that SUMO contains a solver for the Compton electron distribution and does calculate the non-thermal destruction rates for atoms and ions. However, adding molecules to this module is a significant task that lies outside the scope of this paper. Here, we aim to describe the chemical reaction network and CO non-local thermodynamic equilibrium (NLTE) cooling in a test model, so we do not rely on a more detailed treatment of this aspect.

In the case of CO, which is the primary subject of this work’s investigation, the branching fractions leading to dissociation into an atomic ion and electron (the reactions in Eqs. (20) and (21)) are much smaller than the first two reaction types, and are therefore disregarded. For other molecules with no available data for Wi, we assume the same value as for CO.

2.4. Reverse reactions

Between specific internal states, i and j, any reaction has a reverse reaction whose rate obeys the detailed balance: Rreverse = Rforwardgi/gfexp(−ΔE/kT), where g are the statistical weights and ΔE the energy difference. The exponential term is for the assumption of a thermal collision partner. In practice, it is often only the total forward rate, to a sum of target states, that is known. This prevents a direct calculation of state-specific reverse rates. Because of this, in this work, we do not attempt to estimate unknown reverse rates for any forward rate. When both forward and reverse rates are presented in the literature, they typically do not obey a detailed balance relation for the same reason (reactants are in the ground state whereas products are in multiple excited states).

While it may not serve as a very accurate description (per the discussion above), we assume, for practical purposes within the network, that products are formed in their ground states. For molecules that are not explicitly treated as multi-level NLTE elements, this is, in fact, the only state. We further assume, as per the standard in SUMO, that only elements in their ground multiplet state contribute to a chemical reaction (but this fraction is typically close to unity).

2.5. Photoionization

No work that we are aware of to date has calculated self-consistent photoionization rates for molecules in supernova ejecta. While this is something we hope to do in subsequent work, we do not attempt it here (whereas photoionization for atoms and ions is considered as is customary in SUMO). We note that thanks to the use of approximative rates, it is estimated that the UV destruction of molecules is a sub-dominant process (Cherchneff & Dwek 2009).

2.6. CO line emission cooling

The formation of CO is expected to affect the temperature evolution of the SN ejecta, as CO emission in the infrared adds a new efficient cooling channel. A previous work investigating CO cooling in SNe has indicated that the effect is substantial; Liu & Dalgarno (1995) found that CO emission could lower the temperature by several thousand degrees in the O/C zone.

The NLTE cooling by CO ro-vibrational transitions was added into SUMO, adopting energy levels and transition probabilities from Li et al. (2015), including levels up to the vibrational quantum number ν = 6 and rotational quantum number J = 322. Tests showed that this was sufficient, whereas the inclusion of more states did not significantly influence the cooling from CO. Ro-vibrational radiative transitions for the fundamental, first, and second overtones were included. Pure rotational transitions were ignored, as they are not important for the cooling due to the low transition probabilities and small energies involved. Optical depth for the molecular line emission is considered, as it is for the atomic lines, by the Sobolev formalism.

Collision rates were treated with the standard default approximation in SUMO, which uses (for transitions with no measured or calculated values) a collision strength of Υ = 0.004g1g2 for forbidden transitions (and van Regemorter’s formula for allowed). Recombination was ignored in calculating the level populations.

3. Test case: CO formation in SN 1987A

To validate the modeling methods and to perform sensitivity tests, we explored CO formation in a simplified one-zone model of O/C-composition for SN 1987A, devoid of any elements other than oxygen and carbon. Since most CO likely forms in the O/C-zone (Liu & Dalgarno 1995), our calculated CO masses for this test model should be comparable to observationally inferred yields from SN 1987A. The model is based on the O/C zone of a MZAMS = 19 M⊙ core-collapse explosion model of SN 1987A, originally derived from models presented in Woosley & Heger (2007) and previously investigated using SUMO in Jerkstrand et al. (2011, 2012).

This zone has a mass of 0.58 M⊙ and mass fractions of 0.62 O, 0.28 C, 0.078 He, 0.018 Ne, and 4.9 × 10−3 Mg, plus traces of other elements. We note that the He abundance is typically negligible in the C/O zone, thus, the large value here of 8% is an outlier among the 12, 15, 19, and 25 M⊙ models. As we only include C and O here, we normalize their mass fractions to 0.69 and 0.31, respectively. The zone was given a spherical shape, with Vin = 0 km s−1 and Vout = 540 km s−1. This gives it the same density as in Jerkstrand et al. (2012); in that more realistic multi-zone model, the material resides as multiple clumps over a larger velocity span. While ignoring the other SN zones would change the radiation field in the zone, the (ionizing) UV radiation is mostly emitted locally. In addition, as mentioned previously, in this first test model we treat photoionization only for atoms and ions but not for molecules. We adopt the gamma deposition calculation from Jerkstrand et al. (2012) (see Appendix C here). The species included in this step are listed in Table 2 and the reaction network used can be seen in Appendix B. We model the supernova between 150 and 600 days, with time intervals of 50 days.

Table 2.

Neutral atoms, atomic ions, neutral molecules, and molecular ions included in the simulations.

3.1. Results

The resulting temperature and ionization fractions of C and O in the model are shown in Fig. 1. Over the time we investigated, the temperature drops from around 5000 K down to 800 K. The temperature drops significantly when the amount of CO starts to become significant, which is at around 200d (Fig. 2). There is initially around 10% singly ionized carbon and around 1% singly ionized O in these models, which is decreasing with time. We note that the kinks seen at 200 days, in the time evolution of O+ and O2+, are due to a temperature-dependent reaction involving O+ (Eq. (29)) that turns off at lower temperatures (≲2000 K), thus changing the trends seen (see also discussion in Sect. 3.2). The fractions of doubly ionized C and O are orders of magnitude smaller and these play no role in the physical conditions. For CO, the degree of ionization is much smaller than for the atoms, with CO+/CO ≲10−4.

thumbnail Fig. 1.

Black line shows the temperature (left axis) of the SN 1987A model and the colored lines show the ionization degrees (right axis) of carbon (teal) and oxygen (blue).

thumbnail Fig. 2.

Lines represent molecular masses from our model over time and the scatter points are observational estimates of the CO masses. The observational estimates using NLTE models with optical depth consideration (white squares, Liu et al. 1992; and grey diamonds, Liu & Dalgarno 1995) generally give values that are close to the model prediction.

In Fig. 2, the molecular masses against time are shown, with the most abundant molecule being CO. At the first epoch of 150d only a relatively small amount of ∼10−4 M⊙ forms. There is then a rapid increase up to about 250 d, where the CO mass levels out at ∼4 × 10−3 M⊙. Other molecular species only form in small amounts compared to CO, with masses several orders of magnitude lower. Species not plotted only formed in minuscule amounts.

In Fig. 2 also observational estimates of CO masses for SN 1987A are plotted (Spyromilio et al. 1988; Liu et al. 1992; Liu & Dalgarno 1995). These estimates are derived from the same infrared observations of CO first overtone bands, however, several different assumptions are used and the reported CO mass varies by 1–2 orders of magnitude between the different works, indicating significant model dependency. In the initial paper by Spyromilio et al. (1988) LTE and optically thin conditions were assumed, and they derived a CO mass of ∼10−5 − 10−4 M⊙, with a tendency of growth from 200 d to 350 d. Later Liu et al. (1992) found a larger CO mass (∼10−3 M⊙) when considering NLTE and optical depth effects; both effects were found to give a too low CO mass if ignored. In their NLTE model no time-dependence of the CO mass was seen, or a slight decrease with time. Liu & Dalgarno (1995) took a somewhat different approach, using simulated model temperatures instead of having temperature as a free fitting parameter as in the previous works. With this method, they found a slightly higher mass than in Liu et al. (1992) (5 × 10−3 M⊙), with significantly cooler temperatures than from the free-fitting approach. New data allowed the range of epochs to be extended to 550 d; also, they found no time evolution in the inferred CO mass. The papers where optical depth is taken into account (Liu et al. 1992; Liu & Dalgarno 1995) generally produce better spectral fits to the data.

Keeping these differences in mind, as they result in a significant spread in the observational estimates of CO masses, there is a general agreement between the observations and our model results. Between the optically thin LTE results and our models, there are between one and two orders of magnitude of difference, however, for the NLTE models, the difference is less than 0.5 orders of magnitude. Our results are notably similar to the estimates in Liu & Dalgarno (1995) after 250 days.

With our standard SN 1987A model, we come to the same conclusion as Gearhart et al. (1999) and Cherchneff & Dwek (2009) concerning the validity of using a steady-state solution for a small network. If we define a characteristic timescale of a process as the ratio of the quantity divided by the rate of change of that quantity, τ = ϕ / ϕ ˙ $ \tau = \phi / \dot{\phi} $, we can estimate the dynamical timescale as τdyn ≈ n/ṅ ≈ 1/3 × t and the radioactivity timescale as the decay time of 56Co, i.e. τrad = τ56Co = 111days. At a few hundred days post explosion, meaning that at times typical to nebular phase supernovae, these two timescales are comparable (e.g., at 300 days, τdyn ≈ τrad ≈ 100days). In contrast, the timescales of the important chemical reactions τchem tend to be significantly shorter. If we take τchem as [X]/∑𝒟i, where ∑𝒟i is the sum over all the destruction rates of X, we obtain typical values of τchem ≈ 10−2days at 200 days for CO, monotonically increasing to τchem ≈ 4days at 600 days. Because τchem ≪ τdyn, τrad the steady state approximation is valid.

3.2. The principal paths for CO

The panels of Fig. 3 show the reaction rates (i.e., either Ru or Rb depending on the reaction, as described in Eqs. (3) and (5)) for the most important processes divided by the total number density at that time-step. This gives rates per particle (unit s−1) which removes the t−3 background evolution of rates per unit volume due to the homologous expansion. We note that destruction rates are not per CO particle but per any particle; to get the result per CO particle, we divide by the number fraction of CO. The processes creating CO are shown in the left panel and processes that destroy CO are shown in the right panel. As seen in the left panel of Fig. 3, important creation pathways for CO at all times are the radiative association between C and O as

C + O → CO + h ν , $$ \begin{aligned} \mathrm{C} + \mathrm{O} \rightarrow \mathrm{CO} + \mathrm{h}\nu , \end{aligned} $$(23)

thumbnail Fig. 3.

Reactions important for the formation and destruction of CO over time for the SN 1987A model. Left: rates of different reactions creating CO divided by the total number density. Right: destruction rates of CO divided by the total number density.

and neutral exchange with C2

C 2 + O → C + CO . $$ \begin{aligned} \mathrm{C}_2 + \mathrm{O} \rightarrow \mathrm{C} + \mathrm{CO.} \end{aligned} $$(24)

The C2 here forms mainly through radiative association as

C + C → C 2 + h ν . $$ \begin{aligned} \mathrm{C} + \mathrm{C} \rightarrow \mathrm{C}_2 + \mathrm{h}\nu . \end{aligned} $$(25)

At early times, up to around 200 days, reactions involving CO+ have among the highest flows of all CO formation channels. Important reactions are the charge transfer between CO+ and O,

CO + + O → CO + O + , $$ \begin{aligned} \mathrm{CO}^+ + \mathrm{O} \rightarrow \mathrm{CO} + \mathrm{O}^+, \end{aligned} $$(26)

the analogous charge transfer between CO+ and C

CO + + C → CO + C + $$ \begin{aligned} \mathrm{CO}^+ + \mathrm{C} \rightarrow \mathrm{CO} + \mathrm{C}^+ \end{aligned} $$(27)

and the recombination reaction

CO + + e − → CO . $$ \begin{aligned} \mathrm{CO}^+ + \mathrm{e}^- \rightarrow \mathrm{CO.} \end{aligned} $$(28)

At these times, CO+ is in turn primarily created by radiative association between C and O+, and by charge exchange between CO and O+. These channels are, however, effectively quenched at later times as the latter reaction

CO + O + → CO + + O $$ \begin{aligned} \mathrm{CO} + \mathrm{O}^+ \rightarrow \mathrm{CO}^+ + \mathrm{O} \end{aligned} $$(29)

is endothermic with a few tenths of an eV and turns off when the gas gets too cool (the rate is derived and discussed in Petuchowski et al. 1989 based on cross section measurements from Murad 1973). This has two consequences; as seen in Fig. 2, the amount of CO+ decreases significantly between 200 and 250 days, as one of the major creation channels is turned off. The CO creation pathways involving CO+ are suppressed as a result, as seen in the left panel of Fig. 3. Secondly, the dominating destruction channel for CO at early epochs turns off. The combined effect is a rapid increase in CO mass. While there are two other CO creation pathways of a similar efficiency as the ones involving CO+, the rate in Eq. (29) strongly dominates the other destruction processes. When it turns off, the destruction rate, therefore, is reduced by a larger factor than the creation rate, and the CO mass increases. The evolution of CO mass is consequently intimately related to both the formation and destruction of other molecules and to the temperature in a complex way. We note that, in general, the full impact of any atom or molecule in the reaction network for the CO abundance (or any other abundance) cannot be fully established from plots like these. It would be necessary to run the network with that particle extracted. However, if a rate is large, the reaction is potentially important and should be included.

The two dominant destruction pathways from 250 days onwards are the destruction by Compton electrons as

CO + e C − → CO + or C + O $$ \begin{aligned} \mathrm{CO} + \mathrm{e}^-_{\rm C} \rightarrow \mathrm{CO}^+\ \mathrm{or}\ \mathrm{C} + \mathrm{O} \end{aligned} $$(30)

and the neutral exchange with O as

CO + O → C + O 2 . $$ \begin{aligned} \mathrm{CO} + \mathrm{O} \rightarrow \mathrm{C} + \mathrm{O}_2. \end{aligned} $$(31)

These results broadly agree with the findings of previous works concerning CO formation in supernovae. There is a general consensus that the radiative association between C and O is the most important creation pathway for CO in this type of supernova (see e.g., Lepp et al. 1990; Liu et al. 1992; Gearhart et al. 1999; Sarangi & Cherchneff 2013). The importance of other formation processes discussed here also tend to be mentioned. It should be noted, however, that these results are dependent on assumptions about the physical conditions as well as the rate coefficients used in the network. The self-consistent temperature calculation shown here represents one major improvement in reducing free parameters in the modeling.

3.3. Thermal evolution

There is a complex interplay between the temperature and the CO mass. The CO effectively cools the gas which, in turn, affects some of the temperature-sensitive chemical reaction rates that govern the creation and destruction of CO. As mentioned in the previous section, the relationship between CO mass and temperature is non-trivial, however, a cooler environment generally leads to more CO being created. The main destruction channel at higher temperatures, described in Eq. (29), is very temperature-dependent and quickly becomes inefficient when the temperature drops. Once CO starts to form, there is positive feedback; CO then helps to cool the gas by its rovibrational emission and with a cooler temperature, more CO will form. This sustains until critical temperature-dependent reactions, such as charge transfer destruction with O+, turn off.

This can be seen in Fig. 4, where the temperature of our SN 1987A O/C-zone model has been calculated with and without CO cooling included. The effect is substantial; when CO cooling is included the temperature is several thousand degrees lower and the amount of CO formed about an order of magnitude larger. This indicates that even a small amount of CO has large consequences for the thermal evolution. These results broadly agree with previous findings by Liu & Dalgarno (1995).

thumbnail Fig. 4.

Upper: temperature evolution of the SN 1987A model, with (black) and without (pink) CO cooling. Lower: mass of CO formed over time for the SN 1987A model, with (black) and without (pink) CO cooling.

3.4. Comparison to other works

Several theoretical investigations into the CO production in SNe have been carried out since the first detection in SN 1987A. Table D.1 in Appendix D shows an overview of the CO mass reported in selected works that specifically model CO masses which are comparable to those in SN 1987A, much as in this paper, at 300 days.

There is a large spread in results; early works typically found a CO mass of ∼10−5 − 10−6 M⊙ (e.g., Petuchowski et al. 1989; Lepp et al. 1990). Modeling in Liu & Dalgarno (1995), where the cooling effects of CO were taken into account, gave a higher CO mass of ∼10−3 M⊙, while more recent works (e.g., Sarangi & Cherchneff 2013; Sluder et al. 2018) find a CO mass of ∼10−2 M⊙. Our results (4 × 10−3 M⊙ at 300 days) generally agree well with Liu & Dalgarno (1995) and are within an order of magnitude to the findings in Sarangi & Cherchneff (2013) and Sluder et al. (2018).

Differences in networks and reaction rates used and the inclusion of processes such as photoionization, cooling by CO, and Compton destruction in some works are possible explanations for the differences. There are also other possible reasons for the disparity, such as different modeling parameters, for example, the mass and density of the O/C zone, the fractions of O and C, etc. Table D.1 also compares some key modeling parameters and, as we can see, there are significant differences between the mentioned quantities; there is up to a factor five difference between the C and O masses, an order of magnitude difference in number densities, and up to 2500 K difference in temperatures (mostly parameterized rather than calculated).

Earlier works (Petuchowski et al. 1989; Lepp et al. 1990; Liu et al. 1992) typically also assumed a fully mixed model, with several solar masses of helium (He), making the CO destruction by collision with He+ very efficient. For in-depth discussions on the differences between earlier theoretical works, see Gearhart et al. (1999) and Sluder et al. (2018).

4. Sensitivity tests

Here, we investigate how sensitive the resulting molecular masses are to modeling assumptions and uncertainties, which roughly can be divided into two categories; assumptions about the physical properties of the supernova (e.g., deposition energy, density, composition,…) and uncertainties of the chemical network, that is, the rate coefficients.

4.1. Varying the physical parameters

We tested the sensitivity of CO formation on the input parameters density and deposition energy, with Fig. 5 showing the time evolution of the CO mass for different densities in the left panel, and different deposition energies in the right panel.

thumbnail Fig. 5.

Results of varying the physical parameters of the SN1987A model. Scatter points are observational estimates. Left: lines show CO mass of models with original density, twice and half the original density. Right: lines show CO mass of models with original deposition energy, twice and half the original deposition energy.

For different densities, the models show a complex behavior between 150 and 250 days. A rapid onset of CO formation, seen for the original model and the low-density model between 150 and 200 days, starts later in the high-density model. This is due to the previously discussed delicate relationship between the temperature and the CO destruction path through charge transfer with O+ (Eq. (29)). The dense model has a higher temperature and it takes a longer time for this model to cool to temperatures where the aforementioned charge transfer reaction turns off.

After 250 days, there is a positive correlation between density and CO mass, which is more or less uniform over time. In a symmetric network (with the same number of reactants of all processes), the increase in number density should not significantly change the mass of the formed molecules. This is, however, not the case for CO. As discussed in Sect. 3.2, several bimolecular processes contribute significantly to the creation of CO, while destruction is dominated by only one process, namely, destruction by Compton electrons, which is treated as a unimolecular process. While unimolecular processes roughly depend linearly on the number density, bimolecular processes have a square relationship. Therefore, if number density is increased, the sum of the creation processes will increase by a larger factor than the sum of the destructive processes in the case of CO, resulting in a higher CO mass. Related to this point is the fact that the number density of Compton electrons will stay roughly constant upon compression.

In contrast to density, increased deposition energy decreases the CO mass, as seen in the right plot of Fig. 5. When increasing the deposition energy in the model, the temperature and the ionization fractions increase. Initially, as the main destruction channel for CO depends on O+ (Eq. (29)), a larger abundance of O+, as well as a higher temperature, will increase the rate at which CO depletes. The main creation channels for CO, on the other hand, do not depend on any atomic ions, so increases in these will have little impact on the total creation rate. At later times the destruction by Compton electrons takes over as the dominant destruction channel. The Compton electron population increases with higher deposition energy and the outcome of increasing deposition energy is, therefore, a more efficient destruction of CO by Compton electrons.

It should be noted that both the positive correlation between model density and CO mass, and the negative correlation between deposition energy and CO mass, are not general results for all molecules, but rather a consequence of the specific creation and destruction channels for CO. It is most likely that such relationships are unique for each molecule and need to be investigated independently for each molecule of interest.

4.2. Impact of uncertain rate coefficients

The rate coefficients used for describing the thermal collisions involving CO come from multiple sources and are typically derived using different methods (experiments, calculations, or estimates) with a range of uncertainties. We want to investigate how changing the thermal collision rate coefficients by some typical uncertainty factor affects the resulting CO masses, for example, to test the degree of robustness with respect to uncertain reaction rates. This is done using a Monte Carlo approach.

4.2.1. Monte Carlo simulations

We adopt the following method to assess the impact of changing the thermal collision rates. Firstly, each thermal collision rate in the network (as seen in Table B.1) is given a log-normal probability distribution, with a mean equal to the literature value used in the standard network and with a standard deviation σrate (units of dex). Each rate coefficient is then randomly sampled from this distribution to make a randomly sampled network. The CO mass is solved for using a full SUMO run at 300 days using the standard SN 1987A single-zone model for this randomly sampled network. These first steps are then are repeated 500 times, re-sampling the thermal collision rates at each instance to investigate the parameter space. This investigation is are repeated for five different σrate = 0.1, 0.3, 0.5, 0.7, 0.9 dex.

From these Monte Carlo simulations, we get a CO mass distribution for each σrate, which will indicate how the uncertainties of individual rates propagate to the molecular mass uncertainty.

4.2.2. Results of MC simulations

The left panel of Fig. 6 shows examples of the resulting CO mass distributions, for Monte Carlo runs with three different σrate. As seen, the distributions of CO mass can be fit well with a Gaussian in log-normal space. A larger σrate leads to a larger variation in CO masses, as expected. Furthermore, the mean of the CO mass distribution shifts with σrate; with a larger σrate the mean becomes larger.

thumbnail Fig. 6.

Left: CO mass distribution, for the MC tests with three different σrate. Right: σCO (standard deviation of the resulting CO distribution) against the σrate (standard deviation of the individual rates), for details see text. The solid black line is a fit to the data (slope = 0.8), the dashed line is the one-to-one relationship. For more details, see the text.

As the resulting distribution of CO masses can be fitted with a Gaussian, we may compare the standard deviation of the CO mass distribution σCO with the standard deviation of the individual rate coefficients σrate, seen in the right panel of Fig. 6. For the values, σrate = 0.1−0.9, there is a linear relationship between σrate and σCO with a slope slightly smaller than 1. This indicates that there is a roughly linear relationship between the uncertainties of the thermal collision rate coefficients and the calculated mass of CO. Consequently, if the typical rate coefficient is known to a factor of two, then the CO mass uncertainty will be approximately a factor of two.

5. Constraining physical parameters from molecule abundances

If we can assume that the chemistry modeling parameters (i.e., abundance of C and O, the mass of the O/C zone, and the rate coefficients) are correct for a given SN within reasonable uncertainties, it is possible to use the modeled CO mass to constrain the deposition energy and the density by making a comparison with observational estimates of the CO mass. To investigate the feasibility of this method, we made a grid of models spanning a range of deposition energies and densities, with values ranging between 0.1 × ρ0 and 10 × ρ0 for density, and 0.1 × L0 and 10 × L0 for deposition energy (ρ0 = 4.88 × 109 cm−3, and L0 = 6.2 × 1039 erg s−1 are the density and deposition energy for the standard SN 1987A model, described in Sect. 3.1, at 284 days. The epoch of 284 days is chosen to coincide with the observations. We produced 100 models, equidistant in log-space for the mentioned densities and deposition energies. For other input parameters, the same values were used as for the standard model (described in Sect. 3).

The resulting CO mass for the model grid can be seen in Fig. 7. There is a clear inverse relationship between the deposition energy and CO mass, as discussed in Sect. 4.1. The density has a more complex relationship at these times, as the O+ charge transfer reaction (Eq. (29)) is still efficient for some densities. As seen in Fig. 5, 284 days is close to the transition point to the regime where increased density leads to higher CO mass, which is also the general trend seen in Fig. 7.

thumbnail Fig. 7.

Variation of CO mass with deposition energy and density of the standard SN 1987A model at 284d. Each dot represents a model.

The panels in Fig. 8 show the norm between the log of CO masses from the model grid and from different observational estimates, (i.e., ϵ = |log(MCO(model)) − log(MCO(observation))|). Consequently, the 0.3 contour delineates where the difference between the model and observed CO masses is 0.3 dex. As shown in Fig. 7, bands through the (L, ρ) plane produce the same CO mass. A density specification would, in general, allow a unique deposition energy to be picked out. However, a deposition energy specification can in some parts of the plane give degenerate density solutions. The more complex dependency on density is discussed in Sect. 4.1. One particular conclusion can be drawn: if we compare to the NLTE CO mass estimates from Liu et al. (1992) and Liu & Dalgarno (1995), which generally produce the best fits to the spectral shapes, an upper limit on the deposition energy in the O/C zone of SN 1987A can be set of ∼2L0 = 1.2 × 1040 erg −1, at 284 days. Such constraints can be used to test modern 3D hydrodynamic models where the energy deposition in any given zone depends on the morphology and mixing, which, in turn, depends on the properties of the explosion and the progenitor.

thumbnail Fig. 8.

Differences between various observational estimates and the model result, over a grid of models with varying deposition energies and densities. The upper two plots compare to LTE CO mass estimates from Spyromilio et al. (1988) and Liu et al. (1992). The lower two plots are comparisons with NLTE CO estimates from Liu et al. (1992) and Liu & Dalgarno (1995). The dark grey shaded areas show the contours of 0.3 dex difference in CO mass and the light grey shaded areas show the 0.6 dex contours. ρ0 = 4.88 × 109 cm−3, and L0 = 6.2 × 1039 erg s−1 are the density and deposition energy for the standard SN 1987A model.

When more molecule mass estimates become available and can be modeled, as the correlation between deposition energy, density, and molecular mass depends on the unique creation and destruction pathways of each molecule, it could be possible to break degeneracies and get very specific constraints from observationally inferred molecule masses.

6. Summary and conclusions

In this work, we implement molecule formation and NLTE cooling into the spectral synthesis code SUMO and explore molecular formation physics as well as sensitivities to uncertain reaction rates. We study CO formation in a simplified model of the O/C zone of SN 1987A for a time interval of 150–600 days after the explosion. To summarize the results:

  • The thermal evolution is closely connected to the CO formation. The CO emission in the infrared adds a new significant cooling channel and after some time cooling by CO dominates in this environment. We find a positive feedback between CO formation and cooling. Consequently, self-consistent inclusion of CO cooling is important when simulating CO production.

  • The CO mass uncertainty scales approximately linearly with the uncertainties of the thermal collision rates. If the rates are known to a factor of two, the CO mass uncertainty will be a factor two or slightly less (our formal relation between CO mass change and reaction rate uncertainty is 0.8).

  • The physical parameters are important for CO production; density has a positive correlation with CO mass while deposition energy has a negative correlation. This could be used to constrain the physical parameters of the supernova if accurate CO mass observational estimates are available. Using this approach we can put an upper limit on the deposition energy at 284 days in the O/C zone of SN 1987A to ∼2 × L0, where L0 = 6.2 × 1039 erg s−1 is the deposition energy of our standard model of SN 1987A.

  • Our test model reproduces the observationally inferred CO mass in SN 1987A to a satisfactory degree. We have especially good agreement with the Liu & Dalgarno (1995) estimate, for which NLTE and optical depth effects are taken into account (for comparison to other works see also Table D.1). The CO mass settles at ∼4 × 10−3 M⊙ at 250d and stays there until at least 600 d. This corresponds to a condensation efficiency of about 1% of the carbon and oxygen. The CO lowers the temperature in the O/C zone by several 1000 degrees already from the beginning of the nebular phase; this means that the carbon and oxygen in this zone will not contribute significantly to optical emission.

It should be noted that the accuracy of the modeled CO mass hinges on both completeness of the network, meaning that most of the important and relevant reactions are accounted for, and on the accuracy of the used rate coefficients. At times before 250 days, a single reaction, namely, the charge transfer between CO and O+ (Eq. (29)), dominates the destruction rate. If this reaction is disregarded, or if there are reactions of similar importance missing, the CO mass could potentially be quite different in that phase. Similarly, the reported rate coefficients that are crucial for determining the molecular masses can vary greatly between different works and different methods. One example of this is the important radiative association reaction between C and O (C + O → CO + hν, Eq. (23)). In this work, the value from Singh et al. (1999) is used (RC37 in Table B.1), which at 1500 K has a value of 1.77  ×  10−17 cm3 s−1. This is almost a factor of two less than the value reported by Gustafsson & Nyman (2015) (2.84  ×  10−17 cm3 s−1). Other works report a range of values, for example, 1.08  ×  10−17 cm3 s−1 (by Dalgarno et al. 1990, fit to data made by Cherchneff & Dwek 2009), and 2.01  ×  10−17 cm3 s−1 by Franz et al. (2011). This is one of the most well-studied reaction rates; the situation for other reactions is more dire and the uncertainties of the reaction rates are an important source of uncertainty for the calculated molecular masses. Continued research into the reaction rates and different processes that are relevant to molecular formation in these environments is therefore of great importance for understanding supernova chemistry.

There are several aspects of the molecule treatment in SUMO that could be improved and made more realistic, with the end goal of obtaining more realistic synthetic SN spectra. A few planned improvements are mentioned here. While initially stated to only be a minor contribution to the destruction of CO, the photoionization of molecules should be implemented into the code. As SUMO models the full radiation field of the ejecta, we can calculate the photoionization rates of molecules consistently, in contrast to previous works where these rates are estimated or inferred. Related to this is implementing self-consistent destruction rates of molecules by Compton electrons, which also potentially could be directly calculated in the code, instead of using previous estimates. Finally, to produce realistic spectra, full multi-zone models are needed. As these models will contain many more species, the chemical network used needs to be expanded to include the formation of more types of molecules. A possible time dependency should also be investigated for modeling at yet later phases.

2

This data is available through the online material of Li et al. (2015).

Acknowledgments

The authors acknowledge funding from the Swedish National Space Board (SNSB/Rymdstyrelsen), Grant 2017-R 95/17.

References

  1. Abellán, F. J., Indebetouw, R., Marcaide, J. M., et al. 2017, ApJ, 842, L24 [NASA ADS] [CrossRef] [Google Scholar]
  2. Adams, N. G., Smith, D., & Grief, D. 1978, Int. J. Mass Spectrom. Ion Process., 26, 405 [NASA ADS] [CrossRef] [Google Scholar]
  3. Aitken, D. K., Smith, C. H., James, S. D., et al. 1988, MNRAS, 235, 19P [CrossRef] [Google Scholar]
  4. Alge, E., Adams, N. G., & Smith, D. 1983, J. Phys. B At. Mol. Phys., 16, 1433 [NASA ADS] [CrossRef] [Google Scholar]
  5. Andreazza, C. M., & Singh, P. D. 1997, MNRAS, 287, 287 [NASA ADS] [CrossRef] [Google Scholar]
  6. Appleton, J. P., Steinberg, M., & Liquornik, D. J. 1970, J. Chem. Phys., 52, 2205 [CrossRef] [Google Scholar]
  7. Banerjee, D. P. K., Joshi, V., Evans, A., et al. 2018, MNRAS, 481, 806 [NASA ADS] [CrossRef] [Google Scholar]
  8. Biscaro, C., & Cherchneff, I. 2014, A&A, 564, A25 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Biscaro, C., & Cherchneff, I. 2016, A&A, 589, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Brian, J., & Mitchell, A. 1990, Phys. Rep., 186, 215 [NASA ADS] [CrossRef] [Google Scholar]
  11. Catchpole, R. M., Whitelock, P. A., Feast, M. W., et al. 1988, MNRAS, 231, 75 [NASA ADS] [CrossRef] [Google Scholar]
  12. Cherchneff, I., & Dwek, E. 2009, ApJ, 703, 642 [NASA ADS] [CrossRef] [Google Scholar]
  13. Cherchneff, I., & Dwek, E. 2010, ApJ, 713, 1 [NASA ADS] [CrossRef] [Google Scholar]
  14. Cigan, P., Matsuura, M., Gomez, H. L., et al. 2019, ApJ, 886, 51 [NASA ADS] [CrossRef] [Google Scholar]
  15. Clayton, D. D., Liu, W., & Dalgarno, A. 1999, Science, 283, 1290 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  16. Clayton, D. D., Deneault, E. A.-N., & Meyer, B. S. 2001, ApJ, 562, 480 [NASA ADS] [CrossRef] [Google Scholar]
  17. Copp, N. W., Hamdan, M., Jones, J. D. C., Birkinshaw, K., & Twiddy, N. D. 1982, Chem. Phys. Lett., 88, 508 [CrossRef] [Google Scholar]
  18. Dalgarno, A., Du, M. L., & You, J. H. 1990, ApJ, 349, 675 [NASA ADS] [CrossRef] [Google Scholar]
  19. Dessart, L., Hillier, D. J., Waldman, R., & Livne, E. 2013, MNRAS, 433, 1745 [NASA ADS] [CrossRef] [Google Scholar]
  20. Drout, M. R., Milisavljevic, D., Parrent, J., et al. 2016, ApJ, 821, 57 [NASA ADS] [CrossRef] [Google Scholar]
  21. Dwek, E., & Cherchneff, I. 2011, ApJ, 727, 63 [NASA ADS] [CrossRef] [Google Scholar]
  22. Ergon, M., Jerkstrand, A., Sollerman, J., et al. 2015, A&A, 580, A142 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  23. Ertl, T., Janka, H. T., Woosley, S. E., Sukhbold, T., & Ugliano, M. 2016, ApJ, 818, 124 [NASA ADS] [CrossRef] [Google Scholar]
  24. Fahey, D. W., Fehsenfeld, F. C., & Ferguson, E. E. 1981, Geophys. Res. Lett., 8, 1115 [CrossRef] [Google Scholar]
  25. Fehsenfeld, F. C., & Ferguson, E. E. 1972, J. Chem. Phys., 56, 3066 [CrossRef] [Google Scholar]
  26. Fransson, C., Larsson, J., Spyromilio, J., et al. 2016, ApJ, 821, L5 [CrossRef] [Google Scholar]
  27. Franz, J., Gustafsson, M., & Nyman, G. 2011, MNRAS, 414, 3547 [CrossRef] [Google Scholar]
  28. Gearhart, R. A., Wheeler, J. C., & Swartz, D. A. 1999, ApJ, 510, 944 [NASA ADS] [CrossRef] [Google Scholar]
  29. Gerardy, C. L., Fesen, R. A., Nomoto, K., et al. 2002, PASJ, 54, 905 [NASA ADS] [CrossRef] [Google Scholar]
  30. Gustafsson, M., & Nyman, G. 2015, MNRAS, 448, 2562 [CrossRef] [Google Scholar]
  31. Hunter, D. J., Valenti, S., Kotak, R., et al. 2009, A&A, 508, 371 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Jerkstrand, A., Fransson, C., & Kozma, C. 2011, A&A, 530, A45 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Jerkstrand, A., Fransson, C., Maguire, K., et al. 2012, A&A, 546, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  34. Jerkstrand, A., Smartt, S. J., Fraser, M., et al. 2014, MNRAS, 439, 3694 [NASA ADS] [CrossRef] [Google Scholar]
  35. Kotak, R., Meikle, P., van Dyk, S. D., Höflich, P. A., & Mattila, S. 2005, ApJ, 628, L123 [NASA ADS] [CrossRef] [Google Scholar]
  36. Kotak, R., Meikle, P., Pozzo, M., et al. 2006, ApJ, 651, L117 [NASA ADS] [CrossRef] [Google Scholar]
  37. Kotak, R., Meikle, W. P. S., Farrah, D., et al. 2009, ApJ, 704, 306 [NASA ADS] [CrossRef] [Google Scholar]
  38. Larsson, M., Geppert, W. D., & Nyman, G. 2012, Rep. Prog. Phys., 75, 066901 [NASA ADS] [CrossRef] [Google Scholar]
  39. Lepp, S., Dalgarno, A., & McCray, R. 1990, ApJ, 358, 262 [NASA ADS] [CrossRef] [Google Scholar]
  40. Li, G., Gordon, I. E., Rothman, L. S., et al. 2015, ApJS, 216, 15 [NASA ADS] [CrossRef] [Google Scholar]
  41. Liu, W., & Dalgarno, A. 1995, ApJ, 454, 472 [NASA ADS] [CrossRef] [Google Scholar]
  42. Liu, W., & Dalgarno, A. 1996, ApJ, 471, 480 [NASA ADS] [CrossRef] [Google Scholar]
  43. Liu, W., Dalgarno, A., & Lepp, S. 1992, ApJ, 396, 679 [NASA ADS] [CrossRef] [Google Scholar]
  44. Matsuura, M., Indebetouw, R., Woosley, S., et al. 2017, MNRAS, 469, 3347 [NASA ADS] [CrossRef] [Google Scholar]
  45. McElroy, D., Walsh, C., Markwick, A. J., et al. 2013, A&A, 550, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  46. Meikle, W. P. S., Allen, D. A., Spyromilio, J., & Varani, G. F. 1989, MNRAS, 238, 193 [NASA ADS] [CrossRef] [Google Scholar]
  47. Mitchell, G. F. 1984, ApJS, 54, 81 [NASA ADS] [CrossRef] [Google Scholar]
  48. Murad, E. 1973, J. Chem. Phys., 58, 4374 [CrossRef] [Google Scholar]
  49. Nomoto, K., Shigeyama, T., Kumagai, S., & Yamaoka, H. 1991, Supernovae, 176 [CrossRef] [Google Scholar]
  50. O’Connor, E., & Ott, C. D. 2011, ApJ, 730, 70 [NASA ADS] [CrossRef] [Google Scholar]
  51. Oscar Martinez, J., Betts, N. B., Villano, S. M., et al. 2008, ApJ, 686, 1486 [NASA ADS] [CrossRef] [Google Scholar]
  52. Petuchowski, S. J., Dwek, E., Allen, J. E., Jr, & Nuth, J. A., III 1989, ApJ, 342, 406 [NASA ADS] [CrossRef] [Google Scholar]
  53. Pinto, P. A., & Woosley, S. E. 1988, Nature, 333, 534 [NASA ADS] [CrossRef] [Google Scholar]
  54. Pozzo, M., Meikle, W. P. S., Rayner, J. T., et al. 2006, MNRAS, 368, 1169 [NASA ADS] [CrossRef] [Google Scholar]
  55. Prasad, S. S., & Huntress, W. T., Jr 1980, ApJS, 43, 1 [NASA ADS] [CrossRef] [Google Scholar]
  56. Rauscher, T., Heger, A., Hoffman, R. D., & Woosley, S. E. 2002, ApJ, 576, 323 [NASA ADS] [CrossRef] [Google Scholar]
  57. Roche, P. F., Aitken, D. K., & Smith, C. H. 1991, MNRAS, 252, 39P [CrossRef] [Google Scholar]
  58. Sarangi, A., & Cherchneff, I. 2013, ApJ, 776, 107 [NASA ADS] [CrossRef] [Google Scholar]
  59. Singh, P. D., Sanzovo, G. C., Borin, A. C., & Ornellas, F. R. 1999, MNRAS, 303, 235 [NASA ADS] [CrossRef] [Google Scholar]
  60. Sluder, A., Milosavljevic, M., & Montgomery, M. H. 2018, MNRAS, 480, 5580 [Google Scholar]
  61. Smartt, S. J. 2015, PASA, 32, e016 [NASA ADS] [CrossRef] [Google Scholar]
  62. Smartt, S. J., Eldridge, J. J., Crockett, R. M., & Maund, J. R. 2009, MNRAS, 395, 1409 [NASA ADS] [CrossRef] [Google Scholar]
  63. Smith, I. W. M. 2011, ARA&A, 49, 29 [NASA ADS] [CrossRef] [Google Scholar]
  64. Smith, I. W. M., Herbst, E., & Chang, Q. 2004, MNRAS, 350, 323 [NASA ADS] [CrossRef] [Google Scholar]
  65. Spyromilio, J., & Leibundgut, B. 1996, MNRAS, 283, L89 [NASA ADS] [CrossRef] [Google Scholar]
  66. Spyromilio, J., Meikle, W. P. S., Learner, R. C. M., & Allen, D. A. 1988, Nature, 334, 327 [NASA ADS] [CrossRef] [Google Scholar]
  67. Tinyanont, S., Kasliwal, M. M., Krafton, K., et al. 2019, ApJ, 873, 127 [CrossRef] [Google Scholar]
  68. Triani, D. P., Sinha, M., Croton, D. J., Pacifici, C., & Dwek, E. 2020, MNRAS, 493, 2490 [CrossRef] [Google Scholar]
  69. Valiante, R., Schneider, R., Bianchi, S., & Andersen, A. C. 2009, MNRAS, 397, 1661 [NASA ADS] [CrossRef] [Google Scholar]
  70. Viggiano, A. A., Albritton, D. L., Fehsenfeld, F. C., et al. 1980, ApJ, 236, 492 [NASA ADS] [CrossRef] [Google Scholar]
  71. Wakelam, V., Herbst, E., Loison, J.-C., et al. 2012, ApJS, 199, 21 [NASA ADS] [CrossRef] [Google Scholar]
  72. Woosley, S. E. 1988, ApJ, 330, 218 [NASA ADS] [CrossRef] [Google Scholar]
  73. Woosley, S. E., & Heger, A. 2007, Phys. Rep., 442, 269 [NASA ADS] [CrossRef] [Google Scholar]
  74. Yuan, F., Jerkstrand, A., Valenti, S., et al. 2016, MNRAS, 461, 2003 [NASA ADS] [CrossRef] [Google Scholar]
  75. Zhukovska, S., Gail, H.-P., & Trieloff, M. 2008, A&A, 479, 453 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

Appendix A: Wi vs fractional ionization interpolation

Figure 2 in Liu & Dalgarno (1995) shows Wi, the mean energy per ionization or dissociation, as a function of the fractional ionization xe, for different products of the reaction between CO and Compton electrons e C − $ ^{-}_{\mathrm{C}} $. A cubic polynomial is fitted to this data in log-log space, to describe the variation of Wi against xe. The fit takes the form

log ( W i ) = a × log ( x e ) 3 + b × log ( x e ) 2 + c × log ( x e ) + d $$ \begin{aligned} \log (W_i) = a \times \log (x_{\rm e})^3 + b \times \log (x_{\rm e})^2 + c \times \log (x_{\rm e}) + d \end{aligned} $$(A.1)

and can be seen in Fig. A.1. The coefficients a, b, c, and d for each fit is shown in Table A.1.

thumbnail Fig. A.1.

Original data and fits to the results on mean energy per ion pair for CO ionization presented in Liu & Dalgarno (1995). The original data is shown as solid lines and the fit to each relation is shown as dashed lines.

Table A.1.

Coefficients a, b, c, and d, used in Eq. (A.1), for different products from the reaction CO + e C − $ ^-_{\rm C} $.

Appendix B: Model reaction rates

The reaction rates used in the SN 1987A test model. As this model only contains carbon and oxygen, and ions of carbon and oxygen; all included molecules are composites of these.

Table B.1.

Overview of the reaction coefficients for thermal collision reactions.

Table B.2.

Overview of the coefficients for the molecular recombination reactions.

Appendix C: Deposition energy

Figure C.1 shows the deposition energy used for the standard SN 1987A model, which is adopted from Jerkstrand et al. (2012).

thumbnail Fig. C.1.

Deposition energy used for the standard SN 1987A model, taken from the calculations of Jerkstrand et al. (2012).

Appendix D: Overview of other work

Overview of the physical parameters in models for CO formation SN 1987A, from different works.

Table D.1.

Comparison of key properties of different chemical models describing CO formation in SN 1987A.

All Tables

Table 1.

Reaction types included which pertain to molecules.

In the text
Table 2.

Neutral atoms, atomic ions, neutral molecules, and molecular ions included in the simulations.

In the text
Table A.1.

Coefficients a, b, c, and d, used in Eq. (A.1), for different products from the reaction CO + e C − $ ^-_{\rm C} $.

In the text
Table B.1.

Overview of the reaction coefficients for thermal collision reactions.

In the text
Table B.2.

Overview of the coefficients for the molecular recombination reactions.

In the text
Table D.1.

Comparison of key properties of different chemical models describing CO formation in SN 1987A.

In the text

All Figures

thumbnail Fig. 1.

Black line shows the temperature (left axis) of the SN 1987A model and the colored lines show the ionization degrees (right axis) of carbon (teal) and oxygen (blue).
In the text
thumbnail Fig. 2.

Lines represent molecular masses from our model over time and the scatter points are observational estimates of the CO masses. The observational estimates using NLTE models with optical depth consideration (white squares, Liu et al. 1992; and grey diamonds, Liu & Dalgarno 1995) generally give values that are close to the model prediction.
In the text
thumbnail Fig. 3.

Reactions important for the formation and destruction of CO over time for the SN 1987A model. Left: rates of different reactions creating CO divided by the total number density. Right: destruction rates of CO divided by the total number density.
In the text
thumbnail Fig. 4.

Upper: temperature evolution of the SN 1987A model, with (black) and without (pink) CO cooling. Lower: mass of CO formed over time for the SN 1987A model, with (black) and without (pink) CO cooling.
In the text
thumbnail Fig. 5.

Results of varying the physical parameters of the SN1987A model. Scatter points are observational estimates. Left: lines show CO mass of models with original density, twice and half the original density. Right: lines show CO mass of models with original deposition energy, twice and half the original deposition energy.
In the text
thumbnail Fig. 6.

Left: CO mass distribution, for the MC tests with three different σrate. Right: σCO (standard deviation of the resulting CO distribution) against the σrate (standard deviation of the individual rates), for details see text. The solid black line is a fit to the data (slope = 0.8), the dashed line is the one-to-one relationship. For more details, see the text.
In the text
thumbnail Fig. 7.

Variation of CO mass with deposition energy and density of the standard SN 1987A model at 284d. Each dot represents a model.
In the text
thumbnail Fig. 8.

Differences between various observational estimates and the model result, over a grid of models with varying deposition energies and densities. The upper two plots compare to LTE CO mass estimates from Spyromilio et al. (1988) and Liu et al. (1992). The lower two plots are comparisons with NLTE CO estimates from Liu et al. (1992) and Liu & Dalgarno (1995). The dark grey shaded areas show the contours of 0.3 dex difference in CO mass and the light grey shaded areas show the 0.6 dex contours. ρ0 = 4.88 × 109 cm−3, and L0 = 6.2 × 1039 erg s−1 are the density and deposition energy for the standard SN 1987A model.
In the text
thumbnail Fig. A.1.

Original data and fits to the results on mean energy per ion pair for CO ionization presented in Liu & Dalgarno (1995). The original data is shown as solid lines and the fit to each relation is shown as dashed lines.
In the text
thumbnail Fig. C.1.

Deposition energy used for the standard SN 1987A model, taken from the calculations of Jerkstrand et al. (2012).
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.