Globally Accurate Full-Dimensional Potential ... - ACS Publications

Article Cite This: J. Phys. Chem. A XXXX, XXX, XXX−XXX

pubs.acs.org/JPCA

Globally Accurate Full-Dimensional Potential Energy Surface for H2 + HCl Inelastic Scattering Published as part of The Journal of Physical Chemistry virtual special issue “F. Javier Aoiz Festschrift”. Qian Yao,† Masato Morita,‡ Changjian Xie,†,# Naduvalath Balakrishnan,*,‡ and Hua Guo*,† †

Department of Chemistry and Chemical Biology, University of New Mexico, Albuquerque, New Mexico 87131, United States Department of Chemistry and Biochemistry, University of Nevada, Las Vegas, Nevada 89154, United States

Downloaded via BUFFALO STATE on July 25, 2019 at 22:12:57 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

‡

S Supporting Information *

ABSTRACT: A globally accurate full-dimensional potential energy surface (PES) for the inelastic scattering between H2 and HCl is developed on the basis of a large number of points calculated at the coupled-cluster singles, doubles, and perturbative triples level of theory. The machine-learned PES is trained with 42 417 ab initio points using the permutation invariant polynomial-neural network method, resulting in a root-mean-square fitting error of 5.6 cm−1. Both full- and reduced-dimensional quantum calculations for rotationally inelastic scattering are performed on this new PES and good agreement is obtained with previous quantum dynamical results on a reduced-dimensional model. Furthermore, strong resonances are identified at collision energies below 100 K, including cold conditions. This new PES provides a reliable platform for future studies of scattering dynamics with vibrationally excited collision partners in a wide range of collision energies extending to cold and ultracold conditions.

1. INTRODUCTION It has long established that all forms of energy are not equal in collisional processes.1−3 For reactive collisions, Polanyi observed that while translational energy is more effective in promoting a chemical reaction with an early barrier, vibrational energy enhances reactivity for a reaction with a late barrier.4 Such mode specificity can also be viewed in the Sudden Vector Projection (SVP) model as the result of differing coupling strengths of various reactant modes with the reaction coordinate at the transition state,5 including polyatomic reactant molecules.6 Experimental studies of mode specificity require the preparation of vibrational states. Until recently, this can only be achieved with a partial population of desired vibrational levels, which complicates the analysis.1,2 Nowadays, laser based coherent techniques, such as Stimulated Raman Adiabatic Passage (STIRAP or SARP)7,8 and Rapid Adiabatic Passage (RAP),9 have enabled near 100% population of a specific rovibrational state in diatomic as well as polyatomic molecules. These new capabilities open new doors for studying the effects of internal energy on both reactive and nonreactive collisional processes, which allows testing of new theory and models. This is particularly relevant in the cold and ultracold regime,10−12 in which the collision energy is vastly smaller than the energy stored in the vibrational modes of colliding partners, because such processes explore totally different regions of the potential energy surface (PES). In the recent experiments on inelastic collisions of SARP prepared vibrationally and rotationally excited HD(v = 1, j = 2) with H2 in a coexpansion beam, for © XXXX American Chemical Society

example, clear stereodynamical control was achieved at very low (∼1 K) collision temperatures.13,14 Theoretically, a proper understanding of mode specificity requires the construction of the Born−Oppenheimer PES from first-principles and in full dimensionality. This has been a great challenge15 but is recently being overcome, thanks to advances in both electronic structure theory and means to represent the PES using analytical functions.16−21 With the PES, it is possible to calculate the dynamical attributes related to mode specificity, using either quantum mechanical or quasi-classical trajectory methods for nuclear motions.22,23 In our recent investigation of Zare’s SARP work on HD(v = 1, j = 2) + H2 collisions,13,14 full-dimensional quantum scattering calculations were performed on an accurate H2−H2 PES of Hinde.24 These highly accurate quantum dynamical calculations not only reproduced the experimentally observed stereodynamics but also revealed the existence of dynamical resonances dominated by single partial waves25 and suggested a way to stereodynamically control dynamical resonances.26 In the current work, we report the first full-dimensional PES for the H2−HCl system. This system is of great interest in astrochemistry as H2 is the most abundant molecule in the interstellar medium and HCl and other chlorine-containing molecules are a tracer for chlorine chemistry in these environments.27−29 There has recently been much interest in Received: June 21, 2019 Revised: July 2, 2019 Published: July 3, 2019 A

DOI: 10.1021/acs.jpca.9b05958 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

small grids with appropriate coordinates in various regions were used to sample relevant configurations near the stationary points. Then, the selected points are calculated at the CCSD(T)-F12b/AVQZ level and a primitive PES was obtained by fitting the calculated ab initio points using the method detailed in the next subsection. Trajectories were launched with various initial conditions on the primitive PES to explore the configuration space and to generate new points. A new point (r) was added if it is below an energy cutoff (2 eV measured from the global minimum) and satisfies the χ(ri) =

