Complex Rotational Modulation of Rapidly Rotating M-Stars Observed with TESS

TL;DR

This study reports the discovery of highly structured, rapidly rotating M-dwarfs observed with TESS, exploring their unusual light-curve features and potential causes like dust trapping, with implications for stellar magnetosphere interactions.

Contribution

It presents the first detailed analysis of highly structured rotational modulation in M-dwarfs from TESS data, proposing new explanations involving dust and magnetospheric phenomena.

Findings

Identified 10 M-dwarfs with rotation periods less than 1 day and complex light-curve structures.

Discovered 371 M-dwarfs with rotation periods below 4 hours, including the shortest at 1.63 hours.

Structured light-curves cannot be explained by star-spot models alone.

Abstract

We have searched for short periodicities in the light curves of stars with $T_{\rm eff}$ cooler than 4000 K made from 2-minute cadence data obtained in TESS sectors 1 and 2. Herein we report the discovery of 10 rapidly rotating M-dwarfs with highly structured rotational modulation patterns among 10 M dwarfs found to have rotation periods less than 1 day. Star-spot models cannot explain the highly structured periodic variations which typically exhibit between 10 and 40 Fourier harmonics. A similar set of objects was previously reported following K2 observations of the Upper Scorpius association (Stauffer et al. 2017). We examine the possibility that the unusual structured light-curves could stem from absorption by charged dust particles that are trapped in or near the stellar magnetosphere. We also briefly explore the possibilities that the sharp structured features in the lightcurves are produced by extinction by coronal gas, by beaming of the radiation emitted from the stellar surface, or by occultations of spots by a dusty ring that surrounds the star. The latter is perhaps the most promising of these scenarios. Most of the structured rotators display flaring activity, and we investigate changes in the modulation pattern following the largest flares. As part of this study, we also report the discovery of 371 rapidly rotating M-dwarfs with rotational periods below 4 hr, of which the shortest period is 1.63 hr.

Tables

Table 1. Table 2: Rapidly Rotating M Stars with Complex Rotation Modulation Profiles

Object Name	38820496	177309964	201789285	206544316	224283342	234295610	289840928	332517282	425933644	425937691
RA J2000a	7.19516	103.45258	33.88868	18.41884	356.35758	357.98325	317.63150	350.87860	3.69942	5.36556
Dec J2000a	$-67.86237$	-75.70396	$-56.45488$	$-59.65974$	$-40.33782$	$-64.79293$	$-27.18311$	$-28.12114$	$-60.06352$	$-63.85226$
$T_{\_}\mathrm{eff}$ [K]b	2987	3289	…	3112	3198	3082	3120	3032	3176	2733
d [pc] a	$44.1\pm 0.1$	$91.3\pm 0.5$	$45.4\pm 0.2$	$43.1\pm 0.1$	$38.2\pm 0.1$	$48.2\pm 0.1$	$40.4\pm 0.2$	$39.1\pm 0.1$	$44.3\pm 0.2$	$43.9\pm 0.2$
pmra [mas/yr] a	$+89.2\pm 0.1$	$+4.7\pm 0.1$	$+95.3\pm 0.1$	$+96.7\pm 0.1$	$+108.4\pm 0.1$	$74.6\pm 0.1$	$+67.8\pm 0.1$	$+109.3\pm 0.1$	$+88.8\pm 0.2$	$+91.6\pm 0.2$
pmdec [mas/yr]a	$-50.9\pm 0.1$	$39.3\pm 0.1$	$-14.7\pm 0.1$	$-41.0\pm 0.1$	$-78.8\pm 0.1$	$-59.7\pm 0.1$	$-75.1\pm 0.1$	$-123.5\pm 0.1$	$-56.9\pm 0.2$	$-53.9\pm 0.2$
Ga	$14.77\pm \left[2\right]$	$14.45\pm \left[2\right]$	$15.63\pm \left[2\right]$	$12.95\pm \left[2\right]$	$13.66\pm \left[1\right]$	$13.89\pm \left[1\right]$	$13.60\pm \left[1\right]$	$14.74\pm \left[3\right]$	$12.71\pm \left[1\right]$	$14.76\pm \left[2\right]$
G${}_{\_}{}^{}\mathrm{BP}$a	$16.80\pm 0.01$	$16.13\pm 0.01$	$18.04\pm 0.01$	$14.61\pm \left[6\right]$	$15.42\pm \left[4\right]$	$15.65\pm \left[7\right]$	$15.56\pm \left[4\right]$	$16.70\pm 0.02$	$14.31\pm \left[6\right]$	$17.12\pm 0.02$
G${}_{\_}{}^{}\mathrm{RP}$a	$13.42\pm \left[4\right]$	$13.19\pm \left[6\right]$	$14.19\pm \left[5\right]$	$11.72\pm \left[5\right]$	$12.39\pm \left[2\right]$	$12.61\pm \left[4\right]$	$12.27\pm \left[2\right]$	$13.41\pm 0.01$	$11.48\pm \left[3\right]$	$13.34\pm \left[7\right]$
Jc	$11.40\pm 0.03$	$11.36\pm 0.03$	$11.93\pm 0.02$	$9.95\pm 0.02$	$10.53\pm 0.02$	$10.72\pm 0.026$	$10.30\pm 0.02$	$11.41\pm 002$	$9.71\pm 0.02$	$11.02\pm 0.02$
Hc	$10.77\pm 0.02$	$10.73\pm 0.03$	$11.32\pm 0.03$	$9.34\pm 0.03$	$9.95\pm 0.02$	$10.17\pm 0.024$	$9.71\pm 0.03$	$10.83\pm 003$	$9.10\pm 0.02$	$10.48\pm 0.03$
Kc	$10.5\pm 0.02$	$10.48\pm 0.03$	$10.95\pm 0.02$	$9.06\pm 0.03$	$9.69\pm 0.02$	$9.88\pm 0.021$	$9.41\pm 0.02$	$10.51\pm 0.02$	$8.83\pm 0.02$	$10.11\pm 0.03$
$M_{\_}K$c,a	7.278	5.680	7.666	5.888	6.784	6.464	6.379	7.545	5.602	6.896
W1d	$10.31\pm 0.02$	$10.35\pm 0.02$	$10.66\pm 0.02$	$8.90\pm 0.02$	$9.50\pm 0.02$	$9.71\pm 0.022$	$9.06\pm 0.02$	$10.32\pm 0.02$	$8.70\pm 0.02$	$9.91\pm 0.02$
W2d	$10.09\pm 0.02$	$10.15\pm 0.02$	$10.40\pm 0.02$	$8.71\pm 0.02$	$9.30\pm 002$	$9.51\pm 0.019$	$8.81\pm 0.02$	$10.13\pm 0.02$	$8.51\pm 0.02$	$9.66\pm 0.02$
W3d	$9.95\pm 0.05$	$10.00\pm 0.04$	$10.16\pm 0.05$	$8.58\pm 0.02$	$9.13\pm 0.03$	$9.27\pm 0.034$	$8.77\pm 0.03$	$9.90\pm 0.05$	$8.36\pm 0.02$	$9.40\pm 0.03$
W4d	$\gtrsim 8.65$	$\gtrsim 9.35$	$\gtrsim 8.5$	$8.46\pm 0.26$	$\gtrsim 8.5$	$\gtrsim 8.46$	$8.41\pm 0.33$	$\gtrsim 8.5$	$8.17\pm 0.25$	$\gtrsim 8.5$
Radius [$R_{\_}\odot$]e	$0.28\pm 0.01$	$0.49\pm 0.02$	$0.24\pm 0.01$	$0.46\pm 0.01$	$0.34\pm 0.01$	$0.38\pm 0.01$	$0.38\pm 0.01$	$0.25\pm 0.01$	$0.51\pm 0.02$	$0.32\pm 0.01$
Mass [$M_{\_}\odot$]e	$0.25\pm 0.01$	$0.49\pm 0.02$	$0.18\pm 0.01$	$0.35\pm 0.01$	$0.25\pm 0.01$	$0.37\pm 0.01$	$0.29\pm 0.01$	$0.19\pm 0.01$	$0.38\pm 0.01$	$0.30\pm 0.01$
Period [hr]f	$15.727\left(6\right)$	$10.881\left(6\right)$	$3.638\left(1\right)$	$7.727\left(3\right)$	$21.298\left(11\right)$	$18.288\left(13\right)$	$4.780\left(3\right)$	$9.660\left(9\right)$	$11.669\left(6\right)$	$4.816\left(2\right)$
Associationf	THA	CAR	THA/CAR	ABDMG	COL	THA	BPMG	ABDMG	THA	THA
Data sectorsg	1, 2	1, 2	2	1, 2	2	1	1	2	1, 2	1, 2

