Rate constants for a reactive system of astrophysical interest: a statistical study of CH₂⁺

Abstract

A systematic investigation on the different processes involving formation, destruction and (de)excitation of specific rovibrational states of CH⁺ has been carried out based on statistical approaches. Thus, reactive collisions between C⁺(²P) and H₂(v, j) and between CH⁺(v, j) and H for a large number of state-to-state transitions have been studied using a statistical quantum method and a statistical adiabatic channel model. The capabilities of such techniques for the study of the title system are discussed with comparisons to previous quantum mechanical results and experimental data. Integral cross sections as a function of the energy and rate constants in terms of the temperature (up to 1500 K) are obtained and numerical data for astrophysical purposes are provided.

---

1 Introduction

Ionized hydrocarbons play a pivotal role in the chemistry of the interstellar medium (ISM), particularly in regions exposed to intense ultraviolet radiation. Among these, the methylidynium ion (CH⁺) has long been recognized as a fundamental molecular ion that drives key gas-phase reactions leading to the formation of more complex carbon-chain and organic species. First detected in absorption in the diffuse ISM by Douglas and Herzberg,¹ CH⁺ has since been observed in a wide range of astrophysical environments, including planetary nebulae,^2,3 massive star-forming regions,⁴ protoplanetary disks,⁵ and external galaxies.⁶

The recent advent of the James Webb Space Telescope has opened new observational windows in the infrared regime, enabling the detection of rovibrational emission lines from reactive hydrocarbon ions with unprecedented sensitivity and spatial resolution. For the first time, CH⁺ along with the methylium ion (CH₃⁺) have been detected toward both the Orion Bar and the d203–506 protoplanetary disk regions that are heavily irradiated by the trapezium star cluster.⁷ Rovibrational transitions from both the v = 1 → 0 and v = 2 → 1 vibrational bands of CH⁺ were identified, and the chemical formation pumping process was beautifully evidenced, as in the young planetary nebula NGC 7027.³

The detection of CH⁺ emission in such photodissociation regions thus provides critical understanding into the physical conditions and ionization mechanisms that regulate molecular synthesis in space. CH⁺ is known to form efficiently in highly energetic regions where radiative or shock-driven processes provide extra energy to overcome the barrier in its formation reaction C⁺(²P) + H₂. However, to use CH⁺ rovibrational emission as a powerful diagnostic tool requires a very good knowledge of rate constants for the dominant formation, destruction and excitation channels, at the state-to-state as possible.

Given the importance of CH⁺ as discussed above, the following reactions

C⁺ + H₂(v, j) → H + CH⁺(v′, j′) (1)

H + CH⁺(v, j) → C⁺ + H₂(v′, j′) (2)

→ H + CH⁺(v′, j′) (3)

have been subject of numerous investigations. Early experiments by means of a guided ion beam technique by Ervin and Armentrout^8,9 allowed to know about the threshold behaviour of reaction (1). Cross sections for reactions involving different isotopic forms of molecular hydrogen (H₂, D₂ and HD) to form CH⁺ or CD⁺ were measured and rate constants were derived. In particular, a value of 1.2 × 10⁻¹⁶ cm³ s⁻¹ was obtained at T = 300 K for the C⁺(²P) + H₂(v = 0) → CH⁺ + H reaction. Rate constants were also measured between 400 and 1300 K using a high-temperature flowing afterglow apparatus.¹⁰

Among the most recent theoretical efforts to study these processes, Jara-Toro et al.¹¹ used a new approach based on the canonical variational transition state theory to obtain rate constants for the formation of C⁺ according to reaction (2). They concluded that the participation of two pathways, a direct one and an indirect involving either a van der Waals complex or an intermediate CH₂⁺ species, has to be taken into account. Non-adiabatic effects were suggested as responsible for the decrease of the H + CH⁺ rate constants observed in the experiment¹² at low temperature (below 50 K). However, the inclusion of electronic couplings was found unable to describe the dramatic reduction of the experimental rate constant observed when the temperature decreases below 50 K in the investigation performed by del Mazo-Sevillano et al.¹³ Li et al.¹⁴ also treated the non-adiabatic dynamics of the C⁺(²P_1/2,3/2) + H₂ reaction by means of a time dependent wave packet (TDWP) calculation based on a global diabatic potential energy surface (PES) and no noticeable qualitative differences were observed between adiabatic and non-adiabatic rate constants, especially at low temperatures.

Some other numerous TDWP investigations^14–17 have treated the dynamics of reactions (1)–(3). Li et al.,¹⁸ for example, calculated probabilities, integral cross sections (ICSs), differential cross sections (DCSs), rovibrational state distributions and rate constants for the C⁺ + H₂ reaction. In particular, in view of the essentially forward–backward symmetry exhibited by the DCSs at specific collision energies, the authors of ref. 18 concluded that the overall dynamics of the C⁺ + H₂ → CH⁺ + H seems to be governed by a complex-forming mechanism, one of the assumptions required for statistical methods to be applicable in this sort of atom–diatom reactive collisions. In fact, and due to the presence of a relatively deep potential well corresponding to the intermediate CH₂⁺ species between reactants and products, some of the previous theoretical studies have been performed by means of statistical approaches. Despite measured DCSs were not completely symmetric with a slight preference for the forward scattering direction to invoke categorically a dynamics mediated by a complex-forming mechanism, the occurrence of an “osculating” intermediate CH₂⁺ species was assumed as appropriate to describe the overall dynamics of the collision.¹⁹ Moreover, it was shown²⁰ that the threshold behaviour of measured cross sections for the C⁺ + D₂ reaction was consistent with statistical considerations from the theory by Pechukas and Light.²¹ Therefore early applications of the phase space theory (PST) were tried in several occasions^22–27 to obtain cross sections and rate constants for the title reaction or isotopic variants with deuterium atom(s) instead of hydrogen. In particular, the rate constants reported by Gerlich et al.²⁷ were found to reproduce nicely the high temperature flowing afterglow measurements for C⁺ + H₂ by Hierl et al.¹⁰ between 500 and 1300 K.

