What Powered the Kilonova-Like Emission After GRB 230307A in the Framework of a Neutron Star-White Dwarf Merger?

Abstract

The second brightest gamma-ray burst, GRB 230307A (with a duration $T_{90}\sim$ 40 s), exhibited characteristics indicative of a magnetar engine during the prompt emission phase. Notably, a suspected kilonova was identified in its follow-up optical and infrared observations. Here we propose that the origin of GRB 230307A is a neutron star-white dwarf (NS-WD) merger, as this could naturally interpret the long duration and the large physical offset from the center of its host galaxy. In the framework of such a NS-WD merger event, the late-time kilonova-like emission is very likely to be powered by the spin-down of the magnetar and the radioactive decay of ${}^{56}$Ni, rather than by the decay of $r$-process elements as these heavy elements may not be easy to be synthesized in a NS-WD merger. It is demonstrated that the above scenario can be supported by our fit to the late-time observational data, where a mass of ${\sim}10^{-3}\ \rm M_{\odot}$ ${}^{56}$Ni is involved in the ejecta of a mass of ${\sim}0.1\ \rm M_{\odot}$. Particularly, the magnetar parameters required by the fit are consistent with those derived from the early X-ray observation.

1 Introduction

In the past decades, it has been widely accepted that long GRBs are produced by the collapse of massive stars and short GRBs originate from compact binary mergers like mergers of double neutron stars (DNSs) and of a neutron star and a black hole (NS-BH; Eichler et al., 1989; Narayan et al., 1992; Meszaros & Rees, 1992). The smoking-gun evidence to support such a common view is (i) the discovery of Type Ib/c supernovae emerging from the afterglow of some long GRBs (Galama et al., 1998; Bloom et al., 1999; Hjorth et al., 2003; Stanek et al., 2003; Woosley & Bloom, 2006) and (ii) the detection of short GRB 170817A (Goldstein et al., 2017; Savchenko et al., 2017) in counterpart with the gravitational wave event GW170817 (Abbott et al., 2017a). In addition, accompanying the GW radiation, it is expected that DNSs and NS-BH mergers can simultaneously eject an amount of $\sim 10^{-4}-10^{-2}M_{\odot}$ neutron-rich material (Hotokezaka et al., 2013; Radice et al., 2018), which can lead to effective nucleosynthesis of heavy elements through rapid neutron-capture processes ($r$-processes). Then, the radioactive decay of these $r$-process elements can heat the merger ejecta and lead to an ultraviolet-optical-infrared transient emission, namely “kilonova” (Li & Paczyński, 1998; Freiburghaus et al., 1999; Metzger et al., 2010; Metzger, 2017). Such a kilonova model could explain the optical transient AT 2017gfo following the GRB 170817A/GW170817 event (Abbott et al., 2017b; Smartt et al., 2017; Arcavi et al., 2017; Coulter et al., 2017; Pian et al., 2017; Goldstein et al., 2017), in despite of some difficulties (Li et al., 2018).

Recently, suspected kilonova candidates have been surprisingly discovered from some GRBs of duration much longer than 2 seconds, which is a traditional boundary between long and short GRBs, e.g., GRB 060614 with $T_{90}\sim 100\ \rm s$, GRB 211211A and GRB 230307A with $T_{90}\sim 40\ \rm s$ (Della Valle et al., 2006; Yang et al., 2015; Rastinejad et al., 2022; Troja et al., 2022; Yang et al., 2022; Zhu et al., 2022; Wang et al., 2023; Levan et al., 2023; Sun et al., 2023; Dichiara et al., 2023; Gillanders et al., 2023; Yang et al., 2023; Meng et al., 2023). These facts challenge the traditional classification of GRBs, indicating that some long GRBs can be produced by compact binary mergers. Such a merger origin of GRB 230307A is further supported by its large offset of the location from its host galaxy center ($\sim$ 36.60 kpc; Levan et al., 2023; Sun et al., 2023; Dai et al., 2023), as a compact binary can often travel far from its host galaxy before the merger due to the high kick velocity of the NS component in the binary (Kalogera et al., 1998). In principle, a long duration is not completely impossible to be produced by a merger event. For example, if a DNS merger has a high mass ratio $q\gtrsim 1.2$ and high total mass, or a NS-BH merger has a low mass ratio $q\lesssim 3$ and a high pre-merger BH spin, then the merger product (i.e., a BH) could be surrounded by a magnetically arrested massive ($\sim 0.1M_{\odot}$) disk, which could extend the energy release of the BH through an accretion tail of $\dot{M}\propto t^{-2}$ (Gottlieb et al., 2023a, b).

