On the Energy Source of Ultrastripped Supernovae

The progenitor system of binary neutron stars (BNSs) is still under debate. The most likely candidate is an ultrastripped supernova (USSN) that has typically 10 times smaller ejecta mass (${ \mathcal O }(0.1){M}_{\odot }$) than a canonical supernova (Tauris et al. 2015). Because close binary interactions tear off the envelope of the secondary star in a binary system, the second supernova in the system would have little ejecta mass. The small ejecta mass also helps to prevent the binary system from breaking up. Recent progress of transient surveys has led to successful observations of USSNe from the early phase just after the explosion. In particular, De et al. (2018) observed iPTF 14gqr in the shock-cooling phase, indicating the presence of the extended envelope with 500 R_⊙ and 0.01 M_⊙. This result implies that USSNe are of massive star origin.

The energy source of a USSN is assumed to be radioactive ⁵⁶Ni as in the case of canonical stripped envelope supernovae (SNe). From the light-curve fitting, the amount of ⁵⁶Ni in a USSN ejecta is estimated to be 0.05 M_⊙ for SN 2005ek (Tauris et al. 2013), 0.05 M_⊙ for iPTF 14gqr (De et al. 2018), and 0.015 M_⊙ for SN 2019dge (Yao et al. 2020). These values are slightly smaller than those inferred for canonical SNe (see, e.g., Meza & Anderson 2020), which may be consistent with the explosion energy of USSNe estimated to be several times to an order of magnitude smaller than that of canonical SNe.

In order to theoretically calculate the ⁵⁶Ni amount loaded on the SN ejecta, we need detailed hydrodynamics simulations and nucleosynthesis calculations. A lot of these studies have been done for canonical SNe (e.g., Thielemann et al. 2018 and references therein). Although not as systematically investigated as canonical SNe, some previous works have performed calculations on nucleosynthesis in USSNe. Based on the first neutrino-radiation hydrodynamics simulation of USSNe by Suwa et al. (2015), Yoshida et al. (2017) conducted the detailed nucleosynthesis calculations for two cases. Müller et al. (2018) followed a similar simulation up to the shock breakout of a USSN progenitor star and investigated the explosive nucleosynthesis. The explosion energy of these simulations (∼10⁵⁰ erg) was consistent with the values obtained from the light-curve fitting, but the ⁵⁶Ni amount (∼0.01 M_⊙) was insufficient. Additionally, Moriya et al. (2017) have calculated the explosive nucleosynthesis, light curve, and spectrum of USSNe with a simplified explosion simulation. Still the synthesized amount of ⁵⁶Ni in the simulated USSNe is approximately ∼0.03 M_⊙, and not sufficient for SN 2005ek and iPTF 14gqr. We recall this ⁵⁶Ni problem here. It should be noted that all the above studies on nucleosynthesis in USSNe were investigated in a limited number of progenitor models. A systematic study is needed to make a detailed comparison between the theory and observation in USSNe.

The ⁵⁶Ni problem has recently been discussed as an inherent problem not only in USSNe but also in canonical SN explosion simulations (Sawada & Maeda 2019; Suwa et al. 2019; Sawada & Suwa 2021). It has long been known that the explosion energies obtained in simulations are significantly lower than the observed typical values. Recently, updated simulations have been reported to reach 10⁵¹ erg (Bollig et al. 2021; Burrows & Vartanyan 2021), and the explosion energy problem is gradually being solved. It should be noted, however, that almost all other simulations have not yet been able to reproduce it (see, e.g., Janka 2012; Takiwaki et al. 2016; Burrows & Vartanyan 2021; and references therein). Additionally, these simulations have yet to produce a sufficient amount of ⁵⁶Ni. The reason for this is that the explosive nucleosynthesis requires a rapid increase in the explosion energy to produce ⁵⁶Ni, whereas the simulation results show a slow increase.

In USSNe, such an explosion energy problem does not exist. Therefore, ⁵⁶Ni production in USSNe can be discussed more robustly than in normal supernovae. In this study, we perform a systematic study of the energy sources of USSN light curves with a wider range of progenitor models, connecting long-term hydrodynamics simulations with nucleosynthesis calculations to light-curve calculations. We first simulate one-dimensional hydrodynamics and nucleosynthesis, with the explosion model reconstructed from the results of neutrino-radiation hydrodynamic simulations (Suwa et al. 2015). We then compare the numerical results by calculating the light curves of the observed transients iPTF 14gqr and SN 2019dge using the analytical model.

This paper is organized as follows: We describe our model to solve the core-collapse explosion of ultrastripped progenitors and calculate the nucleosynthesis and the SN light curves in Section 2. We show the results of our calculation and compare the ⁵⁶Ni-powered light curves with the observed USSNe in Section 3. We discuss the uncertainties in our model and consider alternative energy sources of USSNe in Section 4. We conclude the paper in Section 5.

