Multimessenger Implications of AT2018cow: High-Energy Cosmic Ray and Neutrino Emissions from Magnetar-Powered Super-Luminous Transients

TL;DR

This paper investigates high-energy cosmic ray and neutrino emissions from magnetar-powered super-luminous transients like AT2018cow, assessing their detectability and implications for multi-messenger astronomy.

Contribution

It provides the first detailed modeling of high-energy particle emissions from FBOTs, especially AT2018cow, and evaluates their detectability with current observatories.

Findings

Neutrino fluence from AT2018cow is below IceCube sensitivity.

Cosmically-integrated neutrino flux from FBOTs aligns with IceCube limits.

High-energy gamma rays may be observable at later times.

Abstract

Newly-born, rapidly-spinning magnetars have been invoked as the power sources of super-luminous transients, including the class of "fast-luminous optical transients" (FBOTs). The extensive multi-wavelength analysis of AT2018cow, the first FBOT discovered in real time, is consistent with the magnetar scenario and offers an unprecedented opportunity to comprehend the nature of these sources and assess their broader implications. Using AT2018cow as a prototype, we investigate high-energy neutrino and cosmic ray production from FBOTs and the more general class of superluminous supernovae (SLSNe). By calculating the interaction of cosmic rays and the time-evolving radiation field and baryon background, we find that particles accelerated in the magnetar wind may escape the ejecta at ultrahigh energies (UHE). The predicted high-energy neutrino fluence from AT2018cow is below the sensitivity of the IceCube Observatory, and estimates of the cosmically-integrated neutrino flux from FBOTs are consistent with the extreme-high-energy upper limits posed by IceCube. High-energy $\gamma$ rays exceeding GeV energies are obscured for the first months to years by thermal photons in the magnetar nebula, but are potentially observable at later times. Given also their potentially higher volumetric rate compared to other engine-powered transients (e.g. SLSNe and gamma-ray bursts), we conclude that FBOTs are favorable targets for current and next-generation multi-messenger observatories.

Tables

Table 1. Table 1: Magnetar parameters of fiducial models

	$P_i$ [ms]	$B_d$ [G]
Case I	10	${10}^{15}$
Case II	2	$2\times {10}^{15}$

Equations

$$\frac{d{E}_{\mathrm{nth}}}{dt}$$

$$\frac{d{E}_{\mathrm{th}}}{dt}$$

$$\frac{d{E}_B}{dt}$$

$$M_{\mathrm{ej}}{v}_{\mathrm{ej}}\frac{d{v}_{\mathrm{ej}}}{dt}$$

$$\frac{d{R}_{\mathrm{ej}}}{dt}$$

$$L_{\mathrm{sd}}$$

$$t_{\mathrm{sd}}\equiv \frac{E_{\mathrm{rot}}}{L_{\mathrm{sd},0}}=1.0P_{i,-2}^2{B}_{d,15}^{-2}d.$$

$$v_{\mathrm{ej}}$$

$$t_d^n\approx \frac{R_{\mathrm{n}}}{c}\left(1+{\tau }_{\mathrm{es}}^n\right),$$

$$t_{d,0}^n\approx 0.03B_{15}^2{P}_{i,-2}^{-4}{\beta }_{\mathrm{ej},-1.3}d.$$

$$t_d^{\mathrm{ej}}\approx \frac{R_{\mathrm{ej}}}{c}\left(1+{\tau }_{\mathrm{es}}^{\mathrm{ej}}\right),$$

$$t_{d,0}^{\mathrm{ej}}=2.1M_{\mathrm{ej},-0.5}^{1/2}{\left(\frac{\kappa }{0.01\mathrm{c}{\mathrm{m}}^2{\mathrm{g}}^{-1}}\right)}^{1/2}{\beta }_{-1.3}^{-1/2}d.$$

$$L_{\mathrm{th}}=\frac{E_{\mathrm{th}}}{t_d^{\mathrm{ej}}}$$

$$B_{\mathrm{n}}\approx {\left(\frac{6{ϵ}_B{L}_{\mathrm{sd}}t}{R_{\mathrm{n}}^3}\right)}^{1/2}\simeq 720{ϵ}_{B,-4}^{1/2}{B}_{15}{P}_{i,-2}^{-2}{\beta }_{n,-2}^{-3/2}{t}_5^{-2}\mathrm{G},$$

$$n_{\mathrm{th}}$$

$$n_{\mathrm{nth}}$$

$$n_p\approx \frac{3M_{\mathrm{ej}}}{4\pi R_{\mathrm{ej}}^3{m}_p}\approx 2.5\times 10^{13}{M}_{\mathrm{ej},-0.5}{t}_5^{-3}{\beta }_{\mathrm{ej},-1.3}^{-3}\mathrm{c}{\mathrm{m}}^{-3}.$$

$$E_{\mathrm{CR}}\simeq \eta Ze{\Phi }_{\mathrm{mag}}$$

$$\frac{dN}{dE}=\frac{9}{8}\frac{c^{2}I}{e\mu }\frac{1}{E}=2\times 10^{42}{E}^{-1}{B}_{15}^{-1}.$$

$$t_{p\gamma ,\mathrm{int}}^{-1}=\frac{c}{2{\gamma }_p^2}{\int }_0^{\infty }d\epsilon \frac{n(\epsilon )}{{\epsilon }^2}{\int }_0^{2{\gamma }_p\epsilon }d{\epsilon }^{\prime }{\epsilon }^{\prime }{\sigma }_{p\gamma }({\epsilon }^{\prime })$$

$$t_{pp}^{-1}=n_p{\sigma }_{pp}{\kappa }_{pp}c$$

$$\frac{r_L}{R_{\mathrm{n}}}\approx 1.7{\eta }_{-1}{ϵ}_{B,-4}^{-1/2}{\beta }_{n,-2}^{1/2},t\gg t_{\mathrm{sd}}.$$

$$t_{p,0}^{\mathrm{ej}}\approx 16M_{\mathrm{ej},-0.5}^{1/2}{\beta }_{\mathrm{ej},-1.3}^{-1},d.$$

$$t_{xp}^{-1}=n_p{\sigma }_{xp}{\kappa }_{xp}c.$$

$$t_{\pi p,0}\approx 2.2{\eta }_{-1}^{1/4}{B}_{15}^{-1/4}{\beta }_{\mathrm{ej},-1.3}^{-3/4}{M}_{\mathrm{ej},-0.5}^{1/4}d.$$

$$E_{\nu ,b}\approx 59.7M_{\mathrm{ej},-0.5}^{-1}{t}_5^3{\beta }_{\mathrm{ej},-1.3}^{3}PeV.$$

$${\left(E^2{F}_{\nu }\right)}_{\mathrm{pk}}\approx \frac{1}{2}{E}_p(t_{\pi p,0})\frac{d{N}_p}{dt}(t_{\pi p,0})\frac{1}{4\pi D^2}$$

$${\left(E^2{J}_{\nu }\right)}_{\mathrm{pk}}\sim {\left(E^2{F}_{\nu }\right)}_{\mathrm{pk}}{t}_{\pi p,0}$$

$${\Phi }_{\nu }=\frac{c}{4\pi }{\int }_0^{z_{\max }}dz\rho (z){\frac{dN}{d{E}^{\prime }}|}_{E^{\prime }=(1+z)E}(1+z)\left|\frac{dt}{dz}\right|$$

$${\tau }_{\gamma \gamma }$$