Table 2. Table 3: Data from Gagné et al. ( 2018 )

Name	Abbr.	Ave. Dist. [pc]	Age [Myr]
AB Doradus	ABDMG	${30}_{\_}^{+20}-10$	${149}_{\_}^{+51}-19$
$\beta$ Pictoris	$\beta$PMG	${30}_{\_}^{+20}-10$	$24\pm 3$
Carina	CAR	$60\pm 20$	${45}_{\_}^{+11}-7$
Columba	COL	$50\pm 20$	${42}_{\_}^{+6}-4$
Tucana Horologium	THA	${46}_{\_}^{+8}-6$	$45\pm 4$

Equations

$$R_{\_}\ast$$

$${\log }_{\_}10(M_{\_}\ast )$$

$$\frac{\Omega }{{\Omega }_{\_}\mathrm{crit}}=\frac{2\pi }{P_{\_}\ast }\sqrt{\frac{R_{\_}{\ast }^3}{G{M}_{\_}\ast }}\simeq 2.8\sqrt{\frac{r_{\_}{\ast }^3}{m_{\_}\ast }}\left(\frac{\mathrm{hr}}{P_{\_}\ast }\right)$$

$$R_{\_}g=\frac{m{v}_{\_}\perp c}{qB}(\mathrm{cgs})$$

$$\frac{R_{\_}g}{R_{\_}\ast }$$

$$d\gtrsim R_{\_}\ast \csc (\pi f).$$

$${\mathcal{M}}_{\_}\mathrm{ring}\simeq \pi \left(R_{\_}{\mathrm{out}}^2-R_{\_}{\mathrm{in}}^2\right)\Sigma$$

$$\sigma ={\sigma }_{\_}\mathrm{Th}(1-u{)}^{-2}\left[1+u+(u/2)(u-3){\sin }^2\theta \right]$$

$$\sigma \simeq {\sigma }_{\_}\mathrm{Th}\left[1+3\left(1-{\sin }^2\theta /2\right)u\right]$$

Full text

stars: variables: general - circumstellar matter - stars: flare - stars: low-mass - stars: pre-main sequence - starspots

††journal: ApJ

1 Introduction

Young M stars often are found to rotate rapidly, i.e., with periods $\lesssim 2$ days, and then are seen to be ‘active’ in a number of different respects. The spectra of such stars may contain hydrogen as well as Ca H&K emission lines (Pettersen & Hawley, 1989); flares on such stars are common (e.g. Davenport 2019; Mondrik et al. 2019; Günther et al. 2019); and the stellar surfaces often have spots that result in large-amplitude rotational modulations (e.g., EY Dra, Vida et al. 2010; GJ 65A and 65B, Barnes et al. 2017). Even though the spots may be distributed around the stellar surface, the rotational modulation patterns are usually fairly simple, i.e., they may be adequately represented by just a few Fourier components. A small percentage of rotationally modulated stars, e.g. KIC 7740983 (Rappaport et al., 2014) and HSS 348 (Strassmeier et al., 2017), exhibit more complex rotational modulations; even these can often be ascribed to particular sizes and arrangements of spots. Sometimes two distinct rotation frequency signatures are seen, indicating the presence of two ‘twin-like’ M stars (Rappaport et al., 2014).

Two other quite distinct categories of young rotating M stars have recently been identified. First, there is a class of young late-type stars that exhibit ‘dipping activity’ (see, e.g., Cody et al. 2014; Ansdell et al. 2016). Such stars can exhibit dips with depths that range from $\sim$5% to 50%, and that may recur periodically, quasi-periodically, or without any apparent underlying pattern. These stars typically exhibit strong emission in the WISE 3 and 4 bands, indicating that dusty disks are present. Second, there is a class of stars found by Stauffer et al. (2017, 2018) that may be related to the dippers, that the authors describe as “rapidly rotating weak-line T Tau stars with orbiting clouds of material at the Keplerian co-rotation radius”. These stars were found in the Upper Scorpius and $\rho$ Oph star forming associations, and exhibit highly structured rotational modulations that, according to the authors, cannot be explained by star spots.

In this work we present the discovery with the Transiting Exoplanet Survey Satellite TESS (Ricker et al., 2015) of 10 rapidly rotating M stars in the Southern Hemisphere with quite complex rotational modulation patterns (the ‘10 complex rotators’). Half of these stars are in the star-forming stellar association Tucana Horologium (‘THA’). In Section 2 we discuss the data from Sectors 1 and 2 that we searched for interesting rapid rotators. In Section 3 we describe the Fourier transform survey of all the 2-minute cadence targets in Sectors 1 and 2, and discuss the Fourier amplitude spectra, the rotational modulation patterns, and other characteristics of the 10 complex rotators. Low and medium resolution spectra of a sample of the stars are presented in Section 4. We attempt to model the light curves with a starspot modeling code in Section 5, but conclude that such a model falls short of explaining the observations. In Section 6 we discuss another model involving charged dust particles trapped in the stellar magnetosphere, but find that this scenario also has difficulties. A model, perhaps the most promising one considered here, that involves a spotted host star rotating obliquely inside a dusty ring is presented in Section 7. In Section 8 we raise the possibility that the highly structured rotational modulations might be due to beamed emission. Changes in the structured modulation pattern around the time of a superflare in one of the rapid rotators are discussed in Section 9. In Section 10 we discuss the stellar associations where these 10 stars are found, and their likely ages. We summarize our work and discuss some related issues in Section 11.

[FIGURE:]

2 TESS Observations

In this study we utilized the 2-minute cadence (short cadence) data collected by TESS during Sectors 1 and 2. The targets selected for short cadence are chosen from the candidate target list (CTL) (Stassun et al., 2018). The light curves were produced by the Science Processing Operations Center (SPOC) pipeline (Jenkins et al., 2016), which is operated by NASA Ames Research Center, and are available online111http://archive.stsci.edu/tess/bulkdownloads.html. For Sectors 1 and 2, 2-minute cadence data were obtained for a total of 24,809 unique targets of which 7,074 were observed in both sectors. Although the primary goal of TESS is to search for terrestrial planets that transit nearby bright stars, the large number of targets monitored with 2-minute cadence enables the study of other astrophysical phenomena on time scales down to the Nyquist limit of 360 cycles d*-1* (4-minute period).

