The white-light superflares from cool stars in GWAC triggers

Abstract

M-type stars are the ones that flare most frequently, but how big their maximum flare energy can reach is still unknown. We present 163 flares from 162 individual M2 through L1-type stars that triggered the GWAC, with flare energies ranging from $10^{32.2}$ to $10^{36.4}$ erg . The flare amplitudes range from $\triangle G=0.84$ to $\sim 10$ mag. Flare energy increases with stellar surface temperature ($T_{\rm eff}$) but both $\triangle G$ and equivalent duration $\log_{10}(ED)$ seem to be independent of $T_{\rm eff}$. Combining periods detected from light curves of TESS and K2, spectra from LAMOST, SDSS and the 2.16 m Telescope, and the Gaia DR3 data, we found that these GWAC flare stars are young. For the stars that have spectra, we found that these stars are in or very near to the saturation region, and $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ is lower for M7-L1 stars than for M2-M6 stars. We also studied the relation between GWAC flare bolometric energy $E_{\rm bol}$ and stellar hemispherical area $S$, and found that $\log_{10}E_{\rm bol}$ (in erg) increases with increasing $S$ (in cm²), and the maximum flare energy $\log_{10}E_{\rm bol,max}\geqslant\log_{10}S+14.25$. For M7-L1 stars, there seem to be other factors limiting their maximum flare energies in addition to stellar hemispherical area.

UTF8gbsn

1 Introduction

Solar flares originated from the release of magnetic energy by the magnetic reconnection in the corona (Haisch et al., 1991). The released energy ranges from $10^{28}$ to $10^{32}$ erg (Shibata & Magara, 2011). The biggest solar flare ever detected is the Carrington flare that occurred in 1859 (Carrington, 1859; Hodgson, 1859), which released energy of about $4\times 10^{32}$ erg (Hudson, 2021). The Carrington Event caused a geomagnetic storm (Hudson, 2021), and interrupted telegraph services (Boteler, 2006). The aurorae can be seen even near the Equator (Moreno Cárdenas et al., 2016).

Like the solar flares, stellar flares also originated from magnetic reconnections (Yan et al., 2021; Yang & Liu, 2019). Most stellar flare energies detected are between $10^{32}and10^{35}$ erg and for giants, their flare energies can even be as large as $10^{38}$ erg (Yang & Liu, 2019; Pietras et al., 2022). The white-light flare energies are often thought to be from blackbody radiations with a temperature of 9000 K - 14,000 K (Kowalski et al., 2013), and even as high as 42000 K (Howard et al., 2020), which is much higher than the photospheric temperature of M dwarfs, so flares with energies as low as of $10^{30}$ erg can be detected in near M-type dwarfs in optical bands (Yang et al., 2017).

The stellar flare activity is related to rotation and also the spectral type. The faster the rotation the stronger the activity, which is the stellar activity-rotation relationship. The relationship has been found by X-ray (Wright et al., 2011; Wright & Drake, 2016), Ca II H & K lines (Boudreaux et al., 2022; Lehtinen et al., 2021), H$\alpha$ (Newton et al., 2017; Li et al., 2023a). M dwarfs are more active than earlier type stars (Althukair & Tsiklauri, 2023) and from M0 to M6 the fraction of flare stars increases from about 10% to over 40% (Günther et al., 2020). M-type stars tend to produce more frequent and powerful flares than the Sun, and the habitable zones of M dwarfs are very near the hosts (Kane, 2018). Therefore, considering their flare energies can be thousands of times higher than the Carrington Event of the Sun and proximities of their habitable zones, the impacts of flares on habitable planets can be several magnitudes higher than those of solar flares on the Earth. The electromagnetic radiation from X-ray to radio (Osten et al., 2005) during flares and coronal mass ejections (CMEs) can be released by powerful flares: the more powerful the flare, the more likely the CME is released (Li et al., 2021). The intense flares can release tremendous UV fluxes and CME which can destroy O₃ (Tilley et al., 2019), and the ultraviolet fluxes can sterilize lives on the planet’s surface (Estrela & Valio, 2018). At the same time, the atmosphere of the planet would be heated, expanded, eroded (Linsky, 2019) and finally even disappeare (Atri & Mogan, 2021). On the other hand, the intense flare can trigger the prebiotic chemistry and then life (Rimmer et al., 2018; Xu et al., 2018; Günther et al., 2020; Chen et al., 2021).

