EDP Sciences
Free Access
Issue
A&A
Volume 494, Number 2, February I 2009
Page(s) 623 - 636
Section Interstellar and circumstellar matter
DOI https://doi.org/10.1051/0004-6361:200810587
Published online 20 November 2008

Chemical modeling of L183 (L134N): an estimate of the ortho/para H2 ratio[*]

L. Pagani1 - C. Vastel2 - E. Hugo3 - V. Kokoouline4 - C. H. Greene5 - A. Bacmann6 - E. Bayet7 - C. Ceccarelli6 - R. Peng8 - S. Schlemmer3


1 - LERMA & UMR8112 du CNRS, Observatoire de Paris, 61 Av. de l'Observatoire, 75014 Paris, France
2 - CESR, 9 avenue du colonel Roche, BP 44348, Toulouse Cedex 4, France
3 - I. Physikalisches Institut, Universität zu Köln, Zülpicher Strasse 77, 50937 Köln, Germany
4 - Department of Physics, University of Central Florida, Orlando, 32816, USA
5 - Department of Physics and JILA, University of Colorado, Boulder, Colorado 80309-0440, USA
6 - Laboratoire d'Astrophysique de Grenoble, Université Joseph Fourier, UMR 5571 du CNRS, 414 rue de la Piscine, 38041 Grenoble Cedex 09, France
7 - Department of physics and astronomy, University College London, Gower street, London WC1E 6BT, UK
8 - Caltech Submillimeter Observatory, 111 Nowelo Street, Hilo, HI 96720, USA

Received 14 July 2008 / Accepted 6 October 2008

Abstract
Context. The high degree of deuteration observed in some prestellar cores depends on the ortho-to-para H2 ratio through the H3+ fractionation.
Aims. We want to constrain the ortho/para H2 ratio across the L183 prestellar core. This is required to correctly describe the deuteration amplification phenomenon in depleted cores such as L183 and to relate the total (ortho+para) H2D+ abundance to the sole ortho-H2D+ column density measurement.
Methods. To constrain this ortho/para H2 ratio and derive its profile, we make use of the N2D+/N2H+ ratio and of the ortho-H2D+ observations performed across the prestellar core. We use two simple chemical models limited to an almost totally depleted core description. New dissociative recombination and trihydrogen cation-dihydrogen reaction rates (including all isotopologues) are presented in this paper and included in our models.
Results. We estimate the H2D+ ortho/para ratio in the L183 cloud, and constrain the H2 ortho/para ratio: we show that it varies across the prestellar core by at least an order of magnitude, being still very high ($\approx$0.1) in most of the cloud. Our time-dependent model indicates that the prestellar core is presumably older than 1.5-2 $\times $ 105 years but that it may not be much older. We also show that it has reached its present density only recently and that its contraction from a uniform density cloud can be constrained.
Conclusions. A proper understanding of deuteration chemistry cannot be attained without taking into account the whole ortho/para family of molecular hydrogen and trihydrogen cation isotopologues as their relations are of utmost importance in the global scheme. Tracing the ortho/para H2 ratio should also place useful constraints on the dynamical evolution of prestellar cores.

Key words: ISM: abundances - ISM: clouds - ISM: structure - astrochemistry - ISM: individual objects: L183 - molecular processes

1 Introduction

Studies of the earliest stages of star formation have increased in the last ten years with the advent of new receivers acquiring better spatial and spectral resolution. Prestellar cores are dense and cold cores where gravitational collapse has not yet occured. In the densest regions of the core (where $n_{\rm H_2}$ is greater than a few 104 cm-3) most heavy molecules containing carbon, nitrogen and oxygen seem to deplete onto the dust grains and only light ions remain in the gas phase. There has been extensive observational evidence of CO and CS depletion in the center of prestellar cores (e.g. Tafalla et al. 2004; Pagani et al. 2005; Bacmann et al. 2002; Caselli et al. 1999; Bergin et al. 2002) which seems to be typical of the majority of dense cores. Nitrogen-bearing species like CN, NH3 and N2H+ appear to subsist longer before freezing-out onto the dust grains (Hily-Blant et al. 2008; Pagani et al. 2007; Tafalla et al. 2006).

In conditions under which heavy species are depleted, H+ and H3+ (and its deuterated counterparts) are the most abundant ions subsisting in the gas phase. H2D+ has been widely detected and mapped in protostars and prestellar cores (Caselli et al. 2003; Vastel et al. 2006; Caselli et al. 2008) through its ortho fundamental line, with abundances high enough to be explained by the CO depletion theory/observations. Although difficult to observe from Earth, the D2H+ molecule has been detected with its para line towards the 16293E prestellar core in the L1689N molecular cloud (Vastel et al. 2004) with an abundance similar to the ortho-H2D+ molecule as suggested by Phillips & Vastel (2003).

Consequently many theoretical studies were performed that included all the deuterated forms of the H3+ ion (e.g. Roberts et al. 2003; Walmsley et al. 2004). However all nuclear spin states (ortho, meta, para, corresponding to the spin state of the protons or deuterons) must be taken into account in order to compare with the observational sets. Moreover the thermicity of the forward/backward reactions strongly depends on the symmetric state of the species. Though the influence of the ortho/para H2 ratio on the chemistry of H2D+ (Pagani et al. 1992) had been described several years before the first detection of the ion (Stark et al. 1999), it is only recently that this specific spin-state dependent chemistry has been studied in detail (Flower et al. 2004,2006a).

The motivation for our study is the many deuterated observations performed in the L183 prestellar core (PSC) and the main aim of this paper is to study the ortho/para H2 ratio from some of these deuterated species observations. To this effect, we constrain two chemical models including all the symmetric states of H2, D2, H3+ and its deuterated counterparts with observations of ortho-H2D+, combined with previous N2H+ and N2D+ observations. These models have been set up using recent dissociative recombination rates computed for H3+ and its isotopologues as well as all non negligible reaction rates between H2 and H3+ and their isotopologues (both presented for the first time in this paper).

2 Observations

2.1 Deuterated H3+

We first observed the ortho-H2D+ (and para-D2H+) with the Caltech Submillimeter Observatory (CSO) monopixel receiver and subsequently took advantage of the newly installed Heterodyne ARray Program, 16 channel 350 GHz band (``B-band'') (HARP-B) camera on the James Clerk Maxwell Telescope (JCMT) to fully map its emission.

2.1.1 CSO observations

The H2D+ and D2H+ observations were carried out at the CSO, between August 2004 and April 2005, under good weather conditions (225 GHz zenith opacity always less than 0.06). Scans were taken, using the chopping secondary with a throw of 3$^{\prime}$, using the reference position: $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$.

 \begin{figure}
\par\includegraphics[width=6cm,clip]{0587fig1.eps}
\end{figure} Figure 1:

CSO map of the H2D+ (110-111) line. The position is indicated in arcseconds for each spectrum and the (0, 0) position corresponds to $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. The Y-axis represents the antenna temperature.

Open with DEXTER

 \begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{0587fig2.eps}
\end{figure} Figure 2:

JCMT map of the H2D+ (110-111) line. The dust peak position (Pagani et al. 2003) is indicated by a cross and corresponds to $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. Contour levels are 0.1, 0.2 and 0.3 K km s-1.

Open with DEXTER

Table 1:   Line parameters from the JCMT and CSO observations. The positions are offsets to the dust peak emission at $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. For non-detected lines we give the 3$\sigma $ upper limit. For JCMT and p-D2H+ at CSO, we give both the Monte Carlo (MC) and the LTE column density estimates.

The 345 GHz (respectively 650 GHz) sidecab receiver with a 50 MHz acousto-optical spectrometer backend was used for all observations with an average velocity resolution of 0.04 km s-1 (respectively 0.02 km s-1). At the observed frequencies of 372.421385(10) GHz for the H2D+ (110-111) and 691.660483(20) for the D2H+ (110-101) lines (Amano et al. 2005), the CSO 10.4 m antenna has a HPBW of about 20 $^{\prime \prime }$ and 11 $^{\prime \prime }$ respectively. We performed a cut in declination across the main dust peak in L183 between (0, -40 $^{\prime \prime }$) and (0, 130 $^{\prime \prime }$) for H2D+ and concentrated on the reference position for D2H+. The system temperature was typically 1000 to 2000 K for H2D+ and 2500 to 3500 K for D2H+.

