Quantum and statistical state-to-state studies of cold Ar + H₂⁺ collisions

Abstract

In this work we present new state-to-state integral scattering cross sections and initial-state selected rate coefficients for the ³⁶Ar (¹S) + H₂⁺ (X²Σ_g⁺,v = 0,j) reactive system for collision energies up to 0.1 eV (with respect to the ³⁶Ar (¹S) + H₂⁺ (X²Σ_g⁺,v = 0,j = 0) channel). To the best of our knowledge, these cross sections are the first fully state resolved ones that were obtained by performing time-independent quantum mechanical and quantum statistical calculations. For this purpose a new full-dimensional ground state ²A′ adiabatic electronic potential energy surface was calculated at the MRCI+Q/aug-cc-pVQZ level of theory, which was fitted by means of machine learning methods. We find that a statistical quantum method and a statistical adiabatic channel model reproduce quantum mechanical initial-state selected cross sections fairly well, thus suggesting that complex-forming mechanisms seem to be playing an important role in the reaction dynamics of the reaction that was studied.

---

1 Introduction

The exothermic reactions,

Ar (¹S) + H₂⁺ (X²Σ_g⁺) → ArH⁺ (X¹Σ⁺) + H (²S) ΔE ≈ −1.30 eV, (1)

Ar⁺ (²P) + H₂ (X¹Σ_g⁺) → ArH⁺ (X¹Σ⁺) + H (²S) ΔE ≈ −1.56 eV, (2)

represent the two main formation routes^1–5 of the argonium molecular cation, ArH⁺, and specifically the ³⁶ArH⁺ isotopic variant, in space. Indeed, the latter has been observed via rotational emissions in the interstellar medium (ISM), among other astronomical environments,^5–9 and its (quantum state resolved) formation rates are thus of interest to the astrochemistry and astrophysics communities. In particular, an accurate thermal rate coefficient for reaction (1) would be useful since, apart from an estimate (∼ 10⁻⁹ cm³ s⁻¹) in ref. 5, it is not available at time of writing. Indeed, to the best of our knowledge there are no theoretical, nor experimental results available for the rate coefficient for reaction (1) involving the ³⁶ArH⁺ isotopologue.

The electronic potential energy surfaces (PESs) for the reactions (1) and (2) are strongly coupled in the entrance arrangements. In fact, the PES for the first reaction is the only reactive one (directly correlating to the ArH⁺ (X¹Σ⁺) product), and the Ar⁺ + H₂ reactants only react to form ArH⁺ + H by an electronically non-adiabatic process.

A lot of our understanding of the reactions (1) and (2) comes from several studies^1–4,10,11 that have primarily focussed on the nonadiabatic reaction dynamics involving the ⁴⁰Ar isotope. Such theoretical approaches were usually based on DIMZO (diatomics in molecule with zero overlap) potential energy surfaces, as first introduced by Kuntz and Roach.¹²

In a more recent theoretical study from 2013,¹³ electronically adiabatic time-dependent quantum wave packet calculations were performed on an analytical fit of the ground state ab initio PES to obtain scattering cross sections for the reaction (1) with H₂⁺ in its ground ro-vibrational state (v_H2⁺ = 0,j_H2⁺ = 0), and specifically for the atomic collider ⁴⁰Ar; the maximum collision energy was constrained to 1.0 eV.

Similarly to the study by Hu et al.,¹³ the focus in the present work is on the reaction (1). However, our aim is twofold. (i) On the one hand we introduce the first quantum mechanical fully state resolved (i.e., state-to-state) integral scattering cross sections (ICSs) for collision energies up to 0.1 eV for the competing inelastic and reactive scattering events involving the collision partners ³⁶Ar and H₂⁺. Restriction to the aforementioned collision energies ensures that the collision happens electronically adiabatically since the lowest two excited states of the ³⁶ArH₂⁺ system are not accessible in this energy regime, thereby justifying the use of a single adiabatic PES.¹⁴ As shown in a recent work,¹⁴ above 0.1 eV non-adiabatic effects come into play, thereby strongly complicating the dynamics of reaction (1).

(ii) On the other hand we wish to investigate the statistical nature of reaction (1) and the applicability of statistical approaches, such as the statistical quantum method (SQM)^15,16 and the statistical adiabatic channel model (SACM),^17,18 that are commonly used in molecular collision theory. The reason being that the PES for reaction (1) (and also (2)) in ref. 14 shows the presence of an energy well separating reactants from products, thus justifying the use of these techniques, as seen in previous investigations where different dynamical features have been successfully reproduced by statistical means.^17–19

