State-to-state integral cross sections and rate constants for the N+(3P)+HD→NH+/ND++D/H reaction:Accurate quantum dynamics studies

Hanghang Chen(陳航航), Zijiang Yang(楊紫江), and Maodu Chen(陳茂篤)

Key Laboratory of Materials Modification by Laser,Electron,and Ion Beams(Ministry of Education),School of Physics,Dalian University of Technology,Dalian 116024,China

Keywords: quantum dynamics,integral cross sections,rate coefficients

1. Introduction

Nitrogen hydrides are of great interest in astrophysical environments. First,nitrogen is one of the most abundant elements in the universe,so NH,NH2,NH3or more complex and the corresponding ions can serve as important tracers for the composition of various parts of the interstellar medium.[1,2]Furthermore, nitrogen hydrides and their ionic variants have great significance in the formation of biologically relevant molecular species under interstellar conditions. The simplest nitrogen hydrogen ion of NH+is mainly generated by the N+(3P)+H2→NH++H reaction,which is also the first step of the synthesis of interstellar ammonia.[1-5]On the reaction dynamics aspect,the reactions of N+(3P)with H2and its isotopic variants (D2and HD) are of fundamental interest due to the thermally neutral characteristic and the domination of a deep well, suggesting that the quantum effects, such as resonances and vibrational zero-point energy,play a crucial role in the dynamical mechanisms.[6,7]

A variety of experimental data of the reactions of N+(3P)+H2/D2/HD, such as cross sections,[8,9]rate constants,[10-13]and product angle distributions,[14,15]have been extensively reported. Hansenet al.[14]measured the product velocity vector distributions of the N+(3P)+H2reaction in the initial energy ranges of 0.98 eV to 3.60 eV using crossed beam technique. The scattering peak implies the predominance of the long-lived NH+2complex at energies below 1.9 eV.Sunderlinet al.[9]measured the integral cross sections(ICSs)of the N+(3P)+HD reaction as a function of collision energy by using a guided ion beam tandem mass spectroscopy. The ICSs of the two product molecules reveal that the intramolecular isotope effect has a marked effect on the N+(3P)+HD reaction. Using a low-temperature 22-pole ion trap apparatus, Zymaket al.[12]extracted the state-resolved rate coefficients for the reaction of N+(3Pja)+H2(j)between 10 K and 100 K.The results clearly indicate that the excitation of N+to the3P2state reduces the reactivity.

There is a general consensus that the reactions of N+with H2and its isotopic variants through a deep potential well,lead to the formation of the long-lived complex NH+2at low energies. Thus, various statistical methods[16-19]have been employed to deal with the reactions of N+with H2and its isotopic variants in terms of theory. Recently,Grozdanovet al.[16]formulated a new statistical model that treats the fine structure of N+ions equally to other internal motions. The cross sections of the N+(3P)+HD reaction indicate that the ND+product is in a dominant position at low collision energy. In the following work,[17]the rate coefficient of the N++H2/HD reaction was determined in the temperature range of 10-200 K.However, the accuracy of statistical methods is limited due to the addition of certain assumptions and approximations.Hence,to describe the reaction characteristics,it is necessary to perform the trajectory[20,21]or quantum mechanical calculations on an exact reactive potential energy surface (PES). In 1998, Russellet al.[22]performed the first quantum mechanical calculation for the reactions of N+(3P)+H2on the PES constructed by Wilhelmssonet al.[23]They investigated several primary quantum mechanical effects,including the existence of quantum resonances and the formation of long-lived intermediate complexes. Regrettably, the previous PESs are not accurate enough due to limitations of earlier methods and computational power, thus the dynamics data of this system cannot be accurately determined. Recently, a high-precision global PES (called YWYC PESs.[24]) of the ground state NH+2system was constructed using a large number of high-levelab initioenergy points and a neural network model. Based on this PES,the state-to-state quantum dynamics for N+(3P)+H2and N+(3P)+D2reactions were studied in detail to explore the dynamics mechanism.

Intramolecular isotype effects in the simplest reaction type of X+HD has been a subject of considerable interest both theoretically and experimentally. For the N+(3P)+H2reaction system,the isotope effects play an important role because of its near thermal neutrality. Previous studies[24]have confirmed that when H2is substituted by D2,there are important changes on the rovibrational distribution and state-resolved angular distribution of the product molecule. However, the quantum dynamics results of the N+(3P)+HD reaction were not reported up to data, and it is also necessary to perform rigorous quantum dynamics calculations on the N+(3P)+HD at the state-to-state level,which not only contributes to further understand the intramolecular isotype effects,but also can provide more comprehensive and accurate data for this important interstellar reaction system.