this system, particularly related to rotational and hyperfine excitations in HCl of interest in astrochemistry.30−33 This system features a much stronger interaction than that between two H2 molecules. Previous studies have identified two Tshaped minima for the HCl−H2 complex.30,34,35 The global minimum corresponding to a geometry with the H atom of the HCl molecule pointing to the H2 molecule is about 200 cm−1 deep, while a secondary minimum has the Cl atom pointing to H2. The PES is given in an analytical form based on a feedforward neural network (NN), trained by 42 417 points calculated with an explicitly correlated coupled cluster method with single, double, and perturbative triple excitations (CCSD(T)-F12b). The permutation of the three hydrogen atoms in this system is explicitly enforced using the permutation invariant polynomial-neural network (PIP-NN) method.21 Furthermore, a careful treatment of the long-range interactions permits its usefulness in cold and ultracold regimes. This PES thus allows a detailed investigation of the importance of vibrations in rotationally inelastic scattering.36−38 Using this new PES, full-dimensional scattering calculations have been carried out with an exact timeindependent quantum method with both molecules in their ground vibrational states, but with HCl in its first few excited rotational states. Reduced-dimensional quantum scattering calculations have also been performed within the rigid-rotor approximation and using a vibrationally averaged fourdimensional (4D) PES. Comparisons with results from full six-dimensional (6D) calculations suggest that the 4D model with the vibrationally averaged potential is quite accurate for rotationally inelastic scattering involving low-lying states of the collision partners. In addition, the good agreement obtained with the previous reduced-dimensional quantum scattering results using a 4D PES30 further validates the accuracy of the new full-dimensional PES. This new PES offers a reliable platform for future studies of inelastic scattering with internally excited H2 and/or HCl molecules, which, along with other similar work,39−43 is expected to give accurate information on the vibrational effects on energy transfer, particularly in cold and ultracold regimes.

6 ∑i = 1 |ri − ri′|2 > 0.1 Å condition for all points (r′i ) already in the training set. This process is repeated until convergence. All points were chosen with the distance between the two diatomic centers of mass (R) no larger than 30 Å. 2-B. Fitting. The PES is divided into two segments to cover the short-range and long-range regions. The short-range part of the PES, denoted as VPIP‑NN, was fit to ab initio points using the PIP-NN method.21 A total of 42 417 points was selected. In the PIP-NN approach, the input layer of the NN is not given in terms of the internuclear distances (rij), which are invariant under rotation, translation, and inversion but do not satisfy the permutation symmetry. The invariance of the PES with respect to permutation of identical atoms is enforced by low-order PIPs,47,48 which can be defined as symmetrized monomials:16 N

G = S ̂ ∏ pijlij

(1)

i 5 Å) between two molecular centers of mass and is fitted by a linear least-squares method using the following expression,49 Vint =

∑

Cl1l2 , m , nυl1l2(r1 ,r2)

l1l 2 , m , n

bnm(θ1 ,θ2 ,φ) Rn

(4)

A more detailed description of this formula and the choice of the parameters are given in Supporting Information. A total of 24 870 ab initio points was used in the fitting and these points were selected from the 42 417 ab initio points used in the PIPNN fitting with R larger than 5 Å. The overall PES is switched from the short-range PES to its long-range counterpart as R increases: VPES = sVPIP − NN + (1 − s)VLR

Figure 2. Three-dimensional plot of the PES as a function of R and θ2 with fixed r1 = 0.74 Å, r2 = 1.28 Å, θ1 = 90°, and φ = 0°. The two Tshaped minima are shown as potential wells.

displayed as a function of R and θ2 for fixed values of other four coordinates, r1 = 0.74 Å, r2 = 1.28 Å, θ1 = 90°, and φ = 0 °. The two T-shaped wells are clearly seen in the figure. The fitting error of the long-range potential is 0.7 cm−1. It should be noted that the functions used in describing the longrange PES give the correct asymptotic behaviors, so the fitting errors will not lead to qualitatively incorrect results at low collision energies. In Figure 3, ab initio data are compared

(5)

where the switching function is defined as s=

1 − tanh[3(R − R switch)] 2

(6)

Table 1. Energies (in cm−1), Harmonic Frequencies (in cm−1), and Geometric Parameters (in Å and deg) of the Stationary Points for the H2 + HCl System frequency species H2 + HCl H2···HCl H2···ClH saddle point

level

energy

ab initio PES ab initio PES ab initio PES ab initio

0.0 0.0 −215.5 −216.9 −102.6 −96.9 −56.9

97.3 99.3 61.7 72.5 105.3i (e)a

122.0 118.6 90.7 89.3 73.3

−65.6

80.7i (e)

75.3

PES

1

2

3

141.6 137.5 102.2 102.1

4

334.6 320.8 150.0 125.2

5

6

r1

r2

R

θ1

θ2

φ

2992.0 2991.8 2985.9 2988.7 2990.8 2991.6 2991.1

4404.5 4402.4 4387.3 4388.7 4400.2 4397.2 4403.5

0.7416 0.7416 0.7428 0.7428 0.7419 0.7419 0.7417

1.2767 1.2767 1.2775 1.2773 1.2769 1.2769 1.2769

3.5538 3.5682 3.3161 3.3066 3.6070

90.0 90.0 90.0 90.0 180.0

0.0 0.0 180.0 180.0 180.0

0.0 0.0 0.0 0.0 0.0

2988.4

4398.1

0.7416

1.2769

3.6251

180.0

180.0

0.0

The “e” denotes the double degeneracy of the mode.

a

C


Article


Figure 3. Comparison of ab initio energies and the PIP-NN PES, long-range PES, and overall PES at different orientations.

Figure 4. Two-dimensional contour plots of the PES along the (R, θ1) and (R, θ2) coordinates. The other coordinates are fixed, as indicated by the figure.

