EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 577 (May 2015) A&A, 577 (2015) A102 Full HTML
Highlight
Free Access
Issue
A&A
Volume 577, May 2015
Article Number A102
Number of page(s) 12
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201425365
Published online 12 May 2015

© ESO, 2015

1. Introduction

Stars gain most of their mass while deeply embedded in the molecular cloud core out of which they first formed. The accretion rate of matter from the collapsing core onto the star is probably highly variable in time. On cloud scales, turbulence can induce order-of-magnitude variations in the accretion rates onto individual cores (Padoan et al. 2014). Within a given core, variations of similar magnitude are possible from mass-loading onto an embedded disk or pseudo-disk (Zhu et al. 2009; Vorobyov & Basu 2010). The disk becomes unstable at some point, resulting in a short-lived accretion burst onto the star. Each burst causes a temporary increase in luminosity, heating up the circumstellar envelope and altering the chemical composition of the gas and dust. This paper analyzes such chemical effects and explores how to use molecular lines as a probe of variable accretion rates and timescales.

Variable or episodic accretion is the leading explanation for the wide range of protostellar luminosities observed in large-scale surveys (Kenyon et al. 1990; Evans et al. 2009; Kryukova et al. 2012; Fischer et al. 2013; Dunham et al. 2014) and for FUor and EXor outbursts in more evolved pre-main-sequence stars (Herbig 1977; Audard et al. 2014). Episodic accretion can also affect processes such as disk fragmentation (Stamatellos et al. 2011, 2012) and lithium depletion (Baraffe & Chabrier 2010). It is therefore important to understand the magnitude-frequency distribution and determine how often bursts of a certain intensity occur. Based on large-scale variability surveys, the strongest bursts – with accretion rates of at least 100 times the quiescent value – happen every 5–50 kyr per protostar (Scholz et al. 2013). The spacings between periodic shocks in jets and outflows support quiescent intervals anywhere from 104 yr down to 10 yr, separating bursts with an unknown range of accretion rates (Devine et al. 1997; Raga et al. 2002; Arce et al. 2013). High-cadence photometric surveys show luminosity variations at the 5–50% level on timescales of hours to weeks (Billot et al. 2012; Stauffer et al. 2014), though only some of that may be related to accretion variability (Morales-Calderón et al. 2011; Rebull et al. 2014). Despite the incomplete statistics, it appears that stronger accretion bursts occur less frequently than weaker ones. Such a distribution is supported by numerical simulations (Zhu et al. 2009; Vorobyov & Basu 2010).

Direct observations of luminosity bursts were traditionally limited to the relatively evolved Class II sources, where the lack of an obscuring dusty envelope makes the bursts visible at optical wavelengths. However, bursts have now also been detected in at least five embedded Class I protostars (Audard et al. 2014) and even in one deeply embedded Class 0 source (Safron et al. 2015). The detection of luminosity flares during the earliest phases of star formation lends credence to the notion of episodic accretion as a wide-spread phenomenon.

Accretion bursts have the ability to alter the chemical composition of the circumstellar material. For example, the 2008 burst in EX Lup sparked the production of crystalline silicates (Ábrahám et al. 2009) and boosted the column densities of mid-infrared H2O and OH lines (Banzatti et al. 2012). In very low-luminosity embedded protostars, whose envelopes are too cold to produce pure CO2 ice, the detection of double-peaked 15-μm absorption profiles is attributed to a higher luminosity at some unknown point in the past (Kim et al. 2011, 2012; see also Poteet et al. 2013).

The chemical effects of embedded accretion bursts primarily result from changes in the temperature of the circumstellar gas and dust. When a protostar enters a burst, the increase in luminosity can alter the D/H ratio of water and other species (Owen & Jacquet 2015). The higher temperatures also lead to the evaporation of some of the icy grain mantles (Lee 2007; Visser & Bergin 2012, hereafter Paper I; Vorobyov et al. 2013). After the burst ends, the luminosity returns to the quiescent value and the envelope cools down almost instantaneously (Johnstone et al. 2013). However, at typical envelope densities, it takes 103–105 yr for the gas to freeze back onto the cold dust grains. The abundance profiles are out of equilibrium with the observed luminosity for all that time. This could explain the presence of spatially extended C18O in eight out of 16 low-mass protostars (Jørgensen et al. 2015).

The best example to date of episodic accretion chemistry in an embedded protostar is the central gap discovered in spatially resolved observations of H13CO+J = 4–3 in IRAS 15398 (Jørgensen et al. 2013). This distribution is consistent with a picture where HCO+ and its isotopologs are destroyed by H2O in the inner envelope when the temperature exceeds 100 K. However, the current temperature at the edge of the H13CO+ gap is only ~30 K. The mismatch between the dust temperature and the H13CO+ morphology implies that IRAS 15398 was hotter in the past and that the chemistry is still adjusting to the current cold environment. Based on the freeze-out timescale of H2O, the simplest explanation for that hotter past is an accretion burst that happened 100–1000 yr ago (Jørgensen et al. 2013).

The aim of the current paper is twofold: to explore episodic accretion chemistry in a sample of 25 spherical envelope models, and to address how certain molecular line ratios can be used as a chronometer of when the most recent accretion burst occurred in any given source. Section 2 introduces the sample and summarizes the available observations, including new Herschel spectra of HCO+ and N2H+ 6–5. Section 3 presents the physical and chemical models. Section 4 discusses the chemical effects of an accretion burst and Sect. 5 shows how to use observed line ratios to derive episodic accretion timescales. Lastly, Sect. 6 presents the main conclusions.

2. Source sample and observations

Paper I explored the chemical aspects of episodic accretion with a set of single-point models at representative densities of 105, 106, and 107 cm-3. The current paper analyzes the chemistry in a sample of observationally determined spherical protostellar envelope models, each covering a range of densities and temperatures. We use the predicted abundance profiles to simulate molecular line spectra and identify observable diagnostics.

Table 1

Observed integrated intensities ( in K km s-1) and half-power beam widths (θ in arcsec)a.

The envelope models are based on the 29 low-mass embedded protostars targeted in the key program “Water in star-forming regions with Herschel” (WISH; van Dishoeck et al. 2011). Source coordinates, distances, and other basic properties are listed in Table 1 of Kristensen et al. (2012). In order to identify the best observable diagnostics, we searched the literature for observations of molecular lines and ice column densities. The most widely available tracers are C18O J = 2–1, 3–2, and 5–4 (20, 26, and 15 sources); H13CO+ 1–0, 3–2, and 4–3 (18, 17, and 16); and N2H+ 1–0 (21). In addition, column densities of CO2, CO, and H2O ice are available for 14, 10, and 12 sources. Lastly, we present previously unpublished spectra of HCO+ 6–5 and N2H+ 6–5 (16 and 7 sources) obtained with the Heterodyne Instrument for the Far-Infrared (HIFI; de Graauw et al. 2010) on the Herschel Space Observatory (Pilbratt et al. 2010). The full data set is summarized in Tables 1 and 2. Appendix A offers more details on the observations and the data reduction.

For lack of molecular line observations, we omit NGC 1333 IRAS3, IRAS 12496, R CrA IRS5, and HH100 IRS from the WISH target list. Our final sample therefore contains 25 protostars, including IRAS 15398 (Sect. 1; Jørgensen et al. 2013). These sources have bolometric luminosities from 0.8 to 35.7 L⊙ and span the full evolutionary sequence from deeply embedded Class 0 to late Class I (Lada 1999).