$$t_{\gamma ,0}\simeq 450P_{i,-2}^{-1}{\beta }_{\mathrm{ej},-1.3}^{-1/2}\mathrm{day}$$

Full text

1 Introduction

Time-domain optical surveys have in recent years discovered a growing number of luminous, rare and rapidly-evolving extragalactic transients, with visual light curves that rise and decay on a timescale of only a few days (hereafter “fast blue optical transients”, or FBOTs; e.g. Drout et al. 2013, 2014; Arcavi et al. 2016; Rest et al. 2018). Although supernova light curves are often powered by the radioactive decay of 56Ni, such an energy source is ruled out for FBOTs (even if the entire ejecta mass were composed of nickel, the resulting optical luminosity would be insufficient to explain the observations). Instead, FBOTs must be powered by an additional, centrally-concentrated energy source or sources, either due to shock interaction between the ejecta and a massive circumstellar shell (e.g. Kleiser et al. 2018), or from a central compact object, such as a rapidly spinning magnetar or accreting black hole.

Indeed, newly-formed magnetars with millisecond rotation periods have widely been invoked as engines responsible for powering the broader class of “superluminous supernovae” (SLSNe; e.g. Kasen & Bildsten 2010; Woosley 2010). Within this framework, the faster evolution of FBOTs requires stellar explosions with lower ejecta masses than of typical SLSNe. Such low ejecta masses \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.1-1M_{\odot} could plausibly arise in a number of scenarios: explosions of ultra-stripped progenitor stars (e.g. Tauris et al. 2015); electron capture supernovae (e.g. Moriya & Eldridge 2016); initially “failed” explosions of very massive stars such as blue supergiants (in which most of the stellar mass collapses into the central black hole, with only a tiny ejected fraction; e.g. Fernández et al. 2018); the accretion-induced collapse of a rotating white dwarf (e.g. Metzger et al. 2008; Lyutikov & Toonen 2018); magnetars created as the long-lived stable remnants of neutron star binary mergers (e.g. Yu et al. 2013; Metzger & Piro 2014; Murase et al. 2018). On the other hand, compared to normal SLSNe or neutron star mergers, FBOTs may be relatively common, occurring at rates of up to $\sim 4-7$ % that of the core collapse supernova rate (Drout et al. 2014). FBOTs also seem to occur exclusively in star-forming galaxies, supporting an origin associated with the deaths of massive stars (Drout et al., 2014). Rapidly-spinning pulsars in ultra-stripped supernovae have also been considered to explain some of the FBOTs (Hotokezaka et al., 2017).

Until recently, all FBOTs were discovered at high redshift ($z>0.1$) and were not identified in real-time, thereby precluding key spectroscopic or multi-wavelength follow-up observations. This situation changed dramatically with the recent discovery of the optical transient AT2018cow (Smartt et al., 2018; Prentice et al., 2018), which shared many characteristics with FBOTs but took place at a distance of only $\approx$ 60 Mpc ($z=0.0141$). AT2018cow reached a peak visual luminosity of L_{\rm opt}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}10^{44} erg s*-1* within only a few days, before fading approximately as a power-law $L_{\rm opt}\propto t^{-\alpha}$ with $\alpha\approx 2.5$ (e.g. Perley et al. 2019). The rapid rise time was indicative of an explosion with a low ejecta mass, M_{\rm ej}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.5M_{\odot}, and a high ejecta velocity, v_{\rm ej}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}0.1 c. The latter also revealed itself through the self-absorbed radio synchrotron emission from the source (de Ugarte Postigo et al. 2018; Ho et al. 2019; Margutti et al. 2018), which was likely produced by the same fast ejecta undergoing shock interaction with a dense surrounding medium. The early spectra were nearly featureless (also consistent with a high photospheric expansion speed), but after a few weeks narrow helium and hydrogen emission lines appeared (Perley et al., 2019), providing clues to the nature of the progenitor star (e.g. Margutti et al. 2018). The persistent blue colors of the transient (Perley et al., 2019), again unlike the behavior of normal supernovae, provides additional evidence that a central ionizing energy source is responsible for indirectly powering the optical light through absorption and “reprocessing” (Margutti et al., 2018).

Perhaps the most remarkable thing about AT2018cow is the direct visibility of the engine itself. AT2018cow was accompanied by bright X-ray emission (Rivera Sandoval et al. 2018; Margutti et al. 2018; Kuin et al. 2019) with a power-law spectrum and luminosity comparable to that of the optical emission, but showing rapid variability on timescales at least as short as a few days. The fact that X-rays appear to be visible from the same engine responsible for powering the optical luminosity provides evidence that the ejecta is globally asymmetric, such that the observed X-rays preferentially escape along low-density polar channels (see Margutti et al. 2018 for additional evidence and discussion). Prentice et al. (2018) proposed that the engine behind AT2018cow is a millisecond magnetar with a dipole field strength of $B_{\rm d}\sim 10^{15}$ G and an initial spin period of $P_{i}\sim 10$ ms. However, other magnetar parameters are consistent with the data and alternative central engine models remain viable, such as the collapse of a blue super-giant star to an accreting black hole, or a deeply embedded shock arising from circumstellar interaction (Margutti et al., 2018).

If AT2018cow is powered by a central compact object, particularly a millisecond magnetar, then it could also be a source of high energy charged particles (Blasi et al., 2000; Arons, 2003; Fang et al., 2012), high-energy neutrinos (Murase et al., 2009, 2014; Fang et al., 2014; Fang, 2015), high-energy gamma-rays (Kotera et al., 2013; Murase et al., 2015; Renault-Tinacci et al., 2018), and gravitational waves (Stella et al., 2005; Dall’Osso et al., 2009; Kashiyama et al., 2016). Such a scenario has also been discussed in the context of pulsar-powered supernova, and the similar low-ejecta mass context of binary neutron star mergers (Piro & Kollmeier 2016; Fang & Metzger 2017). Indeed, two IceCube neutrino track events were reported in spatial coincidence with AT2018cow during a 3.5 day period following the optical discovery. The events are consistent with the rate of atmospheric neutrino background (Blaufuss 2018). While we will show that it is unlikely these neutrinos are physically associated with the transient, they nevertheless motivate a more thorough study of FBOTs as sources of ultrahigh energy cosmic rays (UHECRs) and neutrinos, particularly given their high volumetric rate as compared to other previously considered UHECR sites, such as gamma-ray bursts and SLSNe.

In this paper, we follow the procedure outlined in Fang & Metzger (2017) to calculate the neutrino emission from millisecond magnetars embedded in low ejecta mass explosions, focusing on parameters motivated by the watershed event, AT2018cow. In §2 we describe our model for the radiation and hadronic background. In §3 we describe the acceleration and escape of UHECRs. In §4 we discuss neutrino production and address the neutrino detection prospects, both of AT2018cow and the cosmic background of FBOTs. Throughout this paper we adopt the short-hand notation $q_{x}\equiv q/10^{x}$ cgs.

2 Radiation and Hadron Background

The energetic compact object is embedded within a rapidly expanding ejecta shell filled with baryons and photons. Cosmic rays, which are accelerated close to the compact object, will lose energy and produce neutrinos when traveling through this dense medium. We follow Metzger & Piro (2014) to calculate the density of ejecta baryons, thermal and non-thermal photons. This simplified model assumes the ejecta to be homogeneous and spherically symmetric. While this is reasonable for estimating the environment for cosmic-ray interaction, modeling of AT2018cow indicates that the true ejecta structure is aspherical and thus is more complicated in detail (e.g. Margutti et al. 2018; Liu et al. 2018).