along R with the PIP-NN PES, long-range PES and overall PES at different orientations. It can be readily seen that the longrange PES is in excellent agreement with the ab initio energies. In Figure 4, two-dimensional contours of the PES are shown, which indicate that the switching between the short- and longrange potentials is smooth.

reported so far in full-dimensions. Lanza et al. reported pure rotational transitions employing a 4D PES that does not include stretching of the HCl and H2 bonds.30 They identified propensities for near-resonant rotational transitions when both molecules are initially rotationally excited, as in H2(j1 = 2) + HCl(j2 = 7) → H2(j1 = 0) + HCl(j2 = 9). This transition with an energy gap of ΔE = 3.6 cm−1 is found to be the most dominant at low collision energies compared to other rotational transitions out of these initial states. Many other near-resonant transitions with even smaller energy gaps are possible once the vibrational structure of both molecules is accounted for. Lanza and Lique also reported hyperfine-

3. ROTATIONAL INELASTIC COLLISIONS The H2 + HCl system is a prototype for inelastic collisions involving rotational and vibrational transitions and is amenable to full quantum close-coupling (CC) calculations, at least for low collision energies. However, exact CC calculations of rovibrational transitions in the H2 + HCl system have not been D


Article


which the basis for the angular variables correspond to the eigenstates of the total angular momentum, |JM(lj12(j1j2))⟩, where j12 is the quantum number for the total rotational angular momentum of the two molecules (j12 = j1 + j2), l is the quantum number of the orbital angular momentum (l) for the relative motion between molecules 1 and 2, and J is the quantum number denoting the total angular momentum of the collision complex, J = l + j12. Note that the basis is an eigenfunction of the spatial inversion operator, with the eigenvalue of εI = (−1)j1+j2+l. The total wave function for each J, M, and εI is expanded as

structure resolved cross sections for H2 + HCl collisions using a recoupling scheme within the rigid-rotor approximation.31 Here, we report preliminary calculations of rotational transitions in the H2 + HCl system within the v = 0 vibrational manifold using the full-CC approach and compare our results with those of Lanza et al.30 We also examine the accuracy of different 4D potentials for rotational transitions. These 4D potentials are obtained by fixing the H2 and HCl molecules to their equilibrium geometries, by fixing to vibrationally averaged bond distances, and by averaging the potential with respect to the v = 0 vibrational wave functions of the two molecules. Below, we provide a brief outline of the scattering formalism. The time-independent Schrödinger equation for the collision dynamics is solved with the CC formalism,39,50−53 using the full-dimensional PES computed in this study. 3-A. Quantum Scattering Method. The Hamiltonian for the relative motion of two diatomic molecules (after separating out the center-of-mass motion) can conveniently be expressed using a set of Jacobi vectors (r1, r2, R) in space-fixed coordinate frame as38−41 H(r1,r2,R) = T1(r1) + T2(r2) + T (R) + V (r1,r2,R)

Ψ JMεI =

(7)

+

lim Uint(r1,r2,R) → 0

1

2

1 1

2 2

(11)

vjl

∑ Uint vjlJ ,,εv ′ j ′ l′(R) FvJ′εj ′ l′(R) = 0 I

I

(12)

where μ is the reduced mass of the molecule−molecule system and EC is the relative collision energy. The latter is given by EC = E − εv1j1 − εv2j2, where E is the total energy and εviji denote the ro-vibrational energies of the isolated molecules. Note that the CC equations are independent of M, the projection of J on the space-fixed axis. For a given value of E, J, and εI, the CC equations are solved using the log-derivative matrix propagation method54 and asymptotic boundary conditions are applied at R = Rmax to yield the scattering matrix SJ,εI(E), from which the state-to-state integral cross section is evaluated:

(8)

σα → α′(E) =

(9)

π (2j1 + 1)(2j2 + 1)k 2

∑

×

and the total PES V becomes the sum of the diatomic potentials in this limit. Thus, asymptotically, the total molecular wave function is given by a product of rotational and vibrational wave functions of the two molecules. For paraH2 (nuclear spin I = 0) only even rotational levels are present whereas for ortho-H2 (nuclear spin I = 1) only odd rotational levels exists and ortho-para transition is not allowed in nonreactive collisions of H2. In the CC approach, these asymptotic wave functions and the orbital angular momentum eigenfunctions for the relative motion between the molecules form the basis sets, |ϕvj1χvj2 ⟩|JM(lj12 (j1 j2 ))⟩

I

v ′ j ′ l′

where Uint is the interaction potential between the two molecules and Vi (i = 1, 2) are the isolated diatomic potentials that depend on the bond distances ri = |ri|. In the current system, the diatomic potentials correspond to VH2(r1) and VHCl(r2). The interaction potential vanishes in the limit R → ∞, R →∞

∑ FvjlJMε (R)|ϕvj χvj ⟩|JM(lj12 (j1j2 ))⟩

where v ≡ v1v2 and j ≡ j1j2j12 denote the vibrational and rotational quantum numbers collectively. By substituting eq 11 into the Schrödinger equation with the Hamiltonian of eq 7 in the spherical coordinates yields a set of coupled equations for I the radial expansion coefficients FJε vjl (R), ÄÅ 2 2 ÉÑ ÅÅ ℏ d ÑÑ ℏ2l(l + 1) ÅÅ− + − ECÑÑÑÑFvjlJεI(R ) ÅÅ 2 2 ÅÅÇ 2μ dR ÑÑÖ 2μR

