Spectra of supernovae in the nebular phase

TL;DR

This paper reviews the spectral formation processes of supernovae during the nebular phase, linking spectral features to physical conditions, explosion mechanisms, and ejecta morphology, with derivations and applications of models.

Contribution

It provides a comprehensive review of the physical processes, formalism, and modeling techniques for nebular-phase supernova spectra, including new derivations and applications.

Findings

Line profiles encode ejecta morphology.

Luminosity expressions derived for different physical limits.

Physical processes like gamma-ray deposition and ionization are linked to spectral models.

Abstract

When supernovae enter the nebular phase after a few months, they reveal spectral fingerprints of their deep interiors, glowing by radioactivity produced in the explosion. We are given a unique opportunity to see what an exploded star looks like inside. The line profiles and luminosities encode information about physical conditions, explosive and hydrostatic nucleosynthesis, and ejecta morphology, which link to the progenitor properties and the explosion mechanism. Here, the fundamental properties of spectral formation of supernovae in the nebular phase are reviewed. The formalism between ejecta morphology and line profile shapes is derived, including effects of scattering and absorption. Line luminosity expressions are derived in various physical limits, with examples of applications from the literature. The physical processes at work in the supernova ejecta, including gamma-ray deposition, non-thermal electron degradation, ionization and excitation, and radiative transfer are described and linked to the computation and application of advanced spectral models. Some of the results derived so far from nebular-phase supernova analysis are discussed.

Tables

Table 1. Table 1: Commonly observed SN lines that are expected to be in LTE for some part of the nebular phase. For doublets and triplets only one of the components is listed. The line wavelengths (col. 2) are followed by level IDs in parenthesis by energy ordering. The effective collision strengths Υ Υ \Upsilon (column 6), and resulting critical densities (column 7), are for T = 5000 𝑇 5000 T=5000 K and β S = 1 subscript 𝛽 𝑆 1 \beta_{S}=1 (if the line is optically thick t L ​ T ​ E subscript 𝑡 𝐿 𝑇 𝐸 t_{LTE} can be longer). Column eight lists the length of time the line is in LTE, using Eq. 37 ( t L ​ T ​ E subscript 𝑡 𝐿 𝑇 𝐸 t_{LTE} scales with ( M / 1 ​ M ⊙ ) 1 / 3 ​ ( V / 3000 ​ km s − 1 ) − 1 ​ ( x e / 0.1 ) 1 / 3 ​ ( f / 0.1 ) − 1 / 3 superscript 𝑀 1 subscript 𝑀 direct-product 1 3 superscript 𝑉 3000 superscript km s 1 1 superscript subscript 𝑥 𝑒 0.1 1 3 superscript 𝑓 0.1 1 3 \left(M/1~{}M_{\odot}\right)^{1/3}\left(V/3000~{}\mbox{km s}^{-1}\right)^{-1}\left(x_{e}/0.1\right)^{1/3}\left(f/0.1\right)^{-1/3} ). The t t ​ h ​ i ​ n subscript 𝑡 𝑡 ℎ 𝑖 𝑛 t_{thin} values (tenth column) scale with ( V / 3000 ​ km s − 1 ) − 3 / 2 ​ ( f / 0.1 ) − 1 / 2 superscript 𝑉 3000 superscript km s 1 3 2 superscript 𝑓 0.1 1 2 \left(V/3000~{}\mbox{km s}^{-1}\right)^{-3/2}\left(f/0.1\right)^{-1/2} . The mass scalings are based on typical masses in core-collapse supernova models. Estimates for decays to excited states are not attempted, but labelled ’To exc.’. Atomic data sources are given in Jerkstrand et al. ( 2011 ) .

Ion	Line ($\AA$)	$E_u\left(eV\right)$	$g_u$	$A\left(s^{-1}\right)$	${\rm Y}$	$n_e^{crit}$ (cm-3)	$t_{LTE}$ (d)	${ϵ}_n$	$t_{thin}\left(d\right)$
C I	9850 (4-3)	1.26	5	$2.2\times {10}^{-4}$	0.34	$2.6\times {10}^4$	3700	0.04	25${\left(M/0.1M_{\odot }\right)}^{1/2}$
	8727 (5-4)	2.68	1	0.60	0.20	$2.5\times {10}^7$	380		To exc.
N II	6583 (4-3)	1.89	5	$2.9\times {10}^{-3}$	1.4	$8.5\times {10}^4$	2400	0.14	45${\left(M/0.1M_{\odot }\right)}^{1/2}$
O I	6300 (4-1)	1.97	5	$5.6\times {10}^{-3}$	0.06	$3.8\times {10}^6$	640	0.21	170${\left(M/1M_{\odot }\right)}^{1/2}$
	5577 (5-4)	4.18	1	1.3	0.07	$1.5\times {10}^8$	190		To exc.
Si I	1.64 $\mu$m (4-3)	0.78	5	$2.7\times {10}^{-3}$	0.1	$1.0\times {10}^6$	800	1.0	120${\left(M/0.1M_{\odot }\right)}^{1/2}$
	1.10 $\mu$m (5-4)	1.91	1	0.80	0.02	$3.3\times {10}^8$	120		To exc.
S I	1.08 $\mu$m (4-1)	1.15	5	0.028	0.1	$1.1\times {10}^7$	350	2.6	200${\left(M/0.1M_{\odot }\right)}^{1/2}$
	7725 (5-4)	2.75	1	1.8	0.02	$7.4\times {10}^8$	90		To exc.
Ca II	7291 (3-1)	1.70	6	1.3	7.3	$8.8\times {10}^6$	350	162	480${\left(M/0.01M_{\odot }\right)}^{1/2}$
Fe II	7155 (17-6)	1.96	10	0.15	1.0	$1.2\times {10}^7$	290		To exc.
	1.26 $\mu$m (10-1)	0.99	8	$5.0\times {10}^{-3}$	13	$2.6\times {10}^4$	2200	0.21	55${\left(M/0.1M_{\odot }\right)}^{1/2}$
	1.64 $\mu$m (10-6)	0.99	8	$5.1\times {10}^{-3}$	2.2	$1.5\times {10}^5$	1200		To exc.
Co II	1.02 $\mu$m (9-1)	1.22	9	0.054	0.32	$1.2\times {10}^7$	280	1.9	160${\left(M/0.1M_{\odot }\right)}^{1/2}$
	9338 (10-1)	1.32	7	0.023	0.25	$5.3\times {10}^6$	380	0.5	80${\left(M/0.1M_{\odot }\right)}^{1/2}$
Ni II	7378 (7-1)	1.68	8	0.23	1.2	$1.2\times {10}^7$	280	5.5	90${\left(M/0.01M_{\odot }\right)}^{1/2}$

Table 2. Table 2: Decay data for the most common radioactive isotopes. Annihilation radiation from positrons is included in the gamma fraction (contribution in parenthesis). Fractions smaller than 10 − 3 superscript 10 3 10^{-3} have been put to zero. Data from Nuclear Data Sheets 2011.

Decay	$\tau$	$\hat{Q}$ (keV decay-1)	$f_{\gamma }$	$f_{e+}$	$f_{e-}$	$f_{X-ray}$	Source
56Ni $\to$ 56Co	8.77d	1724	0.996	0	0.0025	0.0013	Si burning
56Co $\to$ 56Fe	111d	3732	0.966 (0.053)	0.032	0.001	0
57Ni $\to$ 57Co	35h	2096	0.924 (0.21)	0.074	0.001	0	Si burning
57Co $\to$ 57Fe	391d	143	0.850	0	0.125	0.025
44Ti $\to$ 44Sc	88y	150	0	0	0.070	0.929	$\alpha$-rich freeze-out
44Sc $\to$ 44Ca	6h	2732	0.78(0.35)	0.22	0	0
55Co $\to$ 55Fe	17h	2429	0.822	0.18	0.001	0	Si burning
55Fe $\to$ 55Mn	3.87y	5.6	0	0	0.71	0.29
22Na $\to$ 22Ne	2.60y	2392	0.918(0.38)	0.082	0	0	C burning
60Co $\to$ 60Ni	5.27y	2600	0.963	0	0.037	0	C burning

Equations

$$F_{\nu }=D^{-2}{\int }_{z=-\infty }^{+\infty }{\int }_{y=-\infty }^{+\infty }{I}_{\nu }(x(\nu ),y,z)dydz$$

$$F_{\nu }=2\pi D^{-2}{\int }_0^{p_{max}(x)}{I}_{\nu }(x(\nu ),p)pdp$$

$$F_{\nu }=2\pi D^{-2}{\int }_{r=x(\nu )}^{R_{max}}{I}_{\nu }(r,x)rdr$$

$$I_{\nu }(r,x+\Delta x)={\int }_{x-\Delta x}^{x+\Delta x}{j}_{\nu }(r,x^{\prime })d{x}^{\prime }$$

$$I_{\nu }(r,x+\Delta x)=\frac{ct}{{\nu }_0}{\int }_{-\Delta {\nu }_0/2}^{+\Delta {\nu }_0/2}{j}_0(r)\varphi (r,{\nu }^{\prime }-{\nu }_0)d{\nu }^{\prime }=\frac{ct}{{\nu }_0}{j}_0(r)$$

$$F_{\nu }=2\pi D^{-2}\frac{ct}{{\nu }_0}{\int }_{r=x(\nu )}^{R_{max}}{j}_0(r)rdr$$

$$\boxed{F_{\nu }=2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{\int }_{V(\nu )}^{V_{max}}{j}_0(V)VdV}$$

$$\frac{d{F}_{\nu }(V)}{dV}=2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_0(V)V$$

$$j_0=\frac{1}{4\pi }{n}_{u}A{\beta }_{S}h{\nu }_0$$

$$F_{\nu }=\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{V}_{max}^2{j}_0\left[1-{\left(\frac{V}{V_{max}}\right)}^2\right]$$

$$j_0(V)=j_{max}\exp \left(\frac{-V^2}{2V_0^2}\right)$$

$$F_{\nu }=2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_{max}{\int }_{V(\nu )}^{\infty }\exp \left(\frac{-V^2}{2V_0^2}\right)VdV$$

$$=2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_{max}{V}_0^2\exp \left(\frac{-V^2}{2V_0^2}\right)$$

$$F_{\nu }=2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_0\frac{1}{2}\left[V_{max}^2-{\left(V_{max}-\Delta V\right)}^2\right]\approx$$

$$2\pi t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_0{V}_{max}\Delta V=\text{mbox}constant$$

$$F_{\nu }$$

$$F_{\nu }=2D^{-2}\frac{ct}{{\nu }_0}{j}_0\Delta z{y}_{max}(\nu )=2t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_0{V}_z{V}_{max}\sqrt{1-{\left(\frac{V}{V_{max}}\right)}^2}$$

$$F_{\nu }=D^{-2}\Delta z2{\int }_{ymin(x)}^{ymax(x)}{I}_{\nu }(x(\nu ),y)dy$$

$$F_{\nu }=2t^2{D}^{-2}\frac{ct}{{\nu }_0}{j}_0{V}_z{V}_{max}\times$$

$$\left[\sqrt{1-{\left(\frac{V}{V_{max}}\right)}^2}-\sqrt{h^2-{\left(\frac{V}{V_{max}}\right)}^2}\right],\left(\frac{V}{V_{max}}\right)<h$$

$$\sqrt{1-{\left(\frac{V}{V_{max}}\right)}^2},\left(\frac{V}{V_{max}}\right)>h$$

$$e^{\frac{-t_{delay}}{{\tau }_{decay}}}=e^{\frac{-2V_{max}t/c}{{\tau }_{decay}}}$$

$$F_{\nu }=2\pi D^{-2}{\int }_0^{pmax(x)}{I}_{\nu }(x(\nu ),p)e^{-\tau (x(\nu ),p)}pdp$$

$$\tau (x,p)={\alpha }_{\nu }\times \left(\sqrt{R^2-p^2}-x\right)$$

$$F_{\nu }=2\pi D^{-2}\frac{ct}{{\nu }_0}{\int }_0^{pmax(x)}{j}_0(r(x,p))e^{-{\alpha }_{\nu }\left(-x+\sqrt{R^2-p^2}\right)}pdp$$