More recently, the reverse process H + CH⁺ has been investigated as well with PST by Halvick et al.²⁸ and cross sections and low-temperature rate constants for both the exchange channel forming H + CH⁺ in eqn (3) and the abstraction channel producing C⁺ + H₂ in eqn (2) were calculated. On the other hand, the existence of a deep potential well (with the minimum located more than 4 eV below the C⁺ + H₂ channel) between reactants and products corresponding to the CH₂⁺ species in the most recent PES describing the process manifested, for instance, in the numerous narrow resonances seen in cross sections,¹⁵ supporting the possible validity of this kind of techniques to investigate the title reaction. However, in recent works, observed deviations from a purely statistical description of the dynamics of the process have been pointed out. First, Zanchet et al.¹⁵ concluded that the reaction is not completely statistical due to the overestimation of predictions from a very simple model based on counting energetically accessible levels in comparison with quantum wave packet probabilities for the C⁺ + H₂(v = 0, 1, j = 0) reaction for a total angular momentum J = 0. Total ICSs also from the vibrational ground state of H₂ obtained with a statistical quasiclassical trajectory (SQCT) method^29–31 were also larger than the corresponding results from a strict quasi-classical trajectory (QCT) calculation.³² The authors of ref. 32 then concluded that the reaction dynamics of the process cannot be described according to a purely unbiased statistical model.

Despite these observations, the H + CH⁺ was chosen as one of the test examples to try the improved version³³ of the statistical adiabatic channel model (SACM) by Quack and Troe.^34–36 Rate constants for this process obtained by means of the SACM were compared with results from a coupled-channel calculation and the agreement was found to be good, especially at low temperatures (20–100 K).³⁷ In their investigation on the dynamics of the reaction in eqn (2), del Mazo et al.¹³ compared cross sections obtained using an adiabatic statistical technique with wave packet results for the H + CH⁺(v = 0, j = 0) → C⁺ + H₂ in both an adiabatic and 3 × 3 diabatic model. Moreover, initial-state-selected rate constants for a variety of rovibrational CH⁺(v, j) states obtained with such a statistical approach were compared with QCT values.

Thus, encouraged by the previous use of statistical techniques in this context, in this work we report results of a systematic calculation of rate constants for the state-to-state processes shown in eqn (1)–(3) by means of the statistical quantum method (SQM), as described in ref. 38 and 39. This approach has been employed successfully for the same astrophysical purpose in the case of the H⁺ + H₂^40–42 and H⁺ + HD⁴³ reactions. The rate constants so obtained, were then incorporated to compute the collisional cooling of a gas by molecular hydrogen under conditions appropriate to the primordial and interstellar media⁴⁴ and also to follow the level populations of H₂ and HD during the adiabatic expansion of the (early) Universe.⁴⁵

Here, we have applied the SQM to processes (1)–(3), in order to calculate ICSs as a function of the energy and then rate constants for temperatures up to 1500 K. Comparisons with available proper quantum mechanical (QM) results have been established for low rovibrational states as (v = 0, 1, j = 0–1). SACM results up to T = 500 K have also been incorporated.

The structure of this paper is as follows: in Section 2 we discuss some basic details of the present theoretical approaches, in Section III we present results for reaction eqn (1) in Section 3.1); for reaction (2) in Section 3.2.1 and for reaction (3) in Section 3.2.2.

2 Theoretical methods

2.1 Statistical quantum method

As explained before in previous works^38,39 the SQM is specifically designed for insertion reactions where the whole process can be treated by separating the formation of an intermediate complex from the reactants and its ulterior fragmentation onto the products. Within such an approximation, the state-to-state probability from the initial (v, j) state of the reactants diatom (with v and j being the QM numbers for the vibrational and rotational levels, respectively) and the final rovibrational (v′, j′) state of the products diatom at the energy E and the total angular momentum J can be expressed as:where p^J_vj(E) is the individual intermediate complex formation probability and to be formed from the initial (v, j) state, and , the fraction of fragmentating complexes into the final rovibrational state (v′, j′), where the sum of the denominator runs over all accessible (v″, j″) channels at the energy E.

The capture probabilities, p^J_vj(E) and , are extracted by separate calculations on each arrangement (reactants and products, respectively) solving the corresponding coupled-channel equations by means of a time independent approach with a log-derivative propagation³⁸ between a capture radius of about 2.9 Å and an asymptotic distance of about 57.6 Å for the C⁺ + H₂ arrangement and between 2.1 Å and 42.4 Å for the H + CH⁺ arrangement. This calculation is performed by means of completely QM techniques under the centrifugal sudden approximation and using the ab initio full PES by Werfelli et al.⁴⁶ As an example, one of the computationally most demanding calculation considered here for reaction (1), with C⁺ + H₂(v = 2, j = 1) at E_c = 1 eV typically involves a maximum value of the total angular momentum of J = 82 and a maximum value of the helicity number Ω = 18 for the reactants arrangement and of Ω′ = 24 for the products arrangement.

The ICS are then calculated as follows:

where k_vj² = 2μ(E − E_vj)/ħ², with E_vj the energy of the initial rovibrational state of the reactant diatom and μ the reduced mass of the system.

ICS were calculated for collision energies sampled between 10⁻⁴ eV and 1.0 eV. With the ICS from eqn (5) we calculate the state-to-state rate constants as a function of the temperature as follows:

where β = (k_BT)⁻¹, k_B is the Boltzmann constant and Q_e(T) is the electronic partition function which, according to ref. 15, takes the following value for reaction in eqn (1):where 91.2 K is the spin–orbit splitting between the ²P_1/2 and ²P_3/2 states of C⁺(²P). For the other reactions in eqn (2) and (3), Q_e(T) = 1. With the maximum value of the collision energy chosen for the ICSs, rate constants are fully converged up to T = 1500 K.

2.2 Statistical adiabatic channel model

The SACM, originally introduced by Quack and Troe in the mid-1970s,^47,48 was recently modified and successfully applied to a number of complex-forming inelastic^33,49–52 and competing inelastic and reactive collisions.^37,53,54 The newer version relies on an ab initio PES and on the diagonalization of the interaction Hamiltonian excluding the nuclear kinetic term for the reactants and products arrangements considered separately. This version of SACM is thus a statistical method similar to the SQM. In the SACM approach, the state resolved integral scattering cross sections are computed according to eqn (5). However, the state-to-state transition probabilities in eqn (4) are computed assuming the capture probabilities to be zero (for closed channels) or unity (for open channels). They are thus of the form,where N(E, J) is the total number of open channels at total scattering energy, E, and total angular momentum, J. To count the aforementioned open channels, adiabatic channel potentials are computed; these are radially coupled potential curves that take into account the contribution of the electronic PES, as well as the effect of relative angular momentum through centrifugal barriers. A more detailed exposition of this modified SACM approach can be found in ref. 55 and 56. Despite its simple generalization beyond the strict atom–diatom processes, SACM does not naturally incorporate QM effects such as tunneling (something that has to be included ad hoc⁵⁴). We have therefore relied on the SQM technique for the production of rate constants for a large number of state-to-state processes for reactions in eqn (1)–(3). SACM predictions up to T = 500 K have been nevertheless obtained to further analyze the statistical character of the title reactions in comparison with the SQM.