The beam efficiency at 372 GHz (respectively 692 GHz) was measured on Saturn and Jupiter and found to be $\sim$60% (respectively $\sim$40%) for point sources. Pointing was monitored every 1.5 h and found to be better than 3 $^{\prime \prime }$. If the emission is extended compared to the beam size of CSO, as appears to be the case for H2D+ then the efficiency is about 70% at 372 GHz and 60% at 692 GHz. The data reduction was performed using the CLASS program of the GILDAS software developed at IRAM and the Observatoire de Grenoble and the LTE data analysis with CASSIS developed at CESR (http://www.cassis.cesr.fr).

These CSO observations confirmed that the H2D+ emission was extended in this source as we could hypothesize from several tracers and we carried on the project at the JCMT, recently equipped with a new 16-pixel camera to be able to fully map the source in a reasonable amount of time.

2.1.2 JCMT observations

The JCMT observations were obtained during semester 07A in service mode, using the HARP-B 16 pixel camera (one pixel, located in a corner, was unavailable). A third of the observations was obtained in jiggle-chop[*] mode and two thirds in position switch (PSw) mode. The jiggle-chop mode appeared to be no faster, the displacement of the telescope in PSw mode seeming minor compared to other overheads, and because the jiggle-chop mode works in the Nyquist regime, each pixel receives much less integration time than in PSw mode. As adjacent pixels had the same off spectrum subtracted, the spatial average did not bring much improvement and we subsequently changed to position switch mode because deep integration on weak signal appears more important than Nyquist sampling for this work. In PSw mode, we made 2 $\times $ 2 pointings to fill the gaps in the camera, thus achieving a full beam sampling. Two such sets were performed to cover the main dust peak and its northern extension (Fig. 2) with one pixel row overlap between the two.

Most of the observations were run in band 1 weather ( $\tau_{225~{\rm GHz}} < 0.05$) while a few were done in band 2 weather leading to rapid degradation of the system temperature. The source was observed only above 40$^\circ$ elevation and the band 1 weather system temperature was in the range 500-1000 K depending on the pixels and on the elevation. To observe both H2D+ and N2H+ (J:4-3) (at 372.672509 GHz), we tuned the receiver half way between the two lines and used a frequency resolution of 61 kHz so that the backend could cover both lines at once.

Data pre-reduction was done with the Starlink software (KAPPA, SMURF and STLCONVERT packages) and subsequently translated into CLASS format for final reduction.


2.2 CO depletion and dust content


All other observations used in this paper have been obtained and published previously. The dust content of L183 both in emission and in absorption has been described in Pagani et al. (2004,2003). The source size is half a degree and contains a long filament extending $\sim$6$^{\prime}$ from north to south. Two peaks are clearly visible, one just south of the middle of the filament (which we call the main peak) with an opacity of $\sim$150 $A_{\rm V}$ and a second one, 3$^{\prime}$ north of the first one (the north peak) with an estimated opacity of $\sim$70 $A_{\rm V}$. These peaks have the characteristics of prestellar cores. Most of the filament have an opacity above 40 $A_{\rm V}$.

Two large scale C18O and C17O maps obtained with the Kitt Peak 12-m telescope fail to trace the dense filament (Pagani et al. 2005). It is now well established that this is due to depletion of CO onto grains. Surprisingly, the opacity at which the depletion begins ($\sim$20 $A_{\rm V}$) is two times higher than what is usually observed in other clouds (e.g. Alves et al. 1999; Bergin et al. 2002; Kramer et al. 1999) though it still appears at a density ($\sim$$\times $ 104 cm-3) that is a typical depletion density threshold (Pagani et al. 2005). Possibly, the low density envelope where depletion has not yet occurred is very extended in this cloud (which is confirmed by its large influence on the C18O (J:2-1) line intensity (Pagani et al. 2002). The depletion factor in volume for CO has been estimated to be 43 on average (Pagani et al. 2005) on the line of sight of the main dust peak and is probably much higher in the inner part of this core where density is above 1 $\times $ 105 cm-3.

2.3 N2H+ and N2D+

N2H+ and N2D+ have been mapped at both low (Kitt Peak 12-m, Pagani et al. 2005) and high (IRAM 30-m) resolutions. From the high resolution data, a strip crossing the main dust peak has been extracted and published (Pagani et al. 2007). In that paper, we performed a detailed analysis of the N2H+ and N2D+ emission with a Monte Carlo model treating exactly the hyperfine structure and line overlap of these species. We derived several physical properties, namely a maximum density of 2 $\times $ 106 cm-3, with a radial dependence proportional to r-1 up to 4000 AU and proportional to r-2 beyond, a kinetic temperature of 7 ($\pm$1) K, a slight depletion of N2H+ in the inner 3000 AU of the core and a deuterium fractionation which is non-measurable at 104 AU (<0.03) and reaches $\sim$0.7 ($\pm$0.15) in the center. As far as we know this is the highest fractionation reported yet for a singly deuterated species. However this may not be exceptional when compared to the detection of triply deuterated species, like ND3 (van der Tak et al. 2002; Lis et al. 2002) and CD3OH (Parise et al. 2004) or to the fact that most reported fractionations are line-of-sight averages and are not derived from detailed profiles.

3 Analysis

We present 3 models in this paper: a Monte Carlo model to compute H2D+ and D2H+ line intensities, a chemical steady-state model and a chemical time-dependent model.

Compared to previous works, we benefit here from two new sets of coefficients and a large set of observations in a PSC. The H3+ + H2 (+isotopologues) set of rate coefficients are extracted from the PhD work of Hugo in advance of publication and the H3+ (and deuterated counterparts) recombination rates have been computed for this work by Kokoouline and Greene and are presented in Appendix B. Rate coefficients, as computed by Hugo, describe all possible interactions between H3+ and H2 isotopologues, including reactive and non-reactive, elastic and inelastic collisional rates, while recombination coefficients describe the dissociative recombination (DR) rates of trihydrogen cation isotopologues. Both sets take into account all the ortho, para and (D3+) meta states.

3.1 Line emission

We have analysed the line emission of both ortho-H2D+ and para-D2H+ using a two-level Monte Carlo code (adapted from Bernes 1979) with our new collisional coefficients. Because the temperature of the PSC is around 7 K for both the gas (Pagani et al. 2007) and the dust (Pagani et al. 2004,2003), the possibility of populating the next rotational level, $J_{{\rm KK}''}$ = 212 at 114 K and 75 K above ground level for ortho-H2D+ and para-D2H+ respectively, is so low that it can be safely ignored. With the new coefficients, the critical densities are 1.1 $\times $ 105 and 4.9 $\times $ 105 cm-3 respectively and thus the lines are close to thermal equilibrium in the inner core. In the case of ortho-H2D+, the difference between LTE ($\sim$2.0 $\pm$ 0.25 $\times $ 1013 cm-2) and Monte Carlo ($\sim$2.3 $\pm$ 0.25 $\times $ 1013 cm-2) column density estimates is typically 10-15% in the direction of the dust peak.

The para-D2H+ line has not been detected (see Caselli et al. 2008) and the 3$\sigma $ upper limit corresponds to a total column density of $\sim$2.4 $\times $ 1013 cm-2 using the Monte Carlo code. The LTE estimate yields $\sim$1.5 $\times $ 1013 cm-2.

 \begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{0587fig3.eps}
\vspace*{5mm}
\end{figure} Figure 3:

Ortho-H2D+ spectra across the main dust peak. East and west sides are folded together and fitted with a Monte Carlo model. Density and temperature profiles are taken from Pagani et al. (2007). The spacing between spectra is 15 $^{\prime \prime }$.

Open with DEXTER

3.2 Deuteration

H3+ ions are formed at a rate 0.96 $\zeta $ by cosmic ray ionization of H2 (Walmsley et al. 2004), rapidly followed by reaction with another H2 to form H3+, and are destroyed in reactions with neutral species and in dissociative recombination with free electrons, negatively charged grains and possibly negatively charged polycyclic aromatic hydrocarbons (PAHs-). In prestellar cores, the primary reservoirs of hydrogen and deuterium are H2 and HD, with HD/H2 = 2(D/H) $_{\rm cosmic}$ $\sim$ 3.2 $\times $ 10-5 (Linsky 2007). The proton exchanging reaction of H3+ with H2 partly regulates the H2 ortho-to-para ratio but has no effect on the H3+ abundance. Concurrently, the reaction with HD forms H2D+ and this primary fractionation is then followed by the subsequent fractionations and produces D2H+ and D3+ (Phillips & Vastel 2003; Roberts et al. 2003):


    $\displaystyle {\rm H_3^+ + HD} \longleftrightarrow {\rm H_2D^+ + H_2 + 232~K}$ (1)
    $\displaystyle {\rm H_2D^+ + HD} \longleftrightarrow {\rm D_2H^+ + H_2 + 187~K}$ (2)
    $\displaystyle {\rm D_2H^+ + HD} \longleftrightarrow {\rm D_3^+ + H_2 + 234~K.}$ (3)

The backward reactions are endothermic with an energy barrier of about 200 K (when considering only the ground level for each species) and were thought to be negligible at the low temperatures found in prestellar cores ($\leq$20 K), in which case the abundance ratios $\frac{[{\rm H}_n{\rm D}_{3-n}^+]}{[{\rm H}_{n+1}{\rm D}_{2-n}^+]}_{n=0,1,2}$ would be greatly enhanced. However, such enhancement can be limited by various processes (see Fig. 4):
  • dissociative recombination of H3+ (and its deuterated counterparts) with free electrons or negatively charged grains (and PAHs-?);

  • reactions with ``proton-friendly'' molecules such as CO and N2 which destroy the trihydrogen cations to produce HCO+ and N2H+;

  • ortho-H2 which can react with the deuterated trihydrogen cation and remove the deuterium (see below).

 \begin{figure}
\par\includegraphics[width=8.5cm,clip]{0587fig4.eps}
\vspace*{5mm}
\end{figure} Figure 4:

Main reactions involved in the H3+ chemical network. When CO and N2 are depleted, the reactions with bold arrows are dominant.

Open with DEXTER

In this modeling, we introduce the backward reactions to Eqs. (1)-(3) as we distinguish between ortho, meta and para states of the different species. When these reactions are completely neglected, the deuteration fractionation is considerably enhanced and observations towards pre-stellar cores cannot be reproduced (Roberts et al. 2003). Indeed, the importance of considering ortho and para states of various H/D carriers in the chemistry of the trihydrogen cation and isotopologues was first discussed by Pagani et al. (1992) and subsequently expanded in a series of papers by Flower and coworkers (Flower et al. 2006b,2004,2006a; Walmsley et al. 2004, hereafter collectively referred to as FPdFW). Not only is this important in comparing the chemical model predictions on the abundance of H2D+ to the observations of the ortho-H2D+ species alone but also because some important reactions are much faster with ortho-H2 than with para-H2 hence no longer negligible.. Indeed, the 170 K internal energy of the lowest ortho-H2 level (J=1) is large with respect to the temperatures of interest and can significantly enhance the Boltzmann factor of endothermic reactions. In some cases, reactions which are endothermic with para-H2 can turn out to be exothermic with ortho-H2 i.e. fast and temperature independent. In fact, the ortho-to-para ratio of H2 is found to be a crucial parameter for the entire chemistry of deuterium.

The key reactions involving ortho-H2 are essentially with meta-D3+, para-D2H+ and ortho-H2D+ as well as para-H2D+ (see rates in Appendix A) because the internal energy of the ortho-H2 alone is not enough to overcome the endothermicity of reactions 1 to 3. Thus only those species which also have an internal energy high enough (so that the sum of the two internal energies is higher than the endothermicity of reactions 1 to 3) can react with ortho-H2 at the Langevin rate in cold gas. Thus, the three reactions involving ortho-H2 with meta-D3+, para-D2H+ and ortho-H2D+ present exothermic or thermoneutral channels to rehydrogenate the cations forming ortho-D2H+, para-H2D+ and ortho/para-H3+ while the reaction between ortho-H2 and para-H2D+ can efficiently convert the latter to ortho-H2D+. The ortho-H2 molecule thus opens a path to climb the 4 step energy ladder back from para-D2H+ to H3+ via para-H2D+ and ortho-H2D+ which can be very efficient in the presence of large ortho-H2 fractions. However, this efficient ladder scheme does not include D3+ because the conversion of ortho-D2H+ to para-D2H+ is strictly forbidden in collisions with H2 and very inefficient in collisons with HD. Nevertheless, these reactions can be a strong limit to the isotopic fractionation of H3+ hence of other species. Any chemical model which includes deuterium chemistry must distinguish between ortho and para states of dihydrogen and trihydrogen cation isotopologues and include reactions between the different spin states following Flower et al. (2006b); Pagani et al. (1992); Flower et al. (2006a) and the present work.

3.3 CO and N2 chemistry

The CO and N2 abundances are critical parameters in the deuteration of the H3+ ion. CO is expected to freeze-out onto the grain mantles at high densities (a few 104 cm-3) and low temperatures ($\leq$20 K) (e.g. Tafalla et al. 2002; Pagani et al. 2005; Bacmann et al. 2002; Bergin et al. 2001; Caselli et al. 1999). With an N2 binding energy similar to the CO binding energy (Bisschop et al. 2006; Öberg et al. 2005), these two molecules are expected to behave similarly. However, observations towards prestellar cores prove that N2H+ (produced from N2) remains in the gas phase at higher densities than CO. This can be explained by the fact that N2H+ is mainly destroyed by CO (Aikawa et al. 2005; Pagani et al. 2005; Caselli 2002), so that the CO freeze-out implies a drop in the destruction rate of N2H+. This would partially balance the lower formation rate due to the N2 freeze-out. Consequently, N2H+ is observed to survive in the gas phase at higher densities ($\sim$106 cm-3). In the case of L183, we have shown that N2H+ partially survives but suffers from depletion at densities starting at $\sim$$\times $ 105 cm-3 to reach a factor 6 +13-3 at the core centre ($\sim$$\times $ 106 cm-3). Because of growing deuterium fractionation, N2D+ abundance still increases towards the PSC center until the N2 depletion becomes predominant over the deuterium enhancement, and in turn, the N2D+ abundance slightly decreases in the inner most part of the core (Pagani et al. 2007).

The N2D+ and N2H+ ions can be produced via the following routes:


     
    $\displaystyle {\rm H_3^+ + N_2} \rightarrow{\rm N_2H^+ + H_2}$ (4)
    $\displaystyle {\rm H_2D^+ + N_2} \rightarrow {\rm N_2D^+ + H_2} \quad\quad\mbox{(for~1/3)}$ (5)
    $\displaystyle \qquad\qquad\;\;\; \rightarrow {\rm N_2H^+ + HD} \quad\quad\mbox{(for~2/3)}$ (6)
    $\displaystyle {\rm D_2H^+ + N_2} \rightarrow {\rm N_2D^+ + HD} \quad\quad\mbox{(for~2/3)}$ (7)
    $\displaystyle \qquad\qquad\;\;\; \rightarrow {\rm N_2H^+ + D_2} \quad\quad\mbox{(for~1/3)}$ (8)
    $\displaystyle {\rm D_3^+ + N_2} \rightarrow{\rm N_2D^+ + D_2.}$ (9)

We assumed that all the H3+ isotopologues react at the Langevin rate  $k_{\rm N_2}$ with N2 (which is inversely proportional to the square root of the reduced mass of the colliding system, hence to the mass of the H3+ isotopologue) and that deuterium and hydrogen nuclei are equiprobably transferred. Consequently, H3+, H2D+, D2H+ and D3+ respectively produce an N2D+:N2H+ ratio of 0:3, 1:2, 2:1 and 3:0. The measured ratio of 0.7 $\pm$ 0.15 in the center of L183 thus implies significant fractions of D2H+ and D3+. It has been shown (Walmsley et al. 2004) that in the case of complete depletion of heavy species (C, N, O ...), D3+ would be the dominant trihydrogen cation isotopologue which would imply that N2D+ is more abundant than N2H+. This is not the case here; nevertheless the N2H+ deuterium fractionation is a good constraint to the abundance of the four trihydrogen cation isotopologues in our chemical model.

At steady state (d[N2H+]/dt = 0 and d[N2D+]/dt = 0), reaction 4 and its isotopic variants (5 to 9) being the main path to produce N2D+ and N2H+, we obtain:

\begin{displaymath}%
{\rm\frac{[N_2D^+]}{[N_2H^+]} = \frac{[H_2D^+]+2[D_2H^+]+3[D_3^+]}{3[H_3^+]+2[H_2D^+]+[D_2H^+]}}\cdot
\end{displaymath} (10)

This ratio has been measured in the cut through the L183 main PSC. We describe in the following how our method can provide an estimate of the ortho/para H2D+ ratio using this variable and subsequently, of the ortho/para H2 ratio itself as well as some indication of the cosmic ionization rate and mean grain size.

3.4 Grain distribution

Recombination of ions with electrons on negatively charged grain surfaces is an important process since it can be much faster than in the gas phase, especially in the case of H+ (Draine & Sutin 1987). The negatively charged grain surface area is therefore a crucial parameter (we can safely ignore positively and multiply negatively charged grains, considered to be very rare in cold environments, Draine & Sutin 1987). The grain size distribution in prestellar cores is unknown since it mostly depends on grain condensation and also on ice condensation (e.g. Vastel et al. 2006). We thus treat the grain radius  $a_{\rm gr}$ as a parameter of the model, assuming all the grains to have the same size and the dust to gas mass ratio to be 0.01. Different values could be advocated (for example, in places where depletion is important, the ices increase the dust mass by up to 25%, Walmsley et al. 2004) but the net result is to change only slightly the average grain size which is not well constrained in the PSCs in any case.

The focusing effect of the Coulomb attraction between charged particles and oppositely charged grains has been included using the Draine & Sutin (1987) formalism:

\begin{displaymath}%
\tilde{J}(Z=-1) = \left(1+\frac{1}{\tau}\right)\left(1+\sqrt{\frac{2}{2+\tau}}\right)
\end{displaymath} (11)

where $\tau$ is the reduced temperature ( $\tau=a_{\rm gr}kT/e^2$, e being the electron charge, k the Boltzman constant). Therefore the recombination rate of the H+ ion on a negatively charged grain can be expressed as:

\begin{displaymath}%
k_{\rm gr}=\sqrt{\frac{8~kT}{\pi m_{\rm H}}}\pi a_{\rm gr}^2(S \times \tilde{J}(Z=-1))
\end{displaymath} (12)

where $a_{\rm gr}$ is the grain radius, $m_{\rm H}$ is the hydrogen mass, T the kinetic temperature and S is the sticking coefficient (S $\leq$ 1). The latter represents the probability that a colliding species will stick onto the grain surface. For ions, Draine & Sutin (1987) concluded that the sticking coefficient should be unity. The same computation can be made to estimate the recombination rate of other ions, H3+, H2D+, D2H+, D3+, $\cdots$ by a simple correction of the atomic mass of the ions (respectively $k_{\rm gr}$/$\sqrt{3}$, $k_{\rm gr}$/$\sqrt{4}$, $k_{\rm gr}$/$\sqrt{5}$, $k_{\rm gr}$/$\sqrt{6}$, $\cdots$). In the case of collisions between charged particles and neutral grains, the attraction due to the polarization of the grain by the charged particle can be expressed through:

\begin{displaymath}%
\tilde{J}(Z=0) = 1+\sqrt{\frac{\pi}{2\tau}}\cdot
\end{displaymath} (13)

Therefore the sticking rate of electrons on neutral grains can be expressed as:

\begin{displaymath}%
k_{\rm e} = \sqrt{\frac{8~kT}{\pi m_{\rm e}}}\pi a_{\rm gr}^2 (S \times \tilde{J}(Z=0))
\end{displaymath} (14)

where $m_{\rm e}$ is the electron mass and S is the sticking coefficient. S is about unity (Umebayashi & Nakano 1980) for a planar surface but a curvature of the grain surface will tend to reduce this parameter. However in the following we will use a factor of about unity as this parameter did not seem to have a large influence on the results in our runs.

The grain abundance [gr] can be expressed using:

\begin{displaymath}%
[{\rm gr}] = \frac{m_{\rm H_2}f_{\rm d/g}}{ \frac{4 \pi}{3} a_{\rm gr}^3 \delta}
\end{displaymath} (15)

where $\delta$ is the mean grain density (assumed to be 3 g cm-3, $f_{\rm d/g}$ is the dust-to-gas mass ratio, and $m_{\rm H_2}$ is the mass of molecular hydrogen. Another important parameter in our model is the abundance of the negatively charged grains ( $\rm [gr]=[gr^0]+[gr^-]$). At steady-state, assuming partial depletion of CO and N2 and total depletion of all the other heavy species:

\begin{eqnarray*}{\frac{{\rm d}[{\rm gr}^-]}{{\rm d}t}}& = &\rm [gr^0][e^-]{\it ...
...[gr^-][HCO^+]{\it k}_ {HCO^+} -[gr^-][DCO^+]{\it k}_ {DCO^+} = 0
\end{eqnarray*}


(here we have neglected HD+, D2+, He+, ...).

3.5 Steady-state chemical model

The code we describe in the following is used to calculate the steady-state abundances of the chemical species found in the different layers of the L183 prestellar core as listed in Table 2.

In the steady-state approximation the abundance species are interlinked via their production rates and their destruction rates (production = destruction).

Since H3+ is produced at a rate 0.96$\zeta $, the H+ abundance can be expressed as (including only the main reactions):

\begin{displaymath}%
[\rm H^+] = \frac{0.04\zeta}{{\it k}_ {rec}[e^-]{\it n}_ {H_2}+{\it k}_ {gr}[gr^-]{\it n}_ {H_2}}\cdot
\end{displaymath} (16)

The main production path of H3+ is via cosmic ray ionization of H2 and proceeds in two steps:

\begin{displaymath}%
\rm\zeta + o{-} H_2 \rightarrow o{-} H_2^+ + e^-
\end{displaymath} (17)


\begin{displaymath}%
\rm\zeta + p{-} H_2 \rightarrow p{-} H_2^+ + e^-
\end{displaymath} (18)

and H2+ rapidly reacts with another H2 to form H3+ but the branching ratios between different combinations of spin states are non-trivial (Oka 2004):
$\displaystyle %
\rm p{-} H_2^+ + p{-} H_2 \rightarrow p{-} H_3^+ + H~~~~~~~~~~$     (19)
       
$\displaystyle \rm p{-} H_2^+ + o{-} H_2 \rightarrow p{-} H_3^+ + H~~~~~2/3$     (20)
$\displaystyle \rm\rightarrow o{-} H_3^+ + H~~~~~1/3$      
       
$\displaystyle \rm o{-} H_2^+ + p{-} H_2 \rightarrow p{-} H_3^+ + H~~~~~2/3$     (21)
$\displaystyle \rm\rightarrow o{-} H_3^+ + H~~~~~1/3$      
       
$\displaystyle \rm o{-} H_2^+ + o{-} H_2 \rightarrow p{-} H_3^+ + H~~~~~1/3$     (22)
$\displaystyle \rm\rightarrow o{-} H_3^+ + H~~~~~2/3.$      

These are different from those advocated by FPdFW who took branching ratios of 1/2 for both species. The ortho-H3+ formation rate from cosmic ray ionization $k_{\rm cr{-} o}$ is therefore the sum of several terms:
$\displaystyle %
k\rm _{cr{-} o} = \rm0.96(1/3[p{-} H_2][o{-} H_2^+]+1/3[o{-} H_2][p{-} H_2^+]
+ \rm 2/3[o{-} H_2][o{-} H_2^+]).$     (23)

The production rate for ortho-H3+ can be expressed as (including only the main reactions. The rates are listed in Table A.1):


$\displaystyle k\rm _{cr{-} o}\zeta + \it k\rm\_1_{oood}[o{-} H_2][o{-} H_2D^+]\...
..._3^+][o{-} H_2]\it n\rm _{H_2}+\it k\rm 1_{pdod}[HD][p{-} H_3^+]\it n\rm _{H_2}$     (24)

which represent respectively the formation from cosmic ray ionization, backward destruction of ortho-H2D+ with ortho-H2, spin conversion of para-H3+ with ortho-H2 and finally, spin conversion of para-H3+ with HD.

Table 2:   Source parameters: distance from the core center, H2 density, N2H+ abundance and N2D+/N2H+ ratio. The position is measured away from the PSC center along the RA axis (from Pagani et al. 2007).

The destruction rate for ortho-H3+ can be expressed as (including only the main reactions):


                                 $\displaystyle (o{-} k\rm _{rec1}[e^-]+(\it k\rm0_{oopp}+\it k\rm0_{oopo})[o{-} H_2]$  
    $\displaystyle \qquad +(\it k\rm 1_{odpd}+{\it k}\rm 1_{odpo}+
{\it k}\rm 1_{odop}+{\it k}\rm 1_{odoo})[HD]$  
    $\displaystyle \qquad +{\it k}\rm _{co}[CO]+{\it k}\rm _{N2}[N_2]+k_{gr1}[gr^-]))[o{-} H_3^+]\it n\rm _{H_2}$ (25)

which respectively represents its destruction by dissociative recombination with electrons, spin conversion with ortho-H2, spin conversion and deuteration with HD, proton transfer reactions with CO and N2 and dissociative recombination on grains.

Similarly, the para-H3+ formation from cosmic ray ionization can be expressed as:

                           $\displaystyle %
k\rm _{cr{-} p}$ = $\displaystyle \rm0.96([p{-} H_2][p{-} H_2^+]+2/3[p{-} H_2][o{-} H_2^+]$  
    $\displaystyle \rm + 2/3[o{-} H_2][p{-} H_2^+]+1/3[o{-} H_2][o{-} H_2^+]).$ (26)

The production rate for para-H3+ is:


                                       $\displaystyle k\rm _{cr{-} p}\zeta+({\it k}\rm\_1_{oopd}[{\rm o{-} H_2}][o{-} H_2D^+]$  
    $\displaystyle \qquad + ({\it k}\rm0_{oopp}+{\it k}\rm0_{oopo})[o{-} H_3^+][o{-} H_2]+{\it k}\rm 1_{odpd}[HD][o{-} H_3^+]$  
    $\displaystyle \qquad + {\it k}\rm\_1_{popd}[o{-} H_2][p{-} H_2D^+])\it n\rm _{H_2}$ (27)

and the destruction rate is:


                                    $\displaystyle (p{-} {\it k}\rm _{rec1}[e^-]+({\it k}\rm0_{poop}+{\it k}\rm0_{pooo})[o{-} H_2]$  
    $\displaystyle \qquad +({\it k}\rm 1_{pdod}+{\it k}\rm 1_{pdpo}+
{\it k}\rm 1_{pdop}+{\it k}\rm 1_{pdoo})[HD]$  
    $\displaystyle \qquad +{\it k}\rm _{co}[CO]+{\it k}\rm _{N2}[N_2]+{\it k}\rm _{gr1}[gr^-]))[p{-} H_3^+]\it n\rm _{H_2}.$ (28)

The N2 abundance has been solved numerically to obtain the observed N2H+ abundance. Electronic abundance is adjusted to reach equilibrium.

In our steady-state model, the H2 ortho/para ratio, the average grain radius and the cosmic ionization rate $\zeta $ are the varying input parameters. Within each layer of the PSC model (Table 2), these parameters are adjusted to match the following:

  • the H2 density;

  • the N2D+/N2H+ ratio at 7 K;

  • the observed ortho-H2D+ column density;

  • the upper limit on the p-D2H+ column density.
Although the full range of grain sizes and ortho-to-para H2 ratios have been explored for each H2 density, we have not allowed solutions in which, for example, the grain size would oscillate from one layer to the next. We have searched for solutions throughout the layers in two different scenarios: (i) the grain size and the ortho-to-para ratio were both kept constant; or (ii) the grain size was allowed to increase and the ortho-to-para H2 ratio to decrease with the H2 density. In both scenarios, $\zeta $ was kept constant throughout the layers. We neglected detailed reactions with D2 as Flower et al. (2004) have shown that its role is negligible in general and we have kept the HD abundance constant which is generally a good approximation.

3.6 Time-dependent chemistry

In a second step, we have constructed a pseudo time-dependent model based on NAHOON, a chemical model, a version of which has been made publicly available by Wakelam[*]. We have modified this model in two ways: 1) we have replaced electron (resp. ion) reactions with neutral (resp. charged) grains as provided in the Ohio State University (OSU) reaction file (delivered with NAHOON) by the set of equations described above (Sects. 3.4 and 3.5), which we have directly included in the program, to take into account Coulomb focusing; 2) we have included the formation of HD and D2 on grain surfaces and we have introduced the spin state of H2 and D2 with the usual assumption that they are formed with the statistical ortho/para spin state ratio of 3 and 2 respectively. We have used the formation rate provided in the OSU reaction file (5 $\times $ 10-17 cm3 s-1) for the formation of molecular hydrogen. Because grains are covered by ice in the environments concerned here, we consider that the only interaction between the atoms and the surface is physisorption. In this case, the formation rates of HD and D2 (in cm3 s-1) is half for HD and 105 times lower for D2 with respect to H2 formation (as dicussed in ). In environments where grains are not covered by icy mantles, on the other hand, one would have to consider chemisorption, which strongly changes the efficiencies of the formation of HD and D2 (Cazaux et al. 2008). We have reduced the set of species and reactions to our needs, limiting ourselves to the most important reactions (see below) but differentiating all ortho and para (plus meta-D3+) species as independent species and including all the detailed rates between the trihydrogen cation and dihydrogen isotopologues (see Sect. 3.7.1) including spin state conversions. However we have included more reactions than in the steady-state model, taking into account reactions with D2, H2+, He+, etc. and allowing the ortho/para H2 ratio and the HD abundance to vary.

