Online material
Appendix A: The condition for the nuclear spin ratio of H_{2} to slow down deuteration processes
Here we derive the condition determining how much o  H_{2} is needed to slow down deuterium fractionation driven by Reaction (1). It is well established that the destruction of H_{2}D^{+} by electrons and CO suppresses the fractionation as well as by H_{2}. H_{2}D^{+} also reacts with HD to form the doubly deuterated species D_{2}H^{+}, leading to further fractionation. Considering the balance of formation and destruction, the steadystate H_{2}D^{+}/H ratio can be calculated as (A.1)where n(i) is the number density of species i, k_{CO} is the rate coefficient of H_{2}D^{+} with CO, k_{e} is the rate coefficient of dissociative electron recombination of H_{2}D^{+}, k_{HD} is the rate coefficient of H_{2}D^{+} with HD, and k_{1f} and k_{1b} are the effective rate coefficients for Reaction (1) in the forward direction and in the backward direction, respectively. Since the forward Reaction (1) is exothermic, k_{1f} does not depend on the OPR(H_{2}). We can say that the presence of o  H_{2} slows down the fractionation when k_{1b}n(H_{2}) is greater than k_{CO}n(CO), k_{e}n(e), and k_{HD}n(HD).
The endothermicity of Reaction (1) in the backward direction depends on the nuclear spin states of H_{2} and H_{2}D^{+}. The following set of reactions can occur: The rate coefficients of the reactions can be found in Hugo et al. (2009), who calculated them with the assumption that the reaction proceeds by a scrambling mechanism in which all protons are equivalent; results are different if longrange hopping is the effective mechanism. The effective rate coefficient k_{1b} can be defined as follows (Gerlich et al. 2002; Lee & Bergin 2015): (A.6)From Eq. (A.6), we can express k_{1b} as a function of orthotopara ratios of H_{2} and H_{2}D^{+} (OPR(H_{2}D^{+})), and gas temperature.
The OPR(H_{2}D^{+}) in the dense ISM is mainly determined by the following reactions (Gerlich et al. 2002; Hugo et al. 2009): Then, assuming steady state, the OPR(H_{2}D^{+}) is given as a function of OPR(H_{2}) and gas temperature: (A.11)From Eqs. (A.6) and (A.11), we can express k_{1b} as a function of OPR(H_{2}) and the gas temperature. In Fig. 1, we show threshold abundances of CO and electrons as functions of the OPR(H_{2}), below which Reaction (1) in the backward reaction is more efficient than the destruction by CO and electrons. Roughly speaking, the backward reaction rate exceeds the destruction rate by CO when OPR(H_{2}) ≳ 20n(CO) /n(H_{2}) at ≤20 K, while it exceeds the recombination rate with electrons when OPR(H_{2}) ≳ 3 × 10^{3}n(e) /n(H_{2}) at ≤20 K. Assuming the canonical HD abundance with respect to H_{2}, 3 × 10^{5}, the condition k_{1b}n(H_{2}) >k_{HD}n(HD) corresponds to OPR(H_{2}) ≳ 6 × 10^{4} at ≤20 K.
Appendix B: Analytical treatment of H_{2} orthotopara ratio
In this Appendix, we derive an analytical solution of OPR(H_{2}) when abundances of ionic species and the H_{2} formation rate on grain surfaces are given. Let us consider the following set of differential equations, which describe the temporal variations of the abundances of o  H_{2} and p  H_{2} in the gas and solid phases:
where F, W, D are the adsorption rate, desorption rate, and rate of H_{2} destruction via e.g., photodissociation, respectively. ⟨ A ⟩ in Eqs. (B.3) and (B.4) is the population of species A in the active surface ice layers, and thus ⟨ A ⟩ n_{gr} is the number density of species A in the surface ice layers per unit volume of gas. R_{H2} is the formation rate of H_{2} on grain surfaces, while b_{o} is the branching ratio to form o  H_{2}. In Eqs. (B.1) and (B.2), we neglect H_{2} formation in the gas phase for simplicity. In Eqs. (B.3) and (B.4), we do not consider nuclear spin conversion on the surface, though it is straightforward to include the spin conversion on the surface in the following analysis.
Combining Eqs. (B.1)−(B.4), we get (B.5)where we used n(o  H_{2}) + n(p  H_{2}) = n(H_{2}) and n(H_{2}) ≫ ⟨ H_{2} ⟩ n_{gr}. The latter should be valid in molecular gases; the binding energy of H_{2} on a H_{2} substrate is only 23 K, corresponding to an sublimation timescale of ~10^{10} s at T = 10 K (Vidali et al. 1991; Cuppen & Herbst 2007), while the adsorption timescale is ~10^{9}/n_{H} yr. The H_{2} formation timescale, τ_{H2}, was defined as n(H_{2}) /R_{H2}. We also assumed that the rate coefficients of reactions to destroy o  H_{2} and p  H_{2} are the same, i.e., D(o  H_{2}) /D(p  H_{2}) = n(o  H_{2}) /n(p  H_{2}). Equation (B.5) describes the time evolution of the OPR(H_{2}). The terms τ_{o → p}, τ_{p → o}, and τ_{H2} are timedependent in general. It is straightforward to solve Eq. (B.5) in the astrochemical simulations without nuclear spin state chemistry, and one may obtain a reasonable approximation of the temporal variation of the OPR(H_{2}). In order to solve Eq. (B.5) analytically, we consider τ_{o → p}, τ_{p → o}, and τ_{H2} as constant. This assumption is valid when the orthotopara spin conversion time scale is longer than the chemical (formation and destruction) timescale of hydrogen and light ions. The solution is given as follows: where τ_{opr} gives the characteristic timescale of OPR(H_{2}) evolution, and OPR_{0}(H_{2}) is the initial OPR(H_{2}). The steadystate value of OPR(H_{2}) (OPR_{st}(H_{2}) ≡ OPR(H_{2})(t → ∞)) is given in Eq. (14), which was derived by Le Bourlot (1991) in a different manner.
When OPR_{st}(H_{2}) ≪ OPR_{0}(H_{2}) ≪ 1, Eq. (B.6) can be simply rewritten as OPR(H_{2})(t) ≈ OPR_{st}(H_{2}) + OPR_{0}(H_{2})exp(−t/τ_{opr}). Then, it takes a greater time than ln [ OPR_{0}(H_{2}) / OPR_{st}(H_{2}) ] τ_{opr} to reach the steady state value.
Appendix C: H_{2} orthotopara ratio when the conversion of hydrogen into H_{2} is almost complete
Evolution of molecular abundances have often been investigated via pseudotime dependent models in which hydrogen is assumed to be in H_{2} at the beginning of the calculation. In such models, the molecular D/H ratio depends on the initial OPR(H_{2}), which is treated as a free parameter (e.g., Flower et al. 2006). In this appendix we analytically derive the OPR(H_{2}) when the conversion of hydrogen into H_{2} is almost complete.
During the H/H_{2} transition, orthopara spin conversion occurs through Reaction (17). H^{+} is primarily formed via cosmicray/Xray ionization of atomic hydrogen, while it is mainly destroyed via recombination with electrons, and charge transfer to other species, such as atomic deuterium and oxygen (e.g., Dalgarno et al. 1973). The former is valid when n(H) /n(H_{2}) >ξ_{H}/ (b_{H+}ξ_{H2}) ≈ 0.05, where ξ_{H} and ξ_{H2} are the ionization rates of atomic hydrogen and H_{2}, respectively. b_{H+} is the branching ratio to form H^{+} for the ionization of H_{2}. At steady state, the number density of H^{+} can be given as (C.1)where (C.2)where k_{(X + Y)} is the rate coefficient of reaction X + Y. We assumed that the number density of electrons is equal to that of carbon ions. The carbon ion is the dominant form of carbon in the H/H_{2} transition region, while oxygen and deuterium are predominantly in atomic form. With Eq. (C.1), we can evaluate β_{1} and β_{2} (Eqs. (15) and (16)) by the following: where x(i) is the abundance of species i with respect to hydrogen nuclei, n_{H} is the number density of hydrogen nuclei, v_{th} is the thermal velocity of atomic hydrogen, and k_{17f} and k_{17b} are the rate coefficients of Reaction (17) in the forward direction and the backward direction, respectively. We used ξ_{H2} = 2.2ξ_{H}. We assumed that the H_{2} formation rate is half of the accretion rate of atomic hydrogen onto dust grains. By substituting Eqs. (C.3) and (C.4) into Eq. (14), we can get the OPR(H_{2}) in the H/H_{2} transition regime under the steady state assumption as a function of ξ_{H2}/n_{H} and gas temperature. In Fig. C.1, we show the OPR(H_{2}) when the the conversion of hydrogen into H_{2} is (almost) complete, i.e., x(H_{2}) = 0.5. Above the thick dashed line, where τ_{o → p} is smaller than τ_{H2}, the steadystate assumption is justified.
It is clear that the OPR(H_{2}) is greater than unity only when H_{2} formation occurs at low ionization (ξ_{H2}/n_{H}< 10^{22} cm^{3} s^{1}) and/or warm conditions. In our fiducial model, this occurs at T_{gas}> 50 K, but the exact value depends on ξ_{H2}/n_{H}. In those regions with such warm gas temperatures, the dust temperature would also be warm. At T_{dust} ≳ 20 K, the H_{2} formation rate may drop considerably, depending on characteristics of chemisorption sites on grain surfaces (e.g., Hollenbach & Salpeter 1971; Cazaux & Tielens 2004; Iqbal et al. 2014), which is not considered here.
Fig. C.1
OPR(H_{2}) when the conversion of hydrogen into H_{2} is almost complete under the steadystate assumption. Above the thick dashed line, where τ_{o → p} is smaller than τ_{H2}, the steadystate assumption is justified. Above the thin dashed line, where β_{1} is larger than β_{2}, the OPR(H_{2}) is similar to the lowtemperature thermalized value of 9exp(−170.5 /T_{gas}). See Appendix C. 

