Extremely Rapid X-Ray Flares of TeV Blazars in the RXTE Era

Blazars, including BL Lac objects and flat-spectrum radio quasars (FSRQs), are a special class of radio-loud active galactic nuclei (AGNs) that have one of their relativistic jets pointing very close to our line of sight (e.g., Urry & Padovani 1995). FSRQs have luminous broad emission lines that are weak or absent in BL Lac objects. Due to Doppler boosting, the emission of a blazar is usually dominated by the jet whose spectral energy distribution (SED) shows two broad humps that smoothly extend from radio to γ-rays. The low-energy hump can extend from radio to soft X-rays. According to the frequency of the first hump, BL Lac objects are further divided into low-frequency peaked BL Lac objects (LBLs; ${\nu }_{\mathrm{peak}}\,\lt {10}^{14}$ Hz), intermediate-frequency peaked BL Lac objects (IBLs; ${10}^{14}\lt {\nu }_{\mathrm{peak}}\lt {10}^{15}$ Hz) and high-frequency peaked BL Lac objects (HBLs; ${\nu }_{\mathrm{peak}}\gt {10}^{15}$ Hz, e.g., Padovani & Giommi 1995; Abdo et al. 2010a). The high-energy hump extends from hard X-rays to γ-rays, even sometimes the very high-energy (VHE) TeV band. Such VHE blazars, typically HBLs, are called TeV blazars.

The low-energy hump is attributed to the synchrotron emission of highly relativistic electrons gyrating in a magnetic field in the jet. The origin of the high-energy hump, however, is still debated. A popular explanation is inverse-Compton emission from the same population of relativistic electrons that produce the synchrotron emission. The seed photons of the inverse-Compton scattering process could be local synchrotron photons (usually for BL Lacs) in the jet and/or external photons from the central engine (usually for FSRQs) or the cosmic microwave background. These models are thus called leptonic models. The hadronic models, on the other hand, attribute γ-ray emission to synchrotron emission of protons (Mücke & Protheroe 2001; Mücke et al. 2003; Fraija & Marinelli 2015) or proton-induced cascades (Mannheim 1998).

The observed photon flux of blazars varies significantly across the electromagnetic spectrum on timescales from minutes to years (e.g., Wagner & Witzel 1995; Ulrich et al. 1997). The origin of the variability is not well understood. Generally, the variability is noise-like (e.g., Kataoka et al. 2001; Chatterjee et al. 2012), similar to the variability of radio-quiet AGNs (e.g., Markowitz et al. 2003). However, blazars are also known to have bursts that show flare-like structures (e.g., Marscher et al. 2010), which may have recognizable patterns (Sasada et al. 2017). The outbursts can be explained by internal shocks of the jets when a new relativistic blob of plasma catches up with an old blob and accelerates particles to ultra-relativistic energies (e.g., Spada et al. 2001). Based on several similarities between the jet emission and corona-disk emission, the ultimate origin of the variability may still be accretion-rate fluctuations of the disk (McHardy 2008).

The shortest variability timescale is a crucial parameter because it serves as an independent constraint on the physical scale of the emission region (Tavecchio et al. 1998), which cannot be easily provided by other observational measurements. Blazars are usually most variable at frequencies just above the two SED humps (e.g., Ulrich et al. 1997; Madejski & Sikora 2016), which usually fall in the hard X-ray and TeV bands in the case of TeV blazars (e.g., Aleksić et al. 2015a; Baloković et al. 2016; Bartoli et al. 2016). In particular, an increasing number of TeV blazars show γ-ray flaring activity on timescales from several to a few tens of minutes that are detected by ground-based Cherenkov telescopes, including both BL Lac objects (e.g., Gaidos et al. 1996; Aharonian et al. 2007; Albert et al. 2007; Arlen et al. 2013) and FSRQs (Aleksić et al. 2011). The minute-scale variability in the TeV band (e.g., Aharonian et al. 2007) has strong implications for our understanding of AGN jets (Begelman et al. 2008). X-ray and TeV emission may be directly related to the same high-energy tail of the relativistic electron population. Indeed, the light curves of HBL-type TeV blazars in the X-ray and TeV bands are usually correlated (e.g., Aleksić et al. 2015b; Furniss et al. 2015b). Attempts to search for extremely rapid X-ray variability have been made (e.g., Cui 2004; Xue & Cui 2005; Paliya et al. 2015; Pryal et al. 2015). The same source can have minute-scale variability in both the X-ray band and TeV band (e.g., Mrk501, Xue & Cui 2005; Albert et al. 2007). However, "orphan" TeV flares that have no X-ray counterparts are occasionally reported (e.g., Krawczynski et al. 2004; Błażejowski et al. 2005; Acciari et al. 2009; Fraija et al. 2015).

