© The Authors 2025

1. Introduction

Up until recently, infrared spectroscopic observations of core-collapse supernovae (SNe) have largely been limited to the near infrared, and primarily used for assessing the presence of molecular emission from the CO first overtone or the presence of He I lines at 1.083 and 2.058 μm (see, e.g., Spyromilio et al. 1988; Gerardy et al. 2002; Hunter et al. 2009; Rho et al. 2018; Hsiao et al. 2019; Rho et al. 2021, or Shahbandeh et al. 2022). The longer wavelength range was explored for SN 1987A by Wooden et al. (1993), and later opened by Spitzer for a sample of nearby objects, including Type II SN 2004dj (Kotak et al. 2005; Meikle et al. 2011), SN 2005af (Kotak et al. 2006), SN 2004et (Kotak et al. 2009), showing evidence of CO and SiO emission at nebular times, as well as a wealth of unique emission lines, typically without optical counterparts and arising from (stable) Ni, Ar, or Ne. Ongoing observations with the James Webb Space Telescope (JWST) are now bringing a revolution in the field by providing high-quality observations from the near- to the mid-infrared with an exquisite signal even for not so nearby transients. So far though, such data for core-collapse SNe have been mostly photometric and employed for constraining their dust properties (see, e.g., Shahbandeh et al. 2023; Zsíros et al. 2024; Szalai et al. 2025; Sarangi et al. 2025).

There is, however, much interest in using infrared spectra to infer properties of the ejecta, the explosion, and the pre-SN star, thus employing spectral lines as probes of the emitting gas, as done recently for Type Ia SNe 2021aefx and 2022pul (see, e.g., DerKacy et al. 2023; Kwok et al. 2023; Ashall et al. 2024; Kwok et al. 2024). Doing this early in the nebular phase of a Type II SN, that is prior to dust formation and the partial obscuration of the metal-rich inner ejecta (thus before 500 d if we rely on what occurred in SN 1987A; Lucy et al. 1989) and before any late-time power contribution from ejecta-CSM interaction appears (see, e.g., Weil et al. 2020), would allow for the inference of key diagnostics not easily accessible with optical data. Similar information for stripped-envelope SNe is yet to come.

Previous radiative-transfer modeling of the infrared properties of core-collapse SNe seems to go back to SN 1987A. From the modeling of the evolution of the iron, cobalt, and nickel lines during the first two years, Li et al. (1993) inferred that Fe clumps occupy 30% of the volume in the inner ejecta as a result of the ⁵⁶Ni-bubble effect (Woosley 1988; Basko 1994). Kozma & Fransson (1998a,b) modeled the evolution of the ejecta temperature, ionization, and the associated optical and infrared line emission using physical models of the explosion for SN 1987A, covering epochs from 200 d out to several years, and thereby constrained some of the ejecta composition. Other works focused on CO and SiO molecule formation and cooling (Liu et al. 1992; Liu & Dalgarno 1994, 1995). More recently, Jerkstrand et al. (2012) performed radiative-transfer modeling with the Monte Carlo code SUMO (Jerkstrand et al. 2011) of nebular-phase optical and infrared observations of SN 2004et, with a simplified treatment of molecules and dust. Liljegren et al. (2020) and Liljegren et al. (2023) implemented molecular networks within the radiative-transfer codes SUMO, and more recently McLeod et al. (2024) operated a similar upgrade in the grid-based code CMFGEN (Hillier & Dessart 2012), in order to solve for molecular formation and cooling together with the standard processes that control the physics of the gas and the radiation. These latter studies confirm earlier work that molecules such as CO and SiO are efficient coolants and can thus control the temperature of the material where they form (i.e., from the O/C or O/Si rich regions of the ejecta) and affect the local emission arising from atoms and ions, potentially shifting some of the optical flux to the infrared (e.g., through emission in the CO fundamental). The modeling of molecular formation and its impact on the radiative transfer is, however, challenging and computationally expensive.

In the present study, molecules and dust were therefore ignored in order to focus exclusively on the cooling of the gas by atoms and ions. Such simulations are prerequisites for more detailed calculations that would include molecules and dust, a topic that is left to future work. For simplicity, and also because our earlier work showed a relative insensitivity of our model results to clumping in nebular-phase spectra, we adopted a smooth ejecta throughout this work – including clumping can be done in the future starting from the present calculations. We thus merely extended the calculations of red supergiant star explosions presented in Dessart et al. (2021b) and those of He-star explosions presented in Dessart et al. (2021a, 2023b), to document their infrared properties over a broad time span from about 200 to about 500 d. Unlike Liljegren et al. (2023), who focused on molecular diagnostics for stripped-envelope SNe, we explored the rich information contained in the infrared spectra from 1 to 30 μm for a wide range of H-rich and H-poor ejecta spanning a broad range of progenitor masses and nucleosynthetic yields, studying in detail the emission from atoms and ions in the gas. We searched for specific emission line diagnostics associated with O, Na, Mg, Ne, Si, S, Ar, Ca, Fe, Co, and Ni, and investigated what information they can provide about the progenitor and its explosion. We also confronted the information revealed by the infrared and optical ranges, for example, between [Ne II] 12.810 μm and [O I] 0.632 μm, in order to build a panchromatic view of these ejecta models. The present simulations are confronted to the nebular phase optical-to-infrared spectral observations of SN 2024ggi in Dessart et al. (2025).

In the next section, we present the numerical setup of the simulations performed for this work. These are largely inherited from Dessart et al. (2021b) and Dessart et al. (2021a, 2023b), but some specific adjustments were made for this new focus on the infrared. We then describe in detail two ejecta models, one for a Type II SN (model s15p2; Section 3) and another for a Type Ib SN (model he6p0; Section 4). In Section 5, we discuss the possibility of using [Ne II] 12.810 μm as an alternative to [O I] 0.632 μm, which has been extensively used to constrain the preSN or progenitor mass (see, e.g., Fransson & Chevalier 1989; Maguire et al. 2012; Jerkstrand et al. 2012; Dessart et al. 2021b). In Section 6, we describe the evolution of a number of key infrared lines such as [Ni II] 6.634 μm or [Ar II] 6.983 μm. In Section 7, we study the evolution of Fe and Co line fluxes to find evidence of the ⁵⁶Ni-decay chain as the power source in our ejecta models. We present our conclusions in Section 8. In the appendix, we provide additional illustrations and results (Appendix A) as well as a comparison to previous work (Appendix B).

2. Numerical setup

This is a continuation of previous work on nebular-phase radiative-transfer modeling carried out with CMFGEN (Hillier & Dessart 2012). Starting from explosions of solar-metallicity red supergiant stars computed by Sukhbold et al. (2016), Dessart et al. (2021b) presented the properties of the Type II SN models but only at a single epoch of 350 d after explosion. Dessart et al. (2023a) evolved a standard Type II SN model (i.e., 15.2 M_⊙ progenitor reaching core collapse as a red supergiant star and exploding with an explosion energy of about 10⁵¹ erg) from 350 to 1000 d with or without the treatment of shock power from ejecta interaction with circumstellar material. Such models with CSM interaction radiate this additional shock power preferentially in the ultraviolet (or X-rays) at late times and thus with little or negligible power emerging in the infrared (broad, boxy profiles are nonetheless present in the infrared as in the optical, for example, with Paα resembling the Hα profile morphology). In all of these simulations, only the two-step decay chain associated with ⁵⁶Ni was accounted for.

Similarly, starting from explosions of solar-metallicity He-stars computed by Woosley (2019) and Ertl et al. (2020), Dessart et al. (2021a, 2022, 2023b) presented radiative transfer models of Type Ib and Ic SNe covering from about 100 up until about 500 d but with a focus on optical properties only. Dessart (2024) studied the ultraviolet to infrared properties of He-star explosion models until about 10 yr after explosion but under the influence of power injection from the compact remnant.

Thus, in this work, we wish to revisit these simulations and shift our focus to the infrared properties to address forthcoming data from the JWST. Our models cover the nebular phase from about 200 out to about 500 d. We used the same progenitor and explosion models as in these works, and the radiative transfer modeling was performed essentially with the same procedure. Thus, to avoid duplication, we provide here only a succinct summary with information exclusively on the relevant adjustments.

The nomenclature in this study is the same as in previous works, with model s15p2 corresponding to the explosion of a solar-metallicity red supergiant star that started originally on the zero-age main-sequence (ZAMS) with a mass of 15.2 M_⊙. Similarly, model he6p0 corresponds to the explosion of a star that had originally a mass of 6.0 M_⊙ on the helium ZAMS (and a mass of 23.33 M_⊙ on the hydrogen ZAMS; Woosley 2019). To limit the computational burden, we limited our sample to models s10p0, s12p0, s15p2, s18p5, and s21p5 for the H-rich ejecta and to models he6p0, he8p0, and he12p0 for the H-deficient ejecta.