2.1. Progenitors

We use the seven preexplosion models of carbon–oxygen (CO) stars, the same as Suwa et al. (2015). These models are computed from an initial CO core of mass 1.45–2 M_⊙ (1.45, 1.5, 1.6, 1.7, 1.8, 1.9, and 2.0 M_⊙) at the central C burning using the stellar evolution code described in Suwa et al. (2015) and Yoshida et al. (2017). We will refer to CO core models with a mass of x.yz M_⊙ as COxyz models. These CO cores correspond to the secondary stars in the close binary system and are expected to lose their H- and He-rich envelopes during the binary evolution with the primary NS. See also Yoshida et al. (2017) and Suwa et al. (2018) for details.

Figure 1 shows the enclosed mass profiles of the USSN progenitor model. For comparison, we also show a canonical SN progenitor model by the dotted line with zero-age main-sequence mass of 15.0 M_⊙ obtained by the stellar evolution code MESA (Paxton et al. 2015). This figure shows that USSN progenitors have smaller core masses than the 15 M_⊙ model. For example, at ρ = 10⁶ g cm⁻³, USSN progenitor masses are 1.4–1.6 M_⊙, whereas the 15 M_⊙ model is 1.8 M_⊙. This result is attributed to the difference in core entropy, as shown in Suwa et al. (2018): the USSN progenitor has smaller core entropy than the 15 M_⊙ model due to the absence of entropy inflow from the high-entropy envelope. As a result, the increase in Chandrasekhar mass due to the finite temperature effect is ineffective, resulting in a smaller core mass when gravitational collapse occurs. Such a relatively small core mass has been reproduced in recent stellar evolution calculations including binary interaction effects (Jiang et al. 2021).

Zoom In Zoom Out Reset image size

2.2. Core-collapse Explosion

In order to estimate the ⁵⁶Ni amount loaded on the SN ejecta, one has to solve the long-term evolution including the fallback effect. Since the neutrino-radiation hydrodynamics simulations like Suwa et al. (2015) are not suitable for 1000 s long simulations because of their computational cost and numerical limitation of tabulated equation of state, we adopt a different approach in this work. For a given progenitor, we trigger an explosion by instantaneously injecting an amount of thermal energy (E_inject) at a certain mass–radius (M_inject). In our calculations, the initial proto-NS mass M_PNS,i is defined as the total mass of matter that has never been ejected outside a radius of R_c = 10¹⁰ cm, where the escape velocity becomes small enough compared to the ejecta velocity. We repeat the simulation with varying (E_inject, M_inject) until both the explosion energy E_expl and the initial proto-NS mass M_PNS,i become consistent with those of Suwa et al. (2015). See Table 1.

Table 1. Properties of the Explosion Models and Summary of Simulation Results

	Our Work			Suwa et al. (2015)		Our Result
Model	MCO	Eexpl	MPNS,i	Eexpl	MPNS	Eexpl	Mej	MNi	Mfb	MNS,f
	(M⊙)	(Bethe)	(M⊙)	(Bethe)	(M⊙)	(Bethe)	(M⊙)	(M⊙)	(10−2 M⊙)	(M⊙)
CO145	1.45	0.17	1.35	0.177	1.35	0.17	0.097	1.63 × 10−2	0.30	1.35
CO15	1.5	0.15	1.36	0.153	1.36	0.15	0.136	1.38 × 10−2	0.54	1.37
CO16	1.6	0.12	1.42	0.124	1.42	0.12	0.151	1.20 × 10−2	0.56	1.43
CO17	1.7	0.12	1.45	⋯	⋯	0.12	0.225	1.09 × 10−2	1.04	1.46
CO18	1.8	0.12	1.49	0.120	1.49	0.12	0.277	4.80 × 10−4	1.23	1.50
CO19	1.9	0.12	1.54	⋯	⋯	0.12	0.307	9.19 × 10−5	2.11	1.56
CO20	2.0	0.12	1.60	⋯	1.60	0.12	0.286	7.78 × 10−5	4.02	1.64

Download table as:  ASCIITypeset image

We note that different pairs of (E_inject, M_inject) can give the same (E_expl,M_PNS,i) even for a given progenitor model. In this case, we employ the case with a sufficiently small M_inject as our physical model, for which the properties of the explosion converges. We also note that no calculations were performed in Suwa et al. (2015) for some progenitor models, where we interpolated M_PNS,i from the results in a CO star model with close mass. In addition, we adopted E_expl = 0.12 Bethe (1 Bethe = 1 × 10⁵¹ erg) for all CO star models heavier than the CO16 model for simplicity, since the explosion energy tended to converge at ∼0.12 Bethe for those CO star models in Suwa et al. (2015). We will discuss impacts of these assumptions in Section 4.1.

