Low temperature HD+ortho-/para-H₂ inelastic scattering of astrophysical interest

Elastic and inelastic collisions between atoms and molecules and/or between molecules and molecules are of great scientific interest in the fields of physical chemistry and chemical physics. The reason is that such processes can provide valuable information about fundamental interactions between different chemical species, their chemical properties, their energy transfer quantum dynamics and many other properties. The pioneering studies focussed on light atoms and molecules, because of their simple nature. For three- and four-atomic systems with a small number of electrons, the potential energy surface (PES) can be computed with relatively high accuracy [1–6]. Consequently, for these systems one can then more easily test different dynamical methods, such as classical, semi-classical, quasi-classical trajectory, and pure quantum-mechanical computational formulations, and compare the results with available experimental data in a controlled fashion. The test methods devised could then be applied to more complex many-atomic systems, wherein a controlled comparison is not possible. Among these small systems the four-atomic H₂+H₂ and H₂+HD scattering processes have attracted significant attention not only in chemical physics but also in astrophysics. In astrophysics H₂ and HD play an important role, because of their abundance in the molecular cloud of the universe [7, 8]. Together with the H₂+H₂ collision, the HD+H₂ collision is also of significant importance in the astrophysics of the early universe. Specifically it is important in the modeling of pre-galactic clouds and planetary atmospheres; in the cooling of primordial gas and in the formation of stars [9–15].

In [16] a semiclassical treatment of H₂-H₂ scattering is formulated. In [17] the author developed and applied a rigid rotor model to study rotational excitation in H₂-H₂ by applying quantum close-coupling scattering calculations. In this approach the distance between the hydrogen atoms in both H₂ molecules was fixed at a constant value based on the average. This model was applied to many different atomic and molecular systems; see, for example [13, 18–20]. The main goal of the work [17] was to compute rotational thermal rate coefficients in the H₂+H₂ system at low temperatures of astrophysical interest. Quantum-mechanical close-coupling calculations for three-dimensional collisions of para-H₂ and ortho-H₂ with HD are performed in [11, 13, 18, 19, 21], where the HD-H₂ potential is derived from the H₂-H₂ potential. A quantum dynamical study of H₂-H₂ collisions is reported for both ortho- and para-H₂ in [22–24]. In [25], the authors considered H₂+H₂ and D₂+D₂ rotational inelastic scatterings with the use of the H₂-H₂ potential energy surface (PES) from [1]. In [22] a full six-dimensional scattering computation has been performed taking into account vibrational relaxation in the H₂+H₂ collision. In that study the H₂-H₂ PES from [2] was used and its anisotropy properties have been studied at low temperatures: 20 K ≲ T ≲ 300 K. HD+HD scattering has been studied experimentally in [26] and theoretically in [27] for a wide range of collision energies. A comprehensive computational and experimental study of total cross section in H₂-H₂, D₂-D₂, and HD-HD scattering for both ortho and para H₂ and D₂ has been reported in [28]. Measurements of energy transfer rates in HD+HD [29] and H₂+H₂ [30, 31] collisions have also been performed.

However, realistic theoretical investigations of the low-energy HD+H₂ collision are lacking, although preliminary quantum calculation of this process has been reported in [11, 13]. Schaefer [11] calculated rate coefficients for the excitation of HD by H₂, for the low-lying rotational levels using a modified older potential for HD and H₂. Flower [13] performed an improved calculation of HD-H₂ scattering with Schwenke's H₂-H₂ PES using a larger rotational basis set [3]. In this paper we report an improved calculation of this problem using a realistic HD-H₂ PES derived from Hinde's recent H₂-H₂ PES [5]. In two papers [18, 19] the PES from work [2] has been applied together with a pure quantum-mechanical dynamical approach. The surprising thing is, that the results of works [18, 19] are closer to the results of older work [11] than to [13]. Therefore, there is a need to carry out new computations with newer PES between HD and H₂.

While the exchange symmetry is broken in HD+H₂ it still possesses many similarities with the H₂+H₂ system. The PESs of H₂-H₂ and HD-H₂ should basically be the same six-dimensional function. This fact follows from the general Born–Oppenheimer approach [32]. At the same time the two collisions H₂+H₂ and H₂+HD should have rather different scattering outputs. This is because the H₂ and HD molecules have fairly different rotational constants, internal symmetries and as a result a different rotational-vibrational spectrum. The HD-H₂ PES can be derived from the H₂-H₂ PES by shifting the center of mass (c.m.) of the H₂ molecule to the c.m. of the HD molecule. Once the exchange symmetry is broken in H₂-H₂ by replacing the H with the D atom in one of the H₂'s then one has the new HD-H₂ PES. In this fashion, we constructed the HD-H₂ PES from the H₂-H₂ PES of Hinde [5] employing all parts of the full HD-H₂ interaction including the HD's dipole moment. Using this HD-H₂ PES we carried out pure quantum-mechanical calculations for inelastic collisions of rotationally excited HD and H₂ molecules, i.e. the process:

$$\begin{eqnarray}&&{\rm{HD}}({j}_{1})+{{\rm{H}}}_{2}({j}_{2})\to {\rm{HD}}({j}_{1}^{\prime })+{{\rm{H}}}_{2}({j}_{2}^{\prime }).\end{eqnarray} \tag{ 1 }$$

The scattering cross sections and their corresponding thermal rate coefficients are computed using a non-reactive quantum-mechanical close-coupling approach. The four-atomic system is shown in figure 1.

Zoom In Zoom Out Reset image size

