LONG-DURATION SUPERLUMINOUS SUPERNOVAE AT LATE TIMES

The origin of the newly discovered class of superluminous supernovae (SLSNe) remains unclear. Reaching peak brightness of $\leqslant -21$ mag and emitting ∼10⁵¹ erg of energy (Gal-Yam 2012), these transients provide a challenge to explain within current paradigms in stellar evolution and explosion theory (see Sukhbold & Woosley 2016, for a recent overview). Some SLSNe evolve rapidly on a timescale of weeks, whereas others evolve slowly on a timescale of months. There is diversity in spectral appearance, with some displaying hydrogen (Ofek et al. 2007; Smith et al. 2008; Gezari et al. 2009; Miller et al. 2009; Chatzopoulos et al. 2011; Benetti et al. 2014; Inserra et al. 2016b), but most appearing like Type Ic SNe (Gal-Yam et al. 2009; Pastorello et al. 2010; Chomiuk et al. 2011; Quimby et al. 2011; Chornock et al. 2013; Howell et al. 2013; Inserra et al. 2013; Lunnan et al. 2013, 2016; Nicholl et al. 2013, 2014, 2015; Vreeswijk et al. 2014).

An important avenue for progress is to study and model SLSNe at late times, when their inner ejecta become visible and constraints on the nucleosynthesis and core structure can be placed. The strengths and shapes of observed emission lines may be compared with models, and element masses, distributions, and mixing may be inferred. There is now a small but growing body of data that presents an opportunity to engage in such analysis. Due to the intrinsically high late-time brightness of the slow-evolving (long-duration) transients with rise and decay times of several months, late-time observations are most easily obtained for these SNe.

The first late-time spectra of a Type Ic SLSN were presented by Gal-Yam et al. (2009), for SN 2007bi (z = 0.1279). The spectra were taken at +367d and +471d post-peak rest-frame, and displayed distinct lines of Ca ii HK, Mg i] 4571, [Fe ii] 5250, [O i] 6300, 6364, [Ca ii] 7291, 7323 + [O ii] 7320, 7330, and O i 7774. The expansion velocities of these lines varied, being 10,000 km s⁻¹ for Mg i] 4571 and [Fe ii] 5250, 5000 km s⁻¹ for [O i] 6300, 6364 and [Ca ii] 7291, 7323, and 3000 km s⁻¹ for O i 7774.

A late-time spectrum of PTF12dam (z = 0.107) at +509d post-peak rest frame was presented by Chen et al. (2015). This spectrum had severe contamination by the bright host, and [O i] 6300, 6364 was the only clearly detected feature (although hampered by noise, Mg i] 4571 may also be seen). In contrast to SN 2007bi, the [O i] 6300, 6364 lines appeared to be significantly stronger than the detection limits for [Fe ii] 5250 and [Ca ii] 7291, 7323. Their expansion velocities were also of order 5000 km s⁻¹.

Yan et al. (2015) reported a spectrum of iPTF13ehe (z = 0.3434) at +251d post-peak rest frame. This spectrum was dominated by a broad Hα (4500 km s⁻¹), as well as a prominent He i 5876 or Na i D line. The SN showed a plateau in the r-band from +50d to +200d, suggesting a phase of circumstellar interaction to produce the Hα line.

Lunnan et al. (2016) presented spectra of PS1-14bj (z = 0.5215) at +202d post-peak rest-frame, and of LSQ14an (z = 0.1637) at +205d. They identified broad (2500 km s⁻¹) [O iii] 4363 and [O iii] 4959, 5007 lines in both SNe. Both also showed Mg i] 4571 and hints of an emerging [O i] 6300, 6364.

Models developed to fit the light curves of SLSNe must also be able to reproduce these spectral features seen at late times. The suggested model scenarios for SLSNe are ⁵⁶Ni powering in either pair-instability SNe (PISNe) (Gal-Yam et al. 2009; Kozyreva et al. 2014; Kozyreva & Blinnikov 2015) or core-collapse SNe (CCSNe) (Moriya et al. 2010; Yoshida & Umeda 2011; Yoshida et al. 2014), central engine powering by a magnetar or accreting black hole (Kasen & Bildsten 2010; Woosley 2010; Dexter & Kasen 2013; Metzger et al. 2015), and strong circumstellar interaction (Smith et al. 2007; Chevalier & Irwin 2011; Ginzburg & Balberg 2012; Moriya et al. 2013; Sorokina et al. 2016), possibly between colliding shells ejected in pulsational pair-instability SNe (PPISNe) (Woosley et al. 2007; Chatzopoulos & Wheeler 2012). Many of these model scenarios can produce light curves in good agreement with observations (e.g., Chatzopoulos et al. 2012, 2013; Inserra et al. 2013; Nicholl et al. 2013) but there is often degeneracy due to the many free parameters.

Very few nebular-phase models for SLSNe have so far been calculated. Initial single-zone modeling of SN 2007bi indicated large element masses and the possibility of a pair-instability origin (Gal-Yam et al. 2009). However, calculations of multi-zone PISN explosion models (Dessart et al. 2013; Jerkstrand et al. 2016) showed that PISNe are cold and neutral at nebular times, and produce red spectra with emission lines of species such as Fe i and Si i, in poor agreement with observed spectra of SN 2007bi and PTF12dam. Models for photospheric-phase spectra have been calculated by Dessart et al. (2012, 2013), and Mazzali et al. (2016). These models generally have favored a central-engine powered scenario over a radioactivity-powered one.

In this paper we aim to place further constraints on the nature of long-duration Type Ic SLSNe by both observational and theoretical work on late-time spectra. We present multi-epoch nebular observations of two SLSNe, LSQ14an and SN 2015bn, covering +250d to +410d. We make a new calibration of the nebular spectra of SN 2007bi, carefully considering host galaxy subtraction using post-explosion host photometry. This gives us a high-quality data set of three long-duration SLSNe, which allows an analysis of similarities and differences within the class and with respect to other SNe. We then investigate spectral models of emission by oxygen-zone material to derive constraints from the emission lines seen.

The paper is structured as follows. In Section 2 we present the new data for LSQ14an and SN 2015bn, as well as the new reductions of SN 2007bi. In Section 3 we compare the nebular spectra of these SLSNe with each other and with other non-superluminous SNe Ic. In Section 4 we study spectral formation in oxygen/magnesium material at powering levels relevant for SLSNe, providing a systematic analysis of dependencies on mass, clumping, and energy deposition. In Section 5 we discuss the results and their implications, as well as uncertainties. In Section 6 we present our summary and conclusions.

Our data set comprises two nebular-phase spectra of each of the three supernovae SN 2007bi, LSQ14an, and SN 2015bn. The LSQ14an and SN 2015bn spectra are new observations, whereas the SN 2007bi spectra are re-reductions of those presented in Gal-Yam et al. (2009). Table 1 summarizes the data. Below the observations and reductions are described in detail.

2.1. New Reduction of Nebular SN 2007bi Spectra

SN 2007bi spectra were published by Gal-Yam et al. (2009). The redshift of z = 0.1279 gives a luminosity distance of 583 Mpc with standard cosmology (H₀ = 72 km s⁻¹ Mpc⁻¹, Ω_m = 0.3, Ω_Λ = 0.7), which is used in later conversion to rest-frame spectral luminosities.

2.1.1. SN 2007bi+367d Spectrum

The reduced VLT FORS2 spectrum at +367d post-peak (rest-frame) from Gal-Yam et al. (2009) is available on WISeREP (Yaron & Gal-Yam 2012).¹¹ This spectrum does not have any host galaxy subtracted. Inspection of the 2D spectral images (taken under seeing conditions of 08) and the contemporary ${BV}\,{R}_{{\rm{c}}}{I}_{{\rm{c}}}$ images, shows that the flux of the SN is effectively unresolved from the host, preventing any direct subtraction. Chen et al. (2015) presented uncontaminated, deep images of the host, taken on 2012 May 25, some five years after explosion. Those images show a compact host with a FWHM measured at 125 in the r-band, under seeing conditions of 08.

The griz magnitudes of the host determined by Chen et al. (2015) show that there must be significant host contamination in the VLT spectrum. For example, Chen et al. (2015) measure $g=22.84\pm 0.1$ for the host, which corresponds to $V=22.56\pm 0.15$ (conversion to Vega mags done using the transformations of Jester et al. 2005). Young et al. (2010) measure $V=22.82\pm 0.24$ for SN 2007bi at +367d after the host was subtracted (using a template taken at +669d), demonstrating that the host flux is similar to that of the SN at the VLT epoch in the V-band. In Gal-Yam et al. (2009), the non-host-subtracted spectrum was scaled to the template-subtracted photometry of Young et al. (2010), which gives a too low luminosity in the SN emission lines.

To attempt an improved calibration, we took the FORS2 spectrum which was wavelength and flux calibrated (with telluric lines removed) but had no other corrections applied. It had standard sky background subtraction applied, but no other attempts to remove the host. We obtained this from S. Valenti (2016, private communication). This spectrum is identical to the WISeREP spectrum, apart from the latter having been scaled down by a factor of 1.8. We then scaled the original spectrum to match the ${BV}\,{R}_{{\rm{c}}}{I}_{{\rm{c}}}$ magnitudes in the FORS2 images from the same epoch (containing both SN + host). We calculated these using aperture photometry matched to a number of SDSS DR12 reference stars (converted to Vega mags), finding $B=22.55\pm 0.08$, $V=21.73\pm 0.09$, ${R}_{{\rm{c}}}=21.46\pm 0.17$, ${I}_{{\rm{c}}}=21.26\pm 0.15$, where the errors are the quadrature sum of the statistical error from aperture photometry and the standard deviation of the zeropoints from the reference stars. The photometry should ideally be S-corrected; however this was not attempted here. Compared to the VLT spectrum, the scaling factors for each band were 1.43(B), 1.37(V), 1.35(R_c), and 1.26(I_c) and we used the mean value 1.35 to scale the VLT spectrum to the photometry.

The next step is to subtract the galaxy component. The optimal method would be to use a spectrum of the host galaxy, which is not available. This leaves two options: either to subtract a simple interpolated function between the photometry points, or use a model which is calibrated to the host photometry. We employ here the second method, using a Starburst99 model (Leitherer et al. 1999) with standard settings. The galaxy model was first reddened in the galaxy rest frame with $E(B-V)=0.3$ mag (see more below on this choice) and a Cardelli et al. (1989) law with R_V = 3.1 (which we use throughout the paper), the wavelengths then transformed to observer frame, and SDSS photometry calculated. A scaling factor was determined that minimized the residuals of synthetic SDSS photometry of this spectrum compared to the Chen et al. (2015) host photometry. With just four data points to constrain the galaxy model (the griz magnitudes in Chen et al. 2015), it is not warranted to investigate various parameters in the Starburst99 models. We limit ourselves to the standard model and an age of 30 Myr, and choose $E(B-V)=0.3$ mag which gives a best fit. $E(B-V)\leqslant 0.2$ mag gives a too blue model, whereas $E(B-V)\geqslant 0.4$ mag gives a too red model. The host extinction of the SN line-of-sight is unknown. The Milky Way reddening is $E(B-V)=0.02$ mag (Gal-Yam et al. 2009; Young et al. 2010). The SN only samples one sightline through the host galaxy, and one will therefore in general have two different extinctions for the SN and the "averaged" host. Recently, Leloudas et al. (2015) estimated the host galaxy reddening to be $E(B-V)=0.13$ mag, a factor ∼2 lower than used here. This value gives a too blue model to fit the photometry, possibly as our model is intrinsically too blue for this galaxy. It is however preferential to subtract a model that fits all the photometry than to have the most accurate value for the extinction.

