The white-light superflares from cool stars in GWAC triggers

GWAC has 36 cameras and can cover about $\sim 5000$ deg² of sky (Xin et al., 2021; Wang et al., 2020). There is a 4K $\times$ 4K CCD for each camera and the pixel size is about $11^{\prime\prime}.7$. The cadence of GWAC is 15 s: 10 s for exposure and 5 s for readout. There is no filter on each camera with a limit on magnitude about $G\sim 15$ mag, but in bad weather the limit magnitude would be much brighter. For a transit that suddenly apparent on a CCD, the GWAC software (Han et al., 2021) triggers one of the two Guangxi-NAOC 60 cm optical telescopes (F60A and F60B) to follow up immediately to check the transit.

We obtained 163 big flares from 162 stars in the GWAC trigger database and all were confirmed by an F60. The flare amplitude $\bigtriangleup G$ ranges from 0.83 to $\sim 10$ mag. Some GWAC images were not saved due to the overload of the master computer and some unknown bugs in the communication program, which is responsible for communications between more than 100 control and calculation computers for GWAC cameras. As a result, we obtained images for 147 GWAC flares and made 147 movies for them. These GWAC light curve data, their plots, and movies are available at https://nadc.china-vo.org/res/r101350/.

In Figure 1, the flare GWAC181208B from Star #18 and the flare GWAC190101A from Star #39 are shown, and the entire flare process can been seen in their movies. The upper two panels are for GWAC181208B, and the bottom two panels are for GWAC190101A. The left two panels show their pre-flare statuses and the right two panels show flare peaks. Star #18 and Star #39 cannot be seen in the upper left and bottom left panels, respectively, in their preflare status, while they appeared at peaks in the upper right and bottom right panels, respectively. The flare GWAC181208B has a common flare profile: a rapid rise (lasting only about 1 minute) followed by a slow decay, with a flare amplitude of $\sim 3.3$ mag, but for the GWAC190101A, its impulsive phase is much slower, which lasted for about 10 minutes.

We searched TESS light curves using the TESS-point Web Tool and found that except Star #147 (GWAC191226A), all other stars have been observed by TESS. We searched TESS and K2 light curves by the Python package, Lightkurve (Lightkurve Collaboration et al., 2018), and found available light curves for 109 stars. For the remaining 52 stars, we extracted light curves from their TargetPixelFiles. For all stars but Star #147 (GWAC191226A), we inspected their Full Frame Images (FFIs) of TESS in MAST (STScI, 2022b) and also Aladin (Bonnarel et al., 2000) by eye to ensure all light curves are not contaminated by nearby bright stars. For all TESS light curves of 161 stars, the bad parts of light curves were removed by hand. The K2 light curves were also obtained if available. Finally, we obtained good light curves from 276 Sectors for 124 stars. From these light curves, we tried to obtain their periods and flares using the algorithm from Li et al. (2023a). In short, the algorithm tries to fit the light curve by a B-spline by iteratively removing flares, then obtains the period by the Lomb-Scargle method (Scargle, 1982; Lomb, 1976) using LombScargle in Astropy (Astropy Collaboration et al., 2013, 2018). After minus the fitted B-spline from the original light curve, flares with 3 successive points higher than 5$\sigma$ were detected. Finally, 1478 flares from 117 stars and periods for 105 stars were obtained. An example can be seen in Figure 1 in Li et al. (2023a). All TESS and K2 light curves used in this work, with flares and periods detected from these light curves are available at https://nadc.china-vo.org/res/r101350/.

We calculated the equivalent duration (ED; Gershberg, 1972) of each flare:

. Here, $t_{0}$ and $t_{1}$ are the flare start and end times in second respectively, $f(t)$ is the flare flux in erg s^-1 at the time $t$, and $f_{0}$ is the stellar quiescent flux in erg s^-1. Then the flare energy $E=ED\times f_{0}$.

