FRB 121102: A Starquake-induced Repeater?

Fast radio bursts (FRBs) are mysterious millisecond-duration radio flashes with high flux densities and prominent dispersive features (Lorimer et al. 2007; Keane et al. 2012; Thornton et al. 2013; Masui et al. 2015; Ravi et al. 2016; Bannister et al. 2017; Bhandari et al. 2017; Caleb et al. 2017). The observed large values of dispersion measure (DM) are in the range of ∼100–2600 pc cm⁻³, which indicate that FRBs are probably of extragalactic or even cosmological origins (e.g., Katz 2016a; Scholz et al. 2016). These transient phenomena stimulate the interests of astrophysicists significantly; this is especially true of FRB 121102, which is the only repeater that has been detected so far and has an estimated burst energy ∼10^37–38 erg (Spitler et al. 2014, 2016). The optical counterpart of the repeater has been identified as a host faint star-forming dwarf galaxy, which is at a redshift of z = 0.193 (Chatterjee et al. 2017; Kokubo et al. 2017; Tendulkar et al. 2017).

A persistent radio source that is thought to be associated with the repeater was identified at a distance of ≲40 pc from the FRB location (Chatterjee et al. 2017; Marcote et al. 2017). Ofek (2017) also found 11 source candidates with luminosities of νL_ν > 3 × 10³⁷ erg s⁻¹, which are spatially associated with disks or star-forming regions of galaxies, rather than in the galactic center. The persistent radio source is likely to be a pulsar wind nebula (Beloborodov 2017; Dai et al. 2017; Kashiyama & Murase 2017). With an active pulsar producing bursts repeatedly, ejecta or ultra-relativistic electron/positron pair winds may sweep up and heat the nebula that emits synchrotron radio emissions. Additionally, Waxman (2017) calculated some stringent constraints on the persistent source's age. The local environment of FRB source would be tested by the variation of DM, which has not shown significant evolution (Yang & Zhang 2017).

It is proposed that FRBs are highly likely to be associated with pulsars, and more than a few efforts have been made to understand their origins. For instance, FRB is supposed to result from a pulsar's magnetosphere suddenly combed by a nearby cosmic plasma stream (Zhang 2017). Also, Dai et al. (2016) suggested that the repeater originated from a highly magnetized pulsar traveling through asteroid belts. Alternatively, in a neutron star (NS)-white dwarf binary system, the accreted magnetized materials may trigger magnetic reconnection that accounts for FRBs (Gu et al. 2016). Another possibility is that the radio emission is produced by the interaction between a highly relativistic flow and a nebula powered by a newborn millisecond magnetar (Murase et al. 2016; Beloborodov 2017; Dai et al. 2017). This process might couple with a long gamma-ray burst or an ultraluminous supernova (Metzger et al. 2017). FRB are also interpreted by the model of a supergiant pulse or giant flare from a young pulsar or magnetar (Popov & Postnov 2010; Kulkarni et al. 2014; Cordes & Wasserman 2016; Katz 2016b).

Very recently, 14 bursts above the threshold of 10 sigma in two 30-minute scans were detected by Breakthrough Listen Digital Backend with the C-band receiver at the Green Bank Telescope (GBT; Gajjar et al. 2017). In this paper, we propose that this repeating burst may arise from a pulsar's starquake. The burst energy and waiting time distributions, as well as the time decaying of the seismicity rate, are shown in Section 2. The scenarios of possible origins will be discussed in Section 3. In Section 4, we make our conclusions.

Considering the unity of telescope selection and wavebands (4–8 GHz), and the completeness of the continuous observations (290 minutes), here we adopt the parameters of repeating bursts from the latest continuous monitoring GBT observations Gajjar et al. (2017), rather than including the Arecibo events (Spitler et al. 2016) and previous GBT events (Scholz et al. 2016). Parameters of these 14 bursts are shown in Table 1.

We generate statistics on the observational parameters and it the results demonstrate the consistency of the burst rate as a function of burst energy. Considering that the binning width can affect the fitting result, here we calculate the cumulative distribution to avoid this problem. The burst energy is proportional to the observed energy density. With a power-law distribution for the number distribution of burst energy N(E) ∝ E^−α, the cumulative distribution can be obtained,