The SACM ICSs for the reactions (1)–(3) mentioned in the introduction have been computed using an energy grid of 1000 logarithmically spaced energies in the collision energy interval [10⁻⁵ cm⁻¹ (1.24 × 10⁻⁹ eV), 16 130 cm⁻¹ (2 eV)]. Care was taken to make sure that the ICSs were fully converged (for the aforementioned energies) with respect to the total angular momentum quantum number and the number of rotational and vibrational states in both the CH⁺ + H and H₂ + C⁺ arrangements. Maximum values for such parameters considered in the present calculation are J^max_tot = 100; v^max_tot = 8 and j^max_tot = 40. The rovibrational energies of the diatomics in said arrangements were computed using the ABC program. State-to-state rate constants are then calculated based on these ICSs using eqn (6) in the temperature range 1–500 K.

3 Results and discussion

Both SQM and SACM calculations have been performed using the PES developed by Werfelli et al.⁴⁶ (hereafter WHHKS PES). This surface includes an improved description of the long range interaction in both reactants and products arrangements as compared to previous attempts such as SH PES,⁵⁷ WS PES⁵⁸ and LZH PES,¹⁶ to describe the interactions ruling the title process. This feature is crucial to reproduce the proper behaviour of the rate constants in the low temperature regime. A larger atomic basis set was used in the ab initio calculations compared to previously available PESs and a larger number of energy points were computed in order to get an accurate analytical representation. The obtained endothermicity, 0.402 eV, for reaction in eqn (1) was found in a much better agreement with the experiment, 0.393 eV, than some other previous surfaces. The existence of a global minimum of more than 4.0 eV with respect to the C⁺ + H₂ asymptote supports the importance of complex-forming mechanisms for the overall reaction. It is worth mentioning that some of calculations considered here as benchmarks for comparison with present statistical predictions have been performed using some of the above mentioned previous PESs.^16,57,58

3.1 The C⁺ + H₂ → H + CH⁺ reaction

The transitions considered for this reaction involves the initial H₂(v = 0, j = 0–13), H₂(v = 1, j = 0–10) and H₂(v = 2, j = 0–5) rovibrational states. We have calculated first ICSs for the C⁺ + H₂(v, j) → H + CH⁺(v′,j′) reaction following the expression in eqn (5). Results for the processes starting from H₂(v = 0, 1, j = 0) as a function of the collision energy E_c (for CH⁺ rotationally and vibrationally summed) are shown in Fig. 1 and 2 respectively, and in the SI Part A, similar cross sections for H₂(v = 0, 1, j = 1) can be found. In all cases we have compared both SQM and SACM predictions with the TDWP ICSs reported by Zanchet et al.¹⁵ calculated on the SH PES.⁵⁷

Fig. 1 Integral cross sections for the C⁺ + H₂(v = 0, j = 0) → H + CH⁺ reaction, obtained with the present SQM (dotted line) and SACM (dashed line) as a function of the collision energy compared with TDWP results (solid lines) from ref. 15 and with two sets of experimental results: (i) the numerical data given by Gerlich et al.²⁷ to fit their experimental results and (ii) the fit reported by González and Aguilar⁵⁹ to experimental measurements by Ervin et al.⁹ Notice that theoretical ICSs have been affected here by the value of the Q_e(T) electronic partition function at T = 300 K for the comparison with experiment.

Fig. 2 ICSs for the C⁺ + H₂(v = 1, j = 0) → H + CH⁺(v′ = 0–2) reaction.

The comparison in Fig. 1 certainly reveals, besides differences in the energy threshold for reaction due the different endothermicities of the corresponding PESs, that the statistical predictions overestimate the TDWP ICS beyond 0.5 eV for the C⁺ + H₂(v = 0, j = 0) → CH⁺ + H reaction. In particular, for the largest value of E_c considered here, 1 eV, the SACM cross section is about 1.7 times larger than the TDWP result and the SQM also exceeds such exact value by a 1.4 factor. For an analogous comparison between the present two statistical approaches and the TDWP cross sections for the final vibrationally selected process, C⁺ + H₂(v = 0, j = 0) → CH⁺(v′) + H, see SI Part A. The above mentioned shift of the threshold for reaction obtained here with the WHHKS PES with respect to the experiment is similar to the one observed in the TDWP study by Wu et al.⁶⁰ using the PES reported by Guo et al.⁶¹ However, the QCT calculations performed by Herráez-Aguilar et al.³² on the PES of ref. 58 and the TDWP cross sections by Zanchet et al.¹⁵ on the SH PES seem to bring theory closer to experiment.

Fig. 1 also includes the numerical fit provided by Gerlich et al.²⁷ for their experimental results. All theoretical cross sections shown in Fig. 1 have been multiplied by the value of the Q_e(T) electronic partition function at T = 300 K (0.407) for such comparison with the experiment. Whereas the statistical predictions constitute a fairly good reproduction of the experiment at low energies, SQM and especially SACM start to overestimate results from ref. 27 beyond 0.6 eV.

These observed deviations with respect to exact QM results are consistent with findings reported by³² when using the above mentioned SQCT method. As discussed in ref. 32 one possible reason for the overestimation of a statistical prediction comes for the fact that this kind of approaches assume that once the intermediate region is reached a sufficiently long-lived collision complex is formed. This is not always the case, especially at high energy, even if complex-forming pathways may play a role in the overall dynamics, either because the lifetime of such and intermediate species is not that long or due to the presence of possible barriers (either centrifugal or dynamical) which are actually not as successfully surmounted as a statistical treatment expects. For some other PESs developed for the present reaction, some minor intermediate barriers or transition states have been found which could be interpreted as possible obstacles preventing from a purely statistical behaviour. The WHHKS PES used here is reported to have a global minimum accessible from a reasonably wide range of angular range away from the purely H–C–H and C–H–H collinear arrangements. On the other hand, the SACM here employed considers the indirect path (formation of a complex and then dissociation) whereas the direct path (breakdown of CH bond) starts to be important at higher temperature.¹¹ But, as we will show here it is possible nevertheless to use these ICSs to extract rate constants in good accord with exact QM results, thus validating the statistical description of the reaction. It is worth mentioning here in support of the capabilities of statistical techniques to address the study of the C⁺ + H₂ reaction, that DCSs obtained with the SQM (not shown here) at some of the collision energies studied by Li et al.¹⁸ employing the LZH PES are in a roughly fair agreement with the TDWP results reported in Fig. 6 of their work.¹⁸