Rapid TeV variability has germinated various models to explain the small timescales. Most models involve some very compact regions moving in the rest frame of the jet. These compact regions could be "jets in a jet" that are either produced by magnetic reconnection processes in a Poynting flux-dominated jet (Giannios et al. 2009) or relativistic turbulence in the jet (Narayan & Piran 2012). The minijets-in-a-jet model can consistently produce the statistical properties of blazar flux (Biteau & Giebels 2012). There are also models involving a red giant star being stripped of its envelope by the jet (Barkov et al. 2012) and models involving beams of magneto-centrifugally accelerated electrons occasionally pointing toward us (Ghisellini et al. 2009).

Several well-studied TeV blazars show rich spectral behavior in X-rays, which may represent the general behavior of the synchrotron peak of all AGN jets. The X-ray spectra are usually curved (Massaro et al. 2004) and can only locally be fitted by a power law. The spectral variation with flux can be complex (Zhang et al. 2002; Cui 2004). Generally, the spectrum hardens when the flux increases (e.g., Gliozzi et al. 2006; Xue et al. 2006; Tramacere et al. 2009), but photon indices can saturate at higher fluxes (Xue & Cui 2005; Giebels et al. 2007). The synchrotron peak usually moves to higher frequencies with increasing flux during outbursts (e.g., Pian et al. 1998), but no correlation between the break energy and the flux exists when a broken power law is adopted to fit the X-ray spectra (Xue & Cui 2005; Giebels et al. 2007; Garson et al. 2010). A cooling break in the spectrum of emitting particles cannot explain these features (Wierzcholska & Wagner 2016), and some special particle acceleration processes may be involved (Madejski & Sikora 2016). There are also energy-dependent lags between the variations of different energy bands. In some flares, soft bands lag behind hard bands (e.g., Zhang et al. 2002), while lags in the opposite direction can also happen (e.g., Ravasio et al. 2004; Sato et al. 2008). Hysteresis in the HR (hardness ratio)–flux diagram is often used as a diagnostic of lags. Clockwise loops (e.g., Acciari et al. 2009; Kapanadze et al. 2016) in the HR–flux plane are a sign of soft lags while counterclockwise loops (e.g., Tramacere et al. 2009) are a sign of hard lags. The same source can exhibit both clockwise and counterclockwise loops; the observed patterns are further complicated by the superposition of flares at different timescales (Cui 2004). The above knowledge of TeV blazars in the X-ray regime comes from studies focusing on timescales of hours to weeks. We will extend this kind of analysis to much smaller timescales in this paper.

The main goal of this paper is to search for X-ray flares at sub-hour timescales from TeV blazars in the entire Rossi X-ray Timing Explorer (RXTE) archival database. We use data from the narrow-field pointing instrument Proportional Counter Array (PCA) on board RXTE, covering a nominal energy range of 2–60 keV. The RXTE satellite was launched in 1995 December and ceased science operation in 2012 January. During its lifetime, it accumulated more than ∼16 Ms of exposure time on TeV blazars in hard X-rays, surpassing any other X-ray observatory. We describe data reduction and the searching results in Section 2. Most of our analysis is based on an assumption that the observed photons are from a flaring component and an underlying constant/slowly varying component, possibly from two separated sites. We describe light curve-model fitting and spectral-model fitting in Sections 3 and 4, respectively. We discuss the implications of our findings in Section 5 and summarize them in Section 6. In the following, we use the ΛCDM model, with H₀ = 67.7 km s⁻¹ Mpc⁻¹ and ${{\rm{\Omega }}}_{{\rm{m}}}=0.307$ (Planck Collaboration et al. 2016).

We retrieved all of the archival RXTE/PCA observations⁷ of TeV blazars (see Table 1). The total number of observations is ∼8400, and the total exposure time is ∼16 Ms. We used the Standard 2 mode data, which have a time resolution of 16 s. The data were reduced following the standard procedure.⁸ We created filter files and good time intervals (GTIs) for each observation according to the suggested screening criteria⁹ for faint sources. Background data were then simulated using the appropriate model.¹⁰ We applied the GTIs to both observational data and simulated background and extracted light curves in initial 16 s bins from channels that correspond to ∼2–20 keV. The light curves of net count rates were calculated using lcmath in the HEASoft (v6.19) package. As only PCU2 among the five proportional counter units (PCUs) of PCA is almost always in operation, we extracted light curves from PCU2 for flaring event selection. We visually inspected every light curve to select events in individual observations that contain a complete or nearly complete sub-hour flaring profile. Specifically, we require the flare to have apparent rise and fall; we also require the existence of a plateau either before the rise and/or after the fall to assess the completeness of the flare and the level of the background component.¹¹ Note that RXTE is in low Earth orbit, and thus an uninterrupted light curve is usually less than ∼50 minutes due to the Earth's occultation or passage through the South Atlantic Anomaly, etc., which limits the timescales of the events investigated. After identifying fast flaring events, we extracted light curves and spectra from all the PCUs available during that observation to achieve a higher signal-to-noise ratio (S/N) for further analysis.