$$\begin{eqnarray}&&N(\gt E)\propto {\int }_{E}^{\infty }{E}^{-\alpha }{dE}\propto {E}^{-\alpha +1}.\end{eqnarray} \tag{ 1 }$$

The fluctuations of events for the cumulative distribution are assumed to follow a random statistic, $\sigma (E)=\sqrt{N(\gt E)}$. The cumulative energy distribution of events for each burst energy (energy density) is well fitted by a power law with an index α_E = α − 1 = 1.16 ± 0.24, shown in Figure 1. This power-law distribution is consistent with the Gutenberg–Richter power law (i.e., N(E) ∝ E⁻², Gutenberg & Richter 1956) which is a well-known earthquake law.

Zoom In Zoom Out Reset image size

Furthermore, the statistics of waiting times contain much significant information concerning the occurrences and correlations of events. The waiting time Δt is defined as the interval time between the adjacent detected FRB events in the continuous monitoring observation. For a simple Poisson process, the cumulative distribution of waiting time can be described by a simple exponential function,

$$\begin{eqnarray}&&N(\gt {\rm{\Delta }}t)\propto {e}^{-\lambda {\rm{\Delta }}t},\end{eqnarray} \tag{ 2 }$$

where λ is the burst rate, which is constant. Also, with an assumed Gaussian distribution of the number distribution N(Δt), the cumulative distribution of the waiting time is

$$\begin{eqnarray}&&N(\gt {\rm{\Delta }}t)\propto {\int }_{{\rm{\Delta }}t}^{\infty }{\exp }^{-\tfrac{{({\rm{\Delta }}t-\tau )}^{2}}{{\sigma }^{2}}}d({\rm{\Delta }}t)\propto 1-\mathrm{erf}({\rm{\Delta }}t).\end{eqnarray} \tag{ 3 }$$

As shown in Figure 2, the cumulative distribution of the waiting time is plotted, fitted by the exponential function with λ = (3.05 ± 0.48) × 10⁻³ s⁻¹ and Equation (3) with τ = (1.13 ± 0.20) × 10³ s and σ = (1.03 ± 0.16) × 10³ s. The waiting time distribution can be represented by a simple Poissonian or Gaussian distribution, which are extracted from statistics of earthquakes (Pepke et al. 1994; Leonard et al. 2001). There may be additional burst events with a detection threshold of around 7 over the full bandwidth of 4 GHz. However, they are not listed in Gajjar et al. (2017) because they are relatively weak and need more analysis. These bursts might have a narrow frequency spread and thus do not show high signal to noise ratio. This may not affect our model significantly while affecting the fitting parameters. Earthquakes from different regions are regarded as random processes of independent and uncorrelated events, while aftershocks that occur in shorter time intervals are correlated.

Zoom In Zoom Out Reset image size

Additionally, combined with magnitudes of burst, the burst rates are plotted with different binning widths of time and fitted by a power law, shown in Figure 3. For short-range temporal correlations between earthquakes, the time decay of the seismicity rates of an aftershock sequence can be interpreted by an empirical relationship, i.e., Omori law (Omori 1894; Utsu 1961; Utsu et al. 1995). By the Omori law, the seismicity rates decay with time and can be expressed by a power law,

$$\begin{eqnarray}&&n(t)=\displaystyle \frac{K}{{(C+t)}^{p}},\end{eqnarray} \tag{ 4 }$$

where n(t) is the seismicity rate and the time decay rate of seismicity is controlled by the third constant p, which typically falls in the range of 0.9–1.5 (Utsu et al. 1995). For an aftershock sequence, the variation of p may be controlled by the structural heterogeneity, stress, and temperature (Utsu et al. 1995). We fixed the parameter C at 1 minute, therefore the number rate at t = 0 is close to K, which is the peak value of the number rate curve. Thus only one parameter p, which controls the decay rate, is estimated. With the occurrence time of 14 FRB events, the burst rates can be well fitted by the power law determined by an index p value of 1.42 ± 0.24, which is consistent with the index for earthquakes. Therefore, the bursts are followed by an earthquake-like aftershock sequence.

Zoom In Zoom Out Reset image size