Present statistical predictions have been also tested against results from the TDWP investigation by Zanchet et al.¹⁵ for the cases in which the reaction initiates from the H₂(v = 1, j = 0, 1) rovibrational states. ICSs for these two possible initial rotational states (v = 1, j = 0) and (v = 1, j = 1) are shown in Fig. 2 and in the SI Part A, respectively. SQM and SACM cross sections for the formation of CH⁺(v′ = 0) in its ground vibrational state are found to underestimate (overestimate) the TDWP cross sections at the low (high) energy regime but the agreement turns reasonably good along the range between 10⁻² and 0.5 eV. It is hard to establish a definitive explanation for the differences observed between both types of results, but it is fair to say that both statistical approaches yield the same behavior with respect to the energy and that, in fact, TDWP cross sections deviate somehow from the expected trend in a log-scale representation. ICSs for the final v′ = 1 and v′ = 2 vibrational states display threshold behaviours which are properly described by SQM and SACM results. The statistical values for v′ = 1 start to slightly overestimate the QM ICSs as seen for the C⁺ + H₂(v = 0, j = 0) case beyond 0.4 eV collision energy.

SQM rate constants were obtained through eqn (6). The comparison with the TDWP results from Zanchet et al.¹⁵ for the C⁺ + H₂(v = 0, j = 0) → CH⁺(v′ = 0, j′ = 0–2) + H reaction can be seen in Fig. 3 up to T = 2000 K. In the inset we include the SACM rates up to T = 500 K for the cases j′ = 0 and j′ = 2. A more extended comparison is shown in the SI Part A. Results shown in Fig. 3 reveal that the statistical predictions certainly constitute a valid counterpart, providing slightly larger values than the QM TDWP values at the entire temperature range under study. The comparison between the TDWP rates from Zanchet et al.¹⁵ and present statistical values for T = 300, 400 and 500 K are listed in Table 1.

Fig. 3 State-to-state rate constants for the C⁺ + H₂(v = 0, j = 0) → H + CH⁺(v′ = 0, j′ = 0–2) reaction as a function of the temperature. Results for CH(v′ = 0) (rotationally summed) are shown in black. Present SQM results (dotted lines) are compared with TDWP rates (solid lines) from ref. 15. SACM results (dashed lines) are included up to T = 500 K in the inset for the cases j′ = 0 and 2. Further comparisons are available in SI Part A.

Table 1 Rate constants (in 10^−x cm³ s⁻¹ with x between parenthesis) for the C⁺ + H₂(v = 0, j = 0) → H + CH⁺(v′ = 0, j′ = 0–2)

T (K)	SQM	SACM	TDWP
j′ = 0
200	1.39(−21)	1.39(−21)	7.29(−22)
300	2.87(−18)	2.75(−18)	1.59(−18)
400	1.26(−16)	1.20(−16)	7.39(−17)
500	1.20(−15)	1.14(−15)	7.34(−16)
j′ = 1
200	2.31(−21)	3.31(−21)	1.16(−21)
300	4.93(−18)	6.80(−18)	2.75(−18)
400	2.21(−16)	3.03(−16)	1.32(−16)
500	2.12(−15)	2.92(−15)	1.33(−15)
j′ = 2
200	1.96(−21)	3.18(−21)	1.08(−21)
300	4.94(−18)	7.58(−18)	3.03(−18)
400	2.41(−16)	3.64(−16)	1.59(−16)
500	2.43(−15)	3.67(−15)	1.69(−15)

---

The experimentally derived value of the rate constant at T = 300 K (1.2 × 10⁻¹⁶ cm³ s⁻¹) by Ervin and Armentrout⁹ compares well with the corresponding thermal SQM rate constant at that temperature (8.9 × 10⁻¹⁷ cm³ s⁻¹). A comparison of the SQM rates for the C⁺ + H₂(v = 0, j = 0) → CH⁺(v′ = 0) + H process and the thermal average obtained with the first rotational states of H₂(v = 0, j = 0–12) following the expression:

where and E_0j are the energies of the rovibrational states H₂(v = 0, j), with the rates measured by Hierl et al.¹⁰ in their high temperature flowing afterglow experiment is shown in Fig. 4. The value of the statistical rate at T = 800 K is 4.67 × 10⁻¹³ cm³ s⁻¹ which remains below the value obtained according to the Arrhenius expression proposed in ref. 10: 2.60 × 10⁻¹² cm³ s⁻¹. The same underestimation of the experimental rates was seen by Li et al.¹⁴ whose TDWP results using their own new PES (also shown in Fig. 4) are in good accord with the present statistical predictions. A similar result can be found in the TDWP investigation made by Guo et al.⁶² (see their Fig. 9, where the experimental results appear to be misplotted) employing the LZH PES.

Fig. 4 Comparison between the presently calculated rate constants for the C⁺ + H₂(v = 0, j = 0) → H + CH⁺(v′ = 0) reaction with the SQM (black dashed line) and the thermal average considering the first (v = 0, j = 0–8) 9 rotational states (black solid line) of H₂ with: the experimental results from the analytical fit k(T) = Ae^−E_a/RT with A = (7.4 ± 0.8) × 10⁻¹⁰ cm³ s⁻¹ and E_a = 0.391 ± 0.005 eV proposed in Hierl et al., J. Chem. Phys., 1997, 106, 10145 (blue triangles) and the TDWP calculation by Li et al.¹⁴ (red solid line).

When we extend a similar comparison for transitions from the vibrationally excited H₂(v = 1, j = 0), the situation remains analogously favourable for a statistical description. The comparison with TDWP results reported by Zanchet et al.¹⁵ for temperatures up to 1500 K is shown in Fig. 5 and reveals that a remarkably good accord between SQM and TDWP rates for C⁺ + H₂(v = 1, j = 0) → CH⁺(v′ = 0, 1) + H. Authors of ref. 15 mentioned that their theoretical value for the rate constant for H₂(v = 1) at T = 800 K, 3 × 10⁻¹⁰ cm³ s⁻¹, was 2–3 times smaller than the experimental estimation by Hierl et al.,¹⁰ 1–2 × 10⁻⁹ cm³ s⁻¹. Present SQM prediction, 3.4 × 10⁻¹⁰ cm³ s⁻¹, in extremely good accord with the TDWP, also overestimates the reported measurement. SACM predictions, here shown up to T = 500 K, exhibit some deviations in the case of formation of vibrationless CH⁺(v′ = 0).

