© ESO, 2015

1. Introduction

Rotationally supported disks have been observed extensively among most protostellar and pre-main sequence evolutionary stages (Li et al. 2014). Recent studies have revealed the existence of such disks in the Class 0 deeply embedded phase (Choi et al. 2010; Tobin et al. 2013; Murillo et al. 2013; Codella et al. 2014; Lee et al. 2014). In the early stages of star formation, the envelope is not yet dispersed and contains enough material to influence the evolution of the star-disk system. The boundary between the disk and the envelope, known as the disk-envelope interface, must then play a role in the formation process. This region, however, is largely unexplored owing to limitations in the resolution and sensitivity, as well as the lack of observed rotationally supported disks, until now. Class 0 sources with confirmed rotationally supported disks grant us the opportunity to study the chemical and physical structure of the disk-envelope interface region, which is crucial for the next step of understanding this region’s role in star formation (Sakai et al. 2014).

Fig. 1 DCO+ compared with continuum and C18O toward VLA1623A and B, marked with a square and triangle, respectively. The black ellipse shows the synthesized beam. Intensity integrated DCO+ (color-scale) with a) 1.3 mm continuum and b) C18O. DCO+ peaks at 0.0433 Jy beam-1 km s-1. Contours are in steps of 3, 5, 10, 15, 20, 40, 60 and 78σ with σ = 1 mJy beam-1 for 1.3 mm continuum, and −5, −3, 3, 5, 10, 15, 20 and 25σ with σ = 13 mJy beam-1 km s-1 for C18O. Velocity map (moment 1, color-scale) for c) DCO+ and d) C18O overlaid with the corresponding intensity integrated map (contours). Contours are in steps of −10, −7, −4, 4, 7, 10 and 11σ with σ = 3 mJy beam-1 km s-1 for DCO+, and the same as in b) for C18O.

Whilst CO isotopologues are good tracers of rotationally supported disks, they freeze out onto dust grains below the evaporation temperature T_ev usually between 30−20 K (Jørgensen et al. 2005). These low temperatures are reached at the edge of the embedded disk (Visser et al. 2009). As a result molecular species whose abundance is enhanced at low temperatures are needed to trace the disk-envelope interface. DCO⁺ emission is known to be optically thin and its abundance is enhanced at a narrow range of temperatures below the CO freeze-out temperature tracing the so-called CO snowline (e.g., Wootten 1987; Roberts et al. 2003; Mathews et al. 2013). Thus, DCO⁺ is a good candidate molecule to trace the chemical, physical and kinematic structure of the disk-envelope interface.

VLA1623-2417 (hereafter VLA1623) is a triple non-coeval protostellar system located in ρ Ophiuchus at d ~ 120 pc (Murillo & Lai 2013). VLA1623A is the prototypical Class 0 source and emits predominantly in the (sub)millimeter range (André et al. 1993; Murillo & Lai 2013). From modeling of ALMA Cycle 0 C¹⁸O observations, Murillo et al. (2013) found that VLA1623A supports a Keplerian disk with a radius of at least 150 AU and a central mass M_∗ of 0.2 M_⊙. The ALMA Cycle 0 observations also detected DCO⁺ molecular line emission toward VLA1623A bordering the C¹⁸O disk. This grants us the opportunity to probe the disk-envelope interface of VLA1623A.

In this paper we present the results of our ALMA observations and simple chemical modeling, aiming to understand the physical structure of the boundary between the envelope and the disk in a Class 0 protostar.

Fig. 2 Intensity (contours) and velocity (color-scale) integrated maps of DCO+ ACA observations. Note the 10 times larger scale of this figure compared with Fig. 1. The positions of VLA1623A, B and W are marked with a square, triangle and circle, respectively. The size of the filled square shows the size of the C18O disk. The red and blue lines indicate the outflow direction while the gray line shows the mayor axis of the disk plane. The black ellipse shows the synthesized beam. Contours are in steps of −3, 3, 5, 7, 10, 15 and 17σ with σ = 220 mJy beam-1 km s-1.

2. Observations

Using the Atacama Large Millimeter/submillimeter Array (ALMA) we observed VLA1623 (pointing coordinates α = 16:26:26.419δ = −24:24:29.988 J2000) during the early science Cycle 0 period on 8 April 2012. Observations were carried out using the extended configuration, comprised of 16 antennas with a maximum baseline of ~400 m, in Band 6 (230 GHz). Total observing time was one hour with a 58% duty cycle. The spectral configuration was set to observe DCO⁺ (3−2), C¹⁸O (2−1), and ¹²CO (2−1) in addition to continuum. N₂D⁺ was also observed but no significant detection was made. Data calibration was done with 3C 279, 1733-130, and Titan for bandpass, gain, and flux calibration, respectively.

In this paper, we present the results and analysis of the DCO⁺ (3−2) (rest frequency: 216.11258 GHz) observations. The DCO⁺ data were calibrated jointly with the continuum, C¹⁸O (2−1), and ¹²CO (2−1) data. Further calibration details and results from the other observed lines can be found in Murillo et al. (2013) and Santangelo et al. (2015). The spectral set-up provided a velocity resolution of 0.0847 km s^-1. The synthesized beam size for the DCO⁺ images is 0.85″ × 0.65″ with PA = 96°. The rms noise of the channel map is 12 mJy beam^-1 for a spectral resolution of 0.0847 km s^-1, giving a peak S/N = 7.

Our ALMA observations provide a maximum scale of 4″ and a field of view (FOV) of 24″, with emission between 4″ and 24″ largely filtered out by the interferometer. The FOV together with the beam size of 0.85″ constrain the scale to which any analysis of the data can be done.

In addition to the ALMA Cycle 0 observations, we present the DCO⁺ (3−2) results from our ALMA Cycle 2 Atacama Compact Array (ACA) observations carried out on 7 August 2014 (pointing coordinates α = 16:26:26.390δ = −24:24:30.688). Total observing time was 2 hours. Data calibration was done with J1517-243 and Mars for flux, J1625-2527 and Mars for gain, and J1733-1304 for bandpass. The rms noise of the DCO⁺ channel map is 73 mJy beam^-1 for a spectral resolution of 0.021 km s^-1, with a synthesized beam of 8.6″ × 4.2″ with PA = −76°. These observations, with a mosaicked area of 6′, provide the DCO⁺ emission between 4″ and 18″ scales.