As we see from Section 2, the FRB repeater exhibits several features commonly found in earthquakes. In fact, nonlinear dissipative systems always show self-organized criticality (SOC) behaviors. With a solid crust or stiff equation of state (EOS) of a pulsar, a star can build up stresses that make the crust crack and adjust the stellar shape in order to reduce its deformation. The characters for these processes of statistical independence, nonlinear coherent growth, and random rise times are consistent with an SOC system (Aschwanden 2011). According to the SOC theory, NS quakes are giant catastrophic events like earthquakes, and are probably accompanied by global seismic vibrations or oscillations. Hence, starquakes in a pulsar share similar statistical distributions with earthquakes. Here we present two possible scenarios.

One scenario involves a normal NS with a solid crust and a superfluid core, and with a strong toroidal magnetic field (i.e., magnetar). An NS can form a solid crust quickly after its birth. The stellar shape changes from oblate to spherical configuration, and thermal and dynamic responses will induce stresses in the crusts. When stress buildup in the dense solid crust extends beyond a yield point a starquake occurs, with a sudden energy release that supports FRBs. This process also brings a slight change in the moment of inertia, with an abrupt jump of angular frequency that is the so-called pulsar glitch (Ruderman 1969; Link & Epstein 1996). For a highly magnetized NS, the stellar crust couples with the magnetosphere. In this scenario, starquakes induce the magnetic curl or twist ejected into the magnetosphere from the crust in a few milliseconds (Thompson et al. 2002, 2017). Electrons in the magnetosphere are suddenly accelerated to ultra-relativistic velocity by magnetic reconnection (Zhang & Yan 2011) and move along the magnetic field lines, resulting in the production of curvature radiation. The characteristic frequency of the curvature radiation is

$$\begin{eqnarray}&&{\nu }_{{\rm{c}}}=\displaystyle \frac{3c{\gamma }^{3}}{4\pi {R}_{{\rm{c}}}}=7.16\times {\left(\displaystyle \frac{\gamma }{100}\right)}^{3}\left(\displaystyle \frac{10\,\mathrm{km}}{{R}_{{\rm{c}}}}\right)\,\mathrm{GHz},\end{eqnarray} \tag{ 5 }$$

where R_c is the curvature radius with a typical value of ∼10 km, and γ is the Lorentz factor of electrons. With a detected FRB frequency ν_c = 4–8 GHz, a γ ≈ 50–100 is required. If generated by such stress-induced reconnections (e.g., Lyutikov 2015), the FRB has an estimated duration timescale

$$\begin{eqnarray}&&{t}_{\mathrm{rec}}\sim \displaystyle \frac{L}{{v}_{{\rm{A}}}}\sim 1\mbox{--}10\,\mathrm{ms},\end{eqnarray} \tag{ 6 }$$

where the scale of the reconnection-unstable zone L ∼ 1–10 km, and the Alfvén velocity is v_A ≃ B/(4πρ_c)^1/2 ∼ 0.01 c, in which ρ_c ≃ 10¹⁴ g cm⁻³ is the average mass density of the crust.

The sudden elastic and magnetic energy release in the crustal stress is estimated to be

$$\begin{eqnarray}&&\delta {E}_{\mathrm{cru}}=4\pi {R}^{2}{h}_{{\rm{c}}}\sigma \delta \varepsilon ,\end{eqnarray} \tag{ 7 }$$

where h_c ≃ R/30 ≃ 0.3 km is the crustal thickness (Thompson et al. 2017), ε is the shear strain, and the total stress including crustal shear stress and Maxwell stress in magnetosphere is

$$\begin{eqnarray}&&\sigma =\sqrt{{(\mu \varepsilon )}^{2}+{\left(\displaystyle \frac{{{BB}}_{z}}{4\pi }\right)}^{2}},\end{eqnarray} \tag{ 8 }$$

in which μ is the shear modulus, B is the surface magnetic field, and B_z is the component of magnetic field perpendicular to the direction of plastic flow. Within the crust, the force balance με ≃ BB_z/(4π) implies that