For the spherically symmetric simulation, we use a 1D Euler code based on hydro1d, ⁷ which employs a Godunov-type scheme to integrate the conservation equations with a gravity source term:

$$\begin{eqnarray}&&\displaystyle \frac{\partial \rho }{\partial t}+\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial }{\partial r}\left({r}^{2}\rho {v}_{r}\right)=0\,,\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {v}_{r}}{\partial t}+{v}_{r}\displaystyle \frac{\partial {v}_{r}}{\partial r}=-\displaystyle \frac{1}{\rho }\displaystyle \frac{\partial P}{\partial r}-\displaystyle \frac{{GM}(r,t)}{{r}^{2}}\,,\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial (\rho e)}{\partial t}+\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial [(\rho e+P){v}_{r}{r}^{2}]}{\partial r}=-\rho {v}_{r}\displaystyle \frac{{GM}(r,t)}{{r}^{2}}\,,\end{eqnarray} \tag{ 3 }$$

where ρ, v_r , P, e, and M(r, t) are density, radial velocity, total pressure, specific energy, and enclosed mass, respectively. The enclosed mass is the sum of the mass inside the inner boundary M_G (t) and the mass in the computational domain interior to the radius r, i.e.,

$$\begin{eqnarray}&&M(r,t)={M}_{G}(t)+4\pi {\int }_{{R}_{\mathrm{in}}}^{r}\rho (r^{\prime} ,t)r{{\prime} }^{2}{dr}^{\prime} .\end{eqnarray} \tag{ 4 }$$

The mass flux flowing through the inner boundary of the computational domain is added to M_G (t). We employ the Helmholtz equation of state (Timmes & Swesty 2000) and neglect the effects of weak interaction and nuclear burning on the dynamics.

We employ a logarithmic grid with 300 mesh per decade in radius, i.e., the grid size ratio is fixed to be Δr/r ≈ 0.008. For t ≤ 10 s, we fix the inner boundary at a mass–radius 0.1 M_⊙ inner from the initial proto-NS mass–radius and the outer boundary at twice the progenitor radius. Here t = 0 s corresponds to the onset of the explosion. For t ≥ 10 s, we move both the inner and outer boundary radii outwards as the SN ejecta expands; every time the forward shock arrives at 85% of the outer boundary radii, we reset the inner boundary R_in → R_in + ΔR_in so that ΔM_G ≤ 10⁻³ M_⊙, and ΔR_in ≤ 5 R_in, and the new outer boundary is set to be at R_out = 10⁵ R_in keeping Δr/r ≈ 0.008. All of the primitive variables are remapped to the new mesh. We continue the calculation up to t ∼ 10 days. Note that we maintain the overall mass conservation within the machine precision.

2.3. Nucleosynthesis

2.4. Supernova Light Curve

We calculate the bolometric light curves of USSNe by the one-zone model (Arnett 1980, 1982; Kasen & Bildsten 2010; Dexter & Kasen 2013) using the ⁵⁶Ni mass M_Ni, the ejecta mass M_ej, and the explosion energy E_expl obtained by the core-collapse explosion simulations as the input parameters. The energy conservation of the SN ejecta is described as

$$\begin{eqnarray}&&\displaystyle \frac{\partial {E}_{\mathrm{int}}(t)}{\partial t}=-P\displaystyle \frac{\partial V(t)}{\partial t}+\dot{Q}(t)-L(t)\,,\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {E}_{\mathrm{kin}}(t)}{\partial t}=P\displaystyle \frac{\partial V(t)}{\partial t}\,,\end{eqnarray} \tag{ 6 }$$

where E_int(t) is the internal energy, E_kin(t) is the kinetic energy, $\dot{Q}(t)$ is the heating rate by some energy injection processes, and L(t) is the radiative cooling rate. The first term in the right-hand side of Equation (5) represents the adiabatic loss of the internal energy, which is converted into kinetic energy through Equation (6). The expansion velocity of the ejecta is given as

$$\begin{eqnarray}&&{v}_{\mathrm{ej}}(t)=\sqrt{2{E}_{\mathrm{kin}}(t)/{M}_{\mathrm{ej}}}\,,\end{eqnarray} \tag{ 7 }$$

and the radius evolves as

$$\begin{eqnarray}&&{R}_{\mathrm{ej}}(t)={\int }_{0}^{t}{v}_{\mathrm{ej}}(t^{\prime} ){dt}^{\prime} .\end{eqnarray} \tag{ 8 }$$

Assuming that the internal energy is dominated by radiation, the adiabatic loss term can be rewritten as

$$\begin{eqnarray}&&P\displaystyle \frac{\partial V(t)}{\partial t}=\displaystyle \frac{{E}_{\mathrm{int}}(t)}{{R}_{\mathrm{ej}}(t)}{v}_{\mathrm{ej}}(t)=\displaystyle \frac{{E}_{\mathrm{int}}(t)}{{t}_{\mathrm{dyn}}},\end{eqnarray} \tag{ 9 }$$