The main path to convert ortho-H2 into para-H2 is via the reaction

\begin{displaymath}%
\rm o{-} H_2 + H^+ \rightarrow p{-} H_2 + H^+
\end{displaymath} (29)

which proceeds seven orders of magnitude faster at 7 K than the reverse reaction.

3.7 Rate coefficients

Many groups have made available gas-phase rate coefficients. The University of Manchester Institute of Science and Technology (UMIST) Database for Astrochemistry contains information on 4500 reactions of which 35% have been measured experimentally, some at temperatures down to 20 K (Woodall et al. 2007). The OSU group provides approximately the same database but focuses more on low temperature chemistry. We accordingly use in our modeling some of the reactions in the latter (with the most recent version OSU2007), considering the low temperatures found in L183. Apart from the rates presented below, all reaction rates involving deuterium have been taken from FPdFW (except recombination on grains for which we have used different sticking probabilities). We disregard the odd branching ratio of N2H+ dissociative recombination reported by Geppert et al. (2004) to consider a single possibility, namely the liberation of dinitrogen (Molek et al. 2007).

3.7.1 H3+ + H2 isotopologue reaction rates

Phase space theory (PST) was used to derive thermal state-to-state rate coefficients for the whole H3+ + H $_2 \rightarrow$ H3+ + H2 system and isotopic variants in the temperature range 5-50 K. This statistical method accounts for such quantities as mass, energy, rotational angular momentum, nuclear spin symmetry and their respective conservation laws. The ergodic hypothesis which is a requisite for PST as well as the full-scrambling hypothesis are assumed according to the topology of the PES (Xie et al. 2005; Yamaguchi et al. 1987). Reactant (product) trajectories are treated with the classical Langevin model. The resulting set of state-to-state rate coefficients deviates from the detailed balance principle by a few percent at worst and is consistent with thermodynamical equilibrium constants. Details will be given in a forthcoming publication (Hugo et al., in prep.)

