Free Access
Issue
A&A
Volume 613, May 2018
Article Number A29
Number of page(s) 10
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361/201731656
Published online 28 May 2018

© ESO 2018

1 Introduction

Water is probably the molecule that appeals most to our imagination, as it is essential for life as we know it. In star-forming regions, water plays an important role as a coolant of warm gas, a dominant carrier of oxygen, and a major constituent of icy grain mantles deep inside a planet-forming disk (see, e.g., Melnick 2009; Bergin & van Dishoeck 2012; van Dishoeck et al. 2013, for reviews). The transition from water being frozen out onto dust grains to being predominantly present in the gas phase occurs at the water snowline: the midplane radius at which 50% of the water is in the gas phase and 50% is in ice. Sincethe selective freeze-out of the major oxygen carrier alters the elemental C/O ratio in both gas and ice, the bulk chemical composition of planets depends on their formation location with respect to the freeze-out radius of water (e.g., Öberg et al. 2011; Madhusudhan et al. 2014; Ali-Dib et al. 2014; Ali-Dib 2017; Walsh et al. 2015; Mordasini et al. 2016; Eistrup et al. 2016; Booth et al. 2017). In addition, planetesimal formation is thought to be significantly enhanced in this region (e.g., Stevenson & Lunine 1988; Schoonenberg & Ormel 2017). The exact location and time evolution of the water snowline in protostellar systems is thus a crucial ingredient in planet formation.

Unfortunately, the water snowline is very hard to observe directly in protoplanetary disks. Because of the large binding energy of H2 O, the transition from ice to gas happens a few AU from the young star where the midplane temperature exceeds ~100–200 K (depending on the vapor pressure). For the nearest star-forming regions, this would already require angular resolutions of ≲ 0.′′ 01. Emission originating from cold water (Eup < 100 K) has been detected in the disk around TW Hya with the Herschel Space Observatory (Hogerheijde et al. 2011; Zhang et al. 2013), but the large Herschel beam (10–45′′) could not spatially resolve the snowline. Furthermore, the only thermal water lines that can be observed from the ground (except for the H2 O line at 183 GHz) are lines from the less abundant isotopolog HO. As such, even ALMA will have great difficulty locating the water snowline in protoplanetary disks.

An alternative approach is to use chemical imaging, as has been carried out for the CO snowline. Because CO is highly volatile, its snowline is located tens of AU from the central star (around ~20 K; see, e.g., Burke & Brown 2010). Although this is far enough to be spatially resolved with ALMA, locating it directly remains difficult; since CO line emission is generally optically thick, it does not reveal the cold disk midplane. Reactions with gaseous CO are the main destruction route for N2 H+, therefore N2 H+ is expected to be abundant only when CO is frozen out. Taking simple chemical considerations into account, an upper limit for the CO snowline location can therefore be derived from N2 H+ emission (Aikawa et al. 2015; van ’t Hoff et al. 2017), as has been performed for the disks around TW Hya and HD 163296 (Qi et al. 2013, 2015), and several protostellar envelopes (Jørgensen 2004; Anderl et al. 2016).

For young protostars, the emission of complex organic molecules (COMs) has been used as tracer of the inner region that is warm enough to sublimate water ice, the so-called hot core, because these species are expected to be trapped in water ice. However, the exact binding energy, and thus spatial extent, differs for different molecules. Compact COM emission has been detected toward NGC 1333-IRAS2A (hereafter IRAS2A) extending between ~ 0.′′ 4 and ~ 1.′′ 0 (e.g., Jørgensen et al. 2005; Maret et al. 2014; Maury et al. 2014), similar to the size (~ 0.′′8) of the HO emission (Persson et al. 2012). The exact relationship between COMs and the water snowline remains unclear, however. For example, the extent of methanol emission in IRAS2A is half the size of the HO emission (Maret et al. 2014), while it is twice the size in NGC 1333-IRAS4A (Anderl et al. 2016). In addition, Persson et al. (2012) concluded that C2 H5CN is likely not related to sublimation of water ice. Since the detection of COMs is very difficult in mature protoplanetary disks, where only methyl cyanide and methanol have been observed so far (Öberg et al. 2015; Walsh et al. 2016), COMs are not suited as tracers of the water snowline in these systems.

The best candidate to chemically trace the water snowline is HCO+ because its most abundant destroyer in warm dense gas is gaseous H2O, that is, (1)

A strong decline in HCO+ is thus expected when H2O desorbs off the dust grains (Phillips et al. 1992; Bergin et al. 1998). Atacama Large Millimeter/submillimeter Array (ALMA) observations of the optically thin isotopolog H13CO+ toward the Class 0 protostar IRAS 15398–3359 indeed revealed ring-shaped emission (Jørgensen et al. 2013). The spatial distribution is consistent with destruction by water in the innermost region, but the inner radius of the H13CO+ emission is further out than expected. This can be explained if the temperature has recently been higher, that is, if the source has undergone a luminosity outburst. Follow-up observations by Bjerkeli et al. (2016) did not detect the HO 414 − 321 high excitation transition (Eup = 322 K). However, an HDO transition with lower upper level energy (Eup = 22 K) was clearly detected. Although HDO emission is also present in the outflow lobes, the observations are consistent with the H2 O-HCO+ anticorrelation scenario.