To accomplish said objectives, a new full-dimensional ab initio ground state PES was computed and fitted by means of machine learning (ML) techniques.¹⁴ In particular Gaussian process (GP) regression was chosen for this purpose, which represents a kernel-based non-parametric fitting strategy that has been used for the purposes of fitting electronic potential energy surfaces (PESs) of molecular systems.^14,20–26 Using this new PES, accurate time-independent quantum mechanical close-coupling (CC) calculations were carried out (related to objective (i)), as well as the aforementioned SQM and much less computationally expensive SACM calculations (objective (ii)).

This article is structured as follows: Section 2 contains a discussion about the use of Gaussian processes for fitting PESs, as well as a discussion about the details of the ab initio electronic structure calculations in relation to the ³⁶ArH₂⁺ system. In addition, the quantum mechanical and statistical approaches used for the computation of state-to-state integral cross sections are introduced. Section 3 is devoted to the analysis of the new ground state electronic PES, as well as the CC, SQM and SACM results that were obtained. In Section 4, we summarize our findings and provide an outlook regarding future studies involving the ³⁶ArH₂⁺ system.

2 Theoretical methods

2.1 Potential energy surface

2.1.1 Ab initio calculations. The ab initio electronic structure calculations were performed for ArH₂⁺ structures in the molecular configuration space defined by the internal coordinates (r_ArH′,r_HH′,φ) (see Fig. 1 for the definition of the internal coordinates). The ranges of the molecular coordinates were chosen to be such that r_ArH′ ∈ [1.00 Å, 4.00 Å], r_HH′ ∈ [0.50 Å, 50.00 Å] and φ ∈ [0°, 180°] in the exit arrangement, and r_ArH′ ∈ [1.00 Å, 50.00 Å], r_HH′ ∈ [0.50 Å, 4.00 Å] and φ ∈ [0°, 180°] in the entrance arrangment, in order to ensure that all the relevant parts of the configuration space for electronically adiabatic reaction dynamics of reaction (1) at energies below ∼ 0.3 eV (with respect to the bottom of the well corresponding to the Ar (¹S) + H₂⁺ (X²Σ_g⁺) asymptote) were probed. Furthermore, because of the necessity to include closed vibrational and rotational manifolds in order to ensure convergence of the scattering matrix (S-matrix) (see Section 2.2 for more information on the scattering calculations), care was taken to include geometries with energies up to the dissociation limit (i.e., up to ∼ 3 eV with respect to the bottom of the well corresponding to the Ar (¹S) + H₂⁺ (X²Σ_g⁺) asymptote). The equilibrium internuclear distances of the diatomics (in their respective electronic ground states) involved are 1.056 Å (for H₂⁺) and 1.281 Å (for ArH⁺).

Fig. 1 The system of internal coordinates which were used for the electronic structure calculations and GP regression fit relating to the ArH₂⁺ system.

The energies defining the ground state PES (transforming according to irreducible representation A′ in C_s symmetry), were calculated at the Davidson-corrected, internally contracted multireference configuration interaction (ic-MRCI+Q) level of theory, in conjunction with the aug-cc-pVQZ one-electron basis set for all the atoms. The multiconfigurational reference wavefunction for the configuration interaction calculations was obtained by means of the state-averaged complete active space self-consistent field (SA-CASSCF) approach, the molecular orbital space of which consisting of 5 inactive orbitals (4A′ + 1A′′) and 6 active orbitals (5A′ + 1A′′), and thus 9 active electrons. The molecular orbitals were obtained by minimizing the average of the electronic energy of the three lowest states (2A′ + 1A′′) with equal weights; this because of the involvement of the excited states in the diabatic fit for the small H₂⁺ internuclear distances.

At this stage it is important to mention that these calculations are in principle very similar to those performed in our previous work on the PESs of the ArH₂⁺ system in ref. 14. However, because of the nature of the quantum scattering calculations (as discussed in Section 2.2) we wanted to perform for reaction (1), it proved to be necessary to obtain a new PES for reasons discussed in more detail in Section 3.1. This new PES only describes the Ar + H₂⁺ entrance asymptote, excluding the Ar⁺ + H₂ asymptote, and is therefore free of spin–orbit effects.

All electronic structure calculations were carried out using the MOLPRO package 2020.2.²⁷

2.1.2 Gaussian process regression. In this subsection we discuss the idea of fitting PESs of molecular systems by means of Gaussian processes. The reader interested in a more detailed account can consult other sources.^28,29

In brief, there is a normal distribution of electronic energies, ρ(E) ∼ [scr N, script letter N](μ,σ), at each point in the molecular configuration space of the system, which is conditioned by means of training geometries, {x_i}, and associated energies, E = ({E_i}), such that the mean of this distribution at a geometry, x_*,