In this paper, accurate quantum TDWP calculations are carried out at the state-to-state level for N+(3P0)+HD(v=0,j=0)→NH++D and N+(3P0)+HD(v=0,j=0)→ND++H reactions based on the YWYC PESs.The present work aims to present accurate total and state-resolved rate constants of the title reactions,and the microscopic mechanism in both the channels is also analyzed in detail. The remaining part of this article is structured as follows: Section 2 briefly introduces the theoretical methods. Section 3 analyzes and discusses the dynamic results of the title reaction. Finally, Section 4 gives conclusions.

2. Method and computational details

For the atom-diatom reactive scattering, the TDWP approach is a mighty measure for studying molecular reaction dynamics at the state-to-state level and has been applied to numerous reaction systems.[25-34]A TDWP method based on the reactant-coordinate-based (RCB) method is employed to study the title reaction in the present work. Here,we will discuss only the essentials, and the detailed theoretical methods have been described in the previous literatures.[30,35-39]In order to represent the wave function in different configurations,three sets of body-fixed (BF) Jacobi coordinatesα,βandγare utilized to represent the reactant N+(3P)+HD,the product NH++D and the product ND++H,respectively. For a given total angular momentumJ, the system’s Hamiltonian in the BF reactant Jacobi coordinates can be written as

whereRrepresents the distance from the N+ion to the mass center of the HD molecule andris the bond length of HD,μRandμrare the corresponding reduced masses,respectively.?Jis the total angular momentum operator of the NHD+system and ?jis the rotational angular momentum operator of HD molecular.Vis the potential energy that can be obtained from PES.

For the specified initial state (v0,j0,l0) of the title reaction,the initial wave packet constructed in the space-fixed(SF)representation can be expressed as

whereG(Rα)is a Gaussian wave packet along the coordinateR,φv0j0(rα)is the rovibrational eigenfunction for the reactant molecule HD, and|JM j0l0ε〉is the total angular momentum eigenfunction with parity of systemε=(-1)j0+l0and needs to be converted to the BF representation before evolution. The propagation of the wave packet is carried out in the reactant Jacobi coordinates in the BF representation using the secondorder split-operator method.

During the wave packet propagation, the construction of absorption potential can effectively avoid the reflection of the wave packet from the grid edge. In our calculations, the absorption potential alongrandRdirections is defined as a sectional function

whereCaandCbdefine the strength of absorbing potentials,andxa,xb, andxendare the positions of the absorbing potentials.

Next, we apply the reactant-coordinate-based(RCB)[30,36]method to extract the state-to-stateS-matrix. Regarding this method, the wave packet evolves exclusively in the body-fixed reactant Jacobi coordinates, and theS-matrix of the two reaction channels can be obtained by a single propagation. Finally,the dynamic results can be obtained using the state-to-stateS-matrix obtained in the helicity representation.

The state-to-state reaction probability is calculated by

wherekBis the Boltzmann constant.

Table 1. Numerical parameters used in the TDWP calculations.

In this work, the TDWP calculations on YWYC PESs are performed for the ground spin-orbital state (3P0) of N+ions and the ground rotational state (v0=0,j0=0) of HD molecules. To perform the TDWP calculations, we conduct the convergence tests on the title reaction by theJ= 0 reaction probability. The calculated maximum partial wave isJ=77, which corresponds to an upper collision energy limit of 0.40 eV for the converged ICSs.For different partial waves,we apply different total propagation time to ensure convergence.The primary numerical parameters are listed in Table 1.

All the datasets presented in this paper,including the ICSs and rate coefficients of the two products for the title reaction are compiled in the supplementary materials(SM).

3. Results and discussion

The schematic diagrams of the global minimum energy path (MEP) for the N+(3P) +HD→NH+/ND++D/H reaction are shown in Fig. 1. It can be seen that there is an extremely deep potential well of 6.48 eV on the global MEP, manifesting that there exist long-lived reactive complexes. Moreover, the deepest position of the potential well corresponds to the perpendicular structure. At this time,the bond length of HD molecule and the distance between the N+ion and the center of the molecule are 3.75a0and 0.50a0, respectively. On the other hand, the N+(3P) +HD reaction is an endothermic reaction with a tiny endotherm of 0.085 eV. The zero-point energies (ZPEs) of HD, NH+and ND+molecules calculated on the analytical YWYC PESs are 0.237 eV, 0.176 eV and 0.127 eV, respectively. Therefore,when considering the ZPE effect, the NH++D channel and ND++H channel exhibit endothermic and exothermic characters,respectively. In conclusion,ZPE has a decisive role in the sign of the energy change in the title reaction.

