Pulsed Accretion in the T Tauri Binary TWA 3A

The LCO 1 m network comprises nine 1 m telescopes located across four global sites: Siding Springs Observatory (Australia), SAAO (South Africa), CTIO (Chile), and McDonald Observatory (USA). Spanning ∼12 TWA 3A orbital periods, our observations were made between 2014 May and 2016 April. Observing visits were scheduled 20 times per orbit while the target was visible (airmass < 2), corresponding to a cadence of 42 hr. Each visit consisted of three images in the UBVR filters. All data are reduced by the LCO pipeline applying bad-pixel, bias, dark, and flat-field corrections. The three images per filter per visit are than aligned, combined, and fit with astrometric solutions using standard IRAF tasks.

The SMARTS 1.3 m telescope at CTIO is outfitted with the ANDICAM detector. Our program requested every-other-night visits of TWA 3A while it was visible (airmass < 2) between 2014 December and 2016 July. Each visit consisted of three images in the B and V filters. Data are reduced with the SMARTS pipeline, which applies bias and flat-field corrections. Each set of images per visit are aligned, combined, and fit with an astrometric solution using standard IRAF tasks.

For each telescope (LCO; SMARTS) and filter set (UBVR; BV), SExtractor (Bertin & Arnouts 1996) is used to perform automated source detection and photometry on each image producing time-series instrumental magnitudes. A source catalog for each data set is then created from spatial matching of the astrometric solutions. Due to poor seeing and/or telescope focus, flux from TWA 3A and TWA 3B could not be consistently separated. As such, SExtractor parameters were optimized to photometer the entire TWA 3 system.

Each source catalog is then fed into an ensemble photometry routine following the Honeycutt (1992) formalism. By selecting non-varying comparison stars interactively, variations from airmass and nightly observing conditions are corrected, resulting in relative-magnitude light curves for each star in the field of view (FOV). We note that this correction does not include color information, which leaves some small systematic error, especially in the U-band where atmospheric corrections are most color dependent.

Relative ensemble magnitudes are then transformed to apparent magnitudes using non-varying stars in the LCO FOV for which published empirical or derived photometry exists. Five such stars are present in the LCO FOV (TYC 7213-797-1, 7213-391-1, 7213-1239-1, 7213-933-1, and 7213-829-1). Their V magnitudes range from 11.18 to 11.57 and (B – V) colors span 0.71 to 2.00. Using empirical measurements where available, we draw B and V values from the All-Sky Compiled Catalogue (Kharchenko 2001). For R-band, we use a colorless transformation from the Carlsberg Meridian Catalog 15 (CMC15; Niels Bohr Institute et al. 2014) ${r}^{\prime }$ to R ($R={r}^{\prime }-0.22\pm 0.12$) derived from 3690 overlapping stars between CMC15 and the Lowell Observatory Near-Earth Object Search (Skiff 2007). Without empirical U-band measurements available, we use the fitted apparent U-band magnitudes derived by Pickles & Depagne (2010). From these five stars we compute zero-point and color transformations that are applied to the rest of the field. Systematic errors associated with this calibration procedure are determined from the root-mean-squared deviation between the transformed published values of the calibrating stars in color–magnitude space. They are 0.18, 0.22, 0.25, and 0.07 mag for U, B, V, and R bands, respectively. (Systematic errors are propagated through the mass accretion rate derivation that follows and are presented in Figure 3.)

Not all five of the calibration stars in the LCO FOV are present in the smaller 6' FOV of the ANDICAM CCD. To transform the SMAR TS data to apparent magnitudes, we bootstrap the apparent magnitudes derived for non-varying stars in the LCO FOV that overlap with the SMARTS FOV and use those to determine zero-point and color transformations.

Figure 1 presents the U-, B-, V-, and R-band light curves for TWA 3 plotted against an arbitrary orbital cycle number set to 1 for the first observed periastron passage. In each panel, vertical dashed lines mark the TWA 3A periastron passage and horizontal dotted lines mark the quiescent value (average of orbital phases 0.2 to 0.4). Brightening events near periastron passages are seen consistently, having the largest increase in the U-band. These events very closely match the accretion behavior predicted for eccentric binaries.

Zoom In Zoom Out Reset image size