The final SN spectrum was obtained by subtracting the host photometry-calibrated galaxy model from the SN+host-photometry-calibrated VLT spectrum. The resulting spectrum is shown in Figure 1 (bottom). The narrow host lines result from UV reprocessing, an effect not included in the Starburst99 models. It is clear from the figure that they contribute little to the host photometry. We estimate a 5% impact from synthetic photometry on exctracted lines, and make this small correction to the fit. For LSQ14an (Section 2.2), this effect is stronger, of order 30%.

Zoom In Zoom Out Reset image size

2.1.2. SN 2007bi+471d Spectrum

There is another spectrum of SN 2007bi, taken with Keck 104 days later (rest-frame) and presented in Gal-Yam et al. (2009) (also available on WISeREP). There is no contemporary photometry to do flux calibration against, and Gal-Yam et al. (2009) used spectral flux calibration using a standard taken 15 minutes before. However, synthetic photometry of this host+SN spectrum yields r = 23.3 mag, 0.9 mag fainter than even the host (${r}_{\mathrm{host}}=22.4$ mag). This means the flux calibration is likely off by a factor of several, and the night was indeed non-photometric (M. Sullivan 2016, private communication), possibly explaining the discrepancy.

The image quality in the 2D frames (from the Keck archive) is 12, meaning that no host galaxy subtraction is possible in the 2D extraction. We can apply the model subtraction method of the preceding section to this spectrum as well, albeit with more uncertainty than for the VLT FORS2 spectrum. We assume that SN 2007bi continued to decline linearly in the R_c-band with the slope given in Gal-Yam et al. (2009) and Young et al. (2010), giving an expected magnitude of ${R}_{c}=23.2\pm 0.2$ at the epoch of the Keck spectrum (this is a simple mean of the Gal-Yam et al. and Young et al. extrapolations). Thus the expected combined magnitude of host and SN is ${R}_{c}=21.9\pm 0.2$. To match the Keck spectrum to this, it needs to be scaled by a factor 3.1. We then subtracted off the same Starburst99 model as used for the VLT FORS2 spectrum correction. We note that the factor 3.1 scaling gives a consistent position of the galaxy component with the VLT spectrum, which should dominate at wavelengths with weak SN lines. If we had chosen an approach to ignore the slope of previous R points and just tried to align the galaxy flux, we would have obtained a very similar scaling.

2.1.3. Final Spectra

Figure 1 shows the calibration process and final spectra. The top panel shows the raw VLT and Keck spectra scaled to SN+host photometry. Also shown is the host photometry which demonstrates that host light makes up about half the flux in the VLT spectrum and about 3/4 in the Keck spectrum.

Synthetic photometry on the galaxy-subtracted VLT spectrum agrees well with the template-subtracted photometry, with errors of +0.31, −0.02, and −0.17 mags in V, ${R}_{{\rm{c}}}$, and ${I}_{{\rm{c}}}$ (Young et al. 2010 did not present a template-subtracted value in B), and means we have applied a self-consistent method to remove the host flux and match the image-subtracted SN photometry.

The final spectra are shown in the bottom panel. Compared to the VLT spectrum presented by Gal-Yam et al. (2009), the rising blue flux is removed and the emission lines are significantly stronger. The spectrum appears more similar to typical Type Ib/c nebular SN spectra. The corrected Keck spectrum is similar to the VLT spectrum, but dimmer and with significantly lower signal to noise. Our final estimates for absolute R-band magnitudes are ${M}_{R}=-16.72$ at +367d and ${M}_{R}=-16.02$ at +471d. The optical flux has declined by a factor 3.0 over 104d, close to the ⁵⁶Co decay (2.55). The Mg i] 4571 line and the Mg i 5180 + [Fe ii] 5250 feature are noticeably weaker in the Keck spectrum. The redder emission lines of [O i] 6300, 6364 and [Ca ii] 7291, 7323 have also weakened, but not as much.

2.2. LSQ14an

LSQ14an was discovered by the La Silla Quest Survey and spectroscopically classified as a SN 2007bi-like Type Ic SN at z = 0.163 by the Public ESO Spectroscopic Survey of Transient Objects (PESSTO¹² ). The classification spectrum of Leget et al. (2014) was taken on 2014 January 2 (MJD 56659) and suggested an epoch of ∼+55d post peak (rest-frame) from comparison with SN 2007bi and PTF12dam. This puts the peak around 2013 October 30 (MJD 56595). If we assume a rise time of 70d, the estimated explosion epoch is ∼2013 August 9 (MJD 56513). With our higher signal-to-noise spectra we refined the redshift estimate to $z=0.1637\pm 0.0001$, corresponding to a luminosity distance of 763 Mpc with standard cosmology.

We observed LSQ14an on two occasions in its nebular phase, with VLT X-shooter at +365d post-peak rest frame, and with VLT FORS2 at +410d. The observations and data reductions are detailed in Appendix A. The full PESSTO data set on the early evolution of LSQ14an is presented by C. Inserra et al. (2016, in preparation).

The spectra are shown in Figure 2. There is strong contamination by the host galaxy lines and continuum, and as for SN 2007bi, the SN is unresolved from the host flux. Hence we need to employ a similar method of modeling the host flux and subtracting it off both the X-shooter and FORS2 spectra to leave a clean SN spectrum.

Zoom In Zoom Out Reset image size

The host galaxy was detected in the Pan-STARRS1 3π survey (cataloged object PSO J193.4492-29.5243). Images in each of the filters ${{grizy}}_{{\rm{P}}1}$ were taken from the PS1 Processing Version 2 stacks. All the individual images which make up these stacked images were taken before MJD 56467 (2013 June 24), 1–2 months before the estimated explosion date of early 2013 August. Hence there is unlikely to be any SN flux in the images. Magnitude measurements of the host were carried out using aperture photometry. We let the aperture size vary until we were confident that it encompassed the whole host flux and avoided covering nearby object flux. The aperture radius we adopted was 2 arcsec. The zeropoint was determined with 15–18 PS1 reference stars in the field (Schlafly et al. 2012; Magnier et al. 2013) resulting in values of ${g}_{{\rm{P}}1}=20.93\pm 0.11$, ${r}_{{\rm{P}}1}=20.50\pm 0.12$, ${i}_{{\rm{P}}1}=20.56\pm 0.10$, ${z}_{{\rm{P}}1}=20.42\pm 0.14$, and ${y}_{{\rm{P}}1}=20.43\pm 0.19$ (AB mags in the filter system described in Tonry et al. (2012)). The position of LSQ14an was measured at R.A. = 193.44924, decl. = −29.52434 in the astrometrically calibrated EFOSC2 images (with SExtractor). In comparison the centroid of the host galaxy was measured at R.A. = 193.44921, decl. = −29.52429. The difference between the two is 0.2 arcsec, which illustrates that it is coincident with (and unresolved from) the dwarf galaxy host.

Figure 2 shows also the host photometry. There is clearly significant contamination by galaxy light. However, distinct broad emission lines are seen on top of this. The situation is somewhat better than initial inspection suggests, because the narrow galaxy lines have an impact on the host photometry.

Shown also is the 30 Myr Starburst99 model, reddened with $E(B-V)=0.1$ mag (in agreement with $E(B-V)$ estimates from the Balmer line ratio (C. Inserra et al. 2016, in preparation)). and scaled to give an average flux level of 0.7 times the observed photometry; this is to compensate for the narrow emission lines which are missing in the model. A factor smaller than 0.7 is not warranted as the lines do not add more than ∼1/3 to the photometry. A larger factor fails to give out much positive SN flux as the host photometry is even somewhat brighter than the host+SN photometry, presumably due to observational or calibration issues. One may attempt to estimate the contribution by the narrow lines more precisely, but this would not help too much here because of the discrepancy that host photometry is somewhat brighter than the SN+host photometry.

Figure 3 shows the host-subtracted optical SN spectra. These spectra also have had the narrow galaxy lines clipped out (see Appendix A). The optical flux decreases by a factor 2 over 45d, or 1.6 mag per 100d. The bottom panel of Figure 6 shows the combined optical and near-infrared spectrum. The NIR flux is weak compared to the optical. The NIR spectrum shows no clearly distinguished lines, apart from possibly O i 1.13 μm.

Zoom In Zoom Out Reset image size

2.3. SN 2015bn

SN 2015bn is alternatively known as PS15ae, CSS141223-113342+004332 and MLS150211-113342+004333. The earliest detected point is on 2014 December 23 by the Catalina Sky Survey and the evolution was then followed by PESSTO after the initial classification as a SLSN in Le Guillou et al. (2015). The SN was monitored from −50d to +250d (with respect to optical maximum) showing it to be a slowly evolving SLSN similar in observational characteristics to SN 2007bi and LSQ14an (Nicholl et al. 2016b). Spectropolarimetry was presented by Inserra et al. (2016a). The redshift is z = 0.1136 from the narrow host galaxy emission lines, giving a luminosity distance of 513 Mpc with standard cosmology. The light curve peaked around 2015 March 27 (MJD 56108) so the SN had a rise time of ≳80d rest-frame.

We observed SN 2015bn at two epochs in its nebular phase, at +250d post-peak rest frame with NTT EFOSC2 and at +315d with VLT X-shooter. The observations and data reductions are detailed in Appendix B. Figure 4 shows the spectra, the host photometry (Nicholl et al. 2016b), and the same Starburst99 model as in the previous section (using $E(B-V)=0$) scaled to this photometry. Figure 5 shows the host-subtracted spectra at rest-frame wavelengths (note that the +250d spectrum has been scaled by 0.2). From +250d to +315d, a more distinct [O i] 6300, 6364 develops in the spectrum, as well as higher contrast for Ca ii HK, Mg i] 4571, and [Fe ii] 5250. We see the same feature at 5500 Å as in SN 2007bi, probably [O i] 5577. The optical flux has decreased by a factor 5 over 65d, or 2.2 mag per 100d. At the end of the red arm we see the first detection of Ca ii NIR in SLSNe at nebular times. This line is strong and broad, with an expansion velocity of 8000–10,000 km s⁻¹.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The top panel of Figure 6 shows the combined optical and NIR spectrum. The NIR flux is weak compared to the optical. There are clear detections of O i 9263, O i 1.13 μm, and Mg i 1.50 μm. These lines are commonly seen in stripped-envelope SNe with nebular NIR spectra, e.g., for the Type IIb SN 2008ax (Taubenberger et al. 2011), SN 2011dh (Ergon et al. 2015; Jerkstrand et al. 2015), the Type Ic SN 2007gr (Hunter et al. 2009), and in the broad-lined SN 2012au (Milisavljevic et al. 2013).