$$F_{\nu }=2\pi D^{-2}\frac{ct}{{\nu }_0}{j}_0{\int }_0^{pmax(x)}{e}^{-{\alpha }_{\nu }\left(-x+\sqrt{R^2-p^2}\right)}pdp$$

$$F_{\nu }=2\pi D^{-2}\frac{ct}{{\nu }_0}{j}_0{\int }_{x^{\prime }=R}^x{e}^{-{\alpha }_{\nu }\left(-x+x^{\prime }\right)}\times -\frac{1}{2}{x}^{\prime }d{x}^{\prime }$$

$$=\pi D^{-2}\frac{ct}{{\nu }_0}{j}_0\frac{1}{{\alpha }_{\nu }^2}\left[e^{-{\tau }_0\left(1-\hat{x}\right)}\left(-{\tau }_0-1\right)+{\tau }_0\hat{x}+1\right]$$

$$F_{\nu }=2\pi D^{-2}\frac{ct}{{\nu }_0}{j}_0{\int }_{x^{\prime }=-R}^x{e}^{-{\alpha }_{\nu }\left(-x-x^{\prime }\right)}\times -\frac{1}{2}{x}^{\prime }d{x}^{\prime }$$

$$=\pi D^{-2}\frac{ct}{{\nu }_0}{j}_0\frac{1}{{\alpha }_{\nu }^2}\left[e^{-{\tau }_0\left(1-\hat{x}\right)}\left(-{\tau }_0-1\right)-e^{2{\tau }_0\hat{x}}\left({\tau }_0\hat{x}-1\right)\right]$$

$${\hat{x}}_{peak}=1-\frac{\ln \left(1+{\tau }_0\right)}{{\tau }_0}$$

$${\tau }_S=\frac{1}{8\pi }\frac{g_u}{g_l}A{\lambda }^3{n}_l\left(1-\frac{g_l}{g_u}\frac{n_u}{n_l}\right)t$$

$${\beta }_S=\frac{1-e^{-{\tau }_S}}{{\tau }_S}$$

Full text

1. 1 Introduction

2. 2 Line profiles

3. 2.1 Spherical symmetry

4. 2.1.1 Uniform sphere

5. 2.1.2 Gaussian profile

6. 2.1.3 Thin shell

7. 2.1.4 Thick shell

8. 2.2 Asymmetric distributions

9. 2.2.1 Uniform disk

10. 2.2.2 Disk with hole

11. 2.2.3 Many clumps

12. 2.2.4 Comments

13. 2.3 Radiative transfer effects

14. 2.3.1 Continuous scattering opacity

15. 2.3.2 Continuous absorptive opacity

16. 2.3.3 Line absorption

17. 3 Line luminosities

18. 3.1 Line transfer in the Sobolev approximation

19. 3.2 Local Thermodynamic Equilibrium

20. 3.3 Optical depth

21. 3.3.1 LTE, optically thin case

22. 3.3.2 LTE, optically thick case

23. 3.4 Non-Local Thermodynamic Equilibrium

24. 3.4.1 Thermal collisionally excited lines

25. 3.4.2 Recombination lines

26. 3.5 Discussion

27. 4 Radioactive powering

28. 4.1 Deposition of decay products

29. 4.1.1 Gamma-rays

30. 4.1.2 Leptons

31. 4.2 Degradation of non-thermal electrons

32. 5 Spectral modelling

33. 5.1 Nebula structure

34. 5.2 Level populations

35. 5.3 Temperature

36. 5.4 Radiative transfer

37. 6 Application of spectral models

38. 6.1 Hydrogen-rich SNe

39. 6.1.1 Hydrogen lines

40. 6.1.2 Oxygen lines

41. 6.1.3 Calcium lines

42. 6.2 Stripped-envelope SNe

43. 6.2.1 Magnesium lines

44. 6.2.2 Oxygen lines

45. 6.2.3 Calcium lines

46. 6.3 Thermonuclear supernovae

47. 7 Conclusions

48. 8 Cross-references 11institutetext: Anders Jerkstrand 22institutetext: Queen’s University Belfast, UK / Max-Planck Institute for Astrophysics, Garching, Germany, 22email: anders@mpa-garching.mpg.de

Spectra of supernovae in the nebular phase

Anders Jerkstrand

1 Introduction

As the supernova expands, reduced column densities reduce the optical depths. Recombination removes free electrons, which further reduces the Thomson opacity. Decreasing temperatures lead to lower populations of excited states, reducing the number of optically thick lines and bound-free continua. After a few months, the nebula becomes mostly optically thin and its deep interior becomes visible. It then joins the class of emission line nebulae, which includes also HII regions, planetary nebulae, and Active Galactic Nuclei. The spectrum changes from having a blackbody character with atmospheric absorption lines imposed to an emission line spectrum rich in spectral fingerprints from the newly synthesized elements. The supernova continues to shine due to radioactive decay of isotopes such as 56Ni produced in the explosion. Eventually the supernova enters the “supernova remnant” phase. There is no generally agreed definition of when this occurs, but the term “remnant” usually refers to spatially resolved supernovae of age $\sim 10^{2}-10^{4}$ years where powering occurs by circumstellar interaction or a central pulsar (as in the case of the Crab). By “nebular phase” we usually refer to epochs from a few months to a few years. Most supernovae become too faint to be observed after this time, unless they enter into a phase of strong circumstellar interaction, which brings them towards a remnant phase.

By studying nebular-phase supernovae we can learn about many important properties of the exploded star. The late-time light curves provide constraints on the amount and distribution of radioactive isotopes created in the explosion. The spectral line strengths allows inferrence of ionic masses, emitting volumes, and physical conditions. The line profiles give information on the expansion velocity, the morphology and mixing of the ejecta, and dust formation. Putting all this information together gives us an opportunity to determine the properties of the progenitor stars, test stellar evolution and nucleosynthesis theory, put constraints on the explosion mechanism, and improve our understanding of the formation of black holes and neutron stars.

The modelling and interpretation of nebular-phase supernova spectra are, however, formidable challenges. Complexities include heterogeneous composition throughout the nebula, fast and differential velocity field, non-thermal processes, and a Non Local Thermodynamic Equilibrium (NLTE) gas state. This means that there is a long way to go from having an observed spectrum to inferring physical properties of the ejecta.

This text is written with the aims to explain the basic aspects of nebular-phase line formation, demonstrate the use of simple models and analytic methods, provide guidance to the ingredients and application of advanced models, and to review of some of the results obtained so far. We begin in Section 2 by studying how line profiles are formed in the expanding nebula. In Section 3 we study the connection between physical conditions and line luminosities. Section 4 reviews how powering occurs in the typical scenario of a radioactive energy source, and Section 5 reviews how physical conditions are calculated once the powering situation is known. Section 6 serves to review the availability of advanced models, and to illustrate some output and results of these.

2 Line profiles

By the time the supernova enters the nebular phase, it has reached homologous expansion (purely radial velocities with $V=r/t$). The line broadening due to the expansion (typically a few 1000 km s*-1*) is about three orders of magnitude larger than the line broadening due to thermal motions ($\Delta V_{0}=\sqrt{kT/m}$, a few km s*-1* for the atoms), and the line profiles are therefore determined by the velocity structure of the nebula, but not by its temperature.

Figure 1 illustrates the situation. Here the supernova is represented as a homologously expanding sphere with maximum velocity $V_{max}$. Let the $x$-axis be along the line of sight, and the $y$ and $z$ axes be perpendicular (z into the paper). Define $\Delta\nu=\nu-\nu_{0}$, where $\nu_{0}$ is the line rest frame centre frequency, and $\nu$ is the observed frequency. The observed flux at frequency $\nu$ has contribution from emission in the sheet perpendicular to the line of sight centred at projected velocity $V_{x}=\Delta\nu/\nu_{0}c$ (distance $x=V_{x}t$ from the centre) having thickness $\Delta V_{x}=\Delta x/t$ equal to the intrinsic line width $\Delta V_{0}$ (which we assume is constant for now). Because the thermal line widths are small ($\Delta V_{0}\ll V_{max}$), the observed line profile provides a “scan” through the nebula, each frequency giving a 2D integration of the emission from the sheet at the corresponding resonance depth.

From a given surface segment of the sheet, with area $dA$, the observed flux is $dF_{\nu}=I_{\nu}d\omega$ where $I_{\nu}$ is the specific intensity in the direction of the observer and $d\omega=4\pi dA/(4\pi D^{2})=dA/D^{2}$ is the solid angle subtended by the segment $dA$ as seen by the observer at distance $D$. Integrating over segments, the total flux is

[TABLE]

2.1 Spherical symmetry

If we assume spherical symmetry, it is enough to use a single perpendicular (cylindrical) coordinate, which we denote as $p$. Each annulus has area $dA=2\pi pdp$, so $d\omega=2\pi pdp/D^{2}$, giving

[TABLE]

where $p_{max}=\left(R_{max}^{2}-x^{2}\right)^{1/2}$, where $R_{max}$ is the outer radius of the nebula. Because $p=\left(r^{2}-x^{2}\right)^{1/2}$, $pdp=rdr$ and

[TABLE]

The specific intensity $I_{\nu}$ is obtained by solving the transfer equation $dI_{\nu}=j_{\nu}ds-\alpha_{\nu}I_{\nu}ds$ through the sheet, where $j_{\nu}$ and $\alpha_{\nu}$ are emission and absorption coefficients. For an optically thin line (the optically thick case gives the same solution with a suitable choice of emissivity, see later), using that $\Delta V_{0}/V_{max}\ll 1$ (so $r$ is kept fixed in the integrand)

[TABLE]

In the comoving frame the emissivity is $j_{\nu}(r)=j_{0}(r)\phi(r,\nu-\nu_{0})$, where $\phi(r,\nu-\nu_{0})$ is the intrinsic line profile, normalized so that $\int_{-\Delta\nu_{0}/2}^{+\Delta\nu_{0}/2}\phi(r,\nu-\nu_{0})d\nu=1$ (in our heuristic picture the line profile is a box but the results hold generally). In the observer frame, ignoring $v/c$ effects on the intensity, $j_{\nu}(r,x)=j_{0}(r)\phi(r,\nu^{\prime}-\nu_{0})$, where the comoving frequency $\nu^{\prime}=\nu\left(1-V_{x}/c\right)$. Because $dx=dV_{x}t=\left(d\nu/\nu_{0}\right)ct$, the integral is

[TABLE]

Finally, the line profile is

[TABLE]

Or, using $r=Vt$, we can also write this as

[TABLE]

One may in principle attempt to determine $j_{0}(V)$ by discretizing this equation and fitting a least-squares solution to an observed line profile. Fransson and Chevalier (1989) suggests a variant of this, where $j_{0}(V)$ at a given $V$ is directly estimated from the derivative of the line profile

[TABLE]

So far, this inverse mapping method has not been applied much in the literature. A desirable goal is to obtain the density distribution $\rho(V)$ of the emitting ion. But connecting $j_{0}(V)$ to $\rho(V)$ requires knowing other functions such as temperature $T(V)$ and electron density $n_{e}(V)$.

For an explicit calculation of the line profile from the gas state, one needs to insert the expression for $j_{0}$. It is given by

[TABLE]

where $n_{u}$ is the number density of the upper levels, $A$ is the radiative decay rate, and $\beta_{S}$ is the local escape probability (see Sect. 3.1), which allows also optically thick lines to be treated in this formalism.

Some limiting cases for line profiles are now derived, with illustrations in Fig. 2.

2.1.1 Uniform sphere

With $j_{0}(V)=\mbox{constant}$, Eq. 7 becomes (using interchangeably $\nu-\nu_{0}=\nu_{0}V/c$)

[TABLE]

The line profile is parabolic in shape. The FWHM of this profile is $\sqrt{2}\times V_{max}$. We can easily understand the parabolic shape as arising from the parabolic function describing the area of the resonance sheets.

2.1.2 Gaussian profile

With

[TABLE]

we get

[TABLE]