where

$$\begin{eqnarray}&&{t}_{\mathrm{dyn}}=\displaystyle \frac{{R}_{\mathrm{ej}}(t)}{{v}_{\mathrm{ej}}(t)}\,\end{eqnarray} \tag{ 10 }$$

is the dynamical timescale of the ejecta. On the other hand, the radiative loss term can be written as

$$\begin{eqnarray}&&L(t)=\displaystyle \frac{{E}_{\mathrm{int}}{t}_{\mathrm{dyn}}}{{t}_{\mathrm{diff}}^{2}},\end{eqnarray} \tag{ 11 }$$

where

$$\begin{eqnarray}&&{t}_{\mathrm{diff}}={\left[\displaystyle \frac{3}{4\pi }\displaystyle \frac{\kappa {M}_{\mathrm{ej}}}{{v}_{\mathrm{ej}}c}\displaystyle \frac{1}{\xi }\right]}^{1/2}\end{eqnarray} \tag{ 12 }$$

is the diffusion timescale t_diff through the ejecta, κ = 0.07 cm² g⁻¹ is the opacity (e.g., Pinto & Eastman 2000; Taddia et al. 2018), and ξ = π²/3 represents the geometrical factor (Arnett 1982). Substituting Equations (9) and (11) into Equations (5) and (6),

$$\begin{eqnarray}&&\displaystyle \frac{\partial {E}_{\mathrm{int}}(t)}{\partial t}=-\displaystyle \frac{{E}_{\mathrm{int}}(t)}{{t}_{\mathrm{dyn}}}+\dot{Q}(t)-\displaystyle \frac{{E}_{\mathrm{int}}{t}_{\mathrm{dyn}}}{{t}_{\mathrm{diff}}^{2}}\,,\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {E}_{\mathrm{kin}}(t)}{\partial t}=\displaystyle \frac{{E}_{\mathrm{int}}(t)}{{t}_{\mathrm{dyn}}}\,.\end{eqnarray} \tag{ 14 }$$

For a given heating rate $\dot{Q}(t)$, we solve Equations (13) and (14) with the fourth-order Runge–Kutta method to obtain E_int(t) and E_kin(t), and then the SN light curve L(t).

By default, we consider the radioactive decay of ⁵⁶Ni as the main heating source of the ejecta. In this case,

$$\begin{eqnarray}&&\dot{Q}(t)={f}_{\mathrm{dep}}\cdot \left({M}_{\mathrm{Ni}}\,{q}_{\mathrm{Ni}}(t)\right)\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}&&{q}_{\mathrm{Ni}}(t)={\epsilon }_{\mathrm{Ni}}\cdot {e}^{-t/{\tau }_{\mathrm{Ni}}}+{\epsilon }_{\mathrm{Co}}\cdot {e}^{-t/{\tau }_{\mathrm{Co}}}\,,\end{eqnarray} \tag{ 16 }$$

where _Ni = 3.22 × 10¹⁰ erg g⁻¹ s⁻¹ and _Co = 6.78 × 10⁹ erg g⁻¹ s⁻¹ are the specific decay energies of ⁵⁶Ni and ⁵⁶Co, and τ_Ni = 8.8 days and τ_Co = 113.6 days are the mean lifetimes of ⁵⁶Ni and ⁵⁶Co, respectively. We consider alternative energy sources in Section 4.2.

Figure 2 shows snapshots of the density profile in the CO145 model. For numerical reasons, we fill the outer region of the star with an ambient medium in hydrostatic equilibrium with ρ ∝ r⁻². The ambient density is chosen as $\rho {(r)={10}^{15}(r/1\,\mathrm{cm})}^{-2}$ g cm⁻³, and the total mass in the computational domain is less than 10⁻⁴ M_⊙ to avoid artificially slowing down the ejecta. While the structure of the outside of the star in the low-density ambient medium might influence the forward shock velocity, it does not significantly affect the dynamics of the inner fallback region (see, e.g., Fernández et al. 2018). The inner region has a structure as ρ ∝ r^−3/2 that is derived by the approximately constant mass accretion rate, $\dot{M}\propto \rho {r}^{2}{v}_{\mathrm{ff}}$ with the freefall velocity v_ff ∝ r^−1/2, when the shock reaches 10¹² cm, corresponding to about t ≈ 1000 s after the explosion. Here, we should note that, in the previous core-collapse explosion calculations for progenitors with massive outer layers, the inward reverse shock self-reflecting at the inner boundary was observed (Ertl et al. 2016). The timing of the self-reflection is sensitive to the position of the inner boundary and artificially alters the fallback dynamics (Vigna-Gómez et al. 2021). In cases of USSN progenitors without massive outer layers, there is no reverse shock causing the above self-reflection problem.