In section 2 we briefly describe the quantum-mechanical approach used in this paper. Section 3 includes the computational results. We compare the cross-sections and rates with those of other authors [11, 13], and our previous calculations [18, 19], where a different HD-H₂ PES derived from the well-known Boothroyd–Martin–Keogh–Peterson (BMKP) H₂-H₂ PES [2], was used. Discussion and conclusions are provided in section 4. The corresponding procedure to obtain a modified HD-H₂ PES from the existing H₂-H₂ surface [5] is presented in the appendix A. Atomic units $(e={m}_{e}={\hslash }=1)$ are used throughout this paper.

In this section we provide a brief account of the present quantum-mechanical close-coupling approach following the method in [17]. The HD and H₂ molecules are treated as linear rigid rotors. The model has been applied in few previous works [13, 17]. In all our calculations with this potential the bond length was fixed at 1.449 a.u. or 0.7668 Å for the H₂ molecule and 1.442 a.u. for HD which is 0.7631 Å. The Schrödinger equation for the (12)+(34) collision in the c.m. frame, where (12) and (34) are diatomic molecules formed by atoms 1–4 is [17]

$$\begin{eqnarray}&&\left[\displaystyle \frac{{\hat{P}}_{{\vec{R}}_{3}}^{2}}{2{{\mathcal{M}}}_{12}}+\displaystyle \frac{{\hat{L}}_{{\hat{R}}_{1}}^{2}}{2{\mu }_{1}{R}_{1}^{2}}+\displaystyle \frac{{\hat{L}}_{{\hat{R}}_{2}}^{2}}{2{\mu }_{2}{R}_{2}^{2}}+V({\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3})-E\right]\\ &&\quad \times {\rm{\Psi }}({\hat{R}}_{1},{\hat{R}}_{2},{\vec{R}}_{3})=0.\end{eqnarray} \tag{ 2 }$$

Here ${\hat{P}}_{{\vec{R}}_{3}}$ is the momentum operator of the kinetic energy of the collision, ${\vec{R}}_{3}$ is the collision coordinate, whereas ${\vec{R}}_{1}$ and ${\vec{R}}_{2}$ are relative vectors between atoms in the two diatomic molecules as shown in figure 1, and ${\hat{L}}_{{\hat{R}}_{1(2)}}$ are the quantum-mechanical rotation operators of the rigid rotors, ${\mu }_{1}$ and ${\mu }_{2}$ are the reduced masses of the HD and H₂ molecules and ${{\mathcal{M}}}_{12}$ is the reduced mass of the two molecules. The vectors ${\hat{R}}_{1(2)}$ are the angles of orientation for rotors $(12)$ and $(34),$ respectively; $V({\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3})$ is the PES of the four-atomic system $(1234),$ and E is the total energy in the (c.m.) system. As in our previous papers [18–20, 27], in the current work the VRTP mechanism is used. It allows us to specify $V({\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3})$ explicitly rather than to expand the potential [33]. Therefore, one can say that in these calculations we used all parts of the modified HD+H₂ PES, although the van der Waals long-range part of it would be more important at low energy collisions. The use and modification of the original H₂-H₂ PESs $V({\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3})$ is discussed in appendix A.

The cross sections $\sigma ({j}_{1}^{\prime },{j}_{2}^{\prime };{j}_{1}{j}_{2},\varepsilon )$ for rotational excitation and relaxation can be obtained from the S-matrix: ${S}_{\alpha {\alpha }^{\prime }}^{J}.$ Specifically, for the excitation from HD(j₁, m₁) + ${{\rm{H}}}_{2}({j}_{2},{m}_{2})$ to HD $({j}_{1}^{\prime }{m}_{1}^{\prime })$ + ${{\rm{H}}}_{2}({j}_{2}^{\prime }{m}_{2}^{\prime })$ they are summed over the final angular momentum projections (${m}_{1}^{\prime }{m}_{2}^{\prime }$) and averaged over the initial projections (${m}_{1}{m}_{2}$) of the HD and H₂ molecules of angular momenta j₁ and j₂. The value of $\sigma ({j}_{1}^{\prime },{j}_{2}^{\prime };{j}_{1}{j}_{2},\varepsilon )$ can be given by the following expression:

$$\begin{eqnarray}\sigma ({j}_{1}^{\prime },{j}_{2}^{\prime };{j}_{1}{j}_{2},\varepsilon ) & = & \displaystyle \frac{\pi }{(2{j}_{1}+1)(2{j}_{2}+1){k}_{\alpha {\alpha }^{\prime }}}\\ & & \times \displaystyle \sum _{{{Jj}}_{12}{j}_{12}^{\prime }{{LL}}^{\prime }}(2J+1){| {\delta }_{\alpha {\alpha }^{\prime }}-{S}_{\alpha {\alpha }^{\prime }}^{J}(E)| }^{2}.\end{eqnarray} \tag{ 3 }$$

The kinetic energy is $\varepsilon =E-{B}_{1}{j}_{1}({j}_{1}+1)-{B}_{2}{j}_{2}({j}_{2}+1),$ where ${B}_{1}=44.7\ \ {{\rm{cm}}}^{-1}$ and ${B}_{2}=60.8\ \ {{\rm{cm}}}^{-1}$ are the rotation constants of rigid rotors $(12)$ and $(34)$ respectively; they are shown in figure 1. Next, E is the total energy of the system, J is the total angular momenta of the four-atomic system, $\alpha \equiv ({j}_{1}{j}_{2}{j}_{12}L),$ where ${j}_{1}+{j}_{2}={j}_{12}$ and ${j}_{12}+L=J,$ ${k}_{\alpha {\alpha }^{\prime }}=\sqrt{2{{\mathcal{M}}}_{12}(E+{E}_{\alpha }-{E}_{{\alpha }^{\prime }})}$ is the channel wavenumber and ${E}_{\alpha ({\alpha }^{\prime })}$ are rotational channel energies. Finally, the relationship between the rotational thermal-rate coefficient ${k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T)$ at temperature T and the corresponding cross section ${\sigma }_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(\varepsilon ),$ can be obtained through the following weighted average formula:

$$\begin{eqnarray}{k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T) & = & \sqrt{\displaystyle \frac{8{k}_{{\rm{B}}}T}{\pi {{\mathcal{M}}}_{12}}}\displaystyle \frac{1}{{({k}_{{\rm{B}}}T)}^{2}}\\ & & \times {\displaystyle \int }_{{\varepsilon }_{s}}^{\infty }{\sigma }_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(\varepsilon ){{\rm{e}}}^{-\varepsilon /{k}_{{\rm{B}}}T}\varepsilon {\rm{d}}\varepsilon ,\end{eqnarray} \tag{ 4 }$$

where k_B is the Boltzmann constant and $\varepsilon =E-{E}_{{j}_{1}}-{E}_{{j}_{2}}$ is the pre-collisional translational energy at the translational temperature T, the and ${\varepsilon }_{s}$ is the minimum kinetic energy for the levels j₁ and j₂ to become accessible.

We used the MOLSCAT program [33] to solve the Schrödinger equation (2). Convergence was obtained for the integral cross sections, $\sigma ({j}_{1}^{\prime },{j}_{2}^{\prime };{j}_{1}{j}_{2},\varepsilon ),$ with respect to the variation of the variables utilized in all considered collisions at different collision energies. For the intermolecular distance R₃ we used from ${R}_{3\mathrm{min}}=3.0$ a.u. to ${R}_{3\mathrm{max}}=30.0$ a.u. We also applied a few different propagators included in the MOLSCAT computer program, and our calculations show that D. Manolopoulos's hybrid modified log-derivative propagator technique [34] would be quite numerically stable and a time effective approach. This method is used in the majority of the calculations.

The maximum value of the total angular momentum J was set at 44 while the number of levels ${N}_{\mathrm{lvl}}$ in the basis set of HD + H₂ was set at 42. Specifically, in the case of HD $({j}_{1})$+para-${{\rm{H}}}_{2}({j}_{2})$, j₁ it has values 0, 1, 2, and 3 and j₂ has values 0, 2, and 4. This combination generates the total number of included levels ${N}_{\mathrm{lvl}}$ = 34. In the case of HD+ortho-H₂ the parameter j₁ has values 0, 1, 2, and 3 and j₂ has values 1, 3, and 5. This results in the total number of levels ${N}_{\mathrm{lvl}}$ = 42. A number of test computations with higher values for the j₁ and j₂ parameters have been carried out. For example, in the case of HD+para-H₂, j₁ was taken as 0, 1, 2, 3, 4 and j₂ as 0, 2, 4, 6. This j₁/j₂ combination produced ${N}_{\mathrm{lvl}}=74.$ We obtained similar results in both cases, confirming the convergence of the calculations.

Because the HD+H₂ total rotational energy transfer cross sections (3) have shape resonances at low energies a large number of energy points are needed in order to effectively reproduce them. We used up to 250 energy points in each computation for each specific rotational transition in the HD and H₂ molecules considered. More space discretization points were used at low collision energies and fewer points in the higher-energy sector.

Below we present results for cross sections and thermal rate coefficients for different quantum-state transitions in HD and H₂ molecules. We reproduced shape resonances in the low velocity region, which are very important in the cooling of the astrophysical media at low temperatures. We compared them with the older quantum-dynamical results of Schaefer [11] and Flower [13] . We also present results [18, 20] using a PES obtained from a modification of the BMKP H₂-H₂ PES [2], viz. Appendix A. In figures 2–4 we show four different results for the integral cross sections in the collisions:

$$\begin{eqnarray}&&{\rm{HD}}(1)+{{\rm{H}}}_{2}(0)\to {\rm{HD}}(0)+{{\rm{H}}}_{2}(0),\end{eqnarray} \tag{ 5 }$$

$$\begin{eqnarray}&&{\rm{HD}}(2)+{{\rm{H}}}_{2}(0)\to {\rm{HD}}(1)+{{\rm{H}}}_{2}(0),\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{\rm{HD}}(2)+{{\rm{H}}}_{2}(0)\to {\rm{HD}}(0)+{{\rm{H}}}_{2}(0),\end{eqnarray} \tag{ 7 }$$

respectively. One can see that the new cross sections obtained with the present PES have the same structure and shape, but also have substantially larger (∼60%) values at medium energies ($v\gt 200\;{{\rm{ms}}}^{-1}$) when compared with the results obtained from the modified BMKP PES [2] and the older result from Schaefer's calculations [11]. Processes (6) and (7) are also important collisions from the astrophysical point of view. They represent transitions from the ${j}_{1}=2$ and ${j}_{2}=0$ state of the HD $({j}_{1})$+${{\rm{H}}}_{2}({j}_{2})$ system. It is seen from table 1 that the corresponding rotational energy of the system is 268.2 cm⁻¹ and rotational relaxation processes to states of lower energy should be relevant. The cross sections from [13] are available only for higher collision velocities, i.e. $v\geqslant$ 300 m s⁻¹. The present cross sections have the same behavior as previous results, but they have larger values than the old ones, specially at higher energies. Next in figure 5 we exhibit the thermal rate coefficients corresponding to the processes considered above. Here, we additionally include the older results from [13], where Schwenke's H₂-H₂ PES [3] was used. One can see that the new rates are substantially larger than those of other calculations. In general, it means that the contribution of the HD+H₂ collision to the HD-cooling function can have even larger contributions than previously expected.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 1.  Total rotational energy ${{\mathcal{E}}}_{\mathrm{rot}}$ in the four-atomic system: HD(j₁)+p-${{\rm{H}}}_{2}({j}_{2})$ and HD(j₁)+o-${{\rm{H}}}_{2}({j}_{2}).$ Here, ${{\mathcal{E}}}_{\mathrm{rot}}={B}_{1}{j}_{1}({j}_{1}+1)$ $+{B}_{2}{j}_{2}({j}_{2}+1),$ where ${B}_{1(2)}$ are the rotational constants of rigid rotors (12) and (34) respectively.