However, for GRB 230307A, the broad-band (0.5–6000 keV) observations during its prompt emission revealed that the central engine of this GRB is probably a rapidly rotating and highly-magnetized remnant NS (Sun et al., 2023). First, its soft X-ray emission in the 0.5–4.0 keV energy range monitored by Lobster Eye Imager for Astronomy (LEIA) shows a different spectrum from the energy band above 15 keV monitored by the Gravitational wave high-energy Electromagnetic Counterpart All-sky Monitor (GECAM). Second, the late decay slope of soft X-ray emission does not conform to the curvature effect. Finally, the soft X-ray light curve is basically consistent with the shape of an initial plateau followed by a steep decay, which has been usually regarded as an observational signature of a spinning-down millisecond magnetar. In summary, the NS-BH merger origin of GRB 230307A can be ruled out. On other hand, although such a remnant mangetar can in principle be produced by a DNS merger, this scenario is still challenged by the long duration of the prompt emission. Furthermore, the existence of the remnant magnetar can also make the DNS merger scenario difficult to explain the $>$ 3 $\mu$m spectroscopic features of GRB 230307A (Gillanders et al., 2023), since the long-term neutrino emission from the magnetar can effectively suppress the synthesis of lanthanide elements (Martin et al., 2015).

In comparison, because the density of white dwarfs (WDs) is much lower than that of NSs, the free-fall timescale of a NS-WD merger can be much longer than that of a DNSs or a NS-BH merger, which may provide a more natural explanation for the origin of GRB 230307A. In the NS-WD merger scenario, a stable remnant magnetar is very probably produced because of the relatively low mass of the remnant (Bartos et al., 2013; Zhang, 2018). However, the problem in this case is how the kilonova-like transient emission can be generated after the merger, since the ejecta from the NS-WD merger could not be an ideal place for synthesizing $r$-process elements (Metzger, 2012; Zenati et al., 2019; Fernández et al., 2019; Kaltenborn et al., 2023). Instead, a small amount of ${}^{56}$Ni can be expected to form (Bobrick et al., 2022; Zhong et al., 2023), the radioactive decay of which can contribute to the thermal emission of ejecta. Besides, it needs to be noticed that the spin-down of the remnant magnetar could power the ejecta emission significantly, as previously suggested for other magnetar-driven transient phenomena (Kasen & Bildsten, 2010; Inserra et al., 2013; Yu et al., 2013; Metzger & Piro, 2014; Yu et al., 2015a; Ai et al., 2022).

Therefore, the purpose of this Letter is to test whether the parameters of magnetar constrained from the early X-ray emission of GRB 230307A can further explain the multi-wavelength afterglow emission plus kilonova-like emission in the framework of a NS-WD merger. In Section 2, we introduce the model and present the fitting result. Conclusion and discussion is drawn in Section 3.

2 Joint fit of Kilonova-like Component and Afterglow

2.1 The Model

As discussed above, the optical excess in the afterglow of GRB 230307A is considered to be powered by the combination of the spin-down of a newborn magnetar and the radioactive decay of ${}^{56}$Ni. Therefore, we first invoke the spin-down luminosity of the magnetar as usual as

$$L_{\rm sd}(t)=L_{\rm sd}(0)\left(1+{t\over t_{\rm sd}}\right)^{-2},$$ (1)