Zoom In Zoom Out Reset image size

Figure 3 shows the time evolution of the SN fallback. We estimate the mass accretion rate as

$$\begin{eqnarray}&&{\dot{M}}_{\mathrm{fb}}(t)=4\pi {R}_{{\rm{c}}}^{2}{\rho }_{{\rm{c}}}(t){v}_{r,{\rm{c}}}(t)\,,\end{eqnarray} \tag{ 17 }$$

where ρ_c and v_r,c are the density and radial velocity at radius R_c = 10¹⁰ cm, respectively. We find a larger fallback mass for a larger CO core mass case; the lightest CO145 model has a total fallback mass of M_fb ≲ 0.01 M_⊙ while the heaviest CO20 model has M_fb ∼ 0.04 M_⊙. These values are comparable to those estimated for canonical SN explosions (e.g., Zhang et al. 2008; Janka et al. 2021); the explosion energy of USSNe smaller than canonical SNe is compensated by the less dense core of ultrastripped progenitors to give a comparable fallback mass. In all models, the fallback rate declines as ∝ t^−5/3 at t ≳ 1000 s (Michel 1988; Chevalier 1989). We note that, in the case of USSNe, there is no enhancement of the fallback rate in a later phase that could occur for a progenitor with a helium layer.

Zoom In Zoom Out Reset image size

Figure 4 shows the time evolution of the temperature profile of the CO145 model. The temperature at the shock front can be well presented by an analytical expression (Woosley et al. 2002) as

$$\begin{eqnarray}&&{T}_{\mathrm{shock}}(r)=7.48\times {10}^{9}\,{\rm{K}}\,{\left(\displaystyle \frac{E(r)}{{10}^{50}\mathrm{erg}}\right)}^{1/4}{\left(\displaystyle \frac{r}{{10}^{8}\mathrm{cm}}\right)}^{-3/4}\,,\end{eqnarray} \tag{ 18 }$$

where E(r) = E_expl + ∣E_bind(r)∣ is an energy in the postshock region, and ∣E_bind(r)∣ is the binding energy of the progenitor outside the radius r. In regions where the postshock temperature exceeds 5 × 10⁹ K, the nuclear statistical equilibrium is realized within a dynamical timescale, and the material turns into iron group elements, mainly ⁵⁶Ni (e.g., Thielemann et al. 1996; Nomoto et al. 2013). The critical radius for the nucleosynthesis can be estimated from Equation (18) as

$$\begin{eqnarray}&&{R}_{{T}_{9}=5}=2.17\times {10}^{8}\,\mathrm{cm}\,{\left(\displaystyle \frac{E(r)}{{10}^{50}\mathrm{erg}}\right)}^{1/3}\,,\end{eqnarray} \tag{ 19 }$$

where T₉ = T/10⁹ K.

Zoom In Zoom Out Reset image size

Figures 5 and 6 show results of our explosive nucleosynthesis calculation for the CO145 and CO20 models, respectively. In the lighter case, a few × 10⁻² M_⊙ of ⁵⁶Ni are synthesized during the explosion. A fraction of them fall back and ∼ 0.01 M_⊙ of ⁵⁶Ni are ejected. In the heavier case, although the amounts of ⁵⁶Ni synthesized are comparable, most of them fall back and only ≲10⁻⁴ M_⊙ of ⁵⁶Ni are ejected. Table 1 summarizes the amount of ⁵⁶Ni ejected in all cases.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Our results show that both the synthesized and ejected amounts of ⁵⁶Ni in USSNe are significantly smaller than those estimated for canonical CCSNe. This trend is mainly due to less compactness of USSN progenitors and comparable proto-NS masses to canonical SN explosions. In more detail, this can be understood as follows. The synthesized amount of ⁵⁶Ni can be roughly calculated as the enclosed mass within the ⁵⁶Ni synthesizing radius ${R}_{{T}_{9}=5}$ (see Equation (19)) minus the NS mass. A smaller explosion energy of USSNe gives a smaller ${R}_{{T}_{9}=5}$, and ultrastripped progenitors have a smaller mass–radius for a fixed radius than the canonical SN progenitor (see Figure 1), thus USSNe have a smaller enclosed mass within the ⁵⁶Ni synthesizing radius. On the other hand, the calculated NS mass is roughly comparable between USSNe and canonical SNe. Resultantly, the synthesized amount of ⁵⁶Ni becomes smaller for USSNe. The ejected amounts of ⁵⁶Ni also become smaller for USSNe since the fallback masses and the proto-NS masses are comparable between USSNe and canonical SNe.