Zoom In Zoom Out Reset image size

A comparison with the He100 PISN model of Jerkstrand et al. (2016) in the same figure shows no strong emission in the [Si i] + [S i] blend around 1.08 μm predicted to be strong. Also [Si i] 1.607, 1.645 μm are predicted to be strong, but this feature falls in the telluric band.

In this section we compare the spectra of SN 2007bi, LSQ14an and SN 2015bn with each other as well as with other SNe. For the three SLSNe we attempt no dereddening as the host line-of-sight extinctions are not known, and the Milky Way extinctions are smaller ($E(B-V)=0.03$ mag for SN 2007bi (Gal-Yam et al. 2009), $E(B-V)=0.07$ mag for LSQ14an (NED¹³ ) and $E(B-V)=0.02$ mag for SN 2015bn (Nicholl et al. 2016b)) than the uncertainties in flux calibration and host galaxy subtractions.

Figure 7 compares the spectra of SN 2007bi, LSQ14an, and SN 2015bn at around one year post-peak in the rest frame. All three SNe show qualitatively similar spectra, with Ca ii HK, Mg i] 4571, Mg i 5180 + [Fe ii] 5250, [O i] 6300, 6364, [Ca ii] 7291, 7323 and O i 7774 all distinct. There are also clear detections of [O i] 5577 and Na i D in SN 2015bn and SN 2007bi and hints of [O i] 5577 in LSQ14an. The similarity in luminosity of both [O i] 6300, 6364 and O i 7774 in all three SNe is striking. LSQ14an has, in addition, strong [O iii] 4363, [O iii] 4959, 5007, and the prominent feature usually identified as the forbidden [Ca ii] 7291, 7323 doublet. Given the strength of the O iii lines, we suggest this likely is dominated by [O ii] 7320, 7330 in LSQ14an (see also Section 4.4). Figure 15 in Appendix C shows comparisons of each of the main emission lines in velocity space.

Zoom In Zoom Out Reset image size

In Figure 8 we compare the spectra of SN 2007bi and SN 2015bn with SN 1998bw spectra at +140d, +215d, and +391d post-explosion (Patat et al. 2001), retrieved from WISeREP. For the comparisons we use z = 0.0085 and $E(B-V)=0.06$ mag for SN 1998bw (Patat et al. 2001). SN 1998bw has a representative spectrum for the class of broad-lined SNe Ic, which have relatively similar nebular phase spectra (e.g., Taubenberger et al. 2006; Milisavljevic et al. 2015; Modjaz et al. 2016). At 200–300d after explosion, SN 1998bw, SN 2002ap (Foley et al. 2003), SN 2003jd (Valenti et al. 2008), SN 2004aw (Taubenberger et al. 2006), and SN 2012ap (Milisavljevic et al. 2015) all show [O i] 6300, 6364 as the strongest line, calcium lines with a Ca ii NIR/[Ca ii] 7291, 7323 ratio of <1, and distinct Mg i] 4571. Some show a distinct line at 5250 Å (where all SLSNe have a line), and some do not.

Zoom In Zoom Out Reset image size

The figure shows that there is strong overall similarity between the SLSNe and SN 1998bw, but with distinct phase shifts. The spectra of SN 1998bw at +140d and +215d show similar emission lines and relative strengths as the spectra of SN 2007bi and SN 2015bn at +315d to +367d after peak. Thus, SLSNe appear to behave like broad-lined SNe Ic but to evolve more slowly. It is particularly noteworthy that SN 1998bw also displays distinct Mg i] 4571 and Mg i 5180+ [Fe ii] 5250. The relative strength of these lines compared to the quasi-continuum is somewhat lower in SN 1998bw than in the SLSNe. Comparing with SN 1998bw at +391d shows that the blue quasi-continuum and the 5250 Å feature have faded compared to the red lines, and the similarity at blue wavelengths is reduced. The [Ca ii] 7291, 7323 lines are, however, more similar in shape at this epoch. At earlier epochs, SN 1998bw has a line at ∼7100 Å not seen in the SLSNe. This line may be [Fe ii] 7155 or possibly He i 7065.

At no epoch does SN 1998bw show a particularly strong O i 7774 line, which is different from the SLSNe. Neither do most other broad-lined Ic such as SN 1997ef (Matheson et al. 2001), SN 2002ap, SN 2003jd, SN 2006aj (Maeda et al. 2007), SN 2009bb (Pignata et al. 2011), or SN 2012ap. There are two known exceptions; SN 1997dq, which showed a strong narrow O i 7774 line at +262d post explosion (Matheson et al. 2001; Mazzali et al. 2004), and SN 2012au at +321d (Milisavljevic et al. 2013). Milisavljevic et al. (2013) showed an overall strong similarity between the nebular spectra of SN 2012au and SN 2007bi. Also, O i 9263 is weaker in the broad-lined Ic SNe.

Another difference is [O i] 5577, which is distinct in the SLSNe but not seen in SN 1998bw. This probably indicates that the temperature is higher in the SLSNe, since the [O i] 5577 line has a high excitation energy.

SN 2015bn has additional coverage of Ca ii NIR and O i 9263. The Ca ii NIR/ [Ca ii] ratio in SN 2015bn is higher than in SN 1998bw or any other Type Ic SN at any epoch. We study the implications of this in Section 4.7.

The relative homogeneity in spectral appearance established in Section 3 warrants an investigation on what physical conditions are needed to produce the observed line luminosities in O, Mg, and Ca lines. Here, we pursue single-zone modeling of oxygen-rich compositions, exploring parameter space over mass, density, and powering level. The aim is to understand the physical mechanisms behind the emission of distinct lines (e.g., [O i] 6300, 6364, [O i] 5577, O i 7774, Mg i] 4571, Mg i 5180, Mg i 1.50 μm, [Ca ii] 7291, 7323, and Ca ii NIR) and the diagnostic use of these. Generalized single-zone modeling is a complementary approach to the calculation of multi-zone explosion models, where constraints may be derived independent of the powering mechanism and parameter space can be explored to a larger extent. Indications of ejecta properties obtained from it can be used to select explosion models to test in follow-up multi-zone modeling.

For the modeling we use both the SUMO spectral synthesis code (Jerkstrand et al. 2011, 2012) and semi-analytic approaches. The SUMO code solves the equations for non-thermal energy degradation, non-local-thermodynamic-equilibrium excitation and ionization balance, temperature, and radiative transfer. It provides a self-consistently calculated spectrum of any ejecta structure, without any free parameters.

Because the forward modeling is always limited to particular compositions and ejecta structures, and involves much complex physics, some of it treated in rough approximations, it is beneficial to complement it with backward modeling using analytic formulae, which we also pursue here.

4.1. SUMO Grid Setup

Each model in the SUMO grid was set up as N = 100 randomly distributed spherical clumps within a volume specified by an expansion velocity V and a time t (assuming homology). Together with a specified total mass M_ej, the density of the clumps is determined by a filling factor f, with vacuum in the $1-f$ space. The radiative transport is computed self-consistently throughout the volume. We compute a model grid with

• 1.  
  M_ej = 3, 10, and 30 M_⊙
• 2.   
  Energy deposition (by gamma-rays or any other high-energy source) $2.0\times {10}^{41}$, $2.5\times {10}^{41}$, $5.0\times {10}^{41}$, $10\times {10}^{41}$, and $20\times {10}^{41}$ erg s⁻¹.
• 3.  
  Composition : pure O, O, and Mg (92% and 8% by mass), and full C-burning ashes (see below). The pure O and OMg compositions were computed for selected values of other parameters only.
• 4.   
  Filling factors : $f=0.001,0.01$, and 0.1.
• 5.  
  V = 8000 km s⁻¹(the characteristic observed line width.)
• 6.   
  Time = 400d since explosion

The models are named by composition, mass, filling factor, and energy deposition; for example O30-0.1-2.0 is a pure oxygen gas of mass 30 M_⊙, with f = 0.1 and deposition $2.0\times {10}^{41}$ erg s⁻¹.

By allowing for a free energy deposition we decouple the problem of powering and spectral formation. For energy input by high-energy particles, the model spectra are insensitive to the properties of the input spectrum as long as the primaries have an energy E ≳ 1 keV (Axelrod 1980; Kozma & Fransson 1992). Thus, the results will be valid for any high-energy input such as gamma-rays from radioactive decay or hard X-rays from a central engine or from shocks.

The choice of C-burning ashes to represent the O-zone composition needs some motivation. In SN explosion models there are three principal oxygen zones: He-burning ash (C/O dominated), C-burning ash (O/Ne/Mg dominated), and Ne-burning ash (O/Si/S dominated). The O/Ne/Mg zone is, in most CCSN models, the most massive of these, in particular at high ${M}_{\mathrm{ZAMS}}$ values (Woosley & Heger 2007). Also in PISN models O/Ne/Mg is the most common composition (Heger & Woosley 2002). In addition, the O/C and O/Si/S zones are prone to formation of CO and SiO molecules at nebular times (Liu & Dalgarno 1995), which can lead to strong cooling and blackout of optical thermal emission. Therefore, C-burning composition is most motivated for single-zone modeling. Clearly, the strong Mg i] 4571 line in SLSNe is also stemming from this zone.

The next question is which composition to use for the C-burning ashes. The composition in nucleosynthesis models depends on the particular progenitor, but the relative abundances of O, Ne, Mg, and Na are relatively insensitive to ${M}_{\mathrm{ZAMS}}$. The Mg/O mass ratio is typically around 0.09, which we use in the pure OMg models (compare to the solar ratio of 0.12). For the C-burning models, we use the composition from the ONeMg zone in the ${M}_{\mathrm{ZAMS}}=25$ M_⊙ solar metallicity model of Woosley & Heger (2007), listed in Table 2. Elements heavier than calcium stem mainly from primordial (solar) abundances, although some modification is made by processes such as n and p capture. Because SLSNe have been shown to arise in regions of subsolar metallicity (Chen et al. 2013, 2015, 2016b; Lunnan et al. 2014; Leloudas et al. 2015; Perley et al. 2016), a reduction of iron-group abundance may be warranted. However, we did some experiments comparing models with 10 times lower metallicity and found only moderate differences.

Table 2.  Composition (Mass Fractions) of C-burning Ashes. The Abundances are Taken from the ${M}_{\mathrm{ZAMS}}=25$ M_⊙ Model of Woosley & Heger (2007)

H	0	S	$6.3\times {10}^{-4}$
He	$1.4\times {10}^{-6}$	Ar	$9.9\times {10}^{-5}$
C	$7.9\times {10}^{-3}$	Ca	$3.1\times {10}^{-5}$
N	$3.5\times {10}^{-5}$	Sc	$1.6\times {10}^{-6}$
O	0.74	Ti	$8.0\times {10}^{-6}$
Ne	0.15	V	$6.1\times {10}^{-7}$
Na	$1.8\times {10}^{-3}$	Mn	$2.8\times {10}^{-6}$
Mg	0.072	Fe	$5.6\times {10}^{-4}$
Al	$6.7\times {10}^{-3}$	Ni	$6.6\times {10}^{-4}$
Si	0.014	Co	$1.3\times {10}^{-4}$