3. Model description

3.1. Physical framework

Spherical density and temperature profiles are available from Kristensen et al. (2012) for our entire sample of 25 protostars. The density follows a power law, n(H2) ∝ r− p, where the slope p is one of three free parameters in a grid of models from the one-dimensional continuum radiative transfer program DUSTY (Ivezić et al. 1999; Jørgensen et al. 2002). The other two are the mass and size of the envelope. Kristensen et al. compared the model grid output to the observed spectral energy distributions (SEDs) and submillimeter brightness profiles to find the best-fit parameters for each source (see their Table C.1).

The WISH target list is biased towards bright protostars, which are likely to be in a relatively active phase of accretion. Nonetheless, for the purpose of exploring the chemical effects of episodic accretion, we treat all sources as if they are currently in a quiescent phase. Alternative assumptions of some or all sources currently experiencing a high accretion rate would produce a different set of chemical models, but would not alter the conclusions in Sect. 4.

Within the best-fit models from Kristensen et al. (2012), we simulate an accretion burst by changing the stellar luminosity and rerunning DUSTY to obtain a new temperature profile. All other parameters remain constant. Our goal is to study the chemical effects of a strong burst, such as might happen at most a few times during the embedded phase. Specifically, the burst luminosity for each source is set to be 100 times higher than the current luminosity; this is the burst intensity inferred for IRAS 15398 (Jørgensen et al. 2013). The simulated bursts last for 100 yr (Vorobyov & Basu 2005), although shorter durations down to at least 1 yr do not affect our results. After each burst ends, we follow the quiescent chemistry for a maximum of 105 yr to cover the full range of potential timescales from, e.g., Scholz et al. (2013). The models are static and ignore any dynamical effects of the envelope collapsing and dissipating on typical timescales of a few 105 yr (Evans et al. 2009; Visser et al. 2011).

Table 2

Observed ice column densities (N(X) in 1017 cm-2).

The heating and cooling rates of the dust are fast enough for our purposes that the dust temperature responds instantaneously to a change in luminosity (Johnstone et al. 2013). We set the gas temperature equal to the dust temperature at all times. As noted in Paper I, the temperature at any point in the envelope is proportional to and thus increases by a factor of 3.2 if the star becomes 100 times brighter (see also Johnstone et al. 2013). As an example, the top panel of Fig. 1 in Sect. 4 shows the “quiescent” and “burst” temperature profiles for IRAS 15398.

thumbnail Fig. 1

Top panel: radial density profile in IRAS 15398 (black), along with the temperatures in the quiescent phase (blue) and during the burst (red). Other panels: radial abundance profiles of CO, H2O, N2H+, and HCO+. The red curves are during the accretion burst. The black curves are during the quiescent phase, at increasing amounts of time after the burst: 10 yr (dotted), 150 yr (solid), 800 yr (dotted), 2600 yr (short-dashed), 6200 yr (dash-dotted), and 24 000 yr (long-dashed). The quiescent profiles culminate in the blue curve at 1 × 105 yr after the burst.

3.2. Chemical network

The basis of our chemical network is the Rate12 release of the UMIST Database for Astrochemistry (McElroy et al. 2013). In addition to standard neutral-neutral and ion-molecule chemistry, Rate12 includes photodissociation and photoionization reactions. These are important only at the outer edges of the envelopes, which we assume are irradiated by the mean interstellar radiation field. The cosmic-ray ionization rate is set to a standard value of 5 × 10-17 s-1 (Dalgarno 2006). X-rays are not included; they have been detected in one T Tauri FUor (Liebhart et al. 2014), but X-ray luminosities in embedded protostars remain highly uncertain and are beyond the scope of this work.

We expand the network with freeze-out and evaporation of all neutral molecules. The freeze-out timescale depends inversely on the density and is typically of similar duration (~103−104yr) as the quiescent phase between subsequent accretion bursts, allowing for an abundance pattern out of equilibrium with the stellar luminosity (Lee 2007; Paper I). Indeed, it is the long freeze-out timescale that drives most of the chemical aspects of episodic accretion explored in this work.

Non-thermal desorption via cosmic rays (Hasegawa & Herbst 1993) or UV radiation (Öberg et al. 2007) is included but only plays a minor role in episodic accretion chemistry. Much more important is thermal desorption. According to spectroscopic observations, protostellar CO evaporates at 25−30 K (Jørgensen et al. 2002, 2004, 2015; Yıldız et al. 2013). This is at least 5 K higher than predicted from binding energies measured in the laboratory, especially the “pure ice” value of 855 K from Bisschop et al. (2006) used in Paper I. Hence, we now adopt a binding energy of 1307 K, measured for CO evaporating from amorphous water ice (Noble et al. 2012). Under protostellar conditions, this higher binding energy is consistent with the empirical evaporation temperature of 25–30 K. Likewise, we increase the binding energy of N2 from 800 K to 1200 K. The binding energy for atomic O has long been set to 800 K (Watson & Salpeter 1972), but new experiments on water ice surfaces offer conclusive evidence that this old estimate is too low (Minissale 2014; Minissale et al., in prep.). We adopt the current best estimate of 1420 K. The binding energies for other molecules are unchanged from Paper I, such as 5773 K for H2O and 2300 K for CO2 (Fraser et al. 2001; Noble et al. 2012).

Our network contains two types of grain-surface reactions. The first is the usual formation of H2 (Black & van Dishoeck 1987). The second type is simple hydrogenation of C to CH4, N to NH3, and O to H2O. This happens one H atom at a time at a rate equal to the adsorption rate of H onto the grain surface multiplied by the relative fraction of the reactant molecule in the ice (Visser et al. 2011). Conversion of CO ice into CO2 ice, as proposed by Kim et al. (2011, 2012), is not included; see Sect. 5.1 for details.

The elemental abundances relative to the total number of H atoms are 0.09 for He, 1.4 × 10-4 for C, 7.5 × 10-5 for N, and 3.2 × 10-4 for O (McElroy et al. 2013). At the start of the first accretion burst, the chemical composition is set to typical pre-stellar core conditions (Maret et al. 2006; Whittet et al. 2009): hydrogen in atomic H (0.005%) and H2 (~100%); carbon in CO (37%), CO ice (35%), and CO2 ice (28%); remaining oxygen in H2O ice; and nitrogen in atomic N (55%), N2 (42%), and NH3 ice (3%). The initial ratio of CO2 ice to CO ice is 0.8:1, the average value observed in three dense molecular clouds believed to be representative of the earliest stage of star formation (Whittet et al. 2009).

3.3. Line radiative transfer

In order to compare our model to observations, we simulate molecular line spectra with the one-dimensional radiative transfer code RATRAN (Hogerheijde & van der Tak 2000). RATRAN solves for the level populations as a function of position, accounting for both collisional and radiative excitation, and then performs ray tracing to synthesize a spectral image. Lastly, the images are convolved to the appropriate telescope beam for each source and transition (Table 1).

The low-J lines of CO and HCO+ are optically thick for most sources, so the optically thin isotopologs C18O and H13CO+ are used instead. The isotope ratios are set to 69 for 12C/13C and 557 for 16O/18O (Wilson 1999). Collision rates are taken from Flower (1999), Daniel et al. (2005), and Yang et al. (2010), as compiled in the Leiden Atomic and Molecular Database1 (LAMDA; Schöier et al. 2005).

4. Effects of accretion bursts on abundance profiles and rotational spectra