We opted to use the SAP_FLUX light curve data and do our own detrending. For each light curve, we split the non-detrended SAP_FLUX time series into separate orbits and then normalized the data for each orbit by the mean flux of that orbit. Since we were primarily interested in finding rapidly rotating M dwarfs, we detrended each orbit using a 3-day spline filter to remove any longer-timescale fluctuations in the flux.

We note that many of the M dwarfs in this study were added to the CTL as a part of the Cool Dwarf specially curated catalog (Muirhead et al., 2018).

3 Analysis

3.1 Fourier Analysis

We used a basic fast Fourier transform (FFT) method to search for rapidly rotating M dwarf stars. The FFT method is fast and efficient when used to search for periodic signals that have high-duty cycles, or, equivalently, smoothly varying rotation profiles. The FFT search used here is similar to the one employed by Sanchis-Ojeda et al. (2013, 2014) but is optimized for the TESS data. Given that the data set essentially consists of samples uniformly spaced in time and minimally affected by window functions, we did not see that any particular advantage would be gained by using the Lomb-Scargle transform (Lomb, 1976; Scargle, 1983). Accurate periods and the corresponding uncertainties for any interesting objects are found via a Stellingwerf transform (Stellingwerf, 1978) which utilizes a phase dispersion minimization algorithm.

The preprocessed data for each orbit are stitched together into one data array, the mean flux is subtracted, and the data gaps are filled with zeroes. The data array was extended by padding with zeroes to a power of 2 number of bins which represents approximately 8 times the length of the non-zero part of the data array. This yields Fourier ‘interpolation’ in frequency. In the case where only 1 TESS Sector (2 satellite orbits) of data ($\sim$20,000 data points) is available, the final data array comprises 217 bins. We then compute the FFT using Python with the SciPy FFT package, which takes $\approx 0.1$ second for each object on a computer with a 2.9 GHz Intel Core i7 Processor. From the complex FFT we form the amplitude spectrum.

We apply the data preparation and Fourier transform operations to the data for all of the 2-minute cadence target stars. However, in this paper we focus on the results of our search for rapidly rotating M stars.

3.2 Finding rapidly rotating M stars

To find rapidly rotating M dwarfs, we first select all stars among the 2-minute cadence targets with reported stellar temperatures $T_{\_}{\rm eff}\lesssim$ 4000 K. Rapidly rotating M dwarfs are known to be active flaring stars; flares can add significant noise to the Fourier transform, and thereby adversely affect the search (see Sect. 3). Therefore flares are removed prior to application of the FFT using a 5-$\sigma$ flare detection threshold and flare fitting technique developed and discussed in Günther et al. (2019).

To detect rapidly rotating M-dwarfs, we compute a local, i.e., 128 frequency bin, mean and local rms ($\sigma$) for each bin in each Fourier amplitude spectrum. We consider a detection to be significant if at least one value in the amplitude spectrum exceeds the local mean by at least 5 local $\sigma$ anywhere in the frequency range222As noted earlier, the TESS short-cadence frequency range reaches 360 cycles d*-1*; however, we truncated our search at 48 cycles d*-1* because even an M star of mass 0.1 $M_{\_}\odot$ and radius 0.15 $R_{\_}\odot$ cannot rotate with a period below 1/2 hr without breaking up due to centrifugal forces (see discussion in Sect. 3.3 and Eqns. 1). of 1 - 48 cycles d*-1*, and, is accompanied by at least one higher harmonic or subharmonic that exceeds its local mean by 3 local $\sigma$. We vetted the targets that pass this threshold by eye to confirm detection, and, in the process, eliminated eclipsing M-star binaries and objects with other obvious non-rotational signals. The result was the detection of approximately 500 potentially rapidly rotating M dwarfs with periods shorter than 1 day in either one or both sectors, which accounts for roughly 1/10th of the M dwarf stars in the CTL for Sectors 1 and 2.

3.3 The collection of rapidly rotating M stars

We found a total of 17 M-dwarf stars to be rotating with periods less than 4 hours among 371 M-dwarfs with periods below 1 day and $\sim$1000 with periods below 2 days. With a period of 1.63 hr, TIC 204294226 is the most rapidly rotating M dwarf discovered thus far.

The 17 fastest rotators are listed in Table LABEL:tbl:rotators along with coordinates and effective temperatures. The ratio of rotation frequency to breakup angular frequency is also listed for each star. The estimated breakup frequency of each star is based on the stellar mass, $M_{\_}*$, and radius, $R_{\_}*$, expressed as simple functions of the absolute K magnitude, $M_{\_}K$, following Mann et al. (2016, 2018):

[TABLE]

where $M^{\prime}_{\_}K\equiv M_{\_}K-7.5$, and the coefficients $a-f$ are:

$a=-0.642$, $b=-0.208$, $c=-0.000843$, $d=0.00787$, $e=0.000142$, and $f=-0.000213$. The ratio of $\Omega/\Omega_{\_}{\rm crit}$ is computed from the expression:

[TABLE]

where $m_{\_}*$ and $r_{\_}*$ are the mass and radius in solar units and $P_{\_}*$ is the rotational period in hours.

3.4 Identifying highly structured rotation profiles

In the process of vetting the rapidly rotating M stars, we discovered a number of stars with highly structured rotational modulation patterns. The first clue came from examining the unusual harmonic structure in their amplitude spectra. An extreme example of this is shown in Fig. 1. While the amplitude spectrum of a rotating M star usually contains 2-4 detectable harmonics, TIC 206544316 exhibits 40 higher harmonics. This is indicative of a highly structured rotation profile.

The raw light curves for four of these objects are shown, as examples, in Fig. 2. Note that each of these has sharp structure that occurs over small fractions of a rotation cycle.

In all, we find 10 similarly structured rotation curves among our sample of rapid rotators. Each of these objects exhibits 10 to 40 significant harmonics in the frequency range 1–96 cycles/day. These contrast with most rapidly rotating M dwarfs that have comparatively simple modulation patterns and therefore exhibit only 1–5 significant harmonics. The folded light curves for these ten rotators with complex modulation patterns are presented in Fig. 3 and some of their properties are presented in Table LABEL:tbl:weird.

We see in Figure 3 an array of exotic rotational modulation patterns. To the extent we can generalize, the amplitudes of the modulations are typically 5% to 10% peak to peak. The amplitudes and structures of the profiles seem roughly independent of rotation period between 2.5 and 21 hours. The structures of the profiles are not simple to describe; the term “scalloped shell” (Stauffer et al., 2017), while descriptive, does only partial justice. ‘Highly structured’, ‘comprised of many Fourier components’, and ‘multi- and sharply-peaked’ further describe the profile shapes.

3.5 Variability of the rotation profiles over time

Fig. 4 shows the light curves of four of our highly structured rotators for which data was obtained from both Sectors 1 & 2, i.e., from TESS orbits 9, 10, 11, and 12333Although there is data for TIC 177309964 from both sectors, they are not plotted because the modulation pattern is quite stable over the course of two months.. Each color-coded folded profile represents the data taken in one TESS orbit, i.e., two weeks. While there are clearly variations in the detailed shapes of the rotational modulations, the basic profiles remain recognizably similar over two months.