where r1 = q1 − q2 and r2 = q3 − q4 (qi, i = 1−4 denote the position vectors of atoms 1−4 in a space-fixed coordinate frame) are the relative vectors of molecules 1 (H2) and 2 (HCl), respectively, and R denotes the relative vector connecting the centers of mass of the two molecules, as shown in Figure 1. The kinetic energy terms T1(r1) and T2(r2) correspond to the relative motion of each of the two atoms in molecules 1 and 2, while T(R) refers to the kinetic energy for relative motion in R. The potential energy function V can be expressed as V (r1,r2,R) = Uint(r1,r2,R) + V1(r1) + V2(r2)

1 R

j12 j′12 , ll ′ , J , εI

(2J + 1)|δαlj

12

, α ′ l ′ j′12

− SαJlj, εI , α ′ l ′ j′12(E)|2 12

(13)

where k = 2μEC/ℏ and α ≡ v1j1v2j2 refers to the combined molecular state. 3-B. Computational Details. The scattering calculations are performed using the TwoBC55 and MOLSCAT56 codes. The latter is mainly used for comparisons with the 4D results of Lanza et al.30 who also used the MOLSCAT code. The fulldimensional and other reduced-dimensional calculations are performed using the TwoBC code. The matrix elements of the interaction potential Uint vjlJ ,,εvI ′ j ′ l′(R ) in eq 12 is evaluated by expanding the interaction 2

(10)

2

potential in a triple series of spherical harmonics52,53,57

in terms of which the total wave function is expanded at each value of the radial coordinate R. The functions |ϕ⟩ and |χ⟩ denote the vibrational eigenstates of molecules 1 and 2 with quantum numbers vi and ji denoting the vibrational and rotational levels. The computations are performed in the total angular momentum (a conserved quantity) representation in

Uint(r1,r2,R) =

∑ Aλ(r1,r2 ,R) Yλ(r1̂ ,r2̂ ,R̂) λ

(14)

where E


Article


Table 2. Rotational Relaxation Cross Sections (Å2) of HCl (j2) in Collisions with para-H2 (j1 = 0) at E = 100 and 500 cm−1 Computed Using the Rigid-Rotor Approximationa

Yλ(r1̂ ,r2̂ ,R̂ ) =

∑ ⟨λ1m λ λ 2m λ |λ12m λ 1

2

12

⟩Yλ1m λ (r1̂ ) Yλ2m λ (r2̂ ) Y λ*12m λ (R̂ )

mλ

1

2

12

j1j2 → j′1j′2

(15)

and λ ≡ λ1λ2λ12. This allows us to evaluate Uint vjlJ ,,εvI ′ j ′ l′(R ) as a

0,1 → 0,0 0,2 → 0,0

53

sum of products of a radial part and an angular part, Uint vjlJ ,,εvI ′ j ′ l′(R ) =

∑ Bαλ,α′(R)f jlJ ,;jλ′ l′ λ

0,1 0,2 0,3 0,3 0,3 0,5 0,5

(16)

Explicit expressions for Bλα,α′(R) and f jlJ ,;λj ′ l′ are given elsewhere53 and not reproduced here. Various parameters employed in our scattering calculations depend on the collision energy as well as the initial and final ro-vibrational states. For the angular expansion of the interaction potential in eq. 14 we adopt 0 ≤ λ1 ≤ 6 and 0 ≤ λ2 ≤ 12 for H2 and HCl, respectively, the same values adopted by Lanza et al.30 Only even values of λ1 contribute due the homonuclear symmetry of H2. The log-derivative matrix is propagated typically from R = 0.5 Å to 50 Å with a step size ΔR = 0.05 Å using Johnson’s method.54 The basis-set size in a J-block is specified by the maximum quantum numbers for the ro-vibrational states of the molecules. The minimum basis set employed in the present calculation is specified as v1,max = 0, j1,max = 2, v2,max = 0, and j2,max = 5, while their values depend on the initial and final ro-vibrational states and collision energy. The maximum values employed in the present calculations are j1,max = 4, j2,max = 6 and Jmax = 30. 3-C. Results. As discussed previously, Lanza et al. reported rotational transitions in HCl by ortho- and para-H2 using an ab initio 4D interaction potential computed in their study and the rigid-rotor approximation for the scattering calculations.30 For comparison, we report similar calculations on a reduceddimensional 4D interaction potential obtained from our full 6D PES. To make a one-to-one comparison, we employ the same parameters reported in the paper of Lanza et al.,30 and use the same code (MOLSCAT56). The reduced-dimensional 4D interaction potential is obtained by taking a cut of the full 6D PES at r1 = 0.76730 Å and r2 = 1.2859 Å for H2 and HCl, respectively. We note that, as mentioned in the paper of Lanza et al.,30 these bond distances are not the equilibrium values but the expectation values of the bond distances in the rovibrational ground states of the isolated molecules. The rotational constants of the molecules are BH2 = 60.853 cm−1 and BHCl = 10.5934 cm−1. In Table 2, we compare rotational quenching cross sections of HCl in the j2 = 1 rotational level due to collisions with paraH2 (j1 = 0). Our results are within 10−20% of Lanza et al.’s results, indicating the overall similarity of the reduced 4D interaction potential extracted from our full-dimensional PES. The agreement becomes better with increasing collision energy. The higher discrepancy at E = 100 cm−1 shows the sensitivity of the scattering properties to the fine details of the interaction potential at low energies. While our primary goal is to report full 6D scattering calculations for H2 + HCl collisions, the computations become expensive due to the large number of rotational levels of HCl needed to yield converged results and the strong anisotropy of the interaction potential (large number of terms in the angular expansion of the interaction potential). Therefore, first we compare the accuracy of different reduced-dimensional 4D