As described in Sect. 3.1, we simulate an accretion burst in all 25 sources by increasing the luminosity to 100 times the observed luminosity for a duration of 100 yr. Starting with the pre-stellar core composition from Sect. 3.2, we evolve the chemistry during the burst at 78–89 radial points per source. After 100 yr, the burst ends and the temperature instantaneously returns to the quiescent profile. Continuing with the abundances from the end of the burst, we evolve the chemistry for another 105 yr and extract abundances at a number of intermediate time steps.

All sources have the same qualitative abundance profiles and molecular line spectra, so we choose IRAS 15398 for a quantitative discussion. The results for the rest of the sample are presented in Sect. 5.

thumbnail Fig. 2

Synthetic spectra for the quiescent phase in IRAS 15398, scaled as indicated. The different colors correspond to the same post-burst time steps as in Fig. 1: 10 yr (black), 150 yr (red), 800 yr (gray), 2600 yr (purple), 6200 yr (green), 24 000 yr (orange), and 1 × 105 yr (blue). The arrows indicate how the intensity of each line changes as a function of time; the colored flags mark the peak intensities at the listed times. The N2H+ 1–0 panel shows the isolated F1 = 0−1 hyperfine component, whereas the 6–5 panel shows the sum over all unresolved components.

4.1. Abundance profiles

IRAS 15398 currently has a luminosity of 1.6 L⊙ and an envelope mass of 0.5 M⊙ (Kristensen et al. 2012). The top panel of Fig. 1 plots the temperature profiles for the quiescent phase (Tq) and the burst (Tb) in IRAS 15398, along with the density profile (n). The quiescent temperature decreases from 250 K at the inner edge of the envelope to 10 K at the outer edge. During the burst, the entire curve goes up by a factor of 3.2 (Sect. 3.1).

The other panels of Fig. 1 show the radial abundance profiles of CO, H2O, N2H+, and HCO+ during the burst (red curves) and at a series of time steps after the burst (black and blue). The basic results for CO are the same as reported elsewhere (Lee 2007; Paper I; Vorobyov et al. 2013). The temperature exceeds 30 K everywhere during the burst and more than 99.9% of CO is in the gas phase at all radii2. After the burst, CO freezes out beyond 250 AU, where the quiescent temperature lies below 25 K. Since the freeze-out rate is proportional to the gas density, depletion is fastest right around the 25 K radius. For the specific density profile in the envelope model of IRAS 15398, it takes 700 yr for the CO abundance to drop by a factor of 10 and 1500 yr for a factor of 100.

The column density of CO ice integrated through the envelope is reduced almost to zero during the burst. Once the envelope cools down, the ice column starts to reform. It reaches 20% of the original amount after 100 yr and 90% after 4000 yr. CO2 has a higher binding energy than CO and evaporates around 45 K. This allows 15% of the integrated CO2 ice column to survive during the simulated burst in IRAS 15398. It takes 600 yr in the quiescent phase to rebuild 90% of the original CO2 ice.

H2O evaporates out to 100 AU during the burst. Like CO, the excess H2O refreezes in the quiescent phase, but it does so about ten times faster because of the higher densities in the inner envelope. All H2O that evaporated during the burst disappears within a few 100 yr in the envelope model of IRAS 15398. On timescales of a few 1000 yr, the freeze-out of CO removes a substantial amount of oxygen from the gas phase. This leads to a drop in the H2O profile near the 25 K radius at 250 AU. The H2O abundance increases again towards the outer envelope because of photodesorption. The integrated column density of H2O ice drops by 65% during the burst; afterwards, it takes only ~10 yr to rebuild the column to 90% of the equilibrium value.

The chemistry of N2H+ is strongly tied to that of CO and H2O because both act as important destruction channels. The depletion of CO and H2O in the quiescent phase thus allows N2H+ to increase in abundance at almost all radii for the first few 1000 yr after the burst. This type of anti-correlation between CO and N2H+ is well known from observations of molecular clouds, embedded protostars, and circumstellar disks (Bergin et al. 2002; Jørgensen 2004; Qi et al. 2013) and was also noted by Lee (2007). On timescales of more than a few 1000 yr, freeze-out of N2 near 400 AU causes a decrease in the N2H+ abundance.

Lastly, the abundance of HCO+ undergoes a seesaw motion around the 25 K radius. CO is the primary parent species of HCO+ and H2O is one of the dominant destroyers. Freeze-out of H2O allows more HCO+ to be formed in the inner envelope, while freeze-out of CO reduces the amount of HCO+ in the outer envelope on longer timescales. The correlation between the protostellar abundances of CO and HCO+ is another well-established observational result (Jørgensen et al. 2004).

The chemical evolution in the other 24 envelope models is qualitatively the same as in IRAS 15398, but the timescales on which the aforementioned changes occur differ depending on the density profiles. IRAS 15398 has relatively high densities in the H2O and CO freeze-out zones compared to the rest of the sample. As a consequence, freeze-out in the other sources is typically slower by up to two orders of magnitude and the chemical effects of the accretion burst take longer to subside.

In summary, the chemical composition of a protostellar envelope changes substantially during an accretion burst. These changes last for long after the system returns to quiescence:

  • enhanced CO in the outer envelope for a few 103–104 yr;

  • reduced CO ice and CO2 ice in the outer envelope for a few 103–104 yr;

  • enhanced H2O in the inner envelope for a few 102–103 yr;

  • reduced N2H+ in the outer envelope for a few 103–104 yr;

  • seesaw pattern for HCO+: enhanced in the inner envelope, reduced in the outer envelope for a few 103–104 yr.

The exact timescales in a given source depend on the density profile. The next section explores the consequences of the abundance changes on the rotational emission lines.

thumbnail Fig. 3

Simulated integrated line intensities and ice column densities in IRAS 15398 as a function of time during the quiescent phase. The N2H+J = 1–0 and 6–5 intensities are summed over all HF components. Where available, the gray bars show the observed intensities and columns from Tables 1 and 2 with 20% uncertainties.

4.2. Rotational spectra and ice column densities

The abundance changes caused by an accretion burst persist for long after the burst has ended. To what extent can the chemical effects of a burst still be observed once the protostar has returned to quiescence? Figure 2 shows a set of nine single-dish spectra corresponding to the abundance profiles at the quiescent time steps from Fig. 1. All spectra are first convolved with the appropriate telescope beam (Table 1) and then continuum-subtracted. The arrow in each panel indicates schematically how the peak intensity changes with time: continuous decrease, continuous increase, or an increase followed by a decrease.

For a more quantitative view of the temporal behavior, the spectra from Fig. 2 are integrated and plotted as a function of time in Fig. 3. Also shown here are the simulated column densities of CO2, CO, and H2O ice. Where available, the observed intensities and column densities from Sect. 2 are included in Fig. 3 as gray bars with 20% error margins. We compare the model results to the observations in Sect. 5.1; the remainder of the current section discusses the model trends on their own.

During the quiescent phase, the simulated line intensities change predictably based on the abundance variations. The freeze-out of CO leads to a monotonic decrease in the optically thin C18O lines. H13CO+ 1–0 and 3–2 are dominated by the cold outer envelope and their intensities also decrease monotonically as CO and HCO+ are depleted. H13CO+ 4–3 and HCO+ 6–5 originate closer to the star and pick up on the seesaw abundance pattern from Fig. 1: the line intensities show an initial rise before depletion in the outer envelope sets in and the intensities go down.