The evolution of non-thermal radiation $E_{\rm nth}$, thermal radiation $E_{\rm th}$, magnetic energy $E_{B}$, ejecta radius $R_{\rm ej}$ and velocity $v_{\rm ej}$ are described by a coupled set of differential equations:

[TABLE]

The system is powered by the dipole spin-down of the magnetar (Ostriker & Gunn, 1969), which injects a luminosity given by

[TABLE]

where $L_{\rm sd,0}=4{\mu^{2}\Omega^{4}}/{9\,c^{3}}=2.6\times 10^{45}\,B_{d,15}^{2}\,P_{i,-2}^{-4}\,\rm erg\,s^{-1}$ is the initial spin-down power111We use the vacuum dipole formula, following the precedent of Özel et al. (2010). In reality the spin-down is instead best described as a force-free MHD wind (e.g. Spitkovsky 2006), but the results are similar to within a factor of a few (e.g. Kashiyama et al. 2016)., $\mu=B_{d}R_{*}^{3}$ is the magnetic moment of a pulsar with dipole field $B_{d}$ and stellar radius $R_{*}$, and $P_{i}=10^{-2}\,P_{i,-2}$ s is the initial spin period. The spin-down time, $t_{\rm sd}$, which characterizes the timescale over which most of the initial rotational energy, $E_{\rm rot}=I\Omega_{i}^{2}/2$, is released is given by

[TABLE]

The pulsar moment of inertia is given by $I=2\,M_{*}R_{*}^{2}\,/5$, where $R_{*}=10$ km and $M_{*}=1.4\,M_{\odot}$ are its radius and mass, respectively. At times $t\gg t_{\rm sd}$, $L_{\rm sd}\propto t^{-2}$, corresponding to a pulsar braking index of $3$ (Ostriker & Gunn, 1969). This decay rate is moderately shallower than that of the bolometric luminosity of AT2018cow ($L\propto t^{-2.5}$; Perley et al. 2019; Margutti et al. 2018); however, a smaller braking index, or additional mechanisms such as the interaction of fall-back from the explosion with the magnetosphere (Metzger et al., 2018), can result in a steeper decay. We also consider magnetic energy losses due to the expansion, which can be present, e.g., if the magnetic field in the nebula is turbulent.

Finally, internal magnetic fields at the scale of $10^{16}$ G can be generated in the core of newborn magnetars if their initial spin period is on the order of a few milliseconds (Stella et al., 2005; Shore et al., 2009). The anisotropic pressure from the toroidal B-field leads to an ellipticity of $\varepsilon_{B,\rm GW}\sim 6.4\times 10^{-4}\,B_{t,16.3}^{2}$ for an average toroidal field strength $B_{t}$ (Stella et al., 2005). Such a quadrupole moment will cause a spin down due to gravitational wave emission at the rate $(d\Omega/dt)_{\rm GW}=(32/5)(G/c^{5})I\varepsilon^{2}_{B,\rm GW}\,\Omega^{5}$. In comparison, the spin-down rate due to the dipole field is $(d\Omega/dt)_{\rm d}=\mu^{2}/(6Ic^{3})\Omega^{3}$. The ratio of the two rates is $(d\Omega/dt)_{\rm GW}/(d\Omega/dt)_{\rm d}=5\%\,P_{i,-3}^{-2}\,B_{d,15}^{-2}\,B_{t,16.3}^{2}$. This fraction could be higher depending on the equation of state of the stellar interior, but should not dominate the energy losses. We thus ignore the quadrupole term in this work, and note that it could impact the spin-down at the level of a few to ten percents. We note that an associated gravitational wave flux from such a source should not be detectable with current instruments (Arons, 2003; Kotera, 2011) except that it is extremely close (Kashiyama et al., 2016), but a stochastic gravitational wave signature from a population of magnetars-powered transients could lead to detection with future instruments (Kotera, 2011).

A significant portion of the magnetar’s rotational energy is ultimately used to accelerate the ejecta, as captured by equation 4. The mean velocity of the ejecta after time $t$ is thus approximately given by

[TABLE]

where v_{\rm 0}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.05\,c is the initial ejecta velocity from the supernova explosion (Prentice et al., 2018) and $M_{\rm ej}=0.3\,M_{\rm ej,-0.5}\,M_{\odot}$ is the ejecta mass. The nebular radius $R_{\rm n}$ is taken to be a fixed fraction of the mean ejecta radius $R_{\rm n}=R_{\rm ej}/5$, where $R_{\rm ej}\approx v_{\rm ej}t$, consistent with observations indicating that the photosphere velocity of AT2018cow declines from \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}0.05\,c at $4$ days to $\sim 0.01-0.02\,c$ over its first two weeks of evolution (Prentice et al., 2018).

At the wind termination shock, the majority of the spin-down power is converted into non-thermal radiation (e.g., Tanaka & Takahara, 2010, 2013). A portion of this is then converted to thermal radiation via absorption by the ejecta walls, as described by equations 1 and 2, respectively. The timescales $t_{d}^{n}$ and $t_{d}^{\rm ej}$ are those required for photons to diffuse radially through the nebula or ejecta shell, respectively. Specifically,

[TABLE]

where $\tau_{\rm es}^{n}=\left(4\,Y\,L_{\rm sd}\,\sigma_{T}/\pi\,m_{e}\,c^{3}\,R_{\rm n}\right)^{1/2}\approx 17\,B_{15}\beta_{\rm ej,-1.3}^{-1/2}\\ t_{5}^{-1/2}P_{i,-2}^{-2}$ is the optical depth due to Thomson scattering by electron positron pairs. The latter are generated by a pair cascade process in the nebula due to its high compactness (Metzger & Piro, 2014; Fang & Metzger, 2017), where $Y\approx 0.1$ is the pair multiplicity (fraction of the spin-down power converted into pair rest mass). Once the photon diffusion time becomes less than the nebula expansion time $t_{\rm d}^{n}=R_{\rm n}/\beta_{\rm n}c$, photons are able to travel freely through the nebula; this occurs after a characteristic time

[TABLE]

Similarly, photons diffuse through the ejecta on a timescale

[TABLE]

where $\tau_{\rm es}^{\rm ej}={3M_{\rm ej}\kappa}/{4\pi R_{\rm ej}^{2}}$ is the optical depth and $\kappa\sim 0.01-0.4\,\rm cm^{2}\,g^{-1}$ is the ejecta opacity at UV/visual frequencies. The free escape of photons takes place once $t_{d}^{\rm ej}<R_{\rm ej}/\beta_{\rm ej}c$, as occurs after the time

[TABLE]

Equations 1 and 2 were derived under the approximation that all of the non-thermal radiation is reprocessed into thermal radiation following its absorption by the ejecta (Metzger et al., 2014). This is justified by the observational fact that the X-ray luminosity of AT2018cow (that escaping from the engine) is less than the optical luminosity (that thermalized) over the first several weeks of evolution (e.g. Margutti et al. 2018).

The thermal luminosity (supernova light curve) is given by the loss term in equation (2),

[TABLE]