a

→ → → → → → →

0,0 0,0 0,0 0,1 0,2 0,3 0,4

Lanza et al.30

this study Total Energy (E = 100 cm−1) 9.795 3.628 Total Energy (E = 500 cm−1) 1.841 1.224 0.508 1.263 2.242 0.371 1.192

7.841 3.060 1.637 1.066 0.492 1.199 2.182 0.420 1.364

The number of the rotational basis for H2 is taken to be 2.

models against the full 6D results. In particular, Faure et al. recently showed that a vibrationally averaged 4D interaction potential yields results in excellent agreement with the 6D results for low energy rotational excitations of CO by H2.38 Similar effects have been demonstrated for OH−He and HCl− H collisions as well.36,37 Here, we consider three different reduced-dimensional 4D models. The first one corresponds to the conventional rigidrotor approximation in which the 4D potential is obtained as a cut of the 6D PES at the equilibrium molecular bond distances, r1,e = 0.74163 Å for H2 and r2,e = 1.3512 Å for HCl; thus we use the symbol Uint(r1,e,r2,e) for the conventional rigid-rotor approximation. As a variant of this method, another 4D potential is obtained by the cut at the expectation values of the bond distances in their ro-vibrational ground states, ⟨r1⟩ = ⟨ϕ00| r1|ϕ00⟩ and ⟨r2⟩ = ⟨χ00|r2|χ00⟩. This is the closest approximation to the 4D potential of Lanza et al.30 As mentioned above, this method has been recognized as a variant of the rigid-rotor approximation giving better agreement than the conventional treatment with the equilibrium bond distances.30,38 We use Uint(⟨r1⟩,⟨r2⟩) as the symbol for this method. As described above, the actual expectation values are ⟨r1⟩ = 0.76730 Å and ⟨r2⟩ = 1.2859 Å. Third, the vibrationally averaged potential,38,58 is obtained by averaging the 6D potential over r1 and r2 as ⟨Uint⟩(X ) = ⟨ϕ00χ00 |Uint(r1,r2 ,X )|ϕ00χ00 ⟩

(17)

where X collectively indicates the remaining four coordinates R, θ1, θ2, and φ. We refer this method with the symbol ⟨Uint⟩(X). Figure 5 presents cross sections for the rotational excitation of HCl(v2 = 0, j2 = 0 → v′2 = 0, j′2 = 1) due to collisions with para-H2(v1 = 0, j1 = 0) as a function of the collision energy. We observe that the vibrationally averaged 4D potential ⟨Uint⟩(X) (red dashed curve) yields nearly identical results as that from the full 6D calculations (black curve) in the whole energy region, including the sharp resonances features. Two rigidrotor calculations (green and blue curves) also show generally good agreement with full-dimensional results. In particular, the conventional rigid-rotor approximation with the equilibrium bond distances (blue curve) displays the largest error, consistent with the trend observed in previous studies of CO−H2.38,58 Similar comparisons for rotational quenching of HCl(v2 = 0, j2 = 1 → v′2 = 0, j′2 = 0) are shown in Figure 6, extending to F


Article


the cold regime. The onset of the Wigner threshold regime in which the cross section diverges inversely as the velocity is evident for collision energies below 10−4 cm−1. In the entire energy region, the vibrationally averaged 4D results show excellent agreement with the full 6D calculations, including the resonant region between 10−3 and 0.1 cm−1 where the cross sections vary by an order of magnitude due to a p-wave shape resonance in the incoming channel. Like in Figure 5, we observe significant deviations between the full-dimensional calculations and the rigid-rotor approximations. In Figure 6a, the rigid-rotor approximation with the expectation values of the bond distances (green curve in the figure) yields results that are significantly shifted in the resonance position at around 10−1 cm−1 compared to the 6D results. The conventional rigid-rotor approximation with equilibrium bond distances (blue curve) also shows a shift in the resonance position at around 10−2 cm−1 but a smaller shift compared to its vibrationally averaged variant. These differences are most likely due to the slight changes in the bound-state structure of the PESs in the two 4D approximations. However, as shown in Figure 6b, the conventional rigid-rotor approximation displays the largest error in the higher energy region. The systematic deviations observed with the 4D rigid-rotor approximations indicate the importance of adequately accounting for the effect of stretching motions in computing rotational transitions. The excellent agreement between the full 6D results and those from the vibrationally averaged 4D potential shows that at least within the v = 0 vibrational level it is possible to evaluate accurate rotational transitions with the vibrationally averaged potential within the 4D treatment. This method is effective even for light molecules such as HCl that are more anharmonic than heavier molecules such as CO. In the vibrationally averaged 4D calculations the vibrational basis, |ϕv1j1⟩ and |χv2j2⟩ in eq 10 are approximated as |ϕ00⟩ and |χ00⟩ in the evaluation of the matrix elements of the interaction potential. All of the results presented above are for rotational transitions of HCl within the vibrational ground state (v2 = 0) with para-H2(v1 = 0, j1 = 0) and rotational ground state (j2 = 0) for either the initial or final state of the HCl molecule. Therefore, the most important matrix element of the interaction potential between the initial and final ro-vibrational basis states, includes the approximation only in the vibrational state of the HCl as |χ00⟩ for |χ10⟩. Both |χ00⟩ and |χ10⟩ are vibrational ground states with different j2; thus, ⟨r2|χ00⟩ and ⟨r2|χ10⟩ are very similar, leading to the excellent performance of the vibrationally averaged 4D potential method. Obviously, this approximation may not be good for processes involving transitions between different vibrational states due to the markedly different behavior of the vibrational wave functions. It also would not work well for transitions from/to highly rotationally excited states, and for molecules whose interaction potential between atoms are very weak because the centrifugal term can have significant effect on the behavior of the vibrational wave functions of the molecules. Finally, in Figure 7, we present a comparison between the full 6D calculations and the vibrationally averaged 4D method for cold collisions for (a) (j2 = 3 → j′2 = 0) and (b) (j2 = 3 → j′2 = 2) with para-H2(v1 = 0, j1 = 0). The excellent agreement indicates comparable values of ⟨r2|χ00⟩, ⟨r2|χ20⟩, and ⟨r2|χ30⟩ due to the small effect of the centrifugal terms on the v2 = 0 vibrational wave functions of HCl. These results illustrate that for low-lying rotational transitions within the vibrational