Fig. 5 Final vibrational state resolved rate constants for the C⁺ + H₂(v = 1, j = 0) → H + CH⁺(v′ = 0, 1) reaction. Present SQM (red dotted lines) and SACM (blue dashed lines) results are compared with rates obtained by the TDWP calculation from Zanchet et al.¹⁵ (in black solid lines).

It is also possible to test the capability of the presently employed statistical techniques at a complete state-to-state level, specifying also the final CH⁺(v′, j′) rotational state, as in Fig. 3. In Fig. 6, the SQM rate constants for the C⁺ + H₂(v = 1, j = 0) → CH⁺(v′ = 0, j′ = 1, 3, 4, 6, 7) + H are compared with the corresponding results from the TDWP study by Zanchet et al.¹⁵ up to T = 2000 K. The same comparison for some other rotational states j′ is shown in the SI Part A. SACM results up to 500 K are also included in the same figure.

Fig. 6 Rate constants as a function of the temperature for the C⁺ + H₂(v = 1, j = 0) → CH⁺(v′ = 0, j′) + H reaction for a variety of different final j′ rotational states (see the SI Part A for some other j′ values). Present SQM (dotted lines) and SACM (dashed lines) are compared with those from ref. 15.

Although the agreement between the TDWP rates and the statistical predictions is not perfect, the figure reveals the capability of these approaches to mimic the overall behaviour of the QM results at the state-to-state level, especially below 200 K. At larger temperatures, the SQM and SACM methods yield values for k(T) which slightly overestimate the QM results. The only noticeable differences with respect to the TDWP are seen for j′ = 1, but according with the scale used in Fig. 6, one can conclude that the statistical rate constants constitute a reasonably good alternative.

Values for the rate constants for the C⁺ + H₂(v = 1, j = 0) → CH⁺(v′ = 1, j′) + H reaction, involving the formation of the CH⁺ species at its first vibrationally excited state are noticeably smaller. Whereas at T = 500 K, k_10,0j′(T) remain typically around ∼10⁻¹¹ cm³ s⁻¹, values of k_10,1j′(T), shown in Fig. 7, do not exceed ∼10⁻¹⁴ cm³ s⁻¹. In that Fig. 7, SQM predictions are compared with the TDWP rates from ref. 15 up to 1500 K. The agreement between both sets of data is very good. In the inset of the figure magnifying the region up to 500 K we have also included the SACM rate constants which are in excellent accord with the QM results.

Fig. 7 Rate constants as a function of the temperature for the C⁺ + H₂(v = 1, j = 0) → CH⁺(v′ = 1, j′ = 0, 4, 8) + H reaction. Present SQM (dotted lines) and SACM (dashed lines) are compared with those from ref. 15. The inset shows the temperature range up to 500 K.

It is worth then pointing out that the above comparisons reveal that present statistical predictions agree with previous quantum mechanical results and experimental data to within a factor of 2 and are thus well suited for astrophysical applications.

3.2 The H + CH⁺ reaction

We have considered both the abstraction arrangement shown in eqn (2) leading to the formation of H₂ (see Section 3.2.1) and the exchange process given in eqn (3) which yields H + CH⁺ (see Section 3.2.2).

3.2.1 The abstraction H + CH⁺ → C⁺ + H₂ reaction. Initial state selected rate constants for the abstraction reaction shown in eqn (2) to produce molecular hydrogen in different final vibrational states, H⁺ + CH⁺(v, j) → H₂(v′), have been calculated using the nuclear-spin multiplicity factors for the final H₂(v′, j′) rotational states according to the expression:where when j is even and when j is odd. Numerical values of the rate constants obtained with the SQM for processes initiated from CH⁺(v = 0, j = 0–20), CH⁺(v = 1, j = 0–14) and CH⁺(v = 2, j = 0–5) are shown in the SI Part B.

In Fig. 8 initial state selected rate constants k_v=0,j(T) (summing on all final rovibrational states (v′, j′)) are shown at T = 500 K as a function of the initial rotational state j. The present SQM and SACM rate constants are compared to those from Faure et al.⁶³ obtained by means of the time-independent QM (TIQM) ABC code⁶⁴ for j = 0–7 and QCT calculations for j = 8–13 using the same WHHKS PES as in this work. The comparison reveals that all rate constants agree within a reasonably good range, being the largest differences (less than a factor of 2) between both statistical predictions and the QM values observed at j = 13 for the SQM and at j = 9 for the SACM. Whereas the SACM results remain larger (more than 1.3 times) for all the rotational states under considerations, the SQM values are smaller for the first states (j < 6) to slightly overestimate the QM rates for the rest of rotational states.

Fig. 8 Rate constants for the CH⁺(v = 0, j) → C⁺ + H₂ abstraction reaction. Present SQM results (black circles) are compared with SACM results (empty triangles) and those from Faure et al.⁶³ (empty squares).

In order to establish a comparison with the various measured values of the rate constant existing in the literature for the H + CH⁺ → C⁺ + H₂ reaction, in this work we have obtained the thermal rate by the corresponding Boltzmann average of the rate constants in eqn (10) following the expression in eqn (9) considering the first rovibrational states CH⁺(v = 0, j: 0–14). The comparison of the so extracted SQM thermal rate coefficient and those reported in the experimental works by Federer et al.,⁶⁵ Plasil et al.¹² and Luca et al.⁶⁶ is shown in Fig. 9. A fairly good agreement is seen between 50 and 100 K with measurements from ref. 65 and 66. Below T = 50 K only some of the SQM rate constants for transitions from large rovibrational initial states CH⁺(v = 0, j ≥ 8) (see for instance the case of j = 14 included in the figure) seem to describe the drastic decreasing trend observed in the experiment by Plasil et al.¹² as the temperature decreases. Contributions from vibrationally excited states such as (v = 1, 2, j = 0), not included in the figure, reveals that they also correspond to less reactive processes. In this sense, present results are consistent with findings from the recent study of ref. 13 where vibrational excitation of CH⁺ was found to play a major role in reducing the reactivity at low temperatures. Fig. 9 also includes the TIQM + QCT results up to 1500 K by Faure et al.⁶³ This calculation is shown to agree with our SQM thermal rate constant to better than a factor of 2. Interestingly, while our SACM predictions agree well with rate constants obtained with the adiabatic statistical approach tested by del Mazo-Sevillano et al.,¹³ as expected, the present SQM results look more similar to the QCT calculations including surface hopping from that study.

