The Nickel Mass Distribution of Stripped-envelope Supernovae: Implications for Additional Power Sources

The "light-curve tail" refers to the late-time evolution of the SN light curve once it enters a phase of linear decline in magnitude versus time. For SESNe, the light-curve tail typically begins at an earlier epoch (t ≳ 60 days post-explosion) than for H-rich Type II SNe, for which the tail is observable only after the H-recombination plateau phase ends (t ≳ 90 days post-explosion). It is widely thought that the tail of core-collapse SNe is powered by the ⁵⁶Co → ⁵⁶Fe radioactive decay chain (Colgate & McKee 1969). At this stage, the ejecta becomes transparent to the stored radiative energy; therefore, the observed luminosity traces the instantaneous heating rate. The γ-rays produced by the radioactive decay heat the ejecta, making the tail of the light curve an appropriate probe for measuring the amount of ⁵⁶Ni produced. The observed luminosity of the tail can then be modeled as

$$\begin{eqnarray}&&L\simeq {L}_{\gamma }\left(1-{e}^{-{(t/{T}_{0})}^{-2}}\right)+{L}_{\mathrm{pos}},\end{eqnarray} \tag{ 1 }$$

(Wygoda et al. 2019), where t is time since the explosion, L_γ is the luminosity produced by the radioactive decay of ⁵⁶Co and ⁵⁶Ni, and L_pos is the total energy release rate of positron kinetic energy. The term in parenthesis is a deposition factor, which represents the incomplete trapping of γ-rays, with T₀ denoting the partial trapping timescale of the tail. The deposition factor is proportional to $1-{e}^{-{t}^{-2}}$ for an explosion in homologous expansion (Sutherland & Wheeler 1984; Clocchiatti & Wheeler 1997). The luminosity terms in Equation (1) can be expressed as

$$\begin{eqnarray}&&{L}_{\gamma }={M}_{\mathrm{Ni}}\left(({\epsilon }_{\mathrm{Ni}}-{\epsilon }_{\mathrm{Co}}){e}^{-t/{t}_{\mathrm{Ni}}}+{\epsilon }_{\mathrm{Co}}{e}^{-t/{t}_{\mathrm{Co}}}\right)\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{L}_{\mathrm{pos}}=0.034{M}_{\mathrm{Ni}}{\epsilon }_{\mathrm{Co}}\left({e}^{-t/{t}_{\mathrm{Ni}}}-{e}^{-t/{t}_{\mathrm{Co}}}\right),\end{eqnarray} \tag{ 3 }$$

where _Ni = 3.9 × 10¹⁰ erg g⁻¹ s⁻¹ and _Co = 6.8 × 10⁹ erg g⁻¹ s⁻¹ are the specific heating rates of Ni and Co decay, respectively, and t_Ni = 8.8 days and t_Co = 111.3 days are their corresponding decay timescales. In this formulation, the escape of positrons, which happens on the timescale of several hundred days, is ignored. While the complete trapping of γ-rays is commonly assumed for H-rich Type II SNe due to their large ejecta masses and correspondingly long T₀, it is important to determine T₀ from the slope of the light-curve tail for SESNe since their ejecta masses are smaller and T₀ is usually comparable to the onset time of the radioactive tail. If the bolometric luminosity L can be ascertained observationally, the only unknowns in Equations (1) and (2) are M_Ni and T₀, which can be determined by fitting the slope and overall normalization of the radioactive tail of the light curve.

While robust, the "light-curve tail" method is not always accessible because the late-time radioactive tails are often faint, and thus more difficult to observe. As a result, M_Ni of SESNe is typically obtained with Arnett's rule, which states that the instantaneous heating rate from the radioactive decay of ⁵⁶Ni and ⁵⁶Co is equal to the bolometric luminosity of the SN at the light-curve peak. This can be rewritten as

$$\begin{eqnarray}&&{L}_{{\rm{p}}}={M}_{\mathrm{Ni}}\left(({\epsilon }_{\mathrm{Ni}}-{\epsilon }_{\mathrm{Co}}){e}^{-{t}_{{\rm{p}}}/{t}_{\mathrm{Ni}}}+{\epsilon }_{\mathrm{Co}}{e}^{-{t}_{{\rm{p}}}/{t}_{\mathrm{Co}}}\right),\end{eqnarray} \tag{ 4 }$$

where t_p and L_p are the peak time and peak bolometric luminosity, respectively. Arnett-like models make several assumptions to solve the thermodynamic differential equation, including homologous expansion of the ejecta, radiation-dominated pressure, spherical symmetry, and a self-similar energy density profile, and also adopt a radiation diffusion approximation (Arnett 1980, 1982).

KK19 showed that the assumption of self-similarity for the energy density profile will limit the accuracy of the Arnett-like models, especially for centrally located heating sources, due to the time-dependent evolution of the diffusion wave through the ejecta. Instead, they propose a new relationship between the peak time, t_p, and peak luminosity, L_p, without assuming self-similarity:

$$\begin{eqnarray}&&{L}_{{\rm{p}}}=\displaystyle \frac{2}{{\beta }^{2}{t}_{{\rm{p}}}^{2}}{\int }_{0}^{\beta {t}_{{\rm{p}}}}{t}^{{\prime} }{L}_{\mathrm{heat}}({t}^{{\prime} }){{dt}}^{{\prime} }\end{eqnarray} \tag{ 5 }$$

where L_heat(t) denotes a generic heating function and β is a dimensionless parameter of the order of unity. When L_heat(t) is powered by ⁵⁶Ni → ⁵⁶Co → ⁵⁶Fe radioactive decay, Equation (5) becomes

$$\begin{eqnarray}\begin{array}{rcl}{L}_{{\rm{p}}} & = & \displaystyle \frac{2{\epsilon }_{\mathrm{Ni}}{M}_{\mathrm{Ni}}{t}_{\mathrm{Ni}}^{2}}{{\beta }^{2}{t}_{{\rm{p}}}^{2}}\left[\left(1-\displaystyle \frac{{\epsilon }_{\mathrm{Co}}}{{\epsilon }_{\mathrm{Ni}}}\right)\left(1-(1+\beta {t}_{{\rm{p}}}/{t}_{\mathrm{Ni}}){e}^{-\beta {t}_{{\rm{p}}}/{t}_{\mathrm{Ni}}}\right)\right.\\ & & \left.+\displaystyle \frac{{\epsilon }_{\mathrm{Co}}{t}_{\mathrm{Co}}^{2}}{{\epsilon }_{\mathrm{Ni}}{t}_{\mathrm{Ni}}^{2}}\left(1-(1+\beta {t}_{p}/{t}_{\mathrm{Co}}){e}^{-\beta {t}_{p}/{t}_{\mathrm{Co}}}\right)\right].\end{array}\end{eqnarray} \tag{ 6 }$$

For this specific form of L_heat(t), the following relationships hold: M_Ni required to reproduce a fixed {t_p, L_p} pair will be directly proportional to the β value adopted. In contrast, if M_Ni is known, then the value of L_p or t_p required to reproduce a given {t_p, M_Ni} or {L_p, M_Ni} pair, respectively, will be inversely proportional to the value of β adopted (if the light curve is powered entirely by radioactive decay).

The parameter β incorporates the fact that L_p does not necessarily trace the radioactive heating rate because the stored internal energy of ejecta may lag or lead the observed luminosity, L(t), at the time of peak. The choice of β critically depends on several physical effects such as the spatial distribution of ⁵⁶Ni, the envelope composition, asymmetries in the heating source or ejecta, and all power sources contributing to the observed luminosity along with their exact heating functions. Using Equation (6), β can be derived from the light curves of SESNe with known L_p and t_p, if there is an independent constraint on M_Ni. We can then observationally calibrate the appropriate values of β using a sample of SESNe of various subtypes. With a sample of calibrated β values, one can potentially apply KK19's model to a wide range of SESNe with only photospheric data coverage to constrain their M_Ni. In addition, comparing the β values obtained from observed SESN light curves to those inferred from numerical light-curve models calculated with different input physics can inform us about the explosion details of SESNe.

Our SESN sample should consist of those SESNe for which a measurement of M_Ni can be made from both the light-curve tail and the light-curve peak. This will allow us to conduct an empirical comparison between M_Ni of SESNe obtained from the Arnett model and the radioactive tail and also provides the means to produce a data-driven calibration of the KK19 β values for SESNe. Therefore, we compiled the photometry of well-observed SESNe in the literature and selected SESNe that meet the following criteria:

• 1.  
  well-sampled coverage of the early rise (i.e., multiple observations before 5 days pre-maximum in at least one band) or the observation of an accompanying GRB/X-ray flash, since a constraint on the epoch of explosion is needed to derive M_Ni from Equations (1) and (2);
• 2.  
  multiple photometric measurements of the tail of the light curve, i.e., epochs ≳60 days, to be able to obtain M_Ni from the radioactive tail;
• 3.  
  reasonable coverage around the light-curve peak, which is required for both computing M_Ni using the Arnett model (Equation (4)) and calibrating the β parameter using KK19's model (Equation (6));
• 4.  
  light curves in at least two bands over the light-curve tail and peak, so that the bolometric luminosity and host galaxy reddening can be computed (see Section 3.5 for the method).

We identified 27 SNe from the literature that satisfy the above criteria and downloaded their photometric data from The Open Supernova Catalog (Guillochon et al. 2017). Our sample consists of eight IIb, eight Ib, four Ic, and seven Ic-BL SNe. These SNe are listed in Table 1 along with their basic properties, including SN type, the host galaxy name, distance estimate, extinction, and the epoch of explosion.

The distances we adopt throughout our analyses are listed in Table 1. We adopt up-to-date host galaxy distances reported in the NASA/IPAC Extragalactic Database (NED). ⁶ We prioritize distances obtained by Cepheids and tip-of-red-giant-branch methods when available and otherwise use cosmology-dependent values. For redshift-dependent distances, we adopt the standard ΛCDM cosmology with a Hubble constant H₀ = 73.24 km s⁻¹ Mpc⁻¹, matter density parameter Ω_M = 0.27, and vacuum density parameter Ω_Λ = 0.73 (Riess et al. 2016) and correct for Virgo, Great Attractor, and Shapley Supercluster infall. For one object (SN 2006el, which exploded in the galaxy UGC 12188), no host galaxy redshift is listed in NED. We therefore adopt the distance given in Drout et al. (2011), which is based on a host galaxy redshift reported in ATel 854 (Antilogus et al. 2006). We note that adopting Planck cosmological parameters (i.e., H₀ = 67.4 km s⁻¹ Mpc⁻¹, Ω_M = 0.315, and Ω_Λ = 0.685, Planck Collaboration et al. 2020) would increase redshift-dependent distances (∼80% of our sample) by ∼8%. Possible systematic effects of this choice on our results will be discussed in Section 5 below.

We adopt the value for galactic extinction along the line of sight to each SN reported in NASA/IPAC Infrared Science Archive ⁷ based on the extinction model of Schlafly & Finkbeiner (2011) and assuming an R_V = 3.1 extinction law. The resulting values are listed in Table 1.

To estimate values for host galaxy extinction, we use the intrinsic SN color-curve templates of Stritzinger et al. (2018b) and attribute the difference between the observed and the intrinsic color to the host galaxy reddening. Specifically, the host extinction can be written in terms of the observed minus intrinsic color as E(X − Y)_host = (X − Y)_obs − (X − Y)_int, where X and Y are the measured magnitudes corrected for the Galactic extinction in two different filters. When computing the extinction, we take the average of the color difference between the observed data and the templates of Stritzinger et al. (2018b) from 5 days to 10 days post-maximum. When determining this average, time of maximum is defined based on the observed filter that we adopt as the "X"-band in the above expression.

Whenever available, we use X − Y = V − R/r/i color indices. Since there is no template provided for V − R intrinsic color in Stritzinger et al. (2018b), we convert observed Johnson R-band photometry to Sloan r-band using the color transformation relation of Jordi et al. (2006) when required. E(X − Y) is then converted to the standard reddening E(B − V) using the bandpass coefficients of Schlafly & Finkbeiner (2011) assuming an R_V = 3.1 Milky Way extinction law. For those SNe for which photometric data are not available in any of the R/r/i filters, we adopt the B − V color index instead. Furthermore, we find that our obtained reddening is not robust to the choice of filters X and Y for several Type IIb SNe (i.e., SN 1993J, SN 2011dh, and SN 2013df). For these SNe, we take the average of E(B − V) values derived using different color indices such as V − r, V − i, and B − V. The final host galaxy extinctions used in our analyses are listed in Table 1.

It worth noting that Stritzinger et al. (2018b) constrained R_V directly for a set of eight SESNe with both optical and IR light curves and moderate reddening. They found values spanning a wide range (1.1 ≲ R_V ≲ 4.3), with tentative evidence for Type Ic SNe showing higher R_V values than Type IIb/Ib SNe. We continue to adopt a standard R_V = 3.1 Milky Way extinction law here, due primarily to the small number of SESNe for which this has been directly constrained, as well as for consistency with previous extinction estimates in the literature. However, we note that extinction corrections in the B/V/r/i bands considered here could vary by ∼0.05–0.3 mag for R_V values between 1.1 and 4.3. Possible systematic effects of this choice on our results will be discussed in Section 5 below.

The estimated epochs of explosion for each SN are presented in Table 1. For several of Type IIb SNe with double-peaked light curves, the epoch of explosion is adopted from the reported values obtained by modeling the shock cooling emission. For SN 1998bw, which is a Type Ic-BL SN associated with GRB 980425, we take the GRB epoch as the explosion epoch. Similarly, for SN 2008D, the epoch of the observed X-ray flash XRO 080109 is taken as the explosion epoch. Since the shock velocities of these SNe are high, we can assume that the GRB or X-ray flash occurs shortly after the explosion time; therefore, GRB 980425 and XRO 080109 provide accurate estimates of the explosion epochs. For the rest of the SNe in our sample, we estimate the explosion time by fitting a power law with the form