Figure 7 shows bolometric light curves of USSNe. The ejecta mass M_ej and explosion energy E_expl, among the three input parameters obtained by hydrodynamical simulations, determine the shape of the light curve, especially the peak time (see Equation (11)). When the ejecta mass M_ej is smaller, the peak of the light curve becomes faster. The peak also becomes faster when the explosion energy E_expl is larger. Figure 7 suggests that the differences in ejecta mass and explosion energy within the range of the present hydrodynamical results make little difference in the shape of the light curve, and the other parameter, nickel mass M_Ni, changes the peak luminosity. Since a USSN from a more massive progenitor has a larger ejecta mass but a smaller ⁵⁶Ni mass, the light curve becomes dimmer and slower; the peak luminosity can be as large as 10⁴² erg s⁻¹ for the relatively small CO core mass cases where the ejected amount of ⁵⁶Ni is M_Ni ∼ 0.01 M_⊙ while the peak luminosity is significantly smaller for the three heaviest cases (CO18, CO19, CO20) for which the ejected amount of ⁵⁶Ni is significantly smaller mainly due to the fallback.

Zoom In Zoom Out Reset image size

Figure 7 also shows the observed bolometric light curves of iPTF 14gqr and SN 2019dge (De et al. 2018; Yao et al. 2020). We find that SN 2019dge is broadly consistent with the relatively small CO core mass cases. ⁹ However, none of our models can reproduce the peak luminosity of iPTF 14gqr. The dotted line indicates the best-fit model obtained by De et al. (2018) with M_Ni = 0.05 M_⊙, M_ej = 0.2 M_⊙, and E_expl = 0.2 Bethe. Such a large amount of ⁵⁶Ni cannot be ejected in our theoretical calculations. This inconsistency can be due to that: (i) our model underestimates the ejected amount of ⁵⁶Ni in USSNe, or (ii) there is an additional energy source other than the ⁵⁶Ni decay at least for relatively bright USSNe. We consider both possibilities in the following section.

4.1.  ⁵⁶Ni Problem?

Our long-term explosion simulation is based on the results of neutrino-hydrodynamics simulations in that the explosion energy E_expl and the initial proto-NS mass M_PNS,i are consistent. However, there are several model uncertainties that could potentially lead to underestimating the ejected amount of ⁵⁶Ni. Here we investigate this point.

First, we discuss the degeneracy in the explosion model: for a given progenitor model, the same (E_expl, M_PNS,i) can be obtained for a series of pairs of injection energy and mass–radius (E_inject, M_inject). Among them, we have employed the case with a sufficiently small M_inject that would better mimic the result of neutrino-hydrodynamics simulation. This choice is somewhat arbitrary and gives a conservative estimate on the synthesized ⁵⁶Ni mass. Figure 8 shows the evolution of the shock downstream temperature of the CO145 model with various M_inject. For larger M_inject cases, the region experiencing temperatures above 5 × 10⁹ K is more extended and thus more ⁵⁶Ni is produced. Nevertheless, Figure 8 indicates that the increase in the synthesized amount of ⁵⁶Ni compared with our fiducial case (M_inject = 1.25 M_⊙) is at most ∼ 0.02 M_⊙ for the CO145 model. Since the proto-NS masses from first-principles calculations is 1.35 M_⊙, the synthesizing radius ${R}_{{T}_{9}=5}$ (see Equation (19)) would have to extend to around 1.4 M_⊙ to synthesize ∼0.05 M_⊙ of ⁵⁶Ni, which is very difficult to achieve. We confirm this is also the case for heavier CO core models. A detailed discussion of the maximum amount of ⁵⁶Ni that can be synthesized, up to the nucleosynthesis calculation, is shown in the Appendix.

Zoom In Zoom Out Reset image size

Next, we discuss the uncertainties of the explosion energy in massive progenitor models. E_expl = 0.12 Bethe adopted for models heavier than CO16 is the only assumption, since the final explosion energy is not obtained by the neutrino-radiation hydrodynamics simulations due to their weak explosions and the limited simulation period. Here we investigate how the explosion energy would affect the hydrodynamic behavior and nucleosynthesis. Figure 9 shows the evolution of the peak temperature behind the shock wave for different explosion energies. We confirm that while the nickel synthesis region is extended by about 0.06 M_⊙ due to larger energy, all of them are still entirely contained within the fallback region and are not expected to escape from the central proto-NS.

Zoom In Zoom Out Reset image size

Another possible indeterminacy includes neutrino-driven wind, which is expected to produce additional ⁵⁶Ni after explosive synthesis. Sawada & Suwa (2021) claimed that the ⁵⁶Ni amount by neutrino-driven wind is strongly related to the mass accretion rate onto the proto-NS just after the onset of the explosion. For the case of a USSN that has a small envelope mass and negligible mass accretion during the important phase, the neutrino-driven wind is not expected to fill the gap described above.