Fig. 9 Rate constants for the H + CH⁺(v, j) → C⁺ + H₂ process. Present thermal SQM and SACM results for different CH⁺(v, j) initial rovibrational states and for the thermal rate constant (eqn (9)) are compared with experimental results from Plasil et al.¹² (solid red squares), Federer et al.⁶⁵ (empty blue squares) and Luca et al.⁶⁶ (solid green triangles). The TIQM + QCT thermal rate constant from ref. 63 has been included for comparison.

3.2.2 The exchange H + CH⁺ → H + CH⁺ reaction. The present study also includes the exchange (also called inelastic) reaction H + CH⁺ → H + CH⁺.

Numerical values of the rate constants up to 1500 K obtained with the SQM are presented in the SI Part C for the transitions leading to the formation of CH⁺(v′, j′) when the reaction is initiated from CH⁺(v, j) with E_vj > E_v′j′ for (v = 0, j = 1–20), (v = 1, j = 0–14) and (v = 2, j = 0–5). However, in order to test the performance of the present statistical techniques in the study of the title reaction, in Fig. 10 we also present rate constants for rotational transitions from the rovibrational ground state (v = 0, j = 0). It is worth reminding two things about the TIQM results obtained with the ABC code reported in the investigation by Faure et al.⁶³ used here as a possible benchmark. On the one hand, as mentioned in ref. 46, the reliability of the ABC code for reactions involving a deep potential well has been questioned and, in fact, due to numerical cost TIQM results by Faure et al. had to be completed with QCT calculations. On the other hand, the method includes both the inelastic and exchange channels. Statistical methods here employed only accounts for the reactive or insertion process which leads from H + CH⁺ reactants to CH⁺ + H products: we have not included an extra factor 2 in the SACM approach to describe inelastic processes. Therefore, any deviation between outputs from these two approaches and the results from ref. 63 can be interpreted as indication of the relevance of any pathway which does not involve the formation of a long-living intermediate complex, such as for example purely inelastic collisions. Nevertheless, both SQM and SACM predictions are found to be of the same relative order as the QM rates. Differences which depend both on the specific H + CH⁺(v = 0, j = 0) → H + CH⁺(v′ = 0, j′) transition and on the specific temperature are observed.

Fig. 10 Rate constants for the exchange H + CH⁺(v = 0, j = 0) → H⁺ + CH(v′ = 0, j′ = 0–4) process. Present SQM (dotted lines) and SACM (dashed lines) are compared with TIQM results (solid lines) from Faure et al.⁶³

Thus, for example, the SQM results are almost in perfect agreement with the TIQM rate constants for j′ = 1 up to 200 K. At larger temperatures the statistical values start to smoothly decrease as opposed to the increasing trend observed for the QM rates. In turn, as the rotational excitation of the product CH⁺(v′ = 0, j′) increases SQM rates are found to slightly overestimate the TIQM results. SACM predictions are found to be consistently smaller than the SQM ones which makes their comparison with the TIQM values slightly worse for j′ = 1 and 2 and slightly better for j′ = 3 and 4.

Also at a state-to-state level, but now fixing the final rovibrational state, we have investigated rate constants corresponding to the production of CH⁺(v′ = 0, j′ = 0) and CH⁺(v′ = 0, j′ = 1). Results for both cases are presented in Fig. 11 and 12, respectively. Rate constants for these two types of processes display a different profile with respect to the previously investigated transitions (see Fig. 10) where the corresponding k(T) increase rapidly in the low temperature range. TIQM rates from Fig. 11 and 12, on the contrary, exhibit a more stable and constant trend as a function of the temperature. Despite the variety observed among all cases here considered, it is worth mentioning the noticeably good agreement with the TIQM rate constants observed for the SQM k_02,00 rate constant between 30 and 500 K, and the k_05,01 rate constant for the entire range up to 1500 K. The SACM calculation also yields k_01,00 and k_02,00 rate constants in specially good accord with the QM rates up to 100 K.

Fig. 11 Rate constants for the relaxation to the ground CH⁺(v′ = 0, j′ = 0) state in the exchange H + CH⁺(v = 0, j = 1–4) → CH⁺(v′ = 0, j′ = 0) reaction. Solid lines are TIQM results from Faure et al.,⁶³ dashed lines are present SACM rates and dotted lines are present SQM values.

Fig. 12 Same as Fig. 11 for the exchange H + CH⁺(v = 0, j = 2–5) → CH⁺(v′ = 0, j′ = 1) reaction.

Conclusions

In this work we have shown the capabilities of statistical approaches to provide useful and accurate information about the reactive collisions of the CH₂⁺ system. Cross sections as a function of the energy and rate constants as a function of the temperature have been obtained by means of a statistical quantum method (SQM) within the centrifugal sudden approximation and a statistical adiabatic channel model (SACM) for the C⁺(²P) + H₂ → H + CH⁺, H + CH⁺ → H + CH⁺ and H + CH⁺ → C⁺(²P) + H₂ reactions at the state-to-state level. The comparison with previously reported quantum mechanical results and experimental data reveals a remarkably good agreement, typically within a factor of 2 in the worst cases, especially for C⁺ + H₂, with in general a slight superiority of the SQM over the SACM, as expected. Thus, in practice, SQM has been applied to calculate rate constants for a large number of specific transitions in each case, in an effort to provide a useful database for future astrophysical applications, as in previous investigations on the H₃⁺ system and isotopic variants. The applicability of results presented here is illustrated through work in progress where the SQM rate constants are being used to model the non-thermal spectrum of CH⁺ towards the planetary nebula NGC 7027 (Sil et al.⁶⁷). The present set of data, which covers all CH⁺ states up to v = 2, j = 5, is indeed fully appropriate to model the 6 pure rotational and 14 ro-vibrational (v = 1 → 0) lines detected by Cernicharo et al. (1997)² and Neufeld et al. (2021),³ respectively, in this astrophysical source (a gas shell expelled by a red giant star). The next target will be the d203–506 protoplanetary disk where the v = 2 → 1 vibrational band was also identified (Zannese et al. (2025)⁷).

Author contributions