for which no single point calculations were performed, is a prediction of the electronic energy at said geometry. The kernels, k(x_i,x_j), are the matrix elements of the covariance matrix, K, in eqn (3), and they are assumed to take a certain functional form. For the purpose of fitting the ground state PES of the title system, we expressed the kernels as follows,²⁰

k(x_i,x_j) = D(x_i,x_j,Σ)[scr M, script letter M]_ν(x_i,x_j,l) + W(x_i,x_j,n), (4)

where [scr M, script letter M]_ν(x_i,x_j,l) stands for the Matérn class of kernels, D(x_i,x_j,Σ) is the dot product kernel, and W(x_i,x_j,n) is the white noise kernel. The explicit mathematical forms of these individual kernels are given by,

D(x_i,x_j,Σ) = Σ² + x_i·x_j (6)

W(x_i,x_j,n) = nδ_xi,xj. (7)

In eqn (5), [scr K, script letter K]_ν is a modified Bessel function of the second kind, and its argument is given by , where d(x^(k)_i,x^(k)_j) represents the euclidean distance between the training points along the internal coordinate k. Γ(ν) is the gamma function and the value of the parameter ν was fixed at 3/2, such that the obtained ML representations are twice differentiable. As is apparent from eqn (5)–(7) there are 6 parameters (in this context often referred to as hyperparameters), namely, n, Σ, ν, and (l₁,l₂,l₃). The values of these parameters (with the exception of ν) are determined in an optimisation procedure that involves the maximisation of the log marginal likelihood,as implemented in the scikit-learn python package.³⁰

The ML model was trained in the Morse representation of the internuclear distances (see Fig. 1), exp(−r_AB), thereby facilitating the accurate representation of both the short- and long-range potential.

2.2 Close-coupling method

Time-independent quantum mechanical close-coupling calculations for the ³⁶Ar + H₂⁺ collisional system on the aforementioned ground state adiabatic electronic PES (1²A′) were performed using the ABC quantum reactive scattering program by Skouteris et al.³¹ Such calculations involve approximately solving the time-independent Schrödinger equation for the motion of the atomic nuclei on a single adiabatic electronic PES by means of a multi arrangement expansion of the nuclear wavefunctions in terms of Delves hyperspherical coordinates.³² In such an approach both reactive and non-reactive (elastic and inelastic) processes can be simulated in full state-to-state detail of the reactants and the products. Performing these ABC calculations requires the input of several parameters, the values of which are determined by means of convergence tests. Several convergence tests at total angular momentum, J = 0, were performed in order to find the optimal values of said parameters: E_int,max = 2.2 eV (maximum internal energy, E_int, in any channel), j_max = 33 (maximum value of the rotational quantum number, j, in any channel), ρ_max = 50a₀, 35a₀ and 25a₀ (maximum values of the hyperradius, ρ) and N_sec = 200 (number of log derivative propagation sectors). Before calculating cross sections, the S-matrix was transformed into a helicity adapted representation based on the prescription in ref. 31. We found that at the maximum collision energy considered (E^max_c = 0.1 eV) the state-to-state ICSs were converged at a total angular momentum quantum number, J = 51. At said value of J, convergence tests suggested that k_max = 4 (maximum value for the helicity quantum number, k) presented to optimal balance between accuracy and computation time.

2.3 Statistical methods

For collisional systems that exhibit a (deep) potential energy well on their electronic PES, statistical theories of molecular scattering have been shown to provide reasonably accurate approximate results. Indeed, a collision proceeding through the formation and decay of an intermediate collision complex (i.e., complex-forming collisions, or sometimes, complex-mediated collisions), can behave statistically if the lifetime of the complex is long enough,^33–35 which is typically the case for deep potential energy wells (with a large density of rovibrational states) and low collision energies (such that the number of open channels is small). Several statistical methods have been developed over the years; examples include, but are not limited to: phase space theory (PST)³⁶ and its modified version, mean potential phase space theory (MPPST),^37–39 the above mentioned SQM,^15,16 its wave packet⁴⁰ and quasi-classical trajectory versions,^41–43 and variations of the statistical adiabatic channel model.^17,18,44,45