At early times ($t\ll t_{d}^{n}$), the kinetic term ($PdV$ work) dominates the energy loss, such that $E_{\rm nth}\sim L_{\rm sd,0}t$ and $E_{\rm th}\sim L_{\rm sd,0}t/\tau_{\rm es}^{n}\propto t^{3/2}$. The luminosity thus scales $L_{\rm th}\sim E_{\rm th}/t_{\rm d}^{\rm ej}\propto t^{5/2}$ at $t\ll t_{d,0}^{\rm ej}$. At late times ($t\gg t_{d,0}^{\rm n}$, $t_{d,0}^{\rm ej}$), $E_{\rm nth}\sim L_{\rm sd,0}t_{\rm sd}\left(t/t_{\rm sd}\right)^{-1}$ and $E_{\rm th}\sim E_{\rm nth}t_{d}^{\rm ej}/t_{d}^{n}\propto t^{-1}$. As a result, $L_{\rm th}$ declines $\propto t^{-2}$, i.e. following the magnetar spin-down luminosity.

Figure 1 compares the time-dependent thermal luminosity calculated from equation 13 to the observed optical light curve of AT2018cow (as integrated from the UV to near IR wavelength bands). We consider two models for the properties of the central magnetar, referred to hereafter as Case I and Case II, respectively (see Table 1). Case I corresponds to a magnetar with dipole field $B_{d}=10^{15}$ G and initial spin period $P_{i}=10$ ms. This model reproduces well the observed optical light curve, as shown in Figure 1, and is consistent with previous fits to the magnetar model (Prentice et al., 2018; Margutti et al., 2018). We use parameters from Case I in our analytic estimates below. Case II corresponds to a magnetar with a somewhat stronger field $B_{d}=2\times 10^{15}$ G and larger rotational energy $P_{i}=2$ ms. The light curve in this case, as shown by a dashed line in Figure 1, is also roughly consistent with that of AT2018cow (though decaying a bit too quickly at late times).

The magnetic field in the magnetar-wind nebula can be estimated as

[TABLE]

where a fraction $\epsilon_{\rm B}$ of the spin-down energy $\sim L_{\rm sd}t$ is assumed to be placed into the magnetic energy of the nebula, $E_{B}=(B_{n}^{2}/8\,\pi)V_{\rm n}$ (e.g. at the wind termination shock), where $V_{\rm n}=4\pi R_{\rm n}^{3}/3$ is the nebula volume. Typical values $\epsilon_{B}\sim 10^{-4}-10^{-2}$ are obtained by modeling pulsar wind nebulae (e.g, Kennel & Coroniti 1984; Torres et al. 2014).

The temperature of the ejecta, and thus that of the thermal radiation field, is approximately given by $T_{\rm th}=\left({3E_{\rm th}}/{4\pi aR_{\rm ej}^{3}}\right)^{1/4}$. The number density of thermal photons is then given by

[TABLE]

The thermal emission of the hot magnetar as considered in Kotera et al. (2015) is subdominant here as the acceleration and interaction sites are distant from the star. As detailed in Section 3, the acceleration should happen beyond the light cylinder. The density of thermal photons from the magnetar is a factor of $(R_{\rm star}/R_{\rm acc})^{2}$ times lower than that near the star, where $R_{\rm star}$ is the stellar radius and $R_{\rm acc}$ the radius of the acceleration site. The contribution of these photons to neutrino production is thus negligible.

The number density of non-thermal photons can likewise be estimated to be

[TABLE]

where we have assumed a flat power-law spectrum $n_{\rm nth}(\varepsilon)\propto\varepsilon^{-2}$, extending from the energy of the thermal radiation $\varepsilon_{\rm min}\sim 3\,k_{B}\,T_{\rm th}$ to the pair creation threshold $\varepsilon_{\rm max}\sim 2\,m_{e}c^{2}$ (Svensson, 1987). The observed X-ray spectrum of AT2018cow ($n(\varepsilon)\sim\varepsilon^{-1.5}$; Margutti et al. 2018) is somewhat harder than assumed in the model. However, this difference does not critically affect our conclusions because the density of the non-thermal radiation is lower than that of the thermal photons (as well as higher in energy), and therefore are generally less important targets for neutrino production.

The baryon density of the ejecta is given by

[TABLE]

While $n_{p}$ scales with time as $\propto t^{-3}$, the evolution of $n_{\rm th}$ changes at $t_{d,0}^{\rm n}$ and $t_{d,0}^{\rm ej}$, introducing features to the cosmic-ray interaction, as discussed below.

3 Acceleration and Escape of UHECRs

3.1 Cosmic Ray Injection

Pulsars and magnetars offer promising sites for particle acceleration. In general, charged particles may tap the open field voltage and gain an energy (Arons, 2003)

[TABLE]

where $\eta=0.1\,\eta_{-1}$ is the acceleration efficiency and $Z$ is the particle charge. For simplicity, we assume cosmic rays of proton composition. The effects of a heavier composition on neutrino production are discussed in Section 5.

Assuming that ions follow the Goldreich-Julian charge density, the rate of cosmic ray injection from the magnetic polar cap is given by $\dot{N}=\mu\Omega^{2}/Zec$ (Arons, 2003). Equation 18 then implies the energy spectrum of accelerated particles is given by

[TABLE]

The specific acceleration mechanism of high-energy particles in the pulsar magnetosphere is still debated (e.g. Cerutti & Beloborodov 2017). For UHECRs, Arons (2003) hypothesized that particles achieve their energy in the relativistic wind through surf-riding acceleration at a radius $10^{3}-10^{4}$ times the light cylinder. The possibility of wake-field acceleration has also been discussed (e.g., Murase et al., 2009; Iwamoto et al., 2017). Recent particle-in-cell simulations (Philippov & Spitkovsky, 2018) support a significant fraction of the accelerated particle energy flux being carried by ions extracted from the stellar surface. The ion energy may reach up to $10-30\%$ of the polar cap voltage. As the misalignment of the rotational and magnetic axises impact the current distribution, the spectrum also depends on the inclination angle of the pulsar. Alternatively, if the acceleration happens close to the star, curvature radiation will be non-negligible and limit the maximum energy. However, particles could be picked up and re-accelerated by the pulsar wind. If the pair multiplicity is low, a significant wind power would go into ions and still allow a UHECR production (Kotera et al., 2015).

3.2 Interaction and Cooling Rates

Accelerated cosmic rays interact with the background of baryons or photons, producing charged pions ($p+p(\gamma)\rightarrow\pi^{\pm}$) which decay into muons and neutrinos ($\pi^{\pm}\rightarrow\mu^{\pm}+\nu_{\mu}(\bar{\nu}_{\mu})$). The muons further decay into electrons and neutrinos ($\mu^{\pm}\rightarrow e^{\pm}+\nu_{\mu}+\bar{\nu}_{\mu}+\nu_{e}(\bar{\nu}_{e})$). Each step of this process may be suppressed by particle cooling.

A proton with Lorentz factor $\gamma_{p}$ interacts with the photon field of spectrum $n(\varepsilon)=dn/d\varepsilon$ on a characteristic timescale given by

[TABLE]

where $\sigma_{p\gamma}$ is the cross section of photopion production. The cooling time is $t_{\,p\gamma}=t_{p\gamma,\,\rm int}/\kappa_{p\gamma}\sim\left(n_{\gamma}\,\sigma_{p\gamma}\,\kappa_{p\gamma}\,c\right)^{-1}$, where $\sigma_{p\gamma}\kappa_{p\gamma}\sim 10^{-28}\,\rm cm^{-2}$ is the inelastic component of the $p\gamma$ interaction cross section. The rate for hadronuclear interaction is likewise given by