Fig. 3 Cycle 0 DCO+ (black) and C18O (gray) PV diagram. Cuts are made at PA = 26° for DCO+ and PA = 35° for C18O. Keplerian (solid line) and Infall plus Keplerian (dashed line) rotation PV models with M∗ = 0.2 M⊙ are overplotted (Murillo et al. 2013). Contours are in steps of −3, 3, 5, 10 and 15σ with σ = 19 mJy beam-1 for C18O, −7, −5, −3, −2, 2, 3, 4, 5σ with σ = 12 mJy beam-1 for Cycle 0 DCO+.

3. Results

3.1. ALMA 12-m Array Cycle 0

The detected DCO⁺ emission (Fig. 1) is located between a velocity range of 2.8 to 5.2 km s^-1 and shows a clumpy structure with two main clumps to the north and south of VLA1623A. The southern clump emission is stronger and is offset by about 2.5″ (300 AU) from VLA1623A and in addition it borders the red-shifted emission of the disk traced in C¹⁸O. The northern clump slightly overlaps the blue-shifted emission of the C¹⁸O disk and borders VLA1623B’s continuum emission. A couple of clumps with emission between 3 and 10σ are observed near the continuum peaks of VLA1623A and B. No significant emission was detected toward VLA1623W, separated by 10″ to the west from VLA1623A, possibly either because of the lack of DCO⁺ or because the emission is too weak and filtered out.

The C¹⁸O line emission tracing the disk was found to be influenced by the outer envelope, with the blue-shifted emission being affected more than the red-shifted emission (Murillo et al. 2013). We expect the same to hold true for the DCO⁺ emission, thus potentially explaining why the northern blue-shifted clump emission is weaker than that for the southern red-shifted clump (Figs. 1a and b).

The velocity weighted (moment 1) map of the DCO⁺ emission (Fig. 1c) shows that its velocity gradient is similar to that of C¹⁸O, with the northern clump being blue-shifted and the southern clump being red-shifted, but with a smaller velocity range. The position-velocity (PV) diagrams of DCO⁺ and C¹⁸O emission (Fig. 3) are constructed and over-plotted with the best approximating thin disk models, Keplerian and Infall plus Keplerian out to 150 AU, obtained by Murillo et al. (2013). It appears as though both line emissions are well described by pure Keplerian rotation out to 300 AU. However the DCO⁺ emission is too weak, peaking at 7σ in the channel map, to carry out further kinematical analysis. The similar velocity signatures of the emission from both species may be due to the disk edge dragging along material from the envelope.

The above simple analysis of the observations indicates that the detected DCO⁺ emission may be tracing a ring which borders the C¹⁸O disk. In addition, DCO⁺’s velocity gradient and PV diagram suggest that the ring is undergoing Keplerian rotation, this is even more plausible given that the C¹⁸O (2−1) emission traces a rotationally supported disk with a Keplerian velocity profile (Murillo et al. 2013).

3.2. ALMA ACA cycle 2

The DCO⁺ (3−2) emission traced with the ACA confirms that the emission is concentrated around VLA1623 (Fig. 2), it does not peak on VLA1623A and instead encircles it, forming a shell-like structure. The ACA detects only weak emission at the position of VLA1623W.

The DCO⁺ emission mapped with the ACA shows the same velocity gradient (Fig. 2) as on the small scales, indicating that the kinematic structure is the same throughout.

4. Analysis

DCO⁺ is detected at the disk-envelope interface of VLA1623A. Given that DCO⁺ is optically thin and a good probe of temperature and CO freeze-out regions, we model the observed emission aiming to probe the physical and chemical structure of the disk-envelope interface. In this section we describe the model used and the results obtained from the model.

4.1. DCO⁺ chemical network and model

We model the observed DCO⁺ emission using a simple chemical network. Such a network, while it may not account for every possible reaction, is preferred as it gives insight into how the physical parameters of the temperature and density affect the observed emission. The network used is a steady-state analytic model that only takes the basic reactions that lead to DCO⁺ production and destruction into account (Table 1, Fig. 4). Since DCO⁺ is formed by the reaction of H₂D⁺ and CO, the rate-determining reactions in our network are ${\mathrm{H}}_{\mathrm{3}}^{\mathrm{+}}\mathrm{+}\mathrm{HD}\mathrm{-}\mathrm{\to }{\mathrm{H}}_{\mathrm{2}}{\mathrm{D}}^{\mathrm{+}}\mathrm{+}{\mathrm{H}}_{\mathrm{2}}$(1)with the back reaction ${\mathrm{H}}_{\mathrm{2}}^{\mathrm{D}}\mathrm{+}\mathrm{+}{\mathrm{H}}_{\mathrm{2}}\mathrm{+}\mathrm{\Delta }\mathit{E}\mathrm{-}\mathrm{\to }{\mathrm{H}}_{\mathrm{3}}^{\mathrm{+}}\mathrm{+}\mathrm{HD}$(2)where the activation energy due to the difference in zero point energy is ΔE ~ 220 K. For reaction 5 in our network (Table 1), for simplicity, we adopt the total rate coefficient summing over all three pathways of the reaction. The rate coefficients for a two-body reaction are given by $\mathit{k}\mathrm{=}\mathit{\alpha }{\left(\frac{\mathit{T}}{\mathrm{300}}\right)}^{\mathit{\beta }}\exp \left(\mathrm{-}\frac{\mathit{\gamma }}{\mathit{T}}\right)\mathrm{c}{\mathrm{m}}^{\mathrm{3}}{\mathrm{s}}^{-1}$(3)where T is the gas temperature, while $\mathit{k}\mathrm{=}\mathit{\zeta }{\mathrm{s}}^{-1}$(4)gives the rate coefficient for cosmic-ray ionization (reaction 1 from Table 1).

Fig. 4 Cartoon of the chemical network used to analyze the observed DCO+ emission. Red denotes the target molecule. Green the input molecular concentrations and orange the species involved in the bottleneck reaction.