The line profile is also Gaussian with a FWHM corresponding to the FWHM of the emissivity function (=2.35 $V_{0}$).

2.1.3 Thin shell

With a thin shell, $j_{0}(V)=j_{0}$, between $V_{max}-\Delta V$ and $V_{max}$, we get

[TABLE]

The line profile is a box, bounded by $\pm V_{max}$.

2.1.4 Thick shell

Let the inner edge be $V_{in}$. For $V<V_{in}$, the lower integration limit is set by $V_{in}$ and the flux is therefore independent of $V$ and constant. For $V_{in}<V<V_{max}$, the solution is the same as the uniform sphere case.

[TABLE]

The line profile is flat topped with parabolic wings.

2.2 Asymmetric distributions

We now study a few non-spherically symmetric configurations, starting with disks. For the disks, the discussion is based on an edge-on viewing angle, but the line profiles will keep their shapes, squeezed in width, for other viewing angles.

2.2.1 Uniform disk

For a uniform disk (consider Fig. 1 as now showing a face-on disk with the observer edge-on to the right), the sheet area equals the thickness of the disk $\Delta z$ times the distance between the inner and outer edges ($2\sqrt{R_{max}^{2}-x^{2}}$). Specifically, we get

[TABLE]

This profile is less sharply peaked than the uniform sphere or Gaussian distributions.

2.2.2 Disk with hole

In this case

[TABLE]

where $y_{max}(x)=\sqrt{R_{max}^{2}-x^{2}}$ and $y_{min}(x)=\sqrt{R_{hole}^{2}-x^{2}}$ for $x<h$, where $h=R_{hole}/R_{max}$ is the normalized hole radius, and $y_{min}=0$ otherwise. Then

[TABLE]

The high velocities have the same solution as the uniform disk. Once the hole is crossed ($V<V_{max}h$), the projected area is $\sqrt{1-\left(\frac{V}{V_{max}}\right)^{2}}-\sqrt{h^{2}-\left(\frac{V}{V_{max}}\right)^{2}}$. The second term grows faster with decreasing $V/V_{max}$, and the line profile develops two horns at plus and minus the velocity of the inner edge of the disk.

Figure 2 summarizes the various line profiles discussed.

2.2.3 Many clumps

If the emission comes from a large number of clumps with a more or less random distribution within a spherically symmetric or axisymmetric region, the line profiles will have a global shape determined by the equations above (using the probability distribution as $j_{0}(V)$), but with small-scale structure determined by the random positions of the clumps.

Assume that we have $N$ identical clumps with (comoving) expansion velocities $V_{c}$ (their ’width’ is $V_{c}$), distributed randomly within a sphere of velocity $V_{max}$. The resulting line profile will have squiggles with statistical properties depending on $N$ and $\epsilon=V_{c}/V_{max}$. Figure 3 shows the result of simulating a line profile with $N=10^{3}$ and $\epsilon=0.05$. A quantitative method to use the statistical properties of these fluctuations to infer the statistical properties of the clump distribution was developed by Chugai (1994). Its deployment needs high-resolution spectroscopy, and is therefore suitable for the most nearby supernovae.

2.2.4 Comments

A final few comments are in place. We have derived line profile shapes neglecting effects of time delays and relativistic effects, apart from Doppler shifts. Adressing the first point, photons arriving from the receding side of the SN were emitted a time $t_{delay}\sim 2V_{max}t/c$ before photons from the approaching side. If there is an evolution of emissivity $j_{0}$ on the evolutionary time-scale $t$, there is therefore a damping of the red side of the line compared to the blue of order $2V_{max}/c$, which is a few percent for typical $V_{max}$. But if evolution occurs on a faster radioactive time-scale, the effect is increased by a factor $t_{delay}/\tau_{decay}$.

[TABLE]

For example, if $t=500d$ and $\tau_{decay}=111d$ (as for 56Co), then $t/\tau_{decay}=4.5$, and if $V_{max}=5000$ km s*-1*, the factor becomes 15%. It may be incresed further if gamma-ray leakage further shortens the time-scale over which emissivity changes.

The source movement also increases the blue intensity as the observed time interval is shorter than in the emitting frame, and vice versa decreases the red side. The Lorentz transformations of specific intensity gives a $\left(\nu/\nu^{\prime}\right)^{3}$ factor, or a change of order $\sim 3V_{max}/c$ in our non-relativistic limit. With $V_{max}/c$ of a few percent, the total additional effect is of order 5-10%.

2.3 Radiative transfer effects

So far we have assumed that apart from self-absorption (which can be treated as a modification to the emissivity $j_{0}$ in the Sobolev approximation, see Sect. 3.1), the photons escape freely the nebula. For any given epoch, this is a reasonable approximation beyond some wavelength. At short wavelengths, various opacities remain, however, for years or decades and can alter both line profiles and luminosities.

2.3.1 Continuous scattering opacity

Continuous scattering may occur by free electrons or dust. This leads to one or several “bounces” for the photons in the homologous (Hubble) flow. Because the comoving frame wavelength is always lower than in the original emitting frame, there is on average a net energy loss, and the line profile becomes distorted with an enhanced red tail.

In the nebular phase the electron scattering optical depth is $\tau_{e}\lesssim 1$, and the majority of photons will experience zero or one scattering events. The line profile distortions are thus relatively mild. Figure 4 shows the resulting line profiles for $\tau_{e}=1,2,3$ in a uniform sphere setup, from a Monte Carlo simulation. The scattering opacity produces blueshifts of the peak, although quite weak. For example, at $\tau_{e}=1$, a shift of $\Delta\lambda/\Delta\lambda_{max}=-0.13$ is obtained (or 390 km s*-1* for a 3000 km s*-1* broad line).

If one considers also the frequency redistribution in the comoving frame due to the thermal motions of the scattering particles, a further symmetric broadening on the scale $\sqrt{kT/m}$ occurs as well. Only electrons provide any significant thermal effect due to their low mass. Still, at 5000 K the thermal electron velocity is just 300 km s*-1*, much lower than the supernova expansion velocity, so this effect can be ignored unless very detailed results are needed.

2.3.2 Continuous absorptive opacity

Continuous absorption (photon destruction) may occur by photoionization or dust. At long wavelengths also free-free absorption may occur. Assume that an absorptive opacity is present with an absorption coefficient $\alpha_{\nu}$ (cm*-1*). Consider for example purely absorbing dust, that reemits at mid-infrared (MIR) wavelengths beyond our consideration. The emergent flux is, by extension of Eq. 2:

[TABLE]

The optical depth is

[TABLE]

where we now denote $R=R_{max}$, and so

[TABLE]

For constant emissivity,

[TABLE]

Substitute $x^{\prime}=\sqrt{R^{2}-p^{2}}$, so $p=\sqrt{R^{2}-x^{\prime 2}}$, $x^{\prime}_{min}=R$ and $x^{\prime}_{max}=\sqrt{R^{2}-p_{max}^{2}}=x$. This last equality holds for positive $x$ only. Then $x^{\prime}dx^{\prime}=-2pdp$ and

[TABLE]

where $\tau_{0}=\alpha_{\nu}R$ and we have denoted $\hat{x}=\Delta\nu/\Delta\nu_{max}$.

For negative $x$, use instead $x^{\prime}=-\sqrt{R^{2}-p^{2}}$. Then $p=\sqrt{R^{2}-x^{\prime 2}}$, $x_{min}^{\prime}=-R$ and $x_{max}^{\prime}=-|x|=x$, $x^{\prime}dx^{\prime}=-2pdp$, and

[TABLE]

These expressions can be used to fit line profiles affected by a destructive opacity to estimate $\tau_{0}$. The line profiles for $\tau_{0}=0.5,1$ and 2 are plotted in Fig. 5. There are two main differences to scattering opacities. First, there is no production of a red tail with $\Delta\lambda/\Delta\lambda_{max}>1$. Second, for a given $\tau_{0}$, the peaks are more strongly blueshifted. The location of the peak is found from equating the derivative of Eq. 30 to zero, giving

[TABLE]

2.3.3 Line absorption

If line opacity can be described as a large set of finely spaced lines, its effect can be treated as a continous opacity as described in the sections above. For scattering/absorption by one or a few lines, no generic treatment is possible, and a diverse set of line profiles can be produced. A common scenario is that only the longer wavelength line in a doublet or triplet emerges, as the bluer components are absorbed.

Figure 6 shows six examples of line profiles from a doublet, separated by $\hat{x}=0.5$. Here, equal emissivity in both lines gives a symmetric line in the optically thin case, a distorted line peaking close to the redder line wavelength for optically thick scattering, and again a symmetric but damped and flattened line for optically thick destructive opacity. The dashed, dot-dashed, and dotted curves show the cases when only the first (blue) line emits.

3 Line luminosities

In this section expressions for line luminosities in different physical limits are examined. Physical conditions, such as temperature and electron density, are here parameters. The processes determining these are discussed in the following sections. As a first step we consider how to treat line transfer in the supernova.

3.1 Line transfer in the Sobolev approximation

Supernovae are still dense enough in the nebular phase that line optical depths can be high. It is in general a difficult problem to solve the radiative transfer through optically thick lines. Sobolev (1957) showed that a great simplification can be achieved in the high velocity gradient limit, meaning situations where the velocity gradient of the expanding nebula is large enough that line profiles are traversed on a length scale smaller than the length scale over which physical conditions change. In the case of homologous expansion the optical depth to traverse a line in this limit is given by the Sobolev optical depth

[TABLE]

where $g_{u}$ and $g_{l}$ are the statistical weights of the upper and lower levels, $n_{u}$ and $n_{l}$ are the number densities, $\lambda$ is the wavelength, and $t$ is time. For photons emitted in the line, the average escape probability can be shown to be

[TABLE]

These simplifications mean that we do not have to compute the detailed transfer through each line, but can treat them as infinitely narrow with optical depth and escape probability given by the expressions above.

The total luminosity in the transition is given by the volume integral

[TABLE]

where $\beta$ is the angle-averaged non-local escape probability (with respect to absorption by photoionization, other lines, dust, etc). We will for now proceed by putting $\beta=1$, this factor can easily be added on to the final formulae if relevant.

3.2 Local Thermodynamic Equilibrium

In Local Thermodynamic Equilibrium (LTE), $n_{u}$ and $n_{l}$ are related by the temperature only. LTE requires that both populations and depopulations of the level occur mainly by thermal collisions, or by radiative interaction with a (possibly diluted) blackbody radiation field. As a starting point, one can consider the competition between thermal collisional deexcitation, with a rate $n_{e}Q(T)$, where $Q(T)$ is the collision rate (cm3 s*-1*), and spontaneous radiative deexcitation, with a rate $A\beta_{S}$. The critical density is defined as the electron density above which collisional deexcitations dominate:

[TABLE]

Note the difference from static media in that there may be a dependency on density as well as temperature. The temperature dependence of $Q(T)$ is usually quite weak, so $n_{e}^{crit}$ depends only weakly on $T$. $Q(T)$ is given by (see e.g. (Osterbrock and Ferland, 2006))

[TABLE]

where $\Upsilon(T)$ is the effective collision strength, which depends on the cross section function for the particular transition, but is typically of order unity and with a normally weak temperature dependency.

The electron density in a uniform sphere is

[TABLE]

where $\mu$ is the mean atomic weight, $x_{e}=n_{e}/n_{nuclei}$ is the electron fraction, and $f$ is the filling factor. To have critical densities below typical nebular densities of $n_{e}\sim 10^{9}$ cm*-3* (i.e. LTE), Eq. 35 with a typical value $Q\sim 10^{-7}$ cm3s*-1* shows that transitions need to be forbidden/semi-forbidden ($A\lesssim 10^{2}$ s*-1*) or effectively forbidden/semi-forbidden ($A\beta_{S}\lesssim 10^{2}$ s*-1*).

The second requirement for LTE, population by thermal collisions, requires the upper level (energy $E_{u}$) to be reachable from the thermal pool, $E_{u}\lesssim\mbox{a few}\times kT$. Because $T\lesssim 5000$ K = 0.4 eV in the nebular phase, that means $E_{u}\lesssim$ a few eV.