[TABLE]

where $\sigma_{\rm pp}\sim 10^{-25}\,\rm cm^{2}$ (at energies $\sim 10^{18}$ eV) and $\kappa_{\rm pp}\sim 0.5$ (e.g., Eidelman et al., 2004). The total interaction rate of protons, due to both $pp$ and $p\gamma$ processes, is then $t_{p}^{-1}=t_{pp}^{-1}+t_{p\gamma}^{-1}$.

The gyro radii of protons in the magnetic field of the nebula, $r_{L}$, is comparable or larger than the nebula size,

[TABLE]

Synchrotron or adiabatic losses of cosmic rays crossing the nebula are thus generally only marginally important. Inverse Compton cooling at such high energies is furthermore suppressed and negligible due to the Klein-Nishina effect.

Protons travel freely when their crossing time is shorter than the interaction time, $t_{\rm cross,ej}<t_{p}$, as occurs after a characteristic timescale

[TABLE]

where we have assumed that at late times $pp$ interactions dominate over $p\gamma$ interactions (however, note that $p\gamma$ interactions are included in our full numerical calculations).

Figure 2 shows the proton interaction timescales as a function of time since explosion. In Case I, the dominant process for protons is hadronuclear interaction. In Case II, the photopion production with the thermal photon background becomes important at late times, t\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}t_{d,0}^{n}. In both cases, the ejecta becomes optically thin at roughly one week, after which time the accelerated protons freely escape to infinity.

3.3 Ultrahigh Energy Cosmic Ray Production

Figure 3 shows the maximum energy of the accelerated UHECR protons at the time when they escape the ejecta freely. The shaded area corresponds to an event with a low ejecta mass $M_{\rm ej}=0.3\,M_{\odot}$, similar to that inferred for FBOTs such as AT2018cow. Depending on the magnetic field of the magnetar, an AT2018cow-like event with P_{i}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}10 ms is a promising source for cosmic rays of energy \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}10^{18} eV (see also Piro & Kollmeier 2016). If the magnetar wind is composed of nuclei with charge $Z$ instead of protons ($Z=1$), then their maximum energy is a factor of $Z$ times higher than for protons and thus in the range of the highest energy cosmic rays observed by Auger (The Pierre Auger Collaboration et al., 2017) and TA (Matthews, 2018), as detailed in Fang et al. (2012).

For comparison, millisecond magnetars born in the normally considered class of superluminous supernovae (SLSNe), for which the ejecta mass is typically much higher ($M_{\rm ej}\sim 10\,M_{\odot}$), only allow the escape of UHE protons for magnetar parameters to the left of the blue-dashed curve. This is because massive ejecta shells require longer to become optically thin, delaying the time of cosmic ray escape and thus reducing the pulsar voltage at this epoch. White data points show magnetar parameters fit to the sample of SLSNe in Nicholl et al. (2017), roughly one third of which appear to be promising UHECR sources.

We conclude this section by estimating the UHECR energy budget of FBOTs similar to AT2018cow. Assuming $P_{i}=2-10\,\rm ms$ and that a fraction $\eta=0.1$ of the magnetar’s rotational energy goes into cosmic rays, the UHECR energy yield per event is $\left(0.2-5\right)\times{10}^{50}$ erg. The total rate of FBOTs is $\sim 4-7\%$ of the core-collapse supernovae rate (Drout et al., 2014). However, the rate of the most luminous members of this class like AT2018cow is likely lower (M. Drout, private communication), perhaps only $\sim 0.5-1\%$ of the core collapse supernova (CCSN) rate, corresponding to a estimated volumetric rate of AT2018cow-like events of $600-1200$ Gpc${}^{-3}\,$yr*-1*. Using a FBOT rate with $\sim 1$% of the CCSN rate, the integrated cosmic ray luminosity density from AT2018cow-like FBOTs is thus roughly estimated to be $\sim\left(0.3-6\right)\times 10^{44}\,\rm erg\,Mpc^{-3}\,yr^{-1}$. This is comparable to the UHECR luminosity density, which is in the order of $10^{44}\,\rm erg\,Mpc^{-3}\,yr^{-1}$ (e.g., Murase & Fukugita, 2019).

4 High-Energy Neutrino Production

4.1 IceCube Observation of AT2018cow

Two IceCube neutrino track events with $\sim 2^{\circ}$ angular resolution were found in spatial coincidence with AT2018cow during a $3.5$ day period between the last non-detection and the discovery, corresponding to a $1.8\,\sigma$ chance coincidence. Assuming an $E^{-2}$ spectrum, a time-integrated $\nu_{\mu}$ flux upper limit $\left(E_{\nu}^{2}\,J_{\nu}\right)_{\rm UL}^{\rm IC}=6.1\times 10^{-2}\,\rm GeV\,cm^{-2}$ is found for this observation period at 90% CL222http://www.astronomerstelegram.org/?read=11785.

Neutrinos with TeV-PeV energies lie in the best sensitivity window for the IceCube Observatory, and are predicted to arrive in the first few hours in our model. Using the neutrino effective area of the IceCube Observatory with its complete configuration of 86 string detectors 333https://icecube.wisc.edu/science/data/PS-IC86-2011 (Aartsen et al., 2014), we compute the detector sensitivity at the declination of AT2018cow as a function of energy. It is shown by the grey dashed curve in Figure 6. Point sources with fluxes above the curve, assuming that they follow an $E^{-2}$ spectrum over a half decade in energy, are excluded at 90% C. L.

4.2 Competitive Cooling and Decay of Pions and Muons

Charged pions and muons decay into neutrinos following their lab-frame lifetimes, $\gamma_{x}\,\tau_{x}$, where $x$ denotes either $\pi$ or $\mu$, and $\tau_{\pi^{\pm}}=2.6\times 10^{-8}\,\rm s$ and $\tau_{\mu}=2.2\times 10^{-6}\,\rm s$ are the rest-frame lifetimes. These lifetimes can be long enough to allow cooling processes to occur. Mesons and muons interact with the hadron background at a rate given by

[TABLE]

Mesons and muons can also interact with the photon background, but as in the case of protons the smaller $\pi\gamma$ cross section renders the $\pi p$ interaction more important, especially at early times.

Figure 4 compares the lifetime of a pion (muon), which is assumed to carry away 20% (15%) of the energy of the parent proton, to the cooling time due to $\pi p$ interactions. This shows that most pions interact with the ejecta baryons, rather than decay, in the first day or longer.

Pion-proton interactions produce lower-energy pions and protons, which then undergo further interactions with background particles. Eventually, higher-order pions reach sufficiently low energies to enable their decay into neutrinos. This happens once $t_{\pi p}>\gamma_{\pi}\tau_{\pi}$, as occurs after a time

[TABLE]

This effect introduces a break at time $t$ in the neutrino spectrum at a characteristic energy (Murase et al., 2009)

[TABLE]

In addition to the higher-order products, multiple pions and other types of mesons such as kaons may be produced from each $pp$ and $p\gamma$ interaction. We take these into account in our numerical simulations presented in Section 4.3.

Finally, similar to protons, the gyroradii of pions and muons are comparable to the nebula size. Synchrotron cooling thus becomes important if $r_{L}/R_{\rm n}\ll 1$, as occurs for high nebular magnetization \epsilon_{B}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}10^{-2}.

4.3 Numerical Procedure