Protostellar envelopes are good targets to test the H13CO+-HO anticorrelation. Because of higher luminosity (due to higher accretion rates) and lower vapor pressure, the snowline is located further away from the star than in protoplanetary disks (10 s – 100 AU instead of a few AU; Harsono et al. 2015; Cieza et al. 2016). In addition, HO has already been observed toward four of these objects (Jørgensen & van Dishoeck 2010; Persson et al. 2012, 2013). The only thing lacking is thus high-spatial resolution images of H13CO+. Here, we present NOEMA observations of H13CO+ toward one source, IRAS2A, and compare these to the HO emission presented by Persson et al. (2012). IRAS2A is a deeply embedded Class 0 protostar in the NGC 1333 region of the Perseus molecular cloud. The quadruple outflows suggest that it is a close binary (Jørgensen et al. 2004a), which has recently been confirmed (0.′′ 6 separation; Tobin et al. 2015). Similar to Persson et al. (2012) we adopt a distance of 250 pc (Enoch et al. 2006).

The observations and results are presented in Sects. 2 and 3, and compared with the HO observations from Persson et al. (2012). In Sect. 4, the integrated H13CO+ emission is analyzed using 1D radiative transfer modeling of a spherical envelope with a parametrized H13CO+ abundance profile. In addition, this abundance profile is validated against the outcome of a full chemical network and low-resolution DCO+ observations. The main conclusions, including that we confirm the predicted anticorrelation between H13CO+ and HO, are summarized in Sect. 5.

2 Observations

IRAS2A (α(2000) = 03h28m55.s58; δ(2000) = 31°14′37.′′10) was observed using the NOrthern Extended Millimeter Array (NOEMA) on December 1, 2015 (C configuration), and April 9, 2016 (B configuration) for a total of 2.3 h on source in B configuration and 1.2 h in C configuration. Combining the observations in the different configurations, the data cover baselines from 16.8 to 456.9 m (14.5 to 396.6 kλ). The receivers were tuned to the H13CO+ J = 3 − 2 transition at 260.255 GHz (1.15 mm), and the narrowband correlator was set up with one unit (bandwidth of 40 MHz) centered on this frequency providing a spectral resolution of 0.078 MHz (0.09 km s−1). In addition, the WideX correlator was used, covering a 3.6 GHz window (259.2–262.8 GHz) at a resolution of 1.95 MHz (2.2–2.3 km s−1).

Calibration and imaging were performed using the CLIC and MAPPING packages of the IRAM GILDAS software1. The standard calibration procedure was followed using the quasars 3C454.3 and 3C84 to calibrate the bandpass, 0333+321 to calibrate the complex gains, and the point sources MWC349 and LkHα101 to calibrate theabsolute flux scale. The continuum visibilities were constructed from line-free channels in the WideX spectrum, and the continuum was subtracted from the line data before imaging. Both the line and continuum data were imaged using default robust weighting. For H13CO+ this resulted in a 0.′′ 93×0.′′68 (PA = 36°) beam, which is comparable to the beam of the HO observations (see below), and a rms of 10 mJy beam−1 in 0.09 km s−1 channels. The continuum image has a 0.′′ 95×0.′′65 (PA = 27°) beam. The continuum flux obtained from a Gaussian fit to the visibilities is 347 mJy. This is consistent with previous observations between 0.8 and 1.47 mm (Jørgensen et al. 2007; Persson et al. 2012; Taquet et al. 2015)within the 20% flux calibration error, assuming a power-law spectrum for thermal dust emission (Fννα, with α ~2.0 − 2.4). The (u, v)-fitted continuum peak position (α(2000) = 03h28m55. s57; δ(2000)  =  31°14′36.′′92) is located ~ 0.′′1 south of the bright component of the 0.′′6 binary identified by Tobin et al. (2015).

The HO 31,3 − 22,0 transition at 203.408 GHz (1.47 mm) was observed using the Plateau de Bure Interferometer (PdBI) in December 2009 and March 2010 and was presented by Persson et al. (2012). The data have a spectral resolution of 0.087 MHz (0.115 km s−1) and a resulting beam size of 0.′′87 × 0.′′72 (PA = 63.5°) using natural weighting.

In addition, IRAS2A was observed using the PdBI in July, August, and November 2010 and March 2011 in the C and D configurations as part of a study of complex organic molecules toward low-mass protostars (Taquet et al. 2015). The WideX backends were used at ~145 GHz providing a bandwidth of 3.6 GHz with a spectral resolution of 1.95 MHz (~   3.5 − 4 km s−1). One of the targeted lines was the DCO+ J = 2 − 1 transition at 144.068 GHz. Phase and amplitude were calibrated by performing regular observations of the nearby point sources 3C454.3, 3C84, and 0333+321. Imaging with natural weighting resulted in a 2.′′ 1× 1.′′7 (PA = −155°) beam.

An overview of the observed molecular lines and their parameters is provided in Table 1.

thumbnail Fig. 1

Channel maps for the H13CO+ J = 3 − 2 transition. Black contours are in steps of 10σ, starting from a 3σ level (30 mJy beam−1 in 0.09 km s−1 channels), and the −3σ contours are shown in blue. The peak intensity is 324 mJy beam−1 at 7.68 km s−1 (vlsr + 0.68 km s−1). The continuum peak position is indicated with a cross and the beam is shown in the bottom left corner of the bottom left panel. Channel velocities with respect to the systemic velocity of vlsr ≈ 7.0 km s−1 are listed in the top right corner of each panel. For the adopted distance to IRAS2A, 1′′ corresponds to 250 AU.