In summary, even taking into account the model's uncertainties, our calculations indicate that it is difficult to synthesize ∼ 0.05 M_⊙ of ⁵⁶Ni in a USSN, which may conflict with the observed light curves of SN 2005ek and iPTF 14gqr. First-principles simulations of the entire USSN explosion and more involved calculations of the ultrastripped progenitors including various effects of binary evolution will be needed to further investigate this issue.

4.2. Alternative Energy Sources

Here we consider alternative energy sources of USSNe other than ⁵⁶Ni decay. One possibility is the fallback accretion onto the newborn NS, e.g., a fraction of the gravitational energy of the accreted matter can be released via outflow from an accretion disk formed by the fallback materials (e.g., McKinney et al. 2012; Dexter & Kasen 2013). Another possibility is the pulsar wind, i.e., the rotational energy of the newborn NS is extracted by the dipole radiation and forms a nascent pulsar wind nebula (e.g., Hotokezaka et al. 2017). We calculate USSN light curves powered by these additional energy sources and compare them in particular with iPTF 14gqr.

4.2.1. Fallback Accretion

To include the additional energy injection to the SN ejecta by the fallback accretion, we modify the heating term in Equation (16) as follows (e.g., Dexter & Kasen 2013; Moriya et al. 2018):

$$\begin{eqnarray}&&\dot{Q}(t)={f}_{\mathrm{dep}}\left({M}_{\mathrm{Ni}}\,{q}_{\mathrm{Ni}}(t)\right)+{\dot{Q}}_{\mathrm{acc}}(t),\end{eqnarray} \tag{ 20 }$$

$$\begin{eqnarray}&&{\dot{Q}}_{\mathrm{acc}}(t)=\eta {\dot{M}}_{\mathrm{fb}}(t){c}^{2}\,,\end{eqnarray} \tag{ 21 }$$

where ${\dot{M}}_{\mathrm{fb}}(t)$ is the fallback rate obtained by our simulation (see Equation (17) and Figure 3), and η is the energy conversion efficiency that is largely uncertain. Here we consider a range of η ≤ 10⁻³ following Dexter & Kasen (2013).

Figure 10 shows the bolometric light curves of the CO145, CO 16, and CO20 models including the energy injection by the fallback accretion. We find that the fallback accretion can make the USSN light curve brighter and faster. In particular, the fallback accretion can significantly increase the peak luminosity of the CO20 model where the ejected amount of ⁵⁶Ni is suppressed due to the intense fallback itself. Even in this case, however, the light-curve shape is incompatible with iPTF 14gqr. Since the fallback accretion rate decreases faster than the ⁵⁶Ni-decay rate, it mainly contributes to the very early phase (t ≲ 5 days). This effect is specific to the USSN, which shows a fast-rising light curve due to its extremely low ejecta mass. Thus, we also confirmed that there was no significant effect when adapting similar energy injection to canonical SNe in the range of typical fallback accretion. We conclude that the fallback accretion could be an important energy source for USSNe especially in the very early phase, but could not be the main energy source for relatively bright USSNe like iPTF 14gqr.

Zoom In Zoom Out Reset image size

4.2.2. Rotation-powered Relativistic Wind

In order to incorporate the energy injection by the rotation-powered relativistic wind to the light-curve calculation, we modify the heating term in Equations (5) and (6) as

$$\begin{eqnarray}&&\dot{Q}(t)={f}_{\mathrm{dep}}\left({M}_{\mathrm{Ni}}\,{q}_{\mathrm{Ni}}(t)\right)+{\dot{Q}}_{{\rm{p}}}(t)\end{eqnarray} \tag{ 22 }$$

$$\begin{eqnarray}&&{\dot{Q}}_{{\rm{p}}}(t)=\displaystyle \frac{{E}_{p}}{{t}_{p}}\displaystyle \frac{1}{{\left(1+t/{t}_{p}\right)}^{2}}\,,\end{eqnarray} \tag{ 23 }$$