Our model takes as input a source density and temperature profile as a function of radius and the parameters needed to calculate the rate coefficients. The CO evaporation temperature T_ev, desorption density n_de and CO abundance X_CO are free parameters. The CO abundance is assumed and not calculated. The CO evaporation temperature dictates when CO is in the gas phase (T>T_ev) or freezes onto the dust grains (T<T_ev). In a similar manner, the desorption density sets the boundary when the freeze-out timescales are too long (n < n_de) compared to the lifetime of the core (Jørgensen et al. 2005). We assume the density profile is equal to the H₂ density n_H2 and the abundance of HD X_HD = 10^-5 with respect to the total hydrogen nuclei density n_H = 2n(H₂). The model returns the calculated concentrations as a function of radius. The results can then be input into excitation radiative transfer programs such as RATRAN (Hogerheijde & van der Tak 2000) for further analysis.

Fig. 5 CO abundance profiles used in the model. The vertical dashed lines show the limits for the Drop abundance profile, evaporation temperature Tev and desorption density nde. X0 denotes the abundance in the constant profile. Xin, XD and X0 are the inner, drop and outer region abundances for the drop profile.

Given that CO is one of the parent molecules of DCO⁺, we study the effect of CO abundance through the use of different abundance profiles following the models detailed in Jørgensen et al. (2005) and Yıldız et al. (2010). The possible CO profiles are Constant and Drop abundance (Fig. 5). The former represents a fixed CO abundance throughout the core, whereas the latter, constrained by T_ev and n_de, is used to account for CO freeze-out in the chemical network. For the constant profile, the abundance is denoted by X₀. For the drop profile, the abundance in the inner, drop and outer regions are denoted by X_in, X_D and X₀, respectively. Previous studies using multiline single dish C¹⁸O have found that the abundance at X_in is lower than X₀ for a number of sources (Alonso-Albi et al. 2010; Yıldız et al. 2010, 2013). One explanation is that some fraction of the CO ice is transformed into more complex and less volatile carbonaceous species in the cold phase. We thus take this effect in our model into account. A lower abundance of CO, due to freeze-out, allows an increase in the abundance of H₂D⁺ (Mathews et al. 2013). However, since both molecules, CO and H₂D⁺, are parent molecules of DCO⁺, a balance must be reached before the effective formation of DCO⁺ takes place. This scenario is found to be common for the envelopes of early embedded protostars (Jørgensen et al. 2005) where the outer region shielded from the protostellar and interstellar radiation heating has low enough temperature that CO freezes out onto dust grains.

In molecular clouds, and consequently in protostellar cores, H₂ chemistry plays a major role, hence the ortho-to-para ratio of H₂ influences the chemical reactions, and has been found to be crucial to the deuterium chemistry (Flower et al. 2006; Pagani et al. 2009). The effect of ortho- and para-H₂ (o-H₂ and p-H₂) is studied in our model through the inclusion of the ortho-to-para ratio (o/p) and the distinction of o-H₂ and p-H₂ in the back reaction of the bottleneck (Eq. (2)) in the chemical network. o-H₂ and p-H₂ are only added in the back reaction since it is here where the distinction has a significant effect (Table 1). We set a lower limit on o/p of 10^-3 at low temperatures, as constrained by models and observations (Flower et al. 2006; Faure et al. 2013). The o-H₂ and p-H₂ reactions and their parameters for the rate coefficients are taken from Walmsley et al. (2004).

When o-H₂ and p-H₂ are included in the network, a thermal (LTE), upper- or lower-limit o/p ratio can be selected. In LTE, the ortho-to-para ratio is given by $\mathrm{o}\mathit{/}\mathrm{p}\mathrm{=}\mathrm{9}\exp \left(\frac{-170}{\mathit{T}}\right)$(5)where T is the gas temperature. Selecting the upper-limit ratio produces 3 times more o-H₂ than p-H₂. Since the back reaction with o-H₂ has a lower activation barrier γ than with p-H₂, the o/p upper-limit implies that H₂D⁺ is being destroyed faster than generated, leading to a decreased production of DCO⁺ (Fig. A.3), since H₂D⁺ is a parent molecule of DCO⁺ (Fig. 4). The lower-limit ratio, on the other hand, implies more p-H₂ which has a higher activation barrier γ for the back reaction, thus H₂D⁺ is generated faster than it is destroyed in turn increasing the DCO⁺ production (Fig. A.3).

As a starting point for our analysis of VLA1623A’s DCO⁺ emission, we use the density and temperature profile of VLA1623 obtained by Jørgensen et al. (2002) where 30 K is at ~1.5′′ assuming a distance of 120 pc. The profile was obtained by fitting single dish JCMT continuum images and the spectral energy distribution (SED) with continuum radiative transfer modeling using DUSTY resulting in a power law density profile of the form n ∝ r^-1.4. In the single dish continuum observations VLA1623A and B are unresolved and the density profile extends well beyond the ALMA field of view. It is expected, however, that VLA1623A dominates at 870 μm and that VLA1623B does not contribute much to the 450 μm continuum (Murillo & Lai 2013). Thus, the density and temperature profile obtained by Jørgensen et al. (2002) is representative of VLA1623A since VLA1623B is not significantly contributing to the continuum emission or the SED used to constrain the profiles.

The results of the analytic chemical network are run through the molecular excitation and line radiative transfer program RATRAN to generate line emission maps. Since the structure we are trying to reproduce is ring-like, we calculate the level populations with the 1D version of RATRAN, and then run the level populations with the 2D ray tracing to form the ring structure. The produced spectral image cubes are convolved with a Gaussian beam with the dimensions of the synthesized beam, continuum subtracted and then an intensity integrated map is generated. Radial profiles are extracted from the resulting images and compared with the observed profiles, which are integrated over the extent of the detected emission.

We find no significant difference in the peak position between the method used here and running the models through the ALMA simulator with the actual ALMA configuration of our Cycle 0 observations. However, the ALMA simulator does show that the emission ≥4″ is indeed filtered out by the Cycle 0 observations, producing somewhat narrower radial profiles. Molecular data for DCO⁺ and C¹⁸O are obtained from LAMDA (Schöier et al. 2005; C¹⁸O: Yang et al. 2010; DCO⁺ extrapolated from Flower 1999). As DCO⁺ is optically thin, the comparison mainly focuses on the position of the peak and the integrated intensity profile with respect to radius, the velocity profile has no effect on the resulting model integrated intensity. Thus, we assume a free-fall velocity profile with a central mass of 0.2 M_⊙ (Murillo et al. 2013) and T_dust = T_gas. The RATRAN output maps are convolved with the observed beam (0.85″ × 0.65″, PA 96.24°) and compared with the observed emission through radial cuts. In addition to comparing DCO⁺, we also compare C¹⁸O in order to further constrain the physical structure of the disk-envelope interface. However, for C¹⁸O we assume a pure Keplerian velocity profile, in agreement with the rotationally supported disk it traces. For both emission lines, we use an inclination of 55° (90° = face - on), in accordance with the results obtained by Murillo et al. (2013).