Because of the abundance anticorrelation between CO and N2H+, the line intensities of C18O and N2H+ evolve largely in opposite directions. As CO and H2O freeze out after the burst, N2H+ becomes more abundant and its 1–0 and 6–5 lines gain intensity. The 1–0 intensity decreases again on timescales of more than ~24 000 yr after the burst. This is due to the depletion of atomic N and the resulting loss of N2H+ in the cold outer envelope. The 6–5 line is not very sensitive to the N2H+ abundance around 1000 AU; its intensity levels off at late times, but does not turn over.

The total observable column densities of CO2 and CO ice increase by factors of 3 and 20 on timescales of 103–104 yr after the burst. The column of H2O ice is only marginally affected by changes in luminosity and remains constant at all times.

5. Line ratios as a diagnostic of episodic accretion timescales

One of the challenges in episodic accretion is to understand the magnitude–frequency distribution: how often do burst of a certain magnitude or intensity occur? The chemical effects explored in this work can help constrain the frequency of very strong bursts, where the accretion rate increases by at least a factor of 100. In the quiescent phase following such a burst, several observables vary monotonically with time. In principle, any of these can be used as a chronometer of when the most recent burst occurred. For a statistical sample of protostars, it then becomes possible to determine the average burst frequency.

With the current set of spherical models, inferring when the most recent burst occurred in any given source will be a rather crude method with uncertainties at the order-of-magnitude level. The following discussion should therefore be seen as an attempt to explore certain methodologies, rather than an attempt to derive firm numbers. For an alternative approach based on spatially resolved observations of C18O, see Jørgensen et al. (2015).

5.1. Models vs. observations

Figure 3 highlights both the successes and the shortcomings of our model for IRAS 15398. Six of the observables are reproduced at some point during the quiescent phase, but the times at which the matches occur range over more than two orders of magnitude: from ~30 yr for CO ice to ~104yr for C18O 2–1. Furthermore, H2O ice is overproduced at all times, and HCO+ 6−5 and N2H+ 1–0 are underproduced.

The overproduction of H2O ice may be due to a lack of grain-surface reactions in our chemical network, which tend to drain oxygen out of H2O and convert it to other species (Schmalzl et al. 2014). N2H+ 1–0 probably suffers from our choice of a spherically symmetric envelope model. In spatially resolved observations, the morphology of N2H+ 1–0 often deviates significantly from spherical symmetry (Jørgensen 2004; Tobin et al. 2011). This can result in enhanced emission in one direction, e.g., from interactions between the outflow and the dense envelope, which our spherical models do not reproduce. The underproduction of HCO+ 6–5 appears to be an excitation effect: this line is enhanced along the outflow cavity walls due to direct ultraviolet heating, similar to what is seen in 12CO and 13CO 6–5 (Spaans et al. 1995; Yıldız et al. 2012).

Kim et al. (2011, 2012) argued that some fraction of CO ice needs to be converted into CO2 ice during each quiescent phase in between subsequent bursts. However, that conclusion may be an artifact of the low initial CO2 abundance in their model. Our model starts with the average CO and CO2 ice abundances observed towards dense cores (Whittet et al. 2009). As a result, it does not require any conversion of CO into CO2 ice to reproduce the observed column densities in IRAS 15398.

The comparison between model predictions and observations is expanded to all 25 protostars in Fig. 4. The trends identified in Figs. 2 and 3 for IRAS 15398 generally hold for the full sample: all C18O lines become weaker as a function of time and usually match the observed intensities at some point. H13CO+ 1−0 and 3−2 also decrease monotonically and are in reasonable agreement with the observations. H13CO+ 4–3 and HCO+ 6–5 reflect the seesaw abundance pattern of HCO+ by first increasing and then decreasing in strength. Both lines are too weak in the models by a factor of 2–10. The N2H+ 1–0 line gains intensity during the quiescent phase in all sources and turns over at late times in about a third of the sample. N2H+ 6–5 initially increases and then levels off. Both N2H+ lines tend to be underproduced relative to the observed intensities. The three ice column densities always increase with time after the burst. CO2 ice and CO ice usually match the observations, but H2O ice is consistently too abundant in the models.

In order to construct a consistent set of models, we treated all sources as if they are currently in the quiescent phase (Sect. 3.1). This is undoubtedly a false assumption – especially for the brighter sources in our target list – and adds another source of discrepancies between models and observations. Future studies of individual sources or large samples will benefit from independent constraints on the accretion rates, e.g., from hydrogen Brγ (Connelley & Greene 2010). If the stellar mass is known from measurements of a rotationally supported disk (Tobin et al. 2012), the accretion rate can also be computed from the relationship between Lbol, M∗, and Ṁ (Adams & Shu 1986).

5.2. Line ratios

The absolute line intensities in Fig. 3 are sensitive to uncertainties in the source distance and in various model parameters like the envelope mass and the temperature profile. A better approach is to take the ratio between two rotational lines whose intensities change in opposite directions after the burst. This effectively acts as a self-normalization against distance and source model, in particular for line pairs at similar excitation levels.

Table 3 lists thirteen observable ratios as candidate chronometers, chosen in part for their sensitivity to the burst, and in part for being easily observable with modern submillimeter facilities. The first four ratios are the low-J line pairs of C18O 2–1, C18O 3–2, H13CO+ 1–0, and H13CO+ 3–2 over N2H+ 1–0. The fifth one is the high-J pair of C18O 5–4 over N2H+ 6–5. The remaining eight are the ratios between the intensities of C18O 2–1, C18O 3–2, H13CO+ 1–0, and H13CO+ 3–2 on the one hand and the column densities of CO2 or CO ice on the other hand.

Table 3

Candidate observable ratios to be used as chronometers of episodic accretion timescales.

thumbnail Fig. 4

Comparison between model predictions and observations. This plot is similar to Fig. 3, except that the model values are now divided by the observed quantities. A ratio of unity means the model is in exact agreement with observations. If the ratio is larger than unity, the model overpredicts the observation, and vice versa. The horizontal scale for each source runs logarithmically from 10 to 105 yr after the burst. The observable quantity in the bottom three panels is the column density of ice along a pencil beam to the central star. In the other panels, it is the integrated intensity of various rotational lines. The N2H+J = 1–0 and 6–5 intensities are summed over all HF components. The sources are ordered from left to right by increasing luminosity.

As an example, the left panel of Fig. 5 shows how the ratio between the integrated intensity of C18O 2–1 and the column density of CO ice evolves as a function of time after the burst for the full source sample. The ratio plots for the other 12 pairs of observables in Table 3 look qualitatively the same. Each curve represents one source model. The ratios decreases with time for all sources, but the curves are scattered across more than four orders of magnitude and are of little diagnostic value in this form.

thumbnail Fig. 5

Left: ratio between the integrated intensity of C18O 2–1 and the column density of CO ice as a function of time during the quiescent phase. Each gray curve corresponds to one source. Right: same, but with the time axis normalized by the characteristic freeze-out timescale for each source.

The scatter in Fig. 5 is due to the range of envelope parameters encountered across the sample. CO freezes out at the same temperature in each source (< 25 K), but because of all the different envelope parameters, it does not do so at the same density. In high-density envelopes, the freeze-out timescale is shorter and the CO gas/ice abundance ratio decreases more rapidly than in lower-density sources. Likewise, the C18O/N2H+ and H13CO+/N2H+ line intensity ratios change at different rates for different sources.

