THE INITIAL CONDITIONS OF CLUSTERED STAR FORMATION. III. THE DEUTERIUM FRACTIONATION OF THE OPHIUCHUS B2 CORE

R. K. Friesen; J. Di Francesco; P. C. Myers; A. Belloche; Y. L. Shirley; T. L. Bourke; P. André

doi:10.1088/0004-637X/718/2/666

1. INTRODUCTION

The Ophiuchus molecular cloud is a source of ongoing low-mass clustered star formation relatively close to the Sun (d ∼ 120 pc; Lombardi et al. 2008). Through low-resolution DCO⁺ observations, Loren et al. (1990) discovered a number of dense cores (labeled A through F) in the central Ophiuchus region with masses of a few to several tens of solar masses, which represented a range in evolutionary status given the presence or absence of embedded protostars. Large-scale submillimeter and millimeter observations of Ophiuchus identified several new cores and revealed that the cores were further fragmented into multiple dense clumps, each with masses M ≲ 1 M_☉ (Motte et al. 1998; Johnstone et al. 2000, 2004; Young et al. 2006; Sadavoy et al. 2010). Analysis of the clump kinematics suggests that most of them are in virial equilibrium and, hence, not transient objects (André et al. 2007). Recent Spitzer Space Telescope studies of infrared emission toward the cores have characterized in detail the starless and protostellar clump population (Jørgensen et al. 2008; Enoch et al. 2009). These objects are thus ideal locations to study the physical characteristics of dense cores in the process of forming small stellar clusters.

To characterize the physical properties of the dense, cluster-forming Ophiuchus cores, we have previously studied Oph B, C, and F through high-resolution NH₃ (1,1) and (2,2) observations (Friesen et al. 2009, hereafter Paper I) and Oph B in both single-dish and interferometer N₂H⁺ 1–0 emission (Friesen et al. 2010, hereafter Paper II). The NH₃ and N₂H⁺ observations probe gas at densities n ∼ 10³–10⁴ cm⁻³ and n ∼ 10⁵ cm⁻³, respectively.

In Paper I, we identified (through the three-dimensional clumpfind structure-finding algorithm; Williams et al. 1994) a number of individual clumps in NH₃ emission in the cores, which did not correlate well with those identified in continuum emission. We determined that Oph B and F were warmer on average than typically found for isolated cores (mean T_K ∼ 15 K), and we found no trend with temperatures and H₂ column density calculated from submillimeter continuum emission. In contrast, Oph C is colder, with temperatures decreasing to a minimum of 9 K at the dust continuum emission peak. Line widths in both the Oph B1 and B2 sub-cores were found to be wide and transonic, despite the presence of three embedded protostars in B2 and none in B1. Oph F, with four embedded protostars, has regions of small and large line widths and complicated line structures. Again similar to isolated objects, NH₃ line widths in Oph C are extremely narrow and consistent with being due to purely thermal motions. In an interesting result, the fractional abundance of NH₃ with respect to H₂, X(NH₃), appeared to decrease with increasing H₂ column density, N(H₂), suggesting that NH₃ may be depleted toward locations of high column density and thus may not be tracing the densest core gas. Finally, Paper I revealed that cores with characteristics of both isolated and clustered environments can coexist in a clustered star-forming environment.

In Paper II, we again used clumpfind to identify small-scale structure in N₂H⁺ emission from Oph B and found that the N₂H⁺ clumps matched well NH₃ clumps from Paper Iand were consequently again offset from continuum objects. In both studies, little difference was found in the gas properties toward embedded protostars compared with the general core gas. Line widths in Oph B2 remain large at the densities traced by N₂H⁺. In B1, however, N₂H⁺ line widths are significantly narrower than found in NH₃ emission. Applying the gas temperatures determined in Paper I, we found that motions in Oph B1 were subsonic with a mean ratio of non-thermal line width to sound speed σ_NT/c_s = 0.7 across the core, while in B2 σ_NT/c_s = 1.3. A clear trend of decreasing X(N₂H⁺) with increasing N(H₂) was found toward Oph B2, with emission suggestive of a N₂H⁺ hole toward the continuum emission peak found in high resolution N₂H⁺ 1–0 emission observed with the Australia Telescope Compact Array.

The results of NH₃ and N₂H⁺ observations toward Oph B2 in Paper I and Paper II indicate that the highest density gas may best be probed by molecular emission lines excited at higher critical densities or by deuterated molecules, which are expected to be abundant in cold cores where significant depletion has occurred. In cold molecular cloud cores, the deuteration process is initiated by the reaction

$\begin{equation} \mbox{H$_3^+$} + {\rm HD} \rightleftharpoons \mbox{H$_2$D$^+$} + \mbox{H$_2$} + \Delta E \end{equation} \tag{ 1 }$

which is exothermic in the forward direction with ΔE = 232 K (Millar et al. 1989). At very low temperatures (≲20 K), the forward reaction dominates, making H₂D⁺ a key molecular ion in the enhancement of deuterated species in a molecular core. CO, which can destroy both H₂D⁺ and H⁺₃ via proton-transfer reactions, is additionally expected to deplete from the gas phase through freeze-out onto dust grains. The deuterium fractionation can thus be elevated orders of magnitude above the cosmic [D]/[H] ratio of 10⁻⁵ for H⁺₃, as well as for other species which react with its deuterated forms.

The ortho-H₂D⁺ 1₁₁–1₁₀ transition at 372 GHz was first detected toward the low-mass young stellar object (YSO) NGC 1333 IRAS 4A by Stark et al. (1999). Multiple detections of H₂D⁺ and doubly deuterated H⁺₃, D₂H⁺, have followed in both starless and protostellar cores (Caselli et al. 2003; Stark et al. 2004; Vastel et al. 2004), including a recent survey of H₂D⁺ emission toward 16 objects by Caselli et al. (2008). Walmsley et al. (2004) show that the H₂D⁺ abundance determined for the dark core L1544 by Caselli et al. (2003) is consistent with a "complete depletion" model, where heavy elements such as C, N, and O (and, consequently, molecules containing these elements) are missing from the gas phase at the highest core densities. In this scenario where no heavy elements are left in the gas phase, H₂D⁺ and multiply-deuterated forms of H⁺₃ are the only gas tracers of the physical conditions of the dense core.

Based on our previous work in Paper I and Paper II, Oph B2 is a good target to investigate the distribution of deuterated species in a more complex star-forming core. To this end, we present in this paper N₂D⁺ 3–2, N₂H⁺ 4–3, and H₂D⁺ 1₁₁–1₁₀ emission maps of Oph B2. We describe the observations in Section 2 and discuss the observed distributions of N₂D⁺ 3–2, H₂D⁺ 1₁₁–1₁₀, and N₂H⁺ 4–3 in Section 3. We fit the observed spectra with a single Gaussian (for H₂D⁺) or a multiple-component Gaussian (for N₂D⁺ and N₂H⁺ due to their hyperfine line structure) and discuss the line velocity centroids and widths in Section 4, and also present column density and fractional abundance calculations. In Section 5, we look at the non-thermal line widths in Oph B2 as a function of critical density of the observed emission lines and discuss trends in the deuterium fractionation. Finally, we present a lower limit to the electron abundance in B2 and summarize our results in Section 6.

2. OBSERVATIONS

In the following section, we discuss the observations presented in this paper. We list the species, transitions, and rest frequencies observed in Table 1.

Table 1. Observed Species, Transitions, and Rest Frequencies

Species	Transition	F'₁	F'	F₁	F	Frequency
						(GHz)
N₂D⁺	3–2	4	5	3	4	231.321906^a
N₂H⁺	4–3	5	5	4	4	372.672494^b
H₂D⁺	1₁₁–1₁₀					372.421385^c

Notes. ^aGerin et al. (2001). ^bPagani et al. (2009a). ^cAmano & Hirao (2005).

Download table as: ASCII Typeset image

2.1. N₂D⁺ at IRAM

An On-The-Fly (OTF) map of the N₂D⁺ 3–2 line at 231.321 GHz in Oph B2 was made at the Institut de Radio Astronomie Millimétrique (IRAM) 30 m Telescope during the 2007 and 2008 winter semesters using the HEterodyne Receiver Array (HERA; Schuster et al. 2004). HERA is a multi-beam spectral line receiver with nine dual polarization pixels arranged in a 3 × 3 array, with pixels separated by 24'' (approximately twice the 231 GHz 11'' beam full width half-power). The maps were performed with the array tracking the sky rotation and rotated by 9 fdg 5 with respect to the equatorial system, such that the spacing between OTF lines was 4''. Observations were frequency-switched in-band with a throw of 2 × 6.9 MHz. The Versatile Spectrometer Assembly autocorrelator was used as the backend, with a 20 MHz bandwidth and 40 kHz channel spacing, corresponding to 0.05 km s⁻¹ per channel at 231 GHz.

The data were calibrated to T*_A units at the telescope using the Multichannel Imaging and Calibration Software for Receiver Arrays. Further calibration, including baselining and folding of the frequency-switched spectra, was performed using the CLASS software package.⁸ The data were smoothed by two spectral channels along the frequency axis to improve sensitivity, for a final spectral resolution of 80 kHz or 0.1 km s⁻¹. A second-order polynomial baseline was fitted to the non-line channels and subtracted from the individual spectra before folding the data. The forward (F_eff) and beam (B_eff) efficiencies of the telescope were interpolated from determined values at 210 GHz and 260 GHz, giving F_eff = 0.92 and B_eff = 0.58 at 231 GHz. The data were then calibrated to main beam temperature via the relation T_MB = (F_eff / B_eff) T*_A. For comparison with the N₂H⁺ 1–0 results in Paper II, we further converted the data to units of Jy beam⁻¹ and convolved the map to a final angular resolution of 18'' to match the N₂H⁺ data.

As the observations were performed in a pooled mode, the 4' × 3' map was observed in three 4' × 1' strips to ensure good sensitivity in potentially limited time. Two-thirds of the map was completed in this manner in winter 2007, while the remaining strip was observed in winter 2008. The sensitivity of the map is consequently not entirely uniform, with the rms noise in the most southern 4' × 1' strip greater (0.15 K (T_MB)) than in the rest of the map (0.08 K (T_MB)). Most emission was found in the central and northern map sections, so this noise increase does not affect our results significantly.