In Fig. 5 we show the time evolution of the Fourier amplitude spectrum for TIC 332517282. For these calculations we used the software package TiFrAn (for details see Rappaport et al. 2014). Even though the amplitude of some of the features in the rotational modulation pattern are clearly changing with time, the fundamental frequency and its numerous harmonics remain essentially constant to within the frequency resolution over the course of a month.

3.6 Comparison with highly structured rotators found in Upper Sco and $\rho$ Oph

The rotational modulation profiles of the 10 M stars discussed in this work are very similar to those that Stauffer et al. (2017) describe as “scalloped shell” modulation profiles. Of course, the 2-minute time resolution of the TESS data, when compared to the 30-minute resolution of the K2 data, allows for higher-frequency structure to be recorded. The ‘scalloped shell’ rotators in Upper Sco and $\rho$ Oph have a typical age of 5-10 Myr, a typical distance of 135 pc, and a median $K_{\_}s$ magnitude of 10.5 (median $M_{\_}K\simeq 4.9$). The comparable rotators in this work are $\simeq 45$ pc away (closer), are $\simeq 45$ Myr old (older), and have a similar median $K_{\_}s$ magnitude, but a median $M_{\_}K\simeq 6.4$, which is notably 1.5 magnitudes less luminous than the ‘scalloped shell’ modulators in Upper Sco and $\rho$ Oph. We also note that the range of periods of our 10 systems of 2.5 to 21 hours is comparable to the range of 6.2 to 15.4 hr of the ‘scalloped shell’ modulators in Upper Sco and $\rho$ Oph. Finally, we point out that the WISE W3 fluxes of the stars in one set are similar to those of the other set, and the upper limits in the W4 band for the two sets of stars are similar.

4 Spectroscopic Observations

We obtained low and medium resolution spectra for a select sample of the rapidly rotating M-dwarfs with complex rotational modulations to search for activity signatures and additional signs of rapid rotation. Low resolution observations were obtained for four of the systems using the Wide Field Spectrograph (WiFeS, Dopita et al., 2007) on the ANU 2.3-m telescope at Siding Spring Observatory, Australia, on the nights of 2019 January 18 and 19. The data covered the wavelength region of $5200-7000$ Å at a resolution of $R\equiv\lambda/\Delta\lambda=7000$, and were reduced per the procedure described in Bayliss et al. (2013). The spectra are displayed in Figure 6. The spectra of all of the M-dwarfs sampled by WiFeS exhibit strong H$\alpha$ emission features with equivalent widths in the range $4-6$ Å, and none shows signs of the composite nature expected if a system is composed of multiple stars. These H$\alpha$ equivalent widths are quite typical of rapidly rotating M stars.

In addition to the low resolution spectra from WiFeS, we also obtained a medium resolution spectrum of TIC 206544316 using the echelle spectrograph on the ANU 2.3 m telescope. The echelle yields a spectral resolution of $R=24000$ over the wavelength region $4200-6725$ Å, and the observations were reduced as per Zhou et al. (2014). The echelle spectrum of TIC206544316 revealed that it is rapidly rotating at a rate of $v_{\_}\mathrm{eq}\sin i\sim 77\,\mathrm{km\,s}^{-1}$, consistent with the equatorial rotational velocity of the star derived from the photometric period. As with the WiFeS spectra, no obvious signatures of spectral blending indicative of multiple stars in the system are apparent.

We find no evidence of ongoing accretion in these M-dwarf systems. The H-$\alpha$ emission strengths as measured from our spectra are in-line with that expected from rapidly rotating field M-dwarfs (Newton et al., 2016), and as such are not indicative of accretion. It is difficult to quantitatively constrain the extent of veiling in our low-resolution spectra. We do not have accurate measurements of the rotational broadening, and as such it is more difficult to match the low-resolution spectra against synthetic templates while simultaneously fitting for line dilutions due to veiling. We suggest future observations at high resolution that may be able to address this problem.

To further probe the spectra of our ten highly structured rotators, especially in the NIR, we utilized the observed spectral energy distribution (SED). We made use of photometric magnitudes from Gaia $G_{\_}\mathrm{BP}$ and $G_{\_}\mathrm{RP}$ (Brown et al., 2018), 2MASS $J$, $H$, $Ks$ bands (Cutri et al., 2003), and WISE $W1$, $W2$, $W3$, and $W4$ where available (Cutri et al., 2013) to describe the SED of all ten of our stars listed in Table LABEL:tbl:weird.

We compared the SEDs to a set of un-reddened optical and near-infrared spectra from the TW Hydra or $\beta$ Pic moving groups (10-25 Myr, Bell et al., 2015), following the method from (Mann et al., 2016), which we briefly summarize here. We restricted our comparison to templates with spectral types or $T_{\_}{\rm{eff}}$ estimates consistent with the stars in Table LABEL:tbl:weird. Spectral types and $T_{\_}{\rm{eff}}$ values for the SED-fitting templates were determined in an identical manner, and hence were on the same overall scale. Gaps in the coverage of each template were filled with PHOENIX BT-SETTL models (Allard et al., 2012). We assumed zero reddening, as all targets are within the local bubble (Sfeir et al., 1999).

After fitting all ten SEDs we reach several conclusions. (1) Four of the SEDs show no NIR excess, e.g., TIC 206544316 and TIC 425933644 (see top panel of Fig. 7) have W4 magnitudes, and both are in excellent agreement with the model+template. TIC 38820496 and TIC 177309964 are well fit without any need for an excess. (2) Three of the systems show a clear NIR excess, e.g., TIC 425937691 (see middle panel of Fig. 7) and TIC 234295610 have W3 magnitudes that are $\gtrsim 4\,\sigma$ above the model, even considering the range of possible model fits. TIC 201789285 is also a fairly clear case, with a possible W2 excess and a larger break at W3. (3) The remaining three systems fall in the ambiguous range, e.g., TIC 224283342 appears to show a W3 excess, but it is only barely significant. TIC 332517282’s W3 excess is not significant, just suggestive (see bottom panel of Fig. 7).

These results may be similar to what we would find if we were to randomly sample stars in these young stellar associations. But, the fact that some show clear excesses does lend credence to the notion that dust is involved. In particular, these excesses are consistent with what we would expect to see from a small dust disk near the star whose projected area is less than an order of magnitude larger than that of the host star, and having equilibrium grain temperatures of $\lesssim 500$ K (see Sect. 7).

5 Modeling the Complex Lightcurves with Spots

5.1 Fitting stellar spots with allesfitter

We studied whether the periodic variability of the 10 M stars may be explained by the presence of spots on the stellar surfaces. For this, we use the allesfitter package (Günther & Daylan, in prep.) to fit spot models of different complexities. Using Nested Sampling (Skilling, 2006), we explore the entire physical parameter space, and compute the Bayesian evidence for different models.

allesfitter444https://github.com/MNGuenther/allesfitter, online Jan 14, 2019 is a publicly available, user-friendly software package built on the packages ellc (light curve and radial velocity models; Maxted 2016); aflare (stellar flares; Davenport et al. 2014); dynesty (static and dynamic Nested Sampling555https://github.com/joshspeagle/dynesty, online Jan 14, 2019); emcee (Markov Chain Monte Carlo sampling; Foreman-Mackey et al. 2013); and celerite (GP models; Foreman-Mackey et al. 2017).