The lines fulfilling both of these criteria are low-lying forbidden transitions in atoms and ions of e.g. C, N, O, Si, S, Ca, Fe, Co, and Ni. Table 1 lists some important transitions, all clearly detected in nebular SNe (there are more transitions fulfilling the criteria but which have not clearly identified, because they are either too weak or blended with other lines).

Some lines that are not in the table warrant comment. To populate a level mainly by thermal collisions requires not only low enough excitation energy, but also a reasonably high abundance compared to the next ionization stage; recombinations may otherwise become the dominant population mechanism. A recombination-dominated situation often arises for neutral elements with low ionization potential, such as Na I and Mg I. A second issue is for resonance lines. These may fulfill $A\beta_{S}\lesssim 10^{2}$ s*-1*, but their high optical depth means they can become dominated by scattering. A line like Na I 5890 Å for example, has both of these properties and would be poorly modelled as an LTE line.

3.3 Optical depth

For lines connected to the ground multiplet, the optical depths can be estimated by taking $M_{l}\approx(g_{l}/g_{g})M_{ion}$, where $g_{g}$ is the total ground multiplet statistical weight. This assumes that most atoms are in the ground multiplet, and that this is in LTE, both normally good approximations. Then, Eq. 32 gives, for a uniform sphere and ignoring stimulated emission

[TABLE]

Define $\epsilon=4.3\times 10^{29}\mu^{-1}g_{u}/g_{g}A\lambda^{3}$ (atomic constants only) and $\epsilon_{n}=\epsilon/10^{14}$. Then, the time at which $\tau_{S}=1$ is

[TABLE]

This equation has stronger dependencies on $M_{ion}$ and $V$ than the inversion of the LTE equation (37) has, thus, estimates for the duration of the optically thick phase are more uncertain than estimates for the duration of the LTE phase. Column 10 in Table 1 lists the epochs ($t_{thin}$) at which the lines become optically thin for typical masses, $V=3000$ km s*-1*, and $f=0.1$.

For lines decaying to excited states ([C I ] 8727, [O I] 5577, [Si I] 1.10 $\mu$m, [S I] 7725, [Fe II] 7155, [Fe II] 1.26 $\mu$m), transition to optical thinness will usually occur early as populations of excited states are much lower than in the ground multiplet. However, if the A-value is large, this should be more carefully checked. The transition time has in this case a strong temperature-dependency through the sensitivity of populations of excited states.

In some cases two lines from the same ion can be used to determine the optical depth and thereby the density. As an example, the [O I] 6300 and 6364 Å lines arise from the same upper level (2p41D), going to the first and second levels in the ground multiplet, respectively, with $A_{6300}=5.6\times 10^{-3}$ s*-1* and $A_{6364}=1.8\times 10^{-3}$ s*-1*. The ratio of their emissivities $R=j_{0}(6300)/j_{0}(6364)$ is

[TABLE]

Ignoring stimulated emission, $\tau_{S,6300}/\tau_{S,6364}=\left(6300^{3}/6364^{3}\right)A_{6300}/A_{6364}n_{1}/n_{2}g_{2}/g_{1}=0.97A_{6300}/A_{6364}=3.0$, where we have assumed LTE within the ground multiplet and $T\gg\Delta E_{ground}$ (230 K) for the ground multiplet populations, so $n_{1}/n_{2}=g_{1}/g_{2}$. Then

[TABLE]

In the optically thick limit $R\rightarrow 1$, and in the optically thin limit $R\rightarrow 3$.

By studying how the ratio transitions from the thick to thin regime over time, the density of O I can be estimated. Initial application for SN 1987A resulted in $\rho_{\rm OI}(t)=1.7\times 10^{-12}\left(t/100\mbox{d}\right)^{-3}$ g cm*-3* (Spyromilio and Pinto, 1991; Li and McCray, 1992). The mass can be estimated from $M_{\rm OI}=\rho_{\rm OI}\frac{4\pi}{3}V_{max}^{3}t^{3}f_{O}$, but the filling factor $f_{O}$ needs to be determined by some other method.

3.3.1 LTE, optically thin case

If the whole ion is modelled in LTE, $N_{u}=Ng_{u}Z(T)^{-1}e^{-E_{u}/kT}$, where $N$ is the total number of ions, and $Z(T)$ is the partition function. Optical thinness means $\beta_{S}=1$. Then, with $N=M_{ion}/\left(\mu m_{p}\right)$, Eq. 34 becomes

[TABLE]

One may also consider a variant, where it is not assumed that the whole ion is in LTE, but only that the upper and lower states are in LTE with each other (thermal collisions dominate transitions in both directions). Then, if the lower state is in the ground multiplet, and we approximate $M_{g}\approx M_{ion}$, Eq. 42 is recovered with $Z(T)$ replaced by $g_{g}$.

The mass of the emitting ion $M_{ion}$ can therefore be estimated if the temperature can be determined. The most robust mass inferrances can be made from lines with $E_{u}\ll kT$ (so $e^{-E_{u}/kT}\rightarrow 1$), which means $\lambda\gtrsim 3\ \mu m\left(T/5000\ K\right)^{-1}$. For temperatures of a few thousand K these are MIR lines. However, mass ratios may be robostly determined also in the $E_{u}\gtrsim kT$ regime, if $E_{u,1}-E_{u,2}\ll kT$. For example, Jerkstrand et al. (2015b) used the LTE and optically thin formula for [Ni II] 7378 and [Fe II] 7155 to estimate the Ni/Fe mass ratio in core-collapse supernovae.

3.3.2 LTE, optically thick case

In the optically thick limit, $\beta_{S}\rightarrow 1/\tau_{S}$, and

[TABLE]

For LTE, $\frac{g_{l}n_{u}}{g_{u}n_{l}}=e^{-h\nu/kT}$, and so

[TABLE]

The luminosity thus depends on the volume and the temperature, but not the mass. Because the line width is proportional to $\lambda$, the peak spectral flux is proportional to $B_{\lambda}(T)$. The peak flux values of separated lines (and ignoring further transfer effects) would therefore follow a blackbody function.

We can use optically thick LTE lines to determine the volume of the emitting region if we know the temperature, or the temperature if we know the volume. Volume determinations are most robust if $h\nu\ll kT$, for which we get

[TABLE]

Because the volume span can be directly inferred from the line expansion velocities ($\mathcal{V}_{span}=4\pi/3\left(V_{max}t\right)^{3}$), in practice this means that $\mathcal{V}/\mathcal{V}_{span}$ gives us the filling factor $f$ of the emitting region; what fraction of the volume is effectively responsible for emission of that line.

Li et al. (1993a) performed LTE modelling for MIR lines with optical depth effects to estimate $f$ and $T(t)$ in SN 1987A. Jerkstrand et al. (2012) showed that [Fe II] 17.94 $\mu$m and [Fe II] 25.99 $\mu$m fall in the optically thick LTE regime for many hundred days in Type II models, and determined $f$ in SN 2004et. In Type Ia SNe, [Ni II] 6.634 $\mu$m, [Ni III] 7.350 $\mu$m, [Co III] 11.88 $\mu$m, and [Fe II] 17.93 $\mu$m can similarly be useful probes (Maeda et al., 2010).

From Table 1, there are few optical/near infrared (NIR) lines where one can be confident to be in the optically thick LTE regime. [Ca II] 7291, 7323 is a candidate in SNe where the emission is from the synthesized calcium. In the early nebular phase, [O I] 6300, 6364, [Si I] 1.64 $\mu$m, [S I] 1.08 $\mu$m, and [Co II] 1.02 $\mu$m would also in many situations be in this regime.

3.4 Non-Local Thermodynamic Equilibrium

Outside LTE, the limiting formula depends on which mechanisms are assumed to dominate the populations and depopulations of the upper state. In general, one obtains expression involving the number abundance of the feeding state, and physical quantities such as temperature and electron density. If the goal is to estimate ion masses, it is desireable that this feeding state should be a ground state, as one may often approximate the mass of ions in the ground state to equal to the total element mass. Sometimes the coupling to a ground state can occur in several steps, as in recombination cascades.

The most common populating mechanisms are thermal collisions, non-thermal collisions, photoexcitation, and recombination. We will here consider two cases in particular, thermal collisionally excited lines, and recombination lines.

3.4.1 Thermal collisionally excited lines

The statistical equilibrium is, letting $f$ denote the feeding state (which may not be the same as the lower state $l$), and $N=n\mathcal{V}$

[TABLE]

Let us consider the regime $n_{e}\ll n_{e}^{crit}$. Then

[TABLE]

The line luminosity becomes, using Eq. 34

[TABLE]

To determine the mass $M_{f}$ we would need to know both temperature and electron density. If neither is known, these lines have to be used in conjunction with other lines to make combined constraints. For example, if we have two emission lines from the same ion, being pumped from the same feeder state, then

[TABLE]

Thus, such a line ratio may be used to determine the temperature at late times (when NLTE and optically thin conditions are more likely).

3.4.2 Recombination lines

If the upper level is populated by recombinations (directly and/or through cascades via recombinations to higher levels), the equilibrium is (assuming $u$ is predominantly emptied in the $ul$ transition)

[TABLE]

where $\alpha_{u}^{eff}$ is the effective recombination rate. Then

[TABLE]

Determination of the mass $M_{+}$ of the recombining ion requires knowledge of the electron density and the temperature. The temperature-dependency of the effective recombination rates are, however, moderate, and so a determination of $M_{+}n_{e}$ with relatively small uncertainty is possible. Recombination lines from two different elements residing in the same zone allows an estimate of the ratio of ion masses, because the electron density cancels out, and the ratio of effective recombination rates will often be almost temperature-independent.

Emission lines from levels more than a few eV above the ground state are often powered by recombination. Examples include H Balmer lines, some He lines, O I 7774, and Mg I 1.50 $\mu$m. Techniques to use O I 7774 as a diagnostic line are discussed in Maurer et al. (2010), and the use of O I 7774, O I 9263, O I 1.129+1.130 $\mu$m, O I 1.316 $\mu$m and Mg I 1.50 $\mu$m recombination lines in Jerkstrand et al. (2015a).

3.5 Discussion

As we have seen, the luminosity in different line limits depends on different combinations of mass, volume, temperature, and electron density. In general one needs multiple lines, preferably formed in different limits, to break the degeneracies and determine unique values for these parameters.

Analytic line formation limits are easy to use, provide an understanding of line luminosity evolutions, and can be very useful diagnostics for some lines. The key to their application is to provide convincing arguments for the validity of the physical regime. This in turn can come from three different approaches: 1) Demonstrated good validity for any reasonable physical conditions. 2) Inferrence from observations. 3) By inspection of forward models.

An added difficulty in SNe, compared to for example HII regions, is the high expansion velocities which makes many lines blended with each other. Consideration of possible blending contaminations should always be done.

4 Radioactive powering

Baade (1945) discovered that supernovae decline on exponential tails with a time-scale of about 70d. The evolution of an explosion without further energy input produces neither enough luminosity nor such an exponential behaviour. This led Borst (1950) to suggest that there is a radioactive power source. Following this were 20 years of speculation on what radioisotope this could be. Finally Colgate and McKee (1969) provided the right answer; it is the second stage of the radioactive decay chain 56Ni $\rightarrow$ 56Co $\rightarrow$ 56Fe. 56Co decays on a time-scale of 111d to 56Fe, and at the same time the source of iron in the Universe had been identified. The solution came after the demonstration of 56Ni production in high-temperature silicon burning by Bodansky et al. (1968).

In SNe, the relevant decay processes are electron captures (EC) and $\beta^{+}$ decays. In the first step the nucleus transmutes by the conversion of a proton to a neutron. In electron captures this energy is emitted as a neutrino, and in $\beta^{+}$ decays the energy is shared between a neutrino and the positron. The nucleus is usually left in some excited state, which then cascades to the ground state by emission of gamma-rays or by ejecting an inner-shell electron in an internal conversion. Following both electron captures and internal conversions, further emission of X-rays and/or Auger electrons occur as the inner hole is filled. The positrons annihilate with free electrons when they are slowed down to thermal energies, producing two 511 keV gamma rays for antiparallel spins and three gamma-rays with total energy 1022 keV for parallel spins.