Download table as:  ASCIITypeset image

All models were computed without charge transfer, as these reactions gave numeric instabilities for some high-density models, and to allow comparison we prefer to use the same physics for all models.

4.2. Varying Composition and Density

In Figure 9 we illustrate the influence of composition and density, comparing ${M}_{{ej}}=30$ M_⊙ models with different compositions and filling factors. For each composition and f-value, we plot the model with energy deposition giving the best reproduction of [O i] 6300, 6364. The model spectra are compared with SN 2007bi at +367d. One should note that with our single-velocity approach we cannot directly address the issue of different line widths, but limit our discussion to line luminosities. Because we are only studying the emission from one of presumably many composition zones in the SN, the goal is not to reproduce the full spectrum, but rather obtain an idea of what this zone contributes. Since the strongest and most prominent lines in the spectra are due to oxygen and magnesium (and calcium) significant parts of the spectrum should be recovered by an OMg composition. We review the main physical processes that produce the model spectra in the following sections.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

4.2.1. Pure Oxygen Composition

At f = 0.1 (model O30-0.1-2.0, Panel A in Figure 9), there is little emission apart from [O i] 6300, 6364. The temperature (T = 4400 K) is too low for [O i] 5577. There is some power in O i 1350, O i 1640, O i 7774, O i 8446 and O i 9263, arising mainly from recombination. [O i] 6300, 6364 and [O i] 5577 are optically thin ($\tau \sim 0.1$) but the recombination lines are optically thick.

O i 7774 is weak compared to the observed line. While about 30% of the deposited energy goes into ionizing O i, the maximum energy fraction that can emerge in O i 7774 is $0.30\times 1064/7774=0.04$ (assuming case B recombination so the effective ionization threshold is 1064 Å for 2p(¹D)). This is small compared with the cooling in [O i] 6300, 6364 which is about 50% of the total deposition. Thermal population of the O i 7774 parent multiplet 3p(⁵P) from the ground state is not possible (the energy gap is 10.8 eV), and any such very hot situation would in any case lead to a very strong [O i] 5577 line (excitation energy of 4.2 eV), which is observed to be weaker than O i 7774. Only 2% of the energy goes to non-thermal excitations in O i, ruling out this mechanism for powering O i 7774 or any other strong O i line. The NIR spectra of SN 2015bn and LSQ14an show a strong O i 9263 line (Figure 6), which is expected as this is also a recombination line with similar expected strength as O i 7774. This is also underproduced in the f = 0.1 model.

Lowering f to 0.01 (model O30-0.01-2.0, panel B) leads to stronger [O i] 5577 because the temperature increases to T = 4820 K (from T = 4400 K). Collisional excitation of this line occurs mainly from the first excited state, 2p(¹D) (i.e., in the [O i] 5577 transition itself), and the fact that there are a factor 1.6 more O i atoms in the 2p(¹D) state in this model also boosts [O i] 5577. The recombination lines strengthen somewhat as the higher clumping leads to lower ionization (electron fraction x_e = 0.053 compared to x_e = 0.13), which reduces the heating fraction and increases the ionization fraction from 0.36 to 0.44. In addition, photoionization begins to contribute. The recombination lines are, however, still much weaker than [O i] 6300, 6364. [O i] 6300, 6364, and [O i] 5577 are marginally optically thick (τ = 1–3), whereas O i 7774, O i 8446 and O i 9263 have $\tau \gg 1$.

Lowering f to 0.001 (model O30-0.001-2.0, panel C), [O i] 5577 strengthens further as the temperature continues to increase (T = 6550 K). The optical recombination lines strengthen again as the ionization fraction goes to 0.49 at x_e = 0.025, and photoionization now dominates the non-thermal rate, effectively transferring energy from the UV lines and continua to the optical recombination lines. In addition, collisional excitation now contributes a similar amount as recombinations to the O i 7774 luminosity. These excitations occur mainly from 3s(⁵S) which is effectively meta-stable having optically thick transitions to the ground state. That O i 7774 can become a cooling line at high density was pointed out by Maurer & Mazzali (2010), and the line holds promise as a tracer of high density.

In this model, all lines seen have $\tau \gg 1$. The model has quite significant O i 7774, O i 8446, and O i 9263. The UV lines are now quenched due to absorption by lines and continua in the high-density clumps. Note some emergence of [O ii] 7320, 7330,¹⁴ which requires high temperature. But overall, none of these 30 M_⊙ ejecta, constrained to fit the [O i] 6300, 6364 lines, produces any significant [O ii] or [O iii] lines.

By contrasting the three different f values, we see that increasing the clumping increases [O i] 5577, O i 7774, O i 8446, and O i 9263 at the expense of [O i] 6300, 6364. The O i 7774, O i 8446, and O i 9263 lines gain from two effects at higher density: increased fraction of non-thermal ionization and photoionization, and more efficient cooling from excited states as their relative populations increase. The observed relation between these lines is better reproduced at high clumping.

4.2.2. O/Mg Composition

We now explore models containing magnesium. Magnesium typically becomes much more ionized than oxygen because its low ionization thresholds at 1620 Å (ground state) and 2510 Å (first excited state) make it prone to photoionization (e.g., Fransson & Chevalier 1989). As discussed in Jerkstrand et al. (2015), the exact solution for the neutral fraction determines whether Mg i] 4571 can contribute to the cooling or not. A high luminosity in this line as observed in SLSNe almost certainly requires cooling and not just recombination. How this is achieved is not easily understood, since Mg i is easy to ionize when the energy deposition is high. The neutral fraction is typically increased by higher density, and this holds potential as a constraint to be explored.

Introducing magnesium at f = 0.1 (model OMg30-0.1-2.0, panel D in Figure 9) gives a temperature T = 4330 K (compared to $T=\ 4400$ K for pure O composition), electron fraction ${x}_{e}=0.15(0.13)$, and neutral magnesium fraction $x(\mathrm{Mg}\,{\rm{I}})=0.004$, as expected much lower than the oxygen neutral fraction (0.90).

Mg ii 2795, 2802 and Mg i] 4571 are distinct. Neutral magnesium, although rare, does about 10% of the cooling. This makes a distinct line, although weaker than observed. The recombination emission in Mg i] 4571 is moderate, and the cooling dominates the line. Note that 10% cooling corresponds to dep × f_heating × f_cool,line = ∼2 × 10⁴¹ × 0.6 × 0.10 = 10⁴⁰ erg s⁻¹, which with a typical line width of 100 Å is 10³⁸ erg s⁻¹ Å⁻¹, close to what is seen in the model. Mg ii 2795, 2802 does 12% of the cooling. Some of this emission fluoresces in O and Mg lines, but a large part escapes.

At f = 0.01, (model OMg30-0.01-2.0, panel E), the temperature is T = 4530 K (compare T = 4820 K pure O composition), ${x}_{e}=0.076(0.053)$, and $x(\mathrm{Mg}\,{\rm{I}})=0.072$. Higher clumping boosts Mg i] 4571, because Mg has now crossed the "ionization runaway" threshold (see Appendix E) and is neutral to 7.2%, doing significant cooling (25%). The line is in LTE and optically thick ($\tau =900$), as such its luminosity can be derived from the volume ${ \mathcal V }$ of the emitting region and the temperature (see e.g., A. Jerkstrand 2016, in preparation¹⁵ ):

$$\begin{eqnarray}&&L=\displaystyle \frac{4\pi { \mathcal V }}{{ct}}\lambda {B}_{\lambda }(T).\end{eqnarray} \tag{ 1 }$$

Once in this regime, the line luminosity cannot easily be increased further; the temperature would need to increase enough with compression so ${B}_{\lambda }(T)$ grows faster than $1/{ \mathcal V }$. This will certainly not happen for lines with $\lambda \gtrsim {hc}/(4.96{kT})$ (peak of the blackbody), and also shorter wavelength lines typically show decline or modest growth. From the perspective of the single transition, increasing the Mg i fraction leads to increased optical depths, which offsets the increased number of atoms available. As a rough picture to understand the 25% cooling in Mg i] 4571, most of the 12% cooling done by Mg ii in the previous model has now been transferred to Mg i] 4571. The Mg ii lines are in the LTE and optically thick regime already at f = 0.1, and thus from Equation (1), reducing ${ \mathcal V }$ by a factor of 10 while obtaining only a small change in T, their luminosity drops by a factor ∼10. Mg i] 4571, on the other hand, is in NLTE at f = 0.1, and can increase its luminosity until it saturates in LTE, compensating for the loss of Mg ii cooling.

The Mg i 5180 multiplet has now also emerged, being formed mainly by recombinations but also by some collisional excitation from 3p(³P). Because most Mg is still singly ionized, ${ \mathcal V }\times {n}_{\mathrm{Mg}{\rm{II}}}$ is constant and the Mg recombination line luminosities increase as ${ \mathcal V }{n}_{e}{n}_{\mathrm{Mg}{\rm{II}}}\propto {n}_{e}$. The Mg i 5180 lines are optically thick ($\tau \sim {10}^{5}$), but the parent population is below LTE.

At f = 0.001 and an energy deposition of $2\times {10}^{41}$ erg s⁻¹ (as used in all models up until now), a much too dim optical spectrum (not shown) is produced, as Mg NIR lines dominate the cooling. To fit the [O i] 6300, 6364 lines, we need instead a deposition of $1\times {10}^{42}$ erg s⁻¹(model OMg30-0.001-10, panel F). This model has temperature T = 6470 K, ionization x_e = 0.11, and Mg neutral fraction $x(\mathrm{Mg}\,{\rm{I}})=0.001$. The higher energy deposition has offset the higher clumping to reduce the Mg i abundance again. But the Mg i] 4571 line is still in the LTE, optically thick regime, ${B}_{4571}(T){ \mathcal V }$ is roughly the same as in the f = 0.01 model, and the line luminosity therefore changes little. All the lines of Mg i] 4571, Mg i 5180, [O i] 5577, [O i] 6300, 6364 are now close to the LTE and optically thick regime; their peak flux therefore follows a blackbody at the zone temperature of 6700 K, per Equation (1) (${F}_{\lambda }\propto L/\lambda \propto {B}_{\lambda }$, since ${\rm{\Delta }}\lambda \propto \lambda$). A multiplet with N lines, all in LTE and optically thick, would give a line feature N times above the blackbody curve. The O i 7774, 8446, and 9263 multiplets tend to have departure coefficients still a factor 2–3 below unity, which offsets this multiplicity effect and gives total line luminosities close to the blackbody curve. [O i] 6300, 6364 has more radiative transfer losses than the other lines, and emerges somewhat below its expected flux.

This model has a good overall optical brightness compared to the observed spectra, and reproduces Mg i] 4571, Mg i 5180 (although somewhat too blue compared to the observed line), and [O i] 6300, 6364. It produces also a strong O i 9263, in good agreement with the observed line in SN 2015bn (Figure 7). It has too strong [O i] 5577 and O i 7774. Mg i 3850 and O i 8446 will scatter in Ca ii HK and Ca ii NIR, respectively, in models including calcium, so their overproduction is not a big concern.