Figure 5. Cross section for the rotational excitation of HCl(v2 = 0, j2 = 0 → v′2 = 0, j′2 = 1) due to collisions with para-H2(v1 = 0, j1 = 0). The solid black curve shows results from the 6D calculations. The red dashed curve is for the vibrationally averaged 4D potential. The rigidrotor calculations using the two 4D cuts of the 6D PES at the expectation values of the bond distances and at the equilibrium bond distances are shown by the green and blue curves, respectively.

Figure 6. Cross sections for rotational quenching of HCl(v2 = 0, j2 = 1 → v′2 = 0, j′2 = 0) due to collisions with para-H2(v1 = 0, j1 = 0): (a) cold regime; (b) low energy region. The solid black curve shows results from the 6D calculations. The red dashed curve is for the vibrationally averaged 4D potential. The rigid-rotor calculations using the two 4D cuts of the 6D PES at the expectation values of the bond distances and at the equilibrium bond distances are shown by the green and blue curves, respectively. G


Article


provide valuable insights into the vibrational effects in energy transfer.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpca.9b05958. Details of the long-range interaction potential and Table S1 giving well depths and geometric parameters (PDF)

■

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. ORCID

Hua Guo: 0000-0001-9901-053X Present Address #

Institute of Modern Physics, Northwest University, Xi’an, Shaanxi, People’s Republic of China. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was supported in part by an ARO MURI grant (W911NF-19-1-0283 to H.G. and N.B.) and NSF grant (PHY1806334 to N.B.). The computation at UNM was done at Center for Advanced Research Computing (CARC).

■

REFERENCES

(1) Crim, F. F. Vibrational state control of bimolecular reactions: Discovering and directing the chemistry. Acc. Chem. Res. 1999, 32, 877−884. (2) Liu, K. Vibrational control of bimolecular reactcions with methane with mode-, bond-, and stereo-selectivity. Annu. Rev. Phys. Chem. 2016, 67, 91−111. (3) Chadwick, H.; Beck, R. D. Quantum state−resolved studies of chemisorption reactions. Annu. Rev. Phys. Chem. 2017, 68, 39−61. (4) Polanyi, J. C. Concepts in reaction dynamics. Acc. Chem. Res. 1972, 5, 161−168. (5) Guo, H.; Jiang, B. The sudden vector projection model for reactivity: Mode specificity and bond selectivity made simple. Acc. Chem. Res. 2014, 47, 3679−3685. (6) Guo, H.; Liu, K. Control of chemical reactivity by transition state and beyond. Chem. Sci. 2016, 7, 3992−4003. (7) Gaubatz, U.; Rudecki, P.; Schiemann, S.; Bergmann, K. Population transfer between molecular vibrational levels by stimulated Raman scattering with partially overlapping laser fields. A new concept and experimental results. J. Chem. Phys. 1990, 92, 5363− 5376. (8) Mukherjee, N.; Zare, R. N. Stark-induced adiabatic Raman passage for preparing polarized molecules. J. Chem. Phys. 2011, 135, 024201. (9) Chadwick, H.; Hundt, P. M.; van Reijzen, M. E.; Yoder, B. L.; Beck, R. D. Quantum state specific reactant preparation in a molecular beam by rapid adiabatic passage. J. Chem. Phys. 2014, 140, 034321. (10) Krems, R. V. Cold controlled chemistry. Phys. Chem. Chem. Phys. 2008, 10, 4079−4092. (11) Quéméner, G.; Julienne, P. S. Ultracold molecules under control! Chem. Rev. 2012, 112, 4949−5011. (12) Balakrishnan, N. Perspective: Ultracold molecules and the dawn of cold controlled chemistry. J. Chem. Phys. 2016, 145, 150901. (13) Perreault, W. E.; Mukherjee, N.; Zare, R. N. Quantum control of molecular collisions at 1 K. Science 2017, 358, 356−359.