To ensure variability observed in the TWA 3 system is indeed from the disk-bearing binary and not the tertiary, we perform point-spread-function photometry on a subset of the SMARTS B-band images where the light from each component can be reliably separated. Figure 2 displays the result where TWA 3A is the clear source of variability. The standard deviation of these light curves are 0.14 and 0.02 mag for TWA 3A and TWA 3B, respectively. In the following, we assume all variability in the TWA 3 system results from the spectroscopic binary.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Before assigning all optical variability to changes in the TWA 3A accretion rate, we inspect our light curves for contributions from stellar flares (magnetic reconnection events at stellar surfaces). In high-cadence photometry of DQ Tau, stellar flares with amplitudes greater than ${\rm{\Delta }}U=0.32$ mag were found to have a temporal contribution of ∼3% (Tofflemire et al. 2017). Assuming the same contribution in TWA 3A corresponds to six measurements. In moderate-cadence photometry, however, these events would likely go undetected having similar amplitudes and colors as accretion variability. Only if a measurement were to contain a large flare or a flare peak, where the photometric color is typically bluer than accretion, would it stand out from the underlying accretion variability (Kowalski et al. 2016).

One measurement in our light curves has a color, magnitude, and temporal behavior that suggest it contains a stellar flare. Occurring at orbital cycle ∼1.2 in Figure 1, it is the brightest U-band measurement by 0.7 mag and the bluest in (U – B) by 0.2 mag, well separated from both observed distributions. The associated flux would correspond to a factor of ∼15 increase in the mass accretion rate compared to other measurements at similar orbital phases and a factor of two greater than the next highest measurement. A measurement only 5 hr later, however, falls securely within the remaining spread for that orbital phase. Accretion events of this scale are expected to be rare and to evolve over much longer timescales (e.g., Cody et al. 2014). And critically, the three U-band images combined in this measurement show a rapid ∼0.2 mag decline over ∼3 minutes. Given these characteristics, we conclude this measurement contains a stellar flare and remove it. In the following analysis, we assume the remaining variability is due to changes in the TWA 3A accretion rate.

Flux-calibrated U-band photometry can be used to derive a mass accretion rate with knowledge of the distance and extinction to the source, the photospheric U-band flux in the absence of accretion, and the stellar mass and radius. With values for the distance and stellar parameters, we determine the extinction and photospheric properties following Herczeg & Hillenbrand (2014). First, we fit a spatially resolved Keck LRIS spectrum of TWA 3A (Herczeg et al. 2009) with a library of empirical weak-lined T Tauri star (WTTS) spectra. The spectra are fit with three free parameters: a flux normalization, an additive accretion spectrum, and the extinction. Our results are consistent with those in Herczeg & Hillenbrand (2014), namely, A_V = 0.04 (±0.30 mag) and a combined TWA 3A spectral type of M4.1 (±0.3 subclasses). Second, the best-fitting WTTS spectrum is convolved with a U-band filter to determine the underlying photospheric contribution of the binary.

Because our U-band measurements are for the entire TWA 3 system, we use a spatially resolved Keck LRIS spectrum of TWA 3B (Herczeg et al. 2009) to determine its U-band contribution. Assuming a distance of 30 pc, we extinction correct the U-band measurements and convert them to U-band luminosities. Subtracting the contribution from TWA 3B and the underlying TWA 3A photosphere, we arrive at the TWA 3A U-band accretion luminosity. Using the model-dependent, empirical relation derived in Gullbring et al. (1998), we calculate the total accretion luminosity from the U-band as follows:

$$\begin{eqnarray}&&\mathrm{log}({L}_{\mathrm{Acc}}/{L}_{\odot })=1.09\ \mathrm{log}({L}_{{U}_{\mathrm{excess}}}/{L}_{\odot })+0.98.\end{eqnarray} \tag{ 1 }$$

For a single star, the mass accretion rate can be determined from the accretion luminosity with the following:

$$\begin{eqnarray}&&\dot{M}\simeq \displaystyle \frac{{L}_{\mathrm{Acc}}{R}_{\star }}{{{GM}}_{\star }}{\left(1-\displaystyle \frac{{R}_{\star }}{{R}_{\mathrm{in}}}\right)}^{-1},\end{eqnarray} \tag{ 2 }$$