${{\mathcal{E}}}_{\mathrm{rot}}$ (cm−1)	HD(j1)	p-${{\rm{H}}}_{2}({j}_{2})$	${{\mathcal{E}}}_{\mathrm{rot}}$ (cm−1)	HD(j1)	o-${{\rm{H}}}_{2}({j}_{2})$
0.0	0	0	121.6	0	1
89.4	1	0	211.0	1	1
268.2	2	0	389.8	2	1
364.8	0	2	658.0	3	1
454.2	1	2	729.6	0	3
536.4	3	0	819.0	1	3
633.0	2	2	997.8	2	3

Figures 2–4 show large shape resonances in the rotational energy transfer cross sections $\sigma (v).$ It is seen that these resonances occur at very low velocity collisions between HD and H₂, i.e. $v\;\sim$ 80m s⁻¹. A possible physical explanation of this cross section behavior is the following. It is known that in the four-atomic H₂-H₂ system there are quasi-bound states. The not-too-deep levels are supported by the long range part of the H₂-H₂ surface's van der Waals well [2]. In the case of the current HD-H₂ system the quasi-bound levels also exist. However, their values (positions) would be slightly different (shifted), because of the differences between the H₂ and HD reduced masses and due to the broken symmetry in the HD-H₂ PES. As a consequence, HD+H₂ collision represents a non-reactive scattering process when in the discreet spectrum of the HD-H₂ system there are weekly bound levels ${\varepsilon }^{\prime }.$ Such a situation was already discussed in [35]. It was shown that when the collision energy is close to a quasi-bound level value, i.e. $\varepsilon \approx {\varepsilon }^{\prime },$ the scattering cross section σ is significantly increased [35]. This general result was obtained in our calculations too.

Now we consider the inelastic cross sections in HD+H₂ for higher rotational energies. For the systems HD(0)+${{\rm{H}}}_{2}(2)$ and HD(1)+${{\rm{H}}}_{2}(2),$ one can see from table 1 that the corresponding rotational energies are 364.8 cm⁻¹ and 454.2 cm⁻¹. In figure 6 (upper plot) we show the integral cross sections for the process:

$$\begin{eqnarray}&&{\rm{HD}}(0)+{{\rm{H}}}_{2}(2)\to {\rm{HD}}(2)+{{\rm{H}}}_{2}(0).\end{eqnarray} \tag{ 8 }$$

While we obtain a relatively good qualitative agreement between cross sections calculated with the older modified BMKP PES and Hinde's PES, there is a dramatic difference with the corresponding result from [11]. The present cross sections are larger than those of [11] by a few orders of magnitude. In figure 6 (lower plot) we show the results for the corresponding thermal rate coefficient ${k}_{20\to 02}(T).$ It is seen that the present rates and those obtained by Flower [13] have a flat temperature dependence, whereas the rate of Schaefer [11], although smaller than the other results, increases monotonically with energy. Because the thermal rate of reaction (8) is relatively large, one can conclude that this channel can also make a substantial contribution to the astrophysical HD-cooling function.

Zoom In Zoom Out Reset image size

In figures 7–10 we show results for the total cross section in the inelastic scattering from the state ${\rm{HD}}(1)+{{\rm{H}}}_{2}(2)$ with rotational energy 454.2 cm${}^{-1},$ viz. Table 1. All de-excitation processes have been computed for this state, namely:

$$\begin{eqnarray}&&{\rm{HD}}(1)+{{\rm{H}}}_{2}(2)\to {\rm{HD}}(0)+{{\rm{H}}}_{2}(2),\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{\rm{HD}}(1)+{{\rm{H}}}_{2}(2)\to {\rm{HD}}(2)+{{\rm{H}}}_{2}(0),\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&{\rm{HD}}(1)+{{\rm{H}}}_{2}(2)\to {\rm{HD}}(1)+{{\rm{H}}}_{2}(0),\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{\rm{HD}}(1)+{{\rm{H}}}_{2}(2)\to {\rm{HD}}(0)+{{\rm{H}}}_{2}(0).\end{eqnarray} \tag{ 12 }$$

The present cross sections exhibit again fairly good qualitative agreement with the results computed with the modified BMKP PES [2], and also with the older results from [11], i.e. the general shape and trend of the behavior of these cross sections are the same in all cases. Also, there is a relatively small bump in the cross sections of the processes (10)–(12) at collision velocity ∼1100 m s⁻¹ which is also reproduced by two PESs. One can see that the process (9) can make a significant contribution to the total HD-cooling function because its cross section is rather large relative to the other channels. Finally, in order to carry out new computations of the astrophysical cooling function, table 2 includes the relevant thermal rate coefficients for the HD + para-H₂ case in the temperature range from 2 K to 300 K. Next, as a test, we choose the initial state HD(2)+H₂(2) with a relatively higher total rotational energy: 633 cm⁻¹. The integral de-excitation cross sections from this state to different lower energy rotational states are shown in figure 11. One can see a fairly good agreement in the shape of the curves between various rotational transition results for the integral cross sections computed with the two different PESs.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In figure 12 we show three different rotational transition cross-sections for the ortho-hydrogen case. Here we chose low lying rotational levels of the two molecules: HD(1) + H₂(1) and HD(2) + H₂(1). Some results from older works [11, 13, 19] are also presented in the figure together with results computed with the newer modified PES from [5]. The corresponding thermal rate coefficients are shown in figure 13 where the older results from [11, 13, 19] are also presented for comparison purposes. Figures 14 and 15 include results for thermal rate coefficients for the transitions from higher rotational states, e.g. HD(0) + H₂(3) and HD(1) + H₂(3).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Finally, in appendix B we present our new results for the thermal rate coefficients which can be used in subsequent computation of the astrophysical HD-cooling function. In table 2 we show thermal rates of different de-excitation processes in low-temperature HD + para-H₂ rotational energy transfer collisions and in table 3 the same data for the HD + ortho-H₂ case.