where the spin-down is assumed to be dominated by the magnetic dipole radiation, $L_{\rm sd}(0)$ is the initial spin-down luminosity, and $t_{\rm sd}$ is the spin-down timescale. In the following calculations, we fixed $t_{\rm sd}$ = 80 s according to the fit of the soft X-ray light curve during the prompt emission phase, where a smoothly broken power law is used with a break time of $\sim$ 80 s (Sun et al., 2023). Meanwhile, the radioactive power per unit mass of ${}^{56}$Ni and its daughter nucleus ${}^{57}$Co is given by

$$\dot{q}_{{\rm r}}=\epsilon_{\rm Ni}e^{-t/\tau_{\rm Ni}}+\epsilon_{\rm Co}\frac{e^{-t/\tau_{\rm Co}}-e^{-t/\tau_{\rm Ni}}}{1-\tau_{\rm Ni}/\tau_{\rm Co}},$$ (2)

where $\epsilon_{\rm Ni}$ = $3.9\times 10^{10}$ $\rm erg\ s^{-1}\ g^{-1}$ and $\epsilon_{\rm Co}$ = $6.8\times 10^{9}$ $\rm erg\ s^{-1}\ g^{-1}$ are energy generation rate per unit mass (Sutherland & Wheeler, 1984; Maeda et al., 2003), and $\tau_{\rm Ni}$ = 8.8 days and $\tau_{\rm Co}$ = 111.3 days are the decay time.

A simplified radiation transfer model given by Kasen & Bildsten (2010) and Metzger (2017) is adopted to describe the emission of the merger ejecta, where the ejecta is separated into $N$ layers with different velocities and densities. The density of each layer is determined by the following density profile of the ejecta (Nagakura et al., 2014):

$$\rho_{\rm ej}(R,t)=\frac{(\delta-3)M_{\rm ej}}{4\pi R_{\max}^{3}}\left[\left(\frac{R_{\min}}{R_{\max}}\right)^{3-\delta}-1\right]^{-1}\left({R\over R_{\max}}\right)^{-\delta},$$ (3)

where $M_{\rm ej}$ is the total mass of the ejecta, and $R_{\min}$ = $v_{\min}t$ ($R_{\max}$ = $v_{\max}t$) is the innermost (outermost) layer of the ejecta by introducing the minimum (maximum) velocity. $R=v_{\rm i}t$ is the radii of each layer assuming there is no dynamic evolution of each layer, where $v_{\rm i}$ is a constant velocity. The evolution of the thermal energy of the $i$th layer is determined by (Yu et al., 2018)

$${dE_{i}\over dt}=\xi L_{\rm sd}+m_{{\rm{Ni},i}}\dot{q}_{{\rm r},i}-{E_{i}\over R_{i}}{dR_{i}\over dt}-L_{i},i=1,$$ (4a)

$${dE_{i}\over dt}=m_{{\rm Ni},i}\dot{q}_{{\rm r},i}-{E_{i}\over R_{i}}{dR_{i}\over dt}-L_{i},i>1.$$ (4b)

Here, the spin-down energy of the magnetar is assumed to be primarily absorbed by the ejecta at its bottom with an efficiency of $\xi$. ${}^{56}$Ni is considered to distribute the same as the ejecta does, where $m_{{\rm Ni},i}$ represents the mass of ${}^{56}$Ni of each layer. $L_{\rm i}$ is the observed luminosity of each layer, which can be calculated by

$$L_{i}={E_{i}\over\max(t_{{\rm d},i},t_{{\rm lc},i})},$$ (5)

where

$$t_{{\rm d},i}={\frac{3\kappa}{4\pi R_{i}c}}{\sum\limits_{j=i}^{n}}m_{j}$$ (6)

is the radiation diffusion timescale for the photons from the $i$th layer to escape from the ejecta (Arnett, 1982), $\kappa$ is the opacity, and $t_{{\rm lc},i}=R_{i}/c$ is the light crossing time. Then, the total luminosity of the merger ejecta is obtained by summing $L_{\rm i}$ up:

$$L_{\rm bol}={\sum\limits_{i=1}^{n}}L_{i}.$$ (7)

Finally, the corresponding temperature of this thermal emission can be given by using the Stefan-Boltzmann law, and the specific luminosity for different colors can be derived from the Planck function by assuming a black-body spectrum for this emission.