From the complete RXTE/PCA database, two new fast flaring events were found. They belong to PKS 2005−489 and S5 0716+714. The light curves of these two events are shown in Figure 1. We also plot an event of Mrk 501 (reported in Xue & Cui 2005). A fast-flaring candidate event of 1ES 1101−232 is reported in Appendix A, which has relatively low credibility because of limited S/N. These flaring observations generally lack simultaneous observations in other wavebands.

Zoom In Zoom Out Reset image size

We checked for potential contamination by soft electron flares that were not screened out by the criteria "ELECTRON2.LE.0.1" in data cleaning. The contemporaneous Electron2s of each event were well below 0.1 and did not show any apparent electron flaring activity that may be responsible for the X-ray flares. The longitudes and latitudes of the satellite at the onsets of the X-ray flares did not cluster in the anomalous high background region (cf. Figure 7 of Xue & Cui 2005). Further support for the genuineness of the flaring events comes from their light curve and spectral features explored below. They behave like well-known X-ray flares of TeV blazars, only at much smaller timescales. In conclusion, we did not find any sign of contamination of soft electrons or any other known sources for any of the flaring events. However, the possibility of an unrelated X-ray transient in the field of view still cannot be ruled out entirely, as the PCA lacks the capacity of imaging.

Light curves in additional energy bands (see the second column of Table 2 and Figures 2–4) are extracted according to the energy-channel conversion table.¹² Thanks to the high data quality, variability at timescales down to the time resolution (16 s) of the light curves (see Figure 2) is seen. Because of the variety of data quality and gain epochs (Jahoda et al. 2006), we do not have uniform definitions for different bands for all three sources. We have four different bands for PKS 2005−489, while we have three bands for Mrk 501 and S5 0716+714; the full band of PKS 2005−489 and Mrk 501 is 1.94–20.30 keV, while the full band of S5 0716+714 is 2.06–10.11 keV.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2.  Fitting Results for the Fast Flaring Events

Band	Energy	${F}_{{\rm{c}}}$	F0	t0	${\tau }_{{\rm{r}}}$	${\tau }_{{\rm{d}}}$	${t}_{{\rm{p}}}$	${F}_{{\rm{p}}}/{F}_{{\rm{c}}}$	ξ	${\chi }_{\nu }^{2}/\mathrm{dof}$
	(keV)	(c/s)	(c/s)	(s)	(s)	(s)	(s)
PKS 2005–489
Soft	1.94–5.47	${2.08}_{-0.02}^{+0.02}$	${3.81}_{-0.26}^{+0.28}$	${303}_{-6}^{+8}$	${23}_{-4}^{+5}$	${160}_{-15}^{+16}$	${343}_{-7}^{+8}$	2.3	0.75	0.957/166
Medium	5.47–10.11	${0.98}_{-0.02}^{+0.02}$	${3.45}_{-0.24}^{+0.26}$	${318}_{-7}^{+10}$	${29}_{-4}^{+8}$	${141}_{-15}^{+14}$	${356}_{-6}^{+9}$	3.2	0.66	1.05/166
Hard	10.11–20.30	${0.33}_{-0.02}^{+0.02}$	${2.43}_{-0.34}^{+0.35}$	${324}_{-14}^{+33}$	${29}_{-8}^{+18}$	${105}_{-30}^{+25}$	${353}_{-12}^{+19}$	5.3	0.56	0.895/166
Full	1.94–20.30	${3.4}_{-0.04}^{+0.04}$	${9.49}_{-0.48}^{+0.47}$	${311}_{-5}^{+6}$	${25}_{-3}^{+4}$	${143}_{-10}^{+11}$	${348}_{-5}^{+6}$	2.8	0.70	1.18/166
S5 0716+714
Soft	2.06–5.31	${0.19}_{-0.08}^{+0.05}$	${2.04}_{-0.32}^{+0.28}$	${1344}_{-100}^{+111}$	${222}_{-43}^{+78}$	${825}_{-171}^{+326}$	${1573}_{-52}^{+66}$	7.5	0.58	1.48/23b
Hard	5.31–10.11	${0.00}_{-0.08}^{+0.04}$	${1.67}_{-0.37}^{+0.21}$	${1216}_{-108}^{+64}$	${161}_{-39}^{+65}$	${1588}_{-285}^{+1903}$	${1551}_{-54}^{+107}$	$\gt 17$ a	0.82	1.62/23b
Full	2.06–10.11	${0.17}_{-0.09}^{+0.08}$	${3.95}_{-0.45}^{+0.38}$	${1258}_{-55}^{+57}$	${178}_{-29}^{+40}$	${1158}_{-179}^{+320}$	${1547}_{-39}^{+48}$	17	0.73	1.15/159b
Mrk 501 first substructure
Soft	1.94–5.82	${5.60}_{-0.04}^{+0.03}$	${3.32}_{-0.88}^{+0.02}$	${1918}_{-59}^{+32}$	${99}_{-34}^{+27}$	${110}_{-41}^{+75}$	${1924}_{-26}^{+14}$	1.30	0.05	0.99/151
Hard	5.82–20.30	${5.35}_{-0.04}^{+0.05}$	${3.29}_{-0.58}^{+0.38}$	${2001}_{-80}^{+8}$	${147}_{-54}^{+18}$	${33}_{-6}^{+65}$	${1961}_{-40}^{+8}$	1.38	−0.64	1.09/151
Full	1.94–20.30	${10.96}_{-0.06}^{+0.06}$	${6.9}_{-1.1}^{+0.2}$	${1961}_{-47}^{+18}$	${119}_{-29}^{+21}$	${67}_{-19}^{+41}$	${1936}_{-18}^{+9}$	1.33	−0.28	1.23/151
Mrk 501 second substructure
Soft	⋯	⋯	${1.78}_{-0.58}^{+0.25}$	${2206}_{-23}^{+90}$	${40}_{-10}^{+135}$	${160}_{-65}^{+71}$	⋯	⋯	0.60	⋯
Hard	⋯	⋯	${1.73}_{-0.34}^{+0.35}$	${2135}_{-16}^{+84}$	${30}_{-11}^{+79}$	${302}_{-115}^{+225}$	⋯	⋯	0.82	⋯
Full	⋯	⋯	${3.7}_{-0.9}^{+0.6}$	${2177}_{-18}^{+54}$	${43}_{-13}^{+55}$	${213}_{-63}^{+83}$	⋯	⋯	0.67	⋯