Fig. 6 Chemical network model results. Left and center columns: results for constant and drop CO abundance, respectively. Top, middle and bottom row: results of Test 1, 2 and 4, respectively. Light blue shaded region shows the observed peak position of DCO+. Right column: input source profile, where the original DUSTY source profile is shown with a dashed line so as to compare with the input for the specific test.

4.2. Modeling results

In this section the results of altering the CO abundance, density and temperature profiles in the chemical model are discussed. Figure 6 presents the models, with the light blue region showing the location of the observed DCO⁺ peak emission with respect to VLA1623A’s position. Model naming follows the scheme xYz where x is the test number from Table 3, Y is either C for constant or D for drop abundance, and z is the case from Table 2.

4.2.1. Chemical properties

To probe the chemical conditions, we establish one case for the constant CO abundance with X_CO = 10^-4 and six cases for the drop CO abundance with X_D ranging from 5 × 10^-6 to 10^-8 and varying T_ev and n_de. Parameter ranges were selected based on trends found in previous work (Jørgensen et al. 2005; Yıldız et al. 2010) and adapted to the current observations. Parameters for each case are listed in Table 2. As a zeroth order test, we compare the abundance profiles to the radial position of the observed DCO⁺ peak in the ALMA 12-m array data in order to investigate which chemical conditions best approximate our observations. This comparison neglects the fact that the peak emission radius also depends on the DCO⁺ excitation, which will be taken in Sect. 4.3 into account.

Figure 6 presents the model abundance profiles as functions of radius for different assumed CO abundances and physical structures. In general, the concentrations of CO and HCO⁺ drop with radius whereas those of H₂D⁺ and DCO⁺ increase with radius due to the lower temperature farther away from the source. Constant CO abundance (Fig. 6, left column) does not appear to produce a DCO⁺ peak within the expected region. Since the radial abundance profile of DCO⁺ is not altered by changing X_CO, we focus instead on the drop CO abundance profile.

While the drop CO abundance profile cannot alone alter the position of the peak, it produces several trends interesting to note. A kink in the abundance of DCO⁺ forms at T_ev for all the cases examined, though in most cases it is relatively small. Altering the abundance in the drop X_D changes the shape of the peak but does not significantly alter its position nor its abundance (Fig. A.2). Decreasing the X_D below 10^-7, however, causes the abundance of the peak to drop by several orders of magnitude and become similar in magnitude to the kink at T_ev, producing two peaks which are not observed (Fig, A.2, case f). Varying T_ev and n_de, in order to constrain the width of the drop, generates a very narrow peak which increases and drops quickly.

To test the effect of the o/p ratio on our chemical network, we vary the ratio from the thermalized value to the upper- and lower-limit for the drop CO abundance profile. In all cases, we find that setting o/p to the lower-limit does not change the position of the DCO⁺ emission since the peak and bulk of the DCO⁺ concentration for the thermalized o/p is already located in the o/p = 10^-3 range (Fig. A.3). In fact, the peak of the modeled DCO⁺ emission starts to decrease as the o/p ratio increases. The lower-limit o/p only alters the inner regions (<100 AU), which the data do not constrain. The upper-limit reduces the overall production of DCO⁺ and pushes the peak outward to larger radii. In conclusion, altering the CO abundance and o/p ratio does not produce the observed results.

4.2.2. Physical properties

Altering the chemistry of the model does not reproduce the observed DCO⁺ peak position, regardless of the case used. Thus we alter the density and temperature profile in order to find the conditions necessary to reproduce the observed emission. We set up 4 tests, including the original profile from Jørgensen et al. (2002) referred to as Test 1, by increasing or decreasing by a constant factor either the density or temperature profile. Parameters for each test are listed in Table 3.

Increasing or decreasing the density profile by one order of magnitude, Tests 2 and 3, does not generate a significant radial shift of the DCO⁺ peak compared with the original profile (Fig. 6, middle row). Altering the density, however, affects the concentration of the peak upwards or downwards with an increase or decrease of the density, respectively.

The factor of 1.5 decrease in the temperature profile with cases b and d with X_D ≈ 10^-6−10^-7 shifts the modeled DCO⁺ peak effortlessly to the observed position (Fig. 6). Decreasing the temperature profile between a factor of 1 to 2 moves the modeled DCO⁺ peak inward (Fig. A.4). For the constant CO abundance, any alteration of the temperature either over-estimates the modeled DCO⁺ abundance or produces a peak too far inward. For the drop CO abundance, a temperature drop less than a factor of 1.5 does not move the peak inward enough, whereas a larger factor moves it too far inward. Cases b and d with X_D ≈ 10^-6−10^-7 shifts the modeled DCO⁺ peak effortlessly to the observed position (Fig. 6 bottom row). Changing the CO abundance in the drop (Table 2) does not alter the location of the peak, as expected from the results presented in Sect. 4.2.1. Examining the temperature profile shows that the DCO⁺ emission peaks at a range of 11−16 K. This is lower than the expected 20 K at ~3″ inferred from radiative transfer modeling of the observed continuum data (Jørgensen et al. 2002). We limit the decrease in the temperature profile to not fall below 8 K (Zucconi et al. 2001), however this limit generates no significant change in the outcome of our model since the limit falls near and beyond the edge of the FOV of our observations (Fig. 6, bottom right panel).

From the results obtained from altering the chemical and physical conditions in our model, we can deduce that the DCO⁺ emission is located at its observed position due to physical conditions, namely a lowered temperature, and not to special chemical conditions in the disk-envelope interface.

4.2.3. Comparison with full chemical network