The scatter in the left panel of Fig. 5 can be reduced by expressing the time axis for each source in units of its characteristic freeze-out timescale (Paper I, erratum): (1)The characteristic freeze-out density nfr is defined as the density at which the temperature reaches Tfr in the quiescent phase. We tested Tfr from 10 to 40 K and found the greatest reduction in scatter for 15 K, where the bulk of the low-J C18O emission originates. At 15 K, nfr ranges from 3.5 × 103 cm-3 in Elias 29 to 5.5 × 106 cm-3 in Ser SMM4. The corresponding freeze-out timescales vary from 2.3 Myr to 1500 yr.

The time-normalized C18O 2–1/CO ice ratios appear as the gray curves in the right panel of Fig. 5. The sources are still scattered across an order of magnitude, but all follow the same trend and an average relationship between normalized post-burst time and gas/ice observable ratio can easily be drawn. Our source sample is probably diverse enough (albeit not unbiased) that we can use this average relationship also for other samples where detailed envelope models are not readily available.

5.3. Future prospects

Using all the observable ratios from Table 3, we can estimate when the most recent accretion burst occurred in each of our sources. However, the results are currently too scattered and the uncertainties too large to produce reliable numbers worth tabulating. For example, Jørgensen et al. (2013) concluded from spatially resolved H13CO+ 4–3 that IRAS 15398 experienced a burst between 100 and 1000 yr ago. We find the same timescale from the C18O gas/CO2 ice ratios, but C18O gas/CO ice suggests a longer timescale of ~6000 yr, whereas H13CO+ 3–2/N2H+ 1–0 yields a value as long as 5 × 104 yr. Furthermore, the observed C18O/N2H+ line ratios do not match any timescale from the model. The situation is no better for most of the other sources. This lack of consistency between different tracers was already evident in Figs. 3 and 4, where some curves match the observations at very different post-burst times (Sect. 5.1).

The most likely reason for these discrepancies is the assumption of spherical symmetry for the envelope models (Sect. 5.1). However, even if the spherical models do not produce reliable timescales, they are still useful as a basis for follow-up studies. Spatially resolved observations offer many advantages over single-dish data (as explored for, e.g., C18O 2–1 by Jørgensen et al. 2015) and can be carried out for statistically significant source samples with current facilities. We will expand the models beyond spherical symmetry and calibrate them against existing interferometric data to identify the most reliable tracers of episodic accretion timescales. These tracers can then be used to derive timescales for many tens or hundreds of protostars and compared against results from other methods. Working with such large samples also reduces the problem that any given source may be an outlier.

6. Conclusions

This paper simulates the chemical effects of episodic accretion in a diverse sample of 25 protostellar envelope models. The intensity of the accretion bursts is set to 100 times the quiescent luminosity, enough to evaporate CO ice all the way to the outer edge of the envelope in each source (Sect. 4.1). The H2O snowline expands in size by a factor of ~10 during the burst. The evaporation of CO and H2O leads to the destruction of N2H+ throughout the envelope. HCO+ is enhanced outside of the H2O snowline and destroyed inside of it.

When the burst ends, the temperatures quickly return to normal. The species that evaporated during the burst start to refreeze onto the cold dust, but this is a relatively slow process. H2O remains enhanced for 102–103 yr and CO for an order of magnitude longer. The abundances of N2H+ and HCO+ are reset on similarly long timescales.

The abundance changes resulting from an accretion burst also have a substantial effect on various molecular lines, as quantified by radiative transfer simulations at a series of time steps (Sect. 4.2). C18O 2–1, 3–2, and 5–4 and H13CO+ 1–0 and 3–2 are enhanced during the burst and decay monotonically after the burst ends. H13CO+ 4–3, HCO+ 6–5, and N2H+ 1–0 first increase after the burst and then decrease. The column densities of CO and CO2 ice increase monotonically after the burst.

The simulated line intensities and ice column densities are compared to single-dish observations gathered from the literature and to previously unpublished spectra of HCO+ and N2H+ 6–5 from the Herschel Space Observatory. Several line ratios are identified that can be used, in principle, as a chronometer of when the most recent strong burst occurred in each source (Sect. 5). In practice, however, the spherical source models are not accurate enough to derive reliable timescales from single-dish data. Interferometric observations are therefore recommended for follow-up work. In order to tackle a statistically significant source sample, the models need to be calibrated against spatially resolved observations to identify robust tracers of episodic accretion timescales.

2

Lee (2007) simulated a lower-luminosity protostar and retained some CO ice in the outer envelope during the burst.

Acknowledgments

This work was supported by the National Science Foundation under grant 1008800 and by a Postdoctoral Fellowship from the European Southern Observatory. J.K.J. was supported by a Lundbeck Foundation Junior Group Leader Fellowship and by the Centre for Star and Planet Formation, funded by the Danish National Research Foundation. This research made use of NASA’s Astrophysics Data System and of the SIMBAD database, operated at CDS, Strasbourg, France.