Open with DEXTER
Table 1

Overview of the molecular line observations toward IRAS2A.

thumbnail Fig. 2

Integrated intensity map for the H13CO+ J = 3 − 2 transition (color scale) toward IRAS2A. The 1.2 mm continuum is overlaid in white contours (top panel) and the HO 31,3 − 22,0 transition in blue contours (middle panel). The continuum contours are 1.8 (1σ) × [3, 10, 25, 50, 75] mJy beam−1, and the HO contours are 9.8 (1σ) × [3, 8, 15, 25, 35] mJy beam−1 km s−1. The beams are depicted in the bottom left corner. The position of the continuum peak is indicated by a black cross (the close binaryis unresolved in these data) and the outflow axes by red and blue arrows. The integrated intensities along the dashed white line are shown in the bottom panel, normalized to their maximum value. The zero flux level and source position are indicated by dotted lines.

Open with DEXTER
thumbnail Fig. 3

Integrated intensity maps for a spherically symmetric ring-shaped H13CO+ distribution. Emission simulated with the radiative transfer code Ratran is convolved with the observed beam (left panel), and simulated with the (u, v) coverage of the NOEMA observations (right panel). The southern peak in the simulations is probably not observed owing to presence of water in the outflow (see text for details). The continuum peak is denoted with a cross, and the beam is shown in the bottom left corner.

Open with DEXTER
thumbnail Fig. 4

Spectra extracted in a beam toward the continuum peak position (left panels) and the H13CO+ peak position (~2′′ NE of the continuum peak, right panels). The top panels show the H13CO+ J = 3 − 2 transition and the bottom panels the HO 31,3 − 22,0 transition. The two spectra for each species are shown on the same vertical scale. The systemic velocity is indicated by the vertical dotted line at vlsr ≈ 7 km s−1.

Open with DEXTER

3 Results

The H13CO+ channel maps presented in Fig. 1 show detection of H13CO+ at velocities between 6.2 and 8.4 km s−1 (vlsr = 7.0 km s−1; Persson et al. 2012), although the emission is brightest for the redshifted velocities. In none of the channels does the H13CO+ emission peak at the source position. Instead, emission peaks are observed at offsets (RA, Dec) of approximately (2′′,1.5′′) and (–3′′,–3′′) in the redshifted channels, and at (–1′′,–6′′) in the blueshifted channels. HO, on the other hand, shows compact emission on source (located by Persson et al. 2012, ~ 0.′′1 southwest ofthe continuum peak) and extends along the southern outflow axis. Overlaying the integrated intensity maps for H13CO+ and HO (see Fig. 2) shows that the H13CO+ emission surrounds the water emission and peaks ~ 2′′ northeast of the continuum peak.

The asymmetry in the spatial distribution of the H13CO+ emission is partly due to the sampling of the (u, v) plane; simulating emission from a spherical envelope model with a ring-like H13CO+ abundance distribution using the same (u, v) coverage as the observations results in emission peaks in the northeast and southwest rather than a spherically symmetric emission pattern (see Fig. 3). The absence of a strong southwestern peak in the data may be attributed to the presence of water in the southern outflow. Alternatively, it could be an effect of the 3D structure of the source, which causes the H13CO+ emission from that part to be blocked. Finally, it could perhaps be the result of asymmetries in the initial filamentary structure at large scales (≫10 000 AU).

Spectra for H13CO+ and HO extracted toward the continuum and H13CO+ peak positions are presented in Fig. 4. Although there is still H13CO+ emission present on source, it is reduced by ~40% compared to the region in which water emission is absent. The narrow width (FWHM  ≈ 0.8 km s−1) of the H13CO+ line at its peak position indicates that the emission does not originate in the outflow. In comparison, the HO line width is ~3.5 km s−1. In addition, it has been shown in Persson et al. (2012) that the outflow component of the HO emission can be distinguished spatially from the compact component. The emission around vlsr   =  3 km s−1 in the HO spectrum is due to dimethyl ether (CH3OCH3; Persson et al. 2012).

The absorption feature for H13CO+ around the systemic velocity is partly because emission from the large-scale envelope is resolved out, as can be seen from the channel maps in Fig. 1 where most of the emission is absent around the systemic velocity. In addition, a comparison with JCMT single dish observations from Jørgensen et al. (2004b) shows that only ~7% of the peak flux is recovered in the NOEMA observations, demonstrating that the bulk of the H13CO+ is located on larger scales. Imaging the data before continuum subtraction shows that some absorption also occurs against the continuum in the central channels, which is indicative of cold gas along the line of sight at scales larger than the scale of interest.

4 Analysis and discussion

4.1 Parametrized abundance profile for H13 CO+