where ${R}_{\mathrm{in}}$ is the magnetospheric disk truncation radius from which material free falls along magnetic field lines (typically assumed to be 5${R}_{\star };$ e.g., Johnstone et al. 2014).

Our measurements, however, are of the combined accretion luminosity from two stars with different masses and radii. Without a theoretical consensus for which star should predominantly receive the mass from circumbinary accretion flows, we assume that each star accretes at the same rate. The accretion luminosity emitted from the primary star alone, ${L}_{\mathrm{Acc},1}$, becomes

$$\begin{eqnarray}&&{L}_{\mathrm{Acc},1}={L}_{\mathrm{Acc},\mathrm{Total}}{\left(1+q\displaystyle \frac{{R}_{\star ,1}}{{R}_{\star ,2}}\right)}^{-1},\end{eqnarray} \tag{ 3 }$$

where q is the mass ratio. The total mass accretion rate for the binary is then

$$\begin{eqnarray}&&\dot{M}\simeq 2\displaystyle \frac{{L}_{\mathrm{Acc},1}{R}_{\star ,1}}{{{GM}}_{\star ,1}}{\left(1-\displaystyle \frac{{R}_{\star ,1}}{{R}_{\mathrm{in}}}\right)}^{-1}.\end{eqnarray} \tag{ 4 }$$

For the stellar radii, we use the Dotter et al. (2008) stellar evolution models to compute the average radii for 0.6 and 0.5 ${M}_{\odot }$ stars between 5 and 10 Myr. They are 1.06 and 0.99 ${R}_{\odot }$, respectively. ${R}_{\mathrm{in}}$ is set to the canonical single-star value of $5{R}_{\star }$. The derived accretion rates range between 0.8 × 10⁻¹¹ and 2.4 × 10⁻¹⁰ ${M}_{\odot }$ yr⁻¹, in good agreement with previous measurements (Muzerolle et al. 2000; Herczeg et al. 2009). Since the stars have near equal masses, the choice to split the accretion rate equally between the two corresponds to only a ±∼5% difference from assigning all the accretion to one star.

The top panel of Figure 3 presents the mass accretion rate phase-folded about the orbital period. The repeated enhanced accretion events observed near periastron increase the accretion rate by a factor of ∼3–10 from the quiescent value.

To determine the significance of the periodic accretion behavior, we perform a Lomb–Scargle periodogram (Scargle 1982) on the accretion rate measurements. The bottom panel of Figure 3 presents the power spectrum. A significant peak is observed above the 99% false-alarm probability with a period of 34.67 ± 0.14 days, in good agreement (1.5σ) with the binary orbital period. (Error in the accretion rate period is derived from a 10⁶ iteration Monte Carlo bootstrap simulation using random sampling with replacement of the $\dot{M}$ and HJD measurement pairs; Press et al. 1992.) Additional higher-frequency peaks occurring at two and three times the peak frequency result from the varying and non-sinusoidal morphology of the enhanced accretion events. The small peak to the left of the primary does not remain significant after filtering the data at the primary frequency.

The morphology of enhanced accretion events contains information on the interaction between the binary orbit and the circumbinary mass flows. In order to compare our observations with numerical simulations, we create an average accretion rate profile as a function of orbital phase. Breaking the orbital-phase-folded data into bins of $\phi =0.05$ (our sampling rate), we calculate the median mass accretion rate for each bin and set the bin error as its standard deviation. The result is presented in the top panel of Figure 4 where, on average, the accretion rate is elevated between orbital phases 0.85 and 1.05, reaching a peak near periastron of ∼4 times the average quiescent value.

Zoom In Zoom Out Reset image size

In the bottom panel of Figure 4, we compare the TWA 3A accretion rate profile to a simulation of binary accretion (Muñoz & Lai 2016) and to the DQ Tau accretion rate profile (Tofflemire et al. 2017). Both have been normalized to the TWA 3A average mass accretion rate. The Muñoz & Lai (2016) simulation shown is a 2D hydrodynamical model using the adaptive-mesh-refinement code AREPO (Springel 2010) for an equal-mass binary with an eccentricity of 0.5. Ten consecutive orbital periods of the simulation were used to create the accretion profile.