Tomás González Lezana: conceptualization (lead); data curation (lead); formal analysis (lead): funding acquisition (equal); investigation (lead); methodology (lead); resources (equal); software (lead); validation (equal); writing – original draft (lead); writing – review and editing (equal). Maarten Konings: conceptualization (supporting); data curation (equal); formal analysis (supporting): funding acquisition (equal); investigation (equal); methodology (equal); resources (equal); software (equal); validation (equal); writing – original draft (equal); writing – review and editing (equal). Jérôme Loreau: conceptualization (supporting); data curation (equal); formal analysis (supporting): funding acquisition (equal); investigation (supporting); methodology (equal); resources (equal); software (equal); validation (equal); writing – original draft (equal); writing – review and editing (equal). François Lique: conceptualization (supporting); data curation (supporting); formal analysis (supporting): funding acquisition (lead); investigation (supporting); methodology (supporting); resources (supporting); software (supporting); validation (supporting); writing – original draft (supporting); writing – review and editing (equal). Milan Sil: conceptualization (supporting); data curation (supporting); formal analysis (supporting): funding acquisition (equal); investigation (supporting); methodology (supporting); resources (supporting); software (supporting); validation (supporting); writing – original draft (equal); writing – review and editing (equal). Alexandre Faure: conceptualization (equal); data curation (supporting); formal analysis (supporting): funding acquisition (lead); investigation (supporting); methodology (supporting); resources (supporting); software (supporting); validation (supporting); writing – original draft (supporting); writing – review and editing (lead).

Conflicts of interest

The authors do not have conflicts to declare.

Data availability

This work is accompanied by supplemental information (SI) Part A (with figures and numerical data for reaction rate constants for the C⁺ + H₂ → CH⁺ + H reaction); SI Part B (with numerical data for reaction rate constants for the CH⁺ + H → C⁺ + H₂ reaction); and SI Part C (with numerical data for reaction rate constants for the CH⁺ + H → CH⁺ + H reaction). Supplementary information is available. See DOI: https://doi.org/10.1039/d5cp04254b.

Acknowledgements

Authors want to thank Alex Zanchet and Octavio Roncero for useful discussions and for providing the TDWP results for Section 3.1 and to Wentao Li for the numerical data for Fig. 4. M. S. acknowledges financial support from the European Research Council (consolidated grant COLLEXISM, grant agreement ID: 811363) T. G. L. acknowledges funding from Grant No. MICIU/AIE/10.13039/501100011033 PID2024-155666NB-I00.

References