To establish the origin of the dip in the H13CO+ emission at the central position, 1D spherically symmetric physical-chemical modeling was performed using the radiative transfer code Ratran (Hogerheijde & van der Tak 2000). For the physical structure we adopted the temperature and density profiles for IRAS2A from Kristensen et al. (2012), derived using the 1D spherically symmetric dust radiative transfer code DUSTY (Ivezić & Elitzur 1997). In this procedure the free model parameters were fitted to the spectral energy distribution (SED) and the spatial extent of the sub-mm continuum emission. Based on the chemical consideration that H13CO+ is particularly abundant when its main destructor H2 O is frozen out, the simplest H13CO+ distribution is an abundance profile with a low abundance inside the 100 K radius (i.e., the water snowline) and a higher abundance at larger radii. Both the inner and outer abundances are varied to get the best match to the observed integrated intensity along the radial cut shown in Fig 2. The absence of the southwest H13CO+ peak (discussed in Sect. 3) cannot be modeled with an axisymmetric 1D model. We therefore restrict the comparison to the northeastern part of the radial cut. For a first analysis, we do not adopt a velocity structure. The difference between the blue- and redshifted channels then cannot be reproduced exactly. Since there is almost no emission along the radial cut in the blueshifted channels, only the redshifted channels are used to construct the zeroth moment maps. Including all channels does not significantly change the observed integrated radial profile, but only affects the H13CO+ abundance required to reach the observed intensity (by a factor of ~3). However, deriving the precise H13CO+ abundance is not the goal of this work. Finally, simulated images representing the observed (u, v) coverage are produced from the Ratran images in GILDAS by creating a (u, v) table with the function uv_fmodel.

Various H13CO+ abundance profiles with their corresponding integrated emission profiles are presented in Fig. 5. For the adopted physical structure, the temperature drops below 100 K around 90 AU. An H13CO+ abundance of 1 × 10−10 at radii >90 AU can roughly reproduce the observed integrated intensity (Fig. 5, second row). However, independentof the inner abundance, the emission always peaks on source, as in the case of a constant H13CO+ abundance at all radii (Fig. 5, top row). Although a large percentage (~90%) of the emission is already being resolved out with the (u, v) coverage of the observations, this thus does not result in a dip in the H13CO+ emission toward the source position.

To reproduce the location of the emission peak, the abundance jump has to be shifted outward to a radius of ~360 AU. To move the100 K radius, representing the water snowline, this far out, the temperature profile has to be multiplied by a factor of 2.2 (Fig. 5, third row). However, the observed integrated intensity cannot be reached with this abundance profile, as the H13CO+ emission becomes optically thick for outer abundances higher than 3 × 10−10. Varying the abundance does not affect the position of the emission peak.

The best match to the observed H13CO+ radial profile can be reached by introducing an additional drop in the abundance around ~4000 AU (Fig. 5, fourth row). An increase in abundance of at least a factor of six is required around 360 AU to get the emission peak at the observed location, but the exact radius of the outer drop does not influence the peak position. Such a decrease in the H13CO+ abundance can be expected around the CO snowline because formation of HCO+ is limited in this region because of a reduction in the parent gaseous CO molecules. A radius of 4000 AU corresponds to 17 K for the Kristensen et al. (2012) temperature profile, consistent with the CO freeze-out temperature for a pure CO ice binding energy. In contrast, for the increased temperature profile the temperature is 38 K at this radius; temperatures below 20 K are not reached within the 20 000 AU radial extent of the model. An increased temperature is thus required to set the inner radius of the H13CO+ abundance at the water snowline, while the original temperature places the outer radius around the CO snowline. This discrepancy is further discussed in Sects. 4.2 and 4.3. The main result of these models is that the observed position of the H13CO+ emission peak can only be reproduced by a reduction in the abundance of at least a factor 6 in the inner ~360 AU.

thumbnail Fig. 5

Various H13CO+ abundance profiles (left panels) and the resulting simulated integrated emission along the northeastern part of the radial cut shown in Fig. 2 (right panels). The observed emission is shown in purple. The location of the abundance jumps (in Kelvin and AU) is indicated in the top right corners. Bottom panel: temperature (gray) and density (blue) profiles for the IRAS2A envelope from Kristensen et al. (2012) used in the top two models. The location of the H2 O and CO snowlines (at 100 K and 20 K, resp.) are indicated by the dotted red and blue lines, respectively. The temperature profile increased by a factor of 2.2, used in the bottom two models, is shown in black. The resulting H2 O snowline is indicated by the dashed red line while the CO snowline now falls outside the adopted radial range. For IRAS2A, 1′′ corresponds to ~250 AU.

Open with DEXTER
thumbnail Fig. 6

Left panels: temperature (gray) and density (blue) profiles for the IRAS2A envelope from Kristensen et al. (2012) (top), and the corresponding H2O (light blue), CO (blue), and HCO+ (purple) abundances predicted by the GRAINOBLE model (bottom). The dashed green line shows the empirically inferred abundance profile for H13CO+ (scaled down by a factor 100). Right panels: as left panels, but with the temperature increased by a factor 2.2. The vertical dotted lines indicate the H2 O and CO snowlines.

Open with DEXTER

4.2 Comparison with a full chemical model