2.2. H₂D⁺ and N₂H⁺ at the JCMT

Simultaneous observations of the H₂D⁺ 1₁₀–1₁₁ and N₂H⁺ 4–3 emission lines at 372.421 GHz and 372.673 GHz, respectively, were performed at the James Clerk Maxwell Telescope (JCMT) over the 07A-08B semesters using the 16-receptor Heterodyne Array Receiver Program B-band receiver (HARP-B; Smith et al. 2008). HARP-B is a 4 × 4 pixel array, with array spacing of 30'' and a 2' × 2' footprint. The JCMT beam at 372 GHz is ∼13'' FWHM. We used the ACSIS correlator as the backend, configured to have a 400 MHz effective bandwidth (e.g., wide enough to encompass both lines), tuned halfway in frequency between the two lines, and a frequency resolution of 61 kHz or 0.05 km s⁻¹ at 372.4 GHz. Observations of five separate HARP footprints were performed in position-switching mode to create an undersampled map of H₂D⁺ and N₂H⁺ over Oph B2. The footprints were placed such that the final beam spacing on the sky is 21'' with a coverage of ∼4' × 2'. The fifth footprint provided additional coverage toward the location of peak continuum emission in B2 but has lower sensitivity compared with the other footprints. The observations were performed in band 1 weather only (τ_225 GHz < 0.05) and generally while Ophiuchus was above 30° in elevation. The mean rms across the H₂D⁺ map is 0.04 K (T*_A) and 0.05 K (T*_A) across the N₂H⁺ map. The main beam efficiency factors, η_MB, of the HARP receptors at 372 GHz are not currently well known. At 345 GHz, the mean η_MB = 0.60, with an rms variation of ∼5% between detectors (Buckle et al. 2009). We therefore estimate η_MB ∼ 0.6 with an uncertainty of ∼10%.

Data reduction was performed using the Starlink software package. Each night's data were first checked for noisy or malfunctioning pixels, and these pixels were flagged individually for each integration. A linear baseline was then removed from each integration, and the data were combined into a final cube with 15'' pixels. To increase sensitivity, the data were averaged along the velocity axis by two channels, giving a final velocity resolution of 0.1 km s⁻¹.

3. RESULTS

We first discuss the integrated intensity of N₂D⁺ 3–2, H₂D⁺ 1₁₁–1₁₀, and N₂H⁺ 4–3 over Oph B2 and compare the species distributions with submillimeter continuum emission and the locations of embedded protostars.

In Figure 1(a), we show the integrated N₂D⁺ 3–2 intensity (convolved to 18'' resolution as described in Section 2.1) in Oph B2 against contours of 850 μm continuum emission observed at the JCMT at 15'' FWHM resolution (original map by Johnstone et al. 2000, combined with all other SCUBA archive data by Jørgensen et al. 2008). Also identified are the locations of 850 μm continuum clumps (Jørgensen et al. 2008) and Class I protostars (classifications and locations taken from Enoch et al. 2009, two of which were previously identified as Elias 32 and 33 (E32 and E33) by Elias 1978; also VSSG 17 and 18 by Vrba et al. 1975). The N₂D⁺ emission is confined to the northern edge of B2 and generally follows the continuum contours in this region, but a very close correspondence is only seen toward the northeast tip of B2. Three maxima of integrated intensity are present, of which two are co-located with continuum clumps while the third is seen to the northwest of the continuum peak. The N₂D⁺ emission is weak toward the continuum emission peak and avoids entirely areas near the protostars.

In Figure 1(b), we plot the N₂D⁺ 3–2 integrated intensity contours over a map of N₂H⁺ 1–0 integrated intensity, observed with the Nobeyama 45 m Telescope at 18'' resolution. The N₂H⁺ observations and analysis were presented in Paper II. We also show the locations of N₂H⁺ clumps, identified through clumpfind in Paper II. Overall, N₂D⁺ emission is found within the extent of N₂H⁺ 1–0 emission. In a similar fashion as the continuum emission, however, we find offsets between the locations of peak N₂D⁺ and N₂H⁺ integrated intensity. In particular, the N₂D⁺ integrated intensity maximum is offset from the integrated N₂H⁺ 1–0 maximum in B2 by ∼20'' to the northeast, ∼1 beam width. The N₂H⁺ 1–0 integrated intensity maximum is itself offset to the east from the continuum peak by ∼22'' (Paper II). The northeastern N₂D⁺ integrated intensity peak, which corresponds well with continuum contours, is also offset from the integrated N₂H⁺ intensity peak by ∼20'' to the east.