The decay power is

[TABLE]

where $N(t)$ is the number of isotopes, $\hat{Q}$ is the average decay energy excluding neutrinos, and $\tau$ is the decay time scale. Note that $\hat{Q}$ is different from the normal $Q$-value, because a significant part of the decay energy will be emitted as neutrinos that escape the remnant.

For a primary isotope (like 56Ni)

[TABLE]

whereas for a secondary isotope (like 56Co)

[TABLE]

where $\tau_{p}$ is the parent decay time scale (e.g. 8.8d for 56Ni). We can write

[TABLE]

where $M_{0}$ is the mass of the primary species (weight $\mu$) synthesised in the explosion. Also other radionuclides are made in the explosion and may provide important input once most of the 56Co has decayed after a few hundred days. These include 57Ni/57Co, 44Ti/44Sc, 55Co/55Fe, 22Na, and 60Co. Decay data for the most common isotopes are listed in Table 2.

4.1 Deposition of decay products

4.1.1 Gamma-rays

The gamma-rays from radioactive decays are at MeV energies. The most important degradation process is Compton scattering. Because the gamma-ray energy is much higher than the excitation and ionization potentials of bound electrons ($\lesssim 1$ keV), the opacity is almost independent of the physical state of the matter - the cross sections for incoherent scattering on various bound electrons are similar to those for free electrons, given by the Klein-Nishina formula. The local absorption thus depends only on the total density of electrons (free + bound).

For $E_{\gamma}>1.022$ MeV pair-production can occur, and for $E_{\gamma}\leq 0.1$ MeV photoelectric absorption becomes important as well. Photoelectric absorption will introduce some dependency on composition, but a relatively small fraction of the initial decay energy will be absorbed by this process.

When the gamma-ray Compton scatters, it loses some fraction of its energy and changes its direction. A detailed solution to the energy deposition requires computation of the multiple scattering processes. The first such calculation in the supernova context was carried out by Colgate et al. (1980). The transfer process can be quite well described with a gray opacity $\kappa_{\gamma}=0.03(Y_{e}/0.5)$ cm2 g*-1* for the case of 56Co gamma-rays, where $Y_{e}$ is the ratio of the total number of electrons to the total number of nucleons. With this opacity, a uniform sphere becomes optically thin to the gamma-rays at

[TABLE]

For Type Ia and some Type Ib/c SNe gamma-ray escape needs to be considered already during the diffusion phase, whereas Type II SNe ($M\sim 10$) enter the tail phase (100-200d) still well before any significant escape occurs; this offers an opportunity to determine the 56Co mass by measuring the bolometric luminosity. Once gamma-ray escape has begun, fitting the observed bolometric light curve in the tail phase to a function such as $P(t)\times\left(1-e^{-t/\tau_{\gamma}}\right)$, can give an estimate for $\tau_{\gamma}$, and from that $ME^{-1/2}$.

If both the radioactive source and ejecta are distributed as a uniform sphere, the mean intensity is, with $\hat{r}=r/R_{max}$

[TABLE]

where $S_{\gamma}$ is the source function. This integral has no analytic solution, but for the $\tau_{\gamma}\ll 1$ limit

[TABLE]

where $j_{\gamma}$ is the gamma-ray emissivity. This function has its maximum at $\hat{r}=0$ and a value of half the maximum at $\hat{r}=1$ (see also Kozma and Fransson (1992)). For $\hat{r}>1$, the intensity falls off as $\hat{r}^{-2}$. Figure 7 shows the $J(\hat{r})$ function in the optically thin limit, as well as the numerical solution for $\tau=10$.

The deposition per mass is $4\pi J_{\gamma}\kappa_{\gamma}$. Thus, Fig. 7 helps to envision how the absorption in a shell depends on its location with respect to the radioactive source (in the optically thin case). Core-collapse explosion models predict density profiles that roughly follow $\rho\propto r^{-2}$ in the inner layers steepening to $r^{-10}$ in the outer. Thus, the bulk of the mass is relatively well described by a set of equal-width shells of similar mass ($4\pi r^{2}\Delta rr^{-2}$ = constant).

Because the gamma field intensity depends sensitively on location, it is important for modelling to use realistic ejecta models capturing the outcome of the mixing processes occurring in the explosion. For Type IIP SNe the situation is relatively satisfactory, because simulations have shown that the metal regions approach the limit of complete macroscopic mixing, and it should be a good approximation to have the 56Ni, Si, and O zones occupy the same volume (macroscopically but not microscopically mixed). For stripped envelope SNe the situation is less clear and modelling becomes more uncertain.

4.1.2 Leptons

The leptons deposit their energy by colliding with the bound and free electrons in the SN ejecta. The effective opacity is around $\kappa_{e+}=$10 cm2 g*-1* (Colgate et al., 1980; Axelrod, 1980) for positrons, and similar for electrons.

Because the leptons are trapped for much longer than the gamma-rays, there is a phase where they take over as the dominant power source. This occurs when $\tau_{\gamma}\approx f_{lep}/f_{\gamma}\ (=0.033\ \mbox{for}\ ^{56}\mbox{Co})$ (Table 2), which for a uniform sphere is at

[TABLE]

A successful identification of this transition phase would allow an estimate of $ME^{-1/2}$. Because 57Co takes over powering at $\sim$1000d, this transition never occurs in Type II SNe for 56Co positrons ($M\sim 10~{}M_{\odot}\rightarrow t_{lep}=1800d$), but in Type Ia and Type Ib/c SNe the predicted transition occurs somewhere between 250d ($M=1.4~{}M_{\odot},E=1$) and 900d ($M=5~{}M_{\odot},E=1$). Unless other effects come into play, the bolometric light curve will then flatten onto the 56Co decay rate, having been steeper before due to $\gamma$-ray escape. In practice, it has proven difficult to demonstrate this transition for several reasons. SNe are dim at late times and estimating the photometry is difficult, with crowded fields and possible contamination by a remnant or companion. The temperatures are low and the true bolometric luminosity is difficult to determine lacking NIR and MIR observations (with molecule and dust formation adding to the problem). Finally, time-dependent effects (freeze-out) come into play at similar epochs and also lead to flattening of the light curve.

The leptons may also escape eventually. If there is no magnetic field, or if the magnetic field is radially combed, their trapping time is about 20 times longer than the gamma-ray trapping time ($t_{lep}^{trap}=620d\left(M/1~{}M_{\odot}\right)\left(E/10^{51}~{}\mbox{erg}\right)^{-1/2}$). But a non-ordered magnetic field even at weak levels will lock the positrons in Larmor orbits and keep them trapped on very small scales:

[TABLE]

where $R_{L}$ is the Larmor radius.

4.2 Degradation of non-thermal electrons

After a few scatterings, the gamma-rays have been downgraded to $\lesssim 0.1$ MeV and are quickly photoabsorbed. In their wake we have a set of mildly relativistic electrons with typical energies between 0.1-1 MeV. These electrons in turn lose their energy by collisions with free and bound electrons leading to ionization, excitation, and heating of the gas. This process is familiar to us from a closer-to-home environment of Aurorae. The collisional ionizations lead to the creation of further high-energy electrons, generally referred to as secondaries. The secondaries have, however, much lower energies than the primaries, $\sim 10-100$ eV.

The degradation process can be modelled with Monte Carlo methods (Shull, 1979; Fransson and Chevalier, 1989), or formal solutions (Lucy, 1991; Xu and McCray, 1991; Kozma and Fransson, 1992). A simplified formal solution is the continuous-slowing-down approximation (Axelrod, 1980), which is quite accurate for the primary electrons, but not for the secondaries.

An important property of this process is that the solution has no strong dependence on the energy of incoming primary particles, as long as its over $\sim$ 1 keV. This allows for generic solvers that do not need information about which high-energy particle source is involved. There is also only a weak dependency on density, leaving the relative abundances of ions and electrons as the main parameters of the problem.

Ionization. Non-thermal ionization determines the ionization balance for the major species and thereby which ions form the spectrum. To solve for the fractions of energy going into different channels, one needs to know the differential cross sections for the collisional ionization processes, i. e. $\sigma(E_{in},E_{out})$. For many elements these differential cross sections are not known, neither experimentally nor observationally, and an approximate treatment is needed. The method of Kozma and Fransson (1992), for example, applies a differential form measured for O I to all ions.

For the integrated cross sections, for electron energies much higher than the excitation/ionization potential $I$, one may use the Bethe approximation for both ionization and excitation cross sections. The cross section for ion $i$, transition $j$, is

[TABLE]

where $E$ is the energy of the electron, $a_{0}$ is the Bohr radius, and $c_{1,ij}$ and $c_{2,ij}$ are constants depending only on level energies and transition strengths. For lower energies, specific calculations or experimental data for the collision cross sections are needed.

For even moderately ionized plasmas with $x_{e}\gtrsim 0.1$, only a moderate fraction, $\sim 1-20\%$ of the energy goes to ionization. This increases towards a plateau value of $\sim$50% for more neutral gas.

Excitation. The coefficient $c_{1,ij}$ is, for bound-bound transitions, proportional to the oscillator strength $f_{ij}$. Allowed transitions ($f_{ij}\sim 1$) will therefore in general be more important than forbidden ones ($f_{ij}\ll 1$). Non-thermal excitations tend to populate high-lying states connected to the ground state by allowed transitions, which leads to UV emissivity. Because of both photoelectric and line opacity, this UV emission does not escape directly, but scatters and undergoes fluorescence. For this reason, models that include non-thermal excitation should preferably also include radiative transfer.

Consideration of non-thermal excitation becomes more important for later epochs; the declining electron fraction gives less heating and more energy going into excitation (and ionization), and the declining temperature makes the resulting fluorescence into the optical and NIR more prominent as thermal emission shifts into the MIR (the so-called “infrared catastrophe”). However, these processes also affect the spectrum in indirect ways. The increased population of excited states can lead to photoionization from these levels and a higher electron fraction. For example, H$\alpha$ in Type II SNe arises mainly as a non-thermal excitation first populates $n=2$, the H I atom is then photoionized from this state, and in the ensuing recombination process H$\alpha$ is emitted (see also Sec. 6.1.1).

Heating. As long as the plasma is ionized to $x_{e}\gtrsim 0.1$, most of the non-thermal energy goes to heating of the gas. The fraction grows monotonically with increasing ionization state. The SN spectrum will therefore be dominated by cooling emission, which is dominated by lines. In Type Ia SNe, the ionization state is high compared to core-collapse SNe, and Axelrod (1980) obtained solutions of the heating fraction close to unity. The ionization balance depends of course critically on the exact value of the ionization fraction, even if it is small.

Figure 8 shows example solutions for pure oxygen and iron plasmas in the low-ionization limit (only neutral and first ionization stages present). The distribution of non-thermal electrons, ionization, temperature, and excitation are all connected equation systems. Thus, one needs to iterate - the non-thermal rates are calculated given an ionization and excitation structure. These are then updated (also following updates of temperature, radiative rates).

5 Spectral modelling

A spectral model defines a scope and a set of physical approximations and assumptions to compute the emergent flux. Spectral models can be divided into two categories 1) Models for individual line luminosities (or sets of lines) and 2) Models for the full spectrum of the SN (or some range of the spectrum). The modeller should ideally not only set up and compute the model, but also use physical reasoning to assess which predictions are robust and which depend more sensitively on ill-constrained assumptions.

There is a large number of physical processes at play in supernova ejecta that each can be treated in several different approximations. The basic ingredients in a nebular model, with different levels of approximation, are listed below together with examples of models implementing them (for a series of papers by the same first-author using the same technique, only the first paper is listed).

• •

Nebula structure

• –

1D, 1-zone (Axelrod, 1980; Ruiz-Lapuente and Lucy, 1992; Mazzali et al., 2001)

• –

1D, multi-zone, no mixing (Fransson and Chevalier, 1989; Ruiz-Lapuente et al., 1995; Liu et al., 1997a; Sollerman et al., 2000; Mazzali et al., 2007)

• –