We note two intrinsic approximations in the following fits using allesfitter. First, broad-band photometric observations, such as those represented by TESS lightcurves, do not allow us to discern any wavelength-dependency of individual stellar surface features (e.g., warm plages vs. cool spots). Second, our approach describes a mean model, averaging over the lifetime and variability time-scales of individual features (which can range from days to weeks).

5.2 Verification of our model on HSS 348

We verify our ability to fit spot models with allesfitter using the example of HSS 348, which shows one of the more complex spot configurations modeled in the literature (Strassmeier et al., 2017). HSS 348 is a B8 + B8.5 double-lined spectroscopic binary system. For comparison, we adopt the authors’ assumption that all spots are on one star. We mask out the eclipse region, and phase-fold the CoRot data on the period reported in Strassmeier et al. (2017).

We ran two models, each containing four dark spots, on the phase-folded data. In the first model, the spin axis of the spotted star is fixed, according to the inclination of the binary system which is known from eclipse observations, to be normal to the binary orbital plane. The dilution is fixed to $50\%$ to account for the light of the other star. We fit for 18 free parameters and adopt the following priors: four parameters for each of the four dark spots (longitude $\mathcal{U}(0\arcdeg,360\arcdeg)$666A uniform probability distribution from $a$ to $b$ is denoted by $\mathcal{U}(a,b)$., latitude $\mathcal{U}(-90\arcdeg,90\arcdeg)$, size $\mathcal{U}(0\arcdeg,30\arcdeg)$, relative brightness $\mathcal{U}(0,1)$), a linear limb darkening term $\mathcal{U}(0,1)$, and the logarithm of a parameter for scaling the photometric error bars $\mathcal{N}(-8,1)$ truncated to $(-23,0)$777i.e. a normal distribution with mean of -8 and standard deviation of 1, which is truncated to the range from -23 to 0.. A small baseline offset taken to be the mean of the residuals is subtracted as part of each step in the sampling.

In the second model we keep the previous settings except that we relax the assumption of a fixed spin axis and fit for the additional parameter, the cosine (prior $\mathcal{U}(0,1)$) of the angle between the axis of the spotted star and the plane of the sky.

Both spot models provide an excellent fit to the variability patterns of HSS 348. We compute the Bayesian evidence for each model in order to compare them quantitatively. We find significant Bayesian evidence for the model with the free spin axis. The best fit of the phase-folded profile and spot configuration for this model are shown in Fig. 8.

A casual comparison of our distribution of spots on the stellar surface with that of Strassmeier et al. (2017) indicates some differences. This is not surprising, as those authors state that their spot solution is not mathematically unique. However, the difference lies basically in the fact that our modeling allows a different temperature for each spot while Strassmeier et al. (2017) assigned a temperature 1000 K below that of the stellar effective temperature to all four spots. A single spot in each quadrant of stellar longitude is sufficient to reproduce the light curve fairly well in both solutions, including the four observed minima. An advantage of using Nested Sampling is that it allows one to search for multiple local and global maxima in the likelihood function. The solution shown in Fig. 8 is in fact the only global likelihood maximum for our model.

5.3 Spots alone cannot describe the rapid rotators

Using models of 2, 5 or 8 spots, which can here be either bright or dark, we try to describe the complex periodic variations of the rapid rotators using the allesfitter code. We apply this to the light curve data extracted by the SPOC pipeline. We do not apply any detrending procedure other than correcting for the dilution of known background stars. Each light curve is phase-folded on the prominent period identified for the respective object. We fit for the following parameters with priors: four parameters for each of the 2, 5, or 8 spots (longitude $\mathcal{U}(0\arcdeg,360\arcdeg)$, latitude $\mathcal{U}(-90\arcdeg,90\arcdeg)$, size $\mathcal{U}(0\arcdeg,30\arcdeg)$, relative brightness $\mathcal{U}(0,2)$), the cosine of the spin axis inclination $\mathcal{U}(0,1)$, a linear limb darkening term $\mathcal{U}(0,1)$, and the logarithm of the parameter for scaling the photometric error bars $\mathcal{N}(-8,1)$ truncated to $(-23,0)$. At each step in the sampling, the mean value of the residuals is computed and subtracted to remove any small baseline offset.

Fig. 9 shows the results of the fit of eight spots to the rotational profile for TIC 201789285. The Bayesian evidence favors the 8-spot model over the 2-spot and 5-spot models. This is likely due to the fact that none of the models yields an entirely adequate fit - the eight spots merely give the most freedom to “bend” the model curve to fit the sharp, high-frequency features. We note that the spot configuration shown in Fig. 9 is locally non-physical; in the model, overlying spots lead to a local stellar surface area with negative flux.

In conclusion, even complex spot configurations cannot describe the sharp features found in the phase-folded light curves of these rapid rotators. While spots might contribute to the overall variability, additional physical effects close to the stars must be involved in the formation of the rotational profiles.

6 Magnetically Constrained Dust

Stauffer et al. (2017, 2018) have suggested that the highly structured rotational modulation patterns in their sample of Upper Sco and $\rho$ Oph stars are produced by “orbiting clouds of material at the Keplerian corotation radius”. Whatever produces the high frequency structure in the phase-folded profiles of those stars is likely to also apply to the stars we report on in this work. Stauffer et al. (2017, 2018) go on to say that “we attribute their flux dips as most probably arising from eclipses of warm coronal gas clouds”. In this section we re-examine the possibility that absorption by dust rather than scattering by warm gas may be responsible for the dips.

In order for dust to be the cause of the sharp rotational modulations, three requirements must be met: (1) there must be a sufficient amount of dust to account for the observed attenuation in flux; (2) the dust must be able to survive for at least days to weeks before sublimating; and (3) the dust must located where it will corotate with the star, perhaps by virtue of being electrically charged and effectively locked onto magnetic field lines of the star (see e.g., Farihi et al. 2017).

Given that these stars are typically 40 Myr old, it is not obvious how much dust might remain close to each of them. A dust sheet comprising \sim$$10^{17} g in submicron dust particles is required to block $\sim$5-10% of a star’s light.

Regarding the third point, the gyroradii of the charged dust particles must be much smaller than orbital distance of the grains, and no larger than the stellar radius. The radius of gyration is:

[TABLE]

where $m$ and $q$ are the mass and charge of the dust grain, $v_{\_}\perp$ is the velocity perpendicular to the magnetic field of strength $B$, and $c$ is the speed of light. If expressed as a ratio with respect to the stellar radius, the requirement becomes

[TABLE]

where, in the top equation, $a$ is the grain radius, $\rho$ the bulk density of the grains, $N$ the number of electron charges on each grain, $e$ the electron charge, $M_{\_}*$ the mass of the host star, $\xi$ is $v_{\_}\perp$ in units of the free-fall speed, $P_{\_}{\rm crit}$ the rotational breakup period of the host star, $B_{\_}s$ the surface magnetic field (assumed dipolar), and $d$ the distance of the grain from the host star. In the bottom equation the expression has been evaluated for plausible values, e.g., the grain size in microns, number of charges per grain in units of 1000, $B_{\_}s$ the surface field in units of 3 kG, $P_{\_}{\rm crit}$ in units of hours, $\xi$ normalized to 0.1, and $\rho$ is taken to be 3 g/cc. For the somewhat arbitrary normalizations used for these parameters, and for $d/R_{\_}*\sim 1$ (i.e., near the stellar surface), it is likely that the dust grains can be shepherded.