Notes.

^aNote our fitting cannot constrain the hard-band quiescent flux level of S5 0716+714, so the hard-band amplitude is unmeasurable. We thus set the full-band amplitude as the lower limit of that of the hard band. ^bWe fitted the model to the soft- and hard-band light curves in bins of 96 s, while the full-band light curve is in bins of 16 s.

Download table as:  ASCIITypeset image

3.1. Method of Fitting

We fitted the light curves with a constant flux plus an exponentially rising and decaying flare following Abdo et al. (2010b)¹³ :

$$\begin{eqnarray}&&F(t)={F}_{{\rm{c}}}+{F}_{0}{({e}^{\tfrac{{t}_{0}-t}{{\tau }_{{\rm{r}}}}}+{e}^{\tfrac{t-{t}_{0}}{{\tau }_{{\rm{d}}}}})}^{-1},\end{eqnarray} \tag{ 1 }$$

where ${F}_{{\rm{c}}}$ represents the constant flux level underlying the flare, and ${\tau }_{{\rm{r}}}$ and ${\tau }_{{\rm{d}}}$ are the characteristic rising and decaying timescales¹⁴ of the flare. t₀ indicates the transition from rising to falling, and the count rate actually peaks at

$$\begin{eqnarray}&&{t}_{{\rm{p}}}={t}_{0}+\displaystyle \frac{{\tau }_{{\rm{r}}}{\tau }_{{\rm{d}}}}{{\tau }_{{\rm{r}}}+{\tau }_{{\rm{d}}}}\mathrm{ln}\left(\displaystyle \frac{{\tau }_{{\rm{d}}}}{{\tau }_{{\rm{r}}}}\right),\end{eqnarray} \tag{ 2 }$$

which equals t₀ only when the flare is symmetrical (${\tau }_{{\rm{r}}}={\tau }_{{\rm{d}}}$). Therefore, we define the amplitude of the flare as ${F}_{{\rm{p}}}/{F}_{{\rm{c}}}$, i.e., the count rate at ${t}_{{\rm{p}}}$ (${F}_{{\rm{p}}}$) over the constant level, instead of ${F}_{0}/{F}_{{\rm{c}}}$. The symmetry of a flare is described by

$$\begin{eqnarray}&&\xi =\displaystyle \frac{{\tau }_{{\rm{d}}}-{\tau }_{{\rm{r}}}}{{\tau }_{{\rm{d}}}+{\tau }_{{\rm{r}}}},\end{eqnarray} \tag{ 3 }$$

whose value is in the range of [−1, 1]. $\xi =-1$ ($=1$) represents completely right (left) asymmetric profiles with a zero falling (rising) timescale; $\xi =0$ indicates a symmetric flare.

The uncertainties of ${t}_{{\rm{p}}}$ and ${F}_{{\rm{p}}}$ have to be propagated from the errors of other parameters. We adopted an MCMC (Markov Chain Monte Carlo) algorithm to fit the light curves, which returns reliable probability intervals of timescales and amplitudes by sampling from their posterior distributions. We first performed weighted least-squares fitting using a numerical minimizer to obtain the best estimates of ${F}_{{\rm{c}}}$, F₀, t₀, ${\tau }_{{\rm{r}}}$, and ${\tau }_{{\rm{d}}}$. Starting from these initial values, we took 1000 random walk steps in parameter space. The samples of ${t}_{{\rm{p}}}$ and ${F}_{{\rm{p}}}$ were calculated according to Equations (1) and (2). Note that we added a second flaring component to the model in the fitting of Mrk 501 (see the bottom panel of Figure 1).