To inspect the legitimacy of the simple abundance profile, we compared it to the outcome of the astrochemical code GRAINOBLE (Taquet et al. 2012a,b, 2013, 2014). In short, GRAINOBLE uses the rate equations approach (Hasegawa et al. 1992) considering three types of reactions (a so-called three-phase model based on Hasegawa & Herbst 1993): reactions for gas-phase species, reactions for species on the ice surface, and reactions for bulk ice species. The gas-phase chemical network is taken from the KIDA database (Wakelam et al. 2015) supplemented by ion-neutral reactions described in Taquet et al. (2016). A cosmic ray (CR) ionization rate of 5 × 10−17 s−1 is used. In addition to thermal desorption, UV photodesorption and chemical desorption upon formation of an exothermic surface reaction have been included as desorption processes. The surface chemical network follows the results of recent laboratory experiments to form the main ice species (see Taquet et al. 2013, for references). For CO and H2 O binding energies have been adopted for an amorphous water ice substrate: 1150 K and 5775 K, respectively (Collings et al. 2004; Fraser et al. 2001). The formation of interstellar ices is first followed during the formation of a dense starless core with the method described in Taquet et al. (2014). This ice chemical composition is then usedas initial conditions for the gas-grain chemical evolution around IRAS2A, adopting the static 1D physical profile derived by Kristensen et al. (2012). Surface reactions have a negligible impact on the HCO+ abundance during the latter protostellar phase. They have a significant importance during the former prestellar phase since they govern the transformation of CO ice into other species, thereby setting the initial ice abundance for the protostellar phase. The initial gas-phase abundances and ice abundances at the start of the protostellar phase are listed in Appendix A.

The left panel of Fig. 6 presents the H2 O, CO and HCO+ abundances for the IRAS2A envelope structure after 104 yr of chemical evolution. The abundance profiles do not change significantly for evolution up to ~105 yr, that is, the duration of the protostellar phase. Going from the outer envelope inward, the temperature rises above 20 K at a radius of ~2700 AU, causing CO to desorb from the ice surface into the gas phase. In the three-phase model, a significant part of the icy CO remains trapped in the ice matrix at temperatures higher than the CO desorption temperature. These CO molecules evaporate together with water when the temperature exceeds ~100 K, resulting in a double abundance jump profile for CO. The water snowline, that is, the radius where 50% of the water is in the gas phase and 50% is frozen out, is located at a distance of ~65 AU from the protostar (corresponding to 127 K). As expected, the HCO+ abundance increases sharply at radii beyond the water snowline, where its main destructor is frozen out, and decreases again outside the CO snowline, where the availability of its parent molecule is reduced. Increasing the temperature by a factor of 2.2, as was required to match the observed H13CO+ emission with a drop in the abundance at the 100 K radius, results in a similar chemical profile, only with the snowlines shifted to larger radii (Fig. 6, right panel). The water snowline is now located at ~225 AU (127 K), while the CO snowline falls outside the modeled radial range. The water snowline location matches the observed radial extent of the HO emission (~1′′ ≈ 250 AU).

For comparison, the best H13CO+ abundance profile derived from the data is added to Fig. 6 (dashed green line). This inferred H13CO+ abundance is higher than expected based on the GRAINOBLE results: similar to the predicted HCO+ abundance instead of a factor ~70 lower. This can most likely be explained by the relatively low CO abundance between the CO and H2 O snowlines in the chemical model, which directly results in a low HCO+ abundance. The CO abundance is low because a large percentage of the CO ice is converted into methanol during the prestellar phase (Table A.2). In addition, in the three-phase model, not all CO evaporates inside its snowline, but remains trapped in water ice (Fig. 6). However, it is the overall trend that the HCO+ abundance sharply rises at the H2 O snowline and drops at the CO snowline that is important.

The radius where H13CO+ has to increase (360 AU) matches the inner radius of the HCO+ distribution in the chemical network with the increased temperature profile (right panel). The water snowline, however, is located ~135 AU closer to the star, at 225 AU. Similar behavior is seen in models for N2 H+; the 50% reduction in CO, or in this case H2O, at the snowline is not yet sufficient to diminish the destruction of N2 H+ or HCO+, respectively, moving the emission ring outward (Aikawa et al. 2015; van ’t Hoff et al. 2017). There are thus two effects at play here: first, the HCO+ abundance is expected to peak outside rather than directly at the water snowline (360 vs. 225 AU). Second, this location for the water snowline is further away from the protostar than expected based on the Kristensen et al. (2012) temperature profile (225 vs. 65 AU).

The outer drop in H13CO+ abundance matches less well with the predicted HCO+ distributions. The best agreement is for the original temperature profile (for the adopted CO binding energy of 1150 K), although HCO+ decreases ~2000 AU (2.4 K) further in than the radiative transfer modeling derived from the H13CO+ data. Owingto the shallow gradient of the temperature profile in the outer envelope, the exact location of the CO snowline is very sensitive to small changes in the thermal structure. The mismatch could therefore be due to uncertainties in the derived temperature. Another important parameter for the snowline location is the assumed binding energy for CO, which is related to the assumed ice structure, that is, whether CO ice is mixed with other molecules or pure. Adopting a lower binding energy, corresponding to a CO ice substrate rather than a H2 O ice substrate, would shift the snowline to larger radii. Anderl et al. (2016) investigated the CO binding energy using NOEMA observations of C18O and N2 H+ toward four young protostars and concluded that the extent of the C18 O emission was best reproduced using a CO binding energy of 1200 K, corresponding to a water ice substrate. In addition, these authors showed that the freeze-out temperature is ~6 K lower assuming a binding energy of 855 K for pure CO ice instead of 1200 K. Since the temperature difference between the CO snowline and the location of the H13CO+ drop is only 2.4 K, a high CO binding energy also seems best for IRAS2A. The small discrepancy between the CO snowline location predicted by the chemical model and the location inferred from the H13CO+ observations may then be the result of asymmetries in IRAS2A, as hinted at by the asymmetric C18 O emission shown by Anderl et al. (2016) for other protostars.