2.3.1 Statistical quantum mechanical method. Specially designed for reactions which take place by means of the formation of an intermediate complex, the SQM approach of ref. 15, 16 and 46 assumes that the state-to-state reaction probability at the total angular momentum J and energy E can be expressed according with the following approximation:where and p^J_vj(E) are the capture probabilities, respectively, for the formation of the collision complex from the initial H₂⁺(v,j) rovibrational state and for its fragmentation into the final ³⁶ArH⁺(v′,j′) state, and the indices in the sum of the denominator, v′′j′′, run over all the energetically open rovibrational states both in reactants and products. As explained in ref. 16, those individual probabilities are calculated by means of a time-independent method which involves a log derivative propagation between a capture radius R_c and R_max, the asymptotic distance, which in this case have been set to ∼3.7 Å and ∼44 Å, respectively. The final value of the capture radius is chosen after some tests in which no noticeable differences are found when this is varied. The reaction probability from eqn (9) is then employed for the calculation of ICSs:where k_vj² = 2μ(E − E_vj)/ħ², with E_vj as the energy of the initial rovibrational state of the reactant diatom H₂⁺ and μ the reduced mass of the system.

The SQM calculations have been performed within the CC framework with no further approximations regarding the helicity Ω coupling. However, given the large number of rovibrational states in the product channel only those (v′,j′) channels with energies beyond a certain cutoff value, E_min, have been included for a proper calculation of the corresponding probabilities. For those other states with energies (including the corresponding centrifugal barrier contributions for each partial wave J) below such a E_min, a constant unity value, has been assumed. Test calculations have been performed in order to ensure the choice of E_min does not affect the state-to-state ICSs. Recent examples of applications of the SQM approach to reactive collisions involving protonated rare gas systems can be found in the literature.^47,48

2.3.2 Statistical adiabatic channel model. Another one of such statistical approaches is the SACM originally introduced by Quack and Troe in the mid-1970s.^44,45 More recently a modified version of their original model was applied to a number of complex-forming inelastic^17,49,50 and competing inelastic and reactive collisions.^18,51 From these studies it was found that this new SACM-inspired approach is capable of reproducing high-level-of-theory state-resolved rate coefficients within a factor of 2, with the error decreasing when the temperature decreases and/or the depth of the energy well increases; and all this at a fraction of the computational cost. This approach therefore provides a great alternative to full quantum calculations, since the latter tend to be prohibitive, if not impossible, for systems with a deep potential energy well, due to the excessive number of closed channels that need to be taken into account.

In the SACM approach, the state resolved integral cross sections are computed according to eqn (10). However, the state-to-state transition probabilities in eqn (9) 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 energy, E, and total angular momentum, J. To provide a count of these 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. More detailed information about this modified SACM approach can be found in ref. 17 and 18.

2.4 Rate coefficients

Based on the state-to-state ICSs, and assuming a Maxwell–Boltzmann distribution of collision energies, temperature-dependent rate coefficients, k_if(T), are computed as follows,where k_B is the Boltzmann constant and E_c is the collision energy, the latter being equal to E − E_vj.

3 Results and discussion

3.1 Potential energy surface

To construct the ground state adiabatic electronic PES (1²A′) of the ArH₂⁺ system, we have performed single point calculations (utilising the electronic structure methods discussed in Section 2.1.1) at a total of 2004 symmetry unique training geometries. The ML fit is characterised by a root-mean-square-error (RMSE) of ∼ 5.97 × 10⁻⁴ eV, based on comparison with the test set data (the RMSE on the training data is ∼ 9.22 × 10⁻¹⁰ eV), which suggests that we have obtained a highly accurate representation of the PES appropriate for chemical dynamics simulations. The test data is comprised of 224 ab initio energies that were not part of the training set for machine learning, and the corresponding molecular structures were generated randomly within the ranges mentioned in Section 2.1.1. This potential is similar to the one reported in ref. 14, the major differences being: (i) the ML fit of this PES ignores the avoided crossing that exists in the entrance arrangement between the two lowest electronic states and pretends that the Ar⁺ + H₂ channel doesn't exist by evolving diabatically for ( being the minimum in the PES along the H₂⁺ internuclear distance for given r_ArH and φ), and consequently (ii) this PES does not include spin–orbit coupling effects since they are of importance only in the Ar⁺ + H₂ channel. As already touched upon in the introduction (Section 1), this approach to the simulation of the reaction dynamics is valid as long as the total scattering energy is lower than the top of the barrier created by the aforementioned avoided crossing because (i) the excited states are not accessible energetically, and (ii) the zero-point vibrational energy of H₂ does not allow for its formation under said energy constraint. Furthermore, such a modified fit is in fact a necessity since the CC calculations (as discussed in Section 2.2) require the construction of a multi-arrangement rovibrational basis, and the vibrational eigenstates in the entrance arrangement can only be computed properly in this way. So, with the exception of the crossing seam created by the interaction between the ground and lowest excited states, the general features of this ground state PES and the one recently published¹⁴ are the same; a more detailed account of the latter can be found in ref. 14. We only show here two-dimensional contour plots (Fig. 2) that provide some information regarding the topography of the PES for fixed bond angles, primarily because it is of interest when discussing the results from the quantum scattering calculations and the statistical calculations. From these plots in Fig. 2 the presence of a potential energy well (∼ −1.8 eV with respect to the Ar + H₂⁺ asymptote for the H–H–Ar linear structure corresponding to the global minimum), as well as a barrier for small bond angles becomes apparent.