$$\begin{eqnarray}&&\delta {E}_{\mathrm{cru}}\simeq 4.2\times {10}^{46}\,\left(\displaystyle \frac{{{BB}}_{z}}{{10}^{30}\,{{\rm{G}}}^{2}}\right){\left(\displaystyle \frac{R}{10\mathrm{km}}\right)}^{2}\delta \varepsilon \,\mathrm{erg}.\end{eqnarray} \tag{ 9 }$$

$\delta \varepsilon$ is smaller than ∼10⁻² (e.g., Hoffman & Heyl 2012). Here, the energy release can meet the energy requirements of FRBs, while not for soft gamma repeaters (SGRs). However, a plastic flow can be initiated when the elastic crustal deformation exceeds a critical value, launching a thermo-plastic wave that dissipates the magnetic energy inside the crust (Beloborodov & Levin 2014). This mechanism might bring much more energy from the inner crust, in which it stores an SGR-required magnetic energy of ≳10⁴⁷ erg with a interior magnetic field of ∼10¹⁶ G (Lander 2016). The Ohmic dissipation in this process can be neglected because of the long timescale (Fujisawa & Kisaka 2014). The timescale for the local energy release is t_tw ∼ 4πη/BB_z, where η is the viscosity. If the energy releases quickly enough for an FRB, a viscosity of ∼10²⁶ erg s cm⁻³ is required. In addition, the transition to hydromagnetic instability of the magnetar core may offer larger energy-supporting SGRs (Thompson et al. 2017).

The other scenario involves a newborn strangeon star (SS), which has a stiff EOS (Lai & Xu 2017) and could release more elastic energy than that of the solid crust of a normal NS. At an early age, SS may shrink in volume abruptly due to solidification-induced stress, and a bulk-variable starquake happens. Basically, there are two kinds of quakes in a solid star: bulk-invariable (type I) and bulk-variable (type II) starquakes (Zhou et al. 2004, 2014). Here, a type II starquake is more likely to be dominant because the elastic energy accumulation of a type I quake is not sufficient enough to produce such short time-interval quakes in this quake sequence (see Equation (39) in Zhou et al. 2014 for a test of t = 100 s). The elastic and gravitational energy release during a type II starquake is

$$\begin{eqnarray}&&\delta {E}_{{\rm{g}}}=\displaystyle \frac{3{{GM}}^{2}}{5R}\displaystyle \frac{\delta R}{R}=\displaystyle \frac{3{{GM}}^{2}}{10R}\displaystyle \frac{\delta {\rm{\Omega }}}{{\rm{\Omega }}}\\ &&\simeq \,{10}^{53}{\left(\displaystyle \frac{M}{1.4{M}_{\odot }}\right)}^{2}\left(\displaystyle \frac{10\,\mathrm{km}}{R}\right)\left(\displaystyle \frac{\delta {\rm{\Omega }}}{{\rm{\Omega }}}\right)\mathrm{erg},\end{eqnarray} \tag{ 10 }$$

where G is the gravitational constant, M is the stellar mass, δR/R is the strain, and δΩ/Ω is the amplitude of a glitch (e.g., 10⁻⁹–10⁻⁶; Alpar & Baykal 1994; Alpar et al. 1996). This energy is large enough to support an FRB and possibly associate with an SGR. Bulk-variable starquakes are accompanied by a change of electrostatic energy (Katz 2017a) and some electrodynamic activities in magnetosphere. A giant quake can power energetic relativistic outflow in order to produce the observed prompt emission of short-hard GRBs, and some aftershocks result following observed X-ray flares (Xu et al. 2006). Starquakes may also lead to the magnetic reconnection that accelerates electrons, and these charges move along the magnetic field lines, emitting curvature radiation.

The duration timescale of the magnetic reconnection in this scenario can be obtained from Equation (6). In addition, these short time-interval quakes might be motivated by an initial shock, which is type I quake dominated. The waiting time of the initial shock can be obtained (Zhou et al. 2004),

$$\begin{eqnarray}&&{t}_{{\rm{i}}}=\displaystyle \frac{{\sigma }_{{\rm{c}}}}{\dot{\sigma }},\end{eqnarray} \tag{ 11 }$$

where σ_c is the critical stress and $\dot{\sigma }$ can be denoted as