$$\begin{eqnarray}f=\left\{\begin{array}{ll}{f}_{0}{(t-{t}_{0})}^{n} & t\gt {t}_{0}\\ 0 & t\leqslant {t}_{0}\end{array}\right.,\end{eqnarray} \tag{ 7 }$$

to the observed fluxes, f, in the band with the best early light-curve coverage. We carry out least-square regression to fit for the power index n, scaling coefficient f₀, and epoch of explosion t₀. For the fitting, we only consider epochs that are pre-maximum-light and within 5 days of the first reported detection as well as any publicly available upper limits from nondetection. The uncertainties on the explosion epoch quoted in Table 1 come directly from the fitting process. The consequences of a possible "dark phase" (Piro & Nakar 2013) between the explosion epoch and the epoch of first light will be discussed in Section 5.1.2.

The relative paucity of SESNe with extensive coverage in UVOIR bands from early to late times is a challenge for computing bolometric light curves that are needed to obtain M_Ni. Here, we focus on obtaining the bolometric luminosities for epochs around the light-curve peak and the late-time tail. In order to leverage the multiband photometric data available for our SESNe sample, we adopt the bolometric correction (BC) coefficients of Lyman et al. (2014, 2016). These color-dependent coefficients were measured by fitting the spectral energy distributions of a sample of SESNe that have coverage in ultraviolet, optical, and infrared wavelengths, and can be utilized as long as light curves for an SN of interest are available in a minimum of two bands.

We first compute the absolute magnitude light curves for all SNe in our sample in the pair of bands indicated in Table 2 (column "BC bands"). We correct for the distances and total line-of-sight extinction described above. Next, we fit the multiband absolute magnitude light curves around the peak and over the radioactive tail with a spline-smoothing function and linear function, respectively. We perform a Monte Carlo (MC) analysis to propagate the uncertainties in the measured magnitudes and distance estimate. The fitted absolute magnitude light curves, which together also provide intrinsic color as a function of time, are then used to calculate a bolometric magnitude light curve by applying the color-dependent BC polynomials of Lyman et al. (2014, 2016). Specifically, M_bol = M_X + BC, where BC is computed using the color indices listed in column "BC bands" of Table 2, and X denotes the first indicated band listed in that column. Finally, we convert bolometric magnitudes to uminosities assuming M_bol,⊙ = 4.74 mag and L_bol,⊙ = 3.83 × 10³³ erg s⁻¹.

Table 2. The Fit Parameters for the SESN Sample

SN Name	Type	BC Bands	$\mathrm{log}\ {L}_{{\rm{p}}}$ (erg s−1)	tp (days)	Tail MNi (M⊙)	T0 (days)	Arnett MNi (M⊙)	Calibrated β	KK19 MNi (M⊙) a	f
SN 1993J	IIb	V − I	42.41 (41.56)	22.0 (0.0)	0.081 (0.005)	110.9 (9.9)	0.15 (0.02)	0.80 (0.16)	0.80 (0.16)	0.23
SN 1994I	Ic	V − I	42.62 (41.56)	8.2 (0.3)	0.048 (0.015)	57.0 (10.3)	0.11 (0.01)	0.00 (0.00)	0.080 (0.007)	0.48
SN 1996cb	IIb	V − R	41.87 (41.12)	21.7 (2.0)	0.030 (0.003)	96.8 (14.6)	0.04 (0.01)	1.17 (0.27)	0.023 (0.004)	−0.01
SN 1998bw	Ic-BL	V − I	43.16 (42.38)	16.6 (0.0)	0.300 (0.013)	108.9 (14.5)	0.67 (0.11)	0.50 (0.19)	0.312 (0.052)	0.38
SN 2002ap	Ic-BL	V − I	42.47 (41.60)	13.1 (2.5)	0.062 (0.003)	113.4 (12.3)	0.11 (0.02)	0.64 (0.20)	0.058 (0.008)	0.27
SN 2003jd	Ic-BL	V − R	42.85 (42.11)	15.0 (2.1)	0.117 (0.014)	101.4 (16.1)	0.29 (0.06)	0.29 (0.21)	0.145 (0.026)	0.47
SN 2004aw	Ic	V − I	42.71 (41.93)	15.6 (6.0)	0.151 (0.030)	150.1 (92.0)	0.22 (0.07)	1.03 (0.26)	0.133 (0.022)	0.08
SN 2004gq	Ib	B − V	42.36 (41.67)	12.3 (4.0)	0.044 (0.008)	144.6 (74.2)	0.08 (0.03)	0.59 (0.32)	0.051 (0.011)	0.29
SN 2005hg	Ib	V − R	43.09 (42.34)	18.8 (2.0)	0.353 (0.039)	143.2 (19.2)	0.63 (0.13)	0.83 (0.23)	0.348 (0.063)	0.21
SN 2006T	IIb	B − V	42.63 (42.00)	15.5 (1.5)	0.082 (0.012)	143.6 (29.5)	0.18 (0.05)	0.49 (0.30)	0.105 (0.025)	0.39
SN 2006el	IIb	V − r2	42.27 (41.63)	21.4 (0.7)	0.052 (0.005)	126.4 (21.2)	0.11 (0.02)	0.70 (0.25)	0.057 (0.013)	0.31
SN 2006ep	Ib	B − V	42.55 (41.83)	15.0 (7.0)	0.058 (0.015)	112.2 (31.1)	0.16 (0.07)	0.29 (0.22)	0.088 (0.017)	0.47
SN 2007gr	Ic	V − I	42.38 (41.60)	9.8 (2.5)	0.047 (0.004)	126.9 (15.7)	0.07 (0.02)	0.77 (0.33)	0.049 (0.008)	0.18
SN 2007ru	Ic-BL	V − I	42.90 (42.13)	11.1(3.0)	0.103 (0.009)	107.7 (14.9)	0.27 (0.07)	0.09 (0.13)	0.147 (0.025)	0.50
SN 2007uy	Ib	B − V	42.80 (42.16)	16.1 (3.5)	0.143 (0.024)	119.9 (25.4)	0.29 (0.09)	0.64 (0.29)	0.166 (0.037)	0.31
SN 2008D	Ib	V − r2	42.08 (41.32)	19.9 (0.2)	0.036 (0.003)	104.6 (14.1)	0.07 (0.01)	0.82 (0.21)	0.036 (0.006)	0.22
SN 2008ax	IIb	V − R	42.31 (41.58)	21.8 (1.2)	0.060 (0.009)	99.2 (17.2)	0.12 (0.02)	0.73 (0.21)	0.063 (0.012)	0.29
SN 2009bb	Ic-BL	V − I	42.78 (42.01)	11.0 (1.1)	0.158 (0.055)	47.8 (14.5)	0.19 (0.04)	1.23 (0.35)	0.109 (0.019)	−0.04
SN 2009jf	Ib	V − I	42.68 (41.91)	22.3 (4.2)	0.164 (0.022)	158.3 (22.7)	0.29 (0.07)	0.87 (0.20)	0.155 (0.026)	0.19
SN 2011bm	Ic	V − I	42.86 (42.09)	34.6 (0.7)	0.574 (0.038)	119.3 (90.3)	0.63 (0.11)	1.71 (0.41)	0.336 (0.057)	−0.33
SN 2011dh	IIb	V − R	42.49 (41.48)	20.1 (0.0)	0.093 (0.002)	102.0 (8.6)	0.17 (0.02)	0.81 (0.12)	0.090 (0.009)	0.22
iPTF13bvn	Ib	V − I	42.31 (41.54)	16.7 (0.8)	0.040 (0.006)	134.9 (15.6)	0.10 (0.02)	0.43 (0.19)	0.055 (0.009)	0.42
SN 2013df	IIb	V − I	42.37 (41.53)	22.6 (1.6)	0.063 (0.004)	128.2 (12.6)	0.14 (0.02)	0.61 (0.14)	0.074 (0.011)	0.37
SN 2013ge	Ib	B − V	42.39 (41.78)	16.8 (4.7)	0.063 (0.012)	122.8 (31.3)	0.11 (0.04)	0.79 (0.33)	0.065 (0.016)	0.24
SN 2014ad	Ic-BL	V − I	42.84 (42.07)	16.4 (3.0)	0.147 (0.023)	96.7 (15.6)	0.32 (0.08)	0.54 (0.20)	0.148 (0.025)	0.36
SN 2016coi	Ic-BL	V − I	42.84 (42.06)	17.6 (2.1)	0.170 (0.016)	135.7 (18.3)	0.34 (0.07)	0.67 (0.20)	0.153 (0.026)	0.30
SN 2016gkg	IIb	V − I	42.29 (41.51)	16.5 (4.0)	0.055 (0.004)	124.8 (15.5)	0.09 (0.02)	0.91 (0.24)	0.050 (0.008)	0.15

Notes.

^a M_Ni values obtained from KK19's model assuming the median values of calibrated β reported in Table 4 (see Section 4.4). ^b Since the BCs provided by Lyman et al. (2014) are given in Johnson bands, magnitudes in Sloan r are first converted to Johnson R using the transformations of Jordi et al. (2006).

Download table as:  ASCIITypeset image

Figure 1 illustrates this process. It presents the absolute magnitude (top panel) and the resulting bolometric luminosity (bottom panel) light curves for two SNe in our sample: SN 2008ax and SN 2007ru. These two objects were specifically chosen to span the range of light-curve coverage available during the early rise and late-time tail for SNe in our sample. The spline and linear fits to the absolute magnitude light curves are shown in the top panel. The x-axis represents time since the inferred epoch of explosion derived in Section 3.4. Gray regions indicate the uncertainties in the bolometric luminosity obtained at each epoch. These uncertainties stem from (in decreasing order of importance) error in the distance estimate, BC error, and photometric error. We also mark the epoch of explosion t₀, the peak luminosity L_p, and peak time t_p of the bolometric light curve in the bottom panel, with shaded regions representing the uncertainty on each parameter.

Zoom In Zoom Out Reset image size

Derived values of L_p and t_p for each SN are listed in Table 2 and plotted in Figure 2. Note that the final error in t_p is a combination of the error in t₀ and the bolometric light curve. Overall, the peak luminosities for our sample span 10^41.87–10^43.09 erg s⁻¹ and rise times span 8.2–34.6 days (SN 1994I and SN 2011bm, respectively). As seen in Figure 2, no correlation is apparent between L_p and t_p for the SESNe in our sample.

Zoom In Zoom Out Reset image size

Comparing our L_p estimates with previous studies, we find that our derived peak luminosities are a factor ∼2 larger than those found by Meza & Anderson (2020), who integrate luminosity over the BVRIJH bands, but they are in good agreement with those of Prentice et al. (2016), who find peak luminosities by summing over the UBVRIJHK bands for a fraction of their SN sample. In addition, our L_p estimates are consistent within our quoted errors with those of Lyman et al. (2016), who also adopt the BC polynomial fits of Lyman et al. (2014).

To constrain M_Ni of our SESNe sample, we model their radioactive light-curve tails using the analytic model of Wygoda et al. (2019) discussed in Section 2. This model is similar to those of Valenti et al. (2008) and Drout et al. (2013) with one minor modification: the positrons' escape is neglected. Since positrons escape on a timescale of a few thousand days, it should not affect our results.

By fitting the bolometric luminosity and slope of the radioactive tail for the SESN sample, we can constrain the two unknown parameters in Equation (1): the nickel mass, M_Ni, and partial trapping timescale of the tail, T₀. The fitting is done in an MC fashion: we run 1000 trials drawing from the distribution of possible luminosities and epochs of explosion. This allows us to propagate the uncertainties in these quantities when obtaining M_Ni and T₀. We only consider epochs of ≥60 days post-explosion, when the ejecta of SESNe are expected to be optically thin, such that the bolometric luminosity will be set by the instantaneous heating rate. In Figure 1, we display the best-fit radioactive tail models for SN 2008ax and SN 2007ru (dotted–dashed red curves; bottom panels). As shown, the model closely matches the evolution of the bolometric radioactive tail. The best-fit parameters, "tail M_Ni" and T₀, for each SN are listed in Table 2. Our best-fit tail M_Ni values range from ∼0.030 to ∼0.574 M_⊙ with a median value of 0.08 M_⊙. T₀ is in the range ∼47.8–158.3 days with a median value of 116.6 days. The points in Figure 2 are color-coded based on these derived tail M_Ni values. Objects with larger M_Ni values also exhibit brighter peak luminosities, as expected for SNe powered predominantly by radioactive decay.

We report the basic statistics of our results (mean, median, and standard deviation) both for the full sample and separated by SESN subtype in Table 3. However, we note that our sample size is relatively small when SNe are categorized by subtypes, especially normal Type Ic SNe, for which only four events met all of our sample criteria outlined in Section 3.1 and whose distribution may be skewed by the extreme event SN 2011bm. Therefore, we conduct an analysis of variance (ANOVA) test to check whether the reported differences between the mean M_Ni of SN subtypes are statistically significant. The result of the ANOVA test indicates that the pairwise comparison of M_Ni between SN subtypes is not generally statistically significant. One exception is Type Ic-BL SNe, for which the reported mean M_Ni was found to be higher than that of the combined sample of all other SESN subtypes with p-value <0.05. This is consistent with previous studies, which have typically found systematically higher M_Ni for Type Ic-BL events (e.g., Drout et al. 2011; Lyman et al. 2016; Prentice et al. 2016).

Figure 3 displays cumulative distribution functions (CDFs) of the derived M_Ni for our sample of SESNe. In the top panel of Figure 3, the M_Ni CDFs are provided for the entire SESN sample obtained using multiple methods: radioactive tail modeling (blue), Arnett's rule (green), and the KK19 model (pink). (See Sections 4.2 and 4.5 for Arnett and KK19 methods, respectively.) In order to account for the errors in individual M_Ni measurements when plotting the CDFs, we run 1000 MC trials in which we sample each M_Ni value based on the distribution defined by its uncertainty and construct a new CDF. These sampled CDFs are overplotted in Figure 3, forming hatched regions that represent the uncertainties associated with the obtained CDFs.

Zoom In Zoom Out Reset image size

In the second panel of Figure 3 we present the CDFs of the tail M_Ni values for each SESN subtype separately. We conduct Kolmogorov–Smirnov (K-S) tests on the CDFs of tail M_Ni estimates for each subtype. A K-S test on the CDF of Type Ic-BL and Ib/c SNe rejects the null hypothesis that these SN types are drawn from the same groups of explosions with a p-value = 0.02. In contrast, we find that Type Ic and Ib SNe are likely to be drawn from the same distributions with p-value = 0.24. Similarly, a K-S test on Type IIb and Ib/c SNe supports the null hypothesis that these samples originate from the same distribution with p-value = 0.10. The p-value further increases to 0.45 when we compare the CDF of Type IIb and Ib SNe.

Radioactive tail nickel masses have previously been estimated for a number of SNe in our sample, and results are consistent. In particular, the tail M_Ni values reported in the recent work of Meza & Anderson (2020) are lower limits because they do not take the partial trapping of γ-rays into account. Our tail M_Ni estimates provided in Table 2 are consistent with these lower limits for SNe that are common to the two samples. Similarly, our M_Ni and T₀ estimates are consistent (within the margin of error) with those obtained from the Katz integral method (Sharon & Kushnir 2020) for eight SESNe in common between the samples.

For comparison, we also measure M_Ni using Arnett's rule, described in Section 2, which has been extensively employed in the literature. For each SN in our sample, we fit for M_Ni in Equation (4) assuming the t_p and L_p values listed in Table 2. Similar to the procedure of deriving tail M_Ni, we run 1000 MC trials to take into account the uncertainties in t_p and L_p when obtaining Arnett M_Ni values. Results for individual SNe are provided in Table 2, while basic statistics of both the full Arnett distribution and SESN subtypes are reported in Table 3. The Arnett M_Ni values span 0.04 M_⊙ to 0.67 M_⊙ with a mean value of 0.22 M_⊙. The third panel of Figure 3 displays Arnett M_Ni CDFs for different SESN subtypes.

Figure 4 (top panel) presents a comparison between the tail and Arnett M_Ni for our sample of SESNe. The results highlight the systematic discrepancy between the two methods. The M_Ni values obtained from the radioactive tail modeling are, on average, a factor of ∼2 smaller than those derived using Arnett's rule. The dashed black line indicates the equality condition between the two models. Despite the scatter in the severity of this discrepancy for different SNe, the Arnett model consistently overestimates M_Ni for every SESN in our sample. The overestimation of M_Ni by the Arnett model is also illustrated in Figure 3 (top panel), where the CDF of the Arnett M_Ni distribution is below that of the tail distribution by a relatively large margin.

Zoom In Zoom Out Reset image size

The means and standard deviations of our Arnett-derived M_Ni values for different SN types closely match the values reported in Lyman et al. (2016), but our median values are lower than those of Prentice et al. (2016) for Type Ib, Ic, and Ic-BL SNe by ∼0.04 M_⊙. This discrepancy is primarily due to different approaches in deriving t_p, which is estimated in Prentice et al. (2016) by measuring the rise time from the half-maximum luminosity to L_p (denoted by t_−1/2) and using a linear empirical correlation for translating t_−1/2 to t_p. We also note that the Arnett M_Ni values of Meza & Anderson (2020) are, on average, 50% lower than our Arnett estimates. This difference can be traced to their anomalously lower peak luminosities as discussed in Section 3.5 above.

The inaccuracy of M_Ni values obtained from Arnett models has also been shown in several radiative transfer numerical simulations of SESNe. For example, Dessart et al. (2015, 2016) found that Arnett's rule overestimated M_Ni of SESNe by 50% and attributed this discrepancy to the assumption of fixed electron scattering opacity in Arnett's models. Similarly, Sukhbold et al. (2016) pointed out that Arnett's rule does not hold for their simulations of Type Ib/c SNe evolved from massive single-star progenitors. We note that the discrepancy we find between our Arnett and tail M_Ni values is approximately twice as large as that quoted by Dessart et al. (2016).

Despite these limitations, Arnett models have been widely used for deriving M_Ni as well as ejecta masses and kinetic energies of SESNe (Drout et al. 2011; Lyman et al. 2016; Prentice et al. 2016, 2018). A few observational studies have also indicated contrasts between results from modeling the early and late-time light curves of SESNe. For example, Valenti et al. (2008) evoke a "two-zone" model to try to resolve an inconsistency between the explosion parameters derived from early and late-time light curves of SN 2003jd, which they fit with Arnett and radioactive models, respectively. In another study, Wheeler et al. (2015) estimate T₀ values for dozens of SESNe by an analytical relation that depends on M_ej and kinetic energy, which were rewritten in terms of the observed rise time and photospheric velocity, assuming Arnett's model. The estimated T₀ values were found to be in tension with the values measured directly from the light-curve tails, which likely is due to the limitations of Arnett's model. More recently, Meza & Anderson (2020) measured M_Ni values for a sample of SESNe using a variety of methods. While their tail M_Ni values are lower limits, they confirm that Arnett values are consistently higher than those derived via other methods.

For the rest of our analyses, we assume that our tail M_Ni estimates are more realistic than those of Arnett's model. This is because the ejecta is expected to be transparent to optical photons over the tail. Therefore, the bolometric luminosity traces the instantaneous heating rate without any further assumption regarding the self-similarity of the energy density profile, which is a fundamental assumption in the Arnett-like models (KK19).

As described Section 1, by comparing measurements of M_Ni for 115 H-rich Type II SNe and 145 SESNe previously published in the literature, Anderson (2019) identified a discrepancy in their distributions, with SESNe displaying a mean M_Ni that was a factor of ∼6 larger than that of their H-rich counterparts. Subsequently, Meza & Anderson (2020) computed M_Ni directly for a smaller sample of SESNe, using a number of methods (Arnett's rule, KK19, and radioactive tail modeling) in a uniform manner, demonstrating that a statistically significant discrepancy remains. However, the tail M_Ni values presented in Meza & Anderson (2020) are strictly lower limits to the true M_Ni, because they assume complete trapping of γ-rays, while the Arnett and KK19-based values may each contain systematic biases (in the latter case because they adopt β values that have yet to be observationally calibrated; see Section 4.4). Thus, while sufficient to robustly demonstrate that SESNe have a different M_Ni distribution than Type II SNe, the magnitude of this discrepancy remains somewhat uncertain.

Here, we compare the CDF of M_Ni for SESNe derived from our radioactive tail measurements to that for the 115 H-rich Type II SNe from Anderson (2019). M_Ni for most of this sample of Type II SNe was calculated using the bolometric luminosity of the radioactive tail of the light curve assuming full trapping of the γ-rays. Figure 5 illustrates the M_Ni CDF of our sample of SESNe (orange curve) and the H-rich Type II SNe (green curve). For reference, we also show the original sample M_Ni values of SESNe compiled from the literature by Anderson (2019, pink curve)—most of which were obtained using Arnett's rule—and the lower limits of M_Ni from Meza & Anderson (2020, light blue curve). We see that our distribution of M_Ni for SESNe measured from the radioactive tail lies between that of the Arnett-based values of Anderson (2019) and the tail upper limits of Meza & Anderson (2020), as expected.

Zoom In Zoom Out Reset image size

Overall we find that our sample of SESNe have a mean value (∼0.12 M_⊙; Table 3) that is a factor of ∼3 larger than that of H-rich Type II SNe (∼0.044 M_⊙). This is a factor of ∼2 smaller than the initial discrepancy reported by Anderson (2019) based on Arnett measurements. We conduct K-S tests on the M_Ni CDFs of H-rich Type II SNe (Anderson 2019) and SESNe in this work. The test gives D-value = 0.52 and p-value = 10⁻⁷, meaning that the CDFs are inconsistent with being drawn from the same distribution. By excluding Type Ic-BL SNe from the test, we find D-value = 0.49 and p-value = 5 × 10⁻⁵, which similarly confirms that Type IIb/Ib/c SNe and H-rich Type II SNe are inconsistent with being drawn from the same distributions.

As in Meza & Anderson (2020) we find that a majority of this discrepancy come from the lack of SESNe in our sample with low M_Ni values. The lowest tail M_Ni in our sample is 0.03 M_⊙, while an incredible ∼48% of Type II SNe have M_Ni lower than this value. If we recompute the K-S tests described above, but consider only Type II SNe with M_Ni > 0.03 M_⊙, we find a p-value = 0.008 for the full sample of SESNe and p-value = 0.06 when Type Ic-BL SNe are excluded. This indicates that the sample of IIb/Ib/Ic SNe are marginally consistent with being drawn from the same population as the high-M_Ni Type II SNe.

In Section 4.2 we demonstrated that Arnett-based models yield M_Ni values that are a factor of 2 larger than those found from modeling the radioactive tail. While obtaining tail-based M_Ni measurements for all SESNe would be ideal, in practice the requisite photometric data exist for only a subset of events. Thus, another means to estimate M_Ni from photospheric data alone would be beneficial.

As discussed in Section 2, KK19 proposed an analytical model that relates the peak luminosity and its epoch to a general heating function without relying on some of the simplifying assumptions adopted by Arnett's models. This new model, described in Equation (6), depends on a dimensionless parameter β in addition to M_Ni, t_p, and L_p. KK19 suggested β = 9/8 for Type Ib/c SNe based on the radiative transfer simulations of SESN light curves from Dessart et al. (2016) and β = 0.82 for Type IIb/pec SNe based on the observed light curve of SN 1987A. However, numerical simulations may not fully represent the behavior of real SESNe, and the light curve of SN 1987A is very different than that of Type IIb SNe. More reliable constraints on β can be obtained from the observed sample of SESNe with independent M_Ni values measured from their radioactive tails.

We use Equation (6) to calculate the value of β inferred for each SESN in our sample given its tail M_Ni, t_p, and L_p provided in Table 2. The uncertainty in the derived β values has contributions from the error in tail M_Ni, L_p, and t_p. In Figure 6 we plot the derived β values versus tail nickel mass for individual SNe. The results exhibit a significant scatter, with β ranging from ≈0.0 to 1.7 and a mean value of 0.70. This is lower, on average, than the β = 9/8 suggested by KK19 for SESNe based on the models of Dessart et al. (2016). This comparison will be examined in more detail in Section 4.6. The horizontal lines in Figure 6 indicate the median β value for each SESN subtype. Table 4 provides summary information on the mean, median, and standard deviation of the derived β values for each SN subtype. The median values of β for Type IIb, Ib, and Ic SNe are roughly similar, but the standard deviation of Type Ic SNe is a factor ∼3 larger due to the small sample size of objects of this type. Type Ic-BL SNe have the smallest median β value of 0.54, which is ≈30% smaller than that of other SESN types.

Zoom In Zoom Out Reset image size

Two SNe in our sample, SN 1994I and SN 2007ru, have β close to zero, which may suggest that the derived tail M_Ni is inadequate to produce the observed peak luminosity, L_p. Recall that for a fixed M_Ni and t_p, a lower β value will yield a higher peak luminosity (see Section 2.3). Both of these SNe are fast declining and will be discussed further in Section 5.

Using the results from Section 4.4 we now assess whether, in practice, the model of KK19 can be used to obtain more reliable M_Ni estimates than Arnett for SESNe from photospheric data alone. In Table 2 we list M_Ni values for each SN that have been calculated using their observed L_p and t_p in conjunction with the median β value for each SN subtype listed in Table 4. Errors listed in Table 2 account only for the errors in L_p and t_p and do not include the impact of the standard deviation in the distribution of β values. Despite the significant scatter in β values found for individual SNe, we find that the procedure of calculating M_Ni assuming the KK19 model and the median β value for each SN subtype offers a significant improvement over Arnett-based measurements, both for the overall distribution of M_Ni values for SESNe and for individual objects. These two effects are demonstrated in Figures 3 and 4, respectively.

In the top panel of Figure 3 we present the CDF of M_Ni calculated using KK19 with the median calibrated β in comparison to those of the radioactive tail and Arnett methods. As shown, not only is the mean value of the KK19 M_Ni distribution the same as that from the M_Ni measured with (shown by a vertical lines; see also Table 3), but the overall CDF of KK19-measured M_Ni values closely approximates that of the M_Ni measurements from the radioactive tail. In the bottom panel of Figure 3 we also display the CDFs of KK19 M_Ni estimates separated by SN subtype. As listed in Table 3, the median values of these distributions are all within ≈10% of those calculated from the radioactive tail.

In Figure 4 we demonstrate that this agreement extends to individual objects. The KK19-measured M_Ni values for each SN (bottom panel) are much closer to their tail counterparts than the Arnett-derived ones (top panel). We find that, on average, the KK19 M_Ni values are within ∼17% of the tail-derived values. The largest variations occur for the Type Ic SN 1994I and SN 2011bm, for which the nickel mass is overestimated by ∼65% and underestimated by ∼40%, respectively. In contrast, the Arnett-based models systematically overpredict M_Ni by a factor of ∼2 (100%) compared to the tail-derived values. We therefore conclude that in cases where the radioactive tail is not observed, the KK19 model with median calibrated β values listed in Table 4 should be used to calculate M_Ni for SESNe.

In addition to providing improved estimates for M_Ni from photospheric data alone, the β values calculated for our observed SESNe encode information on the explosion properties and progenitors of the population. Features such as composition, asymmetry, and additional power sources will impact the degree to which the internal energy of the ejecta lags or leads the observed luminosity at the time of peak. In order to assess whether the observed population of SESNe matches expectations from theory, and to gain insight into the physical processes that dictate β in observed events, we calculate the β values for a set of analytical and numerical light-curve models available in the literature.

4.6.1. Arnett Models

First, as a baseline, we derive β for a grid of light curves calculated using an analytic Arnett model. We take t_p and M_Ni in the ranges 5–40 days and 0.02–0.5 M_⊙, respectively, which correspond to the approximate ranges found for our observed sample. Given a pair of t_p and M_Ni, Arnett's rule gives L_p; we then compute β from Equation (6) using t_p, L_p, and M_Ni. The resulting values of β for this grid are plotted as a gray shaded region in Figure 7. As expected given the inconsistency between Arnett and tail-derived nickel masses shown above, these β values are inconsistent with those of our observed population. In particular, the Arnett models occupy a parameter space with high β values in the range 1.55–1.95, while all but one observed SESNe (SN 2011bm) has a calibrated β < 1.25.

Zoom In Zoom Out Reset image size

4.6.2. Large Grids of Numerical SESNe

4.6.3. Specialized SESN Models

In Section 4.6, we presented the β values inferred for different numerical models of SESNe. The results show that most numerical models give higher β values than observations. In order to disentangle the origin of this discrepancy, in Figure 8 we present the pairwise dependence and histograms of the three quantities that determine β, i.e., L_p, t_p, and M_Ni, for both the numerical models detailed in Section 4.6 and observed SESNe. This figure illustrates that:

• 1.  
  while both the observed sample and model SESNe show a correlation between M_Ni and L_p, the observed objects exhibit considerably (∼0.3–0.4 dex) larger peak luminosities for a given M_Ni;
• 2.  
  the observed SESNe tend to have shorter rise times than the models. The median rise time for the entire sample is 16.6 days, as compared to 19.5 and 26.8 days for the Ertl et al. (2020) and Dessart et al. (2016) models, respectively.

Zoom In Zoom Out Reset image size

As described in Section 2.3 both rise time and peak luminosity are inversely proportional to β. Thus, shorter rise times would primarily act to increase β. Indeed, it appears that the main driver between the different β values of the Ertl et al. (2020) and Dessart et al. (2016) models is that Ertl et al. (2020) find shorter rise times for a given M_Ni. This effect was discussed in Ertl et al. (2020) and is primarily due to their adoption of a constant line opacity. In contrast, the higher luminosity for a given M_Ni would act to lower β, and this is therefore likely the origin of the discrepancy in β displayed in Figure 7. Here we investigate the possible physical origins of this discrepancy, as well as the scatter in effective β displayed by the observed sample.

5.1.1. Systematic Effects due to Distance and Reddening

5.1.2. Dark Period

5.1.3. Composition, Opacity, and Recombination

5.1.4. Mixing of Radioactive Material

5.1.5. Asymmetry

There is growing evidence from a combination of spectropolarimetry, nebular spectroscopy, and resolved SN remnants that some SESNe may be asymmetric (e.g., Valenti et al. 2011; Milisavljevic & Fesen 2015; Tanaka 2017). As shown in Figures 7 and 8, the β values of the mildly asymmetric simulations of Barnes et al. (2018) depend on the observer's viewing angle. This is primarily due to a variation in the observed peak luminosity, with viewing angles with larger L_p leading to lower β values. Thus, depending on their nature, asymmetries in ⁵⁶Ni mixing or ejecta distribution can be reflected in both the mean value and scatter in the β parameters observed for a population of SESNe.

To test the degree to which asymmetry can modify observed β values for a population, we run a set of light-curve simulations with varying degrees of ejecta asymmetry using the multidimensional radiative transfer code Sedona (Kasen et al. 2006). We assume an axisymmetric homologously expanding ejecta profile consisting of a broken power law in density following Chevalier & Soker (1989) and Kasen et al. (2016). To account for deviations from spherical symmetry, we vary the semimajor axis as in Darbha & Kasen (2020) parameterized by A ≡ v_r /v_z , where v_r and v_z are the outer ejecta velocities at the equator and pole, respectively. In total, we run four radiative transfer simulations with A ≤ 1, i.e., prolate ejecta configurations.

We choose fiducial values of M_ej = 2 M_⊙ and v _z = 10⁴ km s⁻¹ for our simulations. To account for the heating, we set the innermost 0.1 M_⊙ to consist of ⁵⁶Ni. Finally, we assume a constant gray opacity of κ = 0.2 cm² g⁻¹. The resulting light curve is then calculated at 10 different viewing angles $\mu =\cos \theta$, in the range μ ∈ {–1, 1}, where μ = 0 and μ = ±1 view the ejecta along the equator and poles, respectively.

We measure the different values of L_p(μ, A) and t_p(μ, A) for the output bolometric light curves, which we then map onto an inferred β based on Equation (6) and the model parameter M_Ni = 0.1 M_⊙. In Figure 9, we show the inferred values of β for different sets of asymmetric ejecta configurations and viewing angles. As expected, β does not vary with viewing angle for the case A = 1, i.e., no ejecta asymmetry. However, increasing the degree of asymmetry results in lower inferred β when viewed along the equator, and higher β when viewed at the poles. This is due to L_p being larger when viewed along the equator, where the projected surface area is largest; similarly, L_p is decreased along the poles for asymmetric ejecta due to a smaller projected surface area (Darbha & Kasen 2020). This is in agreement with the results found in Barnes et al. (2018).

Zoom In Zoom Out Reset image size

While the total spread in β values observed for our asymmetric simulations is ≳1—well matched to the scatter in our observed population—when averaged over all viewing angles, we find that asymmetry acts to increase the mean inferred value of β. This is opposite to the direction that has been observed, where the average β is systematically lower than the spherically symmetric models of Dessart et al. (2016) and Ertl et al. (2020) (as well as the A = 1 model run in this work). Thus, insofar as asymmetry can be represented by an expanding broken power-law ellipsoid, we conclude that asymmetry, although possibly explaining some of the scatter seen in Figure 7, cannot explain the systematically lower inferred values of β. However, we note that when asymmetry is strong other effects such as the development of nonradial flows (Matzner et al. 2013; Afsariardchi & Matzner 2018) and the ejection of nickel-rich clumps to high velocities (Drout et al. 2016) could further influence the light-curve morphology.

5.1.6. Additional Power Sources

In Section 5.1 we demonstrated that the main driver of the low β values for our observed SESNe was their high L_p for a given M_Ni (Figure 8). Thus, another plausible explanation for the origin of the discrepancy between the β distribution of the models and the observed SESNe is that additional power sources beyond the radioactive decay of ⁵⁶Ni contribute to the peak luminosity of SESNe. This was previously proposed by Ertl et al. (2020), who noted that their numerical simulations were unable to reproduce the brighter half of observed Type Ib/c SNe luminosity function (Figure 8). In addition, when modeling a sample of SESNe with the luminosity integral method of Katz et al. (2013), Sharon & Kushnir (2020) required an additional model parameter, which they interpret as a non-negligible amount of emission produced by power sources beyond ⁵⁶Ni. The presence of an additional luminosity source near peak could also explain why the observed SESNe show an even larger discrepancy between Arnett and tail-measured M_Ni than the theoretical models of Dessart et al. (2016).

Within our sample, the need for additional power sources is particularly conceivable for SN 1994I and SN 2007ru, for which we derived a negative β and β close to zero, respectively. In practice, a negative β is not physical in the general definition given by KK19 (see Equation (5)). Rather, this implies that the version of this equation that assumes pure ⁵⁶Ni heating (Equation (6)) is incapable of producing such a bright peak luminosity in the observed rise time when coupled with the M_Ni measured from the radioactive tail. While these are the most extreme cases, other observed SESNe may also have extra power sources contributing to the observed luminosity. In this case, attributing the observed peak luminosity solely to the radioactive decay of ⁵⁶Ni will result in smaller β values than expected based on their composition and ⁵⁶Ni distribution.

Besides the radioactive decay of ⁵⁶Ni, power sources that can contribute to the light curves of CCSNe include shock cooling emission, ejecta interaction with CSM, and energy from a central engine. Here, we assess the viability of each of these sources in explaining the discrepancy between observed and model SESNe, and implications thereof.

Luminosity Required. We begin by evaluating how much excess luminosity would be required to reconcile the β values we derive in Section 4.4 with those of the theoretical models in Section 4.6. To do so, we assume that the average value of β = 1.125 found by the models of Dessart et al. (2016) is accurate for SESNe powered only by radioactive decay. Given their full treatment of radiative transfer/opacity, this is equivalent to assuming that the composition, energetics, and mixing of radioactive material included in Dessart et al. (2016) reflect reality. Then, using β = 1.125, the M_Ni and t_p listed in Table 2, and Equation (6), we calculate the luminosity that ⁵⁶Ni decay can produce at the peak time of the SN. From this, we define an "excess luminosity factor," f, as the fraction of peak luminosity that is in excess over what is expected from a β = 1.125 explosion with a nickel mass as defined by the light-curve tail.

The resulting values of f for each SN are listed in Table 2 and the mean, median, and standard deviation in Table 4. In Figure 10 we plot the excess luminosity factor, f, versus the excess luminosity, f × L_p. We emphasize that these "excess luminosity" values should be taken as order-of-magnitude estimates only, because they (a) neglect any variations from β = 1.125 that can be caused by effects such as enhanced nickel mixing and asymmetry and (b) inherently assume that the radioactive component peaks at the same time as the observed light curve; depending on the nature of the excess luminosity source, this need not be the case. Nevertheless, they provide a useful diagnostic.

Zoom In Zoom Out Reset image size

Overall, we find f values in the range from −0.33 to 0.5. Three events (SN 1996cb, SN 2009bb, and SN 2010bm) have negative f values, reflecting the fact that they had measured β > 1.125, and thus do not require additional power sources. For the remaining 24 events, we find that if the models of Dessart et al. (2016) and our tail M_Ni are accurate, we require that ∼7%–50% of the their peak luminosities come from sources other than ⁵⁶Ni. This translates into considerable excess luminosity in the range from 2.5 × 10⁴¹ to 5.5 × 10⁴² erg s⁻¹ or, equivalently, peak bolometric magnitudes for the excess power source of −14.9 mag > M_bol > −18.2 mag. Out of the six events that would require the most excess luminosity, five are Type Ic-BL, while the events that require the highest fraction of the peak luminosity to come from other sources are mixed between subtypes.

Shock Cooling and CSM Interaction. Both shock cooling and CSM interaction face challenges in explaining this required excess emission. First, shock cooling emission is predicted to be both faint and short-lived for the compact progenitors (R ∼ R_⊙) typically evoked for H-poor SESNe (Nakar & Sari 2010; Rabinak & Waxman 2011). Second, the substantial luminosity required and lack of narrow emission lines in the majority of SESNe spectra limit any potential CSM to dense and confined shell or disk-like configurations (Chevalier & Irwin 2011; Moriya & Tominaga 2012; Smith et al. 2015), which are not predicted from standard models of stellar evolution (Smith 2014). However, recent progress may change this picture: models show that particularly low-mass He stars (≲ 3 M_⊙) can undergo substantial inflation (achieving radii ≳ 100 R_⊙; e.g., Yoon et al. 2010; Kleiser et al. 2018a; Woosley 2019), and there is increasing observational evidence that many CCSNe—of all varieties—undergo a period of enhanced mass loss shortly before core-collapse (e.g Margutti et al. 2015; Drout et al. 2016; Bruch et al. 2021). Thus, a significant fraction of SESN progenitors may have large effective radii (≳20 R_⊙), in which case shock-deposited energy can contribute substantially to their observed luminosity.

As described in Section 4.6.3, Kleiser et al. (2018a, 2018b) model the light curves that would result solely from the diffusion of shock-deposited energy in both of these scenarios. The set of O/He-rich CSM shells modeled in Kleiser et al. (2018b) have masses of 1−4 M_⊙ and are located at radii of ∼15–60 R_⊙, designed to mimic an intense final mass-loss episode. They find transients that rise to peak magnitudes of −15 mag < M_r < −18 mag on timescales of ∼8–25 days—parameters that are exceptionally well matched to our estimates above for the excess luminosity required. Indeed, as shown in Figure 7, when a small amount of ⁵⁶Ni is added, the models of Kleiser et al. (2018a, 2018b) are able to reproduce the low β values of the observed population.

While the particular theoretical models plotted in Figure 7 show slight double-peaked morphology in their optical light curves—which are not typically observed—they were created by adding ⁵⁶Ni to a shock cooling curve that reaches −17 mag and contributes ≳50% of the luminosity near peak. In contrast, we expect the SESNe in our sample to still be dominated by the radioactive decay of ⁵⁶Ni, with a median excess luminosity factor of f = 0.29. In addition, a number of Type Ib/c SNe show evidence for multiple emission components when observed early in the u-band/UV, with SN 2013ge also exhibiting narrow high-velocity absorption features, which may be indicative of a shell-like CSM structure (Drout et al. 2016; Kleiser et al. 2018b). However, the lack of prominent double-peaked structure in many Type Ib/c SNe combined with the need for excess luminosity near maximum light implies that: (i) the He/O-rich material that leads to the large effective radius cannot be too tenuous, in which case it would become optically thin on a timescale of days (Dessart et al. 2018; Woosley 2019), and (ii) the ⁵⁶Ni synthesized in the explosion must be well mixed to avoid a significant dark period and subsequently large offset in time between the two emission components (Kleiser et al. 2018b).

We therefore conclude that shock cooling emission is a viable source for the excess luminosity required to reconcile SESN observations and models. While low-mass He stars can reach the radii required through "natural" stellar evolution processes, such events are predicted to eject very little radioactive ⁵⁶Ni (Kleiser et al. 2018a; Woosley 2019). Thus, reproducing the full observed population via this mechanism likely requires a significant number of SESNe to undergo intense late-stage mass loss due to instabilities (Smith & Arnett 2014; Woosley 2019), internal gravity waves (Quataert & Shiode 2012; Fuller & Ro 2018), or some other physical process during the final nuclear burning stages. We emphasize that this shock cooling emission will rapidly fade once the ejecta cools to the recombination temperature of O/He-rich material (t < 40 days), and thus M_Ni measured at >60 days post-explosion should be unaffected.

Central Engines. On the other hand, a central engine could provide the additional energy source required to power the light curve around peak. We note, in particular, that the we found a lower mean β value and require a higher mean excess luminosity for Type Ic-BL SNe (Table 4), for which central engines are commonly evoked. Both magnetars and collapsars can provide a natural source of extra luminosity, in the form of spin-down energy and fallback accretion (Dexter & Kasen 2013), respectively. Ertl et al. (2020) investigate the former as a power source for SESNe in general, arguing that even the formation of a more moderate millisecond pulsar (e.g., Yoon et al. 2010) may be sufficient for rotational energy to impact the SN light curve without impacting the overall energetics.

We use the magnetar spin-down energy and timescale of Kasen (2017) together with SN timescale and ejecta velocity of Kasen et al. (2016) to estimate the general phase space of magnetar properties required to provide the excess luminosity found above. In Figure 11 we plot rise time versus excess luminosity for our sample of SESNe. Also shown are lines that represent the peak luminosity and time achieved by magnetar models with a fixed period but varying magnetic field (solid) and fixed magnetic field but varying period (dotted). All models assume an opacity of κ = 0.1 cm² g⁻¹, explosion energy of 10⁵¹ erg, and were calculated for both M_ej = 2 M_⊙ (top panel) and M_ej = 5 M_⊙ (bottom panel)—chosen to span the range of SESNe ejecta masses from Lyman et al. (2016). For M_ej = 2 M_⊙, P values of 10–116 ms and B₁₄ (= B/10¹⁴ G) values of 7–59 yield luminosities consistent with requirements. For M_ej = 5 M_⊙, similar magnetic field strengths but shorter periods are necessary with P = 2–32 ms and B₁₄ = 10–60 spanning the phase space occupied by most SESNe in our sample. The ranges of P and B₁₄ found here overlap with—but extend to longer periods and higher magnetic field strengths than—those of magnetar fits to superluminous SN light curves by Nicholl et al. (2017).

Zoom In Zoom Out Reset image size

Figure 12 displays two representative magnetar models for the excess peak luminosity of SN 2008ax. The excess emission factor is 0.29, corresponding to a luminosity of 5.90 × 10⁴¹ erg s⁻¹ (dashed line). The blue curve illustrates a model with B₁₄ = 36, P = 32 ms, and M_ej = 5 M_⊙, while the orange curve represents a magnetar with B₁₄ = 19, P = 116 ms, and M_ej = 2 M_⊙. Both models peak at the level of the excess luminosity calculated above, but they have different magnetar properties and light-curve morphologies. Figure 12 highlights that, unlike shock cooling emission, if a magnetar contributes to the light curve near peak it can also power a non-negligible portion of the tail luminosity. In this case, tail-based measurement of M_Ni would be overestimated, and as a result the fraction of the peak luminosity that must come from sources other than radioactive decay would be underestimated. While we find that fits over 60–120 days post-explosion, as performed in Section 3.5, are not sufficient to distinguish between the L(t) ∝ t⁻² power-law decline of the magnetar model and the exponential form of the radioactive decay of ⁵⁶Ni (Equation (1)) future modeling of a subset of SESNe with data covering a longer baseline (≳ 200 days) could potentially break this degeneracy. If the slope of the tail is shallower than the ⁵⁶Co → ⁵⁶Fe conversion rate this could also be an indication that other power sources are contributing to the light-curve tail (e.g., Ertl et al. 2020). However, none of the SESNe in our sample have such a shallow tail.

Zoom In Zoom Out Reset image size

5.1.7. Possible Limitations of Tail M_Ni Values

The discrepancy in ⁵⁶Ni mass distributions for SESNe and Type II SNe was first identified using Arnett-based measurements of M_Ni for SESNe (Anderson 2019). In Section 4.2 we demonstrated that Arnett's rule overpredicts M_Ni for SESNe by roughly a factor of 2—likely due both to limitations in Arnett's model (KK19) and to possible contributions from additional power sources to the peak luminosity of SESNe (Section 5.1.6). However, in Section 4.3 we find that even when tail-based M_Ni values are used, our population of SESNe have a mean M_Ni value (0.12 M_⊙) that is a factor of ∼3 higher than Type II SNe (0.044 M_⊙). Therefore a critical question remains: what makes the M_Ni distribution of SESNe skew to larger values than that of Type II SNe? Are the M_Ni values for the observed population of SESNe biased? If not, does the discrepancy with Type II SNe come about because SESNe originate from stars of higher ZAMS mass? Or are the ZAMS masses of SESNe and H-rich Type II SNe somewhat similar, but other physical mechanisms are responsible for the observed difference in the M_Ni distribution? Here, we consider each of these questions in turn.

Does the distribution of M_Ni derived accurately represent the true distribution of M_Ni synthesized in the core-collapse of H-poor stars? As highlighted by Meza & Anderson (2020), the discrepancy between SESNe and Type II SNe is primarily due to a lack of observed SESNe with low M_Ni. While the light curves of Type II SNe—powered predominantly by H recombination—remain luminous for ≳100 days regardless of M_Ni, low-M_Ni SESNe may be faint and/or rapidly evolving. In particular, as described in Section 5.1.6, low-mass stripped He stars (M ≲ 3–4 M_⊙; corresponding to ZAMS ≲12–15 M_⊙) can inflate to large radii prior to core-collapse. Stars with these initial masses likely dominate the population of Type II SNe, due to the initial mass function (IMF). However, the explosion of their stripped counterparts may be dominated by bright and rapidly fading cooling envelope emission (Dessart et al. 2018; Kleiser et al. 2018a; Woosley 2019), with minimal contributions from ⁵⁶Ni. Observationally, these could manifest as rapidly fading Type I SNe (Kasliwal et al. 2010; Drout et al. 2013), Type Ibn SNe (Hosseinzadeh et al. 2017), or the broader class of fast blue optical transients (Drout et al. 2014) rather than "traditional" SESNe. Such events were not explicitly included in our sample and have rarely been followed to late enough times to constrain the M_Ni ejected. While the rates of rapidly evolving transients (∼4%–7% of the CCSN rate; Drout et al. 2014) are not sufficient to fully resolve the discrepancy, they could reduce its significance.

Alternatively, the M_Ni distribution shown in Figure 5 may not accurately represent reality if the light curves of SESNe are not solely powered by the radioactive decay of ⁵⁶Ni and the timescale of the additional power source(s) is ≳60 days (see Section 5.1.6). In this case, attributing the full tail luminosity to ⁵⁶Ni would lead to an overestimate of the true values. Within this context, it is notable that while simulations of neutrino-driven core-collapse SNe (Sukhbold et al. 2016; Ertl et al. 2020) are able to self-consistently produce the range of M_Ni values observed in H-rich Type II SNe (∼0.004–0.13 M_⊙; Müller et al. 2017; Afsariardchi et al. 2019), they are unable to produce M_Ni values as high as those derived for the upper ∼30% of our SESN sample (see Figures 7 and 8). It was for this reason that Ertl et al. (2020) invoked magnetars to explain the observed light curves of SESNe. While this may resolve the M_Ni discrepancy, it subsequently requires that magnetars influence SESNe at a higher rate than Type II SNe. This may be a natural consequence of stripped stars retaining a larger fraction of their angular momentum, which in H-rich stars is lost primarily due to rotational braking when the star expands to the RSG phase (e.g., Ertl et al. 2020).

Do SESNe originate from stars of higher ZAMS mass than Type II SNe? If observational biases cannot explain the relative lack of low-M_Ni SESNe, the discrepancy could indicate that SESNe preferentially form from stars of higher ZAMS mass—which are expected to synthesize higher amounts of M_Ni—than Type II SNe. While at face value this may favor the single-star progenitor channel, this is in tension with the observed occurrence rate of SESNe, which is much higher than the explosion rate of stars with ZAMS mass ≳25 M_⊙ that eventually become W-R stars (Smith et al. 2011). In addition, the mass-loss rate of massive stars is still uncertain, and it is not clear whether all massive stars with ZAMS mass ≳25 M_⊙ can strip their H envelope (Smith 2014). Another possible explanation is that SESNe do form stripped binaries and span the same overall ZAMS masses as Type II SNe. However, while the relative rates of Type II SNe are dominated by the IMF, SESNe are skewed to higher relative masses. This may be a natural consequence of the close binary fraction being larger for high-mass stars (Moe & Di Stefano 2017; Moe et al. 2019), although detailed population synthesis calculations are necessary to test this hypothesis. We note that the population of SESNe coming from a distribution skewed to higher ZAMS masses, and hence shorter lifetimes, would be consistent with the fact that they are found closer, on average, to sites of active star formation than Type II SNe (Anderson et al. 2012).

Do additional physical mechanisms modify the M_Ni synthesized in SESNe? Alternatively, if the M_Ni distribution for SESNe is accurate, and SESNe originate from similar ZAMS masses to H-rich Type II SNe, then there must be some physical processes that make the M_Ni of SESNe larger. While it is often assumed stripping of the H and even He envelope via binary interaction should leave the inner core structure of the primary star intact (e.g., Fryer & Kalogera 2001)—in which case the M_Ni distributions of SESNe and Type II SNe should be indistinguishable—this picture may not be complete. In particular:

• 1.  
  In close binaries, fast orbital rotation, tides, magnetic braking, and angular momentum transport can influence the convective core sizes profoundly (e.g., Song et al. 2018). Therefore, changes in the core structure may impact the M_Ni production.
• 2.  
  A fraction of SESNe may originate from the merger of binary stars (Zapartas et al. 2017). In this case, a more massive core capable of producing a larger amount of ⁵⁶Ni may result. Although the rate of this merger channel seems to be relatively small (∼12% of SESNe; Zapartas et al. 2017), it will boost the overall M_Ni of SESNe produced via the binary channel.

Additional analysis is required to distinguish the contributions of each of the above scenarios toward explaining the discrepancy in observed M_Ni for H-poor and H-rich core-collapse SNe.