State-to-state close-coupling quantum-mechanical calculations for rotational de-excitation cross-sections and corresponding thermal rate coefficients of the HD+o-/p-H₂ collisions are presented using a linear rigid rotor model for the HD and H₂ molecules. The symmetrical H₂-H₂ PES of [5] has been appropriately adopted for the current non-symmetrical HD+H₂ system by appropriate translation and rotation. These geometrical operations lead to a new set of angle variables ${\theta }_{1}^{\prime },{\theta }_{2}^{\prime }$ and ${\varphi }_{2}^{\prime }$ for the Jacobi few-body coordinates, a new length of the intermolecular distance ${\vec{R}}_{3}^{\prime }$ and, as a result, to a new HD-H₂ PES. For comparison purposes in this paper we carried out a few calculations with the use of the older BMKP PES[2] which was also modified for HD-H₂. A test of convergence and the results for the cross-sections with the two PESs are obtained for a wide range of values of different parameters, such as the maximum value of the total angular momentum J, the intermolecular distance R₃ and others. The convergence issues were discussed in the beginning of the section 3 of this paper. It is seen from figures 2–4, that in most cases for the lower number quantum transition states the rotational energy transfer cross sections obtained with the use of the modified Hinde's PES have higher values than the cross sections computed with the use of the modified BMKP PES [18, 19]. The same situation occurs in comparisons within the older work [11]. In figure 5 the corresponding thermal rate coefficients are presented. In this case we include results from the paper [13] too. As can be seen the new thermal rates (solid lines) have substantially higher values than the other results. For example, in the rotational energy transfer, process (5), new thermal rates are ∼2 times higher than the other corresponding rates. Because (5) is considered as one of the main contributors to the resulting HD cooling process, one can say that the new HD-cooling function may have substantially greater values than the previous calculations [12]. Next, figure 6 represents results for a very interesting process (8). The interest in this channel lies in the fact that both molecules change their rotational quantum numbers by ${\rm{\Delta }}j$ = 2; therewith HD becomes excited and H₂ de-excited. In this case we obtained a significant deviation from the results of work [11]. However, the cross sections obtained with the modified versions of the Hinde and BMKP PESs are fairly close to each other and have fairly large values. Therefore, this process probably can make a contribution to the HD cooling process. In addition we would like to note that in work [12] the process has also been computed and discussed, and significant differences from [11] were also found.

Further, the following four figures 7–10 represent our integral cross sections and corresponding results from work[11] for the following de-excitation collisions: (9)–(12). In all of these processes the initial state of HD has the rotational state j₁ = 1 and H₂ has the rotational state j₂ = 2. The corresponding total rotational energy of the molecules can be found in table 1. The results in figure 7 are in fairly good agreement with each other and have relatively large values; therefore, this specific channel could also make a substantial contribution to the cooling function. It was found that in figures 8–10 both surfaces, namely, Hinde's and BMKP, both agree on a shape resonance at collision velocity $v\sim 1300\;{{\rm{ms}}}^{-1}$. Figure 11 shows the resulting de-excitation cross sections from the highly located rotational quantum level. Finally, our analysis in this paper would not be complete if we did not undertake computations for the ortho hydrogen case as well. Figure 12 shows results for the lower lying rotational quantum numbers of HD, namely, j₁ = 0 and 1, and H₂: j₂ = 1. It is seen that in this case the cross sections obtained with the Hinde PES are significantly larger than other results. The corresponding thermal rate coefficients are presented in figure 13. Again, as in previous para-hydrogen cases the new rates obtained with Hinde's potential are larger than previous results.

The total cooling function is represented as follows [12, 36, 37]:

$$\begin{eqnarray}&&{{\rm{\Lambda }}}_{{\rm{HD}}}(T)={n}_{{\rm{HD}}}{W}_{{\rm{HD}}}(T),\end{eqnarray} \tag{ 13 }$$

where ${n}_{{\rm{HD}}}$ is the number representing density of HD and ${W}_{{\rm{HD}}}$ is the HD cooling function in units of [erg×s⁻¹]. From work[36] it follows that the astrophysical line-cooling coefficient has the following expression:

$$\begin{eqnarray}&&L(u\to l;n,\tau )=\displaystyle \frac{{n}_{u}{E}_{u\to l}{\omega }_{{ul}}}{{{nn}}_{\nu }^{\prime }}{A}_{u\to l},\end{eqnarray} \tag{ 14 }$$

where ${E}_{{ul}}={E}_{u}-{E}_{l}$ is the energy of the radiated photon, ${\omega }_{{ul}}$ is the probability that the photon escapes absorption [36], E_u and E_l are the energies of the upper and lower states, n_u is the population of the upper level with ${n}_{u}=\mathrm{exp}({-E}_{u}/({k}_{{\rm{B}}}T)),$ ${A}_{u\to l}$ are Einstein's coefficients for the $u\to l$ transition, n and ${n}_{\nu }^{\prime }$ are the numbers of density of the participating species, for example H₂ and HD. For the low densities in the case of optically thin primordial gas one has to have ${\omega }_{{ul}}\approx$ 1. The line-cooling coefficient $L(u\to l;n,\tau )$ is chosen so that the value ${{nn}}_{\nu }^{\prime }L(u\to l;n,\tau )$ is the cooling rate per unit of volume due to the transition $u\to l$ in an atom or molecule ν at an optical depth τ [36]. Therefore, it follows from equations (13) and (14) that if one takes $u\equiv {j}_{1}{j}_{2}$ and $l\equiv {j}_{1}^{\prime }{j}_{2}^{\prime }$ the following formula for ${W}_{{\rm{HD}}}(T)$ can be derived:

$$\begin{eqnarray}&&{W}_{{\rm{HD}}}(T)=\displaystyle \frac{1}{{n}_{{\rm{HD}}}}\displaystyle \sum _{{j}_{1}{j}_{2},{j}_{1}^{\prime }{j}_{2}^{\prime }}{n}_{{j}_{1}{j}_{2}}({E}_{{j}_{1}{j}_{2}}-{E}_{{j}_{1}^{\prime }{j}_{2}^{\prime }}){A}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }},\end{eqnarray} \tag{ 15 }$$

where ${E}_{{j}_{1}{j}_{2}}\gt {E}_{{j}_{1}^{\prime }{j}_{2}^{\prime }}$ and the value of ${n}_{{\rm{HD}}}{W}_{{\rm{HD}}}(T)$ is equal to ${{nn}}_{\nu }^{\prime }L({ul};n,\tau ).$ Now, at the low density limit and taking into account the critical density concept, the total cooling function can be computed with the use of the following formula:

$$\begin{eqnarray}&&{{\rm{\Lambda }}}_{{\rm{HD}}}(T)=\displaystyle \sum _{{j}_{1}{j}_{2},{j}_{1}^{\prime }{j}_{2}^{\prime }}{n}_{{\rm{HD}}}({j}_{1}){n}_{{{\rm{H}}}_{2}}({j}_{2}){k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T)h{\nu }_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }},\end{eqnarray} \tag{ 16 }$$

which is in units of [erg×cm${}^{-3}\times$s⁻¹]. Here, $h{\nu }_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}$ is the emitted photon energy, ${k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T)$ is the thermal rate coefficient (4) corresponding to the rotational transitions ${j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }.$ Therefore, increasing the knowledge of rotational and possibly vibrational excitation and de-excitation rate constants, ${k}_{{jv}\to {j}^{\prime }{v}^{\prime }}(T),$ in atomic and molecular hydrogen-hydrogen collisions, such as HD/H₂+H₂, HD/H₂+H etc, is important in order to understand and be able to model the energy balance in the interstellar medium. For comparison purposes it would be very useful and interesting to carry out new computations of the rotational–vibrational integral cross sections and corresponding thermal rate coefficients for a low-energy HD+H collision. In this case a different H₃ PESs from papers [38, 39] could be applied.

This paper was supported by the Office of Research and Sponsored Programs of St. Cloud State University, USA and CNPq and FAPESP of Brazil.

A few important modifications to the Hinde H₂-H₂ PES [5] were needed for the current non-symmetrical four-atomic collision (1). The application and modification of the original H₂-H₂ BMKP PES were published in [18]. Below in this paragraph we briefly describe the procedure. To compute the distances between the four atoms the BMKP PES uses Cartesian coordinates. Consequently, it was necessary to convert spherical coordinates used in the close-coupling method to the corresponding Cartesian coordinates and compute the distances between the four atoms followed by calculation of the PES [18, 19]. This procedure used a specifically oriented coordinate system OXYZ. As a first step we needed to introduce the Jacobi coordinates $\{{\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3}\}$ and the radius-vectors of all four atoms in the space-fixed coordinate system OXYZ: $\{{\vec{r}}_{1},{\vec{r}}_{2},{\vec{r}}_{3},{\vec{r}}_{4}\}.$ Then the center of mass of the HD molecule has been relocated at the origin of the coordinate system OXYZ, and the ${\vec{R}}_{3}$ was directed to center of mass of the H₂ molecule along the OZ axis. Thus, one could obtain the following coordinate relationships: ${\vec{R}}_{3}=\{{R}_{3},{{\rm{\Theta }}}_{3}=0,{{\rm{\Phi }}}_{3}=0\},$ with ${{\rm{\Theta }}}_{3}$ and ${{\rm{\Phi }}}_{3}$ the polar and azimuthal angles, ${\vec{R}}_{1}={\vec{r}}_{1}-{\vec{r}}_{2},$ ${\vec{R}}_{2}={\vec{r}}_{4}-{\vec{r}}_{3},$ ${\vec{r}}_{1}=\xi {\vec{R}}_{1}$ and ${\vec{r}}_{2}=(1-\xi ){\vec{R}}_{1},$ where $\xi ={m}_{2}/({m}_{1}+{m}_{2})$ [18]. Further, we adopted the OXYZ system in such a way that the HD inter-atomic vector ${\vec{R}}_{1}$ lies on the XOZ plane. Then the angle variables of ${\vec{R}}_{1}$ and ${\vec{R}}_{2}$ are: ${\hat{R}}_{1}=\{{{\rm{\Theta }}}_{1},{{\rm{\Phi }}}_{1}=\pi \}$ and ${\hat{R}}_{2}=\{{{\rm{\Theta }}}_{2},{{\rm{\Phi }}}_{2}\}$ respectively. One can see, that the Cartesian coordinates of the atoms of the HD molecule are [18]: ${\vec{r}}_{1}=\left\{{x}_{1}=\xi {R}_{1}\mathrm{sin}{{\rm{\Theta }}}_{1},{y}_{1}=0,\right.$ $\left.{z}_{1}=\xi {R}_{1}\mathrm{cos}{{\rm{\Theta }}}_{1}\right\}$, ${\vec{r}}_{2}=\left\{{x}_{2}=-(1-\xi )\right.$ ${R}_{1}\mathrm{sin}{{\rm{\Theta }}}_{1}$, ${y}_{2}=0,{z}_{2}=-(1-\xi )$ $\left.{R}_{1}\mathrm{cos}{{\rm{\Theta }}}_{1}\right\}.$ Defining $\zeta ={m}_{4}/({m}_{3}+{m}_{4}),$ we have ${\vec{r}}_{3}={\vec{R}}_{3}-(1-\zeta ){\vec{R}}_{2},$ ${\vec{r}}_{4}={\vec{R}}_{3}+\zeta {\vec{R}}_{2},$ and the corresponding Cartesian coordinates are: ${\vec{r}}_{3}=\left\{{x}_{3}=-(1-\zeta ){R}_{2}\right.$ $\mathrm{sin}{{\rm{\Theta }}}_{2}\mathrm{cos}{{\rm{\Phi }}}_{2},$ ${y}_{3}=-(1-\zeta ){R}_{2}\mathrm{sin}{{\rm{\Theta }}}_{2}\mathrm{sin}{{\rm{\Phi }}}_{2},$ ${z}_{3}={R}_{3}-(1\left.-\zeta ){R}_{2}\mathrm{cos}{{\rm{\Theta }}}_{2}\right\},$ ${\vec{r}}_{4}=\left\{{x}_{4}=\zeta {R}_{2}\mathrm{sin}{{\rm{\Theta }}}_{2}\mathrm{cos}{{\rm{\Phi }}}_{2},\right.$ ${y}_{4}=\zeta {R}_{2}\;\mathrm{sin}{{\rm{\Theta }}}_{2}\mathrm{sin}{{\rm{\Phi }}}_{2},$ $\left.{z}_{4}={R}_{3}+\zeta {R}_{2}\mathrm{cos}{{\rm{\Theta }}}_{2}\right\}.$ In such a manner the cartesian and the Jacobi coordinates are represented together for the four-atomic system HD+H₂ [18].