$$\begin{eqnarray}&&\dot{\sigma }=\displaystyle \frac{3\pi I\dot{P}}{{R}^{3}{P}^{3}}\approx 9.42\times {10}^{27}\,\left(\displaystyle \frac{\dot{P}}{1\,{\rm{s}}\,{{\rm{s}}}^{-1}}\right){\left(\displaystyle \frac{1{\rm{s}}}{P}\right)}^{3},\end{eqnarray} \tag{ 12 }$$

in which I ≈ 10⁴⁵ erg s² is the moment of inertia, P is the rotation period, and $\dot{P}$ is period derivative of the star. Then the waiting time of the initial shock can be written as

$$\begin{eqnarray}\mathrm{log}\left(\displaystyle \frac{{t}_{{\rm{i}}}}{1\,{\rm{s}}}\right) & = & \mathrm{log}\left(\displaystyle \frac{{\sigma }_{{\rm{c}}}}{1\,\mathrm{erg}\,{\mathrm{cm}}^{-3}}\right)+3\mathrm{log}\left(\displaystyle \frac{P}{1\,{\rm{s}}}\right)\\ & & -\mathrm{log}\left(\displaystyle \frac{\dot{P}}{1\,{\rm{s}}\,{{\rm{s}}}^{-1}}\right)-28.\end{eqnarray} \tag{ 13 }$$

The rotation period P ∼ 10 ms for a newborn rapidly rotating SS. From Equation (13), the critical stress is estimated to be 10^19–22 erg cm⁻³, which is consistent with SSs (Zhou et al. 2004).

Starquakes are magnetically powered in a NS, while elastically and gravitationally powered in a SS. In these scenarios, the toroidal oscillation, which might be derived from starquakes (Bastrukov et al. 2007), propagates into the magnetosphere and changes its charge density that brings an induced electric potential (Lin et al. 2015). The electric potential (Ruderman & Sutherland 1975; Chen & Ruderman 1993) is estimated to be

$$\begin{eqnarray}&&{\rm{\Delta }}V\simeq 2.1\times {10}^{12}\,{\left(\displaystyle \frac{{{\rm{\Omega }}}_{\mathrm{osc}}}{10\mathrm{kHz}}\right)}^{1/7}{\left(\displaystyle \frac{B}{{10}^{14}{\rm{G}}}\right)}^{-1/7}{\left(\displaystyle \frac{R}{10\mathrm{km}}\right)}^{4/7}\,V,\end{eqnarray} \tag{ 14 }$$

where the stellar oscillation frequency is estimated to be Ω_osc ∼ c/R ∼ 30 kHz, which enlarges the size of radio beam. Within this picture, a magnetar is most likely to produce electron/positron pair plasma. The electron/positron pair plasma production due to the electric potential is the necessary condition for coherent radio emission. This potential enhances voltage along the gap, which accelerates electrons to higher Lorentz factors, emitting curvature radiation. Then, the pulsar becomes radio loud (i.e., beyond the pulsar death line) until oscillations damp out and the magnetosphere becomes inactive and radio emissions evaporate. Therefore, an FRB may be the "oscillation" of a dead pulsar at near pulsar death line.

We have found that the behaviors of the repeating FRB 121102 are earthquake-like. The distribution of burst energy exhibits a Gutenberg–Richter power-law form which is a well-known earthquake distribution. The distribution of waiting time can be characterized as a Poissonian or Gaussian distribution, which is consistent with earthquakes as well as the local correlated aftershock sequence. The possible origins of the repeater are discussed, including the crustal activity of a magnetar and the solidification-induced stress of a newborn SS. Both possible origins might be associated with SGRs, which are difficult to detect at cosmological distance. Statistic distributions of burst energy and duration time show that FRB 121102 is very similar to SGR 1806–20 (Wang & Yu 2017). In addition, SGR 1806–20 shares some distinctive properties with earthquakes, indicating that SGRs are indeed powered by starquakes (Cheng et al. 1996) and suggests that the giant flares of SGRs are quake-induced (Xu et al. 2006).