3.2. Correcting Light Curve Error Bars

3.3. Light Curve Fitting Results

3.4. Spectral Variability of PKS 2005−489

The flare of PKS 2005−489 has the shortest rising timescales and its data are of the highest quality, so we investigated its spectral variability in detail. We plot the hardness ratios¹⁶ (HRs) in the bottom two panels of Figure 2, which show hardening that corresponds to the flare. We also plot an HR–flux diagram of this event in panel (a) of Figure 5, which shows a "harder when brighter" trend and hysteresis. The loop begins with clockwise motion and then follows a counterclockwise direction. Below, we argue that the clockwise trend at the beginning is due to the superposition of two spectral components.

Zoom In Zoom Out Reset image size

The apparent two-component nature of the light curves suggests that the spectral variability is partially caused by a change of the relative fraction of the two components. In other words, if the flaring component has a different hardness ratio from that of the constant component, then even if neither of the spectra changes over time, the observed overall hardness ratio will still change due to the flux variation of the flaring component (see Sun et al. 2014; Ramolla et al. 2015). However, the hysteresis loop in panel (a) suggests the spectrum of the flaring component is intrinsically variable; otherwise, the track in the HR–flux plane while the flux is rising will be identical to the track while the flux is declining, instead of forming a loop. To discriminate the effects caused by the mixing of different components and the flare's intrinsic spectral variability, we further subtract the constant component ${F}_{{\rm{c}}}$ from the light curve of each band, where ${F}_{{\rm{c}}}$ is from Table 2. We calculated the hardness ratios from the resulting flare-only light curves in panels (c), (d), and (e) and an HR–flux diagram in panel (b). The plots suggest that the flare emerges with a hard spectrum, softens gradually as the flux rises, then hardens near the highest flux, and finally softens as the flare fades away. Only a counterclockwise loop is apparent in the HR–flux plane for the flare-only light curve. Based on the evidence above, we know that the clockwise loop in panel (a) could be caused by the sudden emergence of a hard flaring component and the spectral variability soon follows the spectral variability of this flaring component due to its increasing dominance. The flaring component itself has complex spectral variability and the counterclockwise loop in panel (b) suggests the existence of a hard time lag, albeit being subject to large uncertainties. Note that the overall oblique "8" shape before subtracting the constant component in panel (a) is reminiscent of a flare of Mrk 421 that lasted ∼60 ks and was reported by Garson et al. (2010) using Suzaku data. We reanalyzed the Suzaku data of Mrk 421 and performed the same spectral variability analysis as for PKS 2005−489 above in Appendix B. We reproduced the same results as for PKS 2005−489 in the flare of Mrk 421 and thereby strengthened our conclusions. The similar behavior of spectral variability, seen both in the extremely rapid flare of PKS 2005-489 and in the long-duration flare of Mrk 421, is reminiscent of the scale invariant nature of X-ray flares from TeV blazars (Cui 2004; Xue & Cui 2005).

We also calculated the CCF (Edelson & Krolik 1988; Welsh 1999) of the soft-band and hard-band flare-only light curves but did not find an obvious time lag. Any time lag in this fast flaring event of PKS 2005−489 could be intrinsically small compared with the time resolution of the observation (16 s). The bandpass of our data is narrow, spanning about one order of magnitude, so any energy-dependent lags may not be significant. The size of lags may also be positively correlated to the duration of the flares (Zhang et al. 2002), which means that flares of shorter duration have smaller lags. Note that in Table 2, the differences of t₀, ${t}_{{\rm{p}}}$, and ${\tau }_{{\rm{r}}}$ between different bands are small compared with relatively large error bars.

Motivated by the two-component model in light curve fitting, we carried out spectral fitting in a similar way. By dividing the observation into flaring and quiescent phases according to the light curves, we separately extracted and jointly fitted the spectra of the two phases. The spectra were extracted from the top layers of the operating PCUs. The spectra and models are shown in Figure 6 in the EF_E representation and the photon indices are tabulated in Table 3. We used channels that correspond to 3–20 keV for PKS 2005−489 and Mrk 501 and 3–10 keV for S5 0716+714, respectively. In addition to Galactic absorption, we used a power law and sum of two power laws in XSPEC (v12.9.0) to fit the spectra of the quiescent and flaring phases, respectively. We tied the power law in the quiescent phase to one of two power laws in the flaring phase, for both normalization and photon index (see the legends of Figure 6). We also calculated the X-ray luminosity of the constant component and the total energy of the flaring component in the X-ray band (namely average flare luminosity times the length of the defined flare phase) in Table 3. Note that the spectra should suffer from overestimation of error bars due to the same reason discussed in Section 3.2, but we are not able to correct the error estimation of the spectra as we did in the light curve fitting. One consequence of the overestimation of error bars is that the parameters in the spectral fitting have large confidence intervals. Another consequence is that we might be able to fit the spectra with many models. For example, we can fit the spectra of the flaring phase, quiescent phase, and whole observation with simple power laws, and the Chi-square statistics can still be acceptable. We think that describing the spectra as a combination of different power laws is more physically appropriate and is consistent with the light curve fitting.