1D, multi-zone, artificial mixing (Houck and Fransson, 1996; Kozma and Fransson, 1998a; Sollerman et al., 2000; Dessart and Hillier, 2011; Jerkstrand et al., 2011)

• –

2D (Maeda et al., 2006)

• –

3D (Kozma et al., 2005)

• •

Gamma-ray transfer

• –

Local deposition (Dessart and Hillier, 2011)

• –

Gray transfer (Axelrod, 1980; Houck and Fransson, 1996; Liu et al., 1997a; Kozma and Fransson, 1998a; Mazzali et al., 2001; Jerkstrand et al., 2011)

• –

Compton scattering (Fransson and Chevalier, 1989; Eastman and Pinto, 1993)

• –

Compton scattering + photoelectric (Dessart et al., 2013)

• –

Compton scattering + photoelectric + pair-production (Maurer et al., 2011)

• •

Non-thermal processes

• –

Pure heating (Dessart and Hillier, 2011)

• –

Fixed heating and ionization (0.97-0.03) (Ruiz-Lapuente and Lucy, 1992; Mazzali et al., 2001)

• –

Heating, ionization, excitation - continuous slowing down (Axelrod, 1980; Eastman and Pinto, 1993)

• –

Heating, ionization, excitation - full solution(Kozma and Fransson, 1992; Liu et al., 1997a; Jerkstrand et al., 2011; Dessart et al., 2013)

• •

Level populations

• –

LTE

• –

NLTE with collisional excitation/deexcitation and radiative decay terms (Axelrod, 1980; Ruiz-Lapuente and Lucy, 1992; Mazzali et al., 2001)

• –

NLTE with comprehensive set of processes (Fransson and Chevalier, 1989; Houck and Fransson, 1996; Maurer et al., 2011; Jerkstrand et al., 2012)

• –

NLTE with comprehensive set of processes, time-dependent (Kozma and Fransson, 1998a; Dessart et al., 2013)

• •

Temperature

• –

Non-thermal heating, line cooling (Axelrod, 1980; Ruiz-Lapuente and Lucy, 1992; Mazzali et al., 2001)

• –

Non-thermal and photoelectric heating, line, recombination, and free-free cooling (Houck and Fransson, 1996; Jerkstrand et al., 2011; Maurer et al., 2011)

• –

Non-thermal and photoelectric heating, line, recombination, and free-free cooling, time-dependent (Kozma and Fransson, 1998a)

• •

Radiative transfer

• –

Sobolev approximation locally, no global transport (Liu et al., 1997a; Kozma and Fransson, 1998a; Mazzali et al., 2001; Maeda et al., 2006)

• –

Sobolev approximation locally, global transport (Jerkstrand et al., 2011; Maurer et al., 2011)

• –

Full transport locally and globally, time-dependent (Dessart and Hillier, 2011)

The number of atoms and levels modelled, and the quality of the atomic data library have large impact on the accuracy of the model. Before comparing a model to data, it is important to understand the set-up and limitations, and assess the ability to predict any given observable.

Gamma-ray deposition and non-thermal processes have been covered already, and we discuss the other components in more detail here.

5.1 Nebula structure

One of the biggest challenges to SN spectral modelling is to capture the complex structures evidenced by observations of SN remnants, and also obtained in multi-D explosion simulations. Most models are set up in 1D, and some consideration of these mixing effects is needed. It is important to distinguish between microscopic and macroscopic mixing.

Microscopic mixing.

Each nuclear burning stage gives unique nucleosynthesis products, and these do not readily become mixed with each other on atomic scales in the SN ejecta because the supersonic flow freezes out composition on velocity scales $v_{therm}(t)$, which is $\lesssim$ 10 km s*-1*already after a few days. Diffusion is also inefficient in the very early phases when temperatures are still high. The cross section for atomic collisions is of order $\sigma=\pi a_{0}^{2}\sim 10^{-16}$ cm2, where $a_{0}$ is the Bohr radius. For ion collisions Coulomb interactions give somewhat larger values for typical temperatures. The “optical depth” is $\tau=\sigma nL$ = $\sigma nVt$. Of order $\tau^{2}$ scatterings are needed to travel a distance $\tau$. The time between scatterings is $\Delta t=\lambda/v_{therm}=1/\sigma nv_{therm}$. The total time to diffuse through the nebula is

[TABLE]

For any reasonable values of $M,V_{max},T$, the diffusion time is longer than the age of the Universe, and neglegible mixing occurs.

This idea is supported by several lines of evidence. For instance, models for molecule formation show that microscopically mixed ejecta fail to produce the observed amounts of molecules (e.g. Liu and Dalgarno, 1996; Gearhart et al., 1999). Spectral models with full microscopic mixing also generally fare worse than models without such mixing (Fransson and Chevalier, 1989).

Macroscopic mixing. Macroscopic mixing, on the other hand, is known to occur vigorously during the explosion of many progenitor structures. The dominant mechanism is typically Rayleigh-Taylor mixing, which arises behind the reverse shocks created at composition interfaces where $\rho r^{3}$ shows a positive derivative. The mixing effects are strong in explosion simulations of H-rich SNe, with red supergiants becoming more mixed than blue supergiants (e.g. Herant and Woosley, 1994). Simulations of He core explosions show strong mixing for low-mass progenitors, but weaker for high-mass ones (e.g. Shigeyama et al., 1990). The lack of O/He and He/H interfaces in Type Ic SNe makes mixing harder to obtain.

The implication is that hydrogen-rich SNe are poorly represented by 1D explosion models, and either 2D/3D models or 1D models with artificial macroscopic mixing are preferred. For stripped-envelope SNe there is more freedom within current uncertainties in the mixing processes. In Type Ia SNe, Rayleigh-Taylor instabilities also occur, with buoyant hot 56Ni rising up between downflows of C and O, in particular in deflagration models (e.g. Röpke et al., 2007).

5.2 Level populations

At nebular times LTE can in general not be assumed. The level populations have to be determined by solving a set of NLTE equations describing the various populating and depopulating mechanisms. Letting $n_{ijk}$ denote the number of ions of element $i$ ionization state $j$ and excitation state $k$, each equation has the form

[TABLE]

The rates $\mathcal{R}$ depend in general on $n_{e}$, $T$, radiation field $J_{\nu}$, and non-thermal electron distribution $n_{nt}(E)$. In the Sobolev approximation, they also depend on the level populations $n_{ijk}$ themselves. Charge transfer (ct) also implies dependencies on number densities of colliding partners. All normal processes have $i^{\prime}=i$, but radioactive decay and molecule and dust formation allows $i^{\prime}\neq i$. If we ignore this for now, so dropping the $i$ index, we can break up all processes populating/depopulating a level into the groups of ionization, recombination, excitation, and deexcitation. Letting $R$ denote a radiative process, and $C$ a collisional, these group can be further broken up into

[TABLE]

which gives an overview of the physical processes involved. Let $N$ be the typical number of excited states modelled in each ion. Any given level can transition to any other level in that ion, in the ion below, or above, i.e. to $\sim 3N$ other levels. The number of transition pairs is $N_{ion-stages}3N^{2}/2\sim 10^{4}N_{ion-stages}$ per chemical species for $N=100$. The modelling of just a few species would require specification of $\sim 10^{5}$ transition rates $\mathcal{R}$, each in turn containing 15 individual rates (Eq. LABEL:eq:apa), some of which are temperature dependent. The only rates that are available in batch reading on this scale are spontaneous radiative decay rates. Photoionization and non-thermal ionization rates requires cross sections as function of energy, and simplifying treatments, or inclusion of only some transitions, is necessary. The following considerations help:

• •

For excited non-metastable states, $\mathcal{R}_{exc}=\mathcal{R}_{ion}=\mathcal{R}_{rec}=0$, and $\mathcal{R}_{de-exc}=R_{spont.\ decay}$ are accurate approximations. These states, with $A\gg 1$, are empitied on fractions of a second ($\tau\sim 1/A$) by spontaneous decay, much shorter than any other process can operate.

• •

For ground states, $C_{thermal.\ ion}=0$ and $C_{3-body\ rec.}=0$ are accurate approximations because nebular temperatures are too low for these processes to be competitive.

Steady-state corresponds to ignoring the time-derivative terms in Eq. 67, so a (non-linear) algebraic equation system follows, which is straightforwardly solved by Newton-Raphson iteration. This is justified as long as the reaction time-scales are short compared to dynamic and radioactive time-scales. The slowest reaction is typically radiative recombination, where $t_{rec}\sim 1/\left(\alpha n_{e}\right)\sim 1d\left(10^{7}\mbox{cm}^{-3}/n_{e}\right)$.

The simplest possible scheme for NLTE excitation solutions is to include spontaneous radiative rates, treated in the Sobolev approximation, and thermal collisional rates. This gives a reasonable approximation for plasmas with $x_{e}\gtrsim 0.1$ because heating then accounts for most of the energy deposition (Sect. 4.2), and the reprocessing of this energy is mainly done by collisional cooling. The next natural step is to add recombination, which enables accounting of a non-thermal ionization energy and production of recombination lines. Coupling the radiation field to the NLTE solutions requires addition of line absorption and photoionization rates. The most advanced models include also non-thermal excitations, photoionization from excited states, and charge transfer reactions.

Charge transfer (CT) is an important process at high densities and low ionization; a relatively unique astrophysical environment found in supernovae at late times. By this process electrons jump from one ion to another, for example

[TABLE]

Both species in the outgoing channel may be left in excited states. The important role of this process in governing the ionization balance in SNe was pointed out by Meyerott (1978). The fast rates are of order $\Gamma_{CT}\sim 10^{-9}$ cm3 s*-1*, to be compared with radiative recombination rates of order $\Gamma_{RR}\sim 10^{-12}$ cm3s*-1*. Thus, even an atomic abundance of $10^{-3}$ per electron can lead to recombination by charge transfer dominating the ionization balance.

The process is typically fast between neutral atoms and singly ionized ions, and if the energy defect is small. Ion-ion reactions require high temperatures due to the Coulomb barrier, and are of less importance. Modelling is hampered by many unknown or poorly known rates. Most rates that have been calculated involve H and He, for application in HII regions and planetary nebulae, but SN modelling mainly needs rates between metals.

When the reaction is fast in both directions, the effect of CT is to link the ionization balance of the less abundant element to the more abundant. As an example, for a solar composition gas the primordial O attains the same ionization balance as H by this mechanism.

5.3 Temperature

The temperature evolution is obtained by solving the first law of thermodynamics, which for homologous expansion and pressure from a perfect monoatomic gas is

[TABLE]

where $H(t)$ and $C(t)$ are the heating and cooling rates (which depend on the NLTE solutions $n_{ijk}$). The heating is usually dominated by non-thermal heating, and the cooling by collisional excitation of fine-structure lines (which then decay radiatively).

The thermal equilibrium approximation corresponds to setting $H(T)=C(T)$ and solving the resulting algebraic system for $T$; this is a good approximation if the cooling/heating time-scales are short compared to the dynamic and radioactive decay time scales. This is usually fine for several hundred days into the nebular phase. For example, de Kool et al. (1998) and Kozma and Fransson (1998a) find this approximation to hold for 600-800d for the H zone and several years or decades for the metal zones in models for SN 1987A. Once collisional cooling becomes slow, adiabatic cooling takes over, and the solution to Eq. 70 is $T(t)\propto t^{-2}$ (ignoring the last term involving $x_{e}$ which typically has a small effect).

Once molecules begin to form, the cooling becomes more efficient and the temperature can become significantly lower than models without molecules would suggest. While a few studies have adressed formation and cooling of molecules in single-zone setups, they are typically not included in multi-zone models.

5.4 Radiative transfer

Supernova ejecta remain optically thick to line blocking below $\sim$ 4000-5000 Å for years or even decades. Fig. 9 illustrates this, showing the photon escape probability in a Type IIP SN at 300d. Thus, to model the appearance of the SN at short wavelengths, as well as lines at longer wavelengths influenced by fluorescence, this transfer must be considered. There are two conceptual approaches - solution of the radiative transfer equation and Monte Carlo simulations of the photon propagation.