To explore the limitations of the simple network, we compare the results to a full time-dependent deuterated chemical network, based on the RATE06 version of the UMIST database for Astrochemistry (Woodall et al. 2007) extended with deuterium fractionation reactions (McElroy et al. 2013). Models with only gas phase chemistry and including gas-grain balance (freeze-out, thermal desorption, and cosmic-ray-induced photodesorption) are run. The models start with the same initial abundandances (H₂:HD:CO = 1:3 × 10^-5:1 × 10^-4) and the same physical structure as the simple network. We find that the large network confirms both the general trend found in the simple network and also reproduces the peak position and abundance of DCO⁺ at the observed position for abundances extracted at early times (~10⁵ yr).

Inspection of the full network shows that in addition to the reactions listed in Table 1, the $\mathrm{HC}{\mathrm{O}}^{\mathrm{+}}\mathrm{+}\mathrm{D}\mathrm{-}\mathrm{\to }{\mathrm{H}}_{}^{}\mathrm{+}\mathrm{DC}{\mathrm{O}}^{\mathrm{+}}$ reaction can become important if D/HD is large. The back reaction has a reaction barrier of ~800 K, thus leading to fractionation of DCO⁺ (Adams & Smith 1985). In practice, this reaction is only relevant for much lower densities than encountered in the VLA1623A envelope model. Similarly, dissociative recombination of H$\begin{array}{c}\mathrm{+}\\ \mathrm{3}\end{array}$ with electrons cannot be neglected at low densities. At very early times (<10³ yr, depending on density), reactions with H⁺ and D⁺ become more significant than those with H$\begin{array}{c}\mathrm{+}\\ \mathrm{3}\end{array}$ and H₂D⁺. Finally, in full gas-grain models the bulk of the CO is frozen out in the cold outer envelope resulting in very low DCO⁺ concentrations, unless an efficient non-thermal CO desorption process is included.

4.3. Emission profiles

Aiming to better constrain the physical conditions of the disk-envelope interface we compare the models and observations of both DCO⁺ and C¹⁸O emission computed with RATRAN. For the comparison, the observed southern red-shifted clump in the ALMA 12-m array data for both lines is selected. The reason is two-fold: i) the red-shifted emission of the C¹⁸O disk suffers less absorption from the outer envelope than the blue-shifted lobe as noted above (Murillo et al. 2013); ii) the southern red-shifted clump traced by DCO⁺ is the strongest and most prominent. Similarly, for the model of C¹⁸O only one lobe of the disk is selected for comparison.

Consistent with the abundance plots, the constant CO abundance profiles do not fit the observed DCO⁺ emission in any of the four tests, producing very broad peaks between 3″ and 5″ away from the source for tests 1, 2 and 3; and a peak at 1″ for test 4. In a similar manner, for the drop CO abundance scenario, tests 1, 2 and 3 all produce peaks beyond 3″ with varying broadness, and thus do not approximate the observed DCO⁺ emission. This occurs for all examined cases (Fig. A.5). Models 4Db and 4Dd approximate the observed DCO⁺ emission well (Fig. 7, orange and red lines, respectively). These results are in agreement with the conclusions drawn from the concentrations in the analytics chemical network model (see Sect. 4.2.2).

Comparing the results of the model with the observed C¹⁸O emission (Fig. 8), we find that in all four tests a constant CO abundance over-predicts the extent of C¹⁸O, with the emission well above our 3σ level extending out to 3″ or further, whereas we observe it only out to less than 2″ from the central source. When the drop CO abundance is introduced, the source profile of VLA1623 from Jørgensen et al. (2002) overestimates the C¹⁸O emission in all cases, and generates a second peak at about 3″ at half the intensity of the central peak in almost all cases. The width of the drop does not have any significant effect on the modeled emission. The decreased density and decreased temperature profile tests do not fare better, as they again largely overestimate the C¹⁸O beyond 2″ and even produce secondary peaks. Case a (X_D = 10^-6) for the original profile and decreased density seem to show some promise, however the modeled emission is around 3σ at 1.5″ whereas observations at that radius are closer to 1σ. Surprisingly, we find that the increased density profile produces the best results with case d (X_D = 10^-7), the same as that for DCO⁺ (Fig. 8, red line).

Fig. 7 Comparison between observations and radiative transfer modeling of the chemical network results. Here the best approximations are highlighted, with model 4Dd being the preferred model. The models are referred to by test, abundance profile (Constant “C” or Drop “D”) and case. Thus, 4Dd is the model for test 4 with Drop abundance and case d. Black line shows the DCO+ observations integrated over the southern clump. Gray solid and dashed lines show the 1σ and 4σ levels, respectively.

Fig. 8 Same as in Fig. 7 but for C18O, with 2Dd being the preferred model. Gray solid and dashed lines show the 1σ and 3σ levels, respectively. 2Dd, increased density profile, provides the best approximation to the observed C18O.

Hence we find that the DCO⁺ and C¹⁸O observed emission are not reproduced by the same physical structure (Figs. 7 and 8). DCO⁺ is well modeled by a temperature profile a factor of 1.5 lower than that needed to model C¹⁸O, whereas C¹⁸O is reproduced by a density profile an order of magnitude higher than required for DCO⁺. The abundance in the drop X_D necessary to reproduce the emission is the same for both lines, X_D = 10^-7.

Fig. 9 Radial profiles of the DCO+ ACA observations along the disk plane (top) and along the outflow axis perpendicular to the disk (bottom) as shown in Fig. 2. Radial profile of the DCO+ Cycle 0 observations is shown in light gray. Observed radial profiles are overlaid with models convolved to the resolution of the ACA observations. Gray solid and dashed lines show the 1σ and 3σ levels, respectively, of the ACA emission. These results evidence that the DCO+ emission along the disk plane is best approximated by model 4Dd, whereas along the outflow is best described by model 1C.

4.4. Large scale emission vs. model

We also compare the radial integrated intensity profile of the DCO⁺ ACA observations with the model, both along the plane of the disk and perpendicular to the disk in the outflow direction, since these data probe larger spatial scales than the 12-m array data. For this purpose, the models were run through RATRAN’s ray tracer sky again but instead using the cell size of the ACA observations (0.84″) and then were convolved to match the resolution of the ACA observations. The observed radial profiles are extracted from the red-shifted emission.