In the present astrochemical model, nuclear spin states of the different molecules are treated as distinct species but their rotational states are not considered individually. We thus made the assumption that only the rotational ground states of each nuclear spin species were populated and used the ground state-to-species thermal rate coefficients obtained by summing ground state-to-state thermal rate coefficients over the product channels.

 \begin{figure}
\par\includegraphics[angle=-90,width=16.3cm,clip]{0587fig5.eps}
\end{figure} Figure 5:

N2D+/N2H+ ratio as a function of the ortho/para H2 ratio for all possible electronic abundances and total depletion (the steady-state chemical model). The lower row corresponds to the densest part of the PSC (n(H2) = 2 $\times $ 106 cm-3) and the two horizontal dashed lines the measured range of the N2D+/N2H+ ratio while the upper row corresponds to the external part of the PSC (n(H2) = 7 $\times $ 104 cm-3) with the dashed line representing the N2D+/N2H+ ratio upper limit. The three columns represent different $\zeta $ values as indicated above ($\zeta _0$ = 1 $\times $ 10-17 s-1). The large arrow indicates the direction of increasing electronic abundance.

Open with DEXTER

3.7.2 H3+ isotopologue dissociative recombination rates

The dissociative recombination (DR) rate of H3+ is difficult to measure experimentally. Since the early measurements e.g. by Adams et al. (1984), numerous attempts have been made and are summarized in two papers (Fonseca dos Santos et al. 2007; Bates et al. 1993). In parallel, theoretical work has also been developed with the latest one published by Fonseca dos Santos et al. (2007). Extending upon that work, we present here, in Table B.1, the H3+ updated DR rate (Fonseca dos Santos et al. 2007) along with newly calculated H2D+, D2H+, and D3+ DR rates (see Appendix B for more details). These calculations do not predict the branching ratio of the DR products. We have thus adopted the branching ratios published elsewhere (Datz et al. 1995a; Zhaunerchyk et al. 2008; Strasser et al. 2004; Datz et al. 1995b) which we have applied to the calculated rates. The resulting DR rates at 7 K are listed in Table A.1. Several remarkable effects at low temperature are visible (see Figs. B.1-B.4):

  • the strong departure of the ortho-H3+ DR rate from that of para-H3+. The difference is a factor of 10 at 10 K;

  • D2H+ shows a large DR rate drop, by a factor of 10 at 10 K for both ortho and para species compared to the extrapolated value used by FPdFW;

  • on the contrary, a large increase of the D3+ DR rate is predicted to occur but mostly at temperatures where deuteration is low and therefore the consequences for D3+ abundance are limited.