GWAC has no filter and we used Gaia $G$ to calibrate its photometry. The accuracy is better than 0.1 mag. We used the method in Li et al. (2023a) to calculate the quiescent flux of a star, the zero point flux of $G$ and $T$ passbands are from Li et al. (2023a) and Sullivan et al. (2015), respectively. The parallaxes of stars are from Gaia DR3 (Gaia Collaboration, 2022). Star #9 (GWAC171207A), #95 (GWAC210112A), #124 (GWAC210217A), and #104 (GWAC210202A), have no available parallaxes in Gaia DR3, so their distances are estimated by $M_{Ks}=1.844+1.116(V-J)$ (Raetz et al., 2020), where $Ks$ and $J$ are from 2MASS (Skrutskie et al., 2006) and $V$ is from TIC 8.2 (Paegert et al., 2022). For Star #11 (GWAC181227A), #44 (GWAC190204A), #49 (GWAC211203A), #86 (GWAC210117A), #117 (GWAC190116B) and #120 (GWAC190117A), their flare energies were calculated from their K2 light curves, so the method in Shibayama et al. (2013) was used, with stellar radii and surface temperatures from Huber et al. (2016). Finally, a blackbody with $T=9000$ K was used to estimate the bolometric energy that a flare released. These flares are listed in Table 2.

We searched low-resolution spectra ($R\sim 1800$) in LAMOST DR10 ¹¹1http://www.lamost.org/dr10/ (Deng et al., 2012; Cui et al., 2012) and SDSS DR17 ²²2https://www.sdss4.org/dr17/, and found 68 low-resolution spectra for 46 stars and 17 spectra for 12 stars, respectively. We also obtained 69 spectra for 67 stars with a resolution of 2.34 Å  pixel^-1 with the instrument G5 on the 2.16 m telescope (Zhao et al., 2018). There are five stars observed by LAMOST and also by the 2.16 m Telescope, so 154 spectra for 120 spectra were obtained. Because the resolution of spectra obtained with the 2.16 m Telescope is too low, there is no available radial velocity can be calculated from them.

Except for LAMOST spectra, all other spectra are flux calibrated. Then the photometries of Pan-STARRS1 $g$ and $r$ (Chambers et al., 2016) were, respectively, used to calibrate LAMOST spectra for stars with $r<15$ mag and $r>15$ mag. The zero point fluxes and the filter transmission curves of Pan-STARRS1 $g$ and $r$ are obtained from the Filter Profile Service³³3http://svo2.cab.inta-csic.es/theory/fps/ of the Spanish Virtual Observatory. All spectra cover H$\alpha$, and the equation

is used to fit the H$\alpha$ emission, where $\lambda$ is the wavelength in Å, $F(\lambda)$ is the spectral flux, and $A_{i},i=0,1,2,3,4,5$ is the coefficients to be fitted. Then the apparent H$\alpha$ luminosity $l_{{\rm H}\alpha}=\int_{w_{1}}^{w_{2}}A_{0}\exp(-|(\lambda-A_{1})/A_{2}|^{A_{3}}/A_{3})d\lambda$, where $w_{1}=6544.61$ Å  and $w_{2}=6584.61$ Å. Figure 2 shows one LAMOST spectrum of TIC 336800108 (Star #12; GWAC181129A), with the fitted H$\alpha$ emission shown by the red curve in the insert panel. The function consists of an exponential function $A_{0}\exp(-|(\lambda-A_{1})/A_{2}|^{A_{3}}/A_{3})$ and a linear function $A_{4}\lambda+A_{5}$. The linear function is used to fit the spectral continuum and the exponential function is used to fit the H$\alpha$ emission. The LAMOST instrumental profile is not perfectly gaussian, and this exponential function can fit LAMOST spectral line profiles very well (see Figure 2). In fact, if $A_{3}$ is set to 2, then the exponential function would be the gaussian function.

To obtain the apparent bolometric flux of a star $m_{\rm bol}$, we firstly calculated the bolometric correction of $BC_{J}$ using the equation in Cifuentes et al. (2020):

, where $J$ is from 2MASS and $G$ from Gaia DR3. Then $m_{\rm bol}=J+BC_{J}$. The apparent luminosity $l_{\rm bol}$ in erg cm^-2 s^-1 of a star can be obtained by the definition of apparent bolometric magnitude in Mamajek et al. (2015) by the equation

That is,

Then the fraction of H$\alpha$ luminosity to bolometric luminosity $L_{\rm H\alpha}/L_{\rm bol}=l_{\rm H\alpha}/l_{\rm bol}$. Thus, $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})=\log_{10}(l_{\rm H\alpha})-\log_{10}(l_{\rm bol})$. $l_{\rm H\alpha}$, $l_{\rm bol}$, $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ and radial velocities obtained from spectra with errors are given in Table 3.