The profile along the plane of the disk (peak at ~4″) is well described by the same model as in the small scale, the drop CO abundance profile for case d, with the temperature decreased by a factor of 1.5 (4Dd, Fig. 9). The temperature decreased by a factor <1.5 does not provide a better estimate. The constant CO abundance again overestimates the amount of DCO⁺ produced in the outer regions.

On the other hand, the integrated intensity profile along the outflow direction (peak at distances ≈5″) is very well described by the model 1C (Fig. 9). This indicates that in regions that are not shadowed by the disk, the production of DCO⁺ is as expected from dust continuum models.

5. Discussion

The results of our modeling show that the position of the DCO⁺ emission (~2.5″=~300 AU) is closer to VLA1623A than expected based on the spherically symmetric dust continuum radiative transfer modeling from DUSTY (~4″ = ~480 AU). The results further show that the observed emission peaks at a dust temperature range of 11−16 K. Comparison of the radiative transfer results from RATRAN to the C¹⁸O emission from the disk suggest that the disk has a density higher by one order of magnitude than the emitting structure traced by DCO⁺.

A possible explanation would be that VLA1623A is a very low luminosity object (VeLLO, Young et al. 2004; Dunham et al. 2008) undergoing episodic accretion and just coming out of the quiescent phase. This is highly unlikely, however, since Johnstone et al. (2013) find that the timescale for dust and gas to heat up after an accretion burst is short, on the order of hours to weeks. VLA1623A has been reported of having a bolometric luminosity between 0.4−2 L_⊙ from early observations to more recent work (André et al. 1993; Froebrich 2005; Murillo & Lai 2013; Chen et al. 2013), in contrast to the expected luminosity of 10^-2 L_⊙ for VeLLOs. Hence VLA1623A has not recently come out of the quiescent accretion phase and the location of the observed DCO⁺ emission can not be attributed to the relic of a previous phase of decreased accretion.

Fig. 10 Radiative 2D models based on Harsono et al. (2013). For all models an envelope with Menv = 1 M⊙ irradiated by a 1 L⊙ protostar is assumed, having an outflow cavity with a full opening angle of 30°. The dashed lines indicate the temperature contours for 100, 50 and 20 K. The top and middle panels include a 180 AU and 0.02 M⊙ disk assuming a 0.1 AU scale height, with centrifugal radius Rc = 200 AU (top) and 50 AU (center) orientated in the horizontal direction. The outflow cavity is oriented in the vertical direction. The bottom panel shows the envelope without a disk.

A more plausible explanation to the position and temperature of the region containing the DCO⁺ emission is disk shielding. DCO⁺ is observed to border VLA1623A’s C¹⁸O disk in our present data. This could shield the outer parts of the disk from heating by the central protostar, thus moving in the freeze-out zone of CO, enhancing the production of DCO⁺ closer to the protostar. This scenario is further supported by the result of our simple chemical model that the C¹⁸O disk is more dense than the region of the DCO⁺ emission and than expected from Jørgensen et al. (2002)’s envelope density profile of VLA1623.

We test the possibility of the disk-shielding scenario with radiative 2D disk plus envelope models using radiative transfer methods as in Harsono et al. (2013; Fig. 10). For the models, we assume a central protostellar source of 1 L_⊙ (Murillo & Lai 2013; Chen et al. 2013), an envelope with a mass of 1 M_⊙ (André et al. 1993; Froebrich 2005), with the addition of a disk mass of 0.02 M_⊙ and radius of 180 AU (Murillo et al. 2013). The outflow cavity is assumed to have an opening angle of 30°. A thin disk with a scale height of 0.1 AU is adopted. Disk flaring is not included since we have no information on the flaring of VLA1623A’s disk and the thin disk model approximates the C¹⁸O kinematics well. Two values for the centrifugal radius R_c, 200 AU and 50 AU, are chosen. Within the centrifugal radius the velocity structure of the disk is Keplerian. We find that even for such a thin disk the temperature along the plane of the disk is lower than for the envelope at the same radius (Fig. 10, top and middle panels), thus moving the CO freeze-out zone closer to the protostar along the edge of the disk than in other regions of the core. For either centrifugal radius the temperature beyond 200 AU drops well below 20 K. Finally, we test whether the presence of a disk makes a difference in the location of the CO freeze out region. Figure 10 bottom panel, shows the model with the same conditions as in Fig. 10 middle panel, but without the disk. This shows that omitting the disk causes the CO freeze out region to move outward to ~400 AU, 150 AU further out than the models including a disk.

The temperature conditions required to reproduce the observed DCO⁺ emission with our simple chemical model are therefore in agreement with the results obtained from the 2D radiative transfer disk plus envelope model. These results strongly support the scenario where a disk can shield the regions at its edge from heating by the protostar. This shielding causes the CO freeze-out region to move inward toward the edge of the disk, bringing low-temperature enhanced molecules, such as DCO⁺, closer to the central protostar along the plane of the disk. However, the rest of the envelope is largely unaffected by the disk and a shell of molecules such as DCO⁺ forms at a radius predicted by spherically symmetric models. The ACA observations of DCO⁺ toward VLA1623A provide further evidence for this effect. These results show that the DCO⁺ emission peaks closer to the source along the disk plane than along the outflow direction.

Fig. 11 Cartoon of the disk and envelope of VLA1623A. The ring and disk, traced in DCO+ and C18O respectively, are observed with ALMA 12-m array, while the outer spherically symmetric shell, traced in DCO+, is observed with ALMA ACA.

6. Conclusions

This work presents the results and analysis of ALMA Cycle 0 Early Science Band 6 observations of DCO⁺ (3−2) in the extended configuration toward VLA1623A probing subarcsec scales, as well as Cycle 2 ACA data probing the larger scales. A simple chemical network was setup taking ortho- and para-H₂ into account in the rate-determining reactions. The density and temperature profile of VLA1623A was obtained from fitting the SED and dust continuum data with radiative transfer modeling using DUSTY (Jørgensen et al. 2002). Our simple chemical model coupled with VLA1623A’s physical structure served as input for line radiative transfer calculations with RATRAN with 2D ray-tracing. The CO abundance, density and temperature profiles were altered to study the effect of each parameter on the location of the observed DCO⁺ peak. The results of our observations and analysis can be thus summarized:

• 1.

  DCO⁺ is observed to border the C¹⁸O disk around VLA1623A.Both emission lines show similar velocity gradients(blue-shifted to the north and red-shifted to the south), with DCO⁺emission between 2.8 to5.2 km s^-1. The PV diagrams of C¹⁸O and DCO⁺ suggest that both emission lines are well described by Keplerian rotation. However, the DCO⁺ emission is weak, thus no further kinematical analysis was carried out, and whether the disk extends out to 300 AU cannot be confirmed.

• 2.

  Using a simple chemical network with the inclusion of ortho- and para-H₂ as well as non-LTE line radiative transfer, we model the observed DCO⁺ emission. We find that using a constant CO abundance predicts a DCO⁺ peak at around 4″, twice further out than observed, irrespective of the adopted o/p ratio and density profile. Our model results show that a drop CO abundance with a decreased temperature profile by a factor of 1.5 generates a peak at the same position as the observed emission. Thus, the observed DCO⁺ peak is closer to VLA1623A and with a lower temperature (11−16 K) than that expected from a spherically symmetric physical structure constrained by continuum data and source SED.

• 3.

  The observed DCO⁺ and C¹⁸O emission are not described by the same physical structure. In our model, the C¹⁸O emission is well reproduced by a drop CO abundance with the density profile increased by one order of magnitude at 1″ (≤120 AU) radii, in comparison with the density profile needed to reproduce the observed DCO⁺ emission. A constant CO abundance and a decreased temperature profile over-predict the extent of the C¹⁸O emission.

• 4.

  Disk-shielding is the best possible explanation for the observed DCO⁺ emission toward VLA1623A. Disk-shielding causes the inward shift of the CO freeze-out region along the plane of the disk, lowering the dust temperature to <20 K, generating a ring of molecules whose abundance is enhanced by low temperatures such as DCO⁺. The rest of the envelope is largely unaffected by the disk, thus the CO freeze-out shell predicted by spherically symmetric radiative models is expected to be located further out. This prediction is confirmed by our recent ALMA Cycle 2 ACA observations, which show that the DCO⁺ emission along the outflow axis lies at larger radii, ~5″, consistent with constant CO abundance models without any alteration to Jørgensen et al. (2002)’s temperature and density profile of VLA1623.

• 5.

  The disk-envelope interface in VLA1623A is shown to have a broken transition in density and temperature, with the impact generated by the presence of a disk being observable from small to large scales.

Our observations and modeling results for VLA1623A show the disk-envelope interface to have different physical conditions than other parts of the envelope. Our results also highlight the drastic impact that the disk has on the temperature structure at ~100 AU along the plane of the disk, with the effect being observable even at large scales. We suggest performing further observations to determine whether the unequal physical and chemical conditions observed in the disk-envelope interface of VLA1623A is a common phenomenon in protostellar systems with rotationally supported disks or a special condition of the present source.

Online material

Appendix A: Modeling results: additional figures

Here we present additional figures illustrating details discussed in Sects. 4.2 and 4.3. Figures A.1 and A.2 show the effect of CO abundance X_CO in the drop on the DCO⁺ peak for Test 1 and 4, respectively. Figure A.3 highlights the effect of the upper and lower o/p limits on the DCO⁺ abundance for

the best approximation models to the DCO⁺ observations. The effect of the factor used to decrease the temperature profile on the peak position of the modeled DCO⁺ is explored in Fig. A.4. Figures A.5 and A.6 give a sampling of the observations compared with the radiative transfer modeling of the chemical network results for the different tests, CO abundance profiles and cases.

Fig. A.1 Chemical network results for Drop CO abundance with the Test 1 source profile from Jørgensen et al. (2002) for four cases (Table 2). Case d is the same as in the top row center panel of Fig. 6.

Fig. A.2 Chemical network results for Drop CO abundance with the decreased-temperature (Test 4) source profile for four cases (Table 2). Case d is the same as in the bottom row center panel of Fig. 6.

Fig. A.3 Chemical network results for Drop CO abundance with the decreased-temperature source profile (Test 4) for case d with upper-limit (left), thermal (center) and lower-limit (right) o/p ratio.

Fig. A.4 DCO+ Cycle 0 observations compared with the radial profile of the constant and drop CO abundance models, showing the effect of the factor used to decrease the temperature, for case d. A factor of 1.5 with drop CO abundance provides the best approximation to the observed radial profile.

Fig. A.5 Comparison of the observations and radiative transfer modeling of the chemical network results. The black line shows the DCO+ observations integrated over the southern clump, with the gray solid and dashed lines showing the 1σ and 4σ levels, respectively. Models are referred to by test, abundance profile (Constant “C” or Drop “D”) and case. Color lines show tests 1 (blue), 2 (green), 3 (orange) and 4 (purple) for the top left panel, and cases a (blue), b (green), d (orange) and f (purple) for the other four panels.

Fig. A.6 As in Fig. A.5 but for C18O. Gray solid and dashed lines show the 1σ and 3σ levels respectively. Note that constant CO abundance over-predicts the amount of C18O emission.

Acknowledgments

This paper made use of the following ALMA data: ADS/JAO.ALMA 2011.0.00902.S and 2013.1.01004.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO, and NAOJ. The 2011.0.00902.S data was obtained by N.M.M. while she was a Master student at National Tsing Hua University, Taiwan, under the supervision of S.P.L. The authors thank J. C. Mottram for his help in using RATRAN; and M. Hogerheijde for permission to use the 2D version of RATRAN. Astrochemistry in Leiden is supported by the Netherlands Research School for Astronomy (NOVA), by a Royal Netherlands Academy of Arts and Sciences (KNAW) professor prize, and by the European Union A-ERC grant 291141 CHEMPLAN. This work is supported by grant 639.041.335 from the Netherlands Organisation for Scientific Research (NWO) and by the Netherlands Research School for Astronomy (NOVA) and by the Space Research Organization Netherlands (SRON). S.P.L. acknowledges support from the Ministry of Science and Technology of Taiwan with Grants MOST 102-2119-M-007-004-MY3.

References

All Tables

All Figures