Fig. 2 Two-dimensional contour plots (for fixed bond angles, φ) of the Gaussian process regression representation of the ground state ArH₂⁺ adiabatic PES (energies are in eV and are reported with respect to the ArH⁺ (X¹Σ⁺) + H (²S) asymptote). The energy well, as well as the exothermic character of the H⁺-transfer are clearly visible.

3.2 Scattering and statistical calculations

The present subsection is devoted to an analysis of the state-resolved integral scattering cross sections as a function of collision energy, obtained from quantum mechanical close-coupling calculations (described in Section 2.2) and statistical calculations (described in Section 2.3.1 and Section 2.3.2) all of which were performed on the new ground state adiabatic PES (described in Section 3.1). Quantum and statistical cross sections for inelastic collisions,and for competing reactive scattering events,

³⁶Ar + H₂⁺ (v_H2⁺,j_H2⁺) → ³⁶ArH⁺ (v_ArH⁺,j_ArH⁺) + H, (14)

will be discussed. Due to the total scattering energy being constrained (as discussed in Section 3.1) as well as the values for the E_int,max and j_max parameters for the ABC calculations (Section 2.2), so are the values of the vibrational and rotational quantum numbers in eqn (13) and (14): , and v_ArH⁺ = 0–4,j_ArH⁺ = 0–33.

3.2.1 Inelastic collisions. In Fig. 3 we show the ICSs as of function of collision energy for the inelastic excitations in eqn (13), calculated by means of the CC method and the SQM and SACM approaches. Trivially, cross sections become non-zero at the rotational thresholds; in addition, their maximum magnitudes within the collision energy range studied are similar. Comparing the ICSs in Fig. 3 with reactive ICSs (discussed in the next subsection Section 3.2.2), we notice that inelastic ICSs are generally somewhat larger. This could be related to the existence of a barrier for small bond angles (as apparent in the contour plots in Fig. 2). This barrier hinders the formation of the product ArH⁺, especially for bond angles smaller than 80° for which said barrier is no longer submerged, thereby favouring non-reactive scattering. Regarding the comparison with the SACM and SQM approaches, it is clear from Fig. 3 that both SACM and SQM give very similar results which underestimate CC cross sections by about a factor of 10. This disagreement was not expected and might have several origins. For SACM, presumably, the error is to be found in the use of adiabatic channel potentials to count the number of open channels at a given energy and total angular momentum; specifically, these radial curves are computed by taking into account the angular dependence of the PES in an average way, while fixing the internuclear distance of the diatomic at its equilibrium value. Since said barriers, and primarily the largest values thereof, occur at HH internuclear distances that are larger than the HH internuclear equilibrium distance, their effect is only partly, if at all, accounted for. We therefore assume that the effect of the aforementioned barrier on the statistical cross sections is not properly captured. That would imply that the total number of open channels in the denominator in eqn (11) is overestimated, which then leads to cross sections that are too small.

Fig. 3 Close-coupling, SACM and SQM integral cross sections (in Ų) as a function of collision energy (in eV) for the inelastic excitations involving ³⁶Ar and H₂⁺.

It is difficult to apply the same reasoning in the case of SQM. Given its intrinsically approximate nature, predictions obtained with this approach for the process shown in eqn (13) where rovibrational levels of the same species, H₂⁺, are involved as reactants and products should be understood as result of a collision mediated by the formation and fragmentation of an intermediate complex. Whether this is exactly the same kind of processes described by means of the ABC approach may depend on the precise reaction under study. Therefore, in principle, differences with a purely statistical behaviour can be due to deviations from a purely complex-forming pathway.

3.2.2 Reactive collisions. Regarding the reactive scattering processes in eqn (14), we present in Fig. 4 the ICSs for transitions starting from v_H2⁺ = 0 and j_H2⁺ = 0. We have restricted ourselves to values of j_ArH⁺ = 0–5 in order to keep the figures clear and readable. The CC data show following trends:

Fig. 4 Close-coupling, SACM and SQM cross sections (in Ų) as a function of collision energy (in eV) for the reactive collisions ³⁶Ar + H₂⁺ (v_H2⁺ = 0,j_H2⁺ = 0) → ³⁶ArH⁺ (v_ArH⁺ = 0–4,j_ArH⁺ = 0–5) + H. Last panel: CC ICSs for v_ArH⁺ = 0–4 and j_ArH⁺ = 0.