References

  1. Ábrahám, P., Juhász, A., Dullemond, C. P., et al. 2009, Nature, 459, 224 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  2. Adams, F. C., & Shu, F. H. 1986, ApJ, 308, 836 [NASA ADS] [CrossRef] [Google Scholar]
  3. Aikawa, Y., Kamuro, D., Sakon, I., et al. 2012, A&A, 538, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Arce, H. G., Mardones, D., Corder, S. A., et al. 2013, ApJ, 774, 39 [Google Scholar]
  5. Audard, M., Ábrahám, P., Dunham, M. M., et al. 2014, Protostars and Planets VI, 387 [Google Scholar]
  6. Banzatti, A., Meyer, M. R., Bruderer, S., et al. 2012, ApJ, 745, 90 [NASA ADS] [CrossRef] [Google Scholar]
  7. Baraffe, I., & Chabrier, G. 2010, A&A, 521, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Benz, A. O., Bruderer, S., van Dishoeck, E. F., et al. 2015, A&A, 576, A109 [Google Scholar]
  9. Bergin, E. A., Alves, J., Huard, T., & Lada, C. J. 2002, ApJ, 570, L101 [NASA ADS] [CrossRef] [Google Scholar]
  10. Billot, N., Morales-Calderón, M., Stauffer, J. R., Megeath, S. T., & Whitney, B. 2012, ApJ, 753, L35 [NASA ADS] [CrossRef] [Google Scholar]
  11. Bisschop, S. E., Fraser, H. J., Öberg, K. I., van Dishoeck, E. F., & Schlemmer, S. 2006, A&A, 449, 1297 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Black, J. H., & van Dishoeck, E. F. 1987, ApJ, 322, 412 [NASA ADS] [CrossRef] [Google Scholar]
  13. Blake, G. A., Sandell, G., van Dishoeck, E. F., et al. 1995, ApJ, 441, 689 [Google Scholar]
  14. Boogert, A. C. A., Tielens, A. G. G. M., Ceccarelli, C., et al. 2000, A&A, 360, 683 [NASA ADS] [Google Scholar]
  15. Boogert, A. C. A., Hogerheijde, M. R., Ceccarelli, C., et al. 2002, ApJ, 570, 708 [NASA ADS] [CrossRef] [Google Scholar]
  16. Boogert, A. C. A., Pontoppidan, K. M., Knez, C., et al. 2008, ApJ, 678, 985 [NASA ADS] [CrossRef] [Google Scholar]
  17. Butner, H. M., Lada, E. A., & Loren, R. B. 1995, ApJ, 448, 207 [NASA ADS] [CrossRef] [Google Scholar]
  18. Caselli, P., Myers, P. C., & Thaddeus, P. 1995, ApJ, 455, L77 [NASA ADS] [CrossRef] [Google Scholar]
  19. Connelley, M. S., & Greene, T. P. 2010, AJ, 140, 1214 [NASA ADS] [CrossRef] [Google Scholar]
  20. Cook, A. M., Whittet, D. C. B., Shenoy, S. S., et al. 2011, ApJ, 730, 124 [NASA ADS] [CrossRef] [Google Scholar]
  21. Dalgarno, A. 2006, Proc. Nat. Acad. Sci., 103, 12269 [NASA ADS] [CrossRef] [Google Scholar]
  22. Daniel, F., Dubernet, M.-L., Meuwly, M., Cernicharo, J., & Pagani, L. 2005, MNRAS, 363, 1083 [NASA ADS] [CrossRef] [Google Scholar]
  23. de Graauw, T., Helmich, F. P., Phillips, T. G., et al. 2010, A&A, 518, L6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  24. Devine, D., Bally, J., Reipurth, B., & Heathcote, S. 1997, AJ, 114, 2095 [NASA ADS] [CrossRef] [Google Scholar]
  25. Dunham, M. M., Stutz, A. M., Allen, L. E., et al. 2014, Protostars and Planets VI, 195 [Google Scholar]
  26. Emprechtinger, M., Caselli, P., Volgenau, N. H., Stutzki, J., & Wiedner, M. C. 2009, A&A, 493, 89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Evans, II, N. J., Lee, J.-E., Rawlings, J. M. C., & Choi, M. 2005, ApJ, 626, 919 [NASA ADS] [CrossRef] [Google Scholar]
  28. Evans, N. J., Dunham, M. M., Jørgensen, J. K., et al. 2009, ApJS, 181, 321 [NASA ADS] [CrossRef] [Google Scholar]
  29. Fischer, W. J., Megeath, S. T., Stutz, A. M., et al. 2013, Astron. Nachr., 334, 53 [NASA ADS] [CrossRef] [Google Scholar]
  30. Flower, D. R. 1999, MNRAS, 305, 651 [NASA ADS] [CrossRef] [Google Scholar]
  31. Fraser, H. J., Collings, M. P., McCoustra, M. R. S., & Williams, D. A. 2001, MNRAS, 327, 1165 [NASA ADS] [CrossRef] [Google Scholar]
  32. Gerakines, P. A., Schutte, W. A., Greenberg, J. M., & van Dishoeck, E. F. 1995, A&A, 296, 810 [NASA ADS] [Google Scholar]
  33. Gregersen, E. M., Evans, II, N. J., Zhou, S., & Choi, M. 1997, ApJ, 484, 256 [NASA ADS] [CrossRef] [Google Scholar]
  34. Gregersen, E. M., Evans, II, N. J., Mardones, D., & Myers, P. C. 2000, ApJ, 533, 440 [NASA ADS] [CrossRef] [Google Scholar]
  35. Hasegawa, T. I., & Herbst, E. 1993, MNRAS, 261, 83 [NASA ADS] [CrossRef] [Google Scholar]
  36. Herbig, G. H. 1977, ApJ, 217, 693 [NASA ADS] [CrossRef] [Google Scholar]
  37. Hogerheijde, M. R., & van der Tak, F. F. S. 2000, A&A, 362, 697 [NASA ADS] [Google Scholar]
  38. Hogerheijde, M. R., van Dishoeck, E. F., Salverda, J. M., & Blake, G. A. 1999, ApJ, 513, 350 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  39. Ivezić, Z., Nenkova, M., & Elitzur, M. 1999, User Manual for DUSTY, Univ. of Kentucky Internal Report, http://www.pa.uky.edu/~moshe/dusty [Google Scholar]
  40. Johnstone, D., Hendricks, B., Herczeg, G. J., & Bruderer, S. 2013, ApJ, 765, 133 [NASA ADS] [CrossRef] [Google Scholar]
  41. Jørgensen, J. K. 2004, A&A, 424, 589 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Jørgensen, J. K., Schöier, F. L., & van Dishoeck, E. F. 2002, A&A, 389, 908 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Jørgensen, J. K., Schöier, F. L., & van Dishoeck, E. F. 2004, A&A, 416, 603 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Jørgensen, J. K., Visser, R., Sakai, N., et al. 2013, ApJ, 779, L22 [NASA ADS] [CrossRef] [Google Scholar]
  45. Jørgensen, J. K., Visser, R., Williams, J., & Bergin, E. A. 2015, A&A, accepted [arXiv:1504.02974] [Google Scholar]
  46. Kenyon, S. J., Hartmann, L. W., Strom, K. M., & Strom, S. E. 1990, AJ, 99, 869 [NASA ADS] [CrossRef] [Google Scholar]
  47. Kim, H. J., Evans, II, N. J., Dunham, M. M., et al. 2011, ApJ, 729, 84 [NASA ADS] [CrossRef] [Google Scholar]
  48. Kim, H. J., Evans, II, N. J., Dunham, M. M., Lee, J.-E., & Pontoppidan, K. M. 2012, ApJ, 758, 38 [NASA ADS] [CrossRef] [Google Scholar]
  49. Kristensen, L. E., van Dishoeck, E. F., Bergin, E. A., et al. 2012, A&A, 542, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Kryukova, E., Megeath, S. T., Gutermuth, R. A., et al. 2012, AJ, 144, 31 [NASA ADS] [CrossRef] [Google Scholar]
  51. Lada, C. J. 1999, in The Origin of Stars and Planetary Systems, eds. C. J. Lada, & N. D. Kylafis (Dordrecht: Kluwer), 143 [Google Scholar]
  52. Lee, J.-E. 2007, J. Korean Astron. Soc., 40, 83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Liebhart, A., Güdel, M., Skinner, S. L., & Green, J. 2014, A&A, 570, L11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  54. Mardones, D., Myers, P. C., Tafalla, M., et al. 1997, ApJ, 489, 719 [NASA ADS] [CrossRef] [Google Scholar]
  55. Maret, S., Bergin, E. A., & Lada, C. J. 2006, Nature, 442, 425 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  56. McElroy, D., Walsh, C., Markwick, A. J., et al. 2013, A&A, 550, A36 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  57. Minissale, M. 2014, Ph.D. Thesis, Observatoire de Paris-Meudon & Université de Cergy-Pontoise [Google Scholar]
  58. Morales-Calderón, M., Stauffer, J. R., Hillenbrand, L. A., et al. 2011, ApJ, 733, 50 [NASA ADS] [CrossRef] [Google Scholar]
  59. Noble, J. A., Congiu, E., Dulieu, F., & Fraser, H. J. 2012, MNRAS, 421, 768 [NASA ADS] [Google Scholar]
  60. Öberg, K. I., Fuchs, G. W., Awad, Z., et al. 2007, ApJ, 662, L23 [NASA ADS] [CrossRef] [Google Scholar]
  61. Onishi, T., Mizuno, A., Kawamura, A., Tachihara, K., & Fukui, Y. 2002, ApJ, 575, 950 [NASA ADS] [CrossRef] [Google Scholar]
  62. Ott, S. 2010, in Astronomical Data Analysis Software and Systems XIX, eds. Y. Mizumoto, K.-I. Morita, & M. Ohishi (San Francisco: ASP), ASP Conf. Ser., 434,139 [Google Scholar]
  63. Owen, J. E., & Jacquet, E. 2015, MNRAS, 446, 3285 [NASA ADS] [CrossRef] [Google Scholar]
  64. Padoan, P., Haugbølle, T., & Nordlund, Å. 2014, ApJ, 797, 32 [NASA ADS] [CrossRef] [Google Scholar]
  65. Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1 [CrossRef] [EDP Sciences] [Google Scholar]
  66. Pontoppidan, K. M., Boogert, A. C. A., Fraser, H. J., et al. 2008, ApJ, 678, 1005 [NASA ADS] [CrossRef] [Google Scholar]
  67. Poteet, C. A., Pontoppidan, K. M., Megeath, S. T., et al. 2013, ApJ, 766, 117 [NASA ADS] [CrossRef] [Google Scholar]
  68. Qi, C., Öberg, K. I., Wilner, D. J., et al. 2013, Science, 341, 630 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  69. Raga, A. C., Velázquez, P. F., Cantó, J., & Masciadri, E. 2002, A&A, 395, 647 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  70. Rebull, L. M., Cody, A. M., Covey, K. R., et al. 2014, AJ, 148, 92 [NASA ADS] [CrossRef] [Google Scholar]
  71. Roelfsema, P. R., Helmich, F. P., Teyssier, D., et al. 2012, A&A, 537, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Safron, E. J., Fischer, W. J., Megeath, S. T., et al. 2015, ApJ, 800, L5 [NASA ADS] [CrossRef] [Google Scholar]
  73. San José-García, I., Mottram, J. C., van Dishoeck, E. F., et al. 2015, A&A, submitted [Google Scholar]
  74. Schmalzl, M., Visser, R., Walsh, C., et al. 2014, A&A, 572, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  75. Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F., & Black, J. H. 2005, A&A, 432, 369 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  76. Scholz, A., Froebrich, D., & Wood, K. 2013, MNRAS, 430, 2910 [NASA ADS] [CrossRef] [Google Scholar]
  77. Spaans, M., Hogerheijde, M. R., Mundy, L. G., & van Dishoeck, E. F. 1995, ApJ, 455, L167 [Google Scholar]
  78. Stamatellos, D., Whitworth, A. P., & Hubber, D. A. 2011, ApJ, 730, 32 [NASA ADS] [CrossRef] [Google Scholar]
  79. Stamatellos, D., Whitworth, A. P., & Hubber, D. A. 2012, MNRAS, 427, 1182 [NASA ADS] [CrossRef] [Google Scholar]
  80. Stauffer, J., Cody, A. M., Baglin, A., et al. 2014, AJ, 147, 83 [NASA ADS] [CrossRef] [Google Scholar]
  81. Tielens, A. G. G. M., Tokunaga, A. T., Geballe, T. R., & Baas, F. 1991, ApJ, 381, 181 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  82. Tobin, J. J., Hartmann, L., Chiang, H.-F., et al. 2011, ApJ, 740, 45 [NASA ADS] [CrossRef] [Google Scholar]
  83. Tobin, J. J., Hartmann, L., Chiang, H.-F., et al. 2012, Nature, 492, 83 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  84. van Dishoeck, E. F., Kristensen, L. E., Benz, A. O., et al. 2011, PASP, 123, 138 [NASA ADS] [CrossRef] [Google Scholar]
  85. van Kempen, T. A., van Dishoeck, E. F., Hogerheijde, M. R., & Güsten, R. 2009, A&A, 508, 259 [Google Scholar]
  86. Visser, R., & Bergin, E. A. 2012, ApJ, 754, L18 (Paper I) [NASA ADS] [CrossRef] [Google Scholar]
  87. Visser, R., Doty, S. D., & van Dishoeck, E. F. 2011, A&A, 534, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  88. Vorobyov, E. I., & Basu, S. 2005, ApJ, 633, L137 [NASA ADS] [CrossRef] [Google Scholar]
  89. Vorobyov, E. I., & Basu, S. 2010, ApJ, 719, 1896 [NASA ADS] [CrossRef] [Google Scholar]
  90. Vorobyov, E. I., Baraffe, I., Harries, T., & Chabrier, G. 2013, A&A, 557, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. Watson, W. D., & Salpeter, E. E. 1972, ApJ, 174, 321 [NASA ADS] [CrossRef] [Google Scholar]
  92. Whittet, D. C. B., Cook, A. M., Chiar, J. E., et al. 2009, ApJ, 695, 94 [NASA ADS] [CrossRef] [Google Scholar]
  93. Wilson, T. L. 1999, Rep. Prog. Phys., 62, 143 [Google Scholar]
  94. Yang, B., Stancil, P. C., Balakrishnan, N., & Forrey, R. C. 2010, ApJ, 718, 1062 [NASA ADS] [CrossRef] [Google Scholar]
  95. Yıldız, U. A., Kristensen, L. E., van Dishoeck, E. F., et al. 2012, A&A, 542, A86 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Yıldız, U. A., Kristensen, L. E., van Dishoeck, E. F., et al. 2013, A&A, 556, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  97. Zasowski, G., Kemper, F., Watson, D. M., et al. 2009, ApJ, 694, 459 [NASA ADS] [CrossRef] [Google Scholar]
  98. Zhu, Z., Hartmann, L., & Gammie, C. 2009, ApJ, 694, 1045 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Details of observations