Zoom In Zoom Out Reset image size

Table 3.  Spectral Fitting Results of the Fast Flaring Events

Source	${{\rm{\Gamma }}}_{1}$	${{\rm{\Gamma }}}_{2}$	${F}_{{\rm{quiescent}}}$ a	${F}_{{\rm{flare}}}$ a	${F}_{{\rm{flare}}{\rm{only}}}$ a	${\chi }_{\nu }^{2}/\mathrm{dof}$	${L}_{{\rm{X}}}$ b	Ec
PKS 2005−489	2.73 ± 0.05	1.71 ± 0.18	29.5	42.4	12.9	0.88/38	$3.86\times {10}^{44}$	$2.08\times {10}^{47}$
Mrk 501	2.01 ± 0.02	1.86 ± 0.20	115.7	130.0	14.3	0.76/38	$3.23\times {10}^{44}$	$5.13\times {10}^{46}$
S5 0716+714	2.48 ± 0.27	1.89 ± 0.15	5.7 (pre)	12.3	12.2	0.32/16	$1.90\times {10}^{45}$	$1.06\times {10}^{49}$
			3.8 (post)				$1.25\times {10}^{45}$

Notes.

^aFluxes are in units of 10⁻¹² erg cm⁻² s⁻¹; the energy range is 3–20 keV for PKS 2005−489 and Mrk 501, but 3–10 keV for S5 0716+714. ^bX-ray luminosity of the constant component luminosity in units of erg s⁻¹. ^cEnergy of the flaring component in units of erg.

Download table as:  ASCIITypeset image

Note that due to the limited length of the observation, we have neither a complete flaring nor a quiescent phase of S5 0716+714. Therefore, we extracted the spectra of the two observations just before and after the flare, which are referred to as pre-flare and post-flare, respectively. The gaps between the observations are 3–4 days. The three observations are plotted in the bottom-right panel of Figure 6, where we fitted each of the pre-flare and post-flare light curves with a constant flux. The spectra are jointly fitted with a tied "constant" component photon index, but the normalizations of the "constant" component are free parameters in this case. In other words, the constant component varies on timescales of a few days. If we leave the normalization of the constant component underlying the flare (factor "A" in the legend of top-right panel) free, the value of "A" in the fitting results is negligible but its error bars are large. Additionally, we are not able to constrain the photon index of the flaring component (${{\rm{\Gamma }}}_{2}$) very well in this case. Thus, we fixed $A=0.6{\rm{D}}$, according to the light curve fitting results in the bottom-right panel of Figure 6. The choice of 0.6 or another reasonable value does not affect the result of ${{\rm{\Gamma }}}_{2}\lt {{\rm{\Gamma }}}_{1}$.

The above fitting shows that for all three flares, the flaring components have harder photon indices than the corresponding constant components (see Table 3). The difference is not as apparent for Mrk 501 because of the relatively small variation amplitude (i.e., relatively weak flaring component). The photon index of the constant component ranges from 2.0 to 2.7, indicating that we are observing different declining parts of the synchrotron hump for different sources. The photon indices of the flaring components, on the other hand, lie in a narrower range around 1.8. Every estimated average flux (averaged over the length of the flare phase we defined) of the events is $\sim {10}^{-11}$ erg cm⁻² s⁻¹ (see Table 3). Using the "efficiency limit" of compact sources (e.g., Fabian 1979; Brandt et al. 1999), ${\rm{\Delta }}L/{\rm{\Delta }}t\lesssim 2\times {10}^{41}{\eta }_{0.1}$ erg s⁻¹, the matter to energy conversion efficiency η is greater than 1 for PKS 2005−489 and S5 0716+714, suggesting the existence of significant boosting of the emission (e.g., Remillard et al. 1991). Due to the small flare amplitude, the efficiency of Mrk 501 does not exceed 1.

5.1. The Flaring and Constant Components

5.2. X-Ray Radiation Process of S5 0716+714

5.3. Constraints on the Emission Region

The upper limit on the physical scale of the flaring region can be given by

$$\begin{eqnarray}&&R\approx \displaystyle \frac{{{ct}}_{\mathrm{flare}}\delta }{(1+z)}\approx {10}^{14}\left(\displaystyle \frac{\delta }{30}\right)\left(\displaystyle \frac{{t}_{\mathrm{flare}}}{100\ {\rm{s}}}\right)\,\mathrm{cm},\end{eqnarray} \tag{ 4 }$$