Regarding magnetic field strengths, the rapidly rotating and flaring M dwarf star V374 Peg is similar to the stars investigated in the present paper. Its rotational period is 0.4457 days = 10.70 hours, while its mass and radius are 0.30 $M_{\_}\odot$ and 0.34 $R_{\_}\odot$, all within the parameter ranges of the stars in Table LABEL:tbl:weird. V374 Peg has a low amplitude stable light variation with some minor changes over time, and has exhibited many flares as well as a possible coronal mass ejection (see Vida et al. 2016 for details). More importantly, this star may have a stable global poloidal magnetic field as modeled by Morin et al. (2008). The magnetic field of V374 Peg, in locations, reaches strengths of 3-4 kG, implying an overall 1.5-2 kG field of the star (see Donati et al. 2006). New measurements by Shulyak et al. (2017) reveal a 5.3$\pm$1 kG field on the star. Interestingly, the T-Tauri star MN Lupi appears to have nearly the same rotation period, a large magnetic field of several kG, and an active accretion disk connected to the star via magnetic field lines (Strassmeier et al., 2005).

We next assess whether dust grains could survive for substantial times, e.g., a few days, close enough to the M stars that we are considering in this work. We do this in two steps by (1) computing the equilibrium dust-grain temperatures, $T_{\_}{\rm eq}$, for unshielded grains in circular orbits at each of a range of distances, i.e., orbital radii, from 1 $R_{\_}*$ to 30 $R_{\_}*$, and (2) using $T_{\_}{\rm eq}$ and the bulk properties of the grain material to compute the sublimation time, $\tau_{\_}{\rm sub}$. We do this for grains of three different sizes (0.1, 1, and 10 $\mu$m) and 6 different compositions revolving around a canonical M star with $M_{\_}*=0.3\,M_{\_}\odot$, $R_{\_}*=0.4\,R_{\_}\odot$, and $T_{\_}{\rm eff}=3200$ K. In order to compute $T_{\_}{\rm eq}$ we solve the thermal balance equation by matching the emission of radiation at $T_{\_}{\rm eq}$ to the absorption of radiation from the host star at $T_{\_}{\rm eff}$ (we ignore the cooling effects of sublimation as this effect is rarely important). To solve this equation we utilize the grain absorption cross sections based on Mie scattering theory following Croll et al. (2014) and Xu et al. (2018). In turn, the Mie calculations require, as input, the real and imaginary parts of the indices of refraction ($n$ and $k$, respectively) of the dust grains (taken to be spherically shaped). We study 6 minerals that have been considered for the tails of the so-called “disintegrating exoplanets” (see, e.g., van Lieshout & Rappaport 2018 for a review). These include iron (Querry, 1985), corundum (Begemann et al., 1997), forsterite (Jäger et al., 2003), enstatite (Jäger et al., 2003), and fayalite888https://www.astro.uni-jena.de/Laboratory/Database/jpdoc/f-dbase.html, where the citations in parentheses are original sources of the corresponding $n$ and $k$ values. The sixth mineral is a generic “astronomical silicate” mixture (Draine, 1985). We obtained the indices of refraction from the HITRAN aerosol database (Gordon et al., 2017) and bulk mineral properties, such as the vapor pressure vs. $T_{\_}{\rm eq}$, from Table 5 of van Lieshout et al. (2016).

The results of this grain survival study are shown in Fig. 10. The top panel shows the equilibrium grain temperatures as a function of distance from the host star out to 20 $R_{\_}*$; the corresponding Keplerian orbital periods are labeled along the top axis. For distances of $\lesssim 5\,R_{\_}*$, the temperatures of the grains of iron, corundum, and fayalite (as well as cosmic silicates) are modestly high, while for forsterite and enstatite (the Mg-based silicates) the temperatures can be surprisingly low. This is due to the fact that these latter two minerals have low values of $k$ ($\sim$ 0.001) near the wavelengths they are absorbing and relatively high values of $k$ ($\sim$ 1) at the much longer emission wavelengths. Grain sublimation times for 15 combinations of grain sizes and compositions that we consider are shown as a function of distance in the bottom panel of Fig. 10. (Sublimation times for the astronomical silicate mixture are not given since there is no single set of mineral properties for this collection of materials.) We see that grains with sizes in a relevant range and composed of a number of kinds of minerals can survive for more than a day at distances of $\lesssim 5R_{\_}*$, while grains of 3 of 5 mineral types can survive at distances of $\lesssim 3R_{\_}*$; 2 of 5 can even survive relatively close to the stellar surface.

Perhaps the biggest difficulty with the dust-grain model, or any other model invoking obscuration by intervening material (including coronal gas), lies in the distance from the host star where the material must be located. In order for any attenuating material to produce narrow duty-cycle features in the rotational light curve it must be situated considerably above the stellar surface. If $f$ is the full width of the feature divided by the rotation period, the distance, $d$, from the stellar center, in units of the radius of the host star, $R_{\_}*$, must be

[TABLE]

If, for example, we consider narrow features with duty cycles, $f$, of 1-10% (see, e.g., Figs. 2, 3, and 4), then $d$ should be $\gtrsim 3-30\,R_{\_}*$. (This is basically the reason that starspots, with $d/R_{\_}*=1$ do not work.) Such large factors of $d/R_{\_}*$ in Eqn. (4), raised to the 5/2 power, increase the leading coefficient from 0.06 to values of $\sim$1 to 300 which quickly violates the requirement for meaningful constraint by the magnetic field.

7 A Dusty Ring Model

We now consider a model where the star is orbited by one or more rings composed of dust-size or somewhat larger particles. In this model, there would be no requirement for dust to be magnetically entrained in the stellar magnetic field at distances of many stellar radii where the magnetic field would be relatively weak as in the magnetically constrained dust model discussed in Sect. 6. The ring particles would move in Keplerian orbits at relatively large distances from the star, and therefore the sublimation lifetime would not be an issue even if the particles are dust-like in size. Furthermore, there would be no requirement for small coherent dust structures to form and be maintained.

A sketch of the envisioned geometry is shown in Fig. 11. In terms of radial structure, the disk could be either mostly continuous with some gaps or thin regions, or it could comprise a small number of concentric rings. The ring (or rings) would be viewed at a tilt angle with respect to the edge-on orientation, $\zeta$, that is small enough for the ring to partially occult the star. For a range of orientations of the star’s rotation axis such as that shown in the figure, cool or warm spots will be hidden by the tilted rings for short time intervals as the star rotates and thereby produce rapid photometric changes.

Consider first just a single ring extending from radius $r$ to $r+\Delta r$, or a gap in a more continuous ring of the same dimensions. We assume that the disk is thin, i.e., having a half-thickness $h$ that satisfies $h\ll\Delta r$. With this configuration, there would be a nearly linear feature across the star of extra transmission or attenuation having a thickness $\zeta\Delta r$ and a vertical displacement from the center of the star of $\zeta r$. A requirement for this gap to produce a sharp temporal feature as the star and its spots rotate beneath it, is that $\zeta\Delta r\ll R_{\_}*$. A final requirement is the existence of stellar surface features that can produce photometric modulations of order 1% by this mechanism. In turn, this implies the existence of regions with characteristic sizes $\mathcal{R}$ where $\mathcal{R}\sim\zeta\Delta r\sim 0.1R_{\_}*$, and with surface brightnesses that are sufficiently different from the mean stellar surface brightness to produce $\sim$1% effects.