A hot model like this ubiquitously produces a strong [O i] 5577, which is not distinct in SN 2007bi. There is an unusual "cut-off" in the observed spectra around 5600 Å that may correspond to a line absorption threshold (e.g., in He i 5876 or Na i D), but it appears that a [O i] 5577 as strong as in this model should still reveal itself by its blue edge.

This is an interesting model in the sense that it shows that a pure O/Mg gas can qualitatively reproduce many of the observed spectral features in long-duration SLSNe. However, several key discrepancies remain: blueshift of the 5200 Å feature, and too strong [O i] 5577 and O i 7774. It is unclear whether the observed 5200 Å feature is dominated by Mg i 5180, [Fe ii] 5250, or is a blend. A Mg i 5180 component is plausible because it emerges in the models where Mg i] 4571 approaches its observed strength.

There is another constraint on the emission from magnesium from Mg i 1.50 μm. Figure 11 compares this model to the X-shooter spectrum of SN 2015bn which shows this line. The model reproduces O i 9263 and the Mg i 1.50 μm line, and overall the agreement across the NIR region is reasonably good. The data for LSQ14an are noisier but a feature at the 1.50 μm position may also be present (Figure 6).

Zoom In Zoom Out Reset image size

In summary, adding Mg to the composition (in its typical C-burning abundance) shows that Mg i] 4571 and Mg i 1.50 μm emerge naturally, and consolidates the identification of these lines. Mg i 5180 emerges at low f values, and likely contributes to the 5200 Å feature. The best reproductions of both O and Mg line ratios are at low filling factors. At high f values magnesium tends to be too ionized to make strong Mg i cooling lines, and the Mg i] 4571 luminosity indicates $f\leqslant 0.01$.

4.2.3. C-Burning Composition

Introducing the full composition at f = 0.1 (model C-burn30-0.1-2.0, panel G in Figure 9) gives a temperature T = 4350 K (4330 K for OMg), x_e = 0.14 (0.15), and $x(\mathrm{Mg}\,{\rm{I}})=0.003$ (0.004). Thus, trace elements have a small effect on temperature and ionization at this filling factor. They provide about 10% of the line cooling. They have, however, led to the suppression of O i 1640 and Mg ii 2795, 2802, and the formation of a more pronounced quasi-continuum (made up of many overlapping weak lines). This quasi-continuum is flat and is weak compared to the [O i] and Mg i lines, and does not resemble the 4000–5500 Å component that is often seen in stripped-envelope SNe and may be present in the SLSNe as well. As the OMg model at this f‐value, this model has too dim Mg i lines. The presence of trace metals has weakened Mg i] 4571.

At f = 0.01 (model Cburn30-0.01-2.0, panel H), the temperature is T = 3950 K (4530 K for OMg), x_e = 0.078 (0.076), and $x(\mathrm{Mg}\,{\rm{I}})=0.14$ (0.072). Trace elements have a somewhat larger effect on temperature and ionization at this filling factor. They now provide 70% of the cooling, mainly through Fe i (iron is 75% neutral). Mg i] 4571 is significantly damped. Despite the high Mg i fraction, Mg i] 4571 cannot compete with iron for the cooling, and cools only 5%. The line peaking at 5150 Å is dominated by Fe i, and not Mg i 5180. At this higher density, Na i D has emerged. The neutral fraction of sodium is $x(\mathrm{Na}\,{\rm{I}})=5\times {10}^{-3}$, and Na i provides 5% of the cooling, which makes the luminosity in good agreement with the observed one.

At f = 0.001 (model Cburn30-0.001-10, panel I), the temperature is T = 5490 K (6470 K for OMg), x_e = 0.082 (0.11), and $x(\mathrm{Mg}\,{\rm{I}})=0.002$ (0.001). Trace elements have now a significant effect on the temperature and ionization; they provide 85% of the line cooling. Mg i] 4571 cools only a few percent. At this high density there are also several optically thick iron-group lines in the vicinity of Mg i] 4571, and the line is further quenched by line blocking, as well as blended into the surrounding quasi-continuum. In fact, there are thousands of optically thick lines throughout the optical (due to the primordial iron-group abundances) and a pronounced quasi-continuum is formed; we are approaching a blackbody SED.

In summary, while there are only moderate differences between OMg and C-burn compositions at f = 0.1, the large differences at $f\leqslant 0.01$ shows that at high density, the models are sensitive to the presence of trace elements and thereby to the metallicity of the progenitor. The metals tend to damp the contrast between distinct lines and quasi-continuum, and makes in particular Mg i] 4571 less pronounced.

4.3. Varying Mass

We now consider the models at lower masses, ${M}_{{ej}}=3$ and 10 M_⊙. Each of these models uses the full C-burning composition, but vary in filling factor and energy deposition. Figure 10 shows a subset of nine of these models, where for each M (horisontal variation) and f (vertical variation) combination, the model with energy deposition giving best overall reproduction of O and Mg lines to SN 2007bi is shown.

4.3.1. 3 M_⊙ Models

The 3 M_⊙ models are shown in the left column of Figure 10 (note that the deposition energy is different for each f value). For no deposition value do the 3 M_⊙ models produce strong enough [O i] 6300, 6364 luminosity; increasing the deposition beyond a certain point reduces the [O i] luminosity as O i becomes ionized to O ii (see also Section 4.3.4 and Figure 12).

Zoom In Zoom Out Reset image size

The low mass leads to high temperature and ionization, with the spectra having a strong Mg ii 2795, 2802 and a strong blue quasi-continuum between 1000 and 4000 Å. In a multi-component model this blue flux would probably be reprocessed to longer wavelengths by scattering and fluorescence, and thus may not be relevant to compare with observations. The optical spectrum is at f = 0.1 and f = 0.01 dominated by [O ii] 7320, 7330 and [O iii] 4959, 5007. The O iii fraction is 7%, 2%, and 0.5% in the best-fitting models at $f=0.1,0.01$ and 0.001. Despite the much lower abundance than O ii, O iii provides a similar amount of cooling as O ii because the [O iii] 4959, 5007 transition is connected to the ground state, whereas [O ii] 7320, 7330 connects two excited states. There is also some [O iii] 4363, although a factor ∼10 weaker than [O iii] 4959, 5007, as this line arises from a higher state.

The f = 0.001 model (panel C) fits some observed lines reasonably well, but the [O i] 6300, 6364 line is too weak by a factor 3, and in general the contrast between the emission lines and the quasi-continuum is too low in the 3 M_⊙ models. There is also very little Mg i emission in these models as the ionization state is high (even some Mg III is present).

The 3 M_⊙ models are discussed further in Sections 4.3.4 and 4.4.

4.3.2. 10 M_⊙ Models

The 10 M_⊙ models are shown in the middle panel of Figure 10. These have T = 6000–6500 K, ${x}_{e}\sim 0.5$, and small Mg i abundance. Increasing the mass to 10 M_⊙ allows reproduction of the [O i] 6300, 6364 luminosity as O i is more common. The dramatic difference in [O i] 6300, 6364 strength between 3 M_⊙ and 10 M_⊙ at f = 0.1 illustrates the potentially strong non-linearity of the ionization solutions. In the 3 M_⊙ models at f = 0.1, O i has suffered runaway ionization, whereas in the 10 M_⊙ models this has been avoided.

The [O ii] 7320, 7330 and [O iii] 4959, 5007 lines are less pronounced than in the 3 M_⊙ models as the ionization is lower. [O ii] 7320, 7330 now tends to be too weak to explain the 7300 line, but [O iii] 4959, 5007 makes a reasonable fit to SN 2007bi (while too weak for LSQ14an).

The f = 0.1 and f = 0.01 models have neglegible Mg i] 4571. The f = 0.001 model brings out both Mg i] 4571 and Mg i 5180, and provides an overall reasonable fit to several parts of the spectrum.

4.3.3. 30 M_⊙ Models

To illustrate some dependency on deposition, the models plotted in Figure 10 use a different (higher) energy deposition than those plotted in Figure 9. The f = 0.1 and f = 0.01 models have a deposition $5.0\times {10}^{41}$ erg s⁻¹ instead of $2.0\times {10}^{41}$ erg s⁻¹, and the f = 0.001 model has deposition $20\times {10}^{41}$ erg s⁻¹ instead of $10\times {10}^{41}$ erg s⁻¹.

The best fit is produced for f = 0.01. The higher mass leads to lower temperatures and damping of the blue continuum seen in the 3 M_⊙ and 10 M_⊙ models. The more neutral gas state also leads to easier formation of Mg i lines, although it remains difficult to reproduce the full strength of Mg i] 4571.

4.3.4. A Minimum Oxygen Mass Requirement

Figure 12 shows observed and modeled [O i] 6300, 6364 line luminosities. The "continuum" is subtracted in the following manner. We measure

$$\begin{eqnarray}&&{I}_{1}={\int }_{{\lambda }_{0}-d{\lambda }_{1}}^{{\lambda }_{0}+d{\lambda }_{1}}{F}_{\lambda }d\lambda \end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{I}_{2}={\int }_{{\lambda }_{0}-d{\lambda }_{2}}^{{\lambda }_{0}+d{\lambda }_{2}}{F}_{\lambda }d\lambda \end{eqnarray} \tag{ 3 }$$

where ${\lambda }_{0}$ is the line center wavelength, $d{\lambda }_{i}={\lambda }_{0}{{dV}}_{i}/c$, and we use ${{dV}}_{1}=8000$ km s⁻¹and ${{dV}}_{2}={\rm{10,000}}$ km s⁻¹. We should have, for a linear continuum of average strength F_c:

$$\begin{eqnarray}&&{I}_{1}={L}_{\mathrm{line}}+2{F}_{c}d{\lambda }_{1}\end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray}&&{I}_{2}={L}_{\mathrm{line}}+2{F}_{c}d{\lambda }_{2}\end{eqnarray} \tag{ 5 }$$

then

$$\begin{eqnarray}&&{L}_{\mathrm{line}}={I}_{1}-\displaystyle \frac{{I}_{2}-{I}_{1}}{d{\lambda }_{2}/d{\lambda }_{1}-1}.\end{eqnarray} \tag{ 6 }$$

Figure 12 shows the [O i] 6300, 6364 luminosity in the model grids compared to the observed values. As the top panel shows, the 3 M_⊙ models are always more than a factor 3 too dim, no matter what density and amount of energy deposition we choose. From this a minimum mass of ∼10 M_⊙ is needed. The 10 M_⊙ and 30 M_⊙ model grids (middle and bottom panels) confirm that such O-zone masses can produce sufficient luminosity.

There is little quasi-continuum in [O i] 6300, 6364 in the observed spectra (Figure 7), which makes uncertainty in the continuum subtraction procedure minor. Also note that we have ignored any correction due to extinction, which makes the observed luminosities somewhat underestimated, going in the right direction for the mass limit. Finally, we expect little contamination of this feature by other lines.

4.4. O iii Lines