(1) Within a given vibrational manifold v_ArH⁺, and for the collision energies considered, the ICSs generally increase with increasing rotational excitation of ³⁶ArH⁺.

(2) Cross sections are larger for higher vibrational excitation of ³⁶ArH⁺, suggesting that ³⁶ArH⁺, when formed through reaction (1), is a vibrationally hot product. This has indeed been observed in previous studies on the matter.^1,3 In order to more clearly show the effect of vibrational excitation of ArH⁺, the last panel in Fig. 4 shows the CC data for ArH⁺ (v_ArH⁺,j_ArH⁺ = 0).

Looking more closely at the comparison with the statistical predictions in the same figure (Fig. 4), we notice that the agreement between CC and SACM/SQM is quite good for v_ArH⁺ = 0. However, the agreement becomes progressively worse for higher vibrational excitation of ³⁶ArH⁺. This observation seems to suggest that (v_H2⁺ = 0,j_H2⁺ = 0) → (v_ArH⁺ > 0,j_ArH⁺) transitions in the ³⁶Ar + H₂⁺ reactive system behave less statistically for increasing v_ArH⁺. Presumably, this is related to the large exothermicity of the reaction (∼ 1.32 eV),¹⁴ which would imply a large number of open rovibrational exit channels, and thus a lesser extent of randomisation of the available energy over the rovibrational states of the collision complex. Previous investigations where SQM predictions were compared with exact quantum mechanical results revealed that predictions at the state-to-state level were a much more stringent test for statistical techniques than those processes in which, for instance, only the initial state was selected (see, for example ref. 52 and 53.)

From the results already introduced, we can compute initial-state selected cross sections for the reaction in eqn (14) by summing over all accessible rovibrational states of ³⁶ArH⁺ that are open for a given total scattering energy. In Fig. 5 we compare CC, SQM and SACM cross sections for reactive collisions starting in the rovibrational ground state of H₂⁺. Overall, the results shown in Fig. 5, both quantum mechanical and statistical, agree well. Both statistical approaches (SQM and SACM) make very similar predictions, differing by at most a factor of ∼ 1.2, and they both underestimate the quantum cross section to within a factor of 2. We should also note that the stepwise behaviour of the SACM result smooths out with the inclusion of quantum tunnelling in the entrance channels (see Fig. 5, top panel); tunnelling was accounted for within the Wentzel–Kramers–Brillouin (WKB) approximation,^54,55 implying that the capture probabilities, , which are either 0 or 1 classically (as assumed by eqn (11)), are now replaced by a semiclassical transmission coefficient, the value of which can take any value between 0 and 1.

Fig. 5 Initial-state selected integral cross sections (in Ų) as a function of collision energy (in eV) (top panel), as well as (initial-state selected) rate coefficients (in cm³ s⁻¹) as a function of temperature (in K) (bottom panel), both computed with the close-coupling and statistical methods (SQM and/or SACM) for reactive collisions between ³⁶Ar and H₂⁺ (v_H2⁺ = 0,j_H2⁺ = 0).

As the initial-state selected cross sections for the reaction ³⁶Ar + H₂⁺ (v_H2⁺ = 0,j_H2⁺ = 0) agree very well, so do the initial-state selected rate coefficients (Fig. 5). The CC value varies from ∼ 6.13 × 10⁻¹⁰ cm³ s⁻¹ at 1 K, to ∼ 2.68 × 10⁻⁹ cm³ s⁻¹ at 150 K, which can be compared to the previous estimate ∼ 10⁻⁹ cm³ s⁻¹ of ref. 5 used to estimate the abundance of ArH⁺ in astrophysical environments. In the same Fig. 5 also the close-coupling thermal rate coefficient is shown up to a temperature of 50 K, as well as the Langevin prediction (∼ 2.18 × 10⁻⁹ cm³ s⁻¹) thereof.

4 Conclusions