The spectroscopic observations used in this work are summarized in Tables 1 and 2. The following sections offer more details on how the data were obtained and reduced. The HCO+ and N2H+ 6–5 spectra are published here for the first time; all other observations were gathered from the literature. Typical calibration uncertainties are between 10% and 25%.

Appendix A.1: C18O

Yıldız et al. (2013) compiled all available C18O data from J = 2–1 up to 10–9, with one omission: an integrated intensity of 11.0 K km s-1 in a 22.9′′ beam for the 2–1 line in Ser SMM4 (Hogerheijde et al. 1999). The 6–5 and higher lines are detected in only a handful of sources and are therefore of limited use as an observational diagnostic. The 2–1, 3–2, and 5–4 lines have detection rates of 77%, 100%, and 58% and are well suited for our purpose. Furthermore, with upper-level energies of 16, 32, and 79 K, these three lines together are a good probe of the abundance jump at the CO sublimation front (Yıldız et al. 2013).

Appendix A.2: HCO+

The low-J lines of HCO+ are often optically thick, so we only consider the 6–5 line at 535 GHz (Eu/k = 90 K). This transition was observed with Herschel-HIFI as part of WISH and the related open-time program “Water in low-mass protostars: the William Herschel line legacy” (WILL; proposal code OT2_evandish_4; Mottram et al., in prep.). Five of the spectra were presented by San José-García et al. (2015) and Benz et al. (2015). The other eleven HCO+ lines from this data set are as yet unpublished.

Table A.1 lists the observing dates and identification numbers (ObsIDs), including several duplicate or triplicate observations taken in spectral setups of other lines. Five of these come from the open-time program “Searching for the onset of energetic particle irradiation in Class 0 protostars” (OT2_cceccare_4).

Table A.1

Herschel observations of HCO+ and N2H+.

All observations were taken in double beam switch mode. We use the spectra recorded with the wide-band spectrometer (WBS; 1.1 MHz resolution) and perform basic standard processing in HIPE v10.0.0 (Ott 2010), followed by further analysis in CLASS3. The H and V polarization spectra were averaged after individual inspection and the intensities were converted to main-beam temperature scale through a main-beam efficiency of 0.72 (Roelfsema et al. 2012). Where available, the multiple epochs of data all match to within the calibration uncertainty of 10% and are averaged to reduce the noise. This results in final rms noise levels of 3–7 mK in 0.28 km s-1 bins. The final data reduction step was to subtract a linear baseline.