Lunnan et al. (2016) recently suggested the identification of broad [O iii] 4363 and [O iii] 4959, 5007 lines in PS1-14bj and LSQ14an at +200d. Figure 15 in the Appendix shows a zoom-on in these lines in SN 2007bi, LSQ14an, and SN 2015bn. SN 2007bi shows no clear broad [O iii] lines, even though there may be hints of them. The narrow galaxy lines of Hγ 4341 and [O iii] 4363 may merge together to form a line of some apparent broadening. LSQ14an has strong galaxy lines (clipped out in figures), but there are also distinct broad [O iii] components at +365d and +410d, with expansion velocities of 4000–5000 km s⁻¹. SN 2015bn is similar to SN 2007bi, showing some hints of the [O iii] lines, although they appear somewhat more pronounced in the spectra in Nicholl et al. (2016a).

In the O-zone models, these [O iii] lines emerge for many combinations of mass and density, The [O iii] 4959, 5007 lines have lower excitation energy than [O iii] 4363 and emerge more prominently. The emergence of [O iii] 4959, 5007 is seen for most masses and filling factors, illustrating that [O i] 6300, 6364 and [O iii] 4959, 5007 can be produced in the same physical region. However, [O iii] 4363 is weak in the 10 M_⊙ and 30 M_⊙ models, and only distinct in the 3 M_⊙ models. This line requires high temperature, which is more easily obtained for a low mass.

As an illustration of this point, consider Figure 13, which compares the LSQ14an +410d spectrum with model Cburn3-0.1-10. This model produces [O ii] 7320, 7330, [O iii] 4363, and [O iii] 4959, 5007 in good correspondance with the observed lines. The good agreement suggests that LSQ14an has a low-mass, dilute oxygen region reprocessing a large amount of energy. SN 2007bi and SN 2015bn show no evidence for such a component, but the models show that the O ii and O iii lines are sensitive to deposition, rapidly decreasing in strength if this is lower than around ${10}^{42}$ erg s⁻¹. Thus, it is possible that LSQ14an displays these lines due to its higher energy deposition rather than any fundamental difference in ejecta properties.

Zoom In Zoom Out Reset image size

It should be noted that the bulk of line flux seen in the O ii and O iii lines are at lower velocity (of order 4000 km s⁻¹, see also Figure 15) than the 8000 km s⁻¹ used in the model, which was chosen based mainly on the widths of neutral lines of O and Mg. It is interesting that the O ii/O iii velocity scale is the same as the observed O i 7774 one. Whereas production of O i 7774 is facilitated at high density, production of O ii and O iii lines is facilitated at low density.

The [O iii] lines emerge despite there being a relatively low O iii fraction in the models. They are efficient cooling lines, significantly more so than for example [O ii] 7320, 7330, and can be prominent even with a low O iii abundance. Still, the models with the strongest [O iii] lines also have strong [O ii] 7320, 7330.

High-ionization lines can presumably arise both in an inner pulsar-wind nebula (Chevalier & Fransson 1992; Gaensler & Slane 2006; Metzger et al. 2014), or in an outer interaction region (Chevalier & Fransson 1994; Fransson et al. 2002). Detailed modeling of such scenarios may be able to provide more constraints on this component.

4.5. Oxygen Recombination Lines

All three SNe show clear detection of O i 7774, and SN 2015bn additionally of O i 9263 and O i 1.13 μm. These lines are usually interpreted as recombination lines (Maurer & Mazzali 2010; Jerkstrand et al. 2015). The O i 7774 luminosity is expected to obey, for ${n}_{e}={n}_{{\rm{O}}{\rm{II}}}$ (Jerkstrand et al. 2015, J15a hereafter)

$$\begin{eqnarray}{L}_{7774} & = & { \mathcal V }\times {n}_{e}^{2}{\alpha }^{\mathrm{eff}}(T)h\nu \\ & = & 4.5\times {10}^{41}{{\rm{erg}}{\rm{s}}}^{-1}{\left(\displaystyle \frac{V}{8000{{\rm{km}}{\rm{s}}}^{-1}}\right)}^{3}f\\ & & \times {\left(\displaystyle \frac{t}{400{\rm{d}}}\right)}^{3}{\left(\displaystyle \frac{{n}_{e}}{{10}^{8}{{\rm{cm}}}^{-3}}\right)}^{2}\left(\displaystyle \frac{{\alpha }^{\mathrm{eff}}(T)}{2\times {10}^{-13}\,{{\rm{cm}}}^{3}{{\rm{s}}}^{-1}}\right).\end{eqnarray} \tag{ 7 }$$

Here ${\alpha }^{\mathrm{eff}}(T)$ is the effective recombination rate and ${ \mathcal V }$ is the volume. The observed luminosity in SN 2015bn (and similarly in the others) is ${L}_{7774}=3\times {10}^{40}$ erg s⁻¹ which translates to ${n}_{e}=3\times {10}^{7}{f}^{-1/2}$ cm⁻³ (using V = 8000 km s⁻¹), or ${n}_{e}={10}^{9}$ cm⁻³ for f = 0.001. For V = 3000 km s⁻¹ a factor 4 higher value results. Because we have indications of $f\lt 0.1$ from the analysis in Section 4.2, the electron density becomes constrained to ${n}_{e}={10}^{8}\mbox{--}{10}^{9}$ cm⁻³.

The O i 9263 and O i 1.13 μm recombination lines have theoretical luminosities about a factor 2 and 3, respectively, lower than the O i 7774 line, in good agreement with the observed ratios. This strengthens the case that they are indeed formed by recombination.

One limitation of the analytic treatment is the assumption that ${n}_{e}={n}_{{\rm{O}}{\rm{II}}}$. Charge transfer reactions may cause significant deviation from this equality at high densities (J15a), and this could lead to an underestimate of n_e. As discussed in Section 4.2, at very high density O i 7774 may also obtain contributions by cooling, which would lead to an overestimate of n_e by assuming pure recombination.

The zone mass is

$$\begin{eqnarray}M & = & { \mathcal V }\times {n}_{e}{x}_{e}^{-1}\bar{A}{m}_{p}=\displaystyle \frac{{L}_{7774}}{{n}_{e}{\alpha }^{\mathrm{eff}}h\nu }{x}_{e}^{-1}\bar{A}{m}_{p}\\ & = & 80\ {M}_{\odot }{\left(\displaystyle \frac{{n}_{e}}{{10}^{8}{{\rm{cm}}}^{-3}}\right)}^{-1}{\left(\displaystyle \frac{{x}_{e}}{0.1}\right)}^{-1}\left(\displaystyle \frac{\bar{A}}{16}\right)\\ & & \times {\left(\displaystyle \frac{{\alpha }^{\mathrm{eff}}(T)}{2\times {10}^{-13}{{\rm{cm}}}^{3}{{\rm{s}}}^{-1}}\right)}^{-1}\end{eqnarray} \tag{ 8 }$$

where $\bar{A}$ is the mean atomic weight. The uncertainties in n_e, x_e and ${\alpha }^{\mathrm{eff}}$ are somewhat too large to make a useful constraint, but if these can be constrained further a zone mass can be estimated.

4.6. Mg i 1.50 $\mu {m}$ Recombination Line

J15a demonstrated that the Mg i 1.50 μm recombination line is a useful diagnostic, relatively free from uncertainty in formation and blending. Its luminosity is expected to be given by

$$\begin{eqnarray}{L}_{1.5\mu {\rm{m}}} & = & { \mathcal V }\times {n}_{e}{n}_{\mathrm{Mg}{\rm{II}}}{\alpha }^{\mathrm{eff}}(T)h\nu =6.6\times {10}^{38}\,{\text{erg s}}^{-1}\\ & & \times \left(\displaystyle \frac{{M}_{\mathrm{Mg}}}{1\,{M}_{\odot }}\right)\left(\displaystyle \frac{{n}_{e}}{{10}^{8}\,{{\rm{cm}}}^{-3}}\right)\left(\displaystyle \frac{{\alpha }^{\mathrm{eff}}(T)}{1\times {10}^{-13}\,{{\rm{cm}}}^{3}{{\rm{s}}}^{-1}}\right)\end{eqnarray} \tag{ 9 }$$

where the second step assumes most magnesium to be singly ionized as typically obtained in the models. The observed line in SN 2015bn has an expansion width of 5000 km s⁻¹. The measured luminosity within this velocity is ${L}_{1.5\mu {\rm{m}}}\,=1.7\times {10}^{40}$ erg s⁻¹, but about half is continuum so we take ${L}_{1.5\mu {\rm{m}}}=1\times {10}^{40}$ erg s⁻¹. For comparison, the measured luminosity in the Type IIb SN 2011dh at 200d was $1\times {10}^{38}$ erg s⁻¹, a factor 100 lower (J15a). From Equation (9), ${M}_{\mathrm{Mg}}=15\ {M}_{\odot }{({n}_{e}/{10}^{8}{{\rm{cm}}}^{-3})}^{-1}$. As any proposed model scenario would have ${M}_{\mathrm{Mg}}\lt 15$ M_⊙; this indicates ${n}_{e}\gt {10}^{8}$ cm⁻³. If we take ${n}_{e}={10}^{8}\mbox{--}{10}^{9}$ cm⁻³ as estimated from the O i recombination lines, the Mg mass is 1.5–15 M_⊙. The effective recombination rate depends on temperature as well as on whether case B or case C holds, but the variation between these cases and T = 2500–7500 K is less than a factor 2 (see Table C2 in J15a). As the Mg mass fraction is typically 10% in C-burning ash, the inferred OMg zone mass is very large, 15–150 M_⊙.

Model OMg30-0.001-10 reproduces the Mg i 1.50 μm line reasonably well (Figure 11). This model has ${M}_{\mathrm{Mg}}=3$ M_⊙, ${n}_{e}=1.2\times {10}^{9}$ cm⁻³, and T = 6500 K, giving a predicted luminosity of $2\times {10}^{40}$ erg s⁻¹ using Equation (9). All other models (lower M_ej and/or larger f) produce weaker emission. The Mg i 1.50 μm emission is fainter by a factor of at least 10 in the 3 M_⊙ models, and the 10 M_⊙ models are on the low side as well. The Mg i 1.50 μm line thus provides strong independent support for a large OMg zone mass ($M\gt 10$ M_⊙) and a small filling factor ($f\sim {10}^{-3}$).

4.7. The Ca ii Lines

The observed ratio between Ca ii NIR and [Ca ii] 7291, 7323 in SN 2015bn at +315d is 1.7, unusually high for the nebular phase. There is some uncertainty in the true Ca ii NIR luminosity, because O i 8446 will scatter into the triplet, and the red side may be contaminated by [C i] 8727. In addition, scattering in the triplet and fluorescence from the Ca ii HK lines may contribute to the luminosity. However, O i 8446 is predicted to be only about half as strong as O i 7774, which would make up only 10%–20% of the triplet luminosity. There is also no obvious sign of [C i] 8727 by an asymmetric red side. The absorption trough of Ca ii HK also does not appear to allow much of the observed luminosity in the branching lines. Nevertheless, these contaminations mean an emissivity ratio as low as 1 may be possible, and we take a ratio 1–2 as being admissable.