Fig. 1 DCO+ compared with continuum and C18O toward VLA1623A and B, marked with a square and triangle, respectively. The black ellipse shows the synthesized beam. Intensity integrated DCO+ (color-scale) with a) 1.3 mm continuum and b) C18O. DCO+ peaks at 0.0433 Jy beam-1 km s-1. Contours are in steps of 3, 5, 10, 15, 20, 40, 60 and 78σ with σ = 1 mJy beam-1 for 1.3 mm continuum, and −5, −3, 3, 5, 10, 15, 20 and 25σ with σ = 13 mJy beam-1 km s-1 for C18O. Velocity map (moment 1, color-scale) for c) DCO+ and d) C18O overlaid with the corresponding intensity integrated map (contours). Contours are in steps of −10, −7, −4, 4, 7, 10 and 11σ with σ = 3 mJy beam-1 km s-1 for DCO+, and the same as in b) for C18O.

In the text

Fig. 2 Intensity (contours) and velocity (color-scale) integrated maps of DCO+ ACA observations. Note the 10 times larger scale of this figure compared with Fig. 1. The positions of VLA1623A, B and W are marked with a square, triangle and circle, respectively. The size of the filled square shows the size of the C18O disk. The red and blue lines indicate the outflow direction while the gray line shows the mayor axis of the disk plane. The black ellipse shows the synthesized beam. Contours are in steps of −3, 3, 5, 7, 10, 15 and 17σ with σ = 220 mJy beam-1 km s-1.

In the text

	Fig. 3 Cycle 0 DCO+ (black) and C18O (gray) PV diagram. Cuts are made at PA = 26° for DCO+ and PA = 35° for C18O. Keplerian (solid line) and Infall plus Keplerian (dashed line) rotation PV models with M∗ = 0.2 M⊙ are overplotted (Murillo et al. 2013). Contours are in steps of −3, 3, 5, 10 and 15σ with σ = 19 mJy beam-1 for C18O, −7, −5, −3, −2, 2, 3, 4, 5σ with σ = 12 mJy beam-1 for Cycle 0 DCO+.
In the text

	Fig. 4 Cartoon of the chemical network used to analyze the observed DCO+ emission. Red denotes the target molecule. Green the input molecular concentrations and orange the species involved in the bottleneck reaction.
In the text

	Fig. 5 CO abundance profiles used in the model. The vertical dashed lines show the limits for the Drop abundance profile, evaporation temperature Tev and desorption density nde. X0 denotes the abundance in the constant profile. Xin, XD and X0 are the inner, drop and outer region abundances for the drop profile.
In the text

Fig. 6 Chemical network model results. Left and center columns: results for constant and drop CO abundance, respectively. Top, middle and bottom row: results of Test 1, 2 and 4, respectively. Light blue shaded region shows the observed peak position of DCO+. Right column: input source profile, where the original DUSTY source profile is shown with a dashed line so as to compare with the input for the specific test.

In the text

Fig. 7 Comparison between observations and radiative transfer modeling of the chemical network results. Here the best approximations are highlighted, with model 4Dd being the preferred model. The models are referred to by test, abundance profile (Constant “C” or Drop “D”) and case. Thus, 4Dd is the model for test 4 with Drop abundance and case d. Black line shows the DCO+ observations integrated over the southern clump. Gray solid and dashed lines show the 1σ and 4σ levels, respectively.

In the text

	Fig. 8 Same as in Fig. 7 but for C18O, with 2Dd being the preferred model. Gray solid and dashed lines show the 1σ and 3σ levels, respectively. 2Dd, increased density profile, provides the best approximation to the observed C18O.
In the text

Fig. 9 Radial profiles of the DCO+ ACA observations along the disk plane (top) and along the outflow axis perpendicular to the disk (bottom) as shown in Fig. 2. Radial profile of the DCO+ Cycle 0 observations is shown in light gray. Observed radial profiles are overlaid with models convolved to the resolution of the ACA observations. Gray solid and dashed lines show the 1σ and 3σ levels, respectively, of the ACA emission. These results evidence that the DCO+ emission along the disk plane is best approximated by model 4Dd, whereas along the outflow is best described by model 1C.

In the text

Fig. 10 Radiative 2D models based on Harsono et al. (2013). For all models an envelope with Menv = 1 M⊙ irradiated by a 1 L⊙ protostar is assumed, having an outflow cavity with a full opening angle of 30°. The dashed lines indicate the temperature contours for 100, 50 and 20 K. The top and middle panels include a 180 AU and 0.02 M⊙ disk assuming a 0.1 AU scale height, with centrifugal radius Rc = 200 AU (top) and 50 AU (center) orientated in the horizontal direction. The outflow cavity is oriented in the vertical direction. The bottom panel shows the envelope without a disk.

In the text

	Fig. 11 Cartoon of the disk and envelope of VLA1623A. The ring and disk, traced in DCO+ and C18O respectively, are observed with ALMA 12-m array, while the outer spherically symmetric shell, traced in DCO+, is observed with ALMA ACA.
In the text

	Fig. A.1 Chemical network results for Drop CO abundance with the Test 1 source profile from Jørgensen et al. (2002) for four cases (Table 2). Case d is the same as in the top row center panel of Fig. 6.
In the text

	Fig. A.2 Chemical network results for Drop CO abundance with the decreased-temperature (Test 4) source profile for four cases (Table 2). Case d is the same as in the bottom row center panel of Fig. 6.
In the text

	Fig. A.3 Chemical network results for Drop CO abundance with the decreased-temperature source profile (Test 4) for case d with upper-limit (left), thermal (center) and lower-limit (right) o/p ratio.
In the text

	Fig. A.4 DCO+ Cycle 0 observations compared with the radial profile of the constant and drop CO abundance models, showing the effect of the factor used to decrease the temperature, for case d. A factor of 1.5 with drop CO abundance provides the best approximation to the observed radial profile.
In the text

Fig. A.5 Comparison of the observations and radiative transfer modeling of the chemical network results. The black line shows the DCO+ observations integrated over the southern clump, with the gray solid and dashed lines showing the 1σ and 4σ levels, respectively. Models are referred to by test, abundance profile (Constant “C” or Drop “D”) and case. Color lines show tests 1 (blue), 2 (green), 3 (orange) and 4 (purple) for the top left panel, and cases a (blue), b (green), d (orange) and f (purple) for the other four panels.

In the text

	Fig. A.6 As in Fig. A.5 but for C18O. Gray solid and dashed lines show the 1σ and 3σ levels respectively. Note that constant CO abundance over-predicts the amount of C18O emission.
In the text