In Dessart et al. (2021b) and Dessart et al. (2021a), the ejecta composition was setup slightly differently. For the Type II SN models, we selected the most abundant stable isotope for each species between H and Ni (for Ne, we used both ²⁰Ne and ²²Ne), together with the undecayed ⁵⁶Ni – most species are usually most abundant in the form of one stable isotope (other forms being underabundant or unstable). We then scaled the mass fractions of all selected isotopes except ⁵⁶Ni to ensure a normalization at each depth to unity. This implies that Ar is ³⁶Ar exclusively and that Ni is either ⁵⁶Ni (unstable) or ⁵⁸Ni (stable). For the He-star explosion models (which were computed at a later date), we used a more accurate procedure by accounting for all isotopes of all important species between He and Ni (e.g., we accounted for all Ni isotopes between ⁵¹Ni and ⁷³Ni when considering the Ni abundance) so that no renormalization was needed. In practice, this produces ejecta models with essentially the same yields as with the previous, alternate and more simplistic treatment (e.g., O or Ne), but there are subtle exceptions for some metals.

For example, Ni is present after explosive nucleosynthesis primarily in three forms, namely ⁵⁶Ni, ⁵⁷Ni, and ⁵⁸Ni. In the ejecta models that we took from Sukhbold et al. (2016) and Ertl et al. (2020), the ⁵⁸Ni yield is comparable to that of ⁵⁷Ni (here, it goes from being comparable at low mass to being 40% greater at a high progenitor mass – see Table 1), and both combined typically represent 10% of the original ⁵⁶Ni mass. However, in Dessart et al. (2021b) and Dessart et al. (2021a), we treated only the radioactive decay of ⁵⁶Ni (i.e., ⁵⁷Ni was either ignored or considered stable). This implies that in the Type II SN models, the Ni abundance was essentially accurate (i.e., we ignored ⁵⁷Ni initially together with its decay), whereas in the He-star explosion models, we typically overestimated the mass of stable Ni by O(10%) (i.e., we accounted for ⁵⁷Ni initially but ignored its decay) and correspondingly underestimated the mass of Co (i.e., what should have been ⁵⁷Co at those epochs)¹. This slight inaccuracy in the abundance of stable Ni is not a major concern and is relevant only for the relative strength of Ni and Co lines in the infrared.

To resolve this issue, we recomputed models s15p2, he6p0, he8p0, and he12p0 with a full, detailed account of the isotopic composition of the initial inputs of Sukhbold et al. (2016) and Ertl et al. (2020), and modeled the radiative transfer in such ejecta out to late times with a more complete set of two-step decay chains whose parent isotopes are ⁴⁴Ti, ⁴⁸Cr, ⁵²Fe, ⁵⁶Ni, and ⁵⁷Ni. Most of the models discussed in this work use this latter set of models – we make a special note whenever we use models s10p0, s12p0, s18p5, and s21p5, which were computed with the older treatment (the difference in the results are small and limited to an O(10%) offset in the flux of Ni and Co lines in the infrared). Given the uncertainties surrounding the exact masses of ⁵⁷Ni and ⁵⁸Ni in stellar explosions (see, for example, the discussion in Sukhbold et al. 2016), this is not a major concern. Tests on how composition variations impact the infrared spectra of core-collapse SNe would help when modeling forthcoming JWST observations – this is left to future work.

One important property of these ejecta models is that they were mixed macroscopically by shuffling in mass space the distribution of the main shells in the 1D, unmixed explosion model (for details, see Dessart & Hillier 2020). The advantage of this technique is that it introduces macroscopic mixing without any microscopic mixing. Composition profiles for such shuffled-shell ejecta are shown for completeness in the appendix for models s15p2 (Fig. A.1) and he6p0 (Fig. A.2). We present a summary of the ejecta properties for our model sample in Table 1.

The new calculations presented here were also performed with an up-to-date model atom for all species, and specifically included the updates that were made in Blondin et al. (2023) for the radiative-transfer modeling of optical and infrared spectra of the Type Ia SN 2021aefx. As discussed in that work, these updates were important for the modeling of weaker lines in the infrared. Essentially, the same core model atom was used for metals in both H-rich and H-deficient models, including the treatment of Ni I and a large model atom for Co and Ni. Ne III was included in all calculations because of the need to allow for the potential formation of [Ne III] 15.550 μm – our simulations revealed that this line is in all cases much weaker than [Ne II] 12.810 μm. Thus, in all models, we included (in order of increasing atomic weight and then ionization) He I (40,51), He II (13,30), C I (14,26), C II (14,26), N I (44,104), N II (23,41), O I (21,51), O II (54,123), Ne I (78,155), Ne II (22,91), Ne III (32,80), Na I (22,71), Mg I (39,122), Mg II (31,80), Si I (100,187), Si II (31,59), S I (106,322), S II (56,324), Ar I (56,110), Ar II (134,415), K I (25,44), Ca I (76,98), Ca II (21,77), Sc I (26,72), Sc II (38,85), Sc III (25,45), Ti II (37,152), Ti III (33,206), Cr II (28,196), Cr III (30,145), Fe I (413,1142), Fe II (228,2698), Fe III (96,1001), Co II (112,1005), Co III (88,1075), Ni I (56,301), Ni II (59,1000), Ni III (47,1000), Ba I (17,31), and Ba II (28,61) – the numbers in parenthesis correspond to the number of super-levels and full-levels (for additional information on super levels, see Hillier & Miller 1998). In H-rich ejecta models, we included H I (26,36). Because of the higher ionization in He-star explosion models, we added C III (62,112), N III (25,53), O III (44,86), Mg III (31,99), Si III (33,61), S III (48,98), Ar III (32,346), and Fe IV (48,272).

Figure 1 shows the CMFGEN results for the bolometric light curves of our model sample (i.e., by bolometric we mean all low-energy photons, specifically those falling between 0.1 and 30 μm). Apart from models s10p0 and s12p0, all models have essentially the same initial ⁵⁶Ni mass, and thus the differences in brightness and fading rate reflect the greater γ-ray escape in He-star explosion models owing to their greater expansion rate (see Table 1 as well as Fig. A.4). These bolometric light curves are thus mostly a reflection of the radiative transfer of γ-rays because at such late times the ejecta are optically thin to optical and infrared photons. That is, the decay-power absorbed (rather than emitted – see dashed line in Fig. 1) equals the bolometric luminosity at any time.

Fig. 1. Bolometric light curves for our model sample. The dashed line corresponds to the decay-power emitted in model s21p5.

3. Results for the Type II SN model s15p2

In this section, we discuss the results for the Type II SN model s15p2 (see ejecta properties in Table 1 and composition profiles in Fig. A.1). Figure 2 shows the s15p2 model spectrum from 1 to 30 μm and at an age of 350 d (the corresponding optical spectrum is shown in Dessart et al. 2021b). For a more direct comparison to observations, we show the flux F_ν versus wavelength for an adopted distance of 10 Mpc. At this time, the continuum flux is essentially negligible relative to a myriad of emission lines of varying strengths. Labels indicate the main contributor to each emission feature in Fig. 2 but there are often multiple lines overlapping or the transitions are multiplets (in this case we report the strongest component). The color coding indicates the flux associated with a given species (for example, the color associated with “Ni” includes contributions from Ni I, Ni II, and Ni III).

Fig. 2. Infrared spectrum of Type II SN model s15p2 at 350 d for an assumed distance of 10 Mpc. A logarithmic scale is used for the x-axis. Labels indicate the main contributor to each emission line (i.e., primarily H, Ne, Ar, and iron-group elements). The color coding indicates the associated species. (See Section 3 for discussion.)

Hydrogen lines are predicted up to and beyond the Humphreys series, with the stronger emission stemming from the Paschen series. We predict Paα (1.875 μm), Paβ (1.282 μm), Brα (4.051 μm), Brβ (2.625 μm), Brγ (2.166 μm), Pfα (7.458 μm), Pfβ (4.653 μm), Pfγ (3.740 μm), Huα (12.369 μm), Huβ (7.500 μm), Huγ (5.907 μm), and transitions from the seventh series at 19.057 μm or the isolated 8.758 μm. A few He I lines are present, with the strongest one at 1.083 μm and weaker ones at 2.058 μm and 4.295 μm. The model predicts no sizable emission from C in the infrared, a weak emission from N with [N I] 1.04 μm, and a few weak lines of O with O I at 1.129 μm and 2.763 μm (this is in stark contrast with the strong [O I] 0.632 μm in the optical range; Dessart et al. 2021b) – these lines overlap with extended emission from Fe. A few weak Na I lines are present at 1.138 and 1.140 μm, 2.206 and 2.208 μm, as well as 5.427 and 5.434 μm. Most of the preceding transitions were permitted. With Ne, there is one strong forbidden-line emission at 12.810 μm, whereas [Ne III] 15.550 μm has a negligible strength (about a hundredth that of [Ne II] 12.810 μm). There are a number of Mg I lines but all weak (similar in strength to those of Na I), with transitions at 1.183 μm, 1.502 μm, 1.711 μm, 2.486 μm (many are multiplets). Silicon is mostly Si I and contributes negligible emission. In contrast, sulfur contributes one strong emission at 1.082 μm (which overlaps with He I 1.083 μm) and another at 25.242 μm (both are forbidden transitions). Argon contributes one strong forbidden line with [Ar II] 6.983 μm.