Finally, in regard to a possible occulting ring, we can estimate a lower bound on its mass. The mass of a single uniform ring is given by:

[TABLE]

where $R_{\_}{\rm out}$ and $R_{\_}{\rm in}$ are the outer and inner radii of the ring, and $\Sigma$ is the mass column density of the particles comprising the ring in the direction normal to the ring plane. An illustrative set of parameters that could produce a sharp feature of width a few percent of the modulation period with $\sim$1% amplitude would be: $R_{\_}{\rm in}=10\,R_{\_}*$, $R_{\_}{\rm out}=15\,R_{\_}*$, $\mathcal{R}=0.2\,R_{\_}*$, and $\zeta=0.04\simeq 2^{\circ}$.

A ring consisting of micron-size dust particles with a projected number column density \sim$$1\,\mu{\rm m}^{-2} would involve a minimum total mass to provide the desired light attenuation. In that case, for grains whose bulk material density is $\sim$3 gm/cc, the mass of a grain is \sim$$10^{-11} g, and $\Sigma\approx 0.001\,\zeta\,{\rm g~{}cm}^{-2}$. Plugging this value for $\Sigma$, the above illustrative values of $R_{\_}{\rm in}$, $R_{\_}{\rm out}$, and $\zeta$, and $R_{\_}*\sim 0.3\,R_{\_}\odot$ into Eqn. (6), we find a mass of $\approx 10^{19}$ grams. The mass of the disk ring would be larger for larger particle sizes, $\Delta r$, or $\zeta$.

Finally, we estimate the precession time for an inclined disk with the above parameters. We take the apsidal motion constant, $k_{\_}2$, for the M stars to be that of an $n=3/2$ polytrope with a value of $k_{\_}2\simeq 0.14$ (Brooker & Olle, 1955). From that, we compute the $J_{\_}2$ quadrupole term to be in the range of 0.001 to 0.024 for stellar rotation periods in the range of 1.6 to 4 h. The resultant disk precession times are in the range of 10-70 yr, even for the shortest periods among the highly structured rotators.

8 Radiation Beaming

Beaming of the outgoing radiation is possibly another way of producing sharp features in rotational profiles. This could come about if radiation passing through magnetized plasma near the stellar surface has an anisotropic transmission factor.

Hamada & Kanno (1974) have computed how the Thomson scattering cross section is modified in the presence of a strong magnetic field. For an unpolarized beam, incident at angle $\theta$ with respect to the magnetic field direction, the total scattering cross section is (Hamada & Kanno, 1974):

[TABLE]

where $u$ is the ratio $\nu_{\_}C/\nu_{\_}{\rm rad}$, $\nu_{\_}C$ is the cyclotron frequency, $\nu_{\_}{\rm rad}$ is the frequency of the radiation, and $\sigma_{\_}{\rm Th}$ is the classical Thomson cross section. The cyclotron frequency is: $\nu_{\_}C\simeq 3\,B_{\_}{\rm kG}\,{\rm GHz}$ while the frequency of the radiation in the TESS band is $\approx 3\times 10^{14}$ Hz. Since $u$ is obviously quite small, we can expand Eqn. (7) in a Taylor series to show that:

[TABLE]

Therefore, the magnitude of the difference in the cross section between radiation incident in the direction of the magnetic field vs. perpendicular to it is of the order $u\sim 10$ ppm, and therefore appears to be unable to lead to significant beaming.

Barring some other type of physical effect that can yield beaming that can, in turn, produce a $\sim$5% amplitude with beaming angles smaller than $\sim$20∘, we do not see how beaming can produce the effects we see.

9 Large Flares and Changes to the Rotational Modulation Pattern

The large flares that are sometimes observed in these stars may provide clues to the origin of the complex rotational modulation patterns. In fact, Stauffer et al. (2017, 2018) found that the rotational modulation profile changed significantly after flares were observed in some of the systems in Upper Sco. This prompted us to examine the behavior of our 10 rotators with structured modulation profiles after the occurrences of flares. In the TESS light curves, large flares where the intensity was seen to increase by more than a factor of two were observed on only two of our 10 stars, and, of these, a distinct change in modulation profile was seen only in the light curve of TIC 206544316 after a large 4$\times$ superflare ($E_{\_}\mathrm{bol}=1.7\times 10^{34}$ erg; Günther et al., 2019).

The light curve of this star around the time of this event is shown in Fig. 12. In this figure the average 0.32-day complex modulation pattern is shown in red, while the photometric measurements are shown in black. After the 4$\times$ flare, the three subsequent rotational profiles have larger amplitudes and somewhat different shapes. Two additional smaller flares erupt later, one during the next profile that has a strongly enhanced maximum brightness, and another three cycles later, again at higher maximum light than before or well after the big flare.

About 25 rotations later than the flares shown in Fig. 12 another series of three or four flares with smaller peak brightnesses are observed within an interval of 5 rotations, and they occur at rotational phases similar to those seen in Figure 12.

The two smaller flares seen in Fig. 12 are likely ‘sympathetic flares’, that is, flares that sometimes appear, not merely by coincidence in time, soon after a large flare in analogy to such events that have been observed on the Sun (for the solar case see Török et al. (2011)). Such temporally-correlated solar flares may be produced subsequent to disturbances of the large-scale magnetic field of the Sun from certain flares, as particular types of disturbances of the field may give rise later to other events at other locations on the Sun (see Schrijver & Title (2011)). TIC 206544316 has parameters similar to those of V374 Peg that is inferred to have a strong and stable global magnetic field (see Sect. 6). Thus we can suppose there is a similar strong, global magnetic field on TIC 206544316 which, after being disturbed, could induce additional flaring events.

In light of the above, we interpret the disturbance to the light curve of TIC 206544316 in Fig. 12 in the following way. First, we have fit the large flare by the same model discussed in Günther et al. (2019) and find it to be consistent with the flare having decayed away by the time of the subsequent rotation of the host star. Second, the two prominent dips in the light curve are mostly unaffected by the flare for the next few rotations. However, another smaller feature (the ‘notch’), just after the maximum in the rotational profile, changes permanently after the flare. Third, we believe that the excess amplitude of the profiles in the three rotations following the flare may be attributed to the aftermath of the flare.

Aside from this one dramatic event in TIC 206544316 we find no other case where a flare appears to be associated with a distinct change in either the shape or amplitude of the highly structured rotational profile of any of the other 9 stars of our set.

In the dust model, a large flare could at least temporarily change the modulation profile of the star by, e.g., inducing changes in the configuration of the magnetic field, by increasing the rates at which grains sublimate, or by inducing changes in the amounts of net electrical charge on the grains. Likewise, in the model noted by Stauffer et al. (2017), which invokes corotating coronal gas as a scattering mechanism, a large flare may also be expected to result in temporary changes in the modulation profiles. Our observation of changes in the modulation profile of one star following one flare are not sufficient to give us insight into the physical mechanisms underlying the modulations, and certainly on its own does not confirm or contradict any of the models we have considered.

10 The Stellar Associations

Rapidly rotating M dwarfs are often young and, consequently, close to the stellar associations into which they were born. In such cases, the stellar ages may be estimated from the ages of the stellar associations. To try to figure out which stellar association each of our target stars might belong to, we used the online BANYAN $\Sigma$ Bayesian analysis tool relevant to nearby young associations (Gagné et al., 2018). If provided sufficient information about the 6-dimensional phase-space coordinates of a star (Table LABEL:tbl:weird), BANYAN $\Sigma$ assesses the likelihood that the star belongs to one of the 27 known and well-characterized young associations within 150 pc of the Sun.