Figure 7. Cross sections for rotational quenching of HCl (a) (v2 = 0, j2 = 3 → v′2 = 0, j′2 = 0) and (b) (v2 = 0, j2 = 3 → v′2 = 0, j′2 = 2) induced by collisions with para-H2(v1 = 0, j1 = 0). The solid black curve corresponds to the full 6D calculations and the red dashed curve is from the 4D vibrationally averaged potential.

ground state of HCl induced by para-H2 one can use the 4D vibrationally averaged potential without compromising on the accuracy.

4. SUMMARY AND CONCLUSIONS Extensive ab initio computations are reported for the H2−HCl system using the coupled-cluster singles, doubles, and perturbative triples level of theory. A global 6D PES is constructed using the PIP-NN method in the short-range, but analytical functions in the long-range. The accuracy of the potential is demonstrated by comparing rotational cross sections against recent results of Lanza et al.30 within a 4D rigid-rotor model. Different variants of the 4D potentials are constructed and tested for rotational transitions against the full 6D potential. It is found that, for rotational transitions in the ground vibrational level of HCl, a vibrationally averaged 4D potential yields results in quantitative agreement with the fulldimensional calculations. We plan to report rovibrational calculations involving excited vibrational levels of both molecules and ortho-H2 in a future publication. Such calculations will be numerically challenging, but recent advances in quantum scattering, such as reduction of coupled helicity channels,59 renders such calculations within the reach of the current computational capabilities. Such studies will H


Article


(36) Kalugina, Y.; Lique, F.; Marinakis, S. New ab initio potential energy surfaces for the ro-vibrational excitation of OH(X2Π) by He. Phys. Chem. Chem. Phys. 2014, 16, 13500−13507. (37) Lique, F. Communication: Rotational excitation of HCl by H: Rigid rotor vs. reactive approaches. J. Chem. Phys. 2015, 142, 241102. (38) Faure, A.; Jankowski, P.; Stoecklin, T.; Szalewicz, K. On the importance of full-dimensionality in low-energy molecular scattering calculations. Sci. Rep. 2016, 6, 28449. (39) Yang, B.; Zhang, P.; Wang, X.; Stancil, P. C.; Bowman, J. M.; Balakrishnan, N.; Forrey, R. C. Quantum dynamics of CO−H2 in full dimensionality. Nat. Commun. 2015, 6, 6629. (40) Yang, B.; Balakrishnan, N.; Zhang, P.; Wang, X.; Bowman, J. M.; Forrey, R. C.; Stancil, P. C. Full-dimensional quantum dynamics of CO in collision with H2. J. Chem. Phys. 2016, 145, 034308. (41) Yang, B.; Wang, X. H.; Stancil, P. C.; Bowman, J. M.; Balakrishnan, N.; Forrey, R. C. Full-dimensional quantum dynamics of rovibrationally inelastic scattering between CN and H2. J. Chem. Phys. 2016, 145, 224307. (42) Yang, D.; Huang, J.; Zuo, J.; Hu, X.; Xie, D. A full-dimensional potential energy surface and quantum dynamics of inelastic collision process for H2−HF. J. Chem. Phys. 2018, 148, 184301. (43) Huang, J.; Yang, D.; Zhou, Y.; Xie, D. A new full-dimensional ab initio intermolecular potential energy surface and vibrational states for (HF)2 and (DF)2. J. Chem. Phys. 2019, 150, 154302. (44) Knizia, G.; Adler, T. B.; Werner, H.-J. Simplified CCSD(T)F12 methods: Theory and benchmarks. J. Chem. Phys. 2009, 130, 054104. (45) Dunning, T. H. Gaussian basis sets for use in correlated molecular calculations. I. The atoms boron through neon and hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (46) Werner, H. J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M. Molpro: a general-purpose quantum chemistry program package. WIREs Comput. Mol. Sci. 2012, 2, 242−253. (47) Jiang, B.; Guo, H. Permutation invariant polynomial neural network approach to fitting potential energy surfaces. J. Chem. Phys. 2013, 139, 054112. (48) Li, J.; Jiang, B.; Guo, H. Permutation invariant polynomial neural network approach to fitting potential energy surfaces. II. Fouratom systems. J. Chem. Phys. 2013, 139, 204103. (49) Buckingham, A. D. Permanent and induced molecular moments and long-range intermolecular forces. Adv. Chem. Phys. 1967, 12, 107−142. (50) Takayangi, K. The production of rotational and vibrational transitions in encounters between molecules. Adv. At. Mol. Phys. 1965, 1, 149−194. (51) Green, S. Rotational excitation in H2-H2 collisions: Closecoupling calculations. J. Chem. Phys. 1975, 62, 2271−2277. (52) Alexander, M. H.; DePristo, A. E. Symmetry considerations in the quantum treatment of collisions between two diatomic molecules. J. Chem. Phys. 1977, 66, 2166−2172. (53) Quéméner, G.; Balakrishnan, N. Quantum calculations of H2− H2 collisions: From ultracold to thermal energies. J. Chem. Phys. 2009, 130, 114303. (54) Johnson, B. R. The multichannel log-derivative method for scattering calculations. J. Comput. Phys. 1973, 13, 445−449. (55) Krems, R. V. TwoBC Quantum Scattering Program; University of British Columbia: Vancouver, 2006. (56) Hutson, J. M.; Green, S. MOLSCAT v.14; Engineering and Physical Sciences Research Council: Swindon, 1994. (57) Gray, C. G. On the theory of multipole interactions. Can. J. Phys. 1968, 46, 135−139. (58) Jankowski, P.; Szalewicz, K. A new ab initio interaction energy surface and high-resolution spectra of the H2−CO van der Waals complex. J. Chem. Phys. 2005, 123, 104301. (59) Yang, D.; Hu, X.; Zhang, D. H.; Xie, D. An improved coupledstates approximation including the nearest neighbor Coriolis couplings for diatom-diatom inelastic collision. J. Chem. Phys. 2018, 148, 084101.