Open with DEXTER 
Appendix D: Scaling relation of the HDO/H_{2}O ratio in the surface ice layers with the atomic D/H ratio
Here we derive the scaling relation of the HDO / H_{2}O ratio in the chemically active surface ice layers with the atomic D/H ratio, which is applicable under the UV irradiation conditions where photodesorption regulates the growth of ice mantles. Figure 9 shows the important reactions for H_{2}O ice and HDO ice in our models, and they are considered in the following analysis. Let us consider the following set of differential equations which describe temporal variations of the abundances of H_{2}O, HDO, OH and OD in the active surface ice layers: where k_{ph} and k_{phdes} are the rate coefficients of photodissociation and photodesorption in the active surface ice layers, respectively, of H_{2}O ice and HDO ice. b_{OH} is the branching ratio to form OH + D for HDO ice photodissociation (~0.3, Koning et al. 2013). Combining Eqs. (D.1)−(D.4), the evolutionary equation of the HDO / H_{2}O ratio in the active ice layers (f_{HDO, s} = ⟨ HDO ⟩ / ⟨ H_{2}O ⟩) can be written as follows: (D.5)where we used the inequalities ⟨ OH ⟩ ≪ ⟨ H_{2}O ⟩, ⟨ OD ⟩ ≪ ⟨ HDO ⟩, and f_{HDO, s} ≪ 1, which are verified by our numerical simulations. We defined f_{D, s} = ⟨ D ⟩ / ⟨ H ⟩. We also used the relations, δ_{hop} ≡ k_{hop}(D) /k_{hop}(H) ~ k_{O + D}/k_{O + H} ~ k_{OH + D}/k_{OH + H}, which should be valid because of the much higher hopping rates of atomic H and D compared to those of atomic O and OH radical. The energy barrier difference against hopping between atomic D and H is ~10 K on amorphous water ice (Hama et al. 2012). When the right hand side of Eq. (D.5) is less (more) than zero at given f_{D, s}, the ratio f_{HDO, s} decreases (increases) on a timescale shorter than b_{OH}k_{ph}. Since the photodissociation timescale in gas with low extinction is much shorter than the dynamical timescale (>1 Myr in the case of our cloud formation model), let us assume steady state, which corresponds to the minimum (or maximum) of f_{HDO, s}/f_{D, s}. Then we obtain (D.6)To produce water ice mantles under UV irradiation, the rate of OH formation through the hydrogenation of atomic oxygen (i.e., k_{O + H} ⟨ O ⟩ ⟨ H ⟩) should be larger than the photodesorption rate of H_{2}O (see Fig. 9). When the dynamical timescale is longer than the timescale of the OH formation, the chemistry evolves as the OH formation rate and the H_{2}O photodesorption rate are almost balanced: (D.7)This relation is confirmed in our simulations. The population of OH can be evaluated from the following equation: (D.8)Then (D.9)where we used k_{ph} ⟨ H_{2}O ⟩ ≫ k_{phdes} ⟨ H_{2}O ⟩ ≈ k_{O + H} ⟨ O ⟩ ⟨ H ⟩. From Eqs. (D.6), (D.7), and (D.9), we get the scaling relation, (D.10)where P_{phdes} is k_{phdes}/k_{ph} ~ 0.02 (Andersson et al. 2006; Arasa et al. 2015). The term δ_{hop}f_{D, s}/ (1 + Γ) corresponds to p_{OH → HDO} which is discussed in the main text. In the case with Γ ≪ 1 (or, in other words, OH + H is the dominant formation pathway of H_{2}O ice), we get f_{HDO, s} ≈ f_{D, s}. In another extreme case Γ ≫ 1, we get f_{HDO, s} ≈ (Γ^{1} + 0.02)f_{D, s}, i.e., f_{HDO, s} is much smaller than f_{D, s}. For example, Γ is ~10^{4} at A_{V} = 1 mag in our fiducial simulation, and it further increases with the increase of A_{V}. The minimum of the f_{HDO, s}/f_{D, s} ratio is ~0.02 in our fiducial simulation, which is well reproduced by Eq. (D.10).
Appendix E: Comparisons with observations of sulfurbearing species
As discussed in Sect. 5.1, the evolution of the OPR(H_{2}) depends on the assumed heavy metal abundances. The goal of this appendix is to determine which case, the LM abundances (fiducial case) or the HM abundances, gives the better predictions on the OPR(H_{2}) evolution in the ISM.
There are some estimates of OPR(H_{2}) in cold dense clouds/cores from observations of molecules other than H_{2} in the literature. We do not use them for the constraint here, because how to estimate OPR(H_{2}) from the observations is not wellestablished, and because the evolution of OPR(H_{2}) in the ISM can be in the nonequilibrium manner and may vary among sources depending on their past physical evolution. Instead, we compare observations of sulfurbearing molecules toward diffuse/translucent/molecular clouds in the literature with our model results. We focus on the total abundance of selected sulfurbearing molecules in the gas phase (CS, SO, and H_{2}S) and the HCS^{+}/CS abundance ratio. The former can probe the total abundance of elemental sulfur in the gas phase, though the main form of gaseous sulfur is likely to be the neutral atom S or S^{+} especially in gas with low extinction. The latter, the HCS^{+}/CS ratio, can probe the ionization degree of the gas, which is related to the S^{+} abundance (and abundances of the other heavy metal ions). The HCS^{+}/CS ratio is anticorrelated with the electron abundance in our models, because CS is formed by dissociative electron recombination of HCS^{+} and destroyed by photodissociation. Although the HCS^{+}/CS ratio also depends on the photodissociation timescale of CS, it is common between the fiducial model and model HM at given A_{V}. The rate coefficient and the branching ratios for the electron recombination of HCS^{+} are taken from Montaigne et al. (2005).
Turner (1996) derived the molecular abundances of H_{2}S, CS, and SO in translucent clouds with line of sight visual extinction of up to 5 mag. He found that the total abundance of these three species with respect to hydrogen nuclei is typically 10^{8}–10^{7}. In the dense molecular clouds, TMC1 and L134N, with higher line of sight visual extinction than the samples in Turner (1996), the total abundance is 10^{9}–10^{8} (Ohishi et al. 1992; Dickens et al. 2000). The lower total abundance in the dense molecular clouds implies that depletion of gaseous sulfur is significant in the dense clouds (Joseph et al. 1986). It is not obvious which A_{V} in our model can be compared with the observations towards the translucent and dense molecular clouds. Towards TMC1(CP) and L134N, there is evidence of CO freezeout; the gaseous CO abundance with respect to hydrogen nuclei (4 × 10^{5}) is less than the canonical value of 10^{4} and infrared absorption by CO ice is detected in the line of sight towards the vicinity of them (Whittet 2007, 2013). We assume that our results at A_{V}> 2 mag, where the gaseous CO abundance decreases to less than ~4 × 10^{5} because of the freezeout, can be compared with the observations in TMC1 and L134N, while the results at A_{V}< 2 mag are compared with the observations in clouds with lower line of site extinction.
Figure 11 shows the total abundance of H_{2}S, CS, and SO as a function of visual extinction in our models. Model HM better reproduces the total abundance of the Sbearing molecules observed in the translucent clouds, while the fiducial model better reproduces the total abundance observed in the dense molecular clouds. This result implies that again, the depletion of sulfur from the gas phase occurs during the formation and evolution of molecular clouds, and that model HM underestimates the degree of the sulfur depletion.
The HCS^{+}/CS ratio has been derived toward diffuse clouds (0.08, Lucas & Liszt 2002), clouds with line of sight visual extinction of up to 5 mag (0.01−0.1, Turner 1996), and the dense molecular clouds, TMC1 and L134N (0.06, Ohishi et al. 1992). Roughly speaking, the ratio is in the range of 0.01−0.1 in the diffuse and dense ISM. Figure 12 shows the HCS^{+}/CS ratio and the electron abundance in the fiducial model and model HM. The fiducial model reproduces the HCS^{+}/CS ratio much better than model HM, though the both models tend to underestimate the HCS^{+}/CS ratio compared with observations, i.e., overestimate the ionization degree of the gas.
Considering the HCS^{+}/CS ratio is reproduced by the fiducial model much better than by model HM, we conclude that the fiducial model gives a better prediction for the Sbearing species and thus for the OPR(H_{2}). In particular, model HM seems to overpredict the Sbearing species in the gas phase. The nonsuccess of the model HM probably means that a nonnegligible fraction of sulfur is incorporated into dust grains (Keller et al. 2002) and/or the nonthermal desorption rates of icy sulfur are overestimated in the current model; the partitioning of elemental sulfur between gas and ice in our model depends on the assumed nonthermal desorption rates, especially chemisorption probabilities, which remain uncertain at the current stage. In our models, the dominant sulfur reservoirs in ices are HS and H_{2}S, though there has been no detection of H_{2}S ice in the ISM. HS ice is formed from H_{2}S ice via the hydrogen abstraction reaction, H_{2}S + H, with an activation energy barrier of 860 K (Hasegawa et al. 1992). HS ice is hydrogenated to form H_{2}S ice again, desorbing ~1% of the product H_{2}S through chemisorption. This loop is the main mechanism of the desorption of icy sulfur in our model as in Garrod et al. (2007). If we set the chemisorption probability for the hydrogenation of HS ice to be 0.1%, only 3% of elemental sulfur remains in the gas phase at A_{V} = 3 mag in model HM, though gaseous sulfur is still more abundant than in the fiducial model.
© ESO, 2015