Fig.1. Schematic diagram of the global minimum energy path for the N+(3P)+HD →NH+/ND++D/H reaction.The dotted line represents the reactant or product energy after considering the ZPE effect.The unit of bond length is a0.

The total reaction probabilities of the NH++D and ND++H product channels as a function of collision energy at several different partial waves(J=0,30,60)are plotted in Fig.2. Some significant resonance structures are observed on the probability curves of both channels,especially for low partial waves.This can be attributed to the formation of long-lived intermediate complexes supported by the deep potential wells during the reaction. The resonance amplitudes of both channels are much larger than the other typical complex-forming reactions,such as S+H2,[40-43]F+H2[44-46]and O+H2,[47]due to the deeper potential well. TheJ=0 probability curve of the NH++D channel has a threshold of 0.024 eV, which is consistent with the corresponding endothermicity. On the contrary,the ND++H channel corresponds to an exothermic behavior without a threshold due to the ZPE effect. The probability curve of the ND++H channel generates a threshold and gradually increases with the increase ofJunder the influence of the continuously increasing centrifugal barrier. In addition,regarding large partial waves,the dominant role of the ND++H channel is replaced by the NH++D channel. This is because the enough large centrifugal potential smoothens the structure of the well, leading to the emergence of a direct abstraction process that favors the production of the NH+molecule.

Fig. 2. Reaction probabilities of the NH++D and ND++H product channels as a function of collision energy at several selected partial waves(J=0,30,60).

The collision energy dependence of total ICSs for both the product channels are presented in Fig.3. As shown in the figure,the total ICS of the NH++D channel increases slowly with the collision energy owing to its endothermic properties. When the collision energy is less than 0.15 eV, the total ICS of the ND++H channel decreases significantly with the increase of collision energy, which is a typical feature of exothermic reactions without an energy threshold. However,its ICS witnesses a slowly increasing trend as the collision energy continually increases. The reason for the increase of ICS at high collision energy is that the increasing contribution of large centrifugal barriers leads to an ascending collision energy dependence. The increased centrifugal barrier reduces the effect of the potential well, causing the direct extraction process to come into play,although not yet comparable to the insertion reaction process. On the other hand, intricate resonance exists on the ICSs curve even if all the partial waves have been summed,which is the manifestation of the intermediate complex generated by the deep potential well. To obtain more sufficient information about the intramolecular isotope effect, the ND+/NH+product branching ratio is also shown in Fig.3. As seen,the ND+/NH+product branching ratio increases continuously with the collision energy, but is always less than 1. Therefore, the ND++H channel always plays a major role in the selected collision energies, which indicates that the N+(3P)+HD reaction is dominated by the insertion reaction mechanism.

Fig. 3. Total ICSs of the N+(3P)+HD→NH+/ND++D/H reactions and the product NH+/ND+branching ratio in the collision energy range of 0.02 to 0.4 eV.

Product vibrational state-resolved ICSs and vibrational branching ratio of ICS(v＇= 1)/ICS(v＇= 0) are depicted in Fig. 4. As shown, only the vibrational ground state and the first excited state are populated for the products of both channels in the selected collision energy range.This means that the products of the title reaction are difficult to excite to higher vibrational states. This can be attributed to the fact that the formation of a long-lived complex supported by the deep well in the reaction path,since the product molecules are generated only when the rotational quantum number of the complex is large enough to split the HD bond. The collision energy is mainly converted into the translational energy of the product molecules, whereas the corresponding internal energy is relatively low. Moreover, it is also visible from Fig. 4 that the majority of the NH+or ND+products are distributed on the ground vibrational state. The vibrational branching ratio of ICS(v＇=1)/ICS(v＇=0)for the ND++H and NH++D channels are less than 0.03 and 0.25,respectively. The ICS(v＇=1)presents a threshold and then increases with the collision energy for both channels. Furthermore, the threshold of the NH++D channel is larger.This phenomenon is mainly due to the exothermicity of the ND++H channel caused by the ZPE effect,promoting the opening of itsv＇=1 channel.

Fig.4. Product vibrational state-resolved ICSs and vibrational branching ratios of ICS(v＇ = 1)/ICS(v＇ = 0) for the N+(3P) +HD →NH+/ND++D/H reaction.