Figure 14 (third row) shows the ratio of Ca ii NIR/[Ca ii] emission in on-the-spot NLTE calculations varying $n(\mathrm{Ca}\,{\rm{II}})$, n_e, and T. These models show that a ratio $\gt 1$ requires ${n}_{e}\gtrsim {10}^{8}$ cm⁻³. The curves do not change much for n_e higher than 10⁸ cm⁻³. For ${n}_{e}\geqslant {10}^{8}$ cm⁻³ there are large plateau regions where the ratio is 1–2. As the panels showing departure coefficients and optical depths illustrate, this is the optically thick LTE regime.

Zoom In Zoom Out Reset image size

The regime producing the observed ratio allows for filling factors f = 0.001–1 to give the right luminosity (using the top panels), and temperatures $T\gt 4000$ K, thus no strong constraints on these parameters. The constraint ${n}_{e}\gtrsim {10}^{8}$ cm⁻³ translates to a total number of electrons

$$\begin{eqnarray}&&{N}_{e}\gt { \mathcal V }\times f\,\times \,{10}^{8}\,=\,9\times {10}^{57}f.\end{eqnarray} \tag{ 10 }$$

The corresponding zone mass is

$$\begin{eqnarray}&&M=\displaystyle \frac{{N}_{e}}{{x}_{e}}\bar{A}{m}_{p}\gt 3000{M}_{\odot }{\left(\displaystyle \frac{{x}_{e}}{0.1}\right)}^{-1}\left(\displaystyle \frac{\bar{A}}{40}\right)\,f\end{eqnarray} \tag{ 11 }$$

to reach realistic masses this requires $f\ll 1$, independently confirming a small filling factor as found from Mg and O lines. although the Ca lines likely come from a O-burning zone rather than a C-burning zone, a similar amount of clumping seems reasonable.

Several of the C-burning models in Section 4 produce Ca ii NIR/[Ca ii] ratios larger than unity. In these models the calcium number fraction is small, 10⁻⁵. For example, model C-burn10-0.001-5.0 has $T=6000\,{\rm{K}}$, x_e = 0.1, ${n}_{e}={10}^{8}$ cm⁻³, and $n(\mathrm{Ca}\,{\rm{II}})=7\times {10}^{4}$ cm⁻³. Its Ca ii NIR/[Ca ii] ratio of a few agrees well with the NLTE grid prediction. While the luminosity of the Ca lines is an order of magnitude too faint compared to SN 2015bn, this illustrates the potential of small calcium abundances to make strong lines.

The Ca ii NIR triplet is observed also at two other nebular epochs of +343d and +392d in (Nicholl et al. 2016a). In these spectra the Ca ii NIR/[Ca ii] 7291, 7323 ratio is around 1, somewhat lower than seen here at +315d. The ratio is expected to decrease with time as both electron density and temperature decrease. The Ca ii NIR line is also observed in LSQ14an (Figure 6), but the likely dominance of [O ii] 7320, 7330 over [Ca ii] 7291, 7323 in this SN prevents a determination of the Ca ii NIR/[Ca ii] 7291, 7323 ratio.

The similarity of nebular-phase spectra of SLSNe to those of broad-lined SNe Ic established in Section 3 suggests the possibility of a common origin for these explosions. This idea is also interesting in the context of the newly recognized class of ultra-long GRBs (Levan et al. 2014), the association between one of these and a bright Type Ic SN with early-time spectral appearance similar to SLSNe (Greiner et al. 2015; Kann et al. 2016), as well as progress in modeling rapidly rotating and magnetized cores (e.g., Mösta et al. 2015).

The late-time spectra of SN 1998bw have been modeled as 5–7 M_⊙ of O-rich ejecta powered by 0.6–0.7 M_⊙ of ⁵⁶Ni (Mazzali et al. 2001). While there is some debate whether broad-lined SNe Ic are powered by ⁵⁶Ni or a central engine, existing ⁵⁶Ni models are successful enough to maintain this as the standard picture. Thus, "scaled-up" ⁵⁶Ni powered explosion models must still be pursued for SLSNe, even though the PISN variants of these have not compared well to observations so far (Dessart et al. 2013; Jerkstrand et al. 2016).

The late-time spectra of SN 2015bn and SN 2007bi are shown here to be almost identical (Figure 7). With the multi-epoch coverage of SN 2015bn (and LSQ14an) we see significantly more rapid decline than ⁵⁶Co decay (see also (Nicholl et al. 2016a) and C. Inserra et al. 2016, in preparation). This cannot be achieved in PISN models, which have complete gamma-trapping for at least 500d (Jerkstrand et al. 2016). Some transfer to the near-infrared may occur, but the NIR observations here indicate that this is a weak process. Thus, if these are radioactivity-powered SNe, they must be core-collapse explosions of less massive progenitors where gamma-ray escape can occur earlier.

Analytic light curve models can estimate ejecta masses, but are limited by assumed opacities and explosion energies. The nebular phase can provide independent estimates in a phase when the physics is different. Here we have found ${M}_{{\rm{O}} \mbox{-} \mathrm{zone}}\approx 10\mbox{--}30$ M_⊙ for long-duration SLSNe from a variety of indicators. An immediate consequence is a minimum progenitor mass of ${M}_{\mathrm{ZAMS}}\gtrsim 40$ M_⊙ in a single-star progenitor scenario (e.g., Hirschi et al. 2004; Nomoto et al. 2006; Woosley et al. 2007; Ekström et al. 2012; Chieffi & Limongi 2013).

Several independent observables indicate that this emitting mass occupies a small fraction $f\lesssim 0.01$ of the expansion space given by the line widths, and this puts important constraints on any model scenario. To some extent a value $f\ll 1$ would follow if most luminosity is generated in a dense central region with $V\ll {V}_{\max }$. Inspection of the line profiles indicates, however, that much of the luminosity must come from regions with $V\gtrsim 0.5{V}_{\max }$, suggesting true clumping or shell formation of high-velocity material.

Compression of the ejecta is conceivable in several scenarios. In one-dimensional magnetar models, the ejecta are pushed into a thin shell with $f\sim 0.001$ as they are accelerated by the pulsar wind pressure (Kasen & Bildsten 2010; Chen et al. 2016a). In 2D models, instabilities lead to some degree of fragmentation of this shell, although its main properties persist (Chen et al. 2016a). The instabilities lead to a somewhat reduced clumping level, and the mixing of material toward both lower and higher velocity, thus giving line profiles deviating from the boxy profiles which would be produced by a thin shell. Since none of the observed line profiles are boxy, the instabilities and fragmentation are one way to avoid a uniform shell.

The ⁵⁶Ni bubble effect in radioactivity models (Herant & Benz 1991; Basko 1994) provides also some degree of compression of the O zone, although with more modest factors. As the ⁵⁶Ni zone is hotter than the surrounding substrate due to radioactive decay in the first few days, it expands, an effect ignored in many explosion models. Kozyreva et al. (2016, in preparation) have investigated the ⁵⁶Ni bubble effect for the PISN explosion of a 130 M_⊙ He core, using the radiation hydrodynamics code STELLA. The ⁵⁶Ni zone was diluted by a factor of ∼2, and the expansion created compression in the external Si/S shell, increasing its density by a factor of 2–3 (A. Kozyreva 2016, private communication). Further out also the O zone showed some compression, but by a moderate factor of ∼20% (A. Kozyreva 2016, private communication). In 3D simulations of core-collapse SNe Si and O is somewhat more mixed, possibly allowing a compression of factor 2–3 also of oxygen. This modest effect still leaves the filling factor at values between f = 0.1–1, higher than indicated by the modeling here.

In circumstellar interaction models, formation of a cool dense shell provides a way to compress the material (Chevalier & Fransson 1994). It may be difficult, however, to obtain a high enough mass ($M\gtrsim 10$ M_⊙) in such a shell.

One model which cannot easily meet the requirements derived here is that of colliding shells in subsequent pulsational PISN eruptions. Models for these events predict O-shell masses of one or a few solar masses (Woosley et al. 2007; Yoshida et al. 2016), which would not provide enough total mass to explain the nebular spectra. One should note that a star undergoing pulsational pair-instability eruptions may eventually explode as a core-collapse supernova, in which the mass can be much higher. Thus, we have no strong constraints on such a scenario, only a scenario of collisions between pulsational (non-terminal) eruptions.

The shell collision scenario is problematic also on velocity-grounds. Woosley et al. (2007) lists masses and energies for pulses on a grid of He stars between 48 and 60 M_⊙, and shell velocities are mostly below 1000–2000 km s⁻¹. The same holds for the 54–61M_⊙ CO core models presented by Yoshida et al. (2016). It would clearly be difficult to reproduce the observed high bulk velocities in the SNe studied here. While photospheric spectra can trace velocities in small amounts of fast-moving ejecta, the nebular phase probes the bulk velocity which can be more directly linked to minimum E/M ratios. The light curves in Yoshida et al. (2016) also typically do not reach the SLSN range, unless the final CCSN is involved in the interaction.

Most strong observed nebular emission lines have a secure identification. The largest uncertainty regards the line at 5250 Å. In the models for SN 1998bw, Mazzali et al. (2001) obtains a strong line here due to a blend of [Fe ii] lines. The models in Sollerman et al. (2000) produce a blend of Mg i, Fe i], and [Fe ii] lines, with Mg i 5180 appearing to dominate. Some further discussion of the 5250 Å feature is given in Maeda et al. (2002). Improved understanding of this line is an important goal for further modeling.

Several sources of uncertainty remain for the modeling of these spectra. One of the assumption in the SUMO models studied here is that the ejecta have reached a steady-state phase, when timescales for the radiative transport and atomic processes are short compared to the dynamic time and the timescale for energy input. The biggest concern is whether UV/blue photons are trapped in the transfer for a long enough time that different parts of the spectrum are generated by power input at different times (i.e., energy is still diffusing out). Fluorescence should limit this time-dependence, but the exact magnitude of the effect is hard to judge. Another issue is that the influence of the radiative transport is only qualitatively captured in the simplistic morphology assumed here. While the ejecta are hot and the optical/NIR spectra are expected to be dominated by collisional cooling, the ionization balance of magnesium and similar elements with low ionization potential are governed by the radiation field. The models are also limited by the assumptions of a single density, the number of clumps, and composition. However, these uncertainties are mitigated by analysis with semi-analytic approaches which consistently point in the same direction of massive and dense ejecta.

The metallicity of the progenitor star has an influence on the emission from the OMg zone through several mechanisms including cooling and line blocking. For the solar metallicity carbon burning composition explored here, lines such as Mg i] 4571 tend to be quenched by the trace metals, worsening the fits to observations compared to pure OMg compositions. However, there is dependency on other model parameters such as the clump size, and a separate in-depth study would be needed to see if any firm constraints on metallicity can be placed from nebular spectra. Running a single model at ten times reduced metallicity did not produce any great changes.