As noted before, the exact location of the drop in H13CO+ in the outer envelope does not influence the observed position of the emission peak. The main point is that the chemical model indeed predicts a steep rise in HCO+ abundance just outside the water snowline, and a drop outside the CO snowline. This justifies the abundance profile used in the radiative transfer modeling, and corroborates that the depression in H13CO+ in the center is due to the presence of gas-phase water.

4.3 Tracing the CO snowline: signs of an accretion burst?

Having the snowlines located further out than expected based on the adopted temperature profile derived from SED models, such as found here for the water snowline, could be a sign that the protostar has recently undergone a luminosity outburst. Such an outburst heats up the envelope, causing ices to sublimate off the grains. After the burst, the envelope cools rapidly (Johnstone et al. 2013), while it takes much longer for the molecules to freeze back onto the dust grains (Rodgers & Charnley 2003). As a result, snowlines are shifted away from the star (Visser et al. 2015; Jørgensen et al. 2015; Cieza et al. 2016; Frimann et al. 2017). An increased temperature structure is required to fit the H13CO+ emission, suggesting that the water snowline is further out than expected. In case of an outburst, the CO snowline would also be shifted outward. Observational constraints for the CO snowline can thus help to establish whether IRAS2A has indeed recently undergone an accretion burst.

The deuterated form of HCO+, DCO+, can be used to trace the freeze-out of CO (e.g., Mathews et al. 2013). DCO+ forms via reaction of H2 D+ and CO. Formation of H2 D+ is enhanced at low temperatures (<30 K), but below ~20 K CO starts to freeze-out. This balance results in DCO+ peaking around the CO snowline (Murillo et al. 2015). Contamination of the emission by DCO+ formed in warmer layers from CH2D+ is not likely to be important in protostellar envelopes, unlike in protoplanetary disks (Favre et al. 2015). The DCO+ abundance profile for the IRAS2A envelope structure is predicted with GRAINOBLE by extending the chemical network with a small network that includes reactions between CO and H2 D+, HD and D. The DCO+ abundance indeed peaks around the CO snowline (see Fig. 7), which is located at ~2000 AU for the fiducial temperaturestructure. As such, DCO+ emission can be used to estimate the CO snowline location and constrain the temperature profile in the envelope.

The observed integrated intensity map for DCO+ toward IRAS2A is presented in Fig. 8. Emission is observed south of the continuum peak and peaks ~10′′ off source. Such an asymmetric emission pattern is observed for more protostars (e.g., Murillo et al. 2015). A CO snowline around ~2000 AU, as predictedby the chemical model for the Kristensen et al. (2012) temperature profile, is in agreement with the DCO+ emission. Thus, although the H2O snowline is located at larger radii than predicted based on the IRAS2A temperature profile (225 vs. 65 AU), the location of the CO snowline is as expected. Since the chemical reset after an outburst is faster at high densities, it is unlikely that the CO snowline has already shifted back while the H2 O snowline has not. In addition, Jørgensen et al. (2015) presented C18O images obtained with the Submillimeter Array (SMA) for which the spatial extent matches well with the current luminosity. It is thus unlikely that IRAS2A has recently undergone an accretion burst. Why the temperature profile matches the snowline location for CO but not for H2O is further discussed in Sect. 4.5.

thumbnail Fig. 7

H2O (light blue), CO (blue), HCO+ (purple), and DCO+ (pink) abundances predicted by the GRAINOBLE model for the IRAS2A temperature and density profile from Kristensen et al. (2012). The vertical dotted lines indicate the H2 O and CO snowlines.

Open with DEXTER
thumbnail Fig. 8

Integrated intensity map for the H13CO+ J = 3 − 2 transition (color scale) toward IRAS2A. The DCO+ J = 2 − 1 transition is overlaid in pink contours. Contours are shown in steps of 1σ starting from a 3σ level (15 mJy beam−1 km s−1). The beams are depicted in the bottom left corner. The position of the continuum peak is indicated by a white cross and the outflow axes by red and blue arrows. The length of the arrows is the same as in Fig. 2 to compare the scales. The solid white contour represents the CO snowline radius predicted by the chemical model for the Kristensen et al. (2012) temperature profile (~2000 AU).

Open with DEXTER

4.4 Attenuation of the cosmic ray ionization rate

The chemistry leading to the formation of HCO+ is initiated by the ionization of H2 molecules by cosmic rays (CRs). A reduction in the CR ionization rate therefore results in a lower HCO+ abundance. Numerical simulations performed by Padovani et al. (2013) showed that the CR ionization rate can be reduced by an order of magnitude in the inner region of a collapsing cloud core owing to twisting of the magnetic field lines. However, the exact value of the CR ionization rate in the inner regions of protostellar envelopes is unknown, and more recent work by Padovani et al. (2015, 2016) has shown that jet shocks and protostellar surface shocks can drive particle acceleration, leading to increased CR ionization rates.

Nonetheless, to assess if the effect of CR attenuation is large enough to cause the drop in H13CO+ emission that we observe, we performed an order of magnitude calculation. As such, we approximate the reactions producing H upon CR ionization by (2)

where ζCR is the CR ionization rate. In this formalism we ignore H, as the ionization is the rate limiting step for H. Formation of HCO+ is dominated by the reaction between H and CO, (3)

and to estimate the importance of CRs, only recombination is considered as destruction reaction for HCO+, that is, (4)