The product rotational state-resolved ICSs of both channels at several selected collision energies are shown in Fig.5.It should be noted that the results here include all the vibrational states. As seen from the figure, more rotational states are populated for both the channels as the collision energy increases. Also, the peak of rotational state-resolved ICSs becomes smaller and transfers to higherj＇values. The rotational state distributions of ND+molecule are relatively hotter and wider at a specific collision energy compared with the NH++D channel. The differences in the rovibrational state distributions between the two products are owing to the truth that the vibrational frequency and rotational constant of the ND+molecule are smaller than those of the NH+molecule.

Accurate knowledge of the rate coefficients for the N++HD reaction has important implications for studying the role of HD molecules on the interstellar synthesis of nitrogen hydrides in the early universe. Accordingly,the rate coefficients for the total and two product channels of the N++HD reaction in the temperature range from 20 K to 1000 K are calculated,as shown in Fig.6.As seen,the main population of the rate coefficient is located in the ND++H channel, especially in the low-temperature region, which can be well explained by the thermal properties of these two channels. The rate coefficient of the ND++H channel increases gradually with increasing temperature until 1300 K, and then shows a decreasing trend with further increasing temperature. The difference is that the rate coefficient of the NH++D channel increases consistently with temperature over the studied temperature range.

Fig. 5. Product rotational state-resolved ICSs for the N+(3P)+HD→NH+/ND++D/H reaction at several different collision energies.

Fig. 6. Temperature variation of the rate coefficients for the total and two product channels of the N+(3P)+HD reaction.

As described above,the products of the N++HD reaction are mainly distributed in the ground vibrational state. Therefore,as plotted in Fig.7,the rotational state-resolved rate coefficients of the dominant vibrational states(v＇=0)are studied over the temperature range from 10 K to 2000 K for both the product channels. As seen in Figs.7(a)and 7(c),the state population inversion is observed in the low rotational states of both channels. The maximum population of the NH+and ND+products is located in(v＇=0,j＇=3)and(v＇=0,j＇=4),respectively. The rate coefficient first increases and then decreases with the rotational excitation. Moreover,the rotational state-resolved rate coefficient exhibits a strong temperature dependence. The rate coefficients of ND++H channel for low rotational states(j＇≤7)tend to increase rapidly and then slowly decrease. This indicates that the negative temperature dependence of the total rate coefficient of the ND++H channel at high temperature is primarily caused by the low rotation states. Regarding high rotational states, the rate coefficients increase gradually and then tend to be constant,which is similar to the case for all rotational states of the NH++D channel.

Fig.7. Temperature variation of the rate coefficients for the ND++H and NH++D channels for different product rotational states(v＇=0, j＇).

4. Conclusions and perspectives

In this work, accurate ICSs and rate coefficients for the both channels of the N+(3P)+HD(v0=0,j0=0)reaction are calculated on YWYC PESs using the TDWP method. The collision energy dependence of all dynamics results shows obvious resonance structures, which are caused by the longlived complex supported by the deep well. Under the influence of the ZPE effect, the N+(3P) +HD→NH++D and N+(3P)+HD→ND++H are endothermic and exothermic reactions, respectively. The ND++H channel always maintains a dominant position in the studied collision energy range due to the different thermal properties and the features of the PES.Furthermore,the vibrationally resolved ICSs reflect that both the NH+and ND+products are mainly populated on the ground vibrational state and difficult to excite higher vibrational states.The rotational state distribution of the ND+product is hotter than that of the NH+product within the selected collision energy range. Herein, the rate coefficients of both channels of the N++HD reaction are calculated and compared over the temperature range 10 K≤T ≤2000 K. The results show that the rate coefficients are mainly concentrated in the ND++H channel. In addition,the state population inversion is observed in the rotational state-resolved rate coefficients for the dominant vibrational states(v＇=0)of both channels. The maximum population of product rotational state-resolved rate coefficients is located at(v＇=0,j＇=4)and(v＇=0,j＇=3)of the ND++H and NH++D channels,respectively.

As far as we know, no experiments with such high fine structure states have been performed for the title reaction. We expect that our work will provide an important reference for the further experimental studies at the finer level for this interstellar chemical reaction.

Data availability

The data that support the findings of this study are openly available in Science Data Bank at https://www.doi.org/10.57760/sciencedb.j00113.00034.

Acknowledgement

Project supported by the National Natural Science Foundation of China(Grant No.11774043).