For stars with available spectra, there are 3 stars: Star #94 (GWAC180218A), #107 (GWAC200321A), and #125 (GWAC200317A), have no available H$\alpha$ emission data because of low signal-to-noise fluxes around H$\alpha$ emissions.

The $\bigtriangleup G$ vs. $G_{\rm BP}-G_{\rm RP}$ and $\log_{10}(ED)$ vs. $G_{\rm BP}-G_{\rm RP}$ diagrams are shown in the upper and bottom panels, respectively, in Figure 3. Star #19 (GWAC181229A) has the biggest $\bigtriangleup G\sim 10$ mag, but it is too faint to have available Gaia photometry and had been well studied in Xin et al. (2021), so it is not shown in this figure. From the upper two panels one can see that most $\triangle G$ are between 1 and 2.5 mag, and $\bigtriangleup G$ seems to be independent of $G_{\rm BP}-G_{\rm RP}$. From two lower panels, $\log_{10}(ED)$ also seems to be independent of $G_{\rm BP}-G_{\rm RP}$, or the stellar surface effective temperature.

The bolometric energies of GWAC flares are shown by blue circles in Figure 4. Stars that have available stellar surface effective temperatures in TIC 8.2 (Paegert et al., 2022) are shown in the right panel. From Figure 4 we can see that the flare energy ranges from $10^{32.2}$ to $10^{36.4}$ erg, and decreases with $G_{\rm BP}-G_{\rm RP}$ (increases with the stellar effective temperature).

For comparison, the maximum flare energy of each star in Yang & Liu (2019) from Kepler DR25 is shown in the right panel in gray in Figure 4, from which we can see that the flare energies of GWAC triggers are higher than the maximum ones recorded by Kepler. Yang & Liu (2019) estimated flare energies from the Kepler band, while we estimated flare energies from the $G$ band. However, both assumed that bolometric flare energies are from a blackbody with a temperature of 9000 K. There are 402 Kepler stars with effective temperatures in range of 2700-3700 K for comparison, and each star was monitored by Kepler for 4 yr, which implies that GWAC superflares are at most once every 4 yr in Kepler cool stars. We also compared GWAC superflares with those detected by NGTS shown in Figure 4 in Jackman et al. (2023), and found that GWAC superflares have similar energies as those in Jackman et al. (2023), and significantly higher than those in Jackman et al. (2021), which implies that these superflares may be from the same category of top energetic flares.

The spectral types were assigned by comparing the standard spectra in Kirkpatrick et al. (1991). Using the parallaxes and photometries from Gaia DR3 (Gaia Collaboration, 2022), the H-R diagram is given in Figure 5.

Both Star #40 (GWAC201218A) and #60 (GWAC191030A) have a spectral type of M3.5, but their $G_{\rm BP}-G_{\rm RP}$ are in the range of M2 as shown in Figure 5. Their spectra were obtained with the 2.16m Telescope and are plotted in Figure 6 with the comparing spectrum of the normal M3.5V Star #152 (GWAC180212A), which was also obtained with the 2.16m Telescope. From their spectra we can see that spectra of these 3 stars are similar, but the Na I doublet line around 8190 Å  of Star #40 (GWAC201218A) is significantly weaker than those of Star #60 (GWAC191030A) and #152 (GWAC180212A), which means that Star #60 (GWAC191030A) is a dwarf, while Star #40 (GWAC201218A) is a giant or young star. We checked Star #40 (GWAC201218A) and noticed that its $ruwe=6.958$ in Gaia DR3⁴⁴4https://gea.esac.esa.int/archive/documentation/GDR2/Gaia_archive/chap_datamodel/sec_dm_main_tables/ssec_dm_ruwe.html, which implies its data may be unreliable. As for Star #60 (GWAC191030A), it is still unclear why its $G_{\rm BP}-G_{\rm RP}$ is significantly bluer than other M3.5V stars in Figure 5.

Star #59 (GWAC210109A) is an M5.5 young star in the Taurus star-forming region (Esplin & Luhman, 2019). Star #155 (GWAC211016A) is contaminated by a nearby star with a distance of $<1.^{\prime\prime}2$, so its Gaia data may be unreliable.