The inferred electron densities of ${n}_{e}={10}^{8}\mbox{--}{10}^{9}$ cm⁻³ at ∼400d after explosion are extraordinarily high. As a comparison, the electron density in the Type IIb SN 2011dh at +200d was inferred to be ∼10⁶ cm⁻³ (J15a). Values of over 10⁸ cm⁻³ were inferred for the Type IIn SN 1995n (Fransson et al. 2002). Ben-Ami et al. (2014) inferred ${n}_{e}\sim {10}^{9}$ cm⁻³ for circumstellar matter emitting narrow lines in the (non-SLSN) Type Ic SN 2010mb. Because we still see mostly neutral species, the gas must have ${x}_{e}\lt 1$. The inferred zone masses are large, at least 10 M_⊙, in agreement with modeling of the thermal cooling lines.

We have presented nebular-phase PESSTO and VLT spectra of the long-duration Type Ic SLSNe LSQ14an (+365d and +410d post-peak rest frame) and SN 2015bn (+250d and +315d), respectively, including the first nebular NIR spectra of SLSNe. We have carefully accounted for host galaxy contamination using galaxy models tied to the observed galaxy photometry and also applied this method to recalibrate the nebular phase spectra of SN 2007bi. We show that this makes a large difference to the spectral appearance, and after host correction the spectra of SN 2007bi, LSQ14an, and SN 2015bn show much similarity at about one year post peak, demonstrating homogeneity within the class. These spectra are in turn similar to those of broad-lined SNe Ic such as SN 1998bw at earlier phases (140–200 days), suggesting a possible common origin with GRB-SNe. The major difference seen in these spectral comparisons is that LSQ14an shows strong [O iii] emission lines, and probably [O ii] 7320, 7330, while they are significantly weaker in SN 2015bn and SN 2007bi.

SN 2015bn shows O i 9263, O i 1.13 μm and Mg i 1.50 μm, but no strong emission by [Si i] and [S i] around 1.08 μm as predicted by PISN models. The NIR spectra of LSQ14an and SN 2015bn show strong Ca ii NIR triplet emission at nebular times. NLTE models indicate that an electron density ${n}_{e}\gtrsim {10}^{8}$ cm⁻³ is needed to reproduce the observed Ca ii NIR/[Ca ii] 7291, 7373 ratio of ∼2 at +315d in SN 2015bn. Analytic formulae for the O i recombination lines yield similar values ${n}_{e}={10}^{8}\mbox{--}{10}^{9}$ cm⁻³. These are extraordinarily high electron densities for the nebular phase.

Models for oxygen-zone emission show that only models with $\gtrsim 10$ M_⊙ of O can produce enough [O i] 6300, 6364 luminosity to match observations, irrespective of the powering situation and the density. Too much energy deposition in lower-mass models leads to ionization to O ii before [O i] 6300, 6364 reaches the observed brightness. Typical models produce detectable broad [O iii] 4959, 5007 lines, strengthening the recent identification in observed spectra. The Mg i 1.50 μm recombination line in SN 2015bn indicates a large Mg mass of 1.5–15 M_⊙, which together with a typical mass fraction of 10% gives a similar constraint of $M({\rm{O}} \mbox{-} \mathrm{zone})\gtrsim 15$ M_⊙. These masses require large progenitor ZAMS masses of at least 40 M_⊙.

To produce Mg i] 4571 and several other lines approaching their observed brightness, the O/Mg zone needs to be strongly clumped with a filling factor $f\lesssim 0.01$. A similar constraint is derived for the zone emitting the Ca ii lines.

Overall we can view these results in the context of the three competing physical models which have been advanced to explain SLSNe: ⁵⁶Ni powered luminosity, engine driven explosions, and circumstellar interaction.

• 1.  
  SLSNe Type Ic show little resemblance with PISN spectral models. In addition to optical comparisons published previously, two strong lines of Si i and S i at 1.08 μm are predicted (Jerkstrand et al. 2016) but are not observed. We cannot rule out that the luminosity is actually powered by ⁵⁶Ni, however, and a massive core-collapse scenario is attractive given the strong similarities with SN 1998bw, assuming this was ⁵⁶Ni powered. The requirement for compact regions of emission could possibly be explained by the ⁵⁶Ni bubble effect, but simulations need to demonstrate the low value $f\lesssim 0.01$ for the oxygen material inferred here.
• 2.   
  In the magnetar scenario, pulsar wind compression of the ejecta is one promising mechanism that can explain how the ejecta may be compressed into a dense region with a low filling factor. A thin shell would result in boxy line profiles, which are not observed, but instabilities and fragmentation as demonstrated in recent 2D simulations may produce profiles in better agreement with observed ones. The large ejecta masses and composition similar to type Ic SNe would imply magnetar formation from massive Wolf-Rayet progenitors.
• 3.  
  Dense shell formation in circumstellar interaction models is also a conceivable scenario for forming the compact regions and producing an emitting zone with low filling factor. The fact that we find some evidence of a second, lower-mass and lower-density region producing [O ii] and [O iii] lines might also support an interaction scenario. However a mass of ∼10 M_⊙ in the shell that produces the neutral oxygen and magnesium emission is required. The collision of shells from pulsational pair-instability eruptions is one variant of the CSI scenario, but we find that the ejecta masses and velocities in published models are too low compared to those inferred from observations.

Finally, we note that our analysis is necessarily limited to the group of slowly evolving Type Ic SLSNe. There are no nebular spectra yet available for the faster declining, and more common SN 2005ap-like events. This awaits the discovery of a low-redshift ($z\lesssim 0.05$) event which is close enough to be observed in the nebular phase.

We thank S. Blinnikov, C. Fransson, A. Kozyreva, P Mazzali, and C. Harris for discussion. A.J. thanks the organisers and participants of MIAPP conference "The Physics of Supernovae" and the workshop "Supernovae: The Outliers" for stimulating discussions. A.J. acknowledges funding by the European Union's Framework Programme for Research and Innovation Horizon 2020 under Marie Sklodowska-Curie grant agreement No 702538, and Science and Technology Facilities Council DIRAC computing grants ACSP45 and ACSP74. The research leading to these results has received funding from the European Research Council under the European Union's Seventh Framework Programme (FP7/2007-2013)/ERfC Grant agreement n^o [291222]. Agreements 307260, 320360, and 615929 are also acknowledged. SJS acknowledges funding from STFC grants ST/I001123/1 and ST/L000709/1. This work is based in part on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile as part of PESSTO, (the Public ESO Spectroscopic Survey for Transient Objects Survey) ESO program 188.D-3003, 191.D-0935 and on the VLT ESO Programmes 094.A-0645, 096.D-0191and 296.D-5042. K.M. acknowledges support from the STFC through an Ernest Rutherford Fellowship. S.T. acknowledges support by TRR 33 The Dark universe of the German Research Foundation. T.-W.C. and T.K. acknowledge the support through the Sofia Kovalevskaja Award to P. Schady from the Alexander von Humboldt Foundation of Germany. This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. We have made use of the Weizmann interactive supernova data repository-http://wiserep.weizmann.ac.il. The Pan-STARRS1 Surveys (PS1) have been made possible through contributions of the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, Queen's University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under Grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation under Grant No. AST-1238877, the University of Maryland, and Eotvos Lorand University (ELTE) and the Los Alamos National Laboratory.

A.1. LSQ14an+365d Spectrum (X-shooter)

A.2. LSQ14an+410d Spectrum (FORS2)

B.1. SN 2015bn+250d Spectrum (EFOSC2)

B.2. SN 2015bn+315d Spectrum (X-shooter)

Figure 15 compares the main emission lines in velocity space. There is emissivity in Mg i] 4571, [Fe ii] 5250, and Ca ii NIR out to ∼10,000 km s⁻¹. [O i] 6300, 6364, [Ca ii] 7291, 7323, and O i 7774 appear somewhat narrower, with most flux contained within ∼5000 km s⁻¹. In LSQ14an there is an unusual shape of [Ca ii] 7291, 7323 + [O ii] 7320, 7330.

Zoom In Zoom Out Reset image size

We have added O iii to the spectral synthesis code, using 6 levels and 12 transitions. Energy levels and A-values were taken from NIST¹⁷ (Reader et al. 1980; Moore 1993; Wiese et al. 1996). Collision strengths were taken from Lennon & Burke (1994).

Consider a trace element with total number abundance n, in a gas with electron number density n_e. Let N be the number of ionizing photons emitted per second, and assume only the trace element can absorb these. As long as all these are absorbed, the ionization balance is

$$\begin{eqnarray}&&N=\alpha (T){{Vn}}_{e}n(1-{x}_{n})\end{eqnarray} \tag{ 12 }$$

where x_n is the neutral fraction, and V is the volume. The solution is

$$\begin{eqnarray}&&{x}_{n}=1-\displaystyle \frac{N}{\alpha (T){{Vn}}_{e}n}\end{eqnarray} \tag{ 13 }$$

At full ionization (${x}_{n}\to 0$), the number of recombinations is $\alpha (T){{Vn}}_{e}n$. If N exceeds this, there is no solution. What happens then is that depletion of the neutral abundance continues until escape of ionizing photons occur. The new balance will be given by

$$\begin{eqnarray}&&N(1-{e}^{-\tau })=\alpha (T){{Vn}}_{e}n\end{eqnarray} \tag{ 14 }$$

where $\tau ={{nx}}_{n}\sigma L={{nx}}_{n}\sigma {(3/4\pi )}^{1/3}{V}^{1/3}$. Then, for $\tau \ll 1$, and ignoring the numeric constant

$$\begin{eqnarray}&&{x}_{n}=\alpha (T){V}^{2/3}{{nN}}^{-1}{\sigma }^{-1}\end{eqnarray} \tag{ 15 }$$

Thus, if the gas is still optically thick as ${x}_{n}\ll 1$, x_n will jump down to a value given by this equation; there is catastrophic runaway depletion.

The generic solution is given by solving

$$\begin{eqnarray}&&{\rm{\Gamma }}(1-{e}^{-{\tau }_{0}{x}_{n}})=1-{x}_{n}\end{eqnarray} \tag{ 16 }$$

where ${\rm{\Gamma }}=N/(\alpha {{Vn}}_{e}n)$. Figure 16 shows illustrative solutions, for ${\tau }_{0}={10}^{3}$ and three different Γ values. As Γ increases from zero, through 0.1 and toward 1, there is a gradual reduction in x_n. As ${\rm{\Gamma }}\gt 1$, the solutions jump down to a value ${x}_{n}\lesssim 1/{\tau }_{0}$. Changing the ionization parameter by a factor of 10 from 0.1 to 1 decreases x_n from 0.9 to 0.005, a factor 180. The decrease then slows down and follows $1/{\rm{\Gamma }}$. Thus, the effect of changing Γ by a factor 2 is small for ${\rm{\Gamma }}\ll 1$ and for ${\rm{\Gamma }}\gg 1$. But around ${\rm{\Gamma }}\sim 1$, a factor 2 change in Γ can give huge changes in x_n, of order ${\tau }_{0}$. A similar situation can be shown to hold if n_e is not assumed constant but ${n}_{e}=n(1-{x}_{n})$.

Zoom In Zoom Out Reset image size