We compute the neutrino flux from magnetar-powered supernovae with low-ejecta masses according to the following numerical procedure. At each time $t$ following the explosion, the ejecta radius $R_{\rm ej}$, temperature of the thermal background $T_{\rm th}$, number density of ejecta baryons $n_{b}$, and strength of the nebular magnetic field $B_{n}$, are calculated from equations 1$-$5. Protons of energy $E_{p}$ determined by equation 18 are injected into the ejecta. Their interaction with the baryon ejecta is calculated using a Monte-Carlo approach employing the hadronic interaction model EPOS-LHC (Pierog et al., 2015) as in Fang et al. (2012). Photomeson interactions between cosmic-ray protons and the thermal background is computed based on SOPHIA (Mücke et al., 2000) through CRPropa 3 (Alves Batista et al., 2016).

As we are particularly interested in the earliest phases of the transient (near the time of the observed neutrino coincidence), we take a different approach in treating the pion-proton ($\pi p$) interaction from previous works (e.g. Murase et al. 2009; Fang et al. 2016). The suppression of neutrino production due to $\pi p$ interaction is usually described by a suppression factor $f_{\rm sup}\equiv\min\left(t_{\rm cool}/\gamma_{\pi}\,\tau_{\pi},1\right)$, with $\gamma_{\pi}$ and $\tau_{\pi}$ being the Lorentz factor and the rest-frame lifetime of the charged pion. This analytical approach however misses the secondary and higher-order pions produced by the interaction, which have lower energy and thus decay into neutrinos more easily than their parents. This is a secondary effect, but can play an important role in neutrino production in dense environments.

For reference, Figure 5 shows the distribution of outgoing particles from the interaction of an pion of energy $E=10^{19}$ eV with a proton at rest as calculated using the EPOS-LHC model. Roughly half of the energy of the primary pion is carried away by charged mesons, while $15\%$ is in the form of baryons. These products will ultimately generate neutrinos.

To account for neutrinos produced by secondary pions, we record all charged meson and baryon products from the $pp$ and $p\gamma$ interaction. At any time, each of these intermediate products may undergo one of four processes, depending on the background densities at the current step: 1) cooling by $\pi p$ interaction; 2) cooling by synchrotron radiation; 3) decay into a muon and a muon neutrino (with a 100% branching ratio for $\pi^{\pm}$ and a 63% branching ratio for $K^{\pm}$); and 4) free propagation. The products of $\pi p$ interaction are computed using the EPOS-LHC model (Pierog et al., 2015). All secondaries and their higher-order products are tracked until either their energy falls below TeV (where the atmospheric background dominates over astrophysical sources), or they escape the source without further interaction.

4.4 Neutrino Fluence of AT2018cow

Figure 6 presents the neutrino fluence in the two fiducial models. In Case I, we show results in time intervals normalized to the magnetar spin-down time $t_{\rm sd}$ (roughly one day). Neutrino production is low at both early times $t\ll t_{\rm sd}$, when the bulk of the cosmic rays have yet to be injected, and at late times $t_{\rm pp}>t_{\rm cross,ej}$, when $pp$ interactions become inefficient. Most neutrinos are generated around the time $t_{\rm sd}$, when the magnetar releases most of its rotational energy and a sufficiently dense baryon background still exists for pion production. The turn-over in the neutrino spectrum is determined at early times by the break energy $E_{\nu,b}$ (eq. 26), while at late times the break is determined by the maximum cosmic ray energy $E_{\rm CR}$ (eq. 18). The spectral index before the break is a convolution of the energy distribution of charged pions from the $pp$ and $p\gamma$ interaction (similar to that from $\pi p$ as shown in Figure 5) and the history of particle injection, as described by equation 19.

The bottom panel of Figure 6 shows our results for Case II ($P_{i}=2$ ms, $B_{d}=2\times 10^{15}$ G). The spin-down time in this case is much shorter $t_{\rm sd}=0.01$ d, whereas the neutrino break energy exceeds 1 TeV only after times $t=0.03\,M_{\rm ej,-0.5}^{1/3}\,\beta_{\rm ej,-1.3}^{-1}\,\rm d$. As a result, most cosmic rays are injected too early to generate neutrinos in the energy range of interest. As the magnetar releases most of its energy before the environment becomes optically-thin, significant TeV-PeV neutrinos are produced in the first 0.5 day due to the $\pi p$ interaction. The neutrino spectrum at the earliest epoch features two peaks; the low-energy bump is from pion decay, while the tail at high energies arises from the decay of short-lived mesons other than $\pi^{\pm}$. After about one day, the ejecta becomes sufficiently dilute that meson cooling is no longer severe and the neutrino flux becomes maximal when the suppression is not important. Then, the neutrino spectrum comes to resemble that in Case I. This kind of time evolution of neutrino spectra owing to meson cooling in magnetars was first found in Murase et al. (2009) and Fang et al. (2014).

Figure 7 shows the total neutrino luminosity, $L_{\nu}=\int_{E_{\nu}>\rm TeV}\,E_{\nu}(dN_{\nu}/dE_{\nu})dE_{\nu}dt$, at energies $>1$ TeV as a function of time since explosion. In both Cases I and II the light curve obeys $L_{\nu}\propto t^{-2}$ at times t_{\pi p,0}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}t\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}t_{p,0}^{\rm ej}, similar to the optical and X-ray light curves of AT2018cow (Fig. 1; e.g. Margutti et al. 2018). In the time window when the ejecta is still optically thick to cosmic rays, yet after the pion cooling no longer suppresses neutrino production, the neutrino luminosity tracks the spin-down power of the magnetar, $L_{\rm sd}\propto t^{-2}$ (eq. 6).

The origin of the time dependence of the neutrino luminosity is more complex at earlier and later times. At early times $t\ll t_{\pi p,0}$, both the cosmic ray injection rate $dN/dt$ and maximum energy $E_{p}$ decrease in time as $\left(1+t/t_{\rm sd}\right)^{-1}$, while the suppression factor due to $\pi$ cooling rises as $f_{\rm sup}\propto t^{3}\left(1+t/t_{\rm sd}\right)$. If, as in Case II, the majority of the cosmic ray energy is injected at early times when the system is still optically thick to protons and pions, then the time evolution is also influenced by neutrinos released from $p\gamma$ interaction, $\pi p$ interaction and kaon decay. At late times $t\gg t_{p,0}^{\rm ej}$, the effective optical depth decreases as $\tau_{\nu}\propto t^{-2}$, and thus the neutrino light curve obeys a steeper decay, $L_{\nu}\sim L_{\rm sd}\tau_{\nu}\propto t^{-4}$ (see also Murase et al., 2009; Fang, 2015).

Our predictions for the neutrino fluence from AT2018cow is well below the upper limits placed by the IceCube Observatory for both models, supporting a conclusion that the two detected events are from the background rather than of astrophysical origin. A future event otherwise similar to AT2018cow but occurring $\sim 5$ times closer (at a distance of \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}15 Mpc), would be a promising IceCube source. Prospects may be better for future neutrino telescopes with sensitivity in the EeV energy range, enabling a direct test of FBOTs as cosmic particle accelerators.

4.5 Neutrino Flux of the SLSNe population

The neutrino flux in our models peak around the time $t_{\pi p,0}$, after which pions have sufficient time to decay into neutrinos before cooling. We therefore can estimate the peak flux as

[TABLE]