1. A. E. Douglas and G. Herzberg, Astrophys. J., 1941, 94, 381.
2. J. Cernicharo, X. W. Liu, E. González-Alfonso, P. Cox, M. J. Barlow, T. Lim and B. M. Swinyard, Astrophys. J., 1997, 483, L65–L68.
3. D. A. Neufeld, B. Godard, P. Bryan Changala, A. Faure, T. R. Geballe, R. Güsten, K. M. Menten and H. Wiesemeyer, Astrophys. J., 2021, 917, 15.
4. E. Falgarone, B. Godard, J. Cernicharo, M. de Luca, M. Gerin, T. G. Phillips, J. H. Black, D. C. Lis, T. A. Bell, F. Boulanger, A. Coutens, E. Dartois, P. Encrenaz, T. Giesen, J. R. Goicoechea, P. F. Goldsmith, H. Gupta, C. Gry, P. Hennebelle, E. Herbst, P. Hily-Blant, C. Joblin, M. Kaźmierczak, R. Kołos, J. Krełowski, J. Martin-Pintado, R. Monje, B. Mookerjea, D. A. Neufeld, M. Perault, J. C. Pearson, C. Persson, R. Plume, M. Salez, M. Schmidt, P. Sonnentrucker, J. Stutzki, D. Teyssier, C. Vastel, S. Yu, K. Menten, T. R. Geballe, S. Schlemmer, R. Shipman, A. G. G. M. Tielens, S. Philipp, A. Cros, J. Zmuidzinas, L. A. Samoska, K. Klein, A. Lorenzani, R. Szczerba, I. Péron, P. Cais, P. Gaufre, A. Cros, L. Ravera, P. Morris, S. Lord and P. Planesas, Astron. Astrophys., 2010, 521, L15.
5. W. F. Thi, F. Ménard, G. Meeus, C. Martin-Zadi, P. Woitke, E. Tatulli, M. Benisty, I. Kamp, I. Pascucci, C. Pinte, C. A. Grady, S. Brittain, G. J. White, C. D. Howard, G. Sandell and C. Eiroa, Astron. Astrophys., 2011, 530, L2.
6. N. Rangwala, P. R. Maloney, J. Glenn, C. D. Wilson, J. Kamenetzky, M. R. P. Schirm, L. Spinoglio and M. Pereira Santaella, Astrophys. J., 2014, 788, 147.
7. M. Zannese, B. Tabone, E. Habart, E. Dartois, J. R. Goicoechea, L. Coudert, B. Gans, M. A. Martin-Drumel, U. Jacovella, A. Faure, B. Godard, A. G. G. M. Tielens, R. Le Gal, J. H. Black, S. Vicente, O. Berné, E. Peeters, D. Van De Putte, R. Chown, A. Sidhu, I. Schroetter, A. Canin and O. Kannavou, Astron. Astrophys., 2025, 696, A99.
8. K. M. Ervin and P. B. Armentrout, J. Chem. Phys., 1984, 80, 2978–2980.
9. K. M. Ervin and P. B. Armentrout, J. Chem. Phys., 1986, 84, 6738–6749.
10. P. M. Hierl, R. A. Morris and A. A. Viggiano, J. Chem. Phys., 1997, 106, 10145–10152.
11. R. A. Jara-Toro, O. Roncero and F. Lique, Phys. Chem. Chem. Phys., 2024, 26, 21370–21378.
12. R. Plasil, T. Mehner, P. Dohnal, T. Kotrik, J. Glosik and D. Gerlich, Astrophys. J., 2011, 737, 60.
13. P. del Mazo-Sevillano, R. A. Jara-Toro, O. Roncero and F. Lique, Phys. Chem. Chem. Phys., 2025, 27, 15775.
14. W. Li, B. Dong, X. Niu, M. Wang and Y. Zhang, J. Chem. Phys., 2024, 161, 074302.
15. A. Zanchet, B. Godard, N. Bulut, O. Roncero, P. Halvick and J. Cernicharo, Astrophys. J., 2013, 766, 80.
16. Y. Q. Li, P. Y. Zhang and K. L. Han, J. Chem. Phys., 2015, 142, 124302.
17. P. Sundaram, V. Manivannan and R. Padmanaban, Phys. Chem. Chem. Phys., 2017, 19, 20172–20187.
18. W. Li, L. Liu and Z. Zhu, Chem. Phys., 2022, 562, 111670.
19. B. H. Mahan and T. M. Sloane, J. Chem. Phys., 1973, 59, 5661–5675.
20. I. Maier and B. William, J. Chem. Phys., 1967, 46, 4991–4992.
21. P. Pechukas and J. C. Light, J. Chem. Phys., 1965, 42, 3281–3291.
22. D. G. Truhlar, J. Chem. Phys., 1969, 51, 4617–4623.
23. W. J. Chesnavich and M. T. Bowers, J. Chem. Phys., 1978, 68, 901–910.
24. E. Herbst and S. K. Knudson, Chem. Phys., 1981, 55, 293–297.
25. D. A. Webb and W. J. Chesnavich, J. Phys. Chem., 1983, 87, 3791–3798.
26. K. M. Ervin and P. B. Armentrout, J. Chem. Phys., 1986, 84, 6750–6760.
27. D. Gerlich, R. Disch and S. Scherbarth, J. Chem. Phys., 1987, 87, 350–359.
28. P. Halvick, T. Stoecklin, P. Larrégaray and L. Bonnet, Phys. Chem. Chem. Phys., 2007, 9, 582–590.
29. F. J. Aoiz, V. Sáez Rábanos, T. González-Lezana and D. E. Manolopoulos, J. Chem. Phys., 2007, 126, 161101.
30. F. J. Aoiz, T. González-Lezana and V. Sáez Rábanos, J. Chem. Phys., 2007, 127, 174109.
31. F. J. Aoiz, T. González-Lezana and V. S. Rábanos, J. Chem. Phys., 2008, 129, 094305.
32. D. Herráez-Aguilar, P. G. Jambrina, M. Menéndez, J. Aldegunde, R. Warmbier and F. J. Aoiz, Phys. Chem. Chem. Phys., 2014, 16, 24800–24812.
33. J. Loreau, F. Lique and A. Faure, Astrophys. J., Lett., 2018, 853, L5.
34. M. Quack and J. Troe, Ber. Bunsen-Ges. Phys. Chem., 1974, 78, 240–252.
35. M. Quack and J. Troe, Ber. Bunsen-Ges. Phys. Chem., 1975, 79, 170–183.
36. M. Quack and J. Troe, Ber. Bunsen-Ges. Phys. Chem., 1976, 80, 1140–1149.
37. M. Konings, B. Desrousseaux, F. Lique and J. Loreau, J. Chem. Phys., 2021, 155, 104302.
38. E. J. Rackham, T. González-Lezana and D. E. Manolopoulos, J. Chem. Phys., 2003, 119, 12895.
39. T. González-Lezana, Int. Rev. Phys. Chem., 2007, 26, 29.
40. P. Honvault, M. Jorfi, T. González-Lezana, A. Faure and L. Pagani, Phys. Rev. Lett., 2011, 107, 023201.
41. P. Honvault, M. Jorfi, T. González-Lezana, A. Faure and L. Pagani, Phys. Rev. Lett., 2012, 108, 109903.
42. T. González-Lezana, P. Hily-Blant and A. Faure, J. Chem. Phys., 2021, 154, 054310.
43. T. González-Lezana, P. Hily-Blant and A. Faure, J. Chem. Phys., 2022, 157, 214302.
44. D. R. Flower, G. Pineau des Forêts, P. Hily-Blant, A. Faure, F. Lique and T. González-Lezana, Mon. Not. R. Astron. Soc., 2021, 507, 3564–3571.
45. A. Faure, P. Hily-Blant, G. Pineau des Forêts and D. R. Flower, Mon. Not. R. Astron. Soc., 2024, 531, 340–354.
46. G. Werfelli, P. Halvick, P. Honvault, B. Kerkeni and T. Stoecklin, J. Chem. Phys., 2015, 143, 114304.
47. M. Quack and J. Troe, Ber. Bunsen-Ges. Phys. Chem., 1974, 78, 240–252.
48. M. Quack and J. Troe, Ber. Bunsen-Ges. Phys. Chem., 1975, 79, 170–183.
49. J. Loreau, A. Faure and F. Lique, J. Chem. Phys., 2018, 148, 244308.
50. A. Faure, F. Lique and J. Loreau, Mon. Not. R. Astron. Soc., 2020, 493, 776–782.
51. C. Balança, Y. Scribano, J. Loreau, F. Lique and N. Feautrier, Mon. Not. R. Astron. Soc., 2020, 495, 2524–2530.
52. J. Loreau, A. Faure and F. Lique, Mon. Not. R. Astron. Soc., 2022, 516, 5964–5971.
53. B. Desrousseaux, M. Konings, J. Loreau and F. Lique, Phys. Chem. Chem. Phys., 2021, 23, 19202–19208.
54. M. Konings, T. González-Lezana, S. Camps and J. Loreau, Phys. Chem. Chem. Phys., 2024, 26, 22463–22471.
55. J. Loreau, F. Lique and A. Faure, Astrophys. J., 2018, 853, L5.
56. M. Konings, B. Desrousseaux, F. Lique and J. Loreau, J. Chem. Phys., 2021, 155, 104302.
57. T. Stoecklin and P. Halvick, Phys. Chem. Chem. Phys., 2005, 7, 2446–2452.
58. R. Warmbier and R. Schneider, Phys. Chem. Chem. Phys., 2011, 13, 10285–10294.
59. M. González and A. Aguilar, Chem. Phys., 1989, 132, 443–462.
60. H. Wu, Z. Duan and X. He, Astrophys. J., 2021, 906, 117.
61. L. Guo, H. Ma, L. Zhang, Y. Song and Y. Li, RSC Adv., 2018, 8, 13635–13642.
62. J. Guo, A. J. Zhang, Y. Zhou, J. Y. Liu, J. F. Jia and H. S. Wu, Chem. Phys. Lett., 2017, 689, 121–127.
63. A. Faure, P. Halvick, T. Stoecklin, P. Honvault, M. D. Epée Epée, J. Z. Mezei, O. Motapon, I. F. Schneider, J. Tennyson, O. Roncero, N. Bulut and A. Zanchet, Mon. Not. R. Astron. Soc., 2017, 469, 612–620.
64. D. Skouteris, J. Castillo and D. Manolopoulos, Comput. Phys. Commun., 2000, 133, 128–135.
65. W. Federer, H. Villinger, F. Howorka, W. Lindinger, P. Tosi, D. Bassi and E. Ferguson, Phys. Rev. Lett., 1984, 52, 2084–2086.
66. T. Luca, G. Borodi and D. Gerlich, in Photonic, Electronic and Atomic Collisions, ed. P. D. Fainstein, M. A. P. Lima, J. E. Miraglia, E. C. Montenegro and R. D. Rivarola, World Scientific Publishing Co. Pvt. Ltd., Singapore, 2006, p. 494.
67. M. Sil, A. Faure, H. Wiesemeyer, P. Hily-Blant, T. González-Lezana, J. Forer, J. Loreau and F. Lique, Astron. Astrophys., 2026, submitted.

---