In this study, we have computed the first fully state-resolved quantum mechanical integral scattering cross sections as a function of collision energy for collisions between ³⁶Ar and H₂⁺ (v = 0,j = 0–4) on a new full-dimensional adiabatic ground state potential energy surface (²A′). We have found that the initial-state selected rate coefficients for the reaction ³⁶Ar + H₂⁺ (v = 0,j = 0) → ³⁶ArH⁺ + H, for temperatures up to 150 K, have a magnitude of ∼ 10⁻⁹ cm³ s⁻¹, roughly similar to what was assumed before. Furthermore, we have investigated the statistical nature of the ³⁶ArH₂⁺ collisional system by means of comparison to statistical approaches, such as the statistical quantum method and the statistical adiabatic channel model. Based on the agreement we find between quantum mechanical and statistical initial-state selected cross sections for H₂⁺ in its ground rovibrational state, one could conclude that the collision dynamics of the ³⁶Ar + H₂⁺ reaction has some statistical character, at least in the collision energy range that was studied, and without quantum-state resolution in the product channels. However, a more detailed investigation of said reaction, based on the state-to-state cross sections, suggests that statistical behaviour of the reaction dynamics decreases with vibrational excitation of the ³⁶ArH⁺ product. We also find that SQM and SACM predictions are almost identical. A further comparison has shown that the statistical adiabatic channel model fails to accurately predict ICSs for inelastic collisions involving ³⁶Ar and H₂⁺, presumably because of its inability to account for barriers that hinder the formation of the product ³⁶ArH⁺, as shown by the analysis of the topography of the ground state PES. The next step would be to extend the collision energy range considered in this work. However, as discussed, this would require an electronically non-adiabatic treatment, whereby one should account for the competition with ³⁶Ar⁺ + H₂ channels through charge exchange.

Data availability

Some of the data supporting this article have been included as part of the ESI.†

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

J. L. acknowledges support from KU Leuven through grant no. C14/22/082. The resources and services used in the computations were provided by the VSC (Flemish Supercomputer Center), funded by the Research Foundation-Flanders (FWO) and the Flemish Government. T. G. L. thanks MCIN/AEI/10.13039/501100011033 PID2020-114654GB-I00 for its support.

Notes and references