Assuming steady state and n(e) ≈ n(HCO+), the following relation can then be derived between the HCO+ density, n(HCO+), and the CR ionization rate, ζCR (see Appendix B for more details): (5)

where n(H2) is the molecular density of H2 and k is the reaction rate coefficient for the HCO+ destruction reaction (Eq. (4)). From this it can be seen that an order of magnitude reduction in the CR ionization rate gives a factor of three reduction in HCO+ abundance. Based on our radiative transfer modeling, at least a factor of six is required to fit the observations. Attenuation of the CR ionization rate alone is thus not likely to be the cause of the observed reduced H13CO+ emission toward the center of IRAS2A.

4.5 Physical structure

Finally, a more likely explanation for the discrepancy between the temperatures needed for the CO and H2 O snowlines may be that the physical structure of IRAS2A is more complicated than assumed with the 1D model. For several young protostars, it has been shown that a spherical envelope model cannot reproduce all of the continuum emission on small scales (e.g., Hogerheijde et al. 1999; Looney et al. 2003). This is also the case for IRAS2A (Brown et al. 2000; Jørgensen et al. 2005) and Persson et al. (2016) have derived the presence a compact disk-like structure with a 167 AU radius. The density may therefore change significantly on these scales of disk-envelope transition, that is, the scale at which we observe the HO-H13CO+ anticorrelation. Since the envelope around IRAS2A is relatively optically thick to its own infrared radiation on ~100 AU scales, the exact temperature distribution becomes very dependent on the assumption about the density structure. In addition, the presence of an outflow cavity can shift the 100 K radius outward in specific directions (Bjerkeli et al. 2016), as the temperature along the cavity wall is higher at certain radii than in the envelope midplane (Visser et al. 2012). However, these considerations do not affect the main conclusion that the observations show a clear anticorrelation between HO and H13CO+ emission.

5 Summary and outlook

We have presented subarcsecond resolution observations of the H13CO+ J = 3 − 2 transition at 260.255 GHz toward the Class 0 protostar IRAS2A, and for the first time, compared these directly to observations of HO, which were presented by Persson et al. (2012). While the HO 31,3 − 22,0 emission is centered at the continuum peak within a ~1′′ radius, the H13CO+ emission peaks ~2′′ to the northeast.

Using a 1D envelope model for IRAS2A, this offset in H13CO+ emission can be explained by an abundance decrease of at least a factor of six in the inner ~360 AU, but not by a constant H13CO+ abundance. An increase in H13CO+ abundance just outside the water snowline is consistent with chemical model predictions; the water snowline is then located around 225 AU, whereas the H13CO+ abundance peaks ~135 AU further out at 360 AU. This snowline radius is larger than expected if the 1D spherically symmetric temperature profile from Kristensen et al. (2012) is adopted. DCO+ observations place the CO snowline at the expected location, making an outburst scenario unlikely.

The combined observations of H13CO+ and HO, together with radiative transfer and chemical modeling, thus provide proof of concept for the HCO+-H2O anticorrelation. The next step would be to corroborate this for more protostellar systems. Establishing H13CO+ as a good tracer of the water snowline may allow localization of the snowline in protoplanetary disks for which direct detection is very difficult.

Acknowledgements

We would like to thank the referee for useful comments. In addition, we are grateful to the IRAM staff, in particular Miguel Montargès, for his help with the observations and reduction of the data, and Nadia Murillo and Arthur Bosman for useful discussions. Astrochemistry in Leiden is supported by the European Union A-ERC grant 291141 CHEMPLAN, by the Netherlands Research School for Astronomy (NOVA) and by a Royal Netherlands Academy of Arts and Sciences (KNAW) professor prize. M.L.R.H acknowledges support from a Huygens fellowship from Leiden University. M.V.P. postdoctoral position is funded by the ERC consolidator grant 614264. V.T. has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement no. 664931. Theresearch of J.K.J. is supported by European Union ERC Consolidator Grant “S4F” (grant agreement No 646908) and Centre for Star and Planet Formation funded by the Danish National Research Foundation.

Appendix A Initial abundances in the GRAINOBLE model

Table A.1 lists the initial elemental abundances in the gas phase for the prestellar phase of the GRAINOBLE model, based on the low-metal model from Wakelam & Herbst (2008). The ice composition at the end of this phase (after 1.7 Myr) is listed in Table A.2. These are the initial conditions for the chemical evolution around IRAS2A.

Table A.1

Initial gas-phase elemental abundances.

Table A.2

Ice chemical composition at the end of the prestellar phase.

Appendix B Effect of an attenuated cosmic ray ionization rate

A first order approximation of the effect of CR ionization rate on the HCO+ abundance can be made by considering the main formation and destruction reactions of HCO+. Assuming steady state, that is, assuming that the HCO+ abundance does not change over time, we can write (B.1)

The main formation pathway of HCO+ is through reaction between CO and H, (B.2)

and recombination with electrons is the main destruction route, (B.3)

The reaction of HCO+ with H2O is ignored in this calculation to assess the importance of the CR ionization rate.

The rate of these reactions can be written as (B.4)

where k is the rate coefficient, and n(A) and n(B) are the molecular densities of the reacting species. Equation (B.1) then becomes (B.5)

The H density, n(H), can be calculated by following the same procedure. The formation reactions are (B.6) (B.7)

In order to solve Eq. (B.5) analytically, we summarize the formation of H as follows: (B.8)