Generally speaking, the light curves given by our semi-analytical model can be consistent with the simulation results of Zenati et al. (2020) and Kaltenborn et al. (2023) in some respects, although the specific frequency-dependence of the light curves could still be different. The differences can arise from the different distributions of the ejecta velocity and ${}^{56}$Ni mass and, in particular, the frequency-dependence of the opacity adopted in their simulations. Because of these complex factors, the spectra of the thermal emission of ejecta could deviate from the black body significantly. Nevertheless, due to the overlapping with the bright afterglow emission from the GRB jet, the thermal emission of ejecta can actually emerge only around its peak with a very limited data number. This makes it nearly impossible to constrain the above mentioned detailed properties of the ejecta, whereas the key features of the explosion system can still be inferred from the emission peak. Therefore, it is feasible and convenient to use the semi-analytic model to estimate the masses of the ejecta and ${}^{56}$Ni and the order of magnitude of the opacity for the primary emission band. Such a method has also been widely employed in literature for the fittings of supernova light curves (e.g., Shigeyama et al., 1987; Pinto & Eastman, 2000; Blinnikov et al., 2006; Nicholl et al., 2016, 2017; Liu et al., 2021).

2.2 Fitting Results

Table 1: The fitting parameters, prior bounds, and posterior medians at 1$\sigma$ credible interval.

Parameter	Prior Bounds	Median $\&$ 1$\sigma$ C.I.
Afterglow
log $E_{0}$ (erg)	[50.0, 60.0]	$51.11_{-0.49}^{+0.93}$
log $n_{0}$ (cm${}^{-3}$)	[-6.0, 2.0]	$0.42_{-0.56}^{+1.10}$
log $\theta_{\rm c}$ (rad)	[-2.0, -0.5]	$-0.70_{-0.06}^{+0.05}$
log $\epsilon_{\rm e}$	[-4.0, -0.2]	$-0.92_{-0.76}^{+0.51}$
log $\epsilon_{\rm B}$	[-6.0, -0.2]	$-2.50_{-1.12}^{+0.39}$
$p$	[2.01, 3.0]	$2.03_{-0.02}^{+0.02}$
log $\xi_{\rm N}$	[-5.0, 0.0]	$-1.29_{-0.44}^{+0.74}$
Kilonova-like Component
$v_{\rm min}$ (c)	[0.001, 0.15]	$0.12_{-0.01}^{+0.01}$
$v_{\rm max}$ (c)	[0.18, 0.35]	$0.24_{-0.05}^{+0.06}$
$\delta$	[1.0, 3.0]	$2.58_{-0.65}^{+0.33}$
log $\kappa$ ($\rm cm^{2}$ g${}^{-1}$)	[-1.0, 0.0]	$-0.11_{-0.15}^{+0.08}$
log $M_{\rm ej}$ ($\rm M_{\odot}$)	[-2.0, -0.8]	$-0.93_{-0.19}^{+0.09}$
log $M_{\rm Ni}$ ($\rm M_{\odot}$)	[-5.0, -3.0]	$-3.32_{-0.04}^{+0.03}$
log $\xi L_{\rm sd}(0)$ (erg s${}^{-1}$)	[45.0, 48.5]	$48.09_{-0.07}^{+0.09}$

Figure 1: The corner plot of the posterior probability distributions of 14 parameters listed in Table 1 for the fitting of the multi-wavelength light curves.

There are a total of 14 free parameters in our model. Besides 6 parameters for the merger ejecta and 1 parameter for the magnetar, we also need to invoke another 7 parameters to model the multi-wavelength afterglow light curves of the GRB, which is ascribed to the synchrotron radiation due to the shock interaction of a structured jet and the external medium (Sari et al., 1998; Huang et al., 1999). Here, a Gaussian structure $E(\theta)=E_{0}\rm{exp}(-\theta^{2}/2\theta_{\rm c}^{2})$ is applied for the jet. The truncation angle is set to $\theta_{\rm w}$ = 4$\theta_{\rm c}$, and the viewing angle $\theta_{\rm v}$ is set to 0 rad in our fitting as same as those in Yang et al. (2023). The public python package afterglowpy (Ryan et al., 2020) is used to reproduce the theoretical afterglow emission. The Markov Chain Monte Carlo (MCMC) method with the Python package emcee (Foreman-Mackey et al., 2013) is applied to fit the multi-wavelength light curves by minimizing $\chi^{2}$. We set 5000 walkers in a 56-dimensional Gaussian, and discard the first 2000 steps. The prior bounds, medians, and $1\sigma$ credible intervals for each parameter are given in Table 1 and the corner plot of the posterior probability distributions is shown in Figure 1.