The biggest flares are very rare, especially the flares with energy of $\geqslant 10^{34}$ erg from M dwarfs (Howard & MacGregor, 2022; Jackman et al., 2023), which would have important impacts on the planet’s atmosphere and life (e.g. Miranda-Rosete et al., 2023; Konings et al., 2022; Tilley et al., 2019; Hu et al., 2022). Some biggest flares have been detected by the GWAC, EvryFlare (Howard et al., 2019), ASAS-SN (Schmidt et al., 2019), and NGTS (Jackman et al., 2021). The amplitudes can be $\bigtriangleup R\sim 9.5$ or $\bigtriangleup V\sim 11.2$ (Xin et al., 2023b), $\bigtriangleup R\sim 9.5$ (Xin et al., 2021), $\bigtriangleup V>11$ mag (Schmidt et al., 2016), $\bigtriangleup V\sim 10$ mag (Jackman et al., 2019), and so on.

As the ground instrument of Space-based multi-band astronomical Variable Objects Monitor (SVOM; Wei et al., 2016), the Ground-based Wide Angle Cameras (GWAC) system aims to monitor afterglows of gamma-ray bursts (Xin et al., 2023a), and thus superstellar flares were also detected by the GWAC software (Han et al., 2021). In this work, we collected 163 big flares that occurred between 2017 November and 2023 March from 162 individual stars, and try to explore the mechanism behind them. These flares and properties of host stars are given in Table 1. Several of them have been carefully studied in Xin et al. (2021), Li et al. (2023b), Wang et al. (2022), Wang et al. (2021), and Bai et al. (2023). In Section 2 we will introduce the observations and data; in Section 3, we will present the properties of GWAC flares and their host stars; in Section 4, the rotation-age-activity relationship will be studied; in Section 5 the lower limit of the maximum flare energy that a star can produce is presented; the star age and activity pattern are discussed in Section 6; and at last, the conclusion is given in Section 7.

Table 1: GWAC flares and properties of stars

Star No.	GWAC Name	TIC/K2	Simbad Name	DR3 Name	RAJ2000	DEJ2000	Gmag	…
					(degree)	(degree)	(mag)	…
1	GWAC181208A	TIC 456482672	ATO J000.1128+13.6255	Gaia DR3 2767338005380427264	0.112725	13.62563	14.740	…
2	GWAC220106A	TIC 432551405	2MASS J00013265+3841525	Gaia DR3 2880981530065870720	0.386113	38.69797	18.193	…
3	GWAC180116A	TIC 357411008		Gaia DR3 2769578225960672512	2.798485	15.87654	16.542	…
4	GWAC181206A	TIC 405305098	2MASS J00113451+0659388	Gaia DR3 2742685889532913024	2.893805	6.99415	16.826	…
5	GWAC211124A	TIC 51940383	V* BI Psc	Gaia DR3 2767679884775564928	3.176556	13.13609	15.637	…
…

Note. —

$a$ ’EB’ means eclipse binary, and ’WDMS’ means white dwarf-main sequence binaries.

$b$ The spectral type is from their $G_{\rm BP}-G_{\rm RP}$.