where ${E}_{p}=2\times {10}^{50}\,{P}_{10}^{-2}$ erg is the rotational energy, ${t}_{p}=0.44\,{B}_{14}^{-2}{P}_{10}^{2}$ yr is the spin-down timescale, B₁₄ = (B/10¹⁴ G) is the dipole magnetic field strength, and P₁₀ = (P_i /10 ms) is the spin period of the newborn NS (Kasen & Bildsten (2010). The NS radius is set to be R_ns = 12 km. In Equation (23), we assume the energy conversion efficiency from the rotation-powered wind to the heat as ${ \mathcal O }(1)$, which is appropriate for the bolometric light curve at around the peak (Hotokezaka et al. 2017), and also neglect the spin-down of the NS via gravitational-wave radiation (e.g., Kashiyama et al. 2016). Such wind-driven USSN light curves can have diverse shapes depending on B and P_i (Hotokezaka et al. 2017). Here we focus on searching for the parameter set that can reproduce the light curve of iPTF 14gqr.

Figure 11 shows the bolometric light curves of the CO16 model with the additional energy injection by the rotation-powered wind. We investigate a parameter region of 1.0 × 10⁴⁸ erg < E_p < 4.0 × 10⁵⁰ erg and 1 day < t_p < 1000 days, and show the cases with an fixed rotation energy (E_p = 4.0 × 10⁴⁸ erg) for representative purpose. We find that the bolometric light curve of iPTF 14gqr can be well fitted with E_p = 4.0 × 10⁴⁸ erg and t_p = 19.8 days, namely B ≈ 2.0 × 10¹⁵ G and P_i ≈ 70 ms. A similar set of fitting parameter (B,P_i ) is obtained when other CO core mass models are used; the light-curve shape is still characterized by the diffusion timescale of the USSN ejecta, which is roughly comparable in the original models (see Section 3), while the peak luminosity is mainly determined by ∼ E_p /t_p or (B,P_i ).

Zoom In Zoom Out Reset image size

We note that the above NS parameters for iPTF 14gqr are compatible with the Galactic magnetars (e.g., Enoto et al. 2019), but different from those inferred for superluminous supernovae (B ∼ 10^13–14 G and P_i ∼ 1–10 ms; e.g., Kashiyama et al. 2016). Our result might imply that such a magnetar formation is relatively common in core-collapse explosions in a close binary system.

We have conducted the explosion simulations of ultrastripped progenitors with various CO core masses based on the results of neutrino-radiation hydrodynamics simulations, and consistently calculated the nucleosynthesis and the SN light curves. We found the amount of ⁵⁶Ni synthesized and ejected in USSNe are smaller than those in canonical CCSNe mainly due to less compactness of USSN progenitors and comparable proto-NS masses to canonical CCSNe. Nonetheless, the relatively small CO core mass models can eject sufficient amounts of ⁵⁶Ni to reproduce the light curves of some observed USSNe (e.g., SN 2019dge; Yao et al. 2020). On the other hand, it is difficult for any progenitor model to synthesize ∼ 0.05 M_⊙ of ⁵⁶Ni as inferred for relatively bright USSNe (e.g., SN 2005ek and iPTF 14gqr; Tauris et al. 2013; De et al. 2018). We have investigated the fallback accretion onto and relativistic wind from the newborn NS as alternative energy sources for such USSNe. We found that iPTF 14gqr can be explained by considering an energy injection from a newborn NS with a magnetic field of B_p ∼ 10¹⁵ G and an initial rotation period of P_i ∼ 0.1 s.

We thank T. Yoshida for providing progenitor models. The work has been supported by Japan Society for the Promotion of Science (JSPS) KAKENHI grants 21J00825, 21K13964 (R.S.), 20K04010 (K.K.), 18H05437, 20H00174, 20H01904, and 20H04747 (Y.S.).

Software: MESA (Paxton et al. 2011, 2013, 2015, 2018, 2019), torch (Timmes 1999), hydro1d (https://zingale.github.io/hydro1d/).

Here, as a supplement to the discussion of the degeneracy of the explosion model in Section 4.1, we show the maximum amount of ⁵⁶Ni that can be synthesized by the explosion model with a thermal bomb. We find that the highest amount of synthesized ⁵⁶Ni is given by the model with the outermost energy injection position. It allows us to discuss the most robust maximum limit on the amount of ⁵⁶Ni that can be ejected by USSNe, within the range of reproducing the same (E_expl, M_PNS,i) as neutrino-radiation hydrodynamics simulations. Note, however, that this setup does not allow us to estimate the natural fallback mass accretion rate, so we do not treat this choice as the primary model in this study.

Figure A1 shows the results of our explosive nucleosynthesis calculation with the outermost energy injection position for the CO145 and CO20 models, respectively. CO145 is the model with the highest amount of ⁵⁶Ni synthesis, while CO20 is the model with the lowest amount. The results of other models with similar setups are summarized in Table A1. Moriya et al. (2018) adopted the same thermal bomb injection setup as the current one, so our results are consistent with their results, which are M_Ni ∼ 0.026 M_⊙ at similar explosion energies (E_expl = 0.10 Bethe) and CO star masses (M_CO = 1.50 M_⊙). In order to reproduce the explosion energy of ∼0.2 Bethe and synthesized ⁵⁶Ni mass of ∼0.05 M_⊙ estimated for iPTF 14gqr, the explosive nucleosynthesis would have to proceed up to ∼1.40 M_⊙, for instance, in the CO145 model. Figure A1 shows that it is very difficult to achieve. We emphasize again that this result is based on a model that leads to the maximum amount of ⁵⁶Ni synthesis, and then we can give a more robust conclusion that it is difficult to synthesize ∼ 0.05 M_⊙ of ⁵⁶Ni in a USSN.

Zoom In Zoom Out Reset image size