In all but two sources, the HCO+ 6–5 spectrum shows two emission features: a broad component with a full width at half maximum (FWHM) of 4–14 km s-1 and a narrow component of 1–2 km s-1. These components are also detected in various CO lines and correspond to the outflow (broad) and the quiescent envelope (narrow; Yıldız et al. 2013; San José-García et al. 2015). The two exceptions in our sample are L723 and B335, which only have the narrow component.

Since we are interested exclusively in the quiescent envelope, we fit two Gaussian profiles to each HCO+ spectrum (except L723 and B335) and subtract the outflow emission. Figure A.1 shows the resulting spectra, corrected for the source velocities from Yıldız et al. (2012). Table 1 lists the integrated intensities for the full sample.

thumbnail Fig. A.1

Spectra of HCO+ 6–5, corrected for source velocity and scaled as indicated. If a broad emission component from the outflow was present (in all sources but L723 and B335), it was fitted with a Gaussian profile and subtracted from the spectrum to leave only the narrow emission from the quiescent envelope.

Appendix A.3: H13CO+

Complementing the HCO+ 6–5 line are the 1–0, 3–2, and 4–3 lines of optically thin H13CO+. Their upper-level energies are 4, 25, and 42 K, again offering a good probe of chemical changes induced by the evaporation or freeze-out of CO.

The integrated intensities compiled in Table 1 originate mostly from the ground-based surveys of Hogerheijde et al. (1999) and Jørgensen et al. (2004), with other references as indicated. Some sources have been targeted with multiple telescopes, in which case we choose the observation with the smallest single-dish beam.

Appendix A.4: N2H+

The final two molecular lines are the J = 1–0 and 6–5 transitions of N2H+ (Eu/k = 4 and 94 K), both of which are prone to hyperfine splitting (Caselli et al. 1995). The quantum numbers for the HF levels are F1 and F. The 1–0 transition breaks up into one isolated line (F1 = 0–1) and two clusters of three lines (F1 = 1–1 and 2–1), which are readily resolved with heterodyne spectrographs. The 6–5 transition is dominated by nine HF components within 0.1 km s-1 and therefore appears as a single line in astronomical observations.

Most of the J = 1–0 intensities in Table 1 are taken from Jørgensen et al. (2004) and Emprechtinger et al. (2009), who provided integrated intensities summed over all HF components. Additional data are taken from Mardones et al. (1997), who only showed the spectra of the isolated component and did not tabulate any intensities. We measure the integrated intensities from their Figs. 1 and 2 and multiply by a factor of 9 (based on the Einstein A coefficients and level degeneracies) to approximate the intensity for the full 1–0 band.

The 6–5 lines were observed primarily as part of the Herschel open-time program “The chemistry of nitrogen in dark clouds” (OT1_philybla_1), with two additional sources (L1527 and L1157) covered by OT2_cceccare_4. Both programs targeted N2H+ 6–5 in NGC 1333 IRAS2A and 4A and Ser SMM1. None of these spectra have been published, so we follow the same procedure as outlined for HCO+ 6−5 in Sect. A.2. Table A.1 lists the observing dates and ID codes and Fig. A.2 shows the reduced spectra. The rms noise is 10–14 mK in 0.27 km s-1 bins. N2H+ 6−5 is not detected in L1527. All other sources show a single narrow emission feature with an FWHM of 1–2 km s-1 and no contribution from the outflow. Table 1 lists the integrated intensities for the full sample, including the upper limit for L1527.

thumbnail Fig. A.2

Spectra of N2H+ 6–5, corrected for source velocity and scaled as indicated.

Appendix A.5: Ices

The ice column densities in Table 2 are based on published mid-infrared observations from the ground and from space. The conversion from spectrum to column density requires an absorption band strength measured in the laboratory for the appropriate pure or mixed ices (Gerakines et al. 1995). Most papers on protostellar ices use the same band strengths and we adopt the column densities as published. The two exceptions are the CO2 columns of Zasowski et al. (2009) and the CO column of Tielens et al. (1991), which we scale by a common factor of 1.4 to be consistent with other studies.

All Tables

Table 1

Observed integrated intensities ( in K km s-1) and half-power beam widths (θ in arcsec)a.

In the text
Table 2

Observed ice column densities (N(X) in 1017 cm-2).

In the text
Table 3

Candidate observable ratios to be used as chronometers of episodic accretion timescales.

In the text
Table A.1

Herschel observations of HCO+ and N2H+.

In the text

All Figures

thumbnail Fig. 1

Top panel: radial density profile in IRAS 15398 (black), along with the temperatures in the quiescent phase (blue) and during the burst (red). Other panels: radial abundance profiles of CO, H2O, N2H+, and HCO+. The red curves are during the accretion burst. The black curves are during the quiescent phase, at increasing amounts of time after the burst: 10 yr (dotted), 150 yr (solid), 800 yr (dotted), 2600 yr (short-dashed), 6200 yr (dash-dotted), and 24 000 yr (long-dashed). The quiescent profiles culminate in the blue curve at 1 × 105 yr after the burst.
In the text
thumbnail Fig. 2

Synthetic spectra for the quiescent phase in IRAS 15398, scaled as indicated. The different colors correspond to the same post-burst time steps as in Fig. 1: 10 yr (black), 150 yr (red), 800 yr (gray), 2600 yr (purple), 6200 yr (green), 24 000 yr (orange), and 1 × 105 yr (blue). The arrows indicate how the intensity of each line changes as a function of time; the colored flags mark the peak intensities at the listed times. The N2H+ 1–0 panel shows the isolated F1 = 0−1 hyperfine component, whereas the 6–5 panel shows the sum over all unresolved components.
In the text
thumbnail Fig. 3

Simulated integrated line intensities and ice column densities in IRAS 15398 as a function of time during the quiescent phase. The N2H+J = 1–0 and 6–5 intensities are summed over all HF components. Where available, the gray bars show the observed intensities and columns from Tables 1 and 2 with 20% uncertainties.
In the text
thumbnail Fig. 4

Comparison between model predictions and observations. This plot is similar to Fig. 3, except that the model values are now divided by the observed quantities. A ratio of unity means the model is in exact agreement with observations. If the ratio is larger than unity, the model overpredicts the observation, and vice versa. The horizontal scale for each source runs logarithmically from 10 to 105 yr after the burst. The observable quantity in the bottom three panels is the column density of ice along a pencil beam to the central star. In the other panels, it is the integrated intensity of various rotational lines. The N2H+J = 1–0 and 6–5 intensities are summed over all HF components. The sources are ordered from left to right by increasing luminosity.
In the text
thumbnail Fig. 5

Left: ratio between the integrated intensity of C18O 2–1 and the column density of CO ice as a function of time during the quiescent phase. Each gray curve corresponds to one source. Right: same, but with the time axis normalized by the characteristic freeze-out timescale for each source.
In the text
thumbnail Fig. A.1

Spectra of HCO+ 6–5, corrected for source velocity and scaled as indicated. If a broad emission component from the outflow was present (in all sources but L723 and B335), it was fitted with a Gaussian profile and subtracted from the spectrum to leave only the narrow emission from the quiescent envelope.
In the text
thumbnail Fig. A.2

Spectra of N2H+ 6–5, corrected for source velocity and scaled as indicated.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.