(The full machine-readable form is available at https://nadc.china-vo.org/res/r101350/)

2 Observations and Data

2.1 GWAC

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/.

Figure 1: Two GWAC flares: GWAC181208B and GWAC190101A. The upper two panels are for GWAC181208B (from Star #18), and the bottom two panels are for GWAC190101A (from Star #39). The left two panels are their preflare statuses and the right two panels are the statuses at peaks. In each panel, the upper is the GWAC image, and the bottom is the GWAC light curve with time in minutes on the x-axis and flare amplitude on the y-axis. The black filled circle on the light curve indicates the flare amplitude in the image.

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.

2.2 TESS and K2 Light Curves

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/.

2.3 Flare Energy

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

$$ED=\int^{t1}_{t0}\frac{f(t)-f_{0}}{f_{0}}dt$$

. 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.

Table 2: Flares from TESS or K2 light curves

TIC/K2	BeginTime	EndTime	$\log_{10}E_{bol}$
	(days)	(days)	erg
TIC 156817806	2611.52955	2611.55733	31.78
TIC 156817806	2628.24472	2628.30028	32.35
TIC 156817806	2634.60568	2634.68901	32.26
TIC 156817806	2635.53621	2635.58482	32.05
TIC 318929976	1436.53726	1436.68309	33.38
TIC 13936933	1784.47441	1784.62024	33.03
TIC 118768009	2475.85775	2475.86748	33.17
TIC 118768009	2495.11225	2495.14420	33.88
TIC 118768009	2497.71940	2497.73468	33.63
TIC 303471889	2466.06797	2466.17215	33.85
TIC 303471889	2466.94999	2466.99861	33.55
…

Note. —

’BeginTime’ and ’EndTime’ are Julian days (JD 2547000 days) when a flare starts and ends.

(The full machine-readable form is available at https://nadc.china-vo.org/res/r101350/)

2.4 Spectra

We searched low-resolution spectra ($R\sim 1800$) in LAMOST DR10 (Deng et al., 2012; Cui et al., 2012) and SDSS 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 of the Spanish Virtual Observatory. All spectra cover H$\alpha$, and the equation

$$F(\lambda)=A_{0}\exp(-|(\lambda-A_{1})/A_{2}|^{A_{3}}/A_{3})+A_{4}\lambda+A_{5}$$ (1)

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):

$$BC_{J}=0.0115(G-J)^{3}-0.132(G-J)^{2}+0.735(G-J)+0.576$$ (2)

, 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

$$\log_{10}(l_{\rm bol})=0.4\times(-18.997351-m_{\rm bol})+3$$ (3)

That is,

$$\log_{10}(l_{\rm bol})=0.4\times(-18.997351-J-BC_{J})+3$$ (4)

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.

Table 3: $l_{\rm H\alpha}$, $l_{\rm bol}$, $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ and radial velocities obtained from spectra

Star No.	Telescope	$\log_{10}(l_{\rm bol})$	$\log_{10}(l_{\rm bol})$_err	$\log_{10}(l_{\rm H\alpha})$	$\log_{10}(l_{\rm H\alpha})$_err	$\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$	$\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$_err	rv	rv_err
		[erg cm-2 s-1]	[erg cm-2 s-1]	[erg cm-2 s-1]	[erg cm-2 s-1]			km s-1	km s-1
1	LAMOST	-10.10	0.0256	-13.76	0.0065	-3.65	0.0321	11.9	3.2
1	LAMOST	-10.10	0.0256	-13.72	0.0067	-3.62	0.0323	15.2	3.7
3	216	-10.58	0.0348	-14.07	0.0052	-3.49	0.0400
4	216	-10.81	0.0252	-14.22	0.0108	-3.41	0.0360
6	216	-10.46	0.0260	-14.30	0.0122	-3.83	0.0382
…

Note. — The full machine-readable form is available at https://nadc.china-vo.org/res/r101350/

3 Results

3.1 GWAC flares

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.

Figure 2: A LAMOST spectrum of TIC 336800108. The red line in the insert panel is the fitted H$\alpha$ emission.

Figure 3: The distributions of $\bigtriangleup G$ and $\log_{10}(ED)$ of GWAC flares.

Figure 4: The bolometric energies of GWAC flares. The blue circles are energies from GWAC flares, and the gray circles are the maximum flare energies of individual stars from Yang & Liu (2019).

3.2 Spectral Type

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, 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.

Figure 5: The H-R diagram and spectra type distribution of GWAC flare stars.

Figure 6: The spectra of Star #40 (GWAC201113A), #60 (GWAC191030A) and #152 (GWAC180212A). The spectra are normalized by dividing the fluxes at 7500 Å.

Figure 7: The SED fittings of Star #43 (GWAC190121B) and #76 (GWAC201113A). In each panel, the blue spectrum is the template for a main sequence star (MS), the yellow spectrum is a template for a white dwarf star (WD), the green spectrum is the WD+MS spectrum, and the red circles are the photometries of GALEX FUV, NUV, PS1 g, r, i, z, y, 2MASS J, H, K and WISE W1, W2 bands.

3.3 Flare Frequency Distribution

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

$$\log_{10}({\rm\nu})=\alpha\log_{10}(E_{\rm bol})+\beta$$ (5)

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.

Figure 8: The cumulative flare frequency distributions (FFDs) for 16 stars with more than 20 flares detected in TESS or K2 light curves. The cumulative flare frequencies are denoted by gray circles, and black lines are the fitted functions; The bolometric energies of GWAC flares are shown by red pentagons.

Table 4: Stars with More than 20 Flares

Star No.	TIC/K2	SpT	$\alpha$	$\alpha$_err	$\beta$	$\beta$_err	$\log_{10}\nu_{\rm GWAC}$	$\log_{10}\nu_{\rm GWAC}$_err
			$\left[{\rm day}^{-1}\right]$	$\left[{\rm day}^{-1}\right]$	$\left[{\rm day}^{-1}\right]$	$\left[{\rm day}^{-1}\right]$	$\left[{\rm day}^{-1}\right]$	$\left[{\rm day}^{-1}\right]$
11	EPIC 220404032	M4V	-0.6705	0.0294	21.57	0.98	-2.03	0.06
30	TIC 22944327	M4.5V	-0.4635	0.0246	14.00	0.81	-1.38	0.05
33	TIC 328254412	M5.5V	-0.5478	0.0262	17.12	0.86	-0.96	0.05
44	EPIC 210878635	M4V	-0.8935	0.0530	28.55	1.74	-2.29	0.07
49	EPIC 210417498	M4V	-1.1381	0.1041	36.01	3.38	-3.37	0.07
52	TIC 283866910	M4.5V	-0.9140	0.0745	29.08	2.39	-1.97	0.06
56	TIC 245936201	M3.5V	-0.9368	0.0993	30.56	3.32	-2.98	0.08
77	TIC 461654150	M8V	-0.6492	0.0423	19.80	1.36	-1.89	0.05
86	EPIC 212002525	M4.5V	-0.9269	0.0563	29.61	1.86	-2.79	0.07
87	TIC 175241416	M7V	-0.9632	0.0563	30.13	1.80	-1.99	0.08
93	TIC 471012520	M6V	-1.0777	0.0588	33.61	1.87	-1.66	0.09
96	TIC 251079483	M4.5V	-0.6915	0.0545	21.75	1.78	-1.46	0.09
113	TIC 289534997	M6V	-0.7162	0.0219	22.10	0.72	-2.40	0.04
117	EPIC 248853090	M4.5V	-0.9209	0.0426	29.69	1.42	-2.33	0.06
119	TIC 156151200	M5V	-0.6080	0.0332	18.71	1.08	-2.05	0.06
120	EPIC 201664337	M4V	-1.4847	0.1240	47.52	4.04	-3.95	0.13

Note. — $\alpha$ and $\beta$ are the parameters in Equation 5, and $\alpha$_err and $\beta$_err are respectively their errors. $\log_{10}\nu_{\rm GWAC}$ is the logarithm of the predicted frequency of the GWAC flare energy with the error of $\log_{10}\nu_{\rm GWAC}$_err.

3.4 Kinematics

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.

Figure 9: Left panel: the distributions of tangential velocities of 157 GWAC flare stars; Right panel: Toomre diagram for 52 GWAC flare stars with available radial velocities. The 1$\sigma$ errors are shown in gray lines.

4 Rotation-Age-Activity Relation

The periods of 105/162 stars were obtained from their TESS or K2 light curves by the Lomb-Scargle method. We inspected all of the light curves by eye and found there are no EA and EB eclipses, so the periods should be rotational periods. We also inspected all folded light curves by eye and found that most of them are sinusoidal, so if there are some EW binaries in GWAC flare stars, their real periods would be halved.

The periods are shown in Panel A of Figure 10. Most stars (80/105) have periods shorter than 2 days, and only 4 stars have periods longer than 10 days. In Panel B, GWAC flare stars are shown in red circles in the color-period diagram. The gray and blue circles are, respectively, field and Praesepe stars from Popinchalk et al. (2021). The gray circles show that the most field stars converge to the upper belt with periods greater than several 10 days. Praesepe has an age of $\sim$670 Myr, and its upper periods are about 10-20 days. There is a fast reservoir of $1.2<G-G_{\rm RP}<1.4$ and $0.2<P<2$ days (Popinchalk et al., 2021), where GWAC flare stars with available periods are also clustered. From this panel we can see that GWAC flare stars are still far from the upper period belt where the field stars cluster, and rotate even faster than Praesepe stars, so all these 105/162 GWAC stars should be younger than 670 Myr.

In the works of Wright et al. (2018) and Wright et al. (2011), the convective turnover time $\tau$ is a function of $V-Ks$, but the GWAC flare stars are very red and faint, and then some of them have no available $V$ data. As a result, we firstly obtained $V-Ks$ of some GWAC flare stars from Paegert et al. (2022), then fitted the relation of $V-Ks$ vs. $G_{\rm BP}-G_{\rm RP}$, and finally obtained their predicted $V-Ks$ from $G_{\rm BP}-G_{\rm RP}$:

$$V-Ks=0.15226\times(G_{\rm BP}-G_{\rm RP})^{2}+0.67908\times(G_{\rm BP}-G_{\rm RP})+2.1097$$ (6)

which is shown in Panel C. Then we used the function given by Wright et al. (2018) to calculate the convective turnover time $\tau$:

$$\log_{10}(\tau)=0.64+0.25\times(V-Ks)$$

.

There are 108 spectra with available $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ for 82 GWAC flare stars, which also have available periods. Their $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ vs. Rossby number ($Ro=P/\tau$) diagram is shown in Panel D. For stars that have more than one available $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$, all $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ are shown. We also used the equation presented by Reiners et al. (2014) to calculate the saturation period $P_{\rm sat}$ in days:

$$P_{\rm sat}=(1.1\times 10^{34}/L_{\rm bol})^{1/2}$$ (7)

We found that only Star #51 and #56 are around the saturation $Ro$ and also $P_{\rm sat}$ which are shown in Panel D and E of Figure 10 by black horizontal lines, and all other stars are in the saturation region. There are 117 GWAC flare stars have available $L_{\rm H\alpha}$ in total, and their $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ vs. $G_{\rm BP}-G_{\rm RP}$ diagram is shown in Panel F. In Panel D, E, and F, the blue triangles are stars of M7-L1 and the red circles are stars of M2-M6. For red circles, the filled ones are stars with $P<2$ days and the empty ones are stars with $P>2$ days. The central gray line is the saturation line ($L_{\rm H\alpha}/L=1.49\times 10^{-4}$) and gray dotted lines are the 1$\sigma$ ($\sigma=0.26$) positions given in Newton et al. (2017). From Panel D, E and F, we can see that:

• Stars in Panel D and E are all in or very near to the saturation region, so these stars should be very active.

• Though in the saturation region, $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ are significantly lower for M7 - L1 stars (blue triangles) than for M2 - M6 stars (red circles).

• For M2 - M6 stars (red circles), it seems that $L_{\rm H\alpha}/L_{\rm bol}$ slightly decreases with increasing $Ro$ (as seen in Panel D) and also increasing period (as seen in Panel E), but is independent of $G_{\rm BP}-G_{\rm RP}$ (as seen in Panel F). Especially, the $L_{\rm H\alpha}/L_{\rm bol}$ are lower for stars with $P>2$ days (red empty circles) than for those with $P<2$ days (red filled circles). This declination of $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ with increasing period may imply the increasing stellar age (Kiman et al., 2021).

• As shown in Panel F, $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ seems to decrease with increasing $G_{\rm BP}-G_{\rm RP}$ (i.e. spectral type), which was also presented by Berger et al. (2010) and Schmidt et al. (2015).

In summary, for M2 - M6 stars, their $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ decreases with increasing period and is independent of spectral type, but for M7 - L1 stars, though they are rapid rotators, their $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ are significantly lower than those of M2 - M6 stars.

The functions to calculate the convective overturn time given by Wright et al. (2018, 2011) and Reiners et al. (2014) are frequently used in literature, but these equations were obtained from the stars with masses larger than 0.1 M_⊙ or spectral types no later than M6 (Wright et al., 2018, 2011; Jao et al., 2022). Therefore, we do not know if Ro and saturation periods derived from these functions can apply to stars later than M6. However, in any case, it is true that the $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ is lower for M7 - L0 stars than for M2 - M6 stars as shown in Panel D, E, F of Figure 10, because it is also seen in literature (e.g. Berger et al., 2010; Schmidt et al., 2015; Kiman et al., 2021).

Figure 10: The rotation-age-activity relation. Panel A: The distribution of periods of 105 GWAC flare stars. Panel B: The color-period diagram. The gray and blue circles are, respectively, field stars and Praesepe stars from Popinchalk et al. (2021), and the red circles are GWAC flare stars. Panel C: $V-Ks$ vs. $G_{\rm BP}-G_{\rm RP}$. The fitted relationship (Equation 6) is shown in the red curve. Panel D: $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ vs. Rossby number ($Ro$) for stars with available spectra and $Ro$. Panel E: $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ vs. period for stars with available spectra and periods. Panel F: The $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ vs. $G_{\rm BP}-G_{\rm RP}$ diagram for 117 GWAC stars that have available spectra. The dotted yellow vertical line is $G_{\rm BP}-G_{\rm RP}=4.05$, which is roughly the demarcation for M6 and M7 stars. In Panel D, E, and F, the blue triangles are M7-L1 stars and the red circles are stars M2-M6 stars. For red circles, the filled ones are stars with $P<2$ days and the empty ones are stars with $P>2$ days. For stars that have more than one available $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$, all $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ are shown. The central gray line is the saturation line ($L_{\rm H\alpha}/L=1.49\times 10^{-4}$) and gray dotted lines are the 1$\sigma$ ($\sigma=0.26$) positions given in Newton et al. (2017). The black horizontal short lines are Star #51 and #56, which are around the saturation threshold.

5 Maximum Flare Energy

Some works (e.g. Notsu et al., 2019; Shibata et al., 2013) suggested that for a solar-like star, the upper limit of flare energy is usually determined by the spot size, and thus the stellar hemisphere. Inspired by this, we inspected the relation of the flare bolometric energy $E_{\rm bol}$ vs. the area of stellar hemisphere $S=2\pi R^{2}$ ($R$ is the stellar radius from TIC 8.2 ), as shown in the left panel of Figure 11. From this panel we can see that $E_{\rm bol}$ increases with increasing $S$. The relation of flare energy divided by stellar hemisphere area $\log_{10}(E_{\rm bol}/S)$ vs. stellar surface effective temperature $T_{\rm eff}$ ($T_{\rm eff}$ from TIC 8.2), is shown in the right panel of Figure 11, where the vertical dashed line is $G_{\rm BP}-G_{\rm RP}=3.95$ and roughly the demarcation between M2 - M6 and M7-L1 stars, and the red horizontal line is the median value of $\log_{10}(E_{\rm bol}/S)$ of M2-M6 stars. From this panel we can see that for M2-M6 stars, the distribution of $\log_{10}(E_{\rm bol}/S)$ is roughly same and independent of $T_{\rm eff}$, and significantly higher than those of M7-L1 stars. We found that the maximum $\log_{10}(E_{\rm bol}/S)$ is about 14.25 for GWAC flares in the right panel, which means that for M2-M6 stars, the lower limit of the maximum flare energy is

$$\log_{10}E_{\rm bol}=\log_{10}S+14.25$$ (8)

which is also shown by the yellow dashed line in the left panel of Figure 11.

Figure 11: The $E_{\rm bol}$-$S$-$T_{\rm eff}$ relation. The red circles are M2-M6 stars, and the blue triangles are M7-L1 stars. The yellow dashed line in the right panel is the lower limit of the maximum flare energy, which is $\log_{10}E_{\rm bol}=\log_{10}S+14.25$. In the right panel, the vertical dotted line is $G_{\rm BP}-G_{\rm RP}=3.95$ and roughly the demarcation between M2-M6 and M7-L1 stars, and the red horizontal line is the median value of $\log_{10}(E_{\rm bol}/S)$ of M2-M6 stars.

6 Discussion

6.1 Stellar Age

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.

6.2 The Lower Activity Level of stars later than M6

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.

7 Conclusion

In this work, we presented 163 big flares in GWAC triggers from 162 individual stars with the spectral type from M2 to L1. The flare amplitude $\bigtriangleup G$ ranges from 0.83 to $\sim$10 mag, and the flare energy ranges from $10^{32.2}$ to $10^{36.4}$ erg. From TESS or K2 light curves, we found 1478 flares from 117 stars and calculated periods for 105 stars. Besides, we obtained 154 low-resolution spectra for 120 individual stars with the 2.16 m Telescope in the Xinglong Station, LAMOST and SDSS. Among them, there are 108 available $L_{\rm H\alpha}/L_{\rm bol}$ for 82 individual stars with available periods. We also obtained tangent velocities $V_{\rm T}$ for 157 stars , and velocities relative to LSR for 52 stars. From these data we found that:

• The energy of GWAC flares decreases with increasing $G_{\rm BP}-G_{\rm RP}$ (decreasing stellar surface effective temperature $T_{\rm eff}$) (see Figure 4), but both distributions of flare amplitude $\bigtriangleup G$ and flare equivalent duration $\log_{10}(ED)$ seem to be independent of $T_{\rm eff}$ (see Figure 3).

• If 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}$ show they are the thin disk stars (see Figure 9 and 10). Therefore, they are young stars.

• There are 16 stars with flare numbers greater than 20. From their FFDs 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 years (see Figure 8).

• Stars with spectra and periods are in or very near to the saturation region, so these stars should be very active.

• $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ is higher for M2-M6 stars than for M7-L1 stars. $\log_{10}(L_{\rm H\alpha}/L_{\rm bol})$ of M2-M6 stars is around the saturation level, and decreases with increasing period, while for M7-L1 stars, though they are all fast rotators, their $L_{\rm H\alpha}/L_{\rm bol}$ are significantly lower (see Figure 10).

• The flare energy divided by the stellar hemispherical area $\log_{10}(E_{\rm bol}/S)$ seems to be higher for M2 - M6 stars than for M7-L1 stars (see Figure 11). However, the number of M7-L1 stars is too small to draw a firm conclusion.

• The maximum flare energy of M2 - M6 stars should be larger than $\log_{10}E_{\rm bol}=\log_{10}S+14.25$ (see Figure 11).

The different activity levels between M7-L1 and M2-M6 stars had been reported by Berger et al. (2010) in X-ray, radio and H$\alpha$, but we still do not know if the different activity levels exist in $\log_{10}(E_{\rm bol}/S)$, which should be explored by the sample including a lot of stars later than M6 that release big white-light flare energies. The study would be very interesting because the origin of the activity pattern of cool stars is still unclear.