Star #143 (GWAC210518A) was reported to be a white dwarf + main sequence (WDMS) binary (Rebassa-Mansergas et al., 2012). To find more WDMS stars, for Star #34 (GWAC190901A), #43 (GWAC190121B) and #75 (GWAC190324A), #76 (GWAC201113A), #98 (GWAC180117A), #102 (GWAC190206A), #106 (GWAC201216A), #117 (GWAC190116B), #133 (GWAC200514A) and #135 (GWAC200319A) with abnormally higher FUV and NUV, we fitted their SEDs of GALEX FUV, NUV, PS1 g, r, i, z, y, 2MASS J, H, K and WISE (Wright et al., 2010) W1, W2 bands by a white dwarf and a main sequence star templates (Yuan et al., 2023). Finally, only Star #43 (GWAC190121B) and #76 (GWAC201113A) seem to be well fitted as shown in Figure 7. Certainly, UV excesses of these stars may be from their active chromospheres.

For stars that have no available spectra, their spectral types were assigned by their $G_{\rm BP}-G_{\rm RP}$ in Figure 5. From Figure 5, we can see that spectral types of GWAC flare stars range from M2 to L1, and most stars are M4, which implies that these stars around the convective boundary tend to produce big flares.

The PARSEC isochrones (Chen et al., 2014) of 0.01 Gyr, 0.1 Gyr and 1 Gyr with [Fe/H] $=$ 0 are overplotted in Figure 5, from which we can see that most stars should be younger than 1 Gyr.

There are 16 GWAC flare stars with more than 20 flares in their TESS or K2 light curves, and their cumulative flare frequency distributions (FFDs) are shown in Figure 8, where the cumulative flare frequencies of each star are shown in gray circles with the bolometric flare energies from a blackbody of $T=9000$ K. Each FFD is fitted by the linear function

and shown by a black line. Here, $\nu$ is the cumulative flare frequency in day^-1, and $\alpha$ and $\beta$ are the parameters to be fitted. These 16 stars with their $\alpha$ and $\beta$ and the predicted frequencies of GWAC flares by Equation 5 are shown in Table 4. The GWAC flare energies on the fitted FFD lines are shown by red pentagons in Figure 8. We found that most (13/16) GWAC flares can happen more than once every year, and only three GWAC flares happen once every several or even 20 yr.

To explore kinematic properties of these flare stars, we calculated their tangential velocities ($V_{\rm T}$) for stars with available proper motions and parallaxes in Gaia DR3, and for stars that are not proven to be binaries in this work and have reliable radial velocities in SDSS and LAMOST data, their velocities in the local standard of rest (LSR) ($U_{\rm LSR},V_{\rm LSR},W_{\rm LSR}$) were also calculated by the python package astropy.coordinates.SkyCoord with the velocity of the Sun relative to the LSR of $(U_{\rm LSR},V_{\rm LSR},W_{\rm LSR})_{\odot}=(11.1,12.24,7.25)$ km s^-1 (Schönrich et al., 2010). We notice that there are some stars have several available radial velocities, so the medium velocity was used during calculation. Here, $U_{\rm LSR}$ is directed to the Galactic center, $V_{\rm LSR}$ is in the Galactic rotational direction, and $W_{\rm LSR}$ points to the Galactic North Pole. The $V_{\rm T}$ was calculated by the proper motion $\mu$ and parallax $\varpi$ in Gaia DR3: $V_{\rm T}=4.74\mu/\varpi$, and the total velocity was calculated by $V_{tot}=(U_{\rm LSR}^{2}+V_{\rm LSR}^{2}+W_{\rm LSR}^{2})^{1/2}$. As a result, there are 157 stars have $V_{\rm T}$. After removing spectroscopic binaries, there are 52 stars have velocities relative to the LSR. The results are shown in Figure 9.

We used the certiria in Gaia Collaboration et al. (2018) to select thin disk ($V_{\rm T}<40$ km s^-1 or $V_{\rm tot}<50$ km s^-1), thick disk ($60<V_{T}<150$ km s^-1 or $70<V_{\rm tot}<180$ km s^-1) and halo stars ( $V_{\rm T}>200$ km s^-1 or $V_{\rm tot}>200$ km s^-1 ). From Figure 9 we can see that there is no halo star, and only 8/157 ($V_{\rm T}>60$ km s^-2) or 3/52 ($V_{\rm tot}>70$ km s^-2) are thick disk stars, but speeds of these stars are far from the upper limit of thick disk stars. Therefore, GWAC flare stars are not old and most of them should be young.