In addition, the observed continuous bursts with the modeled occurrence rate λ = 11.0 hr⁻¹ for Poissonian, while τ⁻¹ = 3.2 hr⁻¹ for Gaussian, may originate from some uncorrelated quakes. In that case, it is suggested that these quakes are foreshocks, storing energy and motivating a main quake. Then an aftershock sequence, which may be caused by some local coherent deformations before a new equilibrium sets up, occurs. The motivated shocks are non-Poissonian and not rotation-powered dominated, while the type I starquake may lead to an initial shock that begins when the stresses exceed a certain threshold. Hence, the next repeating FRB might be waiting for ∼10⁶ s because a long time to store elastic energy is needed (Zhou et al. 2014).

The latest FRB volumetric rate, including all of the repeating bursts, is calculated as R_FRB ∼ 10⁻⁵ Mpc⁻³ yr⁻¹ out to redshift of 1 (Law et al. 2017). A pulsar, which has a solid crust or stiff EOS, would naturally have glitches as the result of starquakes. From the statistics of pulsar glitches, the number of glitches per year can be interpreted as (Espinoza et al. 2011),

$$\begin{eqnarray}&&{\dot{N}}_{{\rm{g}}}\simeq 0.003\times {\left(\displaystyle \frac{\dot{\nu }}{{10}^{-15}\mathrm{Hz}{{\rm{s}}}^{-1}}\right)}^{0.47}\,{\mathrm{yr}}^{-1},\end{eqnarray} \tag{ 15 }$$

where $\dot{\nu }$ is the rotational frequency derivative. For a typical millisecond pulsar with P ∼ 10 ms and $\dot{P}\sim {10}^{-21}\,{\rm{s}}\,{{\rm{s}}}^{-1}$, the number of glitches per year can be evaluated to ∼3 × 10⁻⁴ yr⁻¹. A total number of glitches N_g ∼ 3 × 10³ can be estimated during the pulsar lifetime ∼10 Myr. It is hypothesized that the hydrogen-poor superluminous supernovae (SLSNe-I) are powered by millisecond magnetars. The volumetric birth rate of SLSNe-I is R_SLSN = 10⁻⁸ Mpc⁻³ yr⁻¹ (Gal-Yam 2012). Therefore, we estimate the FRB volumetric rate R_FRB ≃ N_gR_SLSN ∼ 3 × 10⁻⁵ Mpc⁻³ yr⁻¹. This inferred event rate from FRB/SLSNe-I-associated events is consistent with the observational FRB events.

Starquakes associating with some X-ray or gamma-ray bursts in a normal NS share similar behaviors with that in a SS, while the X-ray spectra might be different in these scenarios. In an SS atmosphere, thermal X-rays from the lower layer of a normal NS atmosphere are prohibited, and relatively more optical/UV photons and a energy cutoff at X-ray bands are exhibited (Wang et al. 2017). Considering an NS at 1 Gpc with 2–8 keV flux of ∼2 × 10⁻¹⁶ erg cm⁻² s⁻¹, which is consistent with the X-ray limit of Chandra and XMM-Newton, the luminosity is calculated to ∼10³⁹ erg s⁻¹. Such a distant source is too faint to be detected by current X-ray telescopes, except those with a supper-Eddington luminosity. Normal NSs and SSs have different EOSs, which are most likely to be further tested by gravitational wave and electromagnetic radiation from mergers of compact stars (Abbott et al. 2017; Lai et al. 2017).

We expect to detect more repeating events. More constraints on the mysterious origin of FRBs will be given by the statistics' growing samples. The earthquake-like behaviors, including distributions of energy and waiting times for the repeater, are expected to be tested by more continuous monitoring observations.

We are grateful to Stephen Justham at National Astronomical Observatories, Chinese Academy of Sciences (CAS), for discussions. R.X.X. acknowledges the support of National Key R&D Program of China (No. 2017YFA0402602), National Natural Science Foundation of China (NSFC) 11673002 and U1531243. W.Y.W. and X.L.C. acknowledge the support of MoST 2016YFE0100300, NSFC 11473044, 11633004, 11373030, CAS QYZDJ-SSW-SLH017. R.L. and K.J.L. are supported by NSFC U15311243 National Basic Research Program of China, 973 Program, 2015CB857101, XDB23010200, 11690024, 11373011, and funding from the Max-Planck Partner Group.