because the ionization reaction is the rate limiting step. In addition, only the reaction with CO (Eq. (B.3)) is considered as destruction mechanism. We can then write an equation similar to Eq. (B.5), that is, (B.9)

Substituting Eq. (B.9) into Eq. (B.5), and assuming n(e) ≈ n(HCO+), we can solve Eq. (B.5) to get the following expression for the HCO+ density: (B.10)

References


All Tables

Table 1

Overview of the molecular line observations toward IRAS2A.

Table A.1

Initial gas-phase elemental abundances.

Table A.2

Ice chemical composition at the end of the prestellar phase.

All Figures

thumbnail Fig. 1

Channel maps for the H13CO+ J = 3 − 2 transition. Black contours are in steps of 10σ, starting from a 3σ level (30 mJy beam−1 in 0.09 km s−1 channels), and the −3σ contours are shown in blue. The peak intensity is 324 mJy beam−1 at 7.68 km s−1 (vlsr + 0.68 km s−1). The continuum peak position is indicated with a cross and the beam is shown in the bottom left corner of the bottom left panel. Channel velocities with respect to the systemic velocity of vlsr ≈ 7.0 km s−1 are listed in the top right corner of each panel. For the adopted distance to IRAS2A, 1′′ corresponds to 250 AU.

Open with DEXTER
In the text
thumbnail Fig. 2

Integrated intensity map for the H13CO+ J = 3 − 2 transition (color scale) toward IRAS2A. The 1.2 mm continuum is overlaid in white contours (top panel) and the HO 31,3 − 22,0 transition in blue contours (middle panel). The continuum contours are 1.8 (1σ) × [3, 10, 25, 50, 75] mJy beam−1, and the HO contours are 9.8 (1σ) × [3, 8, 15, 25, 35] mJy beam−1 km s−1. The beams are depicted in the bottom left corner. The position of the continuum peak is indicated by a black cross (the close binaryis unresolved in these data) and the outflow axes by red and blue arrows. The integrated intensities along the dashed white line are shown in the bottom panel, normalized to their maximum value. The zero flux level and source position are indicated by dotted lines.

Open with DEXTER
In the text
thumbnail Fig. 3

Integrated intensity maps for a spherically symmetric ring-shaped H13CO+ distribution. Emission simulated with the radiative transfer code Ratran is convolved with the observed beam (left panel), and simulated with the (u, v) coverage of the NOEMA observations (right panel). The southern peak in the simulations is probably not observed owing to presence of water in the outflow (see text for details). The continuum peak is denoted with a cross, and the beam is shown in the bottom left corner.

Open with DEXTER
In the text
thumbnail Fig. 4

Spectra extracted in a beam toward the continuum peak position (left panels) and the H13CO+ peak position (~2′′ NE of the continuum peak, right panels). The top panels show the H13CO+ J = 3 − 2 transition and the bottom panels the HO 31,3 − 22,0 transition. The two spectra for each species are shown on the same vertical scale. The systemic velocity is indicated by the vertical dotted line at vlsr ≈ 7 km s−1.

Open with DEXTER
In the text
thumbnail Fig. 5

Various H13CO+ abundance profiles (left panels) and the resulting simulated integrated emission along the northeastern part of the radial cut shown in Fig. 2 (right panels). The observed emission is shown in purple. The location of the abundance jumps (in Kelvin and AU) is indicated in the top right corners. Bottom panel: temperature (gray) and density (blue) profiles for the IRAS2A envelope from Kristensen et al. (2012) used in the top two models. The location of the H2 O and CO snowlines (at 100 K and 20 K, resp.) are indicated by the dotted red and blue lines, respectively. The temperature profile increased by a factor of 2.2, used in the bottom two models, is shown in black. The resulting H2 O snowline is indicated by the dashed red line while the CO snowline now falls outside the adopted radial range. For IRAS2A, 1′′ corresponds to ~250 AU.

Open with DEXTER
In the text
thumbnail Fig. 6

Left panels: temperature (gray) and density (blue) profiles for the IRAS2A envelope from Kristensen et al. (2012) (top), and the corresponding H2O (light blue), CO (blue), and HCO+ (purple) abundances predicted by the GRAINOBLE model (bottom). The dashed green line shows the empirically inferred abundance profile for H13CO+ (scaled down by a factor 100). Right panels: as left panels, but with the temperature increased by a factor 2.2. The vertical dotted lines indicate the H2 O and CO snowlines.

Open with DEXTER
In the text
thumbnail Fig. 7

H2O (light blue), CO (blue), HCO+ (purple), and DCO+ (pink) abundances predicted by the GRAINOBLE model for the IRAS2A temperature and density profile from Kristensen et al. (2012). The vertical dotted lines indicate the H2 O and CO snowlines.

Open with DEXTER
In the text
thumbnail Fig. 8

Integrated intensity map for the H13CO+ J = 3 − 2 transition (color scale) toward IRAS2A. The DCO+ J = 2 − 1 transition is overlaid in pink contours. Contours are shown in steps of 1σ starting from a 3σ level (15 mJy beam−1 km s−1). The beams are depicted in the bottom left corner. The position of the continuum peak is indicated by a white cross and the outflow axes by red and blue arrows. The length of the arrows is the same as in Fig. 2 to compare the scales. The solid white contour represents the CO snowline radius predicted by the chemical model for the Kristensen et al. (2012) temperature profile (~2000 AU).

Open with DEXTER
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.