Figure 2: The multi-wavelength light curves. The solid line is the afterglow plus kilonova-like emission, while the dashed line is the afterglow-only emission. Both lines are plotted with the median parameters in our fitting. The data are collected from Yang et al. (2023).

Figure 3: Spectral energy distributions. The black points are the unabsorbed spectra at different epochs from Yang et al. (2023). The solid and dashed lines are the afterglow plus kilonova-like emission and afterglow-only emission, respectively.

Figure 4: The bolometric light curves. The dashed and dotted lines represent the absorbed spin-down energy of the ejecta and the radioactive heating power from ${}^{56}$Ni, respectively. The solid line is the total energy of the two energy sources.

With the median values of the model parameters, we plot the multi-wavelength light curves and the temporal evolving spectral energy distributions (SEDs) of the optical emission of GRB 230307A in Figures 2 and 3, respectively. The corresponding bolometric light curve is presented in Figure 4. The bolometric light curve reaches its peak is in about 5 d in our fitting, similar to that in Kaltenborn et al. (2023) (5-8 d), but different from that in Zenati et al. (2020) (10-20 d). The multi-wavelength light curves, spanning from the U-band to I-band, becomes redder over time in Zenati et al. (2020) and Kaltenborn et al. (2023), which can also be found from Figure 2. Furthermore, the obtained mass of ${}^{56}$Ni can be in the reasonable range predicted by the simulations. Figure 4 shows that the radioactive power can dominate the radiation at the late stage, while the magnetar power dominates the first 20 days.

To be specific, the initial absorbed spin-down energy from the magnetar is constrained to be $\xi L_{\rm sd}(0)$ = $1.23\times 10^{48}$ erg s${}^{-1}$ and the isotropic kinetic energy of the afterglow jet is $E_{0}=1.29\times 10^{51}$ erg. The total injected spin-down energy and the kinetic energy of the ejecta can be calculated as $E_{\rm total}\approx\xi L_{\rm sd}(0)$$\times t_{\rm sd}$ = $9.84\times 10^{49}$ erg and $E_{\rm k}=\frac{1}{2}{\sum\limits_{i=1}^{n}}m_{\rm i}v_{\rm i}^{2}$ = $3.56\times 10^{51}$ erg, respectively. Since $E_{\rm k}$ $\gg$ $E_{\rm total}$ and $E_{\rm total}$ $\ll$ $E_{0}$, it is safe to neglect the dynamic evolution of each layer and the energy injection of the magnetar to the afterglow jet. According to Sun et al. (2023), the initial X-ray luminosity is $L_{\rm X}(0)\sim 2\times 10^{48}$ erg s${}^{-1}$ during the prompt emission phase, and is of the same magnitude as $\xi L_{\rm sd}(0)$ in our fitting. Our result suggests that the presence of a magnetar in GRB 230307A is self-consistent in the prompt emission phase and the following multi-wavelength observations.

3 Conclusion And Discussion

GRB 230307A is the first GRB with a magnetar engine signature in the prompt emission phase, the second brightest GRB ever observed (only dwarfed by GRB 221009A), and the third long GRB in association with suspected kilonova emission. The fitting results presented in this Letter demonstrate that the NS-WD merger with a magnetar remnant could be a self-consistent explanation for GRB 230307A both in the prompt emission and the following multi-wavelength afterglow including the optical excesses. The identification of NS-WD merger events from optical transient phenomena usually faces a challenge as their early electromagnetic radiation cannot be recorded initially, making it difficult to identify the nature of the remnant objects. By contrast, GRB 230307A offers a unique opportunity as it was observed by multiple satellites across various wavelengths and, in particular, the early prompt emission of the GRB provides the magnetar parameters, which enables us to test the consistency of the GRB 230307A observations with the NS-WD merger scenario.