We show in Figure 2 H₂D⁺ 1₁₁–1₁₀ spectra across Oph B2 with overlaid 850 μm continuum contours. We find that significant H₂D⁺ emission is present over much of Oph B2, with a larger extent than found in H₂D⁺ observations toward any other cores so far observed (Caselli et al. 2003; Vastel et al. 2006; Pagani et al. 2009a, 2009b). The H₂D⁺ emission again generally follows the continuum emission but avoids the embedded protostars in B2, similar to N₂D⁺. While relatively strong emission is seen toward the continuum peak, the H₂D⁺ integrated intensity maximum (not shown) is also offset the continuum peak by an estimated ∼20''–30'' (due to the map's undersampling, we are unable to exactly determine the location of maximum integrated intensity). In contrast to the N₂D⁺ emission, the H₂D⁺ integrated intensity maximum is found to the northwest rather than the northeast (although N₂D⁺ emission is also found toward the H₂D⁺ peak). The pixel-to-pixel variations seen in Oph B2 may not be significant, however, due to the ∼5% variation determined for the main beam efficiency factors between different HARP receptors (Buckle et al. 2009). We detect some H₂D⁺ emission in northeast B2, where strong N₂D⁺ emission was found, but do not see a significant integrated intensity maximum. (Note, however, that the HARP footprints do not extend past an R.A. of 16:27:35 and are therefore not sensitive to any emission at the tip of the continuum contours in this region.) Overall, while N₂D⁺ emission appears strongest in the north and east of Oph B2, H₂D⁺ emission is the strongest toward B2's western edge.

In Figure 3, we show N₂H⁺ 4–3 spectra across Oph B2 with overlaid 850 μm continuum emission contours. We find bright N₂H⁺ 4–3 emission over most of B2. The greatest integrated N₂H⁺ 4–3 intensity (not shown) is found toward E33, while additional strong N₂H⁺ 4–3 emission is found toward E32 and the 850 μm continuum emission peak. A separate integrated intensity peak is found in northeast B2 (where we also find strong N₂D⁺ 3–2 emission). Overall, while the N₂H⁺ 4–3 emission is generally confined to areas where there is also continuum emission, the strongest N₂H⁺ 4–3 emission is found away from the brightest continuum contours. A comparison between Figure 3 and the N₂H⁺ 1–0 integrated intensity in Figure 1(b) (gray scale) reveals a significant discrepancy between the locations of strong N₂H⁺ 1–0 and 4–3 emission in Oph B. We probe the possible causes of the offset, including the excitation and opacity of the lines, further below.

4. ANALYSIS

4.1. Line Fitting

The rotational transitions of the linear molecular ions N₂H⁺ and N₂D⁺ contain hyperfine structure, which is dominated by the interactions between the molecular electric field gradient and the electric quadrupole moments of the two nitrogen nuclei (Dore et al. 2004). Rotational transitions at higher J contain greater numbers of hyperfine components at smaller frequency intervals, i.e., the satellite components are found more closely spaced and can overlap significantly. In particular, the N₂H⁺ 4–3 and the N₂D⁺ 3–2 transitions contain 21 and 25 hyperfine components, respectively, with relative line strengths ≳0.001 (Pagani et al. 2009a, 2009b; Gerin et al. 2001).

For N₂H⁺ 4–3, we calculate the relative line strengths $s_{\it J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime }$ from the Einstein A coefficients $A_{\it J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime }$ and frequencies $\nu _{\it J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime }$ given by Pagani et al. (2009a, 2009b), using

$\begin{equation} A_{J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime } = \frac{64\pi ^4}{3hc^2} \mu ^2 \nu _{J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime }^3 \frac{J}{[F]} s_{J\,F_1\,F \rightarrow J^\prime \,F_1^\prime \,F^\prime }, \end{equation} \tag{ 2 }$

where μ = 3.4 D is the permanent electric dipole moment, J is the upper rotational quantum number, and [F] = 2F + 1 accounts for the degeneracy of the hyperfine states. For the N₂H⁺ 4–3 transition, 14 components at 8 distinct frequencies make up the main centroid emission and lie within 0.4 km s⁻¹ of each other. At low optical depths the line profile does not differ substantially from a Gaussian profile, while the relative intensities of the satellite lines to the main line component remain small even at large opacities (∼20% at τ ∼ 10). Nevertheless, the satellite components of the line are visible toward several pixels in Oph B (identified as locations A, B, and D on Figure 3, and shown in Figure 4), largely near the embedded protostar E33. At some locations, the line structure is confused with a significant blue emission shoulder, which is not due to any hyperfine structure, but which may be associated with an outflow stemming from either E33 or E32 (Kamazaki et al. 2003; G. White et al. 2010, in preparation). We were thus able to fit well the hyperfine structure toward only one location in the map and constrain the excitation temperature, T_ex, and opacity, τ, toward two others.

**Figure 4.** N₂H⁺ 4–3 spectra toward four locations in Oph B2, labeled A, B, C, and D in Figure 3: 30'' north of E33 (A), toward E33 (B), 21'' northeast of E33 (C), and ∼42'' northeast of E33 (D). Overlaid on the spectra are fits to the full hyperfine structure of the N₂H⁺ 4–3 line (A, B, and D) and a two-component Gaussian fit (C).
Download figure:
Standard image High-resolution image

We fit the hyperfine structure with a custom Gaussian hyperfine structure fitting routine written in IDL, described in detail in Paper I, assuming equal excitation temperatures T_ex for each hyperfine component. The line fitting routine does not take into account the significant overlap between hyperfine components in both the N₂H⁺ and N₂D⁺ transitions, but we find the total N₂H⁺ line opacities are still small where the satellite structure is seen (τ ∼ 2–3) and we therefore do not expect line overlap to introduce substantial errors.

At low opacities, we cannot solve for τ and T_ex independently. Given a reasonable estimate for T_ex, however, we may attempt to fit the line without visible satellite structure. In Paper II, we determined T_ex for N₂H⁺ 1–0 emission across Oph B, and to first order can expect the same T_ex to describe the excitation of all rotational levels of N₂H⁺. This assumption was found to be accurate in a recent study of N₂H⁺ and N₂D⁺ emission in L1544, where Caselli et al. (2002b) were able to determine separately the excitation temperatures for N₂H⁺ 1–0, 3–2 and N₂D⁺ 2–1 and 3–2 and found that they were consistent within uncertainties. We find, however, that while most of the data can be fit reasonably well using the N₂H⁺ 1–0 T_ex values, visual inspection showed that a simple Gaussian fits the data equally well, or significantly better, in many cases, and the line opacities were not well constrained.

Because multiple components of similar line strength overlap to form the main N₂H⁺ 4–3 component, a Gaussian fit will overestimate the line width, while the returned v_LSR will be shifted relative to a full hyperfine structure (HFS) fit. The shift in v_LSR can be corrected (Δv_LSR = 0.036 km s⁻¹), and at low opacities, the increase in Δv is small (≲10% for τ < 1, up to ∼25% for τ ∼ 2). In the following sections, we thus discuss only the results of the Gaussian fitting of the N₂H⁺ 4–3 emission unless otherwise stated.

No hyperfine structure was seen above the rms noise in N₂D⁺ 3–2 emission across Oph B. We are therefore unable to fit well the N₂D⁺ 3–2 hyperfine structure solving for both τ and T_ex. Here, however, equating the N₂D⁺ T_ex with that found for N₂H⁺ 1–0 emission resulted in better constraints on the opacity. The hyperfine fit, using the relative component strengths and velocities from Gerin et al. (2001), is thus more informative than a simple Gaussian fit, and we discuss the results from the full HFS fit in the following sections.

To properly calculate the line width, Δv, we subtracted in quadrature the resolution width, Δv_res (0.1 km s⁻¹ for all species), from the observed line width, Δv_obs, such that $\Delta v = \sqrt{\Delta v_{\rm obs}^2 - \Delta v_{\rm res}^2}$ . The uncertainties reported in the returned parameters are those determined by the fitting routine, not taking the calibration uncertainty into account. The calibration uncertainty does not affect the uncertainties returned for v_LSR or Δv. The excitation temperature, column densities, and fractional abundances discussed further below, however, are dependent on the amplitude of the line emission and thus affected by the absolute calibration uncertainty.

The H₂D⁺ 1₁₁–1₁₀ emission line has no hyperfine structure. We consequently fit the emission at each pixel with a single Gaussian (again in idl) to determine the line amplitude, v_LSR and Δv across Oph B2.

4.2. Centroid Velocity and Line Widths

Table 2 lists the mean, rms, minimum, and maximum values of v_LSR and Δv for N₂H⁺ 4–3, N₂D⁺ 3–2, and H₂D⁺ 1₁₁–1₁₀, along with the same results from NH₃ (1,1) and N₂H⁺ 1–0 presented in Paper I and Paper II for comparison.

Table 2. Line Properties of NH₃, N₂H⁺, N₂D⁺, and H₂D⁺ in Oph B2

Species	Transition	v_LSR				Δv				n_cr
		(km s⁻¹)				(km s⁻¹)				(cm⁻³)
		Mean	rms	Min	Max	Mean	rms	Min	Max
NH₃	(J, K) = (1, 1)	4.14	0.19	3.85	4.58	0.87	0.16	0.42	1.11	2 × 10⁴ ^a
N₂H⁺	J = 1–0	4.04	0.20	3.61	4.43	0.71	0.22	0.38	1.42	1.4 × 10⁵ ^b
H₂D⁺	$J_{K_{-1}K_1}$ = 1₁₁–1₁₀	4.01	0.24	3.65	4.50	0.74	0.19	0.51	1.36	∼10⁵ ^c
N₂D⁺	J = 3–2	4.02	0.19	3.67	4.44	0.53	0.20	0.25	1.36	8 × 10⁵ ^d
N₂H⁺	J = 4–3^e	4.05	0.22	3.39	4.43	0.69	0.18	0.47	1.14	7.7 × 10⁶ ^b

Notes. The velocity resolution of the NH₃ observations was 0.3 km s⁻¹. All other species were observed with 0.1 km s⁻¹ velocity resolution or better. ^aCaselli et al. (2002a). ^bSchöier et al. (2005; n_cr calculated using Einstein A coefficients and collision rates γ_ul at 10 K from the Leiden Atomic and Molecular Database). ^cCaselli et al. (2008). ^dCalculated using Einstein A coefficients and collision rates for N₂H⁺ 3–2 at 10 K (Schöier et al. 2005). Daniel et al. (2007, see their Figure 17) found that despite small differences in Einstein A coefficients and collision rates, critical densities of N₂H⁺ and N₂D⁺ are similar for transitions other than J = 1–0. ^eValues calculated from a Gaussian fit due to a lack of visible hyperfine structure. The reported line widths may be greater by ∼10%–25% than the true values; see Section 4.1.

Download table as: ASCII Typeset image

We show in Figure 5 the N₂D⁺ v_LSR and Δv over Oph B2. We only plot values where S/N of the peak line intensity was ⩾5. The mean v_LSR across B2 is 4.02 km s⁻¹ with an rms variation of 0.19 km s⁻¹. Lower v_LSR values are found toward the western and northeastern regions of the core, with the maximum v_LSR = 4.44 km s⁻¹ in the east toward the NH₃ (1,1) clump B2-A7 (identified in Paper I and labeled on Figure 5; coincident with the N₂H⁺ 1–0 clump B2-N6 from Paper II). The variation in v_LSR does not suggest global rotation of B2. The N₂D⁺ line velocities shift between 3.7 km s⁻¹ and 4.4 km s⁻¹ over ∼40'', or ∼4800 AU, to the west of E32. Nearby, the line v_LSR shifts from 3.7 km s⁻¹ to 4.2 km s⁻¹ on either side of the western embedded protostar.

The mean N₂D⁺ 3–2 line width Δv = 0.53 km s⁻¹ with an rms variation of 0.20 km s⁻¹. In general, the widest line widths are found near E32 and the western protostar (Oph-emb5, Enoch et al. 2009), to a maximum Δv = 1.36 km s⁻¹ in western B2 (but note that the largest Δv values are found only ∼45'' to the northwest of some of the narrowest Δv values). The locations of large line widths correspond well with the locations of dramatic v_LSR shifts described above. Away from the protostars, only small variations in the relatively narrow line widths are seen, notably toward the northeast and southwest, with Δv ∼ 0.3 km s⁻¹. Two objects are associated with localized narrow Δv. An 850 μm continuum clump (162725-24273, Jørgensen et al. 2008), shown on the image, was identified where the southwestern narrow line widths are found. In the east, B2-A7 is also associated with an NH₃, a Δv minimum of 0.36 km s⁻¹. This is similar to the clump line width in N₂H⁺ 1–0 emission at 18'' resolution, but nearly a factor of two larger than found in N₂H⁺ 1–0 interferometer observations of the clump presented in Paper II.

In Figure 6(a) and Figure 6(c), we show the v_LSR determined through Gaussian fitting of the H₂D⁺ 1₁₁–1₁₀ and N₂H⁺ 4–3 lines, respectively. Visual inspection of the v_LSR values shows that the H₂D⁺ emission is similar kinematically to the N₂D⁺ 3–2 values in Figure 5(a), while a small shift in v_LSR of ∼+0.2 km s⁻¹ is seen in N₂H⁺ 4–3 emission toward Elias 32 and the 850 μm continuum emission peak and co-located 850 μm continuum emission peak (co-located with the continuum clump B2-MM8). The shift is small but significant given the otherwise excellent agreement in v_LSR between the different tracers over B2 (see Table 2).

We show in Figure 6(b) and Figure 6(d) the Δv of the H₂D⁺ and N₂H⁺ lines. The H₂D⁺ line widths are large, with a mean Δv = 0.74 km s⁻¹, and do not vary greatly over Oph B2 (Δv rms = 0.19 km s⁻¹). Slightly greater line widths are found in the western half of B2. Gaussian fits to the N₂H⁺ 4–3 data find that the lines are slightly more narrow, with a mean Δv = 0.69 km s⁻¹. As noted in Section 4.1, Gaussian fits likely overestimate the true N₂H⁺ line width by ∼10%–25%. Wider N₂H⁺ lines are also found toward western Oph B2, while narrower lines are seen toward the submillimeter continuum emission peak and in B2's eastern half. N₂H⁺ lines of average width are found toward the embedded protostars but become wider at pixels directly offset from the infrared sources. In particular, directly to the east of E32 and E33, the N₂H⁺ 4–3 line profiles show a significant blue shoulder with emission extending to v_LSR ∼2.5 km s⁻¹. Figure 4(c) shows the result of a two-component Gaussian fit to the line and blue shoulder at location C in Figure 3.

4.3. N₂H⁺ Opacity and Excitation Temperature

Comparison of the N₂H⁺ 1–0 and 4–3 line emission (see Figure 1(b) and Figure 3) shows that while in both cases emission is found over much of Oph B2, the integrated intensity distributions of the two transitions are quite different. It is thus unclear whether the excitation conditions for both transitions would be equal.

Toward several pixels, hyperfine structure is visible. These locations have been labeled A, B, and D in Figure 3, and these spectra are again shown in Figure 4. Using the full hyperfine component fit (overlaid on the spectrum), we were able to determine reasonably well the excitation temperature and opacity of the N₂H⁺ 4–3 line, T_ex,4–3 and τ_4–3 toward location A, with the resulting best-fit total line opacity τ_4–3 = 2.5 ± 1 and T_ex,4–3 = 8 ± 1 K. Both T_ex,4–3 and τ_4–3 are approximately consistent with values found in Paper II for N₂H⁺ 1–0 emission at the same location, where T_ex,1–0 = 7.3 ± 0.1 K and τ_1–0 = 3.1 ± 0.6. The brightest N₂H⁺ 4–3 line emission, however, and additionally the maximum integrated intensity, is located toward E33 (location B in Figure 3, directly south of A). The intensity ratio of the main N₂H⁺ 4–3 component with an apparent red satellite (line strength only 2× the data rms noise level, but at the correct velocity and with a width that matches the main component, see Figure 4(b)) suggests that the N₂H⁺ 4–3 opacity toward E33 remains moderately optically thick, τ_4–3 ∼ 3, and the corresponding T_ex,4–3 ∼ 8. The expected line profile given these values is overlaid on the spectrum in Figure 4. The integrated N₂H⁺ 1–0 emission decreases significantly toward E33 (see Figure 1(b) and Paper II). Additionally, the emission is confused by multiple velocity components which preclude robust fits to the 1–0 spectra, so we do not have a good estimate of the N₂H⁺ 1–0 T_ex at this location. At location D in Figure 3, satellite emission is again seen, but here the hyperfine fits are confused by significant blue shoulder emission. The best-fit results give T_ex ∼ 7.5 K–10 K, with τ ∼ 1 − 3. We overlay the expected line profile for T_ex = 8.5 K and τ = 2 on spectrum D in Figure 4. Elsewhere, no hyperfine structure is visible in the N₂H⁺ 4–3 line.

For opacity ratios between the N₂H⁺ 4–3 and N₂H⁺ 1–0 lines, τ_4–3/τ_1–0 ∼ 0.3 and similar line widths, the derived N₂H⁺ column density, N(N₂H⁺), calculated using either transition following Caselli et al. (2002c, Equation (A1)), is equal if both T_ex,1–0 and T_ex,4–3 ∼ 7 K. At lower opacity ratios, the required T_ex also decreases. On average, the N₂H⁺ 1–0 and 4–3 line widths agree, however we note that since we cannot fit well the hyperfine structure, the N₂H⁺ 4–3 widths are an overestimate of the true value by ∼10%–25%. This uncertainty propagates linearly into the column density calculation. With this caveat, an excitation temperature T_ex ∼ 7 K agrees well with the mean T_ex determined from HFS fitting of the N₂H⁺ 1–0 line in Paper II. We have no column density estimate toward E33 from N₂H⁺ 1–0 emission. Near E32 and E33, however, the greater relative N₂H⁺ 4–3 opacity requires either a higher common T_ex ∼ 10–11 K to obtain similar N(N₂H⁺) estimates, which does not agree with our N₂H⁺ 1–0 T_ex measurement, or the excitation of the N₂H⁺ 4–3 line must be greater than that of the 1–0 line by ∼2–3 K. This moderate increase in excitation temperature of the higher excitation N₂H⁺ line is not ruled out by our data. Detailed modeling of the N₂H⁺ emission near the embedded protostars using radiative transfer codes would help disentangle the excitation and opacity effects.

4.4. Column Density and Fractional Abundance

We next discuss the column densities of N₂D⁺ and H₂D⁺ derived from the observed line emission. The N₂H⁺ column density was calculated in Paper II.

4.4.1. N₂D⁺

We calculate the total N₂D⁺ column density, N_tot(N₂D⁺), from the integrated intensity of the 3–2 transition following Caselli et al. (2002a) for optically thin emission:

$\begin{eqnarray} N_{\rm tot} &=& \frac{8 \pi W}{\lambda ^3 A} \frac{g_u}{g_l} \frac{1}{J_\nu (T_{{\rm ex}}) - J_\nu (T_{{\rm bg}})} \frac{1}{1 - \exp (-h\nu /kT_{{\rm ex}})} \cr &&\times\,\frac{Q_{\rm rot}(T_{{\rm ex}})}{g_l \exp (-E_l/kT_{{\rm ex}}) }, \end{eqnarray} \tag{ 3 }$

where W is the line integrated intensity, A is the Einstein spontaneous emission coefficient, g_u and g_l are the statistical weights of the upper and lower states, respectively, and the equivalent Rayleigh–Jeans excitation and background temperatures are given by J_ν(T_ex) and J_ν(T_bg). The rotational energy E_J = J(J + 1)hB for linear molecules, where the rotation quantum number J = 2 for the lower energy state and B = 38554.717 MHz is the rotational constant (Sastry et al. 1981). The partition function Q_rot(T_ex) = ∑^∞_{J = 0}(2J + 1)exp(−E_J/kT_ex), and we use the integrated value, Q_rot(T_ex) = kT_ex/(hB) + 1/3. We assume that the excitation temperature T_ex is the same for all rotational levels, and use the fitted N₂H⁺ 1–0 T_ex from Paper II as an estimate of the excitation temperature of the N₂D⁺ line at each pixel (as discussed in Section 4.1). All other parameters remaining equal, at the average T_ex of the N₂H⁺ 1–0 line, a decrease in T_ex by 2 K would increase the resulting N₂D⁺ column density by a factor ≲2, while an increase in T_ex by 2 K would decrease the resulting N₂D⁺ column density by ≲0.5. We show N(N₂D⁺) across Oph B2 in Figure 7(a).

**Figure 7.** (a) Total N₂D⁺ column density N(N₂D⁺) (cm⁻²) in Oph B2 determined using Equation (3) (gray scale). Values shown have been divided by 10¹¹. In both plots, only results from pixels where S/N of the peak line intensity is ⩾5 are shown, indicated by the thick black contour. Gray contours show 850 μm continuum emission as in Figure 5. The 18'' (FWHM) beam is shown at lower left. Stars indicate the positions of Class I protostars. (b) Fractional abundance of N₂D⁺ relative to N₂H⁺ and *R_D*, in Oph B2. The scale has been truncated to show the variation in central B2; the maximum *R_D* = 0.16 toward the eastern tip.
Download figure:
Standard image High-resolution image

We then calculate the fractional N₂D⁺ abundance per pixel,

$\begin{equation} X(\mbox{N$_2$D$^+$}) = N(\mbox{N$_2$D$^+$}) / N(\mbox{H$_2$}), \end{equation} \tag{ 4 }$

where the H₂ column density, N(H₂), is calculated from 850 μm continuum emission. The 15'' data were convolved to a final 18'' FWHM resolution to match the convolved N₂D⁺ data. The H₂ column density is then given by N(H₂) = S_ν/[Ω_mμm_Hκ_νB_ν(T_d)], where S_ν is the 850 μm flux density, Ω_m is the main-beam solid angle, μ = 2.33 is the mean molecular weight, m_H is the mass of hydrogen, κ_ν is the dust opacity per unit mass at 850 μm, where we take κ_ν = 0.018 cm² g⁻¹, and B_ν(T_d) is the Planck function at the dust temperature, T_d. As in Paper I and Paper II, we expect the gas and dust to be well coupled in Oph B2, and set T_d = T_K, the gas kinetic temperature determined from NH₃ (1,1) and (2,2) inversion line ratios. We refer the reader to the detailed discussion of the uncertainties in the N(H₂) values presented in Paper II, but note that the derived H₂ column densities have uncertainties of factors of a few. Additionally, due to the negative artifacts introduced into the data from the observational chopping technique (again discussed in detail in Paper II), we limit the analysis to pixels where S_ν ⩾ 0.1 Jy beam⁻¹, though the rms noise level of the continuum map is ∼0.03 Jy beam⁻¹. For a dust temperature T_d = 15 K, this flux level corresponds to N(H₂) ∼ 6 × 10²¹ cm⁻².

The column density distribution is not similar to the integrated intensity map shown in Figure 1(a). Figure 7(a) shows several small peaks of N(N₂D⁺) in B2 which largely correspond with local minima in the excitation temperature (T_ex ∼ 5 K compared with the mean T_ex ∼ 7 K, where the uncertainty in T_ex is ∼0.4–0.6 K). The N(N₂D⁺) peaks remain within the continuum contours but avoid maxima of continuum emission. We find a mean N(N₂D⁺) = 1.8 × 10¹¹ cm⁻² in Oph B2 with an rms variation of 1.3 × 10¹¹ cm⁻². The maximum N(N₂D⁺) = 7.1 × 10¹¹ cm⁻² is found toward the northeast. If we calculate N(N₂D⁺) assuming a constant T_ex = 7 K, the distribution follows the integrated intensity but the mean value and variation remain the same.

The fractional N₂D⁺ abundance (not shown) follows a similar distribution to the N₂D⁺ column density. We find a mean X(N₂D⁺) = 8.2 × 10⁻¹² with an rms variation of 5.1 × 10⁻¹². The maximum X(N₂D⁺) = 2.9 × 10⁻¹¹ is again found toward northeast Oph B2. Within the continuum emission contours, the lowest measured N₂D⁺ abundances (X(N₂D⁺) ∼ 2 × 10⁻¹²) are found toward the continuum emission peak and B2-MM8 clump (emission toward the Elias 33 protostar was not above our S/N limit for the analysis).

4.4.2. H₂D⁺

To calculate the H₂D⁺ column density, N(H₂D⁺), we first estimate the line opacity τ from the H₂D⁺ 1₁₁–1₁₀ line emission using

$\begin{equation} \tau = - \ln \biggr [ 1 - \frac{T_{{\rm MB}}}{J(T_{{\rm ex}}) - J(T_{{\rm bg}})} \biggl ], \end{equation} \tag{ 5 }$

where we must assume a priori an excitation temperature T_ex for the transition. An upper limit is given by the gas kinetic temperature T_K. For the pixels where the H₂D⁺ S/N >5, we found in Paper I a mean T_K = 14.5 K with an rms variation of 1.1 K (slightly less than the mean T_K = 15 K found for the entire Oph B core in Paper I), with a minimum T_K ≲ 13 K toward the highest N(H₂) values. We note, however, that the critical density at which NH₃ is excited (n_cr ∼ 10³–10⁴ cm⁻³) is lower than that of the H₂D⁺ transition by an order of magnitude, and the NH₃ temperatures, or the variation in T_K, may not accurately reflect the conditions where H₂D⁺ is excited.

In their recent H₂D⁺ study of several starless and protostellar cores, Caselli et al. (2008) found that the H₂D⁺ T_ex was similar to or slightly less than the core T_K values, which were estimated through a variety of methods (e.g., continuum observations, NH₃ inversion observations as for B2, and multiple transition observations of molecules such as HCO⁺). Six of the seven protostellar cores in their sample had excitation temperatures which ranged from 9 K to 14 K and similar densities to those calculated for B2 in Paper II. Based on these values, and to better compare our results with Caselli et al., we use a constant excitation temperature T_ex = 12 K in Equation (5). Given T_ex = 12 K, we find a mean τ = 0.13, with a minimum τ = 0.05 and maximum τ = 0.28, which agrees well with the Caselli et al. results in protostellar cores.

The total column density of ortho-H₂D⁺ is then given by Vastel et al. (2006):

$\begin{equation} N(\mbox{ortho-H$_2$D$^+$}) = \frac{8\pi }{\lambda ^3A_{ul}}\frac{Q_{\rm rot}(T_{{\rm ex}})}{g_u} \frac{\exp (E_u/kT_{{\rm ex}})}{\exp (h\nu /kT_{{\rm ex}})-1} \int \tau {\rm d}v, \end{equation} \tag{ 6 }$

where g_u = 9, E_u/k = 17.8 K, and A_ul = 1.08 × 10⁻⁴ s⁻¹ for the 1₁₁–1₁₀ transition (Ramanlal & Tennyson 2004). The integral $\int \tau \mbox{d}v = \frac{1}{2}\sqrt{\pi /(\ln (2)} \tau \Delta v$ . We calculate Q_rot by reducing the H₂D⁺ level structure to a two-level system following Caselli et al. Since the energy of the first excited state above ground is E/k = 17.8 K and that of the second excited state is E/k = 110 K, we expect this approximation to be valid in Oph B2 given the low temperatures determined in Paper I. The partition function depends on T_ex, which we have estimated to be 12 K. In Figure 8, we plot the variation in τ, Q_rot, and N(ortho-H₂D⁺) with T_ex given a line T_MB = 0.5 K and Δv = 0.7 km s⁻¹ (typical values in B2) and have listed the returned values for T_ex between 7 K and 15 K in Table 3. While the variation in returned parameters is fairly large with small T_ex changes at low T_ex values (i.e., factors ∼2 or more around T_ex ∼ 7 K), at higher T_ex the values vary less. For example, an increased T_ex = 14 K rather than T_ex = 12 K would result in a decreased N(ortho-H₂D⁺) by ∼20%, while a decreased T_ex = 10 K would result in an increased N(ortho-H₂D⁺) by ∼30%. We thus assume an uncertainty of ∼25% in the N(H₂D⁺) values reported here.

**Figure 8.** Impact of the assumed excitation temperature, T_ex, on returned H₂D⁺ 1₁₁–1₁₀ line opacity (top), partition function Q_rot (middle), and column density N(ortho-H₂D⁺)(bottom), given a line brightness temperature T_MB and line width Δv similar to observed mean values in Oph B2. At low T_ex, τ and N(ortho-H₂D⁺) vary strongly with T_ex, but the effect of changes in T_ex at higher temperatures becomes small.
Download figure:
Standard image High-resolution image

Table 3. Impact of T_ex on τ, Q_rot, and N(ortho-H₂D⁺)

T_ex	τ	Q_rot	N(ortho-H₂D⁺)
(K)			(10¹³ cm⁻²)
7	0.41	9.7	1.60
8	0.27	10.0	1.11
9	0.20	10.2	0.86
10	0.15	10.5	0.71
11	0.12	10.8	0.61
12	0.10	11.0	0.54
13	0.09	11.3	0.48
14	0.08	11.5	0.45
15	0.07	11.7	0.42

Note. Values were calculated for a line strength T_MB = 0.5 K and width Δv = 0.7 km s⁻¹.

Download table as: ASCII Typeset image

To calculate the total H₂D⁺ column density in B2, an estimate of the ortho- to para-H₂D⁺ (o/p-H₂D⁺) ratio is needed. The ground state para-H₂D⁺ transition occurs at a frequency of 1.37 THz; this line is difficult to observe from the ground and may not be excited at the low temperatures found in star-forming cores. Without a direct measurement of the population of the para state, we look to chemical models to predict the o/p-H₂D⁺ ratio, which show that the o/p-H₂D⁺ ratio is directly dependent on o/p-H₂. Walmsley et al. (2004) find that the o/p-H₂ ratio increases to a steady-state value of ∼5 × 10⁻⁴ at T = 10 K and n = 10⁶ cm⁻³, similar to the density and temperature determined for Oph B2 in Paper I and Paper II. The same models predict o/p-H₂D⁺ ∼0.2–0.3. Using Bonnor–Ebert sphere core models, Sipilä et al. (2010) also predict o/p-H₂D⁺ ratios ≲0.5 at T ∼ 10–15 K in steady state. Their models show that the o/p-H₂D⁺ ratio can dramatically increase to values ∼2–3 at T ≲ 10 K, but these temperatures are unlikely based on our analysis of Oph B2.

It is not clear, however, whether steady-state models accurately represent conditions in star-forming cores like B2. The dynamical timescale of a collapsing core can be much shorter than some of the chemical timescales, leading to large departures of free-fall model predictions from steady-state models (Flower et al. 2006). For example, Flower et al. show that the steady-state model underpredicts the o/p ratio by up to an order of magnitude compared with a free-fall collapse model at gas densities n ≲ 10⁶ cm⁻³ (the models are more consistent at higher densities). Consequently, it may not be unreasonable to expect variations in the o/p-H₂D⁺ ratio across B2 as it evolves. Given this uncertainty in the o/p-H₂D⁺ ratio, we limit discussion to N(ortho-H₂D⁺), but note that the total N(H₂D⁺) is likely greater than N(ortho-H₂D⁺) by a factor of a few.

We show in Figure 9 the resulting N(ortho-H₂D⁺) distribution across Oph B2. Including only pixels where S/N >5, we find the mean N(ortho-H₂D⁺) = 1.4 × 10¹³ cm⁻² with an rms variation of 0.7 × 10¹³ cm⁻² over Oph B2. To this sensitivity limit, we find a minimum N(ortho-H₂D⁺) = 4.1 × 10¹² cm⁻². The maximum N(ortho-H₂D⁺) = 3.3 × 10¹³ cm⁻² is found ∼30'' to the northwest of the 850 μm continuum peak. The largest ortho-H₂D⁺ column densities are found in the western half of Oph B2, with moderate N(ortho-H₂D⁺) values extending to the northeast.

We next calculate X(ortho-H₂D⁺) as described for N₂D⁺ above and show in Figure 9(b) the X(ortho-H₂D⁺) distribution across B2. Similar to the column density distribution, the largest H₂D⁺ abundances, X(ortho-H₂D⁺) ∼ (4–5) × 10⁻¹⁰, are found in the west. A local X(ortho-H₂D⁺) minimum of ∼2 × 10⁻¹⁰ is found toward the central continuum emission peak, with similar abundances toward the northeast.

4.5. Using H₂D⁺, N₂H⁺ , and N₂D⁺ to Determine T_K

The molecular weight of H₂D⁺ (m_mol = 4.01 m_H) is much less than that of N₂H⁺ or N₂D⁺ (29.01 m_H and 30.02 m_H, respectively). If we assume that the H₂D⁺ 1₁₁–1₁₀, N₂H⁺ 4–3, or N₂D⁺ 3–2 transitions are excited in the same environment within Oph B2, we would expect each to trace the same core motions and be described by the same kinetic temperature T_K. This co-location would imply that the non-thermal line widths for the two species should be equal, where σ²_NT = σ²_obs − σ²_th. Given equal T_K, the H₂D⁺ thermal line width, σ_th = (kT_k/m_mol)^1/2, will be measurably larger than that of N₂H⁺ due to its smaller molecular weight. At T_K = 12 K, for example, σ_th = 0.06 km s⁻¹ for N₂D⁺ and 0.16 km s⁻¹ for H₂D⁺. We can use this difference in thermal line widths to determine the kinetic gas temperature by equating the non-thermal line widths and solving for T_K (e.g., in the case of N₂H⁺):

$\begin{equation} T_K = \frac{1}{k_B}\biggl (\frac{1} {m_{\rm H_2D^+}} - \frac{1}{m_{\rm N_2H^+}}\biggr)^{-1} \big[\sigma _{\rm obs}^2\big(\mbox{H$_2$D$^+$}\big) - \sigma _{\rm obs}^2\big(\mbox{N$_2$H$^+$}\big)\big] \end{equation} \tag{ 7 }$

Mean temperatures of 9.2 ± 0.2 K, in good agreement with temperatures derived from NH₃ alone, were derived in this manner toward a number of dense cores by Fuller & Myers (1993) using the molecular species HC₃N and NH₃. The largest source of error in the derived T_K likely arises from the assumptions listed above rather than from the propagation of the small uncertainties on the observed line widths.

Overall, we expect the comparison of H₂D⁺ and N₂D⁺ to be most accurate, given the similar critical densities of the transitions, the fact that N₂D⁺ must form from H₂D⁺, and the similar distribution of integrated intensity (in contrast with the poor correlation between H₂D⁺ and N₂H⁺ 4–3 emission; see Figures 1(a) and 2). The different distributions of H₂D⁺ and N₂H⁺ 4–3 in Oph B2 suggest that these two species are not excited in the same environment in the core. Additionally, since our estimates of the N₂H⁺ line widths are not accurate (given the unresolved hyperfine structure), we will focus only on N₂D⁺ for this analysis.

We show in Figure 10 the temperatures calculated based on Equation (7) for H₂D⁺ and N₂D⁺ 3–2. Only pixels where S/N of both lines was >5 are shown. The derived T_K values span a large range, from unphysical negative values (∼−11 K) to unrealistically large values (∼75 K). Using the mean Δv values of H₂D⁺ and N₂D⁺ from Table 2, we calculate a mean T_K = 16 K, similar to that determined through NH₃ analysis. Multiple pixels with T_K ∼ 12 K are found surrounding the central continuum emission peak. Negative temperatures, however, are found associated with the continuum peak and associated clump B2-MM8.

From Equation (7), we calculate negative T_K values when the observed H₂D⁺ line widths are less than those of N₂D⁺ 3–2. We show in Figure 11 (bottom) Δv derived for H₂D⁺ and N₂D⁺ emission (data values in 15'' pixels plotted), along with the spectra and relative fits toward B2-MM8 (top). Over most of B2, N₂D⁺ line widths are narrower than H₂D⁺ and thus T_K > 0 K. The negative temperatures derive from two adjacent pixels toward B2-MM8 where the N₂D⁺ Δv > Δv of H₂D⁺. Figure 11 shows non-Gaussianity in the spectra of both species toward B2-MM8, which complicates the calculation of a single line width to describe the emission.

Why would we find greater N₂D⁺ Δv than H₂D⁺ Δv only toward the central continuum peak? One possibility is suggested by the spectra in Figure 11, where additional emission in the H₂D⁺ spectrum is clearly seen blueshifted relative to the main Gaussian peak. For the H₂D⁺ line, the best-fit Gaussian does not include the blue line shoulder. Slight blue asymmetry is also visible in the N₂D⁺ line, but the distinction between the main component and the shoulder is less obvious. It is possible that the shoulder emission of N₂D⁺ 3–2 is spread over a larger range in v_LSR due to the contribution of multiple hyperfine components to the observed emission (each of which would contain emission at the central line velocity and at the shoulder velocity), broadening the line overall. The best fit to the N₂D⁺ line thus cannot distinguish well between the main and shoulder emission, resulting in an artificially broad Δv.

5. DISCUSSION

5.1. Line Widths and Density

The mean line widths listed in Table 2 show that on average the observed Δv for a given molecular line decreases with increasing critical line density, suggesting that the gas in B2 becomes more quiescent at higher densities (omitting N₂H⁺ 4–3, where Δv is an overestimate). For an average gas temperature T_K = 15 K, the resulting non-thermal line contribution decreases from σ_NT = 0.36 km s⁻¹ for NH₃ (1,1) emission to σ_NT = 0.22 km s⁻¹ for N₂D⁺ 3–2 emission (from n_cr ∼ 10³–10⁴ cm⁻³ to n_cr ∼ 8 × 10⁵ cm⁻³, respectively). Comparing with the thermal sound speed at 15 K, where $c_s = \sqrt{kT_K/(\mu m_{{\rm H}})}$ , we find an average σ_NT / c_s = 1.6 for NH₃ and σ_NT / c_s = 1.1 for N₂D⁺ 3–2 emission. This comparison shows that despite the decrease in non-thermal motions, the non-thermal motions are equal in magnitude with the sound speed, even at high density. A lower gas temperature would result in even larger magnitudes of σ_NT, while significantly higher temperatures are unlikely.

This decrease in non-thermal motions at higher densities has been observed in starless cores (e.g., Pon et al. 2009; Tafalla et al. 2004; Lada et al. 2003). Frequently in starless cores, however, the narrowest observed line widths are consistent with being due to nearly pure thermal motions, in contrast to the results found here. As σ_NT values decrease, the ability for turbulent pressure to support the core against collapse is reduced, which impacts the stability of the core against its self-gravity. The narrow σ_NT widths found in starless cores are often explained as the result of the dissipation of turbulent motions in the core center (Goodman et al. 1998). If the non-thermal line widths in B2 are due to turbulent motions, the corresponding turbulent pressure (P_NT = mnσ²_NT) is equal to or greater than the thermal pressure (P_T = nkT).

In Paper I, we discussed the source of wide NH₃ line widths in Oph B and concluded that if the relatively large σ_NT values are caused by turbulent motions, they could not be primordial (i.e., inherited from the parent cloud). This is because the damping timescale of such motions (with no external driving force) is approximately the dynamical time (Mac Low & Klessen 2004), and the existence of protostars associated with B2 implies that B2 has likely existed, with high density gas, for at least a dynamical time.

The wide lines in Oph B2 may, however, be due to more organized motions such as infall or outflow. None of the observed molecular lines are optically thick, and we therefore do not see the self-absorbed, asymmetric line profiles frequently used to infer the infall of gas within a core. Infall signatures were observed in CS, H₂CO, or HCO⁺ emission toward several locations in eastern B2 by André et al. (2007), but not in central or western B2. Gurney et al. (2008) observed complicated CO (and CO isotopologues) line structures toward the B2 continuum peak but were unable to distinguish infall conclusively. The analysis of line profiles is complicated by a substantial outflow discovered in CO emission and emanating from one of the two embedded YSOs in B2 (Kamazaki et al. 2003) which extends over 10' (∼0.4 pc) in size.⁹ Upcoming results of CO, ¹³CO, and C¹⁸O 3–2 observations in Ophiuchus (G. White et al. 2010, in preparation) from the JCMT Gould Belt Legacy Survey (GBLS; Ward-Thompson et al. 2007) will be better able to study the infall and outflow motions in Oph B. In conjunction with these data, the process of how the outflow motions affect the dense gas may be studied. There is clear evidence, however, for protostellar influence on the N₂D⁺ line widths shown in Figure 5(b), as discussed in Section 4.2. Additionally, the blue shoulder emission seen in N₂H⁺ and N₂D⁺ east of E32 and E33 coincides with the eastern CO lobe of the outflow and is further evidence for the significant impact of the outflow on the dense gas.

We also noted in Section 4.2 that there is an offset in v_LSR of ∼0.2 km s⁻¹ between N₂H⁺ 4–3 and H₂D⁺ emission toward B2-MM8 and the nearby embedded protostar Elias 32, with the N₂H⁺ emission redshifted relative to the H₂D⁺ emission. Over the rest of Oph B2, the correlation between the v_LSR of the two lines is good, suggesting that the offset is real and localized. Comparison with N₂D⁺ emission at the same location shows that the v_LSR of the H₂D⁺ and N₂D⁺ agree well, and it is consequently the N₂H⁺ emission which is offset in velocity from the rest of the core gas. N₂H⁺ 4–3 has a critical density approximately an order of magnitude greater than the H₂D⁺ and N₂D⁺ lines discussed here. The offset in v_LSR may thus be probing a shift in the gas motions at extremely high densities in the core interior.

5.2. Trends in the Deuterium Fractionation

5.2.1. The N₂D⁺/N₂H⁺ Ratio

The deuterium fractionation, R_D, can be defined as the ratio between the column densities of a deuterated molecule and its hydrogen-bearing counterpart. Here, we define R_D = N(N₂D⁺)/N(N₂H⁺) and calculate R_D in Oph B2 from N(N₂D⁺) calculated in Section 4.4.1 and N(N₂H⁺) calculated in Paper II. Figure 7(b) shows R_D across B2. The distribution is similar to the N₂D⁺ column density in Figure 7(a). The largest values of R_D, to a maximum R_D = 0.16, are found toward northeast B2, while toward the smaller scale maxima, we find R_D ∼ 0.1 or less. We estimate that we are sensitive to R_D ∼ 0.01–0.015 based on the rms noise in the N₂D⁺ integrated intensity map (3σ ∼ 0.2 K km s⁻¹ in the off-B2 pixels), and the mean N₂H⁺ column density over Oph B2, 〈N(N₂H⁺)〉 = 5.9 × 10¹² cm⁻², calculated in Paper II. We find a moderate mean R_D = 0.03 with an rms variation of 0.01 toward B2, comparable to results in Ori B9 toward two N₂H⁺ emission peaks and a protostellar source (R_D = 0.03–0.04, Miettinen et al. 2009). While Pagani et al. (2009a, 2009b) find substantially greater central R_D values through modeling of the starless core L183 (R_D ∼ 0.7 at the core center) which decrease to ∼0.06 at ∼5000 AU, it is not clear how to compare their model of R_D as a function of core radius to our observed column density ratio. Fontani et al. (2006) found an average R_D ∼ 0.01 toward the high-mass star-forming region IRAS 05345+3157, which at high resolution resolved into two N₂D⁺ condensations each with R_D = 0.11 (Fontani et al. 2008). In a study of starless cores, Crapsi et al. (2005) find a range in R_D between a lower limit of 0.02 and a maximum of 0.44 (which was found toward the Oph D core). Crapsi et al. also note that R_D values >0.1 are generally only found for cores with N(N₂H⁺)>10¹³ cm⁻², which are approximately the maximum N(N₂H⁺) in Oph B2. A similar range in R_D was found toward protostellar, Class 0 sources (Emprechtinger et al. 2009). On average, then, our results agree with previous results. Given our map of Oph B2, however, we can also probe the variation of R_D across the core.

In Figure 12, we show R_D as a function of increasing N(H₂) (Figure 12(a)), and also as a function of increasing distance (in projection) from the nearest embedded protostar (Figure 12(b)). Each data point represents an 18'' pixel (to match the beam FWHM), and pixels with S/N <5 in either N₂H⁺ or N₂D⁺ emission are omitted. No clear trend is present in the deuterium fractionation with N(H₂). We find, however, that R_D increases at greater projected distances from the embedded sources, as might be expected from the N₂D⁺ integrated intensity distribution. We fit a linear trend to the data, where R_D = (−0.007 ± 0.002) + (0.37 ± 0.06)d [pc], but note that there is large scatter in the relation. An analysis of the individual N₂H⁺ and N₂D⁺ abundances shows that the increasing R_D trend with projected protostellar distance is dominated by a similar increase in X(N₂D⁺) with projected protostellar distance, while X(N₂H⁺) is approximately constant over the range plotted, shown in Figure 12(d).

5.2.2. H₂D⁺

We next look at the distribution of the ortho-H₂D⁺ abundance over Oph B2. In Figure 12, we show X(ortho-H₂D⁺) as a function of N(H₂) (Figure 12(e)) and also as a function of increasing projected protostellar distance (Figure 12(f)). We show data in 15'' pixels and omit pixels with S/N <5. We find no trend in the fractional ortho-H₂D⁺ abundance with H₂ column density (an increasing trend is found, but not shown, in the ortho-H₂D⁺ column density, N(ortho-H₂D⁺), with N(H₂)). The X(ortho-H₂D⁺) versus N(H₂) relationship is not entirely a scatter plot, as the ortho-H₂D⁺ abundance tends to increase with projected protostellar distance to a maximum X(ortho-H₂D⁺) = 4.9 × 10⁻¹⁰ at d = 0.035 pc, with X(ortho-H₂D⁺) ∼ 1.8 × 10⁻¹⁰ at larger distances. All data points at d > 0.035 pc correspond to the emission found in northeast B2, while the rest show emission closer to the central continuum emission peak. If we fit a linear relationship to the data, omitting the northeast B2 points, we find X(ortho-H₂D⁺) × 10⁻¹⁰ = (1.2 ± 0.5) + (79 ± 23)d [pc].

There is also no clear trend in N(ortho-H₂D⁺) with R_D where both H₂D⁺ and N₂D⁺ are strongly detected across Oph B2, although both species are found with relatively large column densities north and west of the continuum peak. Given the small-scale structure found in the N(N₂D⁺) distribution in B2, the spatial resolution and sampling of our H₂D⁺ observations may not be high enough to discern a correlation. We note, however, that Caselli et al. (2008) also found no clear correlation between R_D and N(ortho-H₂D⁺) in their sample of starless and protostellar cores (with the caveat that some R_D values were derived using column density ratios of deuterated NH₃ and deuterated H₂CO with their undeuterated counterparts rather than N(N₂D⁺)/N(N₂H⁺)). In particular, high values of R_D and relatively low N(ortho-H₂D⁺) were found toward two objects in Ophiuchus (Oph D and 16293E).

5.3. What is Affecting the Deuterium Fractionation in Oph B2?

The largest deuterium fractionation values in Oph B2 are not found toward the continuum emission peak but instead avoid the continuum peak and nearby protostars. The results of the previous section indicate that the distance to the nearest protostar is the dominant parameter impacting the deuterium fractionation in Oph B2. We next discuss several mechanisms by which the embedded sources are most likely to impact the deuterium fractionation distribution in B2, including temperature variation in the gas, liberation of CO from dust grains, the impact of protostellar outflows, and possible increase in the ionization fraction due to X-rays from the protostars.

5.3.1. Temperature and CO Depletion

Temperature can have a substantial impact on the deuterium fractionation in a dense core. As stated in Section 4.4.1, we expect the gas and dust to be well coupled at the densities of Oph B2 such that T_d = T_gas but specify in the following which of the gas or dust temperature is important in the following mechanisms. First, at T_gas ≳ 20 K, the H₂D⁺ reaction in Equation (1) can proceed both forward and backward, resulting in no net increase in H₂D⁺, and consequently other deuterated species formed via reactions with H₂D⁺. Second, the CO which was deposited (adsorbed) as an icy mantle onto dust grains while at low temperatures will evaporate (desorb) back into the gas phase if the dust grains are heated, providing a destruction mechanism for H₂D⁺ and again interrupting the deuteration chain.

The CO depletion factor, f_D(CO), is defined as the canonical CO abundance, X_can(CO), divided by the CO abundance, X(CO) = N(CO)/N(H₂), in the observed region:

$\begin{equation} f_D (\mbox{CO}) = \frac{X_{\rm can}(\mbox{CO})}{X(\mbox{CO})} \end{equation} \tag{ 8 }$

such that f_D values are low where observed CO abundances are high. Low CO depletion factors can dramatically lower the H₂D⁺ abundance at gas temperatures T < 15 K (Caselli et al. 2008). Currently, no estimates of f_D exist in B2 due to a lack of data at spatial resolutions matching the available continuum data. The CO observations of the JCMT GBLS should enable the first calculations of f_D at complementary spatial resolution and allow further analysis of the relationship between CO and the deuterium distribution in B2. The extended distribution of H₂D⁺ and N₂D⁺ across Oph B2 implies, however, that extensive depletion of CO has occurred. At densities of n = 10⁶ cm⁻³, CO adsorption and desorption occur equally frequently when T_d = 18 K (Visser et al. 2009). At T_d > 18 K, CO will thus again become significantly more abundant in the gas phase, while at T_gas > 20 K, enhanced deuterium fractionation via Equation (1) will cease to operate.

An obvious potential source of heating is the presence of embedded protostars within Oph B2. In a one-dimensional dynamic/chemical model of an isolated dense core which collapses to form a protostar, Aikawa et al. (2008, see their Figure 2) find that temperatures greater than 18 K can be found at radii ∼5000 AU (∼0.02 pc) from the protostar at ∼10⁵ yr after its formation. As a function of total luminosity of the core, L_core, however, their model predicts a CO sublimation radius of ≲1000 AU (0.005 pc) for L_core ∼ 0.5 − 1 L_☉. In a more specific example, Stamatellos et al. (2007) model the heating effect of an embedded Class I protostar (L_bol = 10 L_☉) on the surrounding gas in the Oph A core and find only a small temperature increase (≲5 K), limited to the immediate vicinity of the protostar, in agreement with Aikawa et al. (2008).

Both YSOs in central B2 have been classified as embedded (Jørgensen et al. 2008), Class I sources (Enoch et al. 2007). The estimated bolometric luminosities, L_bol ∼ 0.5 L_☉ and ∼1 L_☉ for E32 and E33, respectively (Enoch et al. 2007; van Kempen et al. 2009), are much smaller than used in the heating model by Stamatellos et al. (2007). The bolometric temperatures, T_bol ∼ 300 K and ∼500 K, suggest that both are in the later stages of the Class I phase, as do their negative infrared spectral indices (α_IR = −0.03 and −0.12). Sources with −0.3 < α_IR < 0 are sometimes classified as 'flat-spectrum' sources and are thought to be an intermediate stage between Class I and Class II. Crapsi et al. (2008) suggest that most or all flat-spectrum sources could be Class II objects seen edge-on. Based on these results, E32 and E33 are not likely to be deeply embedded within Oph B and have ages ∼5 × 10⁵ yr or greater based on lifetime estimates of the different protostellar stages (Evans et al. 2009).

The low luminosities of E32 and E33 suggest that they should not have a large impact on the temperature of B2, in agreement with our temperature results in Paper I which showed very little variation in the gas temperature over Oph B. Figure 12 shows that we find no trend in the kinetic temperature T_K, determined using NH₃ inversion transitions in Paper I, with protostellar distance above the ∼1–2 K uncertainties in the data. The derived temperatures are near the required values, however, and the NH₃ temperatures may not be indicative of the conditions of the denser gas where we expect most of the deuterated species to exist. Depletion of H₂D⁺ and N₂D⁺ via CO sublimation would thus require an increasing temperature gradient toward higher densities near the protostars, which is not unreasonable.

Since CO also destroys N₂H⁺, significant evaporation of CO from dust grains should also be observable in the N₂H⁺ distribution. Accordingly, in Figure 1(b) we see a decrease in the integrated N₂H⁺ 1–0 intensity toward E33. No good fit to the N₂H⁺ spectra was found at this location, however, due to complicated line structures, such that we are not able to determine the N₂H⁺ column density or abundance. Figure 12(d) shows no trend on larger scales of X(N₂H⁺) with protostellar proximity. The N₂H⁺ 4–3 distribution shows that N₂H⁺ remains in the gas toward E33, but an analysis of the excitation of the line would be needed to determine the N₂H⁺ column density from these data.

In summary, heating by an embedded protostar can produce slightly greater gas temperatures and evaporate CO from the surface of dust grains, which would decrease the deuterium fractionation in the gas. Since we do not see evidence of a decrease in X(N₂H⁺) toward the protostars, the evaporation of CO from dust grains is likely not the major mechanism for the decrease in abundance of deuterated species. While we also do not see evidence of increased gas temperatures near the protostars in the NH₃ data, these temperatures reflect average values along the line of sight, and do not rule out increased temperatures at higher densities near the protostars. Given the strong dependence of increased deuterium fractionation on temperature, we expect this mechanism is likely important in B2.

5.3.2. The Oph B2 Outflow

The deuterium fractionation may be impacted by the known, extensive CO outflow in B2 described in Section 5.1. It appears that the dense gas in eastern Oph B2 is affected by the outflow, given the larger line widths seen toward the protostars in high density tracers, changes in v_LSR, and the blue shoulder emission which is coincident with the eastern lobe of the outflow (Kamazaki et al. 2003; although large line widths could also indicate the presence of infall motions). We do not find evidence of outflow interaction with the dense gas in northwest Oph B2, suggesting that the outflow may be in front of or behind the bulk of the core gas. Very localized heating of the surrounding gas by the protostar could evaporate CO from the dust grains. The outflow, then containing entrained CO, could mix the CO with the dense core gas, destroying the deuterated species without an increase in the gas temperature on larger scales. Alternatively, shocks associated with the outflow could raise gas temperatures and reintroduce enough CO to the gas phase to affect the deuterium fractionation. Further analysis of the outflow, including observations of shock tracers, is needed to determine whether this scenario is reasonable and can examine the impact of low-mass protostellar outflows on dense gas chemistry. As discussed above, however, we do not find evidence for variation in the N₂H⁺ abundance as would be expected if destruction due to gas-phase CO is the dominant mechanism for decreasing the abundance of deuterated species relative to their non-deuterated counterparts.

5.3.3. Ionization Fraction and Electron Abundance x_e

Variations in the ionization fraction in the dense gas can also affect the deuterium fractionation, in particular in regions of high CO depletion. The electron abundance, $x_e = n_e/n_{\rm H_2 }$ , is the ratio of the electron number density (n_e) with the number density of molecular hydrogen, $n_{\rm H_2}$ , and is equivalent to the ionization fraction. The deuterium fraction can be reduced by an increase in x_e, which in turn increases the rate of dissociative recombination (Caselli et al. 2008). The abundances of deuterated molecules can thus be used as a probe of x_e in dense cores. At high extinction toward starless cores, the creation of ions is expected to be dominated by the cosmic ray flux, with $x_e \propto n_{\rm H_2}^{-1/2}$ (McKee 1989).

Studies investigating x_e in dense cores using the [DCO⁺]/[HCO⁺] ratio have found average line-of-sight x_e ∼ 10⁻⁷ (Bergin & Tafalla 2007, and references therein). In regions of high depletion, the column density ratio of N₂D⁺ to N₂H⁺ is expected to be a better probe of x_e than [DCO⁺]/[HCO⁺]. Since H⁺₃ is expected to be the major molecular ion in cores strongly depleted in heavy elements, we can use observations of H₂D⁺, in conjunction with R_D, to set limits on x_e assuming that the numbers of positive and negative charges are approximately equal (i.e., the core gas is neutral). Given the abundances of N₂H⁺, N₂D⁺, and H₂D⁺ in B2, we can thus derive a lower limit on x_e following Miettinen et al. (2009):

$\begin{equation} x_e > X(\mbox{N$_2$H$^+$}) + X(\mbox{N$_2$D$^+$}) + X(\mbox{H$_3^+$}) + X(\mbox{H$_2$D$^+$}), \end{equation} \tag{ 9 }$

where we neglect multiply deuterated forms of H⁺₃. As discussed in Section 4.4.2, predictions of the o/p-H₂D⁺ ratio vary from ∼0.1 to ≳1, and we use a moderate o/p-H₂D⁺ = 0.5 to estimate X(H₂D⁺) from X(ortho-H₂D⁺). In a simple, steady state analytical model, Crapsi et al. (2004) relate R_D (where R_D = N(N₂D⁺) / N(N₂H⁺)) to the relative H₂D⁺–H⁺₃ abundance, r = [H₂D⁺]/[H⁺₃], and also neglect multiply deuterated forms of H⁺₃ such that R_D ≈ r/(3 + 2r). In an updated result, Miettinen et al. (2009) find R_D ≈ (r + 2r²)/(3 + 2r + r²). Using the Miettinen et al. relation, we then find a mean r = 0.07 across Oph B2 and a lower limit x_e > 1.5 × 10⁻⁸ through Equation (9). We find small variations in the lower limit across B2, with a larger limit x_e > 2.7 × 10⁻⁸ toward Elias 32 and higher values in general near the 850 μm continuum peak. Detailed modeling, in conjunction with a measurement of the CO depletion factor in B2, is needed to improve this analysis of x_e beyond a lower limit only.

Low mass protostars emit X-rays, which could potentially increase the local ionization fraction and consequently decrease locally the deuterium fractionation of gas in Oph B2. In their study of the X-ray emission from young protostars in the central Ophiuchus region with ROSAT, Casanova et al. (1995) detected both Oph B2 protostars in the 1–2.4 keV range, finding luminosities (adjusted for extinction) L_x = 30.7 erg s⁻¹ and L_x = 29.2 erg s⁻¹ for E32 and E33, respectively. We can then estimate the X-ray ionization rate due to the protostars, ζ_x, at a distance r following Silk & Norman (1983, Equation (2)). From Figure 12, we find trends in R_D and X(ortho-H₂D⁺) to projected protostellar distances of order of a few ×0.01 pc. Assuming a constant, uniform density n = 10⁵ cm⁻³, at r = 0.02 pc we find ζ_x ∼ 9 × 10⁻¹⁷ s⁻¹, and at r = 0.04 pc we find ζ_x ∼ 8 × 10⁻¹⁸ s⁻¹. The fractional ionization based on these rates, again following Silk & Norman, is x_e ∼ 9 × 10⁻⁸ and x_e ∼ 3 × 10⁻⁸ at r = 0.02 pc and r = 0.04 pc, consistent with our lower limits calculated above. Clumpy, rather than uniform gas could allow X-rays to penetrate further into the core, potentially increasing the fractional ionization to greater distances. The standard ionization rate due to cosmic rays alone ζ_CR = 1.3 × 10⁻¹⁷ s⁻¹ (Herbst & Klemperer 1973), leading to a fractional ionization x_e ∼ 4 × 10⁻⁸ for n = 10⁵ cm⁻³ (McKee 1989). Our simple calculation shows that at distances comparable to where we find variation in the deuterium fractionation in Oph B2, the protostellar X-ray ionization rate and resulting fractional ionization may be comparable to or greater than expected due to cosmic rays alone. This result suggests that X-ray ionization is a viable mechanism to depress the deuterium fractionation close to embedded protostars. Further work is needed to directly measure the fractional ionization of dense gas near protostars and to quantify the effect of a greater ionization rate on the chemical network of deuterated species in high density gas.

5.3.4. Other Possibilities

While Figure 12 suggests that the distribution of enhanced deuterium fractionation of H⁺₃ and N₂H⁺ in Oph B2 is tied to the presence of the Class I protostars associated with the core, other factors may influence the deuteration chemistry in this complex environment. For completeness, we briefly mention two other possibilities. First, line-of-sight inhomogeneities in Oph B2 would impact our estimate of the deuterium fractionation, as our analysis assumes constant excitation and abundance conditions along the line of sight. Detailed physical and chemical modeling may be useful in probing the effects of possible inhomogeneities but is beyond the scope of this paper due to the complexity of Oph B2. Second, differences in the evolutionary timescale may be important in Oph B2 and in other multiple-star-forming objects. While several protostars have already formed out of the core gas (particularly in the south), there remain multiple dense, starless clumps within B2 which may yet form stars. The timescale on which different parts of the core condense and collapse may impact the dense gas chemistry and produce variations in the deuterium fractionation across Oph B2.

In summary, there a number of different mechanisms by which the deuterium fractionation of N₂H⁺ and H⁺₃ may be decreased near young protostars. We feel that a local increase in the gas temperature is the most likely explanation of this trend, but our simple calculation of the expected increase in the ionization fraction above the standard cosmic ray value is intriguing and worth further study.

6. SUMMARY

In this paper, we have presented the results of N₂D⁺ 3–2, N₂H⁺ 4–3, and ortho-H₂D⁺ 1₁₁–1₁₀ mapping of the cluster-forming Ophiuchus B core. In general, our results show that care must be taken while using deuterated species as a probe of the physical conditions of dense gas conditions in star-forming regions. We summarize in detail our results below.

1.
Significant emission from the deuterated species N₂D⁺ and H₂D⁺ extends over several square arcminutes in B2. In particular, the extent of the H₂D⁺ emission is the largest yet mapped. Both species' integrated intensity distributions are offset from the continuum emission and avoid entirely parts of B2 near embedded protostars and the peak of continuum emission from cold dust.
2.
Through fitting of the spectral lines, we find greater N₂D⁺ line widths near protostars embedded in B2, suggestive of protostellar influence on the dense gas. Complicated line structures are found, including blue line shoulders, in all three tracers. We estimate the non-thermal line widths for the species observed here and compare with NH₃ and N₂H⁺ results from Paper I and Paper II. In B2, non-thermal line widths decrease with increasing critical density of the molecular tracer but remain transonic or only slightly subsonic even at high densities. This result is in contrast with the nearly thermal line widths found in many isolated cores.
3.
We test a method to estimate the gas temperature by equating the non-thermal line widths of N₂D⁺ and H₂D⁺ and solving for T_K. We find both realistic (T_K ∼ 11 K) and unrealistic (T_K < 0 K and T_K ≳ 25 K) values in B2. Complicated line structure in B2 likely affects this analysis, but we expect good results can be found in regions with simple velocity structure.
4.
We find R_D = N(N₂D⁺)/N(N₂H⁺) ∼ 0.03 on average where N₂D⁺ is detected. Small-scale features of enhanced R_D are seen in B2, but most of these do not correlate with previously identified continuum clumps or N₂H⁺ clumps. This R_D is consistent with previous results in protostellar cores and is also within the lower range of values found for starless cores.
5.
The average ortho-H₂D⁺ abundance, X(ortho-H₂D⁺) ∼ 3 × 10⁻¹⁰, given an assumed excitation temperature of 12 K. An accurate measure of the ortho- to para-H₂D⁺ ratio is needed to convert to total H₂D⁺ abundance and has been found in previous studies to vary between 0.1 and >1. Future observations of para-H₂D⁺ transitions, if possible, are needed to constrain this value in dense cores.
6.
An anti-correlation is found with R_D and X(ortho-H₂D⁺) and projected distance to the nearest embedded protostar in B2. This relation indicates that the embedded protostars in B2 are affecting the internal chemistry. Possible mechanisms are through gas heating, evaporation of CO off dust grains, CO mixing caused by the B2 outflow, or by increasing the local ionization fraction. Gas temperatures determined through NH₃ analysis in Paper I are <18–20 K required to evaporate CO and reduce deuterium fractionation, but NH₃ may not be a good tracer of conditions in the high density gas in Oph B2. We calculate that the X-ray flux of the embedded protostars in B2 could increase the fractional ionization of the core gas beyond that expected solely due to cosmic ray ionization to distances ∼0.02–0.03 pc, similar to the distances over which we find variation in the deuterium fractionation. We suggest that a temperature increase is the dominant parameter in this trend, but further research on the effect of young stars and their outflows on both the local ionization fraction and the chemistry of deuterated species in protostellar cores is needed.
7.
We use the observed R_D and X(ortho-H₂D⁺) to estimate a lower limit on the ionization fraction in B2 and find x(e)>1.5 × 10⁻⁸ given an assumed ortho- to para-H₂D⁺ ratio of 0.5.

We thank the anonymous referee and E. Feigelson for comments, which improved the paper. We also thank H. Kirk for providing SCUBA maps of the observed regions and H. Weisemeyer for help with the IRAM observations. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. The James Clerk Maxwell Telescope is operated by the Joint Astronomy Centre on behalf of the Particle Physics and Astronomy Research Council of the United Kingdom, the Netherlands Association for Scientific Research, and the National Research Council of Canada. R.K.F. acknowledges financial support from the University of Victoria and the National Research Council Canada Graduate Student Scholarship Supplement Program. We also acknowledge the support of the National Science and Engineering Research Council of Canada. This research is supported in part by the National Science Foundation under grant number 0708158 (T.L.B.).

Facilities: JCMT - James Clerk Maxwell Telescope, IRAM:30 m - Institute de Radioastronomie Millimetrique 30 meter telescope

THE INITIAL CONDITIONS OF CLUSTERED STAR FORMATION. III. THE DEUTERIUM FRACTIONATION OF THE OPHIUCHUS B2 CORE

Article metrics

Permissions

Author e-mails

Author affiliations

Dates

ABSTRACT

1. INTRODUCTION