The transfer equation can be solved by the method of characteristics. A common approach is to solve high-resolution moment equations (which remove the angular variable) with a variable Eddington factor, which is determined iteratively by a formal solution to the transfer equation (e.g. Hillier and Dessart, 2012).

The basic concept in the Monte Carlo approach is to follow energy packets as they propagate through the ejecta, using random numbers to determine their interactions and trajectory changes. (e.g. Jerkstrand et al., 2011). Two properties of the SN make the machinery relatively straightforward to set up - homology (Hubble flow) makes the expansion isotropic from any point in the ejecta, and the large velocity gradient allows the Sobolev approximation which means that the photon can be transported from line to line, without having to consider line overlap. Different algorithms may be chosen depending on the desired degree of coupling to the gas state.

Iteration between NLTE solutions and radiative transfer have good convergence properties in late-time SN environments, because the dominant excitation and ionization processes are typically collisional (non-thermal and thermal). In practise, only lines and photoionization continua from ground states and (effectively) meta-stable states need to be considered.

What happens to a photon absorbed in a resonance line? For allowed lines, the Regemorter formula gives an estimate of the effective collision strength, $\Upsilon\sim 2\times 10^{6}\lambda^{3}Ag_{u}P_{y}$, where for purposes here we can take $P_{y}=0.1$. Then $Q_{res}\sim 5\times 10^{-16}A\left(\lambda/3000~{}\AA\right)^{3}$ cm3 s*-1*. The ratio of probabilities for collisional deexcitation (i.e. thermalization) and scattering is then

[TABLE]

Because $n_{e}\lesssim 10^{8}$ cm*-3* in the nebular phase only the most optically thick lines ($\tau_{S}\gtrsim 10^{7}$) may suffer thermalization. Note that even if collisional deexcitation is unimportant, multiple scatterings may still lead to enhanced thermalization by other processes such as photoionization as the photons increase their dwell time in the ejecta. Fluorescence gives a capping of this trapping by moving photons out to optically thinner wavelengths after a handful of scatterings.

Modern codes treat the complex global transport through hundreds of thousands of lines. However, simplified approaches are also conceivable. Axelrod (1980) devised an approximate method to treat the global line transfer in a single-zone model. In this method, a position and angle averaged probability to be absorbed in any other transition (using the Sobolev approximation) is calculated (compare with the average continuum escape probability of Osterbrock (1989)). To avoid coupling between the radiation field and the optical depths (the method is in the “no coupling” regime), it is then assumed that absorption occurs only from low-lying states whose populations are calculated without any coupling to the radiation field. The fluorescence is finally obtained by coupling the high-lying state populations through the averaged probabilities, which gives a linear system of equations. The method showed relatively good results for reproducing the blue regions of Type Ia SNe. While the method is outdated for modelling of quite well understood SN classes, it can still find use in rapid exploration of model scenarios for new SN types.

6 Application of spectral models

Spectral models may be used to compare with observations to find the best matches in terms of progenitors, nucleosynthesis, and mixing. Overall, the application of models may be divided into:

• •

Test the viability of a particular explosion model for an observed SN or SN class

• •

Estimate a physical parameter by optimising a model over this parameter

• •

Identify lines

• •

Determine physical conditions and regimes

• •

Estimate bolometic correction factors

Apart from these concrete purposes, a model gives opportunity to improve our general understanding of what is going on in the SN. As in all science, the end use of such fundamental information cannot always be predicted.

6.1 Hydrogen-rich SNe

Nebular multi-zone spectral models of H-rich SNe have been presented by Fransson and Chevalier (1987); Kozma and Fransson (1998a, b); de Kool et al. (1998); Dessart and Hillier (2011); Dessart et al. (2013); Jerkstrand et al. (2011, 2012, 2014, 2015b). This ’first principles’ model set is complemented by a set of more parameterized models studying specific line formation in H (Xu et al., 1992), He (Li and McCray, 1995), O (Li and McCray, 1992), Ca (Li et al., 1993a), and Fe (Li et al., 1993b) in SN 1987A.

The output of the Fransson and Chevalier (1987); Kozma and Fransson (1998a, b) and de Kool et al. (1998) models (all for SN 1987A) are line luminosity tracks, whereas the other papers present spectra. Figure 10 shows two model examples, from Dessart et al. (2013) and Jerkstrand et al. (2014). Current-day models are quite successful at reproducing the main spectral features, such as Mg I] 4571, Na I, [O I] 6300, 6364, H$\alpha$, [Fe II] 7155, [Ca II] 7291, 7323, and Ca II NIR, as well as the underlying quasi-continuum. Three of the most prominent lines in nebular Type IIP spectra are H$\alpha$, [O I] 6300, 6364, and [Ca II] 7291, 7323. The formation of these lines is reviewed in some more detail below.

6.1.1 Hydrogen lines

The models show that the H lines are formed mainly by recombination. The ionization is achieved in a two-step process; first $n=2$ is populated by non-thermal excitations, and then follows a Balmer photoionization. The ionizing photons come partly from H itself in the form of two-photon emission (from $n=2$), and partly from line emission from other elements.

Apart from H$\alpha$, there are no other Balmer lines produced for several hundred days. Kirshner and Kwan (1975), and later Xu et al. (1992) and Kozma and Fransson (1992), demonstrated that this is because the Balmer series is optically thick (“Case C”), and conversion occurs through other series (e.g. H$\beta$ ($n=4$ to $n=2$) converts to Pa$\alpha$ ($n=4$ to $n=3$) + H$\alpha$ ($n=3$ to $n=2$). The Balmer lines are thick because $n=2$ acts like a meta-stable state due to the enormous optical depth in Ly$\alpha$ ($\tau_{S}\sim 10^{10}$).

Modelling of the H lines is also complicated by breakdown of the Sobolev approximation for the Lyman lines - at very high $\tau_{S}$ values one must consider line overlap as well as the possibility that the photons scatter into regions of the ejecta with different composition. Considering Ly$\alpha$ escape, the $n=2$ population may decrease, which in turn decreases the number of photoionizations and the recombination line luminosities throughout. Considering Ly$\beta$ escape, the H$\alpha$ luminosity specifically may decrease while not impacting the other H lines.

An important result from H line modelling in SN 1987A is that H-gas must occupy most of the volume of the central few 1000 km s*-1*, and absorb about half of the gamma rays (Xu et al., 1992). Combined with the lack of flat-topped line profiles, the picture is clear that strong mixing occurs in the explosion and drags hydrogen envelope material down towards low velocities, in line with multi-D hydrodynamic simulations. Kozma and Fransson (1998b) estimate a total H-zone mass of 8 $M_{\odot}$ in SN 1987A (about half of which is H, the rest He). Standard stellar evolution models with a $\sim$8 $M_{\odot}$ H-zone give satisfactory fits to H-lines also in many other Type II SNe (Jerkstrand et al., 2012, 2014).

6.1.2 Oxygen lines

The [O I] 6300, 6364 lines are efficient cooling lines and typically reemit a large fraction of the heating of the oxygen-zone layers. Their strengths are thus indicators of the oxygen mass, which from stellar evolution models is strongly dependent on the progenitor mass. Model luminosity tracks for different $M_{ZAMS}$ are presented in Jerkstrand et al. (2012, 2014) (see also Fig. 11). The method is particularly useful for observations around a year after explosion, when the [O I] lines become optically thin, temperatures and densities are in a regime favorable to [O I] 6300, 6364, and we are confident that $M_{O}\approx M_{OI}$. The current picture points to relatively limited amounts of oxygen produced in Type IIP SNe, and an origin in $M_{ZAMS}\lesssim 20$ $M_{\odot}$ stars (Fig. 11). Another method to determine progenitor masses using the widths of the [O I] lines has been proposed by Dessart et al. (2010).

The luminosity of the O I doublet depends in general on $M_{\rm OI},e^{-\Delta E/T(t)}$, and $n_{e}(t)$ in the optically thin NLTE phase (Sect. 3.4). LTE is valid until quite late times removing the $n_{e}$ dependency, but the exponential dependency on $T(t)$ makes it difficult to make any meaningful estimates of $M_{\rm OI}$ by inverse analytic modelling. Between $T=3000$ K and $T=6000$ K, the inferred O I mass changes by a factor of 45! Thus, one needs strong constraints on the temperature. One approach is to use the [O I] 5577 / [O I] 6300, 6364 ratio as a thermometer. The measurement of [O I] 5577 is intricate as it is weak in Type II SNe, but the method has been demonstrated to be feasible and in good agreement with forward modelling (Jerkstrand et al., 2014).

Few other oxygen lines are distinct in Type II SNe. One exception is O I 1.13 $\mu$m which is often observed to be strong, well above any plausible recombination luminosity. Its high strength likely arises as a fluorescence effect when Ly$\beta$ photons are absorbed in O I 1025, a process that gives indication of mixing of O and H clumps on small scales (Oliva, 1993).

6.1.3 Calcium lines

Kirshner and Kwan (1975), and later Li et al. (1993a), demonstrated that the [Ca II] 7291, 7323 lines in H-rich SNe arise mainly from primordial (solar-abundance) Ca in the H-zone, and not from the synthesized calcium. The reason is the much larger mass of the H zone ($\sim$ 10 $M_{\odot}$) compared to the Si/S/Ca zone ($\sim 0.1$ $M_{\odot}$), which leads to more energy being reprocessed there. In addition, [Ca II] 7291, 7323 is a very efficient cooling channel that emits a large fraction of deposited energy even when the calcium abundance is low (solar). Thus, the [Ca II] emission lines probe mainly the amount of energy processed by the hydrogen gas, and one expects them to be stronger for higher H-zone masses and for closer mixing with 56Ni. Kozma and Fransson (1998b) could confirm the dominance of H-zone emission of [Ca II], but pointed out that a contribution by the O-zone can occur in models using efficient convection and overshooting, where some Ca has been mixed out into the O layers.

Both [Ca II] 7291, 7323 and the Ca II NIR triplet are initially formed by thermal collisional excitation. The Ca II NIR / [Ca II] 7291, 7323 ratio decreases with time as the temperature decreases. At later times fluorescence following UV pumping in the Ca II HK lines becomes important. Li et al. (1993a) estimate an epoch of $\sim$350d for this, while Kozma and Fransson (1998b) finds 500-800d in their models. Fluorescence takes over in the triplet first, roughly when $T\lesssim$ 5000 K. In the limit that both [Ca II] 7291, 7323 and Ca II NIR are driven by HK pumping, their intensity ratios will approach unity.

The high optical depth of the triplet lines means that the Ca II 8498 line can scatter in the Ca II 8542 line (1500 km s*-1* separation), and the Ca II 8542 line can scatter in the Ca II 8662 line (4200 km s*-1* separation). To what extent this happens depends on the distribution of the H-zone gas and its temperature. In Figure 10 (top) it can be seen that full scattering to the 8662 line has occurred in the model, whereas the observed spectrum shows a distinct 8542 line that has not scattered. One should also note that O I 8446 will scatter in the Ca II NIR triplet, and [C I] 8727 can blend with its red wing.

6.2 Stripped-envelope SNe

The first stripped-envelope nebular spectral models were calculated by Fransson and Chevalier (1989), of 4 and 8 $M_{\odot}$ He core explosions (Type Ib SN). The models produced emission lines of Mg I] 4571, [O I] 5577, Na I D, [O I] 6300, 6364, O I 7774, [Ca II] 7300, Ca II NIR + [C I] 8727, as typically observed. The models were important in self-consistently predicting the range of temperatures and ionization expected, $T=2000-8000$ K and $x_{e}\sim 0.1$. It was shown that an assumption of strong microscopic mixing of the SN ejecta gives model spectra discrepant with observations, consistent with the expected inefficiency of atomic diffusive mixing (Section 5.1). The implication is that some caution is needed in interpreting single-zone models for core-collapse SNe, and for highest accuracy multi-zone explosion modelling is needed.