There are 105/162 stars having periods as shown in Panel A of Figure 10. Medina et al. (2022) suggests that for stars with $P_{\rm rot}<10$ days, their ages are from about 600 Myr for 0.2 - 0.3 M_⊙ to about 2 - 3 Gyr for 0.1-0.2 M_⊙. Pass et al. (2022) also suggested that for stars with 0.2 - 0.3 M_⊙ and $P_{\rm rot}<2$ days, their ages are less than 600 Myr, and $2<P_{\rm rot}<10$ days means their ages between 1-3 Gyr. Therefore, these GWAC stars should be younger than 600 Myr - 3 Gyr. By comparing with Praesepe stars in Panel B of Figure 10, we refined their ages as younger than 670 Myr.

There are 38 stars having available $L_{\rm H\alpha}/L_{\rm bol}$ but without available periods. Kiman et al. (2021) proposed (see their Figure 8) that the activity of M-type stars decreases with age. In their Figure 6, for M2 - M6 stars younger than 1 Gyr, their $L_{\rm H\alpha}/L_{\rm bol}$ would significantly greater than $10^{-4}$, and then decreases abruptly to around $10^{-4.25}\approx 5.6\times 10^{-5}$ at 1 - 3 Gyr, while for stars of M7 - M9, $L_{\rm H\alpha}/L_{\rm bol}$ decreases abruptly to significantly $<10^{-4.5}\approx 3.2\times 10^{-5}$ during $\sim$ 30 Myr through 3 Gyr. As a result, for stars in Panel F of Figure 10, M2-M6 and M7-L1 stars should be younger than 1 and 3 Gyr, respectively.

Among the 157 stars with $V_{\rm T}>60$ km s^-1 in the left panel of Figure 9, Star #40, #132, and #145 have neither an available period nor $L_{\rm H\alpha}/L_{\rm bol}$. For the remaining five stars without available $V_{\rm T}$, only Star #19 has neither available period nor $L_{\rm H\alpha}/L_{\rm bol}$. Therefore, there are only 4/162 GWAC flare stars that cannot be determined if they belong to the thin disk. Therefore, except for these four stars, all other 158/162 stars should be young.

In summary, if we think these GWAC flare stars are from the same category, then their periods show they are younger than 670 Myr, $L_{\rm H\alpha}/L_{\rm bol}$ show they are younger than 3 Gyr, and $V_{\rm T}$ shows they are the thin disk stars.

For M and L stars, $L_{\rm H\alpha}/L_{\rm bol}$ decreases with increasing spectral type (Mohanty & Basri, 2003; Berger et al., 2010; Schmidt et al., 2015), which is also shown in Panel F of Figure 10, and the $L_{\rm H\alpha}/L_{\rm bol}$ level is significantly lower for stars later than M6 than for earlier type stars. Kiman et al. (2021) suggested that for M7-M9 stars, they have the same rotation-activity relation as earlier type stars only in their first $\sim$30 Myr, and after that their activity decreases abruptly and keeps at the low level for more than 3-6 Gyr.

Besides $L_{\rm H\alpha}/L_{\rm bol}$, Berger et al. (2010) also found that for the stellar activity there is a breakout around M7 in the quiescent X-ray and radio: for stars ranging from F to M6, they have the similar $L_{\rm X}/L_{\rm bol}$ pattern and $L_{\rm X}/L_{\rm bol}$ is higher than that of later type stars; $L_{\rm radio}/L_{\rm bol}$ and $L_{\rm radio}/L_{\rm X}$ are lower for M1-M6 stars than for later type stars.

In the GWAC sample, for M7-L1 stars, the distributions of their $\bigtriangleup G$ and $ED$ are similar to those of M2-M6 stars as shown in Figure 3, but their $\log_{10}(E_{\rm bol}/S)$ are lower than the median $\log_{10}(E_{\rm bol}/S)$ of M2 - M6 stars as shown in the right panel of Figure 11. Therefore, if the stellar radii given in TIC 8.2 are reliable, then the $\log_{10}(E_{\rm bol}/S)$ are lower for M7 - L1 stars than for M2-M6 stars. However, the number of M7-L1 stars is too small to draw a firm conclusion.

The lower activity level of later type stars has been noticed for a long time. Mohanty et al. (2002) suggested that as the surface temperature decreases, the resistivity of the stellar atmosphere increases, then the charged particles in the magnetic field are carried away due to collisions with neutral gases and thus the magnetic free energy used to produce the flare decreases with the spectral type. However, this theory still needs to be examined by further observations.