4 Results and discussion

4.1 H2D+ linewidth

In order to fit the observed ortho-H2D+ line profiles, we run the Monte Carlo model using the ``best model'' velocity, density and temperature profiles derived from the N2H+ and N2D+ data analyzed in Pagani et al. (2007). However, the linewidth for the two stronger spectra (offsets (0, 0) and (15 $^{\prime \prime }$, 0), Fig. 3) is too wide to be reproduced with the same micro-turbulent width which we have used for N2H+ ( $\Delta v_{\rm turb}$(FWHM$\approx$ 0.14 km s-1). Indeed, the central spectrum linewidth measured by fitting a Gaussian yields 0.45 ($\pm$0.03) km s-1. The thermal contribution is

\begin{displaymath}%
\Delta v_{\rm therm}(FWHM) = 2.336 \times \sqrt{\frac{kT}{m}} = 0.28~{\rm km~s^{-1}}
\end{displaymath} (30)

at 7 K, which implies that $\Delta v_{\rm turb}$(FWHM) contribution is 0.35 km s-1, 2.5 times larger than for N2H+. Here k is the Boltzmann constant and m is the mass of H2D+. If we impose a turbulent velocity similar to the one modeled for N2H+, then the temperature needed to obtain such a wide line would be 16 K which is completely ruled out by N2H+ observations (Pagani et al. 2007). Infall motion limited to the inner core could have explained the H2D+ width if no H2D+ spectrum other than the central one had been wide combined with a large depletion of N2H+. This is not possible because the H2D+ spectrum at (15 $^{\prime \prime }$, 0) has the largest width (0.49 $\pm$ 0.03 km s-1) in a region where N2H+ is hardly depleted. This remains therefore a pending problem.

4.2 N2H+ deuteration

 \begin{figure}
\par\includegraphics[angle=-90,width=16.3cm,clip]{0587fig6.eps}
\end{figure} Figure 6:

Same as Fig. 5 but with a CO/H2 abundance of 10-5 in the outer layer (n(H2) = 7 $\times $ 104 cm-3) and 10-6 in the inner layer (n(H2) = 2 $\times $ 106 cm-3).

Open with DEXTER

4.2.1 Requested conditions

We next discuss the main parameters that control the N2H+ deuteration using the steady-state chemical model.

The models have been run for a temperature of 7 K which prevails in all the layers where N2D+ has been detected in the PSC cut presented in Pagani et al. (2007). We have also run the models for the corresponding density, N2H+ abundance and N2D+/N2H+ ratio of each layer (the parameters are listed in Table 2).

As discussed above, the abundance of ortho-H2 is the main controlling factor of the trihydrogen cation isotopologue abundances and therefore of the N2D+/N2H+ ratio (and similarly of the DCO+/HCO+ ratio, see e.g. Maret & Bergin 2007). We have therefore explored the range of possible solutions for the ortho/para H2 ratio in the two extreme layers of our core profile (n(H2) = 7 $\times $ 104 cm-3 and 2 $\times $ 106 cm-3) for which we have a N2D+/N2H+ ratio of <0.03 and 0.7 $\pm$ 0.15 respectively (Fig. 5). We have done this for three cosmic ray ionization rates (10-17, 10-16, and 10-15 s-1) covering the values generally discussed in the literature (e.g. Maret & Bergin 2007) and for all possible electronic abundances (or average grain sizes as they are linked via the abundance of H+ which is mostly controlled by the grain surface area). In this first run, we have simulated total depletion by adjusting the CO and N2 abundance[*] to 10-8. We have also indicated the range of N2D+/N2H+ ratio measured in both layers. The average grain radius has been varied from 0.01 $\mu $m to 5 $\mu $m and electronic abundance from 10-11 to 10-6 which cover the usually accepted values. We can see that N2D+/N2H+ ratios above 100 are possible in dense gas though they require very low electronic abundances and therefore very small grains which are probably absent from these dense and cold regions due to grain coagulation (see e.g. Stepnik et al. 2003).

In the lower density outer layer where no N2D+ has been detected, the ortho/para H2 ratio must be high enough, i.e. above $\sim$0.01, to prevent any deuteration from occuring whatever the cosmic ray ionization rate. On the contrary, the dense, strongly deuterated layer has solutions only below a maximum ortho/para H2 ratio of 0.01 (or lower for high $\zeta $ rates). Thus the ortho/para H2 ratio across the PSC clearly must vary from above 0.01 to below 0.01. In the case of low cosmic ray ionization rates (10-17 s-1), though the ortho/para H2 ratio of 0.01 seems to be a common solution for both layers, it requires a high electronic abundance (and large grains) in the outer layer and a low electronic abundance (and small grains) in the inner dense part. This is clearly improbable. The temperature being low enough in all the layers, warm layers (above 20 K) cannot be invoked instead of a high ortho/para H2 ratio to limit the deuteration in the outer parts of the PSC. CO total depletion is however questionable and in the model we also used a CO depletion factor of 10 (abundance of 10-5) in the outer layer and a CO depletion factor of 100 (abundance of 10-6) in the inner layer (Fig. 6). This only limits the maximum N2D+/N2H+ ratio which decreases by one order of magnitude. Indeed, the destruction of H3+ by CO dominates over recombination with electrons when their abundance is very low and vice versa. However, the conditions to reach the observed N2D+/N2H+ ratio remain unchanged and therefore only a variable ortho/para H2 ratio can be invoked. Such a variable ortho/para H2 ratio cannot be investigated with a steady-state model because in all layers, the ortho-H2 abundance would eventually decrease to values about of 10-3-10-4 as discussed by FPdFW.

4.2.2 The ortho/para H2 variation

We discuss here the possibilities of making the ortho/para H2 ratio vary across a single PSC.

It is commonly accepted that H2 is formed on grain surfaces with an ortho/para ratio of 3 because of spin statistics and the exothermicity of the reaction H + H  $\rightarrow$ H2. Subsequently, the ortho-H2 is converted into para-H2 following Eq. (29) and to a lesser extent with reactions involving H3+ and its isotopologues. As already discussed by FPdFW, this conversion is slow and has probably not reached steady state in clouds with ages between 105 and 106 years.

We have therefore run the modified Nahoon model to search for a simultaneous solution for all layers which would get the necessary ortho/para H2 gradient across the PSC to achieve both the observed N2D+/N2H+ ratio and ortho-H2D+ abundance profiles. In this model, we have no direct measurement of the CO abundance but for the credibility of the model, we have set the CO abundance to 10-5 in the outer layer and increased the depletion to reach a factor of 100 (i.e. a CO abundance of 10-6) in the core. The results are presented in Figs. 7 and 8. We have searched for a solution where each modeled layer meets the two observational constraints (N2D+/N2H+ ratio and ortho-H2D+ abundance) at the same time, but these solutions are not simultaneous between the different layers. We varied the average grain radius and the cosmic ray ionization rate, $\zeta $. We could find solutions for grains of average radius 0.025 to 0.3 $\mu $m. No solution has been found for grains above 0.3 $\mu $m. Figures 7 and 8 show the case for which the grain average radius is 0.1 $\mu $m and $\zeta $ = 2 $\times $ 10-17 s-1. In this case, the time range inside which all layers meet the requested conditions is 0.6 to 1.7 $\times $ 105 years. Figure 8 shows how the ortho/para H2 ratio evolves for 3 selected layers. We have marked the appropriate time which is the solution for each of these layers as established from Fig. 7. In that figure, we can see that the o/p H2 ratio is below 0.01 for the dense layers and still above 0.01 for the outer layer for which no H2D+ has been detected, as expected from the steady-state model. We can also see that the full o/p H2 relaxation has not yet occurred even for the densest part of the cloud. Smaller grains have a larger interacting surface and therefore lower the abundance of H+ ions which preferentially recombine on negatively charged grains (or PAHs-), consequently slowing down the dominant ortho-H2 relaxation reaction (Eq. (29)). Though smaller grains also imply a lower electronic density, therefore favouring a higher deuteration of N2H+ as shown in Fig. 5, the slower disappearance of ortho-H2 is the dominant process here and finally, smaller grains slow down the deuteration process. For grains of average radius 0.025 $\mu $m, the range of ages matching the range of N2H+/N2D+ observed ratios is 2.7-3.8 $\times $ 105 years, while for grains of average radius 0.3 $\mu $m, the time range is only 3.8-7.2 $\times $ 104  years. Figure 8 suggests that D2 should become a sizeable fraction of available deuterium a short while after the present state (typically 2-3 $\times $ 105  years) and that HD correspondingly should drop slightly.

 \begin{figure}
\par\includegraphics[width=14.5cm,clip]{0587fig7CMYK.eps}
\end{figure} Figure 7:

Time-dependent variation of the N2D+/N2H+ ratio for the 6 layers defined in Table 2. For each colour, the density of the layer is given. In the insert, the ortho-H2D+ abundance is represented with the same color code, zoomed on the epoch of interest. Horizontal dotted lines represent the observed N2D+/N2H+ ratio range for each layer and the observed ortho-H2D+ abundance as derived from the Monte Carlo model applied to the JCMT observations. The part of the chemical solution that fits in both these limits and the common time limits is set in bold. Vertical arrows indicate upper limits for the N2D+/N2H+ ratio. Vertical dashed lines are placed at 0.63 and 1.7 $\times $ 105 years to delimit the period when all layers reach their observed N2D+ enrichment. This case has been computed for an average grain radius of 0.1 $\mu $m and $\zeta $ = 2 $\times $ 10-17 s-1.

Open with DEXTER
 \begin{figure}
\includegraphics[width=14.5cm,clip]{0587fig8CMYK.eps}
\end{figure} Figure 8:

Time-dependent variation of the ortho- and para-H2 and other related species for 3 of the 6 layers defined in Table 2 (layers 1, 3 and 6). Vertical lines are placed at the times corresponding to the observed N2D+/N2H+ ratio in Fig. 7 for the same 3 densities. This case has been computed for an average grain radius of 0.1 $\mu $m and $\zeta $ = 2 $\times $ 10-17 s-1.

Open with DEXTER

4.3 Age of the core and collapse

Though it is normal that dense layers evolve faster than less dense ones, at least to account for a differential ortho/para H2 ratio, the layers are evolving too fast in our model. The densest layers would have reached their present status 2 to 3 times faster than the outer, less dense layers. This could be possible only if the denser layers had reached their steady-state equilibrium. This is not the case here where the densest layer would reach its steady-state equilibrium only after 2 $\times $ 105  years and this would imply a N2D+/N2H+ ratio of 6, almost an order of magnitude higher than observed. The most probable reason for this time discrepancy is that the core has undergone a contraction and therefore all layers were not so dense in the past. While the outer layers have only little evolved in density (the most external one probably started at 0.5-1 $\times $ 104 cm-3 to reach 7 $\times $ 104 cm-3 today), the inner ones have undergone a much greater density increase. As constant density through the core would give no chemical differentiation while a time-frozen density profile as measured here gives too much differentiation, the solution is in between the two. Starting from a uniform gas, the chemical differential evolution of the core should therefore help us to constrain the duration of the contraction and the type of contraction. Of course the model should also include the evolution of depletion which also plays a role in the acceleration of the process.

As the core must have started to contract from a lower density region, typically 104 cm-3, it is clear that all layers have accelerated their chemical evolution while their density was increasing. Therefore, the layer with the longest time to reach the observed N2D+/N2H+ ratio gives a lower limit to the age of the cloud. Depending on the exact average size of the grains, this lower limit is 1.5-2 $\times $ 105 years here. It is even larger because before the cloud underwent contraction, depletion had not yet occured and therefore species like atomic sulfur, S, must have been present in quantities large enough to transfer notable quantities of electric charges from H+ to S+ (H+ + S $\to$ H + S+) and PAHs must have also been abundant enough to help destroy H+ ions (Wakelam & Herbst 2008). All these phenomena contribute to the diminution of the H+ abundance, therefore slowing down the ortho-H2 relaxation process. Indeed, Flower et al. (2006b) show that the relaxation process in some cases takes 3 $\times $ 107  years, typically 50 times slower than in the case presented here (and 15 times slower for similar conditions of grain size and cosmic ray ionization rate but without depletion in their case).

4.4 Para-D2H+

At 7 K, the strongest possible line intensity (LTE case) for the ground transition of para-D2H+ is below 0.3 K because of the Rayleigh-Jeans correction at 691 GHz which becomes large. Moreover the thermalization of the line is difficult to obtain in this source beyond the radius of 3000 AU because the PSC density drops below the para-D2H+ critical density ( $n_{\rm crit}$ = 4.9 $\times $ 105 cm-3) and a slight drop of the excitation temperature turns into an exponential decrease of the brightness temperature. For $T_{\rm ex}$ = 6 K, $T_{\rm bright}$ $\leq$ 0.13 K. Therefore searches for para-D2H+ must reach very low noise levels to have a chance of detection. From the models we present here, we predict an integrated line intensity of 11 mK km s-1 (34 mK peak) with the Monte Carlo model, and an upper limit of 16 mK km s-1 (48 mK peak) in the case of LTE. This is a factor of 3 to 4 below the upper limit we obtained from the observations.

4.5 The chemical profile

 \begin{figure}
\par\includegraphics[width=11.5cm,clip]{0587fig9CMYK.eps}
\end{figure} Figure 9:

PSC profile for the different species. The n(H2), N2H+ abundance, and N2D+/N2H+ ratio profiles are input data. The four trihydrogen cation isotopologue profiles are also grouped together to visualize their relative total abundances and H+ is compared to both e- (upper right box) and to H3+ (2nd upper left box). The profile has been computed for the case presented in Fig. 7, i.e. $\zeta $ = 2 $\times $ 10-17 s-1, $a_{\rm gr}$ = 0.1 $\mu $m. In the two upper boxes, the green curve refers to the green axis on the right. In the top left box, the arrow indicates the N2D+/N2H+ ratio upper limit for that layer.

Open with DEXTER

Finally, we obtain a detailed profile of the PSC which we present in Fig. 9. It represents the solution for the model we presented in Figs. 7 and 8 taking for each layer the values at their respective best-fit time. The large variation of the ortho-H2 species across the core (a factor of 15) makes D3+ change by an even larger amount (2 orders of magnitude) but it does not become the most abundant trihydrogen cation isotopologue in the core center at this stage because the ortho-H2 abundance is not yet low enough. In the present case, for a density of 2 $\times $ 106 cm-3, the inversion between H3+ and D3+ occurs when the abundance of ortho-H2 drops below 3 $\times $ 10-3.

The ortho/para H2 ratio in the outer layer is 0.04 (as expected from Fig. 6 which indicates a lower limit of 0.01). As discussed above, the ortho/para H2 ratio evolution speed is linked to density and grain size, both of which are lower outside the PSC, in its embedding parental cloud. We can thus expect this ratio to be at least 0.05 and probably above 0.1 in the envelope of the cloud.

N2 is an input parameter in our model because we have not included all the nitrogen chemistry. As reactions with H3+ isotopologues are the main path to destroy this molecule (Flower et al. 2006b), we do not introduce a large error in determining its abundance directly from the N2H+ abundance itself and obtain a N2 profile which is probably closer to reality than if we had let the whole N-chemistry freely establish its abundance because of too many unknowns. The abundance profile thus starts at 1.5 $\times $ 10-7 with respect to H2 in the low density layer to diminish to 3 $\times $ 10-8 in the densest layer. The undepleted N2 abundance after attaining the steady state is $\sim$$\times $ 10-5 (Flower et al. 2006b) but this is reached only after $\approx$$\times $ 106 years. As depletion of N2 in the outer layer is possibly still small, we can conclude that N2 has not yet reached its steady state abundance which puts an upper limit on the age of the cloud of about 1 $\times $ 105 years following the estimate of Flower et al. (2006b). However, this depends very much upon several factors, for example the C:O elemental abundance ratio. Consequently this information is only indicative.

Though our model does not deal with the nitrogen chemistry, it seems to indicate that low abundances of N2 are sufficient to explain the observed N2H+ abundance and therefore, the N2 depletion seems not to be a critical factor as long as CO is also depleted. The much debated contradiction between the presence of N2H+ in depleted cores while N2 should deplete like CO would thus be a false problem. Low CO and electronic abundances, limiting the destruction rate of N2H+, seem to be sufficient to compensate for the N2 depletion itself to a large extent. Finally, in the inner core, N2 and N2H+ (+N2D+) follow a similar decreasing trend, suggesting that N2 depletion eventually forces N2H+ to decrease.

While ortho-H2D+ is 83 K above the para ground state, it is more abundant all over the core profile by almost an order of magnitude whereas the thermal equilibrium ratio would be ortho/para H2D+ $\approx$$\times $ 10-5 at 7 K. This demonstrates the efficiency of the ortho-H2 in converting para-H2D+ into ortho-H2D+ and therefore limits the total abundance of H2D+, as the backward channel to H3+ remains open even at 7 K. This ortho/para population inversion does not occur for D2H+ as the species needed to perform this inversion is no longer ortho-H2 but the much rarer ortho-D2. Therefore, the para-D2H+ remains the least abundant of the two spin state species which, combined with the fact that its ground transition is higher in frequency than the one of ortho-H2D+, makes its detection extremely difficult.

5 Conclusions

We have presented a pair of simple chemical models restricted to H-carriers, He plus CO and N2 to account for the observed HCO+, DCO+ (not discussed in this paper), N2H+, and N2D+ ions. We have benefited from new computed reactions rates for both the H3+ + H2 isotopologue combinations and for the H3+ isotopologue dissociative recombination rates which explicitly take into account the nuclear states individually.

With the steady-state model we have shown that the ortho/para ratio of H2 must vary from above 0.01 in the outer parts of the L183 PSC to less than 0.01 in the inner parts to explain the variation of deuteration across the core. Checking the reality of the ortho/para H2 variation with a time-dependent model, we have also found that if the present PSC density profile is static, then the inner layer would have reached its present status 2 to 3 times faster than the outer layers. Because the present status is not in steady-state, the layers should evolve at a similar rate and therefore the density must have been lower in the past. The most probable explanation is that the core has probably evolved from a uniform density cloud to the present centrally condensed PSC. The time-dependent model also suggests that the ortho/para H2 ratio changes by one order of magnitude from $\sim$5% at a density of 7 $\times $ 104 cm-3 down to a few $\times $ 10-3 in the inner dense core. This has two important consequences:

  • It is most probable that most of the cloud, outside the densest regions (i.e. the two PSC and the ridge in between) have an ortho/para H2 ratio also above 5%, and possibly 10%, contrarily to what is usually assumed in models.

  • In principle, it should be possible to fit the PSC profile with this chemical model combined with a dynamical model including depletion, to set an age to this PSC and possibly discriminate between several types of collapse but it is beyond the scope of this paper.
We already have some indication that the age of the PSC is somewhat above 1.5-2 $\times $ 105 years though the N2 abundance suggests a relatively short time (105 years, except if depletion is compensating for its formation) when we compare our adjusted abundance to the formation rate of N2 given by Flower et al. (2006b). The low abundance of N2 needed to explain the observed N2H+ abundance indicates that its depletion is not a real problem, though, obviously, N2H+ would be much more abundant if N2 was not being depleted but CO still was.

Finally, we stress the importance of considering ortho/meta/para chemistry when dealing with the deuteration of the interstellar medium. The importance of the ortho-H2 on the amount of deuteration and the observations limited to only ortho-H2D+ make this inclusion compulsory. Moreover, a complete state-to-state chemical model should be developed to take into account rotational pumping, leading to a higher destruction rate of the deuterated trihydrogen cation and possibly explaining the observed linewidth of ortho-H2D+.

Detecting para-D2H+ would be highly desirable to help constrain the models, but the high frequency and limited transparency of the atmosphere make it a difficult tool to use. Though observable from the ground, because of its weakness in cold dark clouds, which are the only places where it should be found, direct para-D2H+ observations should be made on a large number of Galactic lines of sight using the HIFI receiver on board the Heschel Space Observatory.

Acknowledgements
We would like to thank D. Flower, G. Pineau des Forêts, M. Walmsley, S. Cazaux and V. Wakelam for fruitful discussions and an anonymous referee for her/his careful reading. Part of this work was supported by the National Science Foundation (NSF) grant AST 05-40882 to the CSO. The authors are grateful to the CSO and JCMT staffs for their support at the telescopes. We would like to thank Samantha Santos for the help provided in numerical calculations of thermal DR rate coefficients. This work has been supported by the NSF under Grants No. PHY-0427460 and PHY-0427376, by an allocation of NERSC and NCSA (project # PHY-040022) supercomputing resources. This work has benefited from research funding from the European Community's Sixth Framework Programme under RadioNet contract R113CT 2003 5158187.

References

Online Material

Appendix A: The reaction rate table used in the Nahoon modified chemical model

Table A.1:   Reactions used in the model. The rate coefficients are given for 7 K. Reaction rates less than 10-15 cm-3 s-1 are not taken into account in our models. Reference (1) corresponds to Gerlich (1990), references (2)-(4) correspond to Hugo, OSU 07, and this paper respectively. For OSU 07, branching ratios involving spin states have been infered from quantum mechanical rules. For reactions involving grains, a grain radius of 0.1 $\mu $m and a sticking coefficient of 1 have been considered. (5) Datz et al. (1995b) (6) Datz et al. (1995a), (7) Zhaunerchyk et al. (2008), (8) Larsson et al. (1997), (9) Molek et al. (2007).

Appendix B: On the rate coefficients for dissociation recombination of H3+ and its isotopologues

The dissociative recombination (DR) rate coefficients for ortho- and para-H3+ have been published by Fonseca dos Santos et al. (2007). Here, we present the results obtained for all four H3+ isotopologues. The DR rate coefficients for different species of the nuclear spin are calculated using the approach described in a series of papers devoted to DR theory for triatomic molecular ions. See Fonseca dos Santos et al. (2007); Kokoouline & Greene (2003a,b) for H3+ and D3+ calculations and Kokoouline & Greene (2004,2005) for H2D+ and D2H+. The scope of this paper does not allow us to review the theoretical approach in detail. We only list its main ingredients.

The theoretical approach is fully quantum mechanical and incorporates no adjustable parameters. It relies on ab initio calculations of potential surfaces for the ground electronic state of the H3+ ion and several excited states of the neutral molecule H3, performed by Mistrik et al. (2001).

The total wave function of the system is constructed by an appropriate symmetrization of products of vibrational, rotational, electronic, and nuclear spin factors. Therefore, rovibronic and nuclear spin degrees of freedom are explicitly taken into account.

The electronic Born-Oppenheimer potentials for the four H3+ (and H3) isotopologues have the C3v symmetry group. The C3v symmetry group has a two-dimensional irreducible representation E. The ion has a closed electronic shell. The lowest electronic state of the outer electron in H3 has p-wave character. The p-wave state of the electron also belongs to the E representation. Due to the Jahn-Teller theorem (Landau & Lifshitz 2003), this leads to a strong non-adiabatic coupling between the E-degenerate vibrational modes of the ion and the p-wave states of the incident electron. The coupling is responsible for the fast DR rate (Kokoouline et al. 2001) in H3+. In the present model, only the p-wave electronic states are included because other partial waves have a much smaller effect on the DR probability: the s-wave states have no E-type character and, therefore, are only weakly coupled to the dissociative electronic states of H3; d-wave electronic states are coupled to the E-vibrational modes, but the coupling is rather small because the d-wave of the incident electron does not penetrate sufficiently close to the ionic core owing to the d-wave centrifugal potential barrier.

All three internal vibrational coordinates are taken into account. Vibrational dynamics of the ionic core are described using the hyper-spherical coordinates, which represent the three vibrational degrees of freedom by a hyperradius and two hyperangles. The hyperradius is treated as a dissociation coordinate that represents uniformly the two possible DR channels, the three-body (such as H+H+H) and two-body (such as H2+H). Although the initial vibrational state of the ion is the ground state, after recombination with the electron, other vibrational states of the ionic target molecule can be populated. Therefore, in general, many vibrational states have to be included in the treatment. In particular, the states of the vibrational continuum have to be included, because only such states can lead to the dissociation of the neutral molecule. The vibrational states of the continuum are obtained using a complex absorbing potential placed at a large hyperradius to absorb the flux of the outgoing dissociative wave.

Since the rovibrational symmetry is D3h for H3+ and D3+ and C2v for H2D+ and D2H+, the rovibrational functions are classified according to the irreducible representations of the corresponding symmetry groups, i.e. A1', A1'', A2', A2'', E', and E'' for D3h and A1, A2, B1, and B2 for C2v. We use the rigid rotor approximation, i.e. the vibrational and rotational parts of the total wave function are calculated independently by diagonalizing the corresponding Hamiltonians. In our approach, the rotational wave functions must be obtained separately for the ions and for the neutral molecules. They are constructed in a different way for the D3h and C2v cases. The rotational eigenstates and eigenenergies of the D3h molecules are symmetric top wave functions (see, for example Bunker & Jensen 1998). They can be obtained analytically if the rotational constants are known. The rotational constants are obtained numerically from vibrational wave functions, i.e. they are calculated separately for each vibrational level of the target molecule. The rotational functions for the C2v ions are obtained numerically by diagonalizing the asymmetric top Hamiltonian (Kokoouline & Greene 2005; Bunker & Jensen 1998).

Once the rovibrational wave functions are calculated, we construct the electron-ion scattering matrix (S-matrix). The S-matrix is calculated in the framework of quantum defect theory (QDT) (see, for example, Chang & Fano 1972; Kokoouline & Greene 2005,2003b) using the quantum defect parameters obtained from the ab initio calculation (Mistrik et al. 2001). The constructed scattering matrix accounts for the Jahn-Teller effect and is diagonal with respect to the different irreducible representations $\Gamma$ and the total angular momentum N of the neutral molecule. Thus, the actual calculations are made separately for each $\Gamma$ and N. Elements of the matrix describe the scattering amplitudes for the change of the rovibrational state of the ion after a collision with the electron. However, the S-matrix is not unitary due to the presence of the dissociative vibrational channels (i.e. continuum vibrational states of the ion, discussed above), which are not explicitly listed in the computed S-matrix. The ``defect'' from unitarity of each column of this S-matrix is associated with the dissociation probability of the neutral molecule formed during the scattering process. The dissociation probability per collision is then used to calculate the DR cross-sections and rate coefficients.

The nuclear spin states are characterized by one of the A1, A2, or E irreducible representations of the symmetry group S3 for D3h molecules and by the A or B irreducible representations of the symmetry group S2 for C2v molecules. The irreducible representation  $\Gamma_{\rm ns}$ of a particular nuclear spin state determines its statistical weight and is related to the total nuclear spin $\vec I$ of the state. Here, $\vec I$ is the vector sum of spins $\vec i$ of identical nuclei.

For H3+, the $\Gamma_{\rm ns}=A_1$ states (A2' and A2'' rovibrational states) correspond to I=3/2 (ortho); the $\Gamma_{\rm ns}=E$ states (E' and E'' rovibrational states) correspond to I=1/2 (para). The statistical ortho:para weights are 2:1.

For H2D+, the $\Gamma_{\rm ns}=A$ states (B1 and B2 rovibrational states) correspond to I=1 (ortho); the $\Gamma_{\rm ns}=B$ states (A1 and A2 rovibrational states) correspond to I=0 (para). The statistical ortho:para weights are 3:1.

For D2H+, the $\Gamma_{\rm ns}=A$ states (A1 and A2 rovibrational states) correspond to I=0,2 (ortho); the $\Gamma_{\rm ns}=B$ states (B1 and B2 rovibrational states) correspond to I=1 (para). The statistical ortho:para weights are 2:1.

Finally, for D3+, the $\Gamma_{\rm ns}=A_1$ states (A1' and A1'' rovibrational states) correspond to I=1,3 (ortho); the $\Gamma_{\rm ns}=A_2$ states (A2' and A2'' rovibrational states) correspond to I=0 (para); the $\Gamma_{\rm ns}=E$ states (E' and E'' rovibrational states) correspond to I=1,2 (meta). The statistical ortho:para:meta weights are 10:1:8.

 \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg10.eps}
\end{figure} Figure B.1:

Theoretical DR rate coefficients as functions of temperature for the ortho- and para-species of H3+. The figure also shows the species-averaged rate coefficient. For comparison, we show the rate coefficient obtained from measurements in the TSR storage ring by Kreckel et al. (2005) and the analytical dependence for the coefficient used in earlier models of prestellar cores by FPdFW.

Open with DEXTER

 \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg11.eps}
\end{figure} Figure B.2:

Theoretical DR rate coefficients for the ortho- and para-species of H2D+, as functions of temperature. The rate coefficient averaged over the two species and the analytical expression used in earlier models are also shown.

Open with DEXTER

 \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg12.eps}