Multi-zone spectral models of Type IIb SNe (only small amounts of hydrogen left) have been computed by Houck and Fransson (1996), Maurer et al. (2010), and Jerkstrand et al. (2015a). These models were used to study the behavour of H$\alpha$ from the hydrogen shell, the emission from the helium envelope, and the oxygen line brightness. The SNe modelled so far show nucleosynthesis pointing to an origin in low or intermediate mass stars, where the hydrogen envelope stripping must have occurred by Roche lobe overflow to a companion. Problems to explain apparent late-time H-alpha emission appear to have been resolved by including [N II] 6548, 6583 emission from the helium envelope (Jerkstrand et al., 2015a), which mimics broad H$\alpha$ emission.

The first nebular Type Ic models were calculated by C. Kozma and presented in Sollerman et al. (2000). These were 1D models of energetic explosions of 6 and 14 $M_{\odot}$ CO cores, with the aim to model SN 1998bw - the first SN associated with a gamma-ray burst. As for other SN classes, 1D explosion models gave poor reproduction of observed line profiles, demonstrating that some kind of mixing or asymmetric explosion occurs also in Type Ic SNe. The model set was important in showing that models that give good fits to early-time light curves and spectra can still be rejected from late-time comparisons.

A single-zone approach has been taken in a series of papers by P. Mazzali (e.g. Mazzali et al., 2001, 2004, 2010) (see Figure 12 for an example). The majority of SNe analysed with single-zone models have oxygen mass estimates of 0.5-1.5 $M_{\odot}$, which suggests a low/intermediate progenitor mass range. Some show only a few tenths of solar mass of oxygen (e.g. Sauer et al., 2006), indicating that Nature can produce bare CO cores of as low mass as 2 $M_{\odot}$. The Axelrod method has also been implemented in 2D (Maeda et al., 2006). 2D explosion models offer a way to reproduce narrow lines from intermediate-mass element and broad lines from iron-group elements, as observed in SN 1998bw, if the viewing angle is close to pole-on (Figure 12, right).

Below some further comments are given on the formation of lines of magnesium, oxygen, and calcium, which are typically strong in stripped-envelope SN spectra.

6.2.1 Magnesium lines

The most prominent Mg line is Mg I] 4571, which is the first transition in Mg I. Models show that the neutral fraction is typically small ($x_{\rm MgI}\sim 10^{-3}$) due to efficient photoionization both from ground state and excited states (Jerkstrand et al., 2015a). This puts the line formation in a regime where both cooling and recombination can be important. The luminosity can show a dramatic increase at density thresholds above which cooling takes over. The line may be affected by line blocking to a significant extent, the efficiency of which increases rapidly blueward of 5000 Å. This produces an asymmetric line profile with a blue-shifted peak.

Another distinct line produced that has been observed in many SNe is Mg I 1.504 $\mu$m. This is a pure recombination line, free of line blocking. As models often give $M_{Mg}\approx M_{MgII}$, this line becomes proportional to the magnesium mass, and can be used as a diagnostic of this (Jerkstrand et al., 2015a, 2017).

6.2.2 Oxygen lines

The [O I] 6300, 6364 doublet is, as in Type II SNe, an important coolant of the oxygen layers and therefore a good diagnostic of their mass. Due to the higher expansion velocities in stripped-envelope SNe compared to H-rich SNe, the lines enter the optically thin regime at an earlier epoch, with $t_{\tau 6300}=1$ at 360d for $M=1$, $f=0.1$ and $V=3500$ km s*-1*(using Eq. 39). Because the core expansion is faster than the line separation between the 6300 and 6364 lines (3047 km s*-1*) the lines are blended and the line ratio has to be estimated by fitting the single blended feature.

[O I] 5577 is often distinct early on, and the [O I] 5577 / [O I] 6300, 6364 ratio may be used as a thermometer. This method breaks down quite early as [O I] 5577 falls out of LTE. This occurred after $\sim$150d in a model grid of Type IIb SNe (Jerkstrand et al., 2015a). The [O I] 5577/[O I] 6300, 6364 ratio has been shown to depend on clumping (e.g. Maurer et al., 2010), and holds some promise to be used as diagnotic for this.

When O I 7774, O I 9264, O I 1.129+1.130 $\mu$m and O I 1.316 $\mu$m are formed in the recombination regime, their strengths may be used to estimate the quantity $n_{e}f^{1/2}$, under the assumption that $n_{\rm OII}\sim n_{\rm e}$. Effective recombination rates to be used for this have been calculated by Maurer et al. (2010) and Jerkstrand et al. (2015a). However, models also demonstrate complications and deviation from this regime. The meta-stable behaviour of many excited states in O I leads to significant optical depths for hundreds of days, and this produces a scattering contribution to the lines. At high densities cooling can occur in the O I 7774 transition (Maurer et al., 2010), significantly boosting it over its recombination luminosity. While O I 7774 is relatively free from contaminations, O I 9263, O I 1.13 $\mu$m and O I 1.31 $\mu$m can be significantly blended with other lines (Jerkstrand et al., 2015a).

6.2.3 Calcium lines

[Ca II] 7291, 7323 is typically the main coolant of the explosive oxygen burning ashes, and the calcium line strengths therefore relate to how much energy is reprocessed by these layers. Standard explosion models ejecting $\sim$0.1 $M_{\odot}$ of Si/S/Ca material have been shown to give good agreement with Type IIb SNe (Jerkstrand et al., 2015a).

The Ca II NIR lines, other the other hand, are after $\sim$200d formed mainly by fluorescence following HK absorption, mainly in the 56Ni ashes.

The [O I] 6300, 6364 / [Ca II] 7291, 7323 ratio is often used in the literature as a diagnostic of core mass, but there are many issues with this method. The calcium emission in stripped-envelope SNe comes from an explosively made region, whose size and distribution depend on the explosion energy. As long as the link between progenitor mass and explosion energy is unknown, a diagnotic involving calcium is not robustly linked to core mass.

Figure 13 shows model tracks of the [O I] 6300, 6364 /[Ca II] 7392, 7323 ratio from the models of Fransson and Chevalier (1989), Houck and Fransson (1996) and J15a. The grids indicate increase of this ratio with helium core mass, for other parameters fixed (and keep the caveat mentioned above in mind). The Fransson and Chevalier (1989) grid shows a decrease of the ratio with higher core velocity, and all grids shows strong time dependency. As discussed in Fransson and Chevalier (1989), the very low ratios for the 4.5 $M_{\odot}$ He cores are due to a large fraction of calcium in the O/Mg zone, something seen more rarely in modern explosion models. They are therefore marked in parenthesis.

6.3 Thermonuclear supernovae

The first Type Ia nebular models were computed by Axelrod (1980). These models are single-zone, pure iron-group composition, with relatively simple physics. The model spectra nevertheless showed good resemblance with observed Type Ia spectra, strengthening the exploding white dwarf idea to explain Type Ia SNe. In particular were the strongest observed lines reproduced by Fe II emission at 5200 and 7155 Å, Fe III emission at 4700 Å, and Co III emission at 5900 Å. The cobalt line evolution showed consistency with a declining abundance of cobalt, providing strong support for the 56Co decay model. Axelrod could also show that the “quasi-continuum” between the strong lines was formed by the overlap by a large number of weaker spectral lines, not needing any true continuum emission. Axelrod’s technique, albeit without the global radiative transport component, has been adapted and applied to analyze Type Ia spectra in many subsequent papers. Single-zone models likely work reasonably well due to the homogenous composition of Type Ia SNe. Ruiz-Lapuente and Lucy (1992) developed a method to use Type Ia nebular spectra to determine the extinction towards the SN.

Multi-zone nebular models of the fast deflagration simulation W7 of Nomoto et al. (1984) have been presented in several papers (Ruiz-Lapuente et al., 1995; Liu et al., 1997a; Sollerman et al., 2004; Leloudas et al., 2009; Maeda et al., 2010; Maurer et al., 2011; Mazzali et al., 2011; Fransson and Jerkstrand, 2015) using a variety of codes. Liu et al. (1997a) showed that the very high 58Ni abundance in W7 (Ni/Fe = 4 times solar) gave a much too strong [Ni II] 7378 line, and SNe with MIR data show the same discrepancy for [Ni II] 6 $\mu$m (Leloudas et al., 2009). However, Ruiz-Lapuente and Lucy (1992) estimate a ratio of 3 in SN 1995G, suggesting some variety. Liu et al. (1997a) show how the ionization state and temperature increases with velocity coordinate (to about 8000 km s*-1*), varying from FeI+FeII-dominated mix and $T\sim 3000$ K at the centre to Fe III + Fe IV dominated mix and $T\sim 7000$K at 8000 km s*-1*, at 300d. Liu et al. (1998) discusses the physical ionization mechanisms at late times, finding non-thermal ionization and charge transfer to be more important processes than photoionization. Models generally show continuously decreasing ionization with time, with Fe I becoming the dominant ion after about two years (Sollerman et al., 2004; Fransson and Jerkstrand, 2015). Even before this time, Fe I plays an important role in the optical radiative transfer (Axelrod, 1980).

Ruiz-Lapuente (1996), Liu et al. (1997a), and Mazzali et al. (2015) computed sub-Chandrasekhar models, pointing out that they are hotter and more ionized than W7 in the central regions, producing somewhat higher Fe III/Fe II line ratios. They also have weaker [Ni II] lines due to the smaller amount of 58Ni. Liu et al. (1997a) found a better fit for the sub-Chandra model to a series of observed SNe, whereas Ruiz-Lapuente (1996) favoured the W7 model for SN 1994D. Eastman and Pinto (1993) and Liu et al. (1997b) computed spectra of delayed detonation model DD4 of Woosley and Weaver (1994). Pure deflagration models were studied by Kozma et al. (2005), who found strong O I and C I lines in these to be inconsistent with observations. Two examples of spectral models are shown in Fig. 14.

For the first few hundred days, the ejecta are hot and ionized. Most radioactive energy is converted to heat, and is re-emitted in the optical and NIR. The spectrum is dominated by efficient cooling lines of iron, cobalt, and nickel. At later times, the ejecta pass through a rapid phase of cooling as the cooling switches from excited multiplets to ground multiplet transitions. The bulk of SN emission then moves to the mid-infrared, an effect dubbed the infrared catastrophy. The optical and near-infared regions then become dominated by fluorescence of UV emission (Fransson and Jerkstrand, 2015), and an accurate radiative transfer treatment is necessary. The fluorescence maintains optical output at about 20% of the bolometric luminosity.

7 Conclusions

Supernovae in the nebular phase provide a wealth of information about the interiors of the exploded stars. From the first simple models devised some 35 years ago, we today have sophisticated codes capable of testing stellar evolution, explosion, and nucleosynthesis models in high level of detail. The nucleosynthesis of isotopes such as oxygen, magnesium, and nickel can be estimated, providing constraints on both progenitors and explosion physics. Some signatures need advanced models for interpretation, whereas some can be analyzed with simpler analytic approaches. The vast majority of core-collapse SNe show nucleosynthesis consistent with an origin as $M_{ZAMS}\lesssim 20$ $M_{\odot}$ stars, and there may be a shortage of events from more massive progenitors compared to standard IMF expectations. Line profiles may be used to probe the morphology of the ejecta for the various elements, which puts multi-D hydrodynamic models to the test. Line luminosities provide information about masses, emitting volumes, and physical conditions. The evolution of light curves provides results on masses of radioactive isotopes such as 56Co, 57Co and 44Ti. By continuing efforts at determining the nucleosynthesis and inner ejecta morphology of SNe of different types much remains still to be learned about supernovae and the origin of the elements.

Acknowledgements.

I would like to thank J. Spyromilio, C. Fransson, K. Maeda, R. McCray, S. Taubenberger, J. Sollerman, P. Mazzali, S. J. Smartt, M. Ergon, and P. Ruiz-Lapuente for useful comments on the manuscript.

8 Cross-references

• •

Dust and molecular formation in supernovae

• •

Spectra of supernovae during the photospheric phase

• •

Nucleosynthesis in spherical explosion models of core collapse supernovae

• •

The Multi-Dimensional Character of Nucleosynthesis in Core Collapse Supernovae

• •

Nucleosynthesis in Hypernovae: Gamma Ray Bursts

• •

Nucleosynthesis in Thermonuclear Supernovae