Each of our 10 objects is most likely to be a member of the stellar association given in Table LABEL:tbl:weird. Five out of 10 belong to the association Tucana Horologium. For convenience, we give in Table 3 the names, distances, and ages of the five associations listed in Table LABEL:tbl:weird. These associations all have distances between 30 and 50 pc and the ages of four of the five are between 25 and 50 Myr.

These ages of $\sim$40 Myr are considerably older than the nominal 10 Myr age of the Upper Scorpius association (Pecaut et al., 2012; Feiden, 2016) which houses many of the known dippers as well as the structured rotators found by Stauffer et al. (2017, 2018). It is not clear whether it is reasonable or not to postulate that 40-Myr-old M stars can accrete sufficient amounts of dust into their magnetospheres to produce the observed rotational modulations.

11 Summary and Discussion

In this work we have reported the discovery, using the first two sectors of TESS data, of 10 new rapidly rotating M-dwarfs with highly structured rotational modulation patterns. The overall shapes of the periodic patterns can remain stable for weeks, while individual features can change within days. As part of this study, we also report the discovery of 17 rapidly rotating M-dwarfs with periods $<4$ hr, of which the shortest period is 1.63 hr. These objects are part of a sample of 371 rapidly rotating M-dwarfs also found in the first two sectors of TESS data with $T_{\_}{\rm eff}$ cooler than 4000 K and with rotation periods $<$ 1 day.

We attempt to better understand the nature of the rapidly rotating M-dwarfs with highly structured rotational modulation patterns, which are similar in their properties to a number of stars found in the Upper Sco stellar association (Stauffer et al., 2017). First, we fitted general star-spot models, able to accommodate large numbers of spots, to these periodic variations which exhibit large amounts of structure (Sect. 5). However, we show that spots are not likely to be able to explain the sharp features in the light curves (see also Stauffer et al. 2017). Next we explored the possibility that the highly structured rotational modulation patterns are due to absorption by charged dust particles that are trapped by the stellar magnetosphere (see Sect. 6). We have shown that dust grains may be able to survive for a number of rotation periods in the magnetosphere at the corotation radius. However, such dust would have to be substantially far from the host star (see Eqn. 5) to produce the sharp features observed in the rotation profiles, and the existing magnetic fields at those distances would likely be too weak to trap charged dust particles (see Sect. 6). Beaming of the emitted radiation would be a desirable mechanism, but we show, at least in the case of electron scattering in a magnetic field, that the observed amplitudes of a few percent are not achievable (Sect. 8). We did not explore the scenario, proposed by Stauffer et al. (2017), involving attenuation by scattering from warm coronal material, because it is difficult to conceive of a situation where it would produce significant emission or absorption.

Perhaps the most natural of the scenarios we have explored is that of a spotted host star rotating obliquely inside of a dusty ring of material a few tens of stellar radii away (Sect. 7). While potentially very promising, this model has not been fully fleshed out, and a determination of whether it can quantitatively explain the sharp modulation structures that are observed will require further exploration. However, this is beyond the scope of the present work.

In addition, we explore whether large flares could serve as probes of the highly structured rotational modulation patterns. We found 32 large flaring events with flux increases larger than 2 times the stellar flux. However, we do not see any changes in the intrinsic rotational modulation pattern after most flares. The exception is a large $4\times$ superflare event ($E_{\_}{\rm bol}=1.7\times 10^{34}$ erg) on TIC 206544316, which is the strongest flare observed to date in such a rapidly-rotating M dwarf with highly structured rotational modulation patterns. This superflare apparently changes the subsequent periodic modulation by increasing the maximum brightness for a few rotations, and altering one feature permanently.

Because most of our 10 highly structured rotators are in well-studied stellar associations, their ages are known reasonably well. Therefore, we are able to compare some of the properties of our stars with the ‘scalloped shell’ rotators in Upper Sco and $\rho$ Oph. The latter have ages of 5 - 10 Myr, distances of about 135 pc, and a median $K_{\_}s$ magnitude of 10.5 (median $M_{\_}K\simeq 4.9$). The comparable rotators in this work are some three times closer, and several times older. Both sets of stars have a similar median $K_{\_}s$ magnitude, but a median $M_{\_}K\simeq 6.4$, which is perhaps not surprisingly 1.5 magnitudes less luminous than the ‘scalloped shell’ modulators in upper Sco and $\rho$ Oph which are younger.

A small sample of low mass stars revealing periodic optical variability has also been found to be radio active (see e.g. Harding et al., 2013) – however, we note that their optical modulations seem not as sharp-featured as those of our rapid rotators. For these radio active low mass stars, electron cyclotron maser (ECM) emission is suggested to cause periodic, circularly polarized radio pulses (see e.g. Hallinan et al., 2007). Notably, the ultracool dwarf TVLM 512-46546 shows sharp radio bursts every 1.96 hours (Hallinan et al., 2007) and smooth optical modulations with a similar $\sim$2 h period (Littlefair et al., 2008; Harding et al., 2013). The latter was attributed to inhomogeneous dust clouds in the relatively cool stellar atmosphere (Littlefair et al., 2008). This short period is remarkably similar to the periods of our 10 complex rotators, and suggests that the optical and radio emission are linked by common magnetic phenomena. The modulation signal of TVLM 512-46546, for example, is stable over a time span of five years. Neither spots nor magnetically constrained dust could explain such long-term stability. Whether the modulations of our targets are as long-lived may be studied in the future.

The mechanism that produces the observed highly structured rotational modulation patterns has not yet been confidently identified. In that regard, we have discussed, in Section 7, a promising idea. Continued observation of these 10 complex rotators will also provide opportunities to further understand the relevant physical processes. Ideally, this will encompass observations of time variability in the radio band, in both optical and near-infrared bandpasses, and in H$\alpha$, spectroscopy well beyond that presented in Section 4), and direct imaging. Finding $\sim$100 more of these objects in the upcoming TESS sectors, as we anticipate, will hopefully further elucidate their nature.

This work was carried out as part of a larger search for rapid rotation and high frequency oscillations in the TESS data. The methods described herein can also be applied to other high frequency stellar phenomena.

The authors are grateful to the reviewer, Klaus Strassmeier, for a number of important suggestions which improved the manuscript. We acknowledge the use of public data collected by the TESS mission, which are from pipelines at the TESS Science Office and at the TESS Science Processing Operations Center and made publicly available from the Mikulski Archive for Space Telescopes (MAST). Funding for the TESS mission is provided by NASA’s Science Mission directorate. Resources supporting this work were provided by the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center for the production of the SPOC data products. M. N. G. acknowledges support from MIT’s Kavli Institute as a Torres postdoctoral fellow. We are grateful to Lidia van Driel-Gesztelyi for input about solar analogs to some of the observed features of the systems studied in this work. Software: python (Rossum G., 1995), numpy (van der Walt et al., 2011), scipy (Jones et al., 2001), matplotlib (Hunter et al., 2007), tqdm (doi:10.5281/zenodo.1468033), seaborn (https://seaborn.pydata.org/index.html), allesfitter (Günther & Daylan, in prep.), ellc (Maxted, 2016), aflare (Davenport et al., 2014), dynesty (https://github.com/joshspeagle/dynesty), corner (Foreman-Mackey, 2016).