For a magnetar with $t_{\rm sd}\ll t_{\pi p,0}$ (as in our Case II), the peak flux is somewhat higher than this estimate due to additional contributions from secondaries.

Under the assumption that the peak neutrino flux lasts for a duration $\Delta t\sim t_{\pi p,0}$, the corresponding fluence can be estimated by

[TABLE]

Figure 8 shows the peak fluence for a sample of SLSNe (Nicholl et al., 2017) estimated by equation 28 using the magnetar parameters best-fit to the optical light curves and the distances of the sources. For comparison, we show the discovery potential of a time-dependent search by IceCube for a source at the declination $\delta=16^{\circ}$ (Aartsen et al., 2015) for different source time integrations $\sim t_{\pi p,0}$. The IceCube point-source sensitivity also depends on the declination of the source and is greatest for $\delta\approx 0$ (Aartsen et al., 2017).

4.6 Integrated Neutrino Flux of FBOT Population

The integrated neutrino flux over cosmological distances of magnetar-powered FBOTs with properties similar to AT2018cow is given by

[TABLE]

where $z$ is the source redshift, $|dz/dt|=H_{0}\,(1+z)\,\sqrt{\Omega_{M}\,(1+z)^{3}+\Omega_{\Lambda}}$ taking a flat $\Lambda$CDM with $\Omega_{M}=0.308$ and $H_{0}=67.8\,\rm km\,s^{-1}\,Mpc^{-1}$ (Planck Collaboration et al., 2016). $\rho(z)=\rho_{0}\,g(z)$ is the source birth rate at given redshift. We assume that AT2018cow-like FBOT events occur at $\sim 1\%$ of the core collapse supernova rate ($\rho_{0}=10^{-6}\,\rm Mpc^{-3}\,yr^{-1}$), and track the cosmological star-formation history (Hopkins & Beacom, 2006) with $g(z)\propto(1+z)^{3.4}$ at $0<z<1$, then $g(z)\propto(1+z)^{-0.3}$ up to $1<z<4$, and $(1+z)^{-3.5}$ at $z>4$. Figure 9 shows the integrated neutrino flux, calculated separately using neutrino fluences based on Case I and Case II, respectively. The peak flux is similar to that computed in Murase et al. (2009) after taking into account the difference in the source rates. This is because the peak fluence is determined by the time when pions start to decay rather than interact, and is thus insensitive to the factor $\beta_{\rm ej}^{3/4}\,M_{\rm ej}^{-1/4}$ (see Equation 25 of this work and Equation 4 of Murase et al., 2009). The flux is consistent with the extreme-high-energy upper limits from IceCube at 90% C.L. (Aartsen et al., 2018). Interestingly, if all FBOTs were as luminous as AT2018cow (i.e. if their rate was 4-7% of the CCSN rate; Drout et al. 2014), then the integrated neutrino flux would have overproduced the IceCube limit.

4.7 Gamma-Ray Emission from FBOTs and SLSNe

Neutral pions created by UHE proton interactions produce high-energy $\gamma$ rays via $\pi^{0}\rightarrow 2\,\gamma$. In $pp$ and $p\gamma$ channels, the energy passed by protons to $\gamma$ rays is about $2/3$ and $4/3$ times that given to neutrinos, respectively. High-energy gamma rays will quickly undergo pair production with low-energy photons in the nebula ($\gamma\gamma\rightarrow e^{\pm}$). For GeV to TeV $\gamma$ rays, the peak of the $\gamma\gamma$ pair production cross section lies in the optical to soft X-ray range. The optical depth of the background photon to high-energy $\gamma$ rays is approximately given by

[TABLE]

where $\sigma_{\gamma\gamma}=3/16\,\sigma_{T}$ is the $\gamma\gamma$ pair production cross section and $u_{\rm th}\sim L_{\rm sd}/\left(4\pi\,R_{\rm ej}^{2}\,c\right)$ is the thermal photon density at late times. Thus high-energy $\gamma$ rays can not escape from the ejecta freely until after a time

[TABLE]

The pair production and inverse Compton scattering of the resulting electrons leads to an electromagnetic cascade (e.g. Metzger et al. 2014; Metzger & Piro 2014; Murase et al. 2015). Figure 10 presents the optical depth of the nebula and ejecta to $\gamma$ rays at different energies around day 10. For MeV, sub-GeV, and $>$ GeV $\gamma$ rays, the dominate energy loss process is Compton scattering ($e\gamma\rightarrow e\gamma$), Bethe-Heitler pair production ($p\gamma\rightarrow p\gamma$) and $\gamma\gamma$ pair production respectively. The effective optical depth is defined as $\tau_{\rm eff}=Y_{\rm b}\,\sigma\,\kappa$, where $Y_{\rm b}=3\,M_{\rm ej}/\left(4\pi R_{\rm ej}^{2}\,\mu_{e}\right)$ is the column density of the ejecta, and $\mu_{e}\approx 2$ is the mean molecular weight. The cross section $\sigma$ and inelasticity $\kappa$ of the interactions depend on the $\gamma$-ray energy (equations 40, 46, 48 of Murase et al. 2015, also see Dermer & Menon 2009 and references therein). Figure 10 shows that comparing to high-energy $\gamma$ rays, MeV to sub-GeV $\gamma$ rays have a better chance to be observed at early times. Note that the optical depths due to photon-matter interactions are estimated assuming the spherical geometry. If X-rays from the engine are observed, a more complicated geometry seems necessary (Margutti et al., 2018), and then the escape of gamma rays is likely to be easier in such more realistic setups. Also, we remark that the optical depth to the two-photon annihilation process decreases as energy, so UHE photons can escape from the system.

Murase et al. (2015) computed the $\gamma$-ray and hard X-ray emission from a magnetar-powered CCSN, taking into account details of electromagnetic cascades. Applying the numerical code to AT2018cow-like events to a low-mass ejecta with a magnetar, the GeV gamma-ray flux is estimated to be at the level of $10^{-12}\,\rm erg\,cm^{-2}\,s^{-1}$ at the distance of AT2018cow. If the ejecta is asymmetric, as suggested by X-ray observations, MeV-GeV $\gamma$ rays are likely to leak out from low-density regions in the similar geometry and be observed.

Renault-Tinacci et al. (2018) searched for GeV $\gamma$ rays in the directions of a SLSNe sample with the Fermi-LAT data and found no signals. Assuming that SLSNe population are equally luminous and have an $E^{-2}$ spectrum, Renault-Tinacci et al. (2018) concluded an upper limit at 95% C. L. to the $1-10$ GeV luminosity $L_{\gamma,1-10\,\rm GeV}<9.1\times 10^{41}\,\rm erg\,s^{-1}$. Our scenario is consistent with this limit, as the total spin-down power is below $3\times 10^{41}\,\rm erg\,s^{-1}$ after 100 days.

For AT2018cow, a gamma-ray non-detection was reported by Fermi LAT.444http://www.astronomerstelegram.org/?read=11808 during an one-week interval from day 3 to day 10 since the epoch of detection. The flux upper limit at 10 day is approximately $\sim{10}^{-11}~{}{\rm erg}~{}{\rm cm}^{-2}~{}{\rm s}^{-1}$ (cf. Renault-Tinacci et al., 2018; Murase et al., 2018). HESS reported an 95% C. L. upper limit of $5\times 10^{-13}\,\rm cm^{-2}\,s^{-1}$ based their observation in the third week.555http://www.astronomerstelegram.org/?read=11956 Both are consistent with our model.