where δ is the Doppler boosting factor, and ${t}_{\mathrm{flare}}$ is the observed variability timescale. Because the flares are asymmetric, the size of the emitting region is reflected by the rising timescale ${\tau }_{{\rm{r}}}$ (Zhang et al. 2002). Searching the literature, we find black hole mass estimates for three out of the four sources (except 1ES 1101−232). We plot the rising timescale in the galaxy's frame against black hole mass in Figure 7. Note that the typical Doppler factor of BL Lac objects is 10–20 (20–30 for FSRQs, e.g., Hovatta et al. 2009). The size of the black hole is often thought to be a natural lower limit on the physical scale of the emission region. We also plot this lower limit as a function of black hole mass assuming different Doppler factors. The variability timescales should lie above the lines, if the size of the black hole is a hard lower limit. The rising timescale of PKS 2005−489 is, as far as we know, shorter than any variability timescale of AGNs at any wavelength ever reported (e.g., Remillard et al. 1991; Yaqoob et al. 1997; Gallo et al. 2004; Xue & Cui 2005; Aharonian et al. 2007; Albert et al. 2007; Aleksić et al. 2014; Kara et al. 2016). Figure 7 shows that the Doppler factor of PKS 2005−489 has to be larger than several hundred, which appears unrealistic, in order to support the idea that the black hole sets a lower limit on the physical size of the flaring region. So far, there is no correlation found between the observed minimum variability timescale and the black hole mass (Wagner 2008; Vovk & Babić 2015). The events of PKS 2005−489 and Mrk 421 (Appendix B) display remarkably similar spectral evolution, supporting the same process driving the spectral variability. However, the timescales of the two events are different by two orders of magnitude, reflecting an intrinsic difference between the two jets. Noticeably, the famous Mrk 421 has many more X-ray observations than PKS 2005−489; however, the variability timescale of Mrk 421 has never been found to be as short as the PKS 2005−489 event (Cui 2004; Paliya et al. 2015; Pryal et al. 2015). A lower limit on the variability timescale probably does exist, and probably is not set by the central supermassive black hole.

Zoom In Zoom Out Reset image size

The synchrotron cooling time of emitting electrons is given by ${t}_{\mathrm{cool}}\approx 6\pi {m}_{e}c/{\sigma }_{{\rm{T}}}\gamma {B}^{2}$ (Rybicki & Lightman 1979). The observed photon energy at the synchrotron peak is given by ${E}_{{\rm{p}}}=\delta h\nu \equiv (3{eh}/4\pi {m}_{e}c)\delta {\gamma }^{2}B$ (Rybicki & Lightman 1979). Combining the above two equations, we have

$$\begin{eqnarray}&&{t}_{\mathrm{cool}}=3.04\times {10}^{3}{B}^{-3/2}{\delta }^{-1/2}{E}_{{\rm{p}}}^{-1/2}\ {\rm{s}},\end{eqnarray} \tag{ 5 }$$

where ${E}_{{\rm{p}}}$ is in units of observed keV (Zhang et al. 2002). If we take PKS 2005−489 as an example, and adopt ${t}_{\mathrm{cool}}={\tau }_{{\rm{d}}}\sim 143\,{\rm{s}}$, ${E}_{{\rm{p}}}=10\,\mathrm{keV}$ and $\delta =30$, we can have $B\approx 1.2$ G. The estimated magnetic field of Mrk 501 is similar to that of PKS 2005−489, but the magnetic field of S5 0716+714 is four times weaker. Radio-loud AGNs are potential accelerators of cosmic rays. We estimated the maximum energy of protons, if they can be accelerated in the same region as the electrons, as ${E}_{\max }={eBR}\sim (2\mbox{--}9)\times {10}^{16}$ eV.

5.4. The Rarity of Extremely Rapid Flares of Blazars

5.5. Particle Acceleration

We searched the entire RXTE archival database for rapid X-ray flares of TeV blazars that last less than one hour. We investigated the temporal and spectral properties of the fast flares discovered under a two-component assumption. Our analysis also includes an X-ray flare of Mrk 421 using Suzaku data that has similar spectral variability to that of PKS 2005−489 in Appendix B. Our main findings are as follows:

1. We discovered two new fast X-ray flares from PKS 2005−489 and S5 0716+714 and a candidate flare from 1ES 1101−232. The event of PKS 2005−489 shows, as far as we know, the most rapid variation of AGNs that has been observed at any wavelength. The extremely small timescale (${\tau }_{{\rm{r}}}\lt 30$ s) defies the size that corresponds to the light-crossing time of the supermassive black hole as a lower limit on the size of the flaring region.

2. The flares are usually superimposed on a constant/slowly varying component. The flaring component generally has a harder X-ray spectrum (${\rm{\Gamma }}=1.7\mbox{--}1.9$) than the constant component (${\rm{\Gamma }}=2.0\mbox{--}2.7$).

3. The higher data quality of the X-ray observations can provide more detail than γ-ray observations. The flare-only component shows a counterclockwise pattern in the HR–flux diagram, providing a sign of likely hard lags. This component also has the hardest spectrum right at its appearance.