(14) Perreault, W. E.; Mukherjee, N.; Zare, R. N. Cold quantumcontrolled rotationally inelastic scattering of HD with H2 and D2 reveals collisional partner reorientation. Nat. Chem. 2018, 10, 561− 567. (15) Murrell, J. N.; Carter, S.; Farantos, S. C.; Huxley, P.; Varandas, A. J. C. Molecular Potential Energy Functions. Wiley: Chichester, 1984. (16) Braams, B. J.; Bowman, J. M. Permutationally invariant potential energy surfaces in high dimensionality. Int. Rev. Phys. Chem. 2009, 28, 577−606. (17) Czakó, G.; Bowman, J. M. Reaction dynamics of methane with F, O, Cl, and Br on ab initio potential energy surfaces. J. Phys. Chem. A 2014, 118, 2839−2864. (18) Li, J.; Jiang, B.; Song, H.; Ma, J.; Zhao, B.; Dawes, R.; Guo, H. From ab initio potential energy surfaces to state-resolved reactivities: The X + H2O ↔ HX + OH (X = F, Cl, and O(3P)) reactions. J. Phys. Chem. A 2015, 119, 4667−4687. (19) Behler, J. Constructing high-dimensional neural network potentials: A tutorial review. Int. J. Quantum Chem. 2015, 115, 1032−1050. (20) Dawes, R.; Ndengué, S. A. Single- and multireference electronic structure calculations for constructing potential energy surfaces. Int. Rev. Phys. Chem. 2016, 35, 441−478. (21) Jiang, B.; Li, J.; Guo, H. Potential energy surfaces from high fidelity fitting of ab initio points: The permutation invariant polynomial-neural network approach. Int. Rev. Phys. Chem. 2016, 35, 479−506. (22) Aoiz, F. J.; Bañares, L.; Herrero, V. J. The H + H2 reactive system. Progress in the study of the simplest reaction. Int. Rev. Phys. Chem. 2005, 24, 119−190. (23) Zhang, D. H.; Guo, H. Recent advances in quantum dynamics of bimolecular reactions. Annu. Rev. Phys. Chem. 2016, 67, 135−158. (24) Hinde, R. J. A six-dimensional H2−H2 potential energy surface for bound state spectroscopy. J. Chem. Phys. 2008, 128, 154308. (25) Croft, J. F. E.; Balakrishnan, N.; Huang, M.; Guo, H. Unraveling the stereodynamics of cold controlled HD-H2 collisions. Phys. Rev. Lett. 2018, 121, 113401. (26) Jambrina, P. G.; Croft, J. F. E.; Guo, H.; Brouard, M.; Balakrishnan, N.; Aoiz, F. J. Stereodynamical control of a quantum scattering resonance in cold molecular collisions, Phys. Rev. Lett. 2019, in press. (27) Neufeld, D. A.; Wolfire, M. G. The chemistry of interstellar molecules containing the halogen elements. Astrophys. J. 2009, 706, 1594−1604. (28) Roueff, E.; Lique, F. Molecular excitation in the interstellar medium: Recent advances in collisional, radiative, and chemical processes. Chem. Rev. 2013, 113, 8906−8938. (29) Gerin, M.; Neufeld, D. A.; Goicoechea, J. R. Interstellar hydrides. Annu. Rev. Astron. Astrophys. 2016, 54, 181−225. (30) Lanza, M.; Kalugina, Y.; Wiesenfeld, L.; Lique, F. Nearresonant rotational energy transfer in HCl−H2 inelastic collisions. J. Chem. Phys. 2014, 140, 064316. (31) Lanza, M.; Lique, F. Hyperfine excitation of linear molecules by para- and ortho-H2: Application to the HCl−H2 system. J. Chem. Phys. 2014, 141, 164321. (32) Lanza, M.; Kalugina, Y.; Wiesenfeld, L.; Faure, A.; Lique, F. New insights on the HCl abundance in the interstellar medium. Mon. Not. R. Astron. Soc. 2014, 443, 3351−3358. (33) Kama, M.; Caux, E.; López-Sepulcre, A.; Wakelam, V.; Dominik, C.; Ceccarelli, C.; Lanza, M.; Lique, F.; Ochsendorf, B. B.; Lis, D. C.; et al. Depletion of chlorine into HCl ice in a protostellar core. Astron. Astrophys. 2015, 574, A107. (34) Anderson, D. T.; Schuder, M.; Nesbitt, D. J. Large-amplitude motion in highly quantum clusters: high-resolution infrared absorption studies of jet-cooled H2−HCl and H2−DCl. Chem. Phys. 1998, 239, 253−269. (35) Alkorta, I.; Elguero, J.; Del Bene, J. E. An ab initio investigation of the properties of H2:HX hydrogen-bonded complexes. Chem. Phys. Lett. 2010, 489, 159−163. I


Globally Accurate Full-Dimensional Potential ... - ACS Publications

Recommend Documents