5 Discussion and Conclusions

Rapidly-spinning magnetars have been proposed as the engines responsible for powering the optical light curves of superluminous supernovae. The recent discovery of AT2018cow, the first local example of a superluminous transient with a low ejecta mass (a so-called FBOT), enabled the discovery of coincident time-variable X-ray emission consistent with the presence of a central engine, given that X-rays escape due to the asymmetry of the ejecta (e.g. Prentice et al. 2018; Perley et al. 2019; Margutti et al. 2018; Ho et al. 2019). The engine behind AT2018cow released more than $\sim 10^{50}$ ergs over a timescale of less than a few days, behind a low mass ejecta shell M_{\rm ej}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{<}}}0.3M_{\odot}. If the engine is a millisecond magnetar spinning down in isolation (i.e. neglecting fall-back accretion), then fits to the light curve require an initial spin period of $\sim 2-10$ ms and a dipole magnetic field of $\sim 1-2\times 10^{15}$ G.

For the same reasons they are capable of powerin large luminosities, the extreme magnetars responsible for SLSNe and FBOTs also possess strong electromagnetic potentials, making them potential sites for the acceleration of relativistic particles and even UHECRs. Events with particularly low ejecta masses, such as AT2018cow, provide a way to directly view X-rays from the magnetar nebula (through prompt photo-ionization of the ejecta shell; Metzger et al. 2014) and allow for the timely escape of accelerated cosmic rays relative to the bulk of the normal SLSNe population.

Motivated by the discovery of AT2018cow as a local example of an FBOT which shows direct evidence for a central engine, we have calculated the interaction of high-energy cosmic rays accelerated near the magnetar with background particles of the nebula and ejecta. We take into account the time-evolving thermal and non-thermal radiation field of the nebula, and track all primary and higher-order interaction products down to TeV energies.

Our results are largely insensitive to the modeling of the background radiation fields. The photopion interaction is dominated by the optical/UV emission which is directly observed, and the proton-proton interaction depends on the density of baryons, which is well constrained by the rise time of the explosion and the spectroscopically observed ejecta velocity. Our results do, however, depend on whether the engine behind AT2018cow is truly a magnetar (e.g., as opposed to an internal shock from CSM interaction or an accreting black hole), and whether millisecond magnetars are indeed efficient particle accelerators. Nevertheless, the direct detection of time-variable X-rays from AT2018cow provides greater confidence in the central engine scenario than was available from previous SLSNe samples (which generally show no coincident X-ray emission, likely due to photoelectric attenuation by the larger ejecta shells; Margutti et al. 2018; Margalit et al. 2018).

UHECR protons accelerated by millisecond magnetars in FBOTs and SLSNe escape the source with characteristic energies $\sim$ EeV. In addition, heavy elements may be synthesized efficiently in magnetar winds (Metzger et al., 2011) and destroyed inside nebulae (Horiuchi et al., 2012; Murase et al., 2014). Nuclei heavier than protons could tap $Z$ times more energy from the same electric potential. Their energy losses are dominated by photo-disintegration and hadronuclear interaction. Fang et al. (2012) shows that they may escape from a massive ejecta of CCSNe with $E>10^{20}$ eV if the pulsar has a millisecond spin period and surface dipole field in the range $10^{12}-10^{13}$ G. Neutrino production by heavy nuclei depends on the interaction channel. It could be comparable to that by protons when hadronuclear interactions dominate (Fang, 2015), but lower when photo-disintegration dominates. We will leave a more detailed study to a future work. Engine-driven supernovae, including jet-driven ones, have been discussed as the sources of UHECRs, where particle acceleration sites have been attributed to internal shocks in outflows, external forward and reverse shocks (e.g., Murase et al., 2008; Wang et al., 2007; Chakraborti et al., 2011; Zhang et al., 2018). Our physical model is different in the sense that UHECR acceleration occurs inside pulsar wind nebulae.

For AT2018cow, the cosmic ray energies that can be achieved are lower than for the bulk of the SLSNe population considered previously. On the other hand, if relatively common events similar to AT2018cow provide the dominant UHECR source just above the ankle energy, a tail of rarer more powerful FBOTs or SLSNe (e.g. those born with initial spin periods close to their minimum break-up value $P_{0}\sim 1$ ms instead of P_{0}\mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}}2 ms) could dominate the UHECR budget at the highest energies.

Magnetars with birth properties similar to those required to power AT2018cow may be in some ways optimal for neutrino production. This is because of the apparent coincidence that the magnetar spin-down timescale is comparable to the narrow time window within which the optical depth for cosmic ray interaction is still high, but the secondary pions are no longer efficiently cooled before decaying. FBOTs similar to AT2018cow can therefore be in principle ideal targets for neutrino telescopes.

Depending on the uncertain volumetric rate of FBOTs giving birth to millisecond magnetars with properties similar to those required to explain AT2018cow, we find that such a population could explain $\sim 10-20\%$ of the IceCube astrophysical neutrino background. However, given the sensitivity of current-generation facilities, the detection of individual FBOTs will be challenging without the fortuitous discovery of a source located several times closer than AT2018cow (which, however, was itself already much nearer than the previous cosmological populations of FBOTs and SLSNe). Nevertheless, IceCube neutrino events can be used to optimize the strategies of transient telescopes such as Zwicky Transient Facility for observing FBOTs.

The same interactions giving rise to neutrinos in FBOTs and SLSNe also inevitably give rise to high energy gamma-rays. However, because of the high pair creation optical depth created by the thermal photons of the transient, high-energy gamma-rays require months to years to escape from the magnetar nebula and thus to become detectable from Earth. Nevertheless, well-studied optical light curves of some SLSNe show evidence for “leakage” of engine energy at late times \mathrel{\hbox{\hbox to0.0pt{\hbox{\lower 4.0pt\hbox{\sim}}\hss}\hbox{>}}} months (e.g. Nicholl et al. 2018), which is not observed to be escaping in the soft X-ray band (Margutti et al., 2018; Bhirombhakdi et al., 2018) and therefore may emerge first in the hard X-ray or gamma-ray band.

Our work motivates FBOTs and SLSNe as potentially promising sources for future neutrino and $\gamma$-ray telescopes. The predicted neutrino spectrum peaks at $100$ PeV-EeV, which can be favorably probed by both EeV neutrino detectors such as Askaryan Radio Array (Ara Collaboration et al., 2012), ARIANNA (Barwick et al., 2017), GRAND (GRAND Collaboration et al., 2018) POEMMA (Olinto et al., 2017), as well as TeV-PeV detectors such as KM3NeT (Margiotta, 2014) and IceCube-Gen2 (IceCube-Gen2 Collaboration et al., 2014). Depending on the structure of the ejecta and the distance of the source, the $\gamma$ rays may be observed by existing wide-field telescopes such as Fermi and HAWC, and the next-generation Cherenkov Telescopes Array (CTA Consortium, 2017).

The authors thank Tanguy Pierog for helpful discussions on the pion-proton interaction calculations, Erik Blaufuss and the IceCube Collaboration for useful feedback. B.D.M. is supported by NASA (grants HST-AR-15041.001-A, 80NSSC18K1708, a80NSSC17K0501). K.M. is supported by the Alfred P. Sloan Foundation and NSF grant no. PHY-1620777. I.B. is grateful for the generous support of the University of Florida. K.K. is supported by the APACHE grant (ANR-16-CE31-0001) of the French Agence Nationale de la Recherche.