The Hinde H₂-H₂ PES [5] is a six-dimensional surface which was constructed using recent Raman spectrum data of the (H₂)₂ dimer and which accurately describes the dimer's van der Waals well [5]. It was demonstrated that this PES gives IR and Raman transition energies for the (para-H₂)₂, (ortho-D₂)₂, and (para-H₂)−(ortho-D₂) dimers and is in good agreement with experimental data.

The method to make the Hinde H₂-H₂ PES suitable for the non-symmetric system HD+H₂ is based on a geometrical transformation technique, i.e. a rotation of the three-dimensional space and the corresponding space-fixed coordinate system OXYZ. The new global PES depends on six variables (Jacobi coordinates), $| {\vec{R}}_{1}| ,$ $| {\vec{R}}_{2}| ,$ $| {\vec{R}}_{3}| ,$ ${\theta }_{1},$ ${\theta }_{2},$ and ${\varphi }_{2}$, as shown in figure 1 together with corresponding quantum angular momenta. The initial geometry of the system is designed in such a way that the Jacobi vector ${\vec{R}}_{3}$ connects the c.m.'s of the two H₂ molecules and is directed over the ${OZ}$ axis. We laid out OXYZ in such a manner that the Jacobi vector ${\vec{R}}_{1}$ lies in the XOZ plane. The vector ${\vec{R}}_{2}$ can then be directed anywhere. The spherical coordinates of the Jacobi vectors are: ${\vec{R}}_{1}=\{{R}_{1},{\theta }_{1},0\},$ ${\vec{R}}_{2}=\{{R}_{2},{\theta }_{2},{\varphi }_{2}\},$ and ${\vec{R}}_{3}=\{{R}_{3},0,0\}.$ Because we used the rigid rotor model, the lengths of the HD and H₂ molecules are fixed at equilibrium values, e.g. R₁ = 0.7631 a.u. and R₂ = 0.7668 a.u., thus leaving us with four free variables. We replace one hydrogen atom H with a deuterium atom D, thus shifting the c.m. of one H₂ molecule to another point, that is from O to ${O}^{\prime }.$ The length of the vector $\vec{x}$ is $x=| {\vec{R}}_{1}| /6.$ This is seen in figure A1 . Then we shift the original coordinate system OXYZ along the vector ${\vec{R}}_{1}$ to the new one ${O}^{\prime }{X}^{\prime }{Y}^{\prime }{Z}^{\prime }.$ The origin of the new system, i.e. ${O}^{\prime },$ lies on the c.m. of HD. A new Jacobi vector ${\vec{R}}_{3}^{\prime }$ is defined connecting the c.m.'s of the HD and H₂ molecules. The new intermolecular distance between HD and H₂ is:

$$\begin{eqnarray}&&{R}_{3}^{\prime }=\sqrt{{R}_{3}^{2}+{x}^{2}-2{{xR}}_{3}\mathrm{cos}{\theta }_{1}}.\end{eqnarray} \tag{ 17 }$$

Zoom In Zoom Out Reset image size

Now, if we rotate ${O}^{\prime }{X}^{\prime }{Y}^{\prime }{Z}^{\prime }$ around its ${O}^{\prime }{Y}^{\prime }$ axis in such a way that the ${{OZ}}^{\prime }$ axis is directed over the vector ${\vec{R}}_{3}^{\prime }$ we obtain a new coordinate system ${O}^{\prime }X^{\prime \prime} Y^{\prime \prime} Z^{\prime \prime}$ which should be well designed to carry out computations for the HD+H₂ collision. The rotational angle η satisfies:

$$\begin{eqnarray}&&\mathrm{cos}\eta =\displaystyle \frac{R{{}^{\prime }}_{3}^{2}+{R}_{3}^{2}-{x}^{2}}{2{R}_{3}^{\prime }{R}_{3}}\end{eqnarray} \tag{ 18 }$$

$$\begin{eqnarray}&&\mathrm{sin}\eta =\displaystyle \frac{x}{{R}_{3}^{\prime }\mathrm{sin}{\theta }_{1}}.\end{eqnarray} \tag{ 19 }$$

This transformation converts the initial Jacobi vectors in OXYZ: ${\vec{R}}_{1}=\{{R}_{1},{\theta }_{1},0\},$ ${\vec{R}}_{2}=\{{R}_{2},{\theta }_{2},{\varphi }_{2}\}$ and ${\vec{R}}_{3}=\{{R}_{3},0,0\}$ to the corresponding Jacobi vectors with new coordinates in the new ${O}^{\prime }{X}^{\prime }{Y}^{\prime }{Z}^{\prime }:$ ${\vec{R}}_{1}^{\prime }=\{{R}_{1}^{\prime },{\theta }_{1}^{\prime },{\varphi }_{2}^{\prime }\},$ ${\vec{R}}_{2}^{\prime }=\{{R}_{2}^{\prime },{\theta }_{2}^{\prime },0\}$ and ${\vec{R}}_{3}^{\prime }=\{{R}_{3}^{\prime },0,0\}.$

The coordinate transformations from OXYZ to $O^{\prime \prime} X^{\prime \prime} Y^{\prime \prime} Z^{\prime \prime}$ changes the original Hinde's PES to the new HD-H₂ PES. Rotation of the coordinate system from ${O}^{\prime }{X}^{\prime }{Y}^{\prime }{Z}^{\prime }$ to ${O}^{\prime }X^{\prime \prime} Y^{\prime \prime} Z^{\prime \prime}$ results in a corresponding transformation of the coordinates of the 4-body system as well as changing the distance between the two molecules. One then has the following relations between new and old variables [40]:

$$\begin{eqnarray}&&{\theta }_{1}^{\prime }=\mathrm{arccos}(\mathrm{cos}{\theta }_{1}\mathrm{cos}\eta +\mathrm{sin}{\theta }_{1}\mathrm{sin}\eta ),\end{eqnarray} \tag{ 20 }$$

$$\begin{eqnarray}&&{\theta }_{2}^{\prime }=\mathrm{arccos}(\mathrm{cos}{\theta }_{2}\mathrm{cos}\eta +\mathrm{sin}{\theta }_{2}\mathrm{sin}\eta \mathrm{cos}{{\rm{\Phi }}}_{2}),\end{eqnarray} \tag{ 21 }$$

$$\begin{eqnarray}&&{\varphi }_{2}^{\prime }=\mathrm{arccos}\left(\mathrm{cot}{\phi }_{2}\mathrm{cos}\eta +\displaystyle \frac{\mathrm{cot}{\theta }_{1}\mathrm{sin}\eta }{\mathrm{sin}{\phi }_{2}}\right).\end{eqnarray} \tag{ 22 }$$

In the calculation of HD+H₂ with Hinde's PES one has to use new coordinates ${\theta }_{1}^{\prime },{\theta }_{2}^{\prime },{\varphi }_{2}^{\prime },{R}_{3}^{\prime }.$ However, the original potential has been expressed through the initial H₂-H₂ variables, i.e. ${\theta }_{1},$ ${\theta }_{2},$ ${{\rm{\Phi }}}_{2}$ and R₃. Hence they have to be transformed using (20) and (22). Therefore, in the case of the non-symmetrical HD+H₂ collision one should use the formulas equations (20)–(22) together with equations (18)–(19) and the expression (17) for the new distance ${R}_{3}^{\prime }$ between the center of masses of the H₂ and HD molecules.

In general, any consideration of the HD+H₂, D₂+HD or D₂+D₂ systems should begin with the original H₂-H₂ PES. This six-dimensional function comprises a symmetrical surface over the OZ coordinate axis. This is shown in figures 1 and A1. In spherical coordinates the surface can be described by six variables: R₁, ${\theta }_{1}$, R₂, ${\theta }_{2}$, R₃, and ${\varphi }_{2}.$ The variables are also shown in the figures. The H₂-H₂ PES was obtained in the framework of the Born-Oppenheimer model [32] and can be considered as a symmetrical interaction field. When considering non-reactive scattering problems with participation of hydrogen molecules one needs to solve the Schrödinger equation (2) with the H₂-H₂ potential $V({\vec{R}}_{1},{\vec{R}}_{2},{\vec{R}}_{3}).$ The solution/propagation runs over the ${\vec{R}}_{3}$ Jacobi vector (please see figure A1). Therefore, in the case of the symmetrical H₂+H₂ and D₂+D₂ collisions one can use the original H₂-H₂ PES as it is, i.e. without transformations.

However, in the case of the non-symmetrical (or symmetry-broken) HD+H2/D2 or HD+HD scattering systems one should also apply the original H₂-H₂ interaction field (PES), but the propagation (solution) of the Schrödinger equation runs, in this case, over the corrected Jacobi vector ${\vec{R}}_{3}^{\prime }$ which is directed over the new ${O}^{\prime} {Z}^{\prime \prime}$ axis, as is shown in figure A1.

New data for the thermal rate coefficients, ${k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T),$ equation (4), are listed below. The results are obtained with the modified Hinde H₂-H₂ PES [5]: table 2 includes ${k}_{{j}_{1}{j}_{2}\to {j}_{1}^{\prime }{j}_{2}^{\prime }}(T)$ at low temperatures for the HD + para-H₂ rotational energy transfer collisions and table 3 the same information for the HD + ortho-H₂ case. The astrophysical HD-cooling function can be computed with the use of formula (16).