All other emission lines arise from Fe, Co, and Ni. There is a forest of Fe lines, mostly from Fe II up to two microns, beyond which the model predicts a few isolated, stronger, forbidden lines for Fe I at 24.036 μm, Fe II at 17.931 μm, 24.512 μm, and 25.981 μm, as well as Fe III at 22.919 μm. There are a number of Co lines, with Co II at 10.520 μm, 14.735 μm, 15.455 μm, and 18.799 μm. Co III lines are essentially absent because Co is mostly Co⁺ in our model. Finally, there are a few strong lines of Ni, which arise at 350 d from the sizable abundance of stable Ni isotopes. The strongest Ni line is due to Ni II at 6.634 μm, with weaker transitions (often overlapping with other lines) at 1.939 μm, 10.679 μm, and 12.725 μm. There is also a Ni III transition at 7.347 μm although weaker than [Ni II] 6.634 μm, as well as numerous Ni I transitions at 3.119 μm, 7.505 μm, 11.304 μm and 11.998 μm.

The use of F_ν in Fig. 2 is convenient to show the weak transitions present in the infrared but it can be misleading when evaluating the relative strength of these transitions when widely separated or with respect to optical lines such as [O I] 0.632 μm, Hα, or [Ca II] 0.731 μm. So, to recover a proper perspective on how the flux is distributed across the electromagnetic spectrum, the left panel of Fig. 3 illustrates the normalized cumulative flux distribution, integrated from 0.1 μm in the ultraviolet to 30 μm in the infrared, and shown at multiple epochs spanning 200 d to 500 d. What is first apparent is that the general shape of that distribution evolves little in time, with essentially zero flux falling in the ultraviolet, a fraction rising from 65% to 75% in the optical (below 1 μm), and the rest (35% to 25%) falling in the infrared. The progression of that distribution with wavelength is smooth with a few jumps becoming accentuated as time passes. These jumps correspond to the strongest lines in the spectrum and we can identify [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm, which combined represent over 30% of the total flux emitted at 500 d. In the infrared, the main jumps are due to Paα, [Ni I] 3.119 μm, [Ni II] 6.634 μm, [Ar II] 6.983 μm, the 7.5 μm blend of H I, Ni I and Ni III lines, and finally [Ne II] 12.810 μm.

Fig. 3. Evolution of some spectral properties for model s15p2 over the time span from 200 to 500 d after explosion. We show the evolution of the cumulative flux integrated from 0.1 to 30.0 μm (left) and the evolution of the percentage, fractional luminosity relative to the bolometric luminosity of a few strong lines in the infrared (right). (See Section 3 for discussion.)

The right panel of Fig. 3 shows the fractional flux emitted in a number of strong lines in model s15p2 at epochs between 200 and 500 d. The transitions located at shorter wavelength carry most of the total decay power absorbed in the ejecta (i.e., [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm), whereas the strongest transitions in the infrared carry individually at most 1% (that maximum is reached at 500 d) of the total flux (i.e., [Ni II] 6.634 μm, [Ar II] 6.983 μm, and [Ne II] 12.810 μm). So, in terms of power, the infrared plays a secondary role in model s15p2, although this would be slightly modified in the presence of molecules (see, e.g., Liljegren et al. 2023 and Section 8). Such diagnostic lines of intermediate-mass and Fe-group elements have, however, high scientific content for characterizing the progenitor and explosion properties of Type II SNe. These properties reflect how the decay power emitted is absorbed by the ejecta, in what proportion it is absorbed by the different ejecta regions (i.e., the shells of distinct composition in the preSN model are now the shuffled shells in our ejecta models), and finally what coolants dominate in each of the relevant shells (e.g., the O/Ne/Mg shell versus the H/He shell – see Fig. A.4). In model s15p2, the fractional decay power absorbed by the ejecta drops from 99.5% at 200 d to 66.3% at 500 d. Over the same time span, the fraction of that power absorbed in the H/He shell rises from 38.1% to 55.4%, stays roughly constant at ∼11% in the He-rich shell, drops slightly in the O-rich shell from 21.6% to 18.0% and drops from 29.8% to 15.5% in the Fe/Si-rich shells. These trends reflect the increasing escape of γ rays from the Fe/Si-rich shell which goes mostly to the benefit of the massive H-rich material enshrouding, or macroscopically mixed within the inner, metal-rich ejecta.

At 200 d, the main coolants and processes in model s15p2 are nonthermal excitation of H I (20.3% of the total cooling), He I (9.0%), and O I (3.7%), and collisional excitation of Ca II (17.3%), Fe II (15.7%), Mg II (10.3%), O I (5.6%), Co II (3.1%), C II (2.6%), and Si II (2.0%). Other coolants are typically below 1%. At 350 d, the main coolants and processes in model s15p2 are nonthermal excitation of H I (28.2%), He I (10.5%), and O I (3.5%), and collisional excitation of Mg II (12.7%), Ca II (10.7%), Fe II (10.6%), O I (9.3%), and C II (2.0%). Collisional excitation of metals leading to infrared lines is growing with S I at 1.7%, Ni II at 1.4%, Fe I at 1.7%, Fe III at 0.8%, Co II at 0.9%, but Ar II is still low with only 0.2%. At 500 d, the main coolants and processes in model s15p2 are nonthermal excitation of H I (34.4%), He I (11.6%), and O I (3.8%), and collisional excitation of Fe II (11.8%), O I (9.7%), Ca II (8.4%), Ni II (3.7%), Mg II (3.1%), S I (2.6%), Ni I (1.5%), Ar II (1.2%), and Ne II (1%). These values reflect the fractional line fluxes shown in the right panel of Fig. 3.

Before wrapping up this section, we show in Fig. 4 the morphology in velocity space of two Ni I and Ni II lines predicted in the model s15p2 at 350 d after explosion. To give a more direct connection to the actual flux in each of these lines, we show the quantity λF_λ. For each transition, we show the total flux in each spectral region and draw as a dashed line the flux associated with Ni only. In this figure, we omit [Ni II] 0.738 μm because it is strongly contaminated by [Ca II] 0.731 μm and numerous Fe lines in our ejecta models. Similarly, line overlap affects [Ni II] 1.939 μm and [Ni I] 7.505 μm, whereas [Ni I] 3.119 μm and [Ni II] 6.634 μm are strong and the main contributor to the associated feature. The Ni composition profile shown in the inset confirms that these Ni lines primarily arise from the regions that underwent explosive nucleosynthesis and rich initially in ⁵⁶Ni (and in ⁵⁸Ni) rather than from the H-rich, He-rich, or O-rich material with a primordial abundance of ⁵⁸Ni. There is a small contribution to the Ni II emission from the H-rich material, which is visible in the extended wings of [Ni II] 6.634 μm at ≳2000 km s⁻¹ from line center (Fe II also contributes to cooling the H-rich regions with a primordial composition). The velocity location where the flux jumps in the profile wings corresponds to the ejecta velocity at which the Ni-rich shells reside – these jumps are inherent to the shuffled-shell method (Dessart & Hillier 2020). A less ordered, or even a randomized spatial distribution of the emitting material from the Fe/Si-rich shell would produce a smoother emission profile (i.e., from the spread in projected velocity of that emission). It would, however, have little impact on the actual flux in the lines because γ-rays from ⁵⁶Ni decay permeate throughout the inner ejecta (provided the explosion is not too asymmetric).

Fig. 4. Infrared Ni lines in the models s15p2 at 350 d. We show the quantity λFλ using the total flux (solid) and that due to Ni only (dashed) centered at the rest wavelength of two strong lines of [Ni I] and [Ni II] (see label). A distance of 10 Mpc is assumed. The inset shows the Ni mass fraction versus velocity at 350 d, which thus accounts for all stable Ni isotopes (those that after explosion were neither 56Ni nor 57Ni – these arise here at an 80% level from 58Ni). (See Section 3 for discussion.)

Hence, the profile morphology of Ni lines represents an excellent tracer of the Fe-group elements produced during the explosion and may serve to constrain its distribution in radial or velocity space. Cobalt lines may also be used for that purpose (i.e., as long as the ⁵⁶Co abundance is sizable; by 300 d, Fe is ten times more abundant in the Fe/Si-rich material than Co, and the contrast grows with time), whereas Fe lines may suffer greater contamination from the emission by the abundant, primordial Fe present in H-rich and He-rich material (wherein Fe II acts as an efficient coolant).

4. Results for Type Ib SN model he6p0

The results obtained for the He-star explosion model he6p0 are analogous to those discussed above for model s15p2 (some he6p0 model properties are summarized in Table 1 and illustrated in Fig. A.2). Obvious differences arise from the lack of H I lines from this H-deficient ejecta. Another difference arises from the larger ejecta expansion rate (V_m is 6269 km s⁻¹ compared to 2771 km s⁻¹ in model s15p2, thus about 2.3 times greater) – this arises primarily from the much smaller ejecta mass of 2.82 M_⊙ compared to 10.95 M_⊙, whereas the kinetic energy is comparable (1.1 compared to 0.84 × 10⁵¹ erg). This leads to a weaker continuum flux, stronger and broader emission lines, and a higher ionization of the gas. This higher ionization translates into weaker lines of neutral species (i.e., atoms), and stronger ones from once- and twice-ionized species.

Figure 5 shows the 1–30 μm infrared spectrum of model he6p0 at 200 d after explosion and for an assumed distance of 10 Mpc. Compared to model s15p2, the model predicts (see also labels in figure) stronger and broader He I lines at 1.083 μm, 2.058 μm, and 4.295 μm. No sizable line from C, N, and O is predicted in the infrared, though [O I] 0.632 μm and a weaker [N II] 0.657 μm are present in the optical (see Dessart et al. 2021a, 2023b). As in model s15p2 at 350 d, a few weak Na I lines are present at 1.138 and 1.140 μm, 2.206 and 2.208 μm, as well as 5.427 and 5.434 μm. Weak Mg I lines are also present at 1.183 μm, 1.502 μm, 2.486 μm, and a stronger and more isolated line is predicted at 1.711 μm. The strongest infrared line (as seen in F_ν) is [Ne II] 12.810 μm, with a much weaker [Ne III] 15.550 μm. Because of the higher ionization overall in this model, there are no Si I nor S I lines predicted. Argon contributes one strong emission with [Ar II] 6.983 μm and a weaker one with [Ar III] 8.989 μm. All other emission lines arise from Fe, Co, and Ni. As in model s15p2, there is a forest of overlapping Fe II lines up to two microns but much weaker, likely because of the lower ejecta density (in spite of the earlier date). There is, however, a myriad of forbidden transitions, primarily with [Fe II] at 17.931 μm, 24.512 μm, and 25.981 μm, as well as [Fe III] at 22.919 μm. Cobalt lines are predicted with [Co II] at 1.019 μm, 10.520 μm, 14.735 μm, 15.455 μm, 18.799 μm, and 25.682 μm. There are also lines of [Co III] at 11.885 μm and 16.386 μm. There are a few strong lines of Ni, but with a higher ionization than in model s15p2 (i.e., Ni I lines have a negligible strength). There is the strong [Ni II] 6.634 μm, with weaker transitions of [Ni II] at 1.939 μm, 2.308 μm, 2.911 μm, 10.679 μm, 12.725 μm. Higher in ionization is the strong [Ni III] 7.347 μm, with a weaker component at 10.999 μm.

Fig. 5. Same as Fig. 2 bur for the Type Ib SN model he6p0 at 200 d and a distance of 10 Mpc.

Figure 6 is a counterpart of Fig. 3 but for the model he6p0. The left panel illustrates the normalized cumulative flux distribution integrated from the ultraviolet to the infrared and shown over a time span from 200 d to 520 d. Compared to model s15p2, the distribution is analogous at 200 d, with about 80% of the flux falling in the optical range, but it strongly deviates as time advances. Indeed, the flux progressively shifts to longer wavelength, so that 80% of the flux falls within the infrared at 520 d – the evolution was instead negligible in model s15p2. The vertical jumps in the distribution correspond to lines also present in the s15p2 model but that now have a much greater strength. A major jump occurs across [Ne II] 12.810 μm, and to a lesser extent across [Ni II] 6.634 μm, [Ar II] 6.983 μm, and [Ni III] 7.347 μm.

Fig. 6. Same as Fig. 3 but for model he6p0 over the time span from 200 to 520 d. (See Section 4 for discussion.)

The right panel of Figure 6 shows the evolution of the fractional luminosity relative to the bolometric luminosity of specific lines spanning the optical and the infrared and covering various species, atoms, and ions (see Fig. 3 for the equivalent for model s15p2). As time progresses, the strength of most of these lines, which are forbidden, increases, probably as a result of the decreasing density. But at about 350 d, a shift occurs with a weakening of [O I] 0.632 μm and [Ca II] 0.731 μm, and a continuous strengthening of coolants located in the infrared. By 520 d, [Ni II] 6.634 μm and [Ne II] 12.810 μm are responsible for respectively 10 and 20% of the total ejecta cooling.

This contrast in spectral properties with model s15p2 arises in part from the different distribution of the decay power absorbed within the ejecta. In model he6p0, the fractional decay power absorbed drops from 25.6% at 200 d to only 6.9% at 520 d, time by which the fraction arising from the local absorption of positrons represents 44.5%. This implies that 44.5% of the decay power at 520 d benefits only the regions originally rich in ⁵⁶Ni, that is the regions where explosive burning took place. The decay power absorbed in the Fe/Si-rich regions indeed grows from 25.1% at 200 d to 48.7% at 520 d. Over the same time span, the decay power absorbed in the He-rich material drops from 16.2% to 12.2%, and in the O-rich shell from 58.8% to 39.0%.

At 200 d, the main coolants and processes in model he6p0 are collisional excitation of Mg II (32.0%), Ca II (9.7%), C II (9.2%), O I (8.4%), Co II (5.6%), Fe II (5.4%), Fe III (2.8%), Ne II (2.2%), and Ni II (2.1%), and nonthermal excitation of O II (5.1%) and He I (4.9%). Other coolants are typically at or below 1%. At 350 d, the main coolants and processes in model he6p0 are collisional excitation of O I (15.3%), Ne II (15.2%), Fe II (12.2%), Ca II (9.7%), Ni II (7.9%), Co II (4.4%), Ar II (2.9%), Mg II (2.2%), N II (1.9%), and C II (1.7%) and nonthermal excitation of He I (6.8%), O I (3.8%), and O II (3.2%). At 520 d, the main coolants and processes in model he6p0 are collisional excitation of Ne II (19.6%), Fe II (11.7%), Ni II (10.9%), Ar II (8.3%), O I (5.8%), Ni I (4.9%), Co II (4.8%), Fe I (3.5,%), Si II (1.9%), Fe III (1.8%), and Ni III (1.4%), and nonthermal excitation of He I (9.1%) and O I (5.5%).

Figure 7 is a counterpart of Fig. 4 and illustrates the emission profiles in velocity space of a few [Ni II] and [Ni III] lines in the infrared for model he6p0 at 200 d. Because of the faster ejecta expansion rate, these lines are much broader than in model s15p2, with a maximum extent reaching out to about 5000 km s⁻¹, which represents the velocity of the outermost shell originally rich in ⁵⁶Ni. Because of the large γ ray escape, little Ni emission is predicted outside of these Ni-rich regions, making such line diagnostics excellent tracers of the ashes from explosive burning in the ejecta. The width and strength (i.e., peak flux) of these lines stays essentially constant from 200 to 520 d.

Fig. 7. Same as Fig. 4 but for model he6p0 at 200 d and for [Ni II] and [Ni III] lines ([Ni I] lines have a negligible strength in this model at that time). (See Section 4 for discussion.)

5. Neon as a diagnostic of progenitor mass

In our core-collapse SN ejecta models, the main isotopes of Ne are ²⁰Ne and ²²Ne (other Ne isotopes with atomic mass between 16 and 25 and included in the nuclear network have a negligible yield; for discussion see Sukhbold et al. 2016; Woosley 2019). Neon is essentially unaffected by explosive burning so the Ne composition is set at the time of core collapse. ²²Ne comes from He burning which turns N into Ne, whereas ²⁰Ne comes from carbon burning and thus later in the evolution of the progenitor (see, e.g., Woosley et al. 2002). The abundance of ²²Ne is essentially limited by the primordial abundance of C, N, and O and thus the metallicity, whereas the abundance of ²⁰Ne is much larger and the more so the larger the initial mass or the preSN mass. ²⁰Ne is found in the O/Ne/Mg shell, whereas ²²Ne is primarily located in the He/C/O shell and thus at larger velocities in an unmixed 1D explosion model of the preSN star. In the red supergiant star explosion models of our sample, the ratio M(²⁰Ne)/M(²²Ne) is a few (models s10p0 to s12p0, in which the He/C/O shell is much more massive than the O/Ne/Mg shell) up to a few tens (higher mass models s15p2 to s21p5, where the O/Ne/Mg shell is much more massive than the He/C/O shell). In the models he6p0, he8p0, and he12p0, this mass ratio is ∼30.

Figure 8 shows the Ne yield for the whole set of red supergiant star and He-star explosion models presented in Dessart et al. (2021b) and Dessart et al. (2021a). The Ne yield for a given ZAMS mass is much lower (as for O) for the He-star progenitors of SNe Ibc than for the red supergiant progenitors of SNe II, and the more so for higher adopted Wolf-Rayet wind mass loss rates (see results for the set “x1p0” that uses a nominal mass loss rate and the set “x2p0” that uses a twice higher mass loss rate; for discussion, see Woosley 2019). When the Ne yield is shown versus He-core mass or preSN mass, the models then follow almost a single line, which indicates that the scatter in Ne yield obtained versus ZAMS mass results mainly from preSN mass loss. In any case, the main point shown by Fig. 8 is that the Ne yield increases considerably with He-core or preSN mass, suggesting that Ne could be an alternative to O as a diagnostic of preSN mass, or progenitor mass in H-rich progenitors. For that purpose one would use [Ne II] 12.810 μm in the infrared as an alternative to [O I] 0.632 μm in the optical. One expects [Ne II] 12.810 μm to form under slightly higher temperature or ionization conditions than [O I] 0.632 μm, and in light of the discussion in Sections 3–4, we expect the [Ne II] 12.810 μm emission to thrive in our He-star explosion models.

Fig. 8. Top: total Ne yield versus ZAMS mass for the ejecta arising from He-star progenitor models (black symbols; sample from Dessart et al. 2021a, originally from Ertl et al. 2020; x1p0, x1p5, and x2p0 correspond to different scalings for the adopted wind mass loss rates; see Woosley 2019 for details) and for red supergiant star models (filled red stars; sample from Dessart et al. 2021b, originally from Sukhbold et al. 2016). Bottom: same as top, but now shown versus the preSN mass for He-star models and versus the final He-core mass for the red supergiant star models. (See Section 5 for discussion.)

The top panel of Fig. 9 illustrates the [Ne II] 12.810 μm emission profile in velocity space at 350 d for a subset of our H-rich and H-deficient ejecta models. A distance of 10 Mpc is assumed. Although the different models have a different ⁵⁶Ni mass, which can account in part for the differences shown, the major driver behind the diversity amongst these models are the decay power absorbed by the material from the O/Ne/Mg shell and its ionization. The former sets the limits of how much power may be radiated by lines forming in the O/Ne/Mg shell. The latter tunes how much the cooling of the O/Ne/Mg-rich material may operate through [O I] 0.632 μm (i.e., in case of a lower ionization) or through [Ne II] 12.810 μm (i.e., in case of a higher ionization).

Fig. 9. Properties of the [Ne II] 12.810 μm emission line for a sample of core-collapse SN models. Top: profile morphology versus wavelength at 350 d after explosion for models s12p0, s15p2, s21p5, he6p0, he8p0, and he12p0. The spectra are scaled to a distance of 10 Mpc and an additional scaling is applied to the He-star explosion models for better visibility. Bottom: evolution of the percentage fraction of the bolometric flux that is radiated in [Ne II] 12.810 μm for our full model sample. The color coding differs between the two panels. A similar plot for [O I] 0.632 μm is shown in the left panel of Fig. A.3. (See Section 5 for discussion.)

In all models, the simulated [Ne II] 12.810 μm emission line extends in velocity space out to the outermost layers of the O/Ne/Mg shell, corresponding to ≳2000 km s⁻¹ in our model s15p2 and ∼6000 km s⁻¹ in our model he6p0. In these line profiles, each slanted part of the wing corresponds to a range of velocities in the model ejecta that contributes to emission in the [Ne II] 12.810 μm line. In our model s15p2, the O/Ne/Mg shells are located in regions between 700 and 1000 km s⁻¹, between 1500 and 1700 km s⁻¹, and around 2000 km s⁻¹ (see Fig. A.1). In model he6p0, the contributing shells are between 2000 and 3500 km s⁻¹, between 4000 and 4500 km s⁻¹, and between 5200 and 6000 km s⁻¹ (see Fig. A.2). Each ledge in the profile wing corresponds to regions that contribute no flux (a gap in velocity space with no material from the original O/Ne/Mg shell). Hence, as for Ni lines before (Figs. 4 and 7), the [Ne II] 12.810 μm line can be used to constrain the distribution in velocity space of the ²⁰Ne-rich material from the O/Ne/Mg shell of the preSN star.

The bottom panel of Fig. 9 shows the evolution of the fractional integrated flux radiated in [Ne II] 12.810 μm from about 200 to 500 d and for all models in this study (over that time span, the fractional decay power absorbed by the O-rich material stays roughly constant in all models, except in models s10p0 and s12p0 in which it drops significantly; Fig. A.4). In all cases, the fractional flux radiated in [Ne II] 12.810 μm increases with time as the ejecta density drops. However, in red supergiant star explosion models, the total line flux is weak, rising from about 0.01% at 200 d up to about 1.0% at 500 d. A trend toward greater [Ne II] 12.810 μm fractional flux with main sequence mass is present, but weak and scattered. In He-star explosion models, the [Ne II] 12.810 μm fractional flux is much greater, rising from about 1% at 200 d up to 10–30% at 520 d. Ejecta with a greater expansion rate, and therefore ionization, show stronger [Ne II] 12.810 μm — of all three He-star explosion models, it is the model with the lowest Ne abundance (i.e., he6p0) that has the strongest [Ne II] 12.810 μm line flux. As indicated in Table 1, the Ne yield is comparable between models s15p2 and he6p0, or between models s21p5 and he12p0, and thus the difference in the fractional line flux is dominated by the differences in ejecta structure (i.e., expansion rate, and subsequently ionization) rather than by differences in Ne yield.

These simulations suggest that within an homogeneous set (H-rich progenitors only, or H-deficient progenitors only, but for a similar ionization), the [Ne II] 12.810 μm may be used as a substitute to [O I] 0.632 μm (results for this line are shown in Fig. A.3) for evaluating the progenitor or the preSN mass (the [O I] 0.632 μm diagnostic also suffers from degeneracies in that context). In some models and epochs, the infrared range is essential to constrain ejecta cooling since much power may be channeled into this spectral region, with [Ne II] 12.810 μm acting as a major coolant of the O-rich material (thus weakening dramatically [O I] 0.632 μm as well as the optical luminosity).

6. The diversity of infrared spectral lines in core-collapse SNe

In this section, we document the evolution of a number of infrared lines from about 200 to 500 d in our model sample. This extends what was presented earlier in detail for models s15p2 (Section 3) and he6p0 (Section 4), or for the full model sample in the case of [Ne II] 12.810 μm (Section 5). The evolution of the fractional line flux for the optical lines [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm is shown in the appendix (Fig. A.3) as well as in previous papers (Dessart et al. 2021b,a, 2023b). Because of the numerical procedure, the line strength of Fe-group elements are underestimated at the 10% level in models s10p0, s12p0, s18p5, and s21p5 (see Section 2 for discussion).

Figure 10 shows the evolution of the fractional line flux radiated in [Ni II] 1.939 μm, [Ni I] 3.119 μm, [Ni II] 6.634 μm, [Ar II] 6.983 μm, [Ni III] 7.347 μm, and [Fe II] 17.931 μm. We recover the trends discussed in Sections 3–4, with a general increase with time in all cases except for [Ni II] 1.939 μm in the He-star explosion models (probably because the lower level of that transition lies above the ground state and thus weaker than for example [Ni II] 6.634 μm, which is a transition between the lowest two levels of Ni II). The more massive, slower expanding, denser H-rich ejecta exhibit systematically lower fractional fluxes of at most 1% in those lines compared to H-deficient ejecta. This is not an ionization effect since, for example, the Ni I, Ni II, and Ni III line fluxes are greater in model he6p0 than in any of the H-rich ejecta models at say 350 d. What likely dominates here is the combination of a lower ejecta density and a greater fraction of decay power absorbed by the Fe/Si-rich material in He-star explosion models (see also the bottom-right panel of Fig. A.4).

Fig. 10. Evolution from ∼200 d to ∼500 d of the fractional flux associated with a selection of Ni, Ar, and Fe lines in the infrared for our sample of red supergiant star and He-star explosion models. (See Section 6 for discussion.)

The diversity within each set of H-rich and H-deficient models is more subtle. The decay power absorbed in the Fe/Si-rich material in models s15p2 and s21p5 is up to two times greater than that in the lower mass models, which boosts the power to be radiated in Ni or Ar lines. This is, however, mitigated by the difference in ejecta density (or expansion rate V_m; see Table 1) with V_m being lower in s15p2 and s21p5. In the H-deficient models, the ejecta density seems to drive the difference between the faster expanding model he6p0 compared to the slower expanding model he12p0, the latter developing strong forbidden lines in the infrared at later times.

Because of these complications, it may be difficult to draw firm conclusions on Ar or Ni yields from an observed infrared spectrum of a core-collapse SN. Indeed, the metal yields from models s15p2 and he6p0 are comparable (see Table 1) in spite of the drastic difference between their infrared spectra, largely driven by the dynamical properties of their ejecta. However, the morphology of these metal lines offers an excellent constraint on the spatial or velocity distribution of these metals in the ejecta (as shown for models s15p2 in Fig. 4 and he6p0 in Fig. 7). Stable Ni but also Ar and Fe trace directly the material produced through explosive nucleosynthesis (i.e., ⁵⁶Ni), and could be used to compare to the corresponding distribution of Ne from the O/Ne/Mg material and inferred from [Ne II] 12.810 μm (Section 5).

7. Evidence of ⁵⁶Ni decay inferred from the evolution Fe and Co lines

The evolution of the strength of Fe and Co lines may provide evidence of ⁵⁶Ni decay as the power source in our ejecta models, as an alternative to the evolution of the bolometric luminosity at nebular times. Indeed, the mass of Fe and Co in the ejecta regions originally rich in ⁵⁶Ni varies by about a factor of ten over the time span from 200 to 500 d, and thus for a given power, one would naively expect the Fe lines to strengthen considerably more than the Co lines. The exact mass ratio to consider can be ambiguous since isotopes other than ⁵⁶Fe and ⁵⁶Co are present in the region of formation of Fe and Co lines (i.e., essentially exclusively those regions where explosive nucleosynthesis took place). The decay of ⁵⁷Ni to ⁵⁷Co prevents the Co abundance from continuously dropping at 200–500 d. The challenge is to identify Fe and Co lines for that purpose.

The high ionization of SNe Ia ejecta is conducive to the formation of [Fe III] and [Co III] lines, which have strong emission in the optical at 0.469 μm and 0.589 μm. The evolution of the line strength of these two emission features at nebular times gives unambiguous evidence of ⁵⁶Ni decay as the power source in SNe Ia (Kuchner et al. 1994), or at least in the representative sample of SNe Ia used by these authors. Unfortunately, massive star explosions are typically denser and cooler, so that these Fe III and Co III lines are absent from their optical spectra. However, infrared spectra of massive star explosions contain numerous lines from Fe and Co, as discussed in previous sections and observed in the past (see, for example, Wooden et al. 1993 or Kotak et al. 2009).

Figure 11 shows the evolution of the fractional flux in a variety of Fe and Co lines in the infrared for models s15p2 and he6p0. Compared to the optical lines of Fe III and Co III present in SN Ia spectra, we see that these metal lines (from a mixture of once- and twice-ionized Fe and Co) are weak, representing less than 1% of the total power in model s15p2 and at most a few percent in model he6p0. In most cases and for both models, this fractional line flux increases in time for both Fe and Co lines, and in general this increase is greater for the Fe lines relative to Co lines. This is as expected given the ten fold increase in Fe mass relative to Co mass in the Fe/Si-rich shell of these two models. When estimating the Fe and Co masses, we only account for material within the Fe/Si-rich shells (i.e., where these lines form). Furthermore, we account for all Fe and Co isotopes since they all contribute to radiating the decay power absorbed in these regions.

Fig. 11. Evolution from 200 to about 500 d of the strength of some Fe (solid) and Co (dashed) lines in the infrared and for models s15p2 (top) and he6p0 (bottom).

But we also see great disparity between lines. For example, in model s15p2, [Co II] 10.520 μm and [Fe II] 17.931 μm follow essentially the same trajectory, and thus would have a constant flux ratio from 200 to 500 d, as if there was no change in the relative Fe and Co abundances. But other lines exhibit a flux ratio that follows closely the evolution of the Fe and Co masses in the Fe/Si-rich shell, such as the flux ratio of [Fe II] 17.931 μm and [Co II] 14.735 μm in model he6p0 (Fig. 12). Unfortunately, the same choice of lines for model s15p2 yields a severe mismatch with the evolution of M_Fe/M_Co(Fe-sh). Evidence of ⁵⁶Ni decay as the power source of our ejecta models (and by extension of observed core-collapse SNe) is thus not robustly revealed by the evolution of the Fe and Co infrared emission lines in those core-collapse SN ejecta. One potential reason for this complexity, or potential lack of sensitivity, is the efficient cooling of the Fe/Si-rich material by lines other that Fe and Co, and specifically by Ar II, Ca II, and Ni II lines (stable Ni is as abundant as Co in the Fe/Si-rich material at about 300 d, and it dominates over Co both in abundance and as a coolant at later times). Optical-depth effects might also play a role. Most infrared lines of Fe II and Co II are optically thick but Fe II lines are more optically thick and may be affected for longer.

Fig. 12. Evolution from 200 to 500 d of the ratio of the total flux emitted in [Fe II] 17.931 μm and [Co II] 14.735 μm in model he6p0 (solid). The dashed line shows the corresponding evolution of the total Fe and Co mass within the Fe/Si-rich material, adjusted through an additional power of 1.07.

8. Conclusions

We have presented the infrared, nebular-phase properties of a sample of massive star explosion models that arise from red supergiant stars (models originally from Sukhbold et al. 2016) or from H-deficient, stripped stars (models originally from Woosley 2019 and Ertl et al. 2020). This work complements our previous studies which focused on the photometric and spectroscopic properties in the optical primarily (Dessart et al. 2021b,a, 2022; Dessart et al. 2023b). Because many of the models discussed in those works have comparable properties (say model s14p5 and s15p2 or he5p0 and he6p0) or do not seem to have obvious observational counterparts (e.g., he3p30), we selected a subset of models with greater differences (primarily in initial mass), namely model s10p0, s12p0, s15p2, s18p5, and s21p5 for those with hydrogen at core collapse, and he6p0, he8p0, and he12p0 for those without. To improve the physical consistency of our calculations, we rerun most models with a better account of the full isotopic composition, up-to-date atomic data and model atoms, and treated the two-step decay chains associated with the parent isotopes ⁴⁴Ti, ⁴⁸Cr, ⁵²Fe, ⁵⁶Ni, and ⁵⁷Ni (exceptions are models s10p0, s12p0, s18p5, and s21p5, which accounted only for the ⁵⁶Ni decay chain, as in Dessart et al. 2021b). The time span covers from ∼200 to ∼500 d. Dust is typically expected to form after about 500 d in the inner ejecta of Type II SNe (Sarangi 2022), as also observed in SN 1987A (Lucy et al. 1989), and thus it is largely absent over the time span covered here. Because of space limitation, a comparison to other work is provided in Appendix B.

The focus of this work is on the emission from atoms and ions in the ejecta material, which will serve in a subsequent step as the initial conditions for simulations that account for molecule formation and cooling (McLeod et al. 2024). Molecules such as CO, SiO, SiS impact core collapse SN spectra by potentially modifying the temperature and ionization of the emitting gas and radiating a fraction of the decay-power absorbed in selective ejecta regions (see, e.g., Liljegren et al. 2023). For example, CO forms in the O/C shell (i.e, or in blobs or clumps of C-rich O-rich material that was originally located in the O/C shell of the preSN star) and contributes to cooling that gas, thereby depleting flux that would otherwise have come out in lines such as [O I] 0.632 μm. This molecular cooling is lacking in the present simulations. We thus overestimate some lines fluxes by some amount and do not predict any infrared emission for example from CO, in contradiction with observations such as those of the Type II-peculiar SN 1987A (Spyromilio et al. 1988), the Type II-P SNe 2004et (Kotak et al. 2009) and SN 2017eaw (Rho et al. 2018) or the Type Ic SN 2020oi (Rho et al. 2021).

However, whilst molecular emission can be efficient, its impact is limited at nebular times to the regions where molecules form. For example, in our model s15p2, only 5% of the total decay power absorbed by the ejecta is absorbed by the O/C-rich material. And within those regions, CO molecules may account for some but probably not the entire power absorbed, which means that CO emission is fundamentally limited to be at most just a few percent of the luminosity in that model. In our model he6p0, the decay power absorbed in the O/C-rich shell is 10–15% of the total (in both these models, this fraction represents a quarter of the power absorbed in the O-rich shell), so the impact is potentially greater, but it is not clear molecules can form in great abundance in those lower density, fast expanding ejecta (relative to Type II SN ejecta). In other words, there may be more power for the CO molecules to radiate but their formation is perhaps not as efficient. In the simulations of Liljegren et al. (2023), large CO masses are predicted in the case of a high compression of the O/C shell material, but compression of that material following the ⁵⁶Ni bubble effect (see, e.g., Woosley 1988) is likely less efficient in SNe Ibc than in SNe II². Furthermore, molecular cooling will alter the total flux emitted in lines but much less (and perhaps not at all) their profile shapes, which represent a key diagnostic for example for chemical mixing. Finally, numerous lines are largely unaffected by molecular formation, such as [Ne II] 12.810 μm since neither CO nor SiO forms in the O/Ne/Mg material where ²⁰Ne is abundant but where both C and Si are underabundant (see, e.g., Liljegren et al. 2023). So, there is much to learn from nebular-phase radiative-transfer models even when molecular formation and cooling are ignored.

In the Type II SN models, about 75% of the flux is emitted below 1 μm at all times (this rises to 90% below 2 μm), with little evolution over time, so the infrared contains a small fraction of the total flux (molecular emission would raise this fraction a little – see above). In the infrared, our model spectra exhibit a myriad of H I lines from the Paschen up to the seventh series, He I lines at 1.083 and 2.058 μm, essentially no lines of C, N, or O, weak lines of Na I and Mg I, but the strongest line from intermediate mass elements is [Ne II] 12.810 μm. At higher atomic mass, there are strong emission lines from [S I] at 1.082 μm, as well as 25.242 μm, and from Ar with [Ar II] 6.983 μm.

However, many lines in the infrared arise from Fe, Co, and Ni, which are present throughout the spectrum from 1 to 30 μm. Lines from Ni arise from the presence of stable Ni in the regions where explosive nucleosynthesis took place and include [Ni I] at 3.119 and 7.505 μm, [Ni II] at 1.939, 6.634, 10.679, and 12.725 μm, and [Ni III] 7.347 μm. This wealth of emission lines from iron-group elements, which are predicted to be present in all our models, contrasts with the scarcity of such lines in the optical. One exception is [Ni II] 0.738 μm, which was unambiguously identified in SN 2012ec (Jerkstrand et al. 2015), although this was in a large part due to the unusual weakness of [Ca II] 0.731 μm. The associated emission line profiles can directly reveal the magnitude of ⁵⁶Ni-mixing in the explosion. There are also numerous lines of Fe with [Fe II] at 17.931, 24.512, or 25.981 μm, and Co lines with [Co II] at 10.520 or 14.735 μm. All these Fe, Co, and Ni transitions in the infrared are forbidden and essentially all strengthen with time as the ejecta density drops from 200 to 500 d. Combined with the luminosity drop by a factor of about 20 over that time span, many of these emission lines retain a near-constant flux (or luminosity), which should make them as detectable at 500 d as at 200 d. Results for other Type II SN models in our sample are presented in Figs. A.3–A.4–A.5.

In the He-star explosion models (which would correspond to Type Ib or Ic SNe), the infrared properties are qualitatively analogous but with marked quantitative differences arising from the lower ejecta density (the models have essentially the same explosion energy but a lower ejecta mass relative to the Type II SN models). This leads to a greater escape of γ-rays and a greater fraction of the power arising from the local deposition of positrons within the Fe/Si-rich regions. Consequently, the continuum flux level is weaker, forbidden lines are stronger, and the ionization is higher, leading to strong emission from [Ne II] 12.810 μm, [Ni II] 6.634 μm, [Ar II] 6.983 μm, [Ni III] 7.347 μm, or [Co II] 10.520 μm. The spectral energy distribution exhibits a drastic shift to the infrared with time. In model he6p0, the fraction of the flux emitted beyond 1 μm is 15% at 200 d but nearly 80% at 520 d. This increase in the strength of infrared lines in He-star explosion models is even greater than in model s15p2, such that numerous lines in the infrared evolve essentially at constant luminosity. Results for other He-star explosion models in our sample are presented in Figs. A.3–A.4–A.5.

We investigated whether Ne, whose abundance increases with He-core mass or preSN mass (in a similar fashion to O), could be used as a progenitor mass indicator. There is a tendency for models with a greater Ne mass to exhibit a greater flux in [Ne II] 12.810 μm. However, there is a scatter within each model type (i.e., Type II or Type Ibc SN model), as well as between ejecta models of different initial mass. For example, models s15p2 and he6p0 have a similar Ne yield (and O yield), but the fractional line flux radiated in [Ne II] 12.810 μm is 10–100 times greater in the He-star explosion model at all times (this contrast is reduced when considering the actual line flux since the luminosity is greater in model s15p2 because of the greater γ-ray trapping efficiency; this contrast is also smaller when considering [O I] 0.632 μm – see Fig. A.3). The profile of the [Ne II] 12.810 μm emission can nonetheless serve to map the distribution in velocity space of the Ne-rich material. In He-star explosion models, [Ne II] 12.810 μm and the infrared as a whole is needed to capture the full SN luminosity since the optical represents a decreasing fraction of the total flux as time passes.

We find that Fe and Co lines generally strengthen in time in all models but each line tends to follow its own trajectory. This likely arises from the fact that each line has a different sensitivity to the density, temperature, and ionization evolution of the region where they form – many of these lines correspond to transitions between excited states, which would thus tend to weaken in time for the benefit of transitions tied to the ground state. The other reason is that the Fe/Si-rich material cools through numerous lines of species other than Fe and Co, and in particular Ar, Ca, and Ni, and the role of these other coolants also varies in time. We thus find that the evolution of the Fe and Co line fluxes does not provide a robust evidence of ⁵⁶Ni decay as the ejecta power source of our core-collapse SN models, even though it is the only power source in all the simulations of this work.

The present simulations will now be compared to full optical to mid-infrared observations of nearby core-collapse SNe (see, for example, Dessart et al. 2025). They will also be used as initial conditions for radiative transfer calculations that include molecular formation and cooling (McLeod et al., in prep.).

Acknowledgments

LD acknowledges fruitful discussion with Wynn Jacobson-Galán and Rubina Kotak. This work was granted access to the HPC resources of TGCC under the allocations 2023 – A0150410554 and 2024 – A0170410554 on Irene-Rome made by GENCI, France.

References

Appendix A: Additional information and results for our sample of core-collapse SN models

Figures A.1 and A.2 illustrate the ejecta properties (i.e., composition and velocity versus Lagrangian mass) for the H-rich SN model s15p2 and the H-deficient he6p0 model. The complicated composition results from the inherent onion-shell structure of the original, unmixed 1D ejecta models from Sukhbold et al. (2016) and Ertl et al. (2020), but also from the shuffling in mass space that we applied. As discussed in Dessart & Hillier (2020), the approach in all cases was to split all shells (but only a small fraction of the outermost shell) into three subshells and redistribute them in the same order, starting from the innermost ejecta layer. We thus obtain a pattern of alternating shells of different composition repeated three times, connecting in the outer ejecta with the outermost shell that was largely excluded from the shuffling. The present H-rich (H-poor) ejecta models are the same as those presented in Dessart et al. (2021b) (Dessart et al. 2021a).

To complement the results presented in the main text for infrared lines, we illustrate the evolution of the fractional luminosity emerging in the strong optical lines of [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm in Fig. A.3. Figure A.4 also shows the evolution of the fractional powers absorbed in the ejecta, its regions of distinct composition, as well as the relative contribution of positrons, which are absorbed locally in the ⁵⁶Ni-rich regions, as a function of time. Finally, we show in Fig. A.5 the evolution of the fractional luminosity emerging in a sample of infrared lines in models other than s15p2 and he6p0, for which a similar information is provided in the right panel of Figs 3 and 6.

The explosion energy of our H-rich ejecta models (all taken from Sukhbold et al. 2016) is typically below 10⁵¹ erg and thus probably on the low side (e.g., compared to SN 1987A). So, we tested the impact of increasing the kinetic energy of the ejecta, which for simplicity was achieved by scaling both the radius and the velocity by some factor α (this leaves the SN age unaffected but raises the kinetic energy by a factor α²) and scaling the densities correspondingly by 1/α³ in order the maintain the same ejecta mass and yields. When raising in that manner the ejecta kinetic energy by a factor of two in model s15p2 at 350 d, the decay power absorbed dropped from 87.9% to 68.1%, which led to an overall reduction of the emergent flux throughout the optical and infrared ranges – the decay-power absorbed was reduced following the enhanced escape of γ-rays. A second effect is the different spatial distribution of that absorbed decay power. Lines forming in the H-rich ejecta weakened (e.g., H I lines, including those from the Balmer series in the optical and from the higher series in the infrared), whereas the [O I] 0.632 μm line flux increased relative to Hα. The same effect would have been produced by reducing the ejecta mass. This trend suggests that some diversity in Type II nebular-phase spectra may be driven by variations in the trapping efficiency of their ejecta (i.e., driven by variations in M_ej or E_kin) rather than variation in oxygen yield. However, in this comparison, strong metal lines in the infrared (e.g., [Ni II] 6.634 μm, [Ar II] 6.983 μm, or [Ne II] 12.810 μm) retained a comparable strength or even strengthened, likely because the reduced ejecta density favored the emission from these forbidden lines relative to recombination lines such as those of H I. This complex sensitivity of nebular-phase spectral properties to variations in ejecta properties (not limited to composition variations) is a challenge for radiative-transfer modeling.

Fig. A.1. Composition profile for the H-rich, red supergiant star explosion model s15p2 at 350 d. We show the variation of the mass fraction versus Lagrangian mass in the inner 5 M⊙ for species H, He, C, O, Ne, Mg, S, Ar, Fe, and Ni (with 58Ni shown as a thick dashed line) – these abundances account for radioactive decay where appropriate. The thin dashed line corresponds to the velocity profile (see y-axis at right).

Fig. A.2. Same as Fig. A.1 but for the H-deficient, He-star explosion model he6p0 at 200 d.

Fig. A.3. Same as Fig. 10 but for important line coolants present in the optical range. From left to right, we show the model results for [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm.

Fig. A.4. Evolution of various powers as a function of the total decay power emitted (top-left panel) or absorbed (all other panels) by the ejecta for our set of simulations and over the time span from about 200 to about 500 d after explosion. From left to right and top to bottom, we show the fractional power absorbed in the ejecta, that fraction coming from positrons, and the fractions absorbed in the H-rich (He-star explosion models are excluded since H deficient), He-rich, O-rich ,and Fe/Si-rich material.

Fig. A.5. Evolution of a few strong coolants across the optical and infrared spectral ranges for models other than s15p2 (shown in the right panel of Fig. 3) and he6p0 (shown in the right panel of Fig. 6).

Appendix B: Comparison to earlier work

We can compare the results for our type II SN models to those of Jerkstrand et al. (2012) who employed similar explosion models (12–19 M_⊙ stars from Woosley & Heger 2007 instead of the models from Sukhbold et al. 2016 used here) – similar comparisons have been presented in Dessart et al. (2021b) but mostly for the optical diagnostics. We neglected molecules and dust, whereas Jerkstrand et al. (2012) treated them, although in an approximate manner: molecular cooling was introduced as a prescribed “thermostat” for the material from the O/Si and O/C shells where SiO and CO are expected to form, and molecular and dust emission were both added as a separate infrared emission component to the emergent spectrum (this was later improved by modeling molecules with a molecular network and then treating molecules with the gas in the NLTE solution (Liljegren et al. 2020) or in the full radiative transfer solution (Liljegren et al. 2023)). Chemical mixing was treated with a comparable aim (i.e., by allowing for macroscopic mixing without microscopic mixing) but this is achieved in the Monte-Carlo code SUMO through a virtual grid in which the inner, mixed ejecta is composed of spherical clumps of distinct composition and density and randomly distributed within some prescribed velocity (Jerkstrand et al. 2011), whereas in the grid-based code CMFGEN, this is achieved by means of a shuffling in mass space of spherical shells from the 1D, unmixed ejecta (Dessart & Hillier 2020). Line overlap and blanketing is taken into account in CMFGEN, whereas the Sobolev approximation is used in SUMO. The temperature solution or the sources of atomic data may also deviate between the two codes. The shuffled-shell structure in the CMFGEN input models forces photons emitted within the inner ejecta layers to cross all overlying shells before escaping and may thus overestimate optical-depth effects relative to SUMO models.

One difference is the prediction by SUMO of sizable variations in density amongst clumps of distinct composition, for example when modeling nebular-phase spectra of Type II SNe such as 2004et (Jerkstrand et al. 2012). Similar density variations have been inferred by Li et al. (1993) or Spyromilio & Pinto (1991) in the context of SN 1987A, although Kozma & Fransson (1998b) find little sensitivity of their results to the filling factors given to material of distinct composition. The impact of clumping was studied by Dessart et al. (2021b) but was found to be limited unless it induced an ionization shift. In the largely recombined ejecta of their Type II SN models, Dessart et al. (2021b) found that even for a factor of ten in clumping (i.e., a uniform 10% filling factor), the [O I] 0.632 μm line flux decreased by order 10%, causing a small strengthening of Na I D and Mg I] 4571 Å. Dessart et al. (2021b) also implemented a radial stretching of the Fe-rich shells to mimic the ⁵⁶Ni-bubble effect but this was found to have negligible effect on the resulting Type II SN spectra, in contrast to the findings of Li et al. (1993) or Jerkstrand et al. (2012). Using smooth ejecta models, Dessart et al. (2021a) reproduced satisfactorily the optical observations of a number of well observed Type II SNe one year after explosion, including for example SN 1987A, SN 2012aw, or SN 2013ej. One might argue that infrared metal lines are better diagnostics of ejecta clumping, but the unclumped s15p2 model presented here was found to yield a satisfactory match to infrared metal lines observed in Type II SN 2024ggi (Dessart et al. 2025).

Dessart et al. (2021a) found a greater impact of clumping on Type Ibc model spectra at 200 d because these ejecta were partially ionized when smooth, such that clumping enhanced recombination, for example quenching Fe II emission for the benefit of [O I] 0.632 μm emission from the O-rich material. However, when clumping was raised beyond a factor of a few, such models overpredicted the strength of Na I D and Mg I] 4571 Å. As we surmise in Footnote 2, the density variations arising from the ⁵⁶Ni-bubble effect in SNe Ibc should be weaker than in SNe II. Recently, Ergon et al. (2024) advocated a strong ⁵⁶Ni-bubble effect and clumping in the Type IIb SN 2020acat to explain its nonstandard, overly luminous peak and short rise time, which is perplexing. One alternative to extreme clumping in SN 2020acat may be large scale ejecta asymmetry and high-velocity ⁵⁶Ni along our line of sight.

Different studies also predict distinct emitting regions for a variety of lines. Dessart et al. (2021b) find that [Ca II] 0.731 μm emission arises from the Fe/Si-rich material, whereas Jerkstrand et al. (2012) find that it arises from primordial Ca in the H-rich gas, as found earlier by Li & McCray (1993) and Kozma & Fransson (1998a). When modeling SN Ibc ejecta, which are H deficient, Dessart et al. (2021a) found the same, Fe/Si-rich material at the origin of the [Ca II] 0.731 μm emission. Another difference is the [Ar II] 6.983 μm emission, which was found here to arise from the Si-rich material, whereas Jerkstrand et al. (2012) found that it arises from the He/C shell. [Ne II] 12.810 μm is predicted somewhat stronger in our model s15p2 than in the 15 M_⊙ model of Jerkstrand et al. (2012), with an increasing offset of a factor of 4–5 at 500 d, although this may in part be due to the use of a gray dust in Jerkstrand et al. (2012), something that may overestimate the dust optical-depth in the mid-infrared.

Understanding the origin of these differences is beyond the scope of this paper but the smooth-ejecta model s15p2 presented here yields a satisfactory match to the observed optical-to-infrared nebular-phase spectra of SN 2024ggi (Dessart et al. 2025), which is also a close analog in the optical range and at nebular times of SN 2004et (Dessart et al., in prep). What the present modeling indicates is not that such ejecta are not clumped, but rather than the impact of clumping on nebular-phase spectra is typically small (and function of ejecta ionization) in the simulations with CMFGEN based on shuffled-shell ejecta and thus justifies conducting a study with smooth ejecta. However, it does raise a concern about the robustness of the inferred ejecta clumping from nebular-phase spectral modeling, since if it holds, it should invalidate models that ignore clumping. This requires further study.

All Tables

All Figures

	Fig. 1. Bolometric light curves for our model sample. The dashed line corresponds to the decay-power emitted in model s21p5.
In the text

	Fig. 2. Infrared spectrum of Type II SN model s15p2 at 350 d for an assumed distance of 10 Mpc. A logarithmic scale is used for the x-axis. Labels indicate the main contributor to each emission line (i.e., primarily H, Ne, Ar, and iron-group elements). The color coding indicates the associated species. (See Section 3 for discussion.)
In the text

	Fig. 3. Evolution of some spectral properties for model s15p2 over the time span from 200 to 500 d after explosion. We show the evolution of the cumulative flux integrated from 0.1 to 30.0 μm (left) and the evolution of the percentage, fractional luminosity relative to the bolometric luminosity of a few strong lines in the infrared (right). (See Section 3 for discussion.)
In the text

Fig. 4. Infrared Ni lines in the models s15p2 at 350 d. We show the quantity λFλ using the total flux (solid) and that due to Ni only (dashed) centered at the rest wavelength of two strong lines of [Ni I] and [Ni II] (see label). A distance of 10 Mpc is assumed. The inset shows the Ni mass fraction versus velocity at 350 d, which thus accounts for all stable Ni isotopes (those that after explosion were neither 56Ni nor 57Ni – these arise here at an 80% level from 58Ni). (See Section 3 for discussion.)

In the text

	Fig. 5. Same as Fig. 2 bur for the Type Ib SN model he6p0 at 200 d and a distance of 10 Mpc.
In the text

	Fig. 6. Same as Fig. 3 but for model he6p0 over the time span from 200 to 520 d. (See Section 4 for discussion.)
In the text

	Fig. 7. Same as Fig. 4 but for model he6p0 at 200 d and for [Ni II] and [Ni III] lines ([Ni I] lines have a negligible strength in this model at that time). (See Section 4 for discussion.)
In the text

Fig. 8. Top: total Ne yield versus ZAMS mass for the ejecta arising from He-star progenitor models (black symbols; sample from Dessart et al. 2021a, originally from Ertl et al. 2020; x1p0, x1p5, and x2p0 correspond to different scalings for the adopted wind mass loss rates; see Woosley 2019 for details) and for red supergiant star models (filled red stars; sample from Dessart et al. 2021b, originally from Sukhbold et al. 2016). Bottom: same as top, but now shown versus the preSN mass for He-star models and versus the final He-core mass for the red supergiant star models. (See Section 5 for discussion.)

In the text

Fig. 9. Properties of the [Ne II] 12.810 μm emission line for a sample of core-collapse SN models. Top: profile morphology versus wavelength at 350 d after explosion for models s12p0, s15p2, s21p5, he6p0, he8p0, and he12p0. The spectra are scaled to a distance of 10 Mpc and an additional scaling is applied to the He-star explosion models for better visibility. Bottom: evolution of the percentage fraction of the bolometric flux that is radiated in [Ne II] 12.810 μm for our full model sample. The color coding differs between the two panels. A similar plot for [O I] 0.632 μm is shown in the left panel of Fig. A.3. (See Section 5 for discussion.)

In the text

	Fig. 10. Evolution from ∼200 d to ∼500 d of the fractional flux associated with a selection of Ni, Ar, and Fe lines in the infrared for our sample of red supergiant star and He-star explosion models. (See Section 6 for discussion.)
In the text

	Fig. 11. Evolution from 200 to about 500 d of the strength of some Fe (solid) and Co (dashed) lines in the infrared and for models s15p2 (top) and he6p0 (bottom).
In the text

	Fig. 12. Evolution from 200 to 500 d of the ratio of the total flux emitted in [Fe II] 17.931 μm and [Co II] 14.735 μm in model he6p0 (solid). The dashed line shows the corresponding evolution of the total Fe and Co mass within the Fe/Si-rich material, adjusted through an additional power of 1.07.
In the text

Fig. A.1. Composition profile for the H-rich, red supergiant star explosion model s15p2 at 350 d. We show the variation of the mass fraction versus Lagrangian mass in the inner 5 M⊙ for species H, He, C, O, Ne, Mg, S, Ar, Fe, and Ni (with 58Ni shown as a thick dashed line) – these abundances account for radioactive decay where appropriate. The thin dashed line corresponds to the velocity profile (see y-axis at right).

In the text

	Fig. A.2. Same as Fig. A.1 but for the H-deficient, He-star explosion model he6p0 at 200 d.
In the text

	Fig. A.3. Same as Fig. 10 but for important line coolants present in the optical range. From left to right, we show the model results for [O I] 0.632 μm, Hα, and [Ca II] 0.731 μm.
In the text

Fig. A.4. Evolution of various powers as a function of the total decay power emitted (top-left panel) or absorbed (all other panels) by the ejecta for our set of simulations and over the time span from about 200 to about 500 d after explosion. From left to right and top to bottom, we show the fractional power absorbed in the ejecta, that fraction coming from positrons, and the fractions absorbed in the H-rich (He-star explosion models are excluded since H deficient), He-rich, O-rich ,and Fe/Si-rich material.

In the text

	Fig. A.5. Evolution of a few strong coolants across the optical and infrared spectral ranges for models other than s15p2 (shown in the right panel of Fig. 3) and he6p0 (shown in the right panel of Fig. 6).
In the text