Issue |
A&A
Volume 509, January 2010
|
|
---|---|---|
Article Number | A98 | |
Number of page(s) | 10 | |
Section | Interstellar and circumstellar matter | |
DOI | https://doi.org/10.1051/0004-6361/200913350 | |
Published online | 26 January 2010 |
Modelling line emission of deuterated H3+ from prestellar cores
O. Sipilä1,3 - E. Hugo2 - J. Harju1,3 - O. Asvany2 - M. Juvela1,3 - S. Schlemmer2
1 - Observatory, PO Box 14, 00014 University of Helsinki, Finland
2 - I. Physikalisches Institut, Universität zu Köln, Germany
3 -
Present address: Department of Physics, PO Box 64, 00014 University of Helsinki, Finland
Received 25 September 2009 / Accepted 5 November 2009
Abstract
Context. The depletion of heavy elements in cold cores of
interstellar molecular clouds can lead to a situation where deuterated
forms of
are the most useful spectroscopic probes of the physical conditions.
Aims. The aim is to predict the observability of the rotational lines of
and
from prestellar cores.
Methods. Recently derived rate coefficients for the
isotopic system were applied to the ``complete depletion'' reaction
scheme to calculate abundance profiles in hydrostatic core models. The
ground-state lines of
(372 GHz) and
(692 GHz) arising from these cores were simulated. The excitation
of the rotational levels of these molecules was approximated by using
the state-to-state coefficients for collisions with
.
We also predicted line profiles from cores with a power-law density distribution advocated in some previous studies.
Results. The new rate coefficients introduce some changes to the
complete depletion model, but do not alter the general tendencies. One
of the modifications with respect to the previous results is the
increase of the
abundance at the cost of other isotopologues. Furthermore, the present model predicts a lower H2D+ (o/p) ratio, and a slightly higher D2H+
(p/o) ratio in very cold, dense cores, as compared with previous
modelling results. These nuclear spin ratios affect the detectability
of the submm lines of H2D+(o) and D2H+(p). The previously detected H2D+ and D2H+ lines towards the core I16293E, and the H2D+ line
observed towards Oph D can be reproduced using the present
excitation model and the physical models suggested in the original
papers.
Key words: ISM: clouds - radio lines: ISM - astrochemistry - radiative transfer - ISM: abundances - ISM: molecules
1 Introduction
The rate coefficients for reactive and inelastic collisions between the H3+ ion and H2 with all possible deuterated variants and nuclear spin symmetries have been recently derived and a few of them have been experimentally tested by Hugo et al. (2009). Besides making state-specific astrochemical models viable, the cross-sections derived in this work make it possible to model the intensities of the dipole-allowed rotational transitions of H2D+ and D2H+.
The H2D+ and D2H+ ions are potentially useful probes of pre-protostellar cores (see e.g. Pagani et al. 2009; Caselli et al. 2003; Vastel et al. 2004; van der Tak et al. 2005). Together with the isotopologues H3+ and D3+, they are likely to belong to the most abundant ions in the dense and cold nuclei of prestellar cores, where the usual tracer molecules can have practically disappeared from the gas phase (Walmsley et al. 2004; Roberts et al. 2003) - this state of matters is called ``complete depletion'' in the latter paper. While it has been observed that heavier substances such as CN can exist in the gas phase at densities of the order of 106 cm-3 (Hily-Blant et al. 2008), the assumption of complete depletion remains a valid first approximation at high densities.
The astrophysical significance of these ions is related to the facts that H3+ originates almost directly in the cosmic ray ionization of H2, and that it initiates the ion-molecule reactions in dense clouds (Herbst & Klemperer 1973). Furthermore, the deuteration of H3+ depends on the ortho-para ratio of H2, which is likely to be non-thermal and evolve with time in interstellar clouds (Pagani et al. 2009; Gerlich et al. 2002; Pagani et al. 1992; Flower et al. 2006).
In this paper we apply the newly derived chemical rate coefficients to the ``completely depleted'' case first discussed by Walmsley et al. (2004), and use the state-to-state coefficients to predict the ground-state rotational lines of ortho-H2D+ and para-D2H+ from prestellar cores. The two core models used in these simulations correspond to the observed properties of the cores Oph D and I16293E in Ophiuchus. The rate coefficients of Hugo et al. (2009) for some important deuteration reactions differ from those adopted in Roberts et al. (2003), and in the series of papers by Flower, Pineau des Forêts & Walmsley, and we discuss changes implied to the complete depletion model. The Hugo et al. (2009) rate coefficients have been previously used by Pagani et al. (2009) in the modelling of the prestellar core L183, where also CO, N2 and their derivatives were included in the chemical reaction scheme. The organization of the present paper is as follows: in Sect. 2, we describe the physical and chemical models and, in Sect. 3, we present the modelling results. The results are discussed in Sect. 4, and finally, the main points of this discussion are summed up in Sect. 5.
2 Model
In what follows, we describe the model assumptions. These include the physical core model, the description of the dust grain component, the chemical reaction scheme, and the adopted rate coefficients. We also give a brief account of the methods and programs used in solving the chemical abundances and the emitted molecular line radiation.2.1 Core model
For the core model we used a modified Bonnor-Ebert sphere (BES), which is a non-isothermal cloud in hydrostatic equilibrium (previously discussed by, e.g., Evans et al. 2001), heated externally by the interstellar radiation field, ISRF. The possibility of a hydrostatic or near-equilibrium configuration is suggested by observations towards some prestellar cores which seem to represent an advanced stage of chemical evolution (characterized by a high degree of molecular depletion), and have near-thermal linewidths (e.g., Bergin & Tafalla 2007).
We first calculated a density profile resulting from an isothermal model (Bonnor 1956). The density profile was then fed into a Monte Carlo radiative transfer program (Juvela & Padoan 2003; Juvela 2005)
for the temperature calculation. For this purpose the core model was
divided into concentric shells. We used the grain size distribution of Mathis et al. (1977) (MRN) with two types of dust grains: carbonaceous and silicate (the grain model is further discussed in Sect. 2.2.4). The choice of using the MRN distribution implies that the abundance of very small grains (with radii less than 0.01
m)
is assumed to be negligible, and that grain coagulation is not
considered. The temperature profile was calculated separately for the
two dust types. In each position of the core, the final temperature
profile was taken to be the arithmetic average of the temperatures of
both grain types.
The gas and dust temperatures were assumed to be equal. While this is probably a good approximation in the dense (>
)
central part of the core (Burke & Hollenbach 1983), it may not hold as well in the outer layers of the core (
,
Bergin et al. 2006). In the present
analysis, which concentrates on the very dense nucleus of the core, we
have ignored this phenomenon. The model core is assumed to be embedded
in a molecular cloud with extinction corresponding to
.
This assumption agrees with observations from the vicinity of Oph D (Chapman et al. 2009).
To determine the exact core density profile using the computed
temperature distribution, we slightly modified the method used in Bonnor (1956) by making the following substitutions:
![]() |
(1) |
![]() |
(2) |
where