RXTE ceased science operation in 2012 January, but X-ray observatories like Chandra, XMM-Newton, Swift, and NuSTAR are still accumulating exposures on TeV blazars. For example, XMM-Newton has 9 Ms of exposure on TeV blazars as of 2017 July. A natural follow-up work is to search for additional rapid X-ray flares from TeV blazars in these archival databases. With the increased sample, we may study more details of the spectral variability of the flares, thus shedding light on their origin and particle acceleration processes. We can also study the correlations between flares and the properties of the jets and central engine, as well as the reason for the small flaring timescales.

We thank the anonymous referee for his/her helpful comments that improved the paper. S.F.Z., Y.Q.X., and Y.J.W. acknowledge support from the 973 Program (2015CB857004), NSFC-11473026, NSFC-11421303, the CAS Frontier Science Key Research Program (QYZDJ-SSW-SLH006), and the Fundamental Research Funds for the Central Universities. S.F.Z. and W.N.B. acknowledge support from Chandra X-ray Center grant G04-15093X.

We found one fast X-ray flaring event of 1ES 1101−232, which is an HBL at a relatively high redshift (z = 0.186) compared with other HBLs. The full-band (2–20 keV) light curve in 16 s bins (Figure 8) shows elevated flux in the second half of the observation. When shown in 96 s bins, the light curve clearly manifests an almost complete flare profile and can be fitted using the model of Equation (1). The amplitude (${F}_{{\rm{p}}}/{F}_{{\rm{c}}}=1.93$) of this flare is modest, and the rising timescale (${\tau }_{{\rm{r}}}=60$s) is almost as fast as PKS 2005−489. We only reported it as a candidate because the flux preceding the flare seems to be in continuous declining. The noisy data prevent us from performing data analysis in the same fashion as other sources (sub-band light curve fitting, joint spectral fitting, etc.).

Zoom In Zoom Out Reset image size

The flare-only component of PKS 2005−489 shows an interesting spectral variability pattern, but the noisy data prevent us from making strong conclusions. The oblique "8" pattern in Figure 5(a) is reminiscent of a flare of Mrk 421, which allows us to confirm our finding in another source using data from a different satellite.

We downloaded Suzaku data (ObsID 703043010) of Mrk 421 (Garson et al. 2010). We only used data from the X-ray Imaging Spectrometer (XIS) on board Suzaku. The data were reprocessed and screened using Aepipeline included in the Suzaku FTOOLS. In addition to the standard screening criteria, we also required the cut-off rigidity to be larger than 6 ${\rm{GV}}/c$ following Garson et al. (2010). We extracted light curves separately from the cleaned event files of XIS0 and XIS3 using xselect in initial 16 s bins. The source region has an inner radius of 35 pixels and outer radius of 408 pixels, while the background region is an annulus with an inner radius of 432 pixels and outer radius of 464 pixels. The light curves were extracted in two energy bands: 0.5–2 keV and 2–10 keV, which are referred to as the soft band and hard band, respectively. Light curves of the same energy band but different detectors (XIS0 and XIS3) were then merged. We rebinned the light curves in a size of 5760 s in accordance with the orbital period of Suzaku. Finally, the background light curves were subtracted to obtain the estimation of net count rates from the source.

We clipped the light curves and kept only the segment of Flare 2 that is defined in Garson et al. (2010). We fitted the soft- and hard-band flares using a model that is analogous to Equation (1). Only the constant component underlying the flare was replaced by a slowly varying component using a linear function ${F}_{{\rm{c}}}=m+\mathrm{Slope}\times t$. This component is in long-term decline in both bands. The fitted models were shown in the top panel of Figure 9. The error bars of the light curves are exceedingly small, which renders the fitting statistically unacceptable and suggests ultimate incorrectness of the model due to the existence of sub-structures. However, variation of the smooth model curves in the figure match reasonably well with observational data.

Zoom In Zoom Out Reset image size

The timescales are ${\tau }_{{\rm{r}}}=9630\,{\rm{s}}$ and ${\tau }_{{\rm{d}}}=9261\,{\rm{s}}$ in the soft band and ${\tau }_{{\rm{r}}}=14756\,{\rm{s}}$ and ${\tau }_{{\rm{d}}}=6796\,{\rm{s}}$ in the hard band, which are consistent with the pattern that hard band rises slowly and decays fast. The fitted long-term declining background was subtracted from the observed light curve of each band to leave out the flare-only component. We plotted hardness ratio variation with time and flux in the bottom panel of Figure 9 and the right panel of Figure 10. Figure 10 shows the same patterns as panels (a) and (b) of Figure 5. The effective photon indices (${{\rm{\Gamma }}}_{{\rm{Eff}}}\sim 1.84$) are also close to the photon indices of the flaring components that are listed in Table 3.

Zoom In Zoom Out Reset image size