1. F. Aguillon and M. Sizun, J. Chem. Phys., 1997, 106, 9551–9562.
2. S. Chapman and R. K. Preston, J. Chem. Phys., 1974, 60, 650–659.
3. S. Chapman, J. Chem. Phys., 1985, 82, 4033–4043.
4. M. Sizun, J.-B. Song and E. A. Gislason, J. Chem. Phys., 2002, 116, 2888–2895.
5. P. Schilke, D. A. Neufeld, H. S. P. Müller, C. Comito, E. A. Bergin, D. C. Lis, M. Gerin, J. H. Black, M. Wolfire, N. Indriolo, J. C. Pearson, K. M. Menten, B. Winkel, A. Sánchez-Monge, T. Möller, B. Godard and E. Falgarone, Astron. Astrophys., 2014, 566, A29.
6. P. Dagdigian, Mon. Not. R. Astron. Soc., 2018, 477, 802–807.
7. M. J. Barlow, B. M. Swinyard, P. J. Owen, J. Cernicharo, H. L. Gomez, R. J. Ivison, O. Krause, T. L. Lim, M. Matsuura, S. Miller, G. Olofsson and E. T. Polehampton, Science, 2013, 342, 1343–1345.
8. H. S. P. Müller, S. Muller, P. Schilke, E. A. Bergin, J. H. Black, M. Gerin, D. C. Lis, D. A. Neufeld and S. Suri, Astron. Astrophys., 2015, 582, L4.
9. R. M. García-Vázquez, M. Márquez-Mijares, J. Rubayo-Soneira and O. Denis-Alpizar, Astron. Astrophys., 2019, 631, A86.
10. M. Baer and H. Nakamura, J. Chem. Phys., 1987, 87, 4651–4664.
11. J.-B. Song and E. A. Gislason, Chem. Phys., 2003, 293, 231–237.
12. P. J. Kuntz and A. C. Roach, J. Chem. Soc., Faraday Trans. 2, 1972, 68, 259–280.
13. M. Hu, W. Xu, X. Liu, R. Tan and H. Li, J. Chem. Phys., 2013, 138, 174305.
14. M. Konings, J. N. Harvey and J. Loreau, J. Phys. Chem. A, 2023, 127, 8083–8094.
15. E. J. Rackham, F. Huarte-Larranaga and D. E. Manolopoulos, Chem. Phys. Lett., 2001, 343, 356–364.
16. E. J. Rackham, T. Gonzalez-Lezana and D. E. Manolopoulos, J. Chem. Phys., 2003, 119, 12895–12907.
17. J. Loreau, F. Lique and A. Faure, Astrophys. J., 2018, 853, L5.
18. M. Konings, B. Desrousseaux, F. Lique and J. Loreau, J. Chem. Phys., 2021, 155, 104302.
19. K. M. Hickson, P. Larrégaray, L. Bonnet and T. González-Lezana, Int. Rev. Phys. Chem., 2021, 40, 457–493.
20. A. Christianen, T. Karman, R. A. Vargas-Hernández, G. C. Groenenboom and R. V. Krems, J. Chem. Phys., 2019, 150, 064106.
21. O. T. Unke, D. Koner, S. Patra, S. Käser and M. Meuwly, Mach. Learn.: Sci. Technol., 2020, 1, 013001.
22. H. Sugisawa, T. Ida and R. V. Krems, J. Chem. Phys., 2020, 153, 114101.
23. C. Qu, Q. Yu, B. L. Van Hoozen, J. M. Bowman and R. A. Vargas-Hernández, J. Chem. Theory Comput., 2018, 14, 3381–3396.
24. J. Dai and R. V. Krems, J. Chem. Theory Comput., 2020, 16, 1386–1395.
25. Q. Song, Q. Zhang and Q. Meng, J. Chem. Phys., 2020, 152, 134309.
26. A. Kamath, R. A. Vargas-Hernández, R. V. Krems, T. Carrington and S. Manzhos, J. Chem. Phys., 2018, 148, 241702.
27. H.-J. Werner, P. J. Knowles, F. R. Manby, J. A. Black, K. Doll, A. Heßelmann, D. Kats, A. Köhn, T. Korona, D. A. Kreplin, Q. Ma, T. F. Miller, A. Mitrushchenkov, K. A. Peterson, I. Polyak, G. Rauhut and M. Sibaev, J. Chem. Phys., 2020, 152, 144107.
28. C. E. Rasmussen and C. K. I. Williams, Gaussian Processes for Machine Learning, The MIT Press, Cambridge, 2006.
29. R. V. Krems, Phys. Chem. Chem. Phys., 2019, 21, 13392–13410.
30. F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot and E. Duchesnay, J. Mach. Learn. Res., 2011, 12, 2825–2830.
31. D. Skouteris, J. Castillo and D. Manolopoulos, Comput. Phys. Commun., 2000, 133, 128–135.
32. G. A. Parker and R. T. Pack, J. Chem. Phys., 1993, 98, 6883–6896.
33. R. B. Bernstein, A. Dalgarno, H. S. W. Massey and I. C. Percival, Proc. R. Soc. London, Ser. A, 1963, 274, 427–442.
34. J. C. Light, Discuss. Faraday Soc., 1967, 44, 14–29.
35. W. H. Miller, J. Chem. Phys., 1970, 52, 543–551.
36. P. Pechukas, J. C. Light and C. Rankin, J. Chem. Phys., 1966, 44, 794–805.
37. P. Larrégaray, L. Bonnet and J.-C. Rayez, J. Phys. Chem. A, 2006, 110, 1552–1560.
38. P. Larrégaray, L. Bonnet and J.-C. Rayez, J. Chem. Phys., 2007, 127, 084308.
39. L. Bonnet, P. Larrégaray and J.-C. Rayez, Phys. Chem. Chem. Phys., 2007, 9, 3228–3240.
40. S. Y. Lin and H. Guo, J. Chem. Phys., 2006, 124, 031101.
41. F. J. Aoiz, V. Sáez Rábanos, T. González-Lezana and D. E. Manolopoulos, J. Chem. Phys., 2007, 126, 161101.
42. F. J. Aoiz, T. González-Lezana and V. Sáez Rábanos, J. Chem. Phys., 2008, 129, 094305.
43. P. G. Jambrina, F. J. Aoiz, C. J. Eyles, V. J. Herrero and V. Sáez Rábanos, J. Chem. Phys., 2009, 130, 184303.
44. M. Quack and J. Troe, Ber. Bunsenges. Phys. Chem., 1974, 78, 240–252.
45. M. Quack and J. Troe, Ber. Bunsenges. Phys. Chem., 1975, 79, 170–183.
46. T. González-Lezana, Int. Rev. Phys. Chem., 2007, 26, 29.
47. D. Koner, L. Barrios, T. González-Lezana and A. N. Panda, Molecules, 2021, 26(14), 4206.
48. M. J. Montes de Oca-Estévez, B. Darna, B. García-Ruiz, R. Prosmiti, T. González-Lezana and D. Koner, ChemPhysChem, 2023, 24, e202300450.
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. B. Desrousseaux, M. Konings, J. Loreau and F. Lique, Phys. Chem. Chem. Phys., 2021, 23, 19202–19208.
52. T. González-Lezana, O. Roncero, P. Honvault, J.-M. Launay, N. Bulut, F. Javier Aoiz and L. Bañares, J. Chem. Phys., 2006, 125, 094314.
53. E. Carmona-Novillo, T. González-Lezana, O. Roncero, P. Honvault, J.-M. Launay, N. Bulut, F. Javier Aoiz, L. Bañares, A. Trottier and E. Wrede, J. Chem. Phys., 2008, 128, 014304.
54. H. A. Kramers, Z. Phys., 1926, 39, 828–840.
55. G. Wentzel, Z. Phys., 1926, 38, 518–529.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d4cp02179g

---