Equation (3) was integrated numerically and the resulting density profile was fed back into the radiative transfer program to calculate a new temperature profile. This procedure was repeated a few times until the density profile converged.
2.2 Chemistry model
We adopted the complete depletion model discussed in Walmsley et al. (2004) and Flower et al. (2004). The reaction scheme is fairly simple: the gas-phase reactants include only H, H+, H2, H2+, H3+ and their deuterated forms, He, He+, and free electrons. Because of the low temperature associated with the environment, the zero-point energies become relevant and thus the deuterated forms and different nuclear spin modifications have to be considered explicitly. The reaction set also includes cosmic ray ionization, H2 formation on grains, and some other grain processes. We next discuss the gas-phase and grain-surface reactions separately.
2.2.1 Gas phase chemical reactions
In all our models, we assumed that the core is electrically neutral, i.e. that electron abundance equals the difference between the total abundance of positive ions and the abundance of negatively charged grains.
The chemical reactions and the associated rate coefficients were compiled using four sources: Walmsley et al. (2004), Flower et al. (2004), Hugo et al. (2009), and Pagani et al. (2009). The reaction set of Walmsley et al. (2004) was complemented by Flower et al. (2004) to include the nuclear spin modifications of ,
,
and
.
Hugo et al. (2009) took the analysis of the H3+ + H2
reacting system further by using a microcanonical approach to derive
state-to-state rate coefficients. This approach differs from that of Flower et al. (2004)
in that all internal states of the reacting system are explicitly taken
into account, resulting in detailed state-to-state rate coefficients
for different configurations of the system. We adopted their ground
state-to-species rate coefficients for all reactions relevant to the H3+ + H2 system (including the deuterated forms) and substituted them into the combined data of Walmsley et al. (2004) and Flower et al. (2004). Furthermore, for dissociative recombination (DR) reactions of
and its deuterated forms we used the newly calculated rate coefficients presented in
Appendix B of Pagani et al. (2009).
2.2.2 Grain reactions
The H2 and HD molecules are mainly formed on grain surfaces (Hollenbach & Salpeter 1971; Gould & Salpeter 1963), whereas for the formation of D2, gas-phase reactions are more important. Nevertheless, the formation of all three isotopologues on dust grains, and the corresponding destruction terms for atomic H and D were included in the chemical system. We assumed that ortho and para forms are produced according to their statistical population ratios 3:1 and 2:1 for H2 and D2, respectively.
Grain processes, in particular the formation of ortho-H2 and the attachment and destruction of positive ions on grain surfaces, are important for the overall degree of deuteration of the cloud (see Sect. 4). The details of these processes (grain surface properties, desorption mechanisms, etc.) were not considered in the present model. It was assumed that chemical reactions on grains (along with the associated adsorption and desorption) take place instantaneously.
The rate coefficients for grain reactions (including, e.g.,
electron attachment and recombination of positive ions on grains),
listed in Walmsley et al. (2004), were taken from Flower & Pineau des Forêts (2003), in which an MRN distribution with
m,
m was assumed. These limits correspond to an effective grain radius of
m
. Since the rate coefficient for a grain reaction is proportional to
,
we scaled the Flower & Pineau des Forêts (2003) coefficients by
,
where
is the effective grain radius of the assumed size distribution. The Coulomb factor,
,
which takes into account the electric interaction between dust grains
and gas phase neutrals and/or ions, was included as an additional
factor. Its mathematical form has been discussed in detail in Draine & Sutin (1987), and recently in Pagani et al. (2009).
2.2.3 Molecular line radiation
The abundance distributions of chemical species were determined by a
chemistry program in physical conditions corresponding to those in
different parts of the model core - in practice in each of the
concentric shells used in the temperature calculations. The radial
(o) and
(p)
abundance distributions, together with the density and temperature
profiles were used as input for a Monte Carlo radiative transfer
program (Juvela 1997) to predict observable line emission.
The excitation of the rotational transitions of
(o) and
(p) in collisions with para and ortho H2 were calculated using the state-to-state rate coefficients from Hugo et al. (2009). The data concerning the line frequencies and Einstein A coefficients were obtained from Miller et al. (1989), Ramanlal & Tennyson (2004), and Amano & Hirao (2005).
![]() |
Figure 1: Optical thickness per column density according to the opacity data of Ossenkopf & Henning (1994), denoted in the figure by ``OH 1994'', (solid curve) and Li & Draine (2001), denoted by ``LD 2001'' (dashed curve). The modified extinction curve (see text) is also plotted, denoted by ``LD mod'' (dotted curve). |
Open with DEXTER |
2.2.4 Dust grain model
The MRN model was adopted for the grain size distribution. In the density and temperature profile calculations, we used two types of grains; carbonaceous and silicate. The optical properties of these components were adopted from Li & Draine (2001). The densities of the carbonaceous and silicate grains were taken to be 2.5 g cm-3 and 3.5 g cm-3, respectively. In the chemistry model, no separation between the grain materials was done, and a grain material density of 3.0 g cm-3 was assumed. For the comparison with the Flower et al. (2004) results (Sect. 3.1) we used, however, the same assumptions of the grain material density (2.0 g cm-3) and dust-to-gas mass ratio (0.013) as in Walmsley et al. (2004).
![]() |
Figure 2:
Steady state abundances of electrons, protons, H3+ and its deuterated forms, calculated using the Flower et al. (2004) ( upper panel) and Hugo et al. (2009)
( lower panel) rate coefficients in an isothermal (T = 10 K) model with single size grains (a = 0.1 |
Open with DEXTER |
The model of Li & Draine (2001) describes
dust in diffuse medium. In dense clouds, coagulation and growth of icy
grain mantles are expected to modify the grain properties and,
according to the models of Ossenkopf & Henning (1994),
the extinction curve is expected to become flatter in the far-infrared
and at longer wavelengths. Since in our cores the central densities are
high,
cm-3, we modified the
Li & Draine (2001) dust opacities accordingly.
This leads to temperatures of the order of 6 K in the core center.
The curves are plotted in Fig. 1.
3 Results
In this Section we discuss the differences arising from the use of the modified chemical reaction rate coefficients used in this paper compared to those used in Flower et al. (2004). We then present results of chemical modelling carried out using a hydrostatic core model. We conclude the section with a brief discussion on simulated line emission spectra.
3.1 Comparison with Flower et al. (2004) results using homogeneous models
![]() |
Figure 3:
Steady-state abundance ratios of the nuclear spin modifications of H2, D2, H3+ and its deuterated forms as functions of
|
Open with DEXTER |
![]() |
Figure 4:
Steady-state o/p ratios of |
Open with DEXTER |
We ran a set of homogeneous models corresponding to those presented in Walmsley et al. (2004) and in Flower et al. (2004). The steady state fractional abundances of the principal ions as functions of the gas density are shown in Fig. 2. The upper diagram corresponds to the Flower et al. (2004) model and is a reproduction of their Fig. 1. The lower plot is produced using the Hugo et al. (2009) rate coefficients for the H3+ + H2 isotopic system and the DR rate coefficients from Pagani et al. (2009). The figures were produced assuming an isothermal (T = 10 K) model with single sized grains (a = 0.1 m) and a grain material density of
2.0 g cm-3.
The abundance ratios of the most important nuclear spin variants as functions of the gas density are shown in Fig. 3, and the temperature dependence at low temperatures is illustrated in Fig. 4. The former should be compared with Fig. 2 of Flower et al. (2004) and with Fig. 3 of Walmsley et al. (2004), and the latter with Figs. 5 and 6 in Flower et al. (2004). Following Hugo et al. (2009) we have assigned meta-D3+ with the modification having the lowest ground state energy, corresponding to the A1 representation of the symmetry group S3, and ortho-D3+ with the E representation (see their
Sect. II.A.1). The ortho and meta appelations for
are therefore interchanged with respect to those used in Flower et al. (2004), Pagani et al. (2009), and most other papers. As in Flower et al. (2004), the D3+ (o/m) and D2 (p/o) ratios fall rapidly with increasing
.
However, the
(p/o) ratio is nearly constant unlike in the previous study. The H3+
(o/p) ratio is about three times lower in the present model. There are
also a couple of differences in the temperature dependence
(Fig. 4). The H2D+ (o/p) ratio levels off between 2-3 when approaching very low temperatures. Using the Flower et al. (2004) rate coefficients the corresponding value is
10 in the same conditions. Secondly, the D2H+ (p/o) reaches slightly higher
(0.3 - 0.4) at very low temperatures than in the Flower et al. (2004) model. We note that as our model follows the chemistry until steady state, the H2 (o/p) ratio becomes much lower than in the models of Pagani et al. (2009). As a consequence, we obtain a clearly higher degree of deuteration of H3+
and somewhat different nuclear spin ratios for its isotopologues as
compared with this previous study. We return to this matter briefly in
Sect. 4.1.
All the plots presented so far assume single sized grains with a radius of 0.1 m. As discussed in Sect. 3.3 of Flower et al. (2004), the grain size has a considerable effect on the abundances and nuclear spin ratios. The D2 (o/p), D2H+ (o/p), and D3+ (m/o) ratios peak near the grain radius
m (see Fig. 9 in Flower et al. 2004), whereas the H2 (o/p), H3+ (o/p), and H2D+ (o/p) ratios increase monotonically towards smaller grain radii.
3.2 Hydrostatic models
We studied the influence of the grain size distribution and the
cosmic ionization rate on the chemistry of a depleted, hydrostatic
core. The model BE sphere was assumed to have a central density of
cm-3 and an outer radius of 2400 AU, where the density drops to about 105 cm-3. These parameters correspond roughly to the dense nucleus of the Oph D core according to Harju et al. (2008).
The effective grain radius was varied by changing the lower limit of
the grain size distribution and keeping the upper limit constant at
0.3
m. Modelling was carried out with effective grain radii of 0.05
m, 0.1
m and 0.2
m. The slope of the extinction curve (Fig. 1) was modified to produce a central core temperature of about 6 K (see Sect. 2.2.4). The chemical reaction network was integrated until
years, which was well within steady state (all the abundances settled into steady state at times <106 years).
Figure 5 shows the core density profile corresponding to effective grain radius
m. The temperature profiles for three different grain sizes are superposed. Changing
has a marked effect on the temperature profile, but alters the density
distribution only slightly (the central density is kept fixed in our
models). The core gets colder as the effective grain radius gets
smaller; this is because the total grain surface area (effectively the
total opacity) grows larger as the effective grain radius decreases.
Figure 6
shows the steady state fractional abundance profiles of the principal
ions for the selected values of effective grain radius and for a cosmic
ray ionization rate of
s-1. The grain radius decreases from top to bottom. As shown in Fig. 5,
this decrease is accompanied by a slight drop in the average
temperature. The fractional ion abundances generally decrease towards
the core centre, i.e., towards higher densities, in accordance with the
results from homogeneous models shown in Fig. 2. The density dependence of D3+
is less marked, which leads to an increased degree of deuterium
fractionation in the centre when the effective grain radius decreases.
The fractional H+ abundance decreases drastically with
decreasing grain size. This is caused by accentuated recombinations on
negatively charged grains as explained in Walmsley et al. (2004).
The radial profiles of the
(o/p),
(o/p),
(o/p),
(p/o), and
(o/m) ratios are presented Fig. 7. While most ratios have rather flat distributions, the
(o/m) ratio drops towards the centre, i.e. towards higher densities.
The density dependence is, however, smoother for small grain sizes
(bottom panel). The ortho/para ratios of H2, H3+, and H2D+ increase from top to bottom, along with the decreasing grain size. This is caused by an intensified replenishment of H2(o) when the total grain surface area becomes larger. The slight increase of H2 (o/p) towards the outer parts reflects the decreasing density (see
Fig. 3). The
(o/p) ratio follows the H2 (o/p) ratio closely. The
(o/p) ratio correlates with the H2 (o/p) ratio, but rises towards low temperatures (see Fig. 4). These tendencies result in a slight increase of the
(o/p) ratio towards the core centre. The D2H+ (p/o) ratio changes very little as a function of radius and from model to model. Like the D3+ (o/m) ratio, the D2H+ (o/p) ratio has a flat maximum around
m (Flower et al. 2004), and depends weakly on the density (Fig. 3).
The effect of changes in the cosmic ionization rate, ,
is illustrated in Fig. 8. This shows the radial distributions of the principal ions for the cases where
is made 10 times lower (top) and 10 times higher (bottom) than in the models discussed above. Lowering
decreases
the degree of ionization. An enhanced cosmic ray ionization increases
the ion abundances, but decreases the deuteration through an
intensified electron recombination (Walmsley et al. 2004; Flower et al. 2004). The effect is particularly marked for
near the outer edge of the core.
![]() |
Figure 5:
The radial core density profile corresponding to effective grain radius
|
Open with DEXTER |
![]() |
Figure 6:
Radial fractional (with respect to
|
Open with DEXTER |
3.3 H2D+ and D2H+ spectra
We calculated the spectral line profiles of the ground state transitions of
(
110-111, 372 GHz) and
(
110-101, 692 GHz) from models corresponding to the observed properties of the prestellar cores I16293E and Oph D.
The very narrow 372 GHz line of
detected with APEX towards the density maximum of Oph D has been
suggested to originate in a hydrostatic core (Harju et al. 2008).
In accordance with this suggestion the physical model used for
Oph D is a hydrostatic core heated externally by the interstellar
radiation field (ISRF). Judging from previous maps, in particular the
ISOCAM
m image from Bacmann et al. (2000), we assumed an angular radius
of
.
We adopted the recent distance estimate of 120 pc to the Ophiuchus complex (Lombardi et al. 2008),
which implies that the adopted radius corresponds to 2400 AU. The
obscuration provided by the surrounding molecular cloud was assumed to
correspond to
.
We varied the central H2 density (n0), the effective grain radius (
), and the cosmic ray ionization rate (
),
and calculated the appropriate abundance profiles using the chemistry
model described above. Finally, the populations of the rotational
levels of
,
and the
110-111 line profiles were calculated using the state-to-state rate coefficients for
(para and ortho separately) collisions from Hugo et al. (2009)
and our Monte Carlo radiative transfer program. Calculations with
different values for the model parameters were carried out until a
reasonable agreement with the observed
profile was met. This model was then used to predict the
line at 692 GHz as observed with APEX.
![]() |
Figure 7:
Radial profiles for the ratios of the different spin states of H2, H3+ and its deuterated forms at
|
Open with DEXTER |
![]() |
Figure 8:
Same as the middle panel in Fig. 6, but for cosmic ray ionization rates
|
Open with DEXTER |
The results of our simulations are presented in Fig. 10, which also shows the observed
spectrum. The best fit was obtained using the following parameters:
,
m, and
.
The model core has a central temperature T = 6.36 K. The fractional H2D+(o) abundance is
in the core center, rising to
at the edge. A Gaussian fit to the model spectrum yields an FWHM of 0.35 km s-1; the peak optical depth is
.
A molecular line profile for the D2H+(p) 110-101
transition using the same core model parameters is shown in the lower
panel of the figure. A Gaussian fit to the spectrum yields an FWHM of 0.27 km s-1, the peak optical depth is
.
The D2H+(p) abundance is
at the core center and
at the edge - essentially a flat profile with a minor maximum at around 1200 AU.
Both
and
have been detected with the CSO towards I16293E by Vastel et al. (2004). To our knowledge this work has produced the only rather firm detection of the 692 GHz line so far. Vastel et al. (2004) estimated the
and
column densities from the observed line profiles by assuming line-of-sight homogeneity and an excitation temperature of
K.
We assumed that these column densities are approximately valid, and
adjusted the fractional abundances accordingly. The physical model of
the core was adopted from Stark et al. (2004), allowing for the fact that the cloud is probably nearer than previously thought (120 pc, Lombardi et al. 2008). The model consisted of a compact, homogeneous nucleus with a density
up to a radius of 750 AU, surrounded by an envelope with the density power law
extending to a distance of 6000 AU (
). A dust temperature
K derived by Stark et al. (2004)
was adopted as the temperature of the nucleus. The dust temperature was
gradually increased toward the core edge, reaching a value
K at 6000 AU. It should be noted, however, that the observations of Vastel et al. (2004) suggest a lower value for the region where the
and
emission originate. Assuming thermal broadening only, the line widths listed in their Table 1 imply kinetic temperatures
(
)
and
(
).
These estimates are upper limits as they neglect the possible turbulent
broadening and instrumental effects. In order to reproduce the column
densities derived by Vastel et al. (2004), the fractional
(o) abundance was set to
10-10 in the nucleus, and it was let to fall radially to
10-11 at the outer edge. The
(p) abundance was set at 0.6 times the
(o) abundance, which is within the column density error range reported by Vastel et al. (2004).
The line profiles calculated from this model are shown in Fig. 11. The simulated H2D+(o) line has a peak optical depth
.
A Gaussian fit to the spectrum yields an FWHM of
0.46 km s-1 and a peak antenna temperature
K. The simulated D2H+(p) line yields
,
an
FWHM of
0.41 km s-1 and a peak
K. The antenna temperatures of the simulated H2D+(o) and D2H+(p) spectra agree well with the observations of Vastel et al. (2004),
but the lines are too broad, probably because the assumed kinetic
temperature is too high. We return to this issue in Sect. 4.3.
4 Discussion
We now discuss in more detail the results presented in the last section, starting with the differences in homogeneous models arising from using different rate coefficients. We then discuss the abundance distributions in our hydrostatic core models, and conclude by looking into the simulated line emission spectra.
4.1 Homogeneous models
In the models calculated using the new reaction rate coefficients derived by Hugo et al. (2009), the deuteration of
proceeds faster than in the models of Flower et al. (2004). As a consequence, D3+ becomes the dominant ion at relatively low densities. The total
rate coefficients of Hugo et al. (2009) for the deuteration reactions
are larger than those used by Flower et al. (2004) by an average factor of
4-5. This difference originates in the discrepancy between the laboratory measurements of Gerlich et al. (2002) and Hugo et al. (2009, see their Sect. IV.C)
concerning the forward rate coefficients at very low temperatures. On
the other hand, the total rate coefficients for the ``backward''
reactions (k-1, k-2, k-3)
are similar in the two reaction schemes. This can be traced to the fact
that the ``backward'' rate coefficient for the reaction
(o) +
(o)
(o/p) +
derived by Hugo et al. (2009) agrees with the estimate of Gerlich et al. (2002, Sect. 3.3).
The ratios K=k/k- are thus larger in Hugo et al. (2009) which implies that chemical equilibrium is established at higher concentrations of deuterated species. The equilibrium constants in the present model are, however, smaller than in the models of Roberts et al. (2004) and Caselli et al. (2008). The forward rate coefficients of Hugo et al. (2009) are similar to those used in the two previous studies while the backward rate coefficents are several orders of magnitude larger owing to the non-thermal H2 (o/p) ratio.
As discussed in Flower et al. (2006), the ortho/para ratios of
and
correlate with
through proton exchange reactions. The H3+ ion is primarily produced in the para form via the reaction of
(p) with
(p). The para-ortho and ortho-para conversions through reactions with
(o)
are, however, rapid, and compete only against the deuteration reaction
with HD. In the conditions considered here a good approximation of the
o/p ratio can be obtained from

where












Also for
,
the nuclear spin changing reactions with
are effective and an approximate value of the o/p ratio can be obtained through balancing these according to
Eq. (7) of Gerlich et al. (2002). Below
8 K, however, primary production from
and destruction through the secondary deuteration reaction
start to dominate. At very low temperatures the
(o/p) ratio settles somewhere between 2 and 3 (Fig. 4).
The production and destruction of ortho and para D2H+ occur mainly via the deuteration reactions
and
.
As discussed in Flower et al. (2006), the lower energy ortho state is favoured in the reverse reactions with
(reactions 11 and 12 in their
Sect. 3.4; Table VIII in Hugo et al. 2009). Moreover,
is reduced through para-ortho conversion in the reaction with HD. In the present model the p/o ratio is
0.2-0.4 below 10 K and is not sensitive to
.
We note that while the
(o/p) ratio is similar to that found in the models of Pagani et al. (2009, Fig. 9), the
(p/o)
ratio is larger by a factor of 2-4 in the present model. By performing
some test runs we could establish that the main reason for this
apparent discrepancy is in the different H2 (o/p) ratio. In the time-dependent model of Pagani et al. (2009), the initial H2
(o/p) ratio is larger than unity, and the chemical evolution is
followed until the model matches with the constraints based on
observations of the
ratio toward the L183 core. The resulting H2 (o/p) ratio (
0.005-0.05) is much higher than in our models (less than 10-4). At a high H2 (o/p) ratio,
(p) is mainly destroyed in the reaction
(p) + H2(o)
(p) + HD, and the destruction by HD (including para-ortho conversion) is not significant.
The principal reactions determining the ortho/meta ratio of
are discussed in Flower et al. (2004) (Appendix B, with interchanged meta and ortho appellations with respect to the present work). The (lowest energy) meta form of
is primarily formed via deuteration of
by HD, and via ortho-meta conversion in reaction
.
The destruction of
is overwhelmingly dominated by dissociative recombination with electrons or negatively charged grains. For
,
the reaction
in the backward direction and the ortho-meta
conversion mentioned above are comparable with electron recombination
as destruction mechanisms. The formula presented in the appendix of Flower et al. (2004)
approximates the equilibrium o/m ratio with a high accuracy. The ratio
is proportional to the electron density, and decreases as the gas
density increases. On the other hand, it changes very little as a
function of temperature in the range considered here (see Figs. 3
and 4).
Finally, we note that the
(p/o) ratio correlates with the
(o/m) ratio (see Fig. 3). The higher energy
is principally formed from
via dissociative recombination with free electrons and on dust grains, whereas
forms primarily from
(DR
rate coefficients from Pagani et al. (2009) and branching ratios
from Flower et al. (2004) have been used here). The destruction of
is dominated by deuterium exchange reactions with
and deuteration of
,
while for
the deuteration of
and
are the most important destruction reactions. The resulting equilibrium
(p/o) ratio is roughly an order of magnitude lower than the
(o/m) ratio.
4.2 Hydrostatic models
In the hydrostatic models (Figs. 6 and 7)
combined effects of changing grain size, density and temperature can be
seen. The decrease of the average grain size from the top panel to the
bottom signifies an increase in the total grain surface area. This
results in a diminishing fractional abundance of free electrons and a
growing abundance of negatively charged grains (owing to the sticking
of electrons onto neutral grains), and an increase of the
(o/p) ratio (because of enhanced
formation on grain surfaces). The radial fall off of the density (from
left to right) is accompanied by an elevated temperature at the outer
edge of the core.
The
ion is the most abundant deuterated form of
in the centre of the core. In the following we use the
abundance ratio to quantify the degree of deuteration. The
ratio changes as a function of grain size. For the largest grain size,
m,
is
3.3 in the core center, and
0.08 at the edge. As the grain size is decreased to
m,
increases both in the center (to
8.6) and at the edge (to
0.6). For still smaller grains,
m,
drops slightly (to 8.5) in the center, but increases to
1.4 at the edge.
The changes in the deuterium fractionation as functions of the density and the effective grain radius are plotted in Fig. 9 for isothermal models with
K.
One can see in the plot that for each value of density, there is a
deuteration maximum which shifts to lower values of a as the density is decreased. The sloping down of
on the right, towards larger grains, is caused by an increasing abundance of free electrons. The decrease of
on the left side, towards smaller grains, is caused by the increase in the abundance of H2(o) which hinders deuteration by converting
back to
,
and by the increase of negatively charged grains which replace electrons as the principal destructors of
.
Figure 8 shows radial abundance profiles in a chemical model otherwise similar to the middle panel of Fig. 6, but for cosmic ionization rates
s-1 and
s-1.
It can be immediately seen in Fig. 8, as compared with the middle panel in Fig. 6, that the abundances of
and
follow the change in the cosmic ray ionization rate. The H+ abundance is proportional to
via dissociation of H2 by cosmic rays and via dissociative charge transfer reaction between He+ and H2 (see, e.g. Walmsley et al. 2004, Appendix A). The formation of
depends on
through
,
which is the principal product of the collision between
and a cosmic ray proton (e.g., Pagani et al. 2009; Walmsley et al. 2004). The
and
abundances are dragged upwards with
with an intensified cosmic ray ionization, whereas the
abundance remains roughly constant in the core center and falls off in the outer regions. So, quite reasonably, the
ratio which depends on DR reactions diminishes strongly with an increasing
,
whereas the effects on
and
which depend on HD and
,
are less notable.
![]() |
Figure 9: Dependence of the D fractionation on the average grain size. The curves correspond to different values of the gas density. |
Open with DEXTER |
4.3 Molecular line profiles
The observed
(o)
profile from Oph D can be reproduced reasonably well using a
hydrostatic core model together with the ``complete depletion''
chemistry model with the rate coefficients of Hugo et al. (2009) (see Fig. 10). The model giving the best agreement with the observed
(o) spectrum (see Sect. 3.3) predicts a very weak
(p) line at 692 GHz (
K with APEX). This line has a lower peak opacity (
)
than the
(o) line, whereas the radial
distributions of the two lines are similar. One should note that even an optically thick
(p)
(110-101) line
from a cold source is weak in the brightness temperature scale because
the frequency lies far away from the range where the Rayleigh-Jeans
approximation is valid.
The physical model adopted from Stark et al. (2004) (see Sect. 3.3) for the core I16293E produces easily the bright
(o) and
(p) lines observed by Vastel et al. (2004) (see Fig. 11, solid curves). The relatively high temperature (16-20 K) assumed by Stark et al. (2004)
gives, however, rise to too broad lines as compared with the
observations. Another problem with this assumption is that the
chemistry model predicts a clearly higher D2H+(p)/H2D+(o) abundance ratio around 16 K (see Fig. 7 in Flower et al. 2004).
We carried out another simulation assuming that the kinetic temperature
is 10 K in the center and that it increases gradually to 12 K at
the edge. The fractional H2D+(o) and D2H+(p)
abundances were assumed to be roughly equal, in agreement with the
prediction of the chemistry model. The modelled spectra are shown in
Fig. 11 (dashed curves). Now both intensities and line widths are within the range reported by Vastel et al. (2004). The H2D+(o) and D2H+(p) column densities, both
,
are slightly higher than those derived by Vastel et al. (2004).
![]() |
Figure 10: Predicted line profiles for the H2D+ 110-111 (ortho) and the D2H+ 110-101 (para) transitions as observed with APEX for the best-fit core model (solid curves). Also plotted is the observation toward Oph D by Harju et al. (2008) for the H2D+ 110-111 transition (dotted curve). |
Open with DEXTER |
![]() |
Figure 11: Predicted molecular line profiles for the H2D+ (110-111) and the D2H+ (110-101) transitions for the prestellar core 16293E as observed with the CSO, using the core parameters of Stark et al. (2004) (solid curve). The modified model (see text) is also plotted (dashed curve). |
Open with DEXTER |
5 Conclusions
We studied the chemistry in molecular cloud cores in the complete depletion limit, taking advantage of new state-to-state chemical reaction rate coefficients for the H3+ + H2 reacting system (Hugo et al. 2009). We used the modelling results to predict line emission spectra for the H2D+ 110-111 (372 GHz) and the D2H+ 110-101 (692 GHz) transitions in a core model corresponding to the prestellar core Oph D. We also carried out simulations of the above transitions for the I16293E core. The simulated spectra were compared with observations by Harju et al. (2008) and Vastel et al. (2004) for the respective cores.
In this paper, we used steady-state chemical abundances in the
molecular line simulations. It has been noted in the literature that
the timescale for chemical development can be an order of magnitude
larger than the dynamical timescale, and thus chemical steady-state
might not be reached before collapse begins (see for example Flower et al. 2006).
In general, however, ions reach chemical steady-state earlier than
neutral species, due to the speed of the ion-molecule reactions. In
fact, our chemical models indicate that ions reach steady-state before
years. On the other hand, the relatively high
(o) and
(p)
column densities derived toward Oph D and I16293E imply an appreciable
degree of deuteration, which suggests that the cores are in an advanced
chemical state. Taking into account these arguments, the assumption of
steady-state seems, in the present context, appropriate.
The simulated profile of the H2D+ (110-111)
line from a self-consistent hydrostatic core model (including chemistry
modelling) is in reasonably good agreement with observations toward
Oph D. The model predicts that the D2H+ (110-101)
line at 692 GHz probably cannot be detected in this source. In fact,
this transition is difficult to detect in emission in any cold object
because of its high frequency. The same is true for the H2D+(p) and D2H+(o) ground state lines, which lie at even higher frequencies of 1.4 THz.
The Oph D core seems particularly appropriate for using a hydrostatic
model. However, the chemistry model used here begins to lose its
validity in the outer parts, where the density drops below 106 cm-3. In these conditions one can no longer ignore heavier substances in the gas phase, which modify the abundances of
and its deuterated forms. So, in reality, the
and
abundances which in our models rise toward the core edge are likely to
turn down outside the dense nucleus. This effect can be seen in the
models of Pagani et al. (2009).
It is obvious that chemical modelling should be extended to account for
the presence of heavier species. This would also justify studies of
larger cores, presenting larger density gradients.
Our simulations of molecular line profiles toward the I16293E
core are in good agreement with the observational results presented in Vastel et al. (2004). While line emission simulations utilizing the physical core parameters of Stark et al. (2004) produce enough emission, we found that to produce linewidths comparable to those in Vastel et al. (2004),
the core temperature profile and chemical abundances should be
modified. This modification is also justified by results from chemical
modelling. The H2D+(o) and D2H+(p) column densities are found to become comparable when the temperature drops to 10 K.
The deuteration of H3+ proceeds much faster in the present reaction scheme using the Hugo et al. (2009) results than in the Flower et al. (2004) model where the deuteration rate coefficients were adopted from Gerlich et al. (2002). As a consequence, the steady-state
abundance ratio becomes larger than in the models of Flower, Pineau des
Forêts & Walmsley. However, it should be noted that already in Walmsley et al. (2004), D3+ becomes the most abundant isotopologue of H3+ at very high densities, and that this phenomenon was predicted previously by Roberts et al. (2003).
The low-temperature rate coefficients from Gerlich et al. (2002) are lower than those by Hugo et al. (2009) by a factor of 4. The experimental circumstances have been thoroughly discussed in Hugo et al. (2009),
concluding that further experimental investigations, preferably with a
different setup, are urgently needed. The ortho/para ratios of the
deuterated forms of H3+ are also modified in the new model. Along with the changes in abundances, this has an impact on the observability of H2D+ and D2H+ toward prestellar cores. We find at very low temperatures a lower H2D+ (o/p) ratio, and a slightly higher D2H+ (p/o) ratio than predicted by using the coefficients adopted in Flower et al. (2004).
The state-to-state coefficients calculated by Hugo et al. (2009)
provide an opportunity to refine both the chemistry model and the model
for the excitation of the rotational transitions of interest. Besides
collisions with para and ortho H2, it would seem reasonable to examine the excitation of rotational levels of
and
in collisions with HD. The next step forward is to combine the full
state-to-state reaction network with radiative transfer calculations.
The realization of this
improvement looks particularly feasible for the complete depletion
model.
O.S. and J.H. acknowledge support from the Academy of Finland grant No. 117206. We would also like to thank D. Flower, G. Pineau des Forêts, C.M. Walmsley and L. Pagani for their comments on this paper, and the anonymous referee for his/her helpful report.
References
- Amano, T., & Hirao, T. 2005, J. Mol. Spec., 233, 7 [Google Scholar]
- Bacmann, A., André, P., Puget, J.-L., et al. 2000, A&A, 361, 555 [Google Scholar]
- Bergin, E. A., & Tafalla, M. 2007, ARA&A, 45, 339 [Google Scholar]
- Bergin, E. A., Maret, S., van der Tak, F. F. S., et al. 2006, ApJ, 645, 369 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
- Bonnor, W. B. 1956, MNRAS, 116, 351 [CrossRef] [Google Scholar]
- Burke, J. R., & Hollenbach, D. J. 1983, ApJ, 265, 223 [NASA ADS] [CrossRef] [Google Scholar]
- Caselli, P., van der Tak, F. F. S., Ceccarelli, C., & Bacmann, A. 2003, A&A, 403, L37 [Google Scholar]
- Caselli, P., Vastel, C., Ceccarelli, C., et al. 2008, A&A, 492, 703 [Google Scholar]
- Chapman, N. L., Mundy, L. G., Lai, S.-P., & Evans, N. J. 2009, ApJ, 690, 496 [NASA ADS] [CrossRef] [Google Scholar]
- Draine, B. T., & Sutin, B. 1987, ApJ, 320, 803 [NASA ADS] [CrossRef] [Google Scholar]
- Evans, II, N. J., Rawlings, J. M. C., Shirley, Y. L., & Mundy, L. G. 2001, ApJ, 557, 193 [NASA ADS] [CrossRef] [Google Scholar]
- Flower, D. R., & Pineau des Forêts, G. 2003, MNRAS, 343, 390 [NASA ADS] [CrossRef] [Google Scholar]
- Flower, D. R., Pineau des Forêts, G., & Walmsley, C. M. 2004, A&A, 427, 887 [Google Scholar]
- Flower, D. R., Pineau des Forêts, G., & Walmsley, C. M. 2006, A&A, 449, 621 [Google Scholar]
- Gerlich, D., Herbst, E., & Roueff, E. 2002, Planet. Space Sci., 50, 1275 [NASA ADS] [CrossRef] [Google Scholar]
- Gould, R. J., & Salpeter, E. E. 1963, ApJ, 138, 393 [NASA ADS] [CrossRef] [Google Scholar]
- Harju, J., Juvela, M., Schlemmer, S., et al. 2008, A&A, 482, 535 [Google Scholar]
- Herbst, E., & Klemperer, W. 1973, ApJ, 185, 505 [NASA ADS] [CrossRef] [Google Scholar]
- Hily-Blant, P., Walmsley, M., Pineau Des Forêts, G., & Flower, D. 2008, A&A, 480, L5 [Google Scholar]
- Hollenbach, D., & Salpeter, E. E. 1971, ApJ, 163, 155 [NASA ADS] [CrossRef] [Google Scholar]
- Hugo, E., Asvany, O., & Schlemmer, S. 2009, J. Chem. Phys., 130, 164302 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
- Juvela, M. 1997, A&A, 322, 943 [Google Scholar]
- Juvela, M. 2005, A&A, 440, 531 [Google Scholar]
- Juvela, M., & Padoan, P. 2003, A&A, 397, 201 [Google Scholar]
- Li, A., & Draine, B. T. 2001, ApJ, 554, 778 [Google Scholar]
- Lombardi, M., Lada, C. J., & Alves, J. 2008, A&A, 480, 785 [Google Scholar]
- Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 [NASA ADS] [CrossRef] [Google Scholar]
- Miller, S., Tennyson, J., & Sutcliffe, B. 1989, Mol. Phys., 66, 429 [NASA ADS] [CrossRef] [Google Scholar]
- Ossenkopf, V., & Henning, T. 1994, A&A, 291, 943 [Google Scholar]
- Pagani, L., Salez, M., & Wannier, P. G. 1992, A&A, 258, 479 [Google Scholar]
- Pagani, L., Vastel, C., Hugo, E., et al. 2009, A&A, 494, 623 [Google Scholar]
- Ramanlal, J., & Tennyson, J. 2004, MNRAS, 354, 161 [NASA ADS] [CrossRef] [Google Scholar]
- Roberts, H., Herbst, E., & Millar, T. J. 2003, ApJ, 591, L41 [NASA ADS] [CrossRef] [Google Scholar]
- Roberts, H., Herbst, E., & Millar, T. J. 2004, A&A, 424, 905 [Google Scholar]
- Stark, R., Sandell, G., Beck, S. C., et al. 2004, ApJ, 608, 341 [NASA ADS] [CrossRef] [Google Scholar]
- van der Tak, F. F. S., Caselli, P., & Ceccarelli, C. 2005, A&A, 439, 195 [Google Scholar]
- Vastel, C., Phillips, T. G., & Yoshida, H. 2004, ApJ, 606, L127 [NASA ADS] [CrossRef] [Google Scholar]
- Walmsley, C. M., Flower, D. R., & Pineau des Forêts, G. 2004, A&A, 418, 1035 [Google Scholar]
Footnotes
- ...
m
- The rate coefficients of grain reactions listed in Table
A.1. of Walmsley et al. (2004)
correspond to
m, contrary to what is said in the caption.
- ... total
- By the ``total rate coefficients'' we mean the resultant
forward and backward rate coefficients, k1
and k-1, in the presentation
omitting the nuclear spin variants, e.g.,
. In this presentation, the
formation rate can be written
, where k1 is a combination of seven individual rate coefficients weighted according to the fractions of
or
of the total
abundance (see Table VIII of Hugo et al. 2009). The other deuteration reactions and their reverse reactions are composed of eight reactions. Besides the temperature, the total rate coefficients and the equilibrium constant, K=k/k-, depend on the relative populations of nuclear spin variants, in particular o/p H2 (e.g. Gerlich et al. 2002).
All Figures
![]() |
Figure 1: Optical thickness per column density according to the opacity data of Ossenkopf & Henning (1994), denoted in the figure by ``OH 1994'', (solid curve) and Li & Draine (2001), denoted by ``LD 2001'' (dashed curve). The modified extinction curve (see text) is also plotted, denoted by ``LD mod'' (dotted curve). |
Open with DEXTER | |
In the text |
![]() |
Figure 2:
Steady state abundances of electrons, protons, H3+ and its deuterated forms, calculated using the Flower et al. (2004) ( upper panel) and Hugo et al. (2009)
( lower panel) rate coefficients in an isothermal (T = 10 K) model with single size grains (a = 0.1 |
Open with DEXTER | |
In the text |
![]() |
Figure 3:
Steady-state abundance ratios of the nuclear spin modifications of H2, D2, H3+ and its deuterated forms as functions of
|
Open with DEXTER | |
In the text |
![]() |
Figure 4:
Steady-state o/p ratios of |
Open with DEXTER | |
In the text |
![]() |
Figure 5:
The radial core density profile corresponding to effective grain radius
|
Open with DEXTER | |
In the text |
![]() |
Figure 6:
Radial fractional (with respect to
|
Open with DEXTER | |
In the text |
![]() |
Figure 7:
Radial profiles for the ratios of the different spin states of H2, H3+ and its deuterated forms at
|
Open with DEXTER | |
In the text |
![]() |
Figure 8:
Same as the middle panel in Fig. 6, but for cosmic ray ionization rates
|
Open with DEXTER | |
In the text |
![]() |
Figure 9: Dependence of the D fractionation on the average grain size. The curves correspond to different values of the gas density. |
Open with DEXTER | |
In the text |
![]() |
Figure 10: Predicted line profiles for the H2D+ 110-111 (ortho) and the D2H+ 110-101 (para) transitions as observed with APEX for the best-fit core model (solid curves). Also plotted is the observation toward Oph D by Harju et al. (2008) for the H2D+ 110-111 transition (dotted curve). |
Open with DEXTER | |
In the text |
![]() |
Figure 11: Predicted molecular line profiles for the H2D+ (110-111) and the D2H+ (110-101) transitions for the prestellar core 16293E as observed with the CSO, using the core parameters of Stark et al. (2004) (solid curve). The modified model (see text) is also plotted (dashed curve). |
Open with DEXTER | |
In the text |
Copyright ESO 2010
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.