Received 2008 August 27; accepted by ApJL 2009 February 23 Preprint typeset using LATEX style emulateapj v. 03/07/07 HOW VARIABLE IS ACCRETION IN YOUNG STARS? Duy Cuong Nguyen1 , Alexander Scholz2 , Marten H. van Kerkwijk1 , Ray Jayawardhana1 , Alexis Brandeker3 arXiv:0902.4235v1 [astro-ph.SR] 25 Feb 2009 Received 2008 August 27; accepted by ApJL 2009 February 23 ABSTRACT We analyze the variability in accretion-related emission lines for 40 Classical T Tauri stars to probe the extent of accretion variations in young stellar objects. Our analysis is based on multi-epoch highresolution spectra for young stars in Taurus-Auriga and Chamaeleon I. For all stars, we obtain typically four spectra, covering timescales from hours to months. As proxies for the accretion rate, we use the Hα 10% width and the Ca II-λ8662 line flux. We find that while the two quantities are correlated, their variability amplitude is not. Converted to accretion rates, the Ca II fluxes indicate typical accretion rate changes of 0.35 dex, with 32% exceeding 0.5 dex, while Hα 10% width suggests changes of 0.65 dex, with 66% exceeding 0.5 dex. We conclude that Ca II fluxes are a more robust quantitative indicator of accretion than Hα 10% width, and that intrinsic accretion rate changes typically do not exceed 0.5 dex on timescales of days to months. The maximum extent of the variability is reached after a few days, suggesting that rotation is the dominant cause of variability. We see a decline of the inferred accretion rates towards later spectral types, reflecting the Ṁ vs. M relationship. There is a gap between accretors and non-accretors, pointing to a rapid shutdown of accretion. We conclude that the ∼ 2 orders of magnitude scatter in the Ṁ vs. M relationship is dominated by object-to-object scatter instead of intrinsic source variability. Subject headings: accretion, accretion disks — circumstellar matter — stars: formation — stars: low-mass, brown dwarfs — planetary systems 1. INTRODUCTION Young low-mass stars that are still accreting while contracting towards the main-sequence are identified observationally with T Tauri stars (e.g. Bertout 1989; Appenzeller & Mundt 1989). One of their defining properties is variability (e.g. Knott 1891; Joy 1945; Rydgren et al. 1976), which can be traced to changes in continuum excess flux caused by varying accretion flow from a circumstellar disk (Hartigan et al. 1991); other possible sources of variability include cool spot rotation and extinction events. In addition to the continuum, this variability is seen in a number of emission lines, including the Hydrogen Balmer series, and has been attributed to rotation, variable accretion, and winds (e.g. Johns & Basri 1995; Batalha et al. 2001; Alencar & Batalha 2002; Alencar et al. 2005). Photometric monitoring campaigns for large numbers of objects have revealed that hot spots formed by gas accretion, co-rotating with the objects, are one of the most important sources of variability in young stars (e.g. Herbst et al. 1994; Bouvier et al. 1995; Fernandez & Eiroa 1996). However, it remains unclear to what extent the photometric and spectroscopic variability seen in typical T Tauri stars is directly related to changes in the accretion rate. As we argue below, this issue can be addressed by spectroscopic monitoring, but so far spectroscopic studies have been limited to small Electronic address: nguyen@astro.utoronto.ca 1 Department of Astronomy & Astrophysics, University of Toronto, 50 St. George Street, Toronto, ON M5S 3H4, Canada; nguyen, mhvk, rayjay@astro.utoronto.ca 2 SUPA, School of Physics & Astronomy, University of St. Andrews, North Haugh, St. Andrews, KY16 9SS, United Kingdom; as110@st-andrews.ac.uk 3 Department of Astronomy, Stockholm Observatory, SE-106 91 Stockholm, Sweden; alexis@astro.su.se samples. Inferences about accretion are usually made in the framework of the magnetospheric accretion scenario (see the review by Bouvier et al. 2007). Recent observations have shown the presence of a correlation between inferred mass accretion rate and central object mass, extending over several orders of magnitude: Ṁ ∝ M α , with α ∼ 2 (e.g. Muzerolle et al. 2003; Natta et al. 2004; Mohanty et al. 2005). Apart from this correlation, all Ṁ vs. M plots feature a large scatter: at any given object mass, the accretion rates show a dispersion of about two orders of magnitude. A number of ideas have been put forward to interpret these findings. The accretion rate vs. mass correlation has been attributed to a BondiHoyle flow to the star-disk system (Padoan et al. 2005), to initial rotational velocities of collapsing cores (Dullemond et al. 2006), to a dispersion in disk parameters (Alexander & Armitage 2006), to a complex magnetic field geometry (Gregory et al. 2006), and to a declining disk ionization with stellar mass (Muzerolle et al. 2003; Mohanty et al. 2005). Alternatively, it has been argued that the correlation is not physical, but rather reflects selection effects; therefore, the dominant feature is the scatter itself (Clarke & Pringle 2006). All these scenarios predict specific values for α (0–2) and specific properties of the scatter. A basic question about the Ṁ vs. M relation is whether scatter around it is due to object-to-object variations (i.e. the accretion rate for any individual object remains mostly constant) or due to a variable accretion rate in any given object (Scholz & Jayawardhana 2006). At least for some objects, strong changes in accretion related lines have been reported, e.g. for the brown dwarf 2MASSW J1207334-393254 (Scholz et al. 2005; 2 Nguyen et al. Stelzer et al. 2007), but it is not clear if such variations are common or not. Since most previous studies were based on single-epoch measurements of accretion rates, it was impossible thus far to address this issue. Previous work has established that a number of optical emission lines originate in the accretion flow and are affected by the accretion rate. In particular, the Hα 10% width and the Ca II-λ8662 line flux are found to correlate well with the accretion rate (Muzerolle et al. 1998; Natta et al. 2004; Mohanty et al. 2005). These empirical indicators facilitate studies of accretion, and for the first time, allow us to investigate accretion in multiepoch spectra for large samples of objects. Here, we use both diagnostics to probe the intrinsic accretion variability in young stars in Taurus-Auriga and Chamaeleon I, to provide new observational limits for the aforementioned scenarios. 2. OBSERVATIONS AND DATA ANALYSIS We obtained multi-epoch high-resolution optical spectra of 40 members in the ∼ 2 Myr old Chamaeleon I and Taurus-Auriga star forming regions. The targets span the spectral types from F2 to M5 based on published classifications, and consist of accretors without suspected close companions selected from the list of Nguyen et al. (2009). The data were collected using the echelle spectrograph MIKE (Bernstein et al. 2003) on the Magellan Clay 6.5 meter telescope at the Las Campanas Observatory, Chile on 15 nights during four observing runs between 2006 February and 2006 December. Each target was observed typically at four epochs with baselines of hours, days, and months (only one target was observed just twice; 24, 12 and 3 targets were observed four, five and six times, respectively). For the data reduction, we used customized routines running in the ESO-MIDAS environment (described in detail in Brandeker et al., in preparation). The Hα 10% widths were determined as follows. First, we estimated the continuum level by linearly interpolating flux values in the range of 500 km s−1 to 1000 km s−1 on either side of the line. Next, the maximum flux level of Hα emission was measured with respect to this continuum level. Finally, the crossing points of the Hα emission with the 10% flux level were identified, and the width was measured. Sometimes absorption components in the Hα emission line profile falls below the 10% flux level. In such cases, we ensured that we measured the width consistently for all epochs, i.e. we measured widths to the edges of any blue- or red-shifted absorption features. We did not correct for underlying photospheric absorption in Hα. We derived the Ca II-λ8662 emission fluxes (FCa II ) from the observed emission equivalent widths. To determine the widths, we integrated the emission above the continuum level. For emission profiles attenuated by a broad absorption feature, we used the median flux within 0.2 Å of the absorption minima as an approximate continuum level for integration, similar to what was done by Muzerolle et al. (1998). To infer the emission fluxes from the equivalent widths, we must know the underlying photospheric continuum flux. We used the continuum flux predicted by the PHOENIX synthetic spectra for a specified Teff and surface gravity. We inferred Teff from our spectral types, and assumed a surface gravity of log g = 4.0 (cgs units). In our estimate, we ignore veiling, which could lead to an underestimate of the line flux. Indeed, for five targets shared by Mohanty et al. (2005), our results were lower by 0.05 to 0.41 dex. We tried to measure veiling from our spectra, but we found the relatively poor S/N prevented us from reaching sufficient accuracy (S/N of our spectra was typically ∼ 20 at Hα, whereas literature studies for veiling, e.g. Hartigan et al. (1989), use spectra with S/N & 75). However, for our purposes of studying variability, the bias due to veiling is not important. 3. ACCRETION INDICATORS AND THEIR VARIABILITY Previous studies have shown that the Hα 10% width and the Ca II-λ8662 flux are correlated with accretion rates determined using the traditional method based on optical veiling (Natta et al. 2004; Mohanty et al. 2005), with Ca II showing significantly less scatter around the fit (Herczeg & Hillenbrand 2008). In Fig. 1, we compare these two accretion indicators for our sample, with the points set to the average values from the multi-epoch spectra, and the ‘error bars’ indicating the range of values (the ranges are dominated by variability, not by measurement uncertainty). We find a clear linear correlation between the Hα 10% width and log flux in Ca II-λ8662, with the deviations comparable to the scatter. This provides reassurance that both parameters are mainly determined by the same physical quantity, which we identify as the mass accretion rate based on published findings. For individual objects, the degree of correlation between the two accretion indicators varies, which may be explained by the varying extent to which accretion affects different emission lines. We do not see a significant correlation between the variability in the accretion indicators and stellar mass, Spitzer IR excess, or star forming region in our sample. In Fig. 2, we compare the variability in the two accretion indicators. Here we use (max − min) for the Hα 10% width and (max/min) for the Ca II-λ8662 flux, because these quantities are directly proportional to accretion rate changes (Natta et al. 2004; Mohanty et al. 2005). One sees that the data are not correlated (correlation coefficient r2 ∼ 0.02), indicating that at least one of these parameters does not reflect changes in the accretion rate. Using the correlations log Ṁ = 1.06 log FCa II − 15.40 of Mohanty et al. (2005), and log Ṁ = −12.89 + 0.0097 Hα 10% [km s−1 ] of Natta et al. (2004), we convert the variability seen in the spectroscopic indicators to accretion rate changes. Based on Ca II-λ8662 fluxes, the median accretion rate variability in our spectra is 0.35 dex, with 13/40 or 32% exceeding 0.5 dex and only 1/40 or 2.5% exceeding one order of magnitude. These numbers should be treated as upper limits, as the uncertainties in measuring Ca II fluxes likely contribute ∼ 0.1 dex to the scatter. In contrast, using the Hα 10% width gives a median variability of 0.65 dex, with the majority exceeding 0.5 dex (26/40 or 65%), and still a substantial fraction exceeding one order of magnitude (16/40 or 40%). This result indicates clearly that the Ca II-λ8662 flux is a more robust quantitative diagnostic of accretion rate. Our empirical result is consistent with the models of Classical T Tauri stars of Azevedo et al. (2006), which Accretion Variability Fig. 1.— Ca II-λ8662 fluxes vs. Hα 10% widths in our sample: shown are the averages from our multi-epoch spectra with the error bars indicating the minimum-maximum range. From linear regression analysis, the probability that the data are not correlated is < 10−6 ; the best linear fit is overdrawn. 3 line fluxes or equivalent widths, the 10% width also has the advantage of not being heavily affected by uncertainties in estimating the underlying continuum. This property is beneficial especially for faint objects: the typical variation in 10% width from continuum uncertainty is ∼ 15 km s−1 . Based on our findings, however, the 10% width should not be trusted for quantitative measurements. Nevertheless, it is still useful as a qualitative indicator of accretion. A major problem in the 10% width measurement is its strong dependence on the line profile. The majority of the accreting stars in our sample show absorption components in their Hα profile, caused by relatively cool gas, e.g. in stellar or disk winds, or parts of the accretion flow, seen in projection against the hot shock front. Since these absorption features are variable in intensity and in position relative to the emission, the maximum intensity and therefore the 10% level can vary considerably from epoch to epoch. This effect can cause large changes in 10% width not necessarily associated with accretion rate changes. To probe the timescales of variability in the two line indicators, we used our time series baselines of hours to months: we examined indicator changes between all combinations of epochs for each object, e.g. if an object was observed at epochs A, B & C, then we reviewed changes between epochs A-B, A-C, and B-C. We show the changes graphically in Fig. 3. Both for the Hα 10% width and for the Ca II-λ8662 flux, we obtain consistent results: the extent of the variability increases on timescales ranging from hours to several days. On longer timescales, the amount of variability saturates. Thus, the dominant timescale for accretion-related variability is several days. This idea suggests the variability is determined near the star, where the timescales are sufficiently short. Probably, one major factor is rotation, since typical rotation periods in young stars are in the range of 1–10 d (Herbst et al. 2007). In addition, there may also be a connection to the characteristic infall timescale, typically hours, or magnetic field reconnection events, perhaps several days. Note that we cannot probe long-term variations on timescales of years with our dataset. However, the photometric results of Grankin et al. (2007), with observations secured over more than 20 years, find most Classical T Tauri stars to be fairly stable on long timescales. 4. THE ORIGIN OF THE SCATTER IN ACCRETION RATE VS. MASS Fig. 2.— Variability in spectroscopic diagnostics of accretion measured in our multi-epoch spectra for individual objects: (max/min) in Ca II-λ8662 flux vs. (max − min) of the 10% width. Both quantities are supposed to be directly related to accretion rate changes. We do not find a significant correlation in this dataset (correlation coefficient r 2 ∼ 0.02). find that Ca II broad components are formed predominantly in the accretion flow, and track accretion rate, while the hydrogen lines are much more affected by stellar winds. The main reason to use Hα for accretion rate measurements is that it is much easier to measure, and is sensitive to lower levels of accretion. In comparison with Our new constraints on accretion variability in a large sample allows us to constrain the origin of the distribution of data points in the accretion rate vs. mass diagram. In Fig. 4, we show the two accretion indicators used in this study as a function of spectral type. For the Hα 10% width, we also include the non-accreting objects. In both cases, a clear trend is seen: the spectroscopic indicator drops towards later spectral types, which reflects the Ṁ vs. M correlation reported in previous studies (see §1). As can be seen in the Hα 10% widths, the opposite trend is seen for the non-accretors, which is likely due to increasing levels of magnetic activity towards later spectral types. This leads to a ‘U-shaped’ distribution of data points: accretors and non-accretors are separated by a large gap at F–K spectral types, but are harder to 4 Nguyen et al. Fig. 3.— Change in Ca II-λ8662 flux (upper panel) and Hα 10% width (lower panel) vs. time. Data from Cha I are represented by triangles, and those from Tau-Aur are drawn as squares. The maximum extent of the variability is reached after a few days and does not further increase towards longer timescales. distinguish in M-type objects. As can be seen in both panels, the total amount of scatter in the accretion indicator cannot be explained by intrinsic accretion rate variability, which is shown as error bars. This is particularly obvious from the Ca II-λ8662 fluxes, which we identified as the more robust accretion indicator in §3. Converted to accretion rates, the objects cover 1–2 orders of magnitude at any given spectral type, while the variability accounts only for a small fraction of this scatter (∼ 0.35 dex, see §3). In other words, on timescales of days to months the objects do not move significantly in the diagram. Although we cannot exclude the presence of significant variability on longer timescales, the scatter in the Ṁ vs. M plot is most likely dominated by object-toobject variations. This confirms previous claims by Natta et al. (2004) and is contrary to suggestions by Scholz & Jayawardhana (2006) based on a much smaller sample. Thus, stellar mass or parameters that strictly scale with stellar mass are clearly not the only factors that affect the accretion rate. One possible explanation for the spread of accretion rates is the different evolutionary stages of the objects considered here. Accretion rates are expected to drop with age, either following a power law, if the timescale is determined by the slow viscous evolution of the disk (Clarke & Pringle 2006), or in a more rapid process, consistent with the rapid inner disk clearing inferred from the scarcity of ‘transition’ objects with optically thin inner disks. The presence of the clear gap between the accreting and non-accreting populations in Fig. 4, particularly well-defined for the higher mass objects is Fig. 4.— Ca II-λ8662 flux (upper panel) and Hα 10% width (lower panel) vs. spectral types. Error bars indicate the full range of the variations. In the lower panel, we also show the non-accretors (hollow symbols) with Hα 10% widths < 200 km s−1 (after subtracting rotational broadening). Objects without Hα emission are shown with a 10% width of 0 km s−1 . The right axis shows the accretion rates inferred from the relations from Mohanty et al. (2005) and Natta et al. (2004). For the M and K type stars, one sees a clear rise in Hα 10% width towards earlier spectral type, but at K0 and earlier, it appears to plateau, with no objects showing widths larger than ∼ 600 km s−1 . thus highly interesting. The scarcity of objects in the transition from accretors to non-accretors argues for a rapid evolution between states, as suggested in disk evolution models including mass loss from the disk (e.g. Clarke et al. 2001; Armitage et al. 2003). This may be related to the similarly rapid timescale of transition between objects with IR excess and those without. Apart from the evolutionary stage, a number of other factors have been suggested to influence the accretion rates. If disk ionization is a major factor, our findings imply that the main source of ionization is probably external (e.g. cosmic rays), as the stellar ionising radiation is expected to scale with mass and thus would not allow for a wide range of accretion rates at constant mass. In general, our results favour scenarios where initial or environmental conditions are more important than stellar parameters. For instance, a model where the dispersion in initial disk parameters, e.g. mass and radius, determines the distribution of data points in the Ṁ vs. M plot, as recently described by Alexander & Armitage (2006), would be consistent with our results. To constrain the influence of stellar and environmental parameters on the appearance of the Ṁ vs. M plot, it would be necessary to isolate them from evolutionary effects. It would be of particular interest to investigate accretion rates versus other parameters along lines of constant mass in the Ṁ vs. M diagram. Unfortunately, for these tasks the currently available samples may not be sufficiently large. In principle, however, such analy- Accretion Variability ses have the potential to provide important information on the nature of the accretion process in young stars. Finally, the causal relationships between our accretion indicators and accretion rates are not fully understood, and we use them merely as empirical tools, based on published correlations. Nevertheless, the fact that indicators such as Ca II-λ8662 line flux correlate remarkably well with the current empirical framework of measuring accretion rates gives some confidence that these indicators relate to physical reality. 5 We thank the anonymous referee for a very helpful review and constructive critical comments that greatly improved the clarity of the letter. DCN acknowledges gratefully the hospitality of the astronomy group in St. Andrews during his visit. This work was supported in part by NSERC grants to RJ and MHvK and an Early Researcher Award from Ontario to RJ. REFERENCES Alencar, S. H. P., & Batalha, C. 2002, ApJ, 571, 378 Alencar, S. H. P., Basri, G., Hartmann, L., & Calvet, N. 2005, A&A, 440, 595 Alexander, R. D., & Armitage, P. J. 2006, ApJ, 639, L83 Appenzeller, I., & Mundt, R. 1989, A&A Rev., 1, 291 Armitage, P. J., Clarke, C. J., & Palla, F. 2003, MNRAS, 342, 1139 Azevedo, R., et al. J. 2006, A&A, 456, 225 Batalha, C., Lopes, D. F., & Batalha, N. M. 2001, ApJ, 548, 377 Bernstein, R., et al. 2003, SPIE, 4841, 1694 Bertout, C. 1989, ARA&A, 27, 351 Bouvier, J., et al. 2007, in Protostars and Planets V, ed. B. Reipurth, D. Jewitt, & K. Keil, 479–494 Bouvier, J., et al. 1995, A&A, 299, 89 Clarke, C. J., Gendrin, A., & Sotomayor, M. 2001, MNRAS, 328, 485 Clarke, C. J., & Pringle, J. E. 2006, MNRAS, 370, L10 Dullemond, C. P., Natta, A., & Testi, L. 2006, ApJ, 645, L69 Fernandez, M., & Eiroa, C. 1996, A&A, 310, 143 Grankin, K. N., et al. 2007, A&A, 461, 183 Gregory, S. G., et al. 2006, MNRAS, 371, 999 Hartigan, P., et al. 1989, ApJS, 70, 899 Hartigan, P., et al. 1991, ApJ, 382, 617 Herbst, W., et al. 2007, Protostars and Planets V, 297 Herbst, W., et al. 1994, AJ, 108, 1906 Herczeg, G. J., & Hillenbrand, L. A. 2008, ApJ, 681, 594 Johns, C. M., & Basri, G. 1995, AJ, 109, 2800 Joy, A. H. 1945, ApJ, 102, 168 Knott, G. 1891, The Observatory, 14, 97 Krumholz, M. R., McKee, C. F., & Klein, R. I. 2006, ApJ, 638, 369 Mohanty, S., Jayawardhana, R., & Basri, G. 2005, ApJ, 626, 498 Muzerolle, J., Hartmann, L., & Calvet, N. 1998, AJ, 116, 455 Muzerolle, J., et al. 2003, ApJ, 592, 266 Natta, A., et al. 2004, A&A, 424, 603 Nguyen, D. C., et al. 2009, accepted to ApJ Padoan, P., et al. 2005, ApJ, 622, L61 Rydgren, A. E., Strom, S. E., & Strom, K. M. 1976, ApJS, 30, 307 Scholz, A., & Jayawardhana, R. 2006, ApJ, 638, 1056 Scholz, A., Jayawardhana, R., & Brandeker, A. 2005, ApJ, 629, L41 Stelzer, B., Scholz, A., & Jayawardhana, R. 2007, ApJ, 671, 842