\end{figure} Figure B.3:

Theoretical DR rate coefficients for the ortho- and para-species of D2H+ as a function of temperature. The rate coefficient averaged over the two species and the analytical expression used in earlier models are also shown.

Open with DEXTER

 \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg13.eps}
\end{figure} Figure B.4:

Theoretical DR rate coefficients for the ortho-, para-, and meta-species of D3+. The rate coefficient averaged over the two species and the one used in earlier models are also shown.

Open with DEXTER

Figures (B.1)-(B.4) summarize the obtained DR thermal rate coefficients calculated separately for each nuclear spin species of the four H3+ isotopologues and the numerical values are listed in Table B.1. For comparison, the figures also show the analytical dependences used in previous models of prestellar core chemistry (FPdFW). As one can see, the rates for different nuclear spin species are similar to each other (for a given isotopologue) at high temperatures. However, for lower temperatures, the rates for different ortho/para/meta-nuclear spin species significantly differ from each other. The difference in behavior at low temperatures is explained by different energies of Rydberg resonances present in DR cross-sections at low electron energies. The actual energies of such resonances are important for the thermal average at temperatures below or similar to the energy difference between ground rotational levels of different nuclear spin species. At higher temperatures, the exact energy of the resonances is not important. The averaged rate is determined by the density and the widths of the resonances, which are similar for all nuclear spin species over a large range of collision energies.

Table B.1:   Dissociative recombination rates of H3+, H2D+, D2H+, and D3+ for each individual nuclear spin state species.


Footnotes

... ratio[*]
Appendices A and B are only available in electronic form at http://www.aanda.org
... jiggle-chop[*]
http://www.jach.hawaii.edu/software/jcmtot/het_obsmodes.html
... Wakelam[*]
http://www.obs.u-bordeaux1.fr/radio/VWakelam/Valentine%20 Wakelam/Downloads.html
... abundance[*]
Here, we look for general solutions of the N2D+/N2H+ ratio as a function of several parameters, we thus do not try to fit the N2 abundance to obtain the observed N2H+ abundance.

All Tables

Table 1:   Line parameters from the JCMT and CSO observations. The positions are offsets to the dust peak emission at $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. For non-detected lines we give the 3$\sigma $ upper limit. For JCMT and p-D2H+ at CSO, we give both the Monte Carlo (MC) and the LTE column density estimates.

Table 2:   Source parameters: distance from the core center, H2 density, N2H+ abundance and N2D+/N2H+ ratio. The position is measured away from the PSC center along the RA axis (from Pagani et al. 2007).

Table A.1:   Reactions used in the model. The rate coefficients are given for 7 K. Reaction rates less than 10-15 cm-3 s-1 are not taken into account in our models. Reference (1) corresponds to Gerlich (1990), references (2)-(4) correspond to Hugo, OSU 07, and this paper respectively. For OSU 07, branching ratios involving spin states have been infered from quantum mechanical rules. For reactions involving grains, a grain radius of 0.1 $\mu $m and a sticking coefficient of 1 have been considered. (5) Datz et al. (1995b) (6) Datz et al. (1995a), (7) Zhaunerchyk et al. (2008), (8) Larsson et al. (1997), (9) Molek et al. (2007).

Table B.1:   Dissociative recombination rates of H3+, H2D+, D2H+, and D3+ for each individual nuclear spin state species.

All Figures

  \begin{figure}
\par\includegraphics[width=6cm,clip]{0587fig1.eps}
\end{figure} Figure 1:

CSO map of the H2D+ (110-111) line. The position is indicated in arcseconds for each spectrum and the (0, 0) position corresponds to $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. The Y-axis represents the antenna temperature.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{0587fig2.eps}
\end{figure} Figure 2:

JCMT map of the H2D+ (110-111) line. The dust peak position (Pagani et al. 2003) is indicated by a cross and corresponds to $\alpha _{2000}$ = $\rm 15^h54^m08\hbox{$.\!\!^{\rm s}$ }50$, $\delta _{2000}$ = $-02^\circ 52^\prime 48^{\prime \prime }$. Contour levels are 0.1, 0.2 and 0.3 K km s-1.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{0587fig3.eps}
\vspace*{5mm}
\end{figure} Figure 3:

Ortho-H2D+ spectra across the main dust peak. East and west sides are folded together and fitted with a Monte Carlo model. Density and temperature profiles are taken from Pagani et al. (2007). The spacing between spectra is 15 $^{\prime \prime }$.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8.5cm,clip]{0587fig4.eps}
\vspace*{5mm}
\end{figure} Figure 4:

Main reactions involved in the H3+ chemical network. When CO and N2 are depleted, the reactions with bold arrows are dominant.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=16.3cm,clip]{0587fig5.eps}
\end{figure} Figure 5:

N2D+/N2H+ ratio as a function of the ortho/para H2 ratio for all possible electronic abundances and total depletion (the steady-state chemical model). The lower row corresponds to the densest part of the PSC (n(H2) = 2 $\times $ 106 cm-3) and the two horizontal dashed lines the measured range of the N2D+/N2H+ ratio while the upper row corresponds to the external part of the PSC (n(H2) = 7 $\times $ 104 cm-3) with the dashed line representing the N2D+/N2H+ ratio upper limit. The three columns represent different $\zeta $ values as indicated above ($\zeta _0$ = 1 $\times $ 10-17 s-1). The large arrow indicates the direction of increasing electronic abundance.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=16.3cm,clip]{0587fig6.eps}
\end{figure} Figure 6:

Same as Fig. 5 but with a CO/H2 abundance of 10-5 in the outer layer (n(H2) = 7 $\times $ 104 cm-3) and 10-6 in the inner layer (n(H2) = 2 $\times $ 106 cm-3).

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=14.5cm,clip]{0587fig7CMYK.eps}
\end{figure} Figure 7:

Time-dependent variation of the N2D+/N2H+ ratio for the 6 layers defined in Table 2. For each colour, the density of the layer is given. In the insert, the ortho-H2D+ abundance is represented with the same color code, zoomed on the epoch of interest. Horizontal dotted lines represent the observed N2D+/N2H+ ratio range for each layer and the observed ortho-H2D+ abundance as derived from the Monte Carlo model applied to the JCMT observations. The part of the chemical solution that fits in both these limits and the common time limits is set in bold. Vertical arrows indicate upper limits for the N2D+/N2H+ ratio. Vertical dashed lines are placed at 0.63 and 1.7 $\times $ 105 years to delimit the period when all layers reach their observed N2D+ enrichment. This case has been computed for an average grain radius of 0.1 $\mu $m and $\zeta $ = 2 $\times $ 10-17 s-1.

Open with DEXTER
In the text

  \begin{figure}
\includegraphics[width=14.5cm,clip]{0587fig8CMYK.eps}
\end{figure} Figure 8:

Time-dependent variation of the ortho- and para-H2 and other related species for 3 of the 6 layers defined in Table 2 (layers 1, 3 and 6). Vertical lines are placed at the times corresponding to the observed N2D+/N2H+ ratio in Fig. 7 for the same 3 densities. This case has been computed for an average grain radius of 0.1 $\mu $m and $\zeta $ = 2 $\times $ 10-17 s-1.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=11.5cm,clip]{0587fig9CMYK.eps}
\end{figure} Figure 9:

PSC profile for the different species. The n(H2), N2H+ abundance, and N2D+/N2H+ ratio profiles are input data. The four trihydrogen cation isotopologue profiles are also grouped together to visualize their relative total abundances and H+ is compared to both e- (upper right box) and to H3+ (2nd upper left box). The profile has been computed for the case presented in Fig. 7, i.e. $\zeta $ = 2 $\times $ 10-17 s-1, $a_{\rm gr}$ = 0.1 $\mu $m. In the two upper boxes, the green curve refers to the green axis on the right. In the top left box, the arrow indicates the N2D+/N2H+ ratio upper limit for that layer.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg10.eps}
\end{figure} Figure B.1:

Theoretical DR rate coefficients as functions of temperature for the ortho- and para-species of H3+. The figure also shows the species-averaged rate coefficient. For comparison, we show the rate coefficient obtained from measurements in the TSR storage ring by Kreckel et al. (2005) and the analytical dependence for the coefficient used in earlier models of prestellar cores by FPdFW.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg11.eps}
\end{figure} Figure B.2:

Theoretical DR rate coefficients for the ortho- and para-species of H2D+, as functions of temperature. The rate coefficient averaged over the two species and the analytical expression used in earlier models are also shown.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg12.eps}
\end{figure} Figure B.3:

Theoretical DR rate coefficients for the ortho- and para-species of D2H+ as a function of temperature. The rate coefficient averaged over the two species and the analytical expression used in earlier models are also shown.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8cm,clip]{0587fg13.eps}
\end{figure} Figure B.4:

Theoretical DR rate coefficients for the ortho-, para-, and meta-species of D3+. The rate coefficient averaged over the two species and the one used in earlier models are also shown.

Open with DEXTER
In the text


Copyright ESO 2009

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

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

Initial download of the metrics may take a while.