By comparing with some simulation studies in literature (e.g., Zenati et al., 2020; Kaltenborn et al., 2023), the NS-WD merger origin of GRB 230307A could further be supported by its detailed emission features. For example, Kaltenborn et al. (2023) showed a prominent emission line at 20000 Å at 30.37 d, which can be found in the observation of GRB 230307A at 29 d. In the framework of the NS-WD merger model, the notable emission features of GRB 230307A at $\sim$ 2.1 $\mu$m (at 29 and 61 days) and $\sim$ 4.4 $\mu$m (at 61 days) could be attributed to the Ca ii line as well as the first overtone and second overtone features of the CO molecular line (Levan et al., 2023; Gillanders et al., 2023), although these features were previously ascribed to the $r$-process elements in Gillanders et al. (2020). Anyway, it should be noticed that our knowledge of the spectra of NS-WD mergers is actually very limited. For example, the early optical spectroscopy from 5000 Å to 9200 Å recorded by the Gemini South telescope demonstrated the absence of Si ii absorption lines around 6300 Å at 2.4 d, which is however expected at early stage in NS-WD mergers (Zenati et al., 2020). Such a situation also appeared in another NS-WD merger candidate event, AT 2018kzr, which showed no absorption lines above 5000 Å at 1.9 d and a weak Si ii absorption line at 3.8 d (Gillanders et al., 2020). In order to clarify these detailed emission features, a more elaborated model definitely needs to be considered in future works, by taking into account the complex nuclear heating and cooling, and the multigroup opacity in different energy band even in different region of the ejecta. Furthermore, the potential significant influence of the remnant magnetar on all of these physical properties should also be investigated carefully.

Nevertheless, for GRB 230307A, actually only the peak of the ejecta emission can be resolved from the observational data, the characteristic of which is generally determined by the basic properties of the ejecta and energy sources (cf., the Arnett’s law for supernova light curves (Arnett, 1982)). Therefore, the simplified model can still be workable in our fittings. In future, it is expected that upcoming instruments like the Jonit Science Center for China Space Station Telescope (CSST) would be capable of measuring the spectra more precisely and, especially, the establishment of a collaborative network involving multiple telescopes would hold the potential to reveal more detailed features of the light curves and spectra of such kilonova-like transients. In particular, if the GRB afterglow emission can decay more quickly, then the ejecta emission would be detected at time much later than its peak, where the emission spectra could significantly deviate from the black body and many line features could appear. In that cases, the comparison between the observational data and a more elaborated model would lead to a more reliable constraint on the details of the nucleosynthesis and opacity of the ejecta, as investigated in theory in Kaltenborn et al. (2023) and Zenati et al. (2019, 2020).

Finally, we would like to point out that the formation of a system including a NS and a ${}^{56}$Ni-rich ejecta could also be caused by the accretion-induced collapse (AIC) of a WD in companion with a non-degenerate donor or another WD (Nomoto & Kondo, 1991; Taam & van den Heuvel, 1986; Duncan & Thompson, 1992; Usov, 1992; Yoon et al., 2007), which therefore could also lead to a bright optical transient (Yu et al., 2015b, 2019b, 2019a). Nevertheless, the AIC processes are deemed less likely to produce GRBs; otherwise, the observed rate of GRBs would be too high to be consistent with the observations. This may hint that the magnetic fields of the AIC-formed NSs are relatively lower than those of the NS-WD merger products, according to the statistical inferences given by Yu et al. (2017). Besides the association with the GRB emission, the NS-WD mergers can also be different from the AIC events by their gravitational wave (GW) radiation, which would be detected by future millihertz and decihertz space-borne GW detectors (Yin et al., 2023; Kang et al., 2024). In any case, it is still necessary to investigate the theoretical differences of the optical transients between the AIC events and NS-WD mergers.