Nonadiabatic Dynamics on the Two Coupled Electronic PESs: The H+

J. Phys. Chem. A 2010, 114, 10357–10366

10357

Nonadiabatic Dynamics on the Two Coupled Electronic PESs: The H+ + O2 System F. George D. Xavier* Department of Chemistry, Indian Institute of Technology Madras, Chennai 600036, India ReceiVed: July 29, 2010; ReVised Manuscript ReceiVed: August 13, 2010

Multistate adiabatic and diabatic PESs were computed for the H+ + O2 collision system in Jacobi coordinates, (R,r,γ) using the cc-pVTZ basis set and the ic-MRCI level of theory. In addition, all possible interaction potentials and nonadiabatic coupling matrix elements among those different electronic states were also computed. Comparisons with earlier computed interaction potentials were made wherever possible, and the differences between them is attributed to the multistate diabatization and the chosen level of theory and basis set. Focusing our attention on the ground-state (GS) and the first excited-state (ES) PES, quantum dynamics were performed using the 2 × 2 diabatic potential submatrix obtained from the multistate (four) diabatic potential matrix within the VCC-RIOSA scheme at two experimentally reported collision energies, Ecm ) 9.5 and 23 eV. The scattering quantities were computed for two experimentally observed collision processes, namely, the inelastic vibrational excitation (IVE), H+ + O2 + 3 + 3 2 (X3Σg ,V ) 0) f H + O2 (X Σg ,V′), and the vibrational charge transfer (VCT), H + O2 (X Σg ,V ) 0) f H ( S) 2 + O (X Πg,V′′). Comparisons were made with experimental results and found an overall improvement relative to the earlier computed results, and the discrepancies, if any, could be brought down to minimum by further modification in employed ab initio PESs and the interaction potential. I. Introduction Proton-molecule collisions occur dominantly in the interstellar spaces, leading to the identification of many bound protonated molecular ions like H3+, N2H+, HCO+, HOC+, HCS+, HO2+, and so forth in them.1,2 They are proposed to be formed through various elastic/inelastic or charge-transfer processes. As a result, the proton-molecule interaction has been the subject of many experimental and theoretical studies, and a wealth of varied experimental information on inelastic vibrational excitation (IVE) and vibrational charge transfer (VCT) excitation has become available over the years using the molecular beam and the proton energy-loss spectroscopy techniques.3-9 The focus has been on collision energies in the range of 0-30 eV, where extensive vibrational-rotational excitations of the target molecules occur along with possible VCT processes. There have been several experimental studies on proton collisions with diatomic and polyatomic targets such as H2, N2, O2, CO, NO, HF, HCl, CH4, SF6, and so forth leading to an interesting finding of marked selectivity in vibrational excitations of apparently similar molecules. The dynamics of complex formation via translation-to-vibration energy exchange and types of resonances exhibited by the complex has also been studied.10 In particular, the experiments6,9 were performed on N2, CO, NO, and O2 at Ecm ) 9.5 eV in the range of scattering angles of 5 e θcm e 21°, and the results were compared in terms of transition probability (P0fV′), average vibrational energy transfer (∆Evib), and so forth. Noll and Toennies7 did detailed experiments on the H+ + O2 system for both IVE and the VCT channels at Elab ) 23.7 eV in the range of 0 e θlab e 11°, and the results obtained were discussed in detail. From the theoretical side, the collision systems studied in the past include H+ + CO,11 H+ + NO,12 H+ + H2,13 and H+ + N2.14 There had been few early ab initio studies on the H+ + O215-17 system and then came studies18,19 involving global potential energy surfaces (PESs) of the ground state (GS) and * E-mail: [email protected].

first excited state (ES) using the projected valence bond (PVB) approach and the semiempirical diatomics-in-molecule (DIM) approach. Quantum dynamics on these two sets of PESs qualitatively agree with the experiments.20,21 Recently, the quantum dynamics using global ab initio PESs were performed in our group,22-24 and they were found to explain the experimental results with limited success. The aim of the present study is to bring about an improvement over the previously reported studies with mimimum computational effort by taking out 2 × 2 diabatic potential submatrix from the computed 4 × 4 matrix and subjecting it to dynamical computation. In this way, this work differs from the earlier work in our group,22-24 which was fully based on the 2 × 2 diabatic matrix. This paper is orgainzed as follows. In section II, we provide the details of ab initio computation, in section III, the details of adiabatic and diabatic potentials, in section IV, vibronic couplings, and in section V, quantum dynamics at two collision energies, Ecm ) 9.5 and 23 eV, which includes five subsections on different scattering quantities, followed by summary and conclusion in section VI. II. Electronic Structure Calculations Ab initio calculations involving diabatization were carried out on the bound HO+ 2 molecule for four lowest electronic states b,r in Jacobi coordinates (R b,γ) for the triplet spin state using the r is the molecular bond length of O2, b R MOLPRO package.25 b is the vector connecting the center of mass of O2 and the projectile (H+), and γ is the angle between these two vectors, defined by

cos γ )

( ) b Rb r b |R || b| r

Dunning’s26 correlation-consistent polarized valence triple-ζ (ccpVTZ) basis set, applied to H, C, and O atoms and multiconfiguration self-consistent field (MCSCF)27 followed by the internally contracted multireference configuration interaction

10.1021/jp107098v  2010 American Chemical Society Published on Web 09/01/2010

10358

J. Phys. Chem. A, Vol. 114, No. 38, 2010

Xavier TABLE 1: Number of Reference Configurations (Nref) and Number of Roots (Nroot) Treated in Each Irreducible Representation and the Corresponding Number of Generated (Ntot) and Selected (Nsel) Symmetry-Adapted Functions for a Threshold of 0.32 × 10-6 Hartree at R ) 2aoa γ 3

0° 45° 90°

state

Nref/Nroot

Ntot

Nsel

3

382/4 304/2 686/4 343/4 343/2

10744181 9933007 21477678 10721927 10722093

326572 240744 642952 322036 238874

B1( Π) A 2( 3Σ - ) 3 A′′ 3 B1 3 A2 3

Note that R is the distance between H+ and the O2 center of mass. a

Figure 1. Correlation diagram of HO+ 2.

with single and double excitations (ic-MRD-CI)28 level of theory, which does not involve cluster-corrected energies such as Davidson’s correction or Pople’s correction, were used in the computations. The grid points in the scattering coordinates depend on the approach of H+ toward the molecular target, defined by γ, which is equally spaced from 0 to 90° with the step of 15° interval. For example, when γ ) 0°, R ) 2.0-7.0(0.2),8-15(1.0),whenγ)45°,R)0.8-7.0(0.2),8-15(1.0), and when γ ) 90°, R ) 0.2-7.0(0.2),8-15(1.0). The grid along r varied as 1.5-3.5(0.1). The number in the parentheses denotes the step size in the stated interval and is expressed in units of bohr. This grid pattern was deliberately chosen to include the interaction region fully in the scattering calculations for all of the angular approaches of H+. On average, 6000 ab initio data points were computed for 7 γ values, a computationally expensive task. The equilibrium bond distance (req) of O2 was fixed at 2.293ao. The reference geometry for the diabatic calculation was fixed at R ) 16ao, at which both adiabatic and diabatic potentials remain the same throughout the entire calculations. The calculations at displaced geometries were performed relative to that of the reference geometry. All the calculations were done in the highest abelian subgroup of the point group to which the molecule belongs. For example, when γ ) 0°, the molecule belongs to C∞V, but the calculations are done only in the C2V group, which is the abelian subgroup of C∞V. When 0 < γ < 90°, the molecule belongs to the Cs group, and the calculations are done in this symmetry as it is an abelian group. When γ ) 90°, the molecule belongs to the C2V group, and the calculations are done in the same group. Since there is a change of symmetry as γ tends from 0 to 90°, it is instructional to show the correlation of molecular orbitals as a function of γ in Figure 1. Note that the angular dependence cuts sharp near the 90° end. The chosen basis set generates 74 molecular orbitals, out of which the lowest 9 are occupied in the Hartree-Fock (HF) level. The nine HF orbitals are arranged in the order of increasing energy in the ground-state configuration as follows: (1σ+)2,(2σ+)2,(3σ+)2,(4σ+)2,(5σ+)2,(1π)2xy,(1π)2xz(2π)1xy, 1 for γ ) 0° (linear molecule), (1a′)2,(2a′)2,(3a′)2,(4a′)2, (2π)xz 2 (5a′) ,(6a′)2,(1a′′)2,(7a′)1,(2a′′)1 for γ ) 45°, and (1a1)2,(1b2)2, (2a1)2,(2b2)2,(3a1)2,(4a1)2,(1b1)2,(3b2)1,(1a1)1 for γ ) 90°. The superscript denotes the number of electrons, and the + sign denotes the reflection symmetry of molecular orbitals across the plane containing the three nuclei. The two lowest pairs of electrons (four electrons) were treated as core electrons in the MCSCF level, and the remaining 12 electrons were involved in the excitation. The

convergence is sensitive to the number of wave functions in the primary configuration space (p-space) and can be forced upon by increasing the value of p-space (default is 0.4). The average weightage (0.25) given to each electronic state remains equal. In the CI calculations, the two lowest orbitals are kept frozen, and remaining orbitals are taken into account in the electron correlation process. A small selection threshold29 of 0.32 × 10-6 hartree has been used in the present treatment. More details of the present MRD-CI calculations can be found in Table 1. The radial coupling matrix elements are obtained by a finite-difference scheme with an increment of 0.0002ao.30 III. Adiabatic and Diabatic Potentials The diabatization approach used in the generation of diabatic potential energy curves (PECs) or potential energy surfaces (PESs) is extensively discussed in the literature31 regarding its accuracy and validity. The approach that we follow is the one already used by Werner and co-workers32 in the study of photodissociation of H2S on the coupled PESs to explain the experimental observations. We present in Figure 2 adiabatic as well as diabatic PECs of the two lowest electronic states for three molecular orientations. The diabatic potentials are obtained from the diabatization of the four lowest electronic states, but the potentials belonging to the two lowest electronic states have been shown and will be subjected to the scattering calculations. The left panel shows PECs as a function of R at req, with diabatic PECs shown by the dashed line for γ ) 45°, and the right panel shows the PECs as a function of r at fixed R. The characteristics of the PECs (left panel) are obvious from the diagram. Each PEC of the bound HO+ 2 correlates to the specific electronic state of bound O2 in the asymptotic region. The GS and the first ES in the γ ) 0 and 90° orientations correlate 2 2 to H+ + O2(X3Σg ) (IVE channel) and the H( S) + O(X Πg) (VCT channel), and just the reverse is the correlation pattern in γ ) 45°. This happens because the two lowest states in γ ) 0 and 90° are of different symmetries. Due to this reason, there is no radial coupling between them and hence no diabatization. Therefore, diabatic and adiabatic PECs merge together for these two orientations and differ significantly for other angular approaches, 0 < γ < 90°. The curves are running parallel to each other at large R, and the γ-dependent energy difference at R ) 15ao between the IVE and the VCT channels is ∼1.8 eV, which is comparable with the reported experimental exoergic value of 1.52 eV.9 At R ) ∞, the computed exoergic value (a γ-independent quantity) turns out to be 2.02 eV. The shallow well present in the first ES PECs at extended R supports the weak van der Waals complex of HO2+, and its depth is γ-dependent. At γ ) 15°, the depth is ∼0.67 eV, a value comparable to the value of 1.0 eV computed earlier.15 Both diabatic and adiabatic PECs merge together at large R regions, as

Nonadiabatic Dynamics and the H+ + O2 System

J. Phys. Chem. A, Vol. 114, No. 38, 2010 10359

Figure 2. Adibatic PECs (solid line) of the GS and the first ES for three molecular orientations as a function of R and r at a fixed r and R, respectively. The diabatic PECs (dashed line) are shown for γ ) 45°.

seen for γ ) 45°. The hump-like structure that exists in the PEC of the first ES is characteristic of the interaction with the second ES of the same symmetry. Another distinct feature is the crossing of curves observed in all orientations. At γ ) 0 and 90°, where there is no diabatization, it occurs due to the difference in symmetry, whereas at γ ) 45°, this occurs in diabatic curves and is absent in the adiabatic curves. The point of crossing also changes with respect to γ (it decreases in R and increases in r as γ tends to 90°). In the plots as a function of r at a fixed R (right panel), we can see that there are crossings at two places in γ ) 0° with one at the extended r region; this is not the case with the other two angles. The hump denoting the interaction with the second ES at γ ) 0° is fully absent at γ ) 90°, suggesting the active role of the higher-lying ES in the scattering dynamics in the collinear approach of the projectile. The irregularity seen in the diabatic curve (dashed curve) at γ ) 45° at the extended r may arise from the involvement of higherlying excited states. In Figure 3, we show adiabatic PESs of the GS and the first ES as a function of R and r for the γ ) 0 and 90° orientations scaled relative to the energy of the asymptotic products, H+ + O2(X3Σg-), to be zero. The arguments provided above to PECs

are also applicable to these PESs. The PESs at γ ) 0° undergo crossing along R and r coordinates, whereas the same at γ ) 90° undergo crossing only along the R coordinate alone, as can be seen. The interaction well depth in the GS is deeper and sharper in the case of the former than that in the latter, and they are reported in the present work to be 3.78 and 3.15 eV, respectively. The deepest well occurs at γ ) 42.7° with 4.84 eV, corresponding to the most stable geometry of HO2+, whose bond length and bond angle are reported in Table 2 in bold fonts and compared with previously reported theoretical results. The off-diagonal elements (Vij, i * j) in the diabatic potential matrix (coupling potential or interaction potential) connect any two electronic states. They are computed by evaluating the integral Vij ˆ el|ψdj 〉, where ψi,jd is the diabatic wave function. We ) Vji ) 〈ψdi |H present the interaction potential V12 between the GS and the first ES as a function of R at different r for γ ) 45° and compare with those computed by Sidis et al.20 in Figure 4. We present, in Figure 5, the interaction potential as a function R at req for various γ and compare with Sidis’s work. The difference between the present work and that of Sidis et al. arises from the difference in the diabatization procedure adopted, which, in the present work, is

10360

J. Phys. Chem. A, Vol. 114, No. 38, 2010

Xavier

Figure 4. Comparison of the coupling potential (V12) in the present work with that computed by Sidis et al.20 for γ ) 45°.

Figure 3. Ab initio smoothed adiabatic PESs of the GS (13Σ-/13B1) and the first ES (13Π/13A2) for γ ) 0 and 90° orientations as a function of R and r.

TABLE 2: Computed Equilibrium Geometry Parameters in the Valence Coordinates of the HO2+ Ion in the GS (13A′) theory

r(O-O)/(Å)

r(O-H)/(Å)

R(O-O-H)/(deg)

this work earlier work36 earlier work37 earlier work38 earlier work39 earlier work2

1.239 1.22 1.215 1.237 1.231 1.227

1.01 0.995 0.985 1.007 1.009 1.01

111.52 111.54 111.0 118.8 112.461 112.74

based on the CI approach while in the earlier work, it was based on unitary transformation of selected adiabatic basis sets, named M2 diabatization. We are of the view that the differences existing in the coupling potential will alter the course of the dynamics calculation and will provide us with better results. The first-order nonadiabatic coupling matrix elements (NACME) can be computed using the three-point numerical formula given below

〈| |〉 ψa1

d a 1 ψ ) 〈ψa (R + ∆R)|ψa2(R0 - ∆R)〉 dR 2 2∆R 1 0

(1)

Figure 5. Comparison of the coupling potential (V12) for various γ at req with that of Sidis et al.20

where ψa is the adiabatic wave function, ∆R is the small increment fixed as 0.0002ao, and R0 is the point at which the radial coupling is computed. The computed radial couplings are shown in Figure 6 for γ ) 45° as a function of R and r and then of γ and R at req, and the sharp changes with sign reversal can be seen in both of the plots. This is due to the phase changes of wave functions in CI calculations, apart from the presence of low-lying ESs. The radial coupling is present in all γ’s, less in the middle and higher at both ends along γ, and it goes to 0 at large R for the entire range of r, suggesting that the two states are widely separated. In the region where it is high, irrespective of sign reversal, the electronic states are closely spaced. The scattering calculation requires data of PESs in the longrange region much beyond whatever raw ab initio points computed. Because the task of computing raw points is timeconsuming and tedious, we generate the points by modeling those regions using multipolar expression given below and connect them with existing ab intio raw data, and thus, the R range extends up to 50 bohr, keeping the r range unchanged for each γ.

Nonadiabatic Dynamics and the H+ + O2 System

J. Phys. Chem. A, Vol. 114, No. 38, 2010 10361

Figure 6. The radial coupling between the GS and the first ES for γ ) 45° and req.

Vas(R, r;γ) ≈

R0(r) R2(r) Q(r) P2(cos γ) P2(cos γ) 3 4 R 2R 2R4 (2)

where Vas is the asymptotic long-range potential, Q(r) is the quadrupole moment, R0(r) and R2(r) are the polarizability components, and the Pi’s are Legendre polynomials. We have computed the Q(r), R0(r), and R2(r) as a function of r of O2(X3Σg-) and fitted them with the function N

x(r) )

∑ Cx(i)( 1r )

i

(3)

i)0

where x(r) stands for Q(r), R0(r), or R2(r). Cx are the coefficients used in the fitting. IV. Vibronic Couplings The vibronic couplings are computed to study the magnitude of electronic transitions of HO2+ between the GS and the first ES during the vibrational motion of O2, and this information is needed to compute the vibrational mode selective scattering quantities. The numerical formula to compute this quantity is given below

VV(R;γ) ) 〈χV'(r)|V12(R, r, γ)|χV''(r)〉

Figure 7. Vibronic couplings for γ ) 45°.

Numerov method34 for seven equally spaced orientations of γ between 0 and 90° for each partial wave (l). The maximum partial waves (lmax) required depend on the collision energies and the type of channel. We computed scattering properties in the center of mass (cm) frame of reference and compared them with experimental results, which are available in the laboratory frame of reference (lab) after converting them into the cm frame of reference through the following relation

(4) sin θcm )

where VV(R;γ) is the state-dependent vibronic coupling denoted by the suffix V, V′ refers to the vibrational quantum number of O2 (0, 1, 2, ...), and V′′ refers to the vibrational quantum numer of O+ 2 (0, 1, 2, ...). V12(R,r,γ) is the interaction potential between GS and the first ES of HO2+. χV′(r) is the vibrational wave function of O2, and χV′′(r) is that of O2+. Approximately 20 vibrational states of O2 and O2+ were included in the computation. This integral is summed over r to make it an r-independent quantity and to study the effect as a function of R and γ. We present in Figure 7 some vibronic coupling transitions as a function of R for γ ) 45°. The intensity of transition follows the order V01 > V02 > V13 > V03 > V12 > V23, and for other orientations too, roughly the same pattern of order is followed. V. Dynamics The quantum dynamics was performed within the VCCRIOSA33 framework making use of vibronic couplings as input at two collision energies, Ecm ) 9.5 and 23 eV. The closecoupled (CC) Schro¨dinger equation corresponding to a 20 × 20 matrix of vibronic couplings was solved by the sixth-order

sin θlab m1 cos θlab + √m22 - m21 sin2 θlab m2

(

) (5)

where θ is the scattering angle and m1, m2, and m3 are the masses of the particles in the triatomic system. Opacity Function. The angle-dependent opacity is computed for the IVE and VCT processes through the relation

σl(jV f V';γ) )

()

|

π jl (2l + 1) TVV' (γ) 2 kjV

|

2

(6)

where l is the orbital momentum quantum number of the projectile, j is the rotational momentum of the target molecule, k is the wave vector, and T is the γ dependent T-matrix element. We have assumed j ) 0 throughout the calculation, indicating that the target molecule does not rotate during collisions. It is defined as the relative probability that a given process occurs as a function of the impact parameter (related to each contributing l value) of the initial vibrational state (V) and of the final

10362

J. Phys. Chem. A, Vol. 114, No. 38, 2010

Xavier

Figure 8. The opacity function as a function of partial waves (l) in the units of p at 23 eV, with the actual opacities multiplied by a factor of 5 for all angles except 15° for clarity in the VCT process in the second row.

vibrational state (V′). We present in Figure 8 the computed opacities as a function of the number of partial waves (l) for IVE and VCT processes at 23 eV. The lmax at 9.5 eV is 500 and 200, and the same at 23 eV is 800 and 300 for IVE and VCT channels, respectively. The relative probability of various process soon goes to 0 once lmax is reached. The first row shows the relative probability of elastic and inelastic transitions, which follow the order 0 f 0 > 0 f 1 > 0 f 2 > 0 f 3 and holds good for all γ. The same order is followed by the VCT channel in the case of γ ) 45°, while for γ ) 15 and 75°, the order is little different, 0 f 1 > 0 f 2 > 0 f 0 > 0 f 3. Since the angular dependence plays crucial role in the relative probability (known as steric effect) of a particular transition, we show in the second row the influence of γ on the elastic process (V ) 0 f V′ ) 0) and the first VCT (V ) 0 f V′′ ) 0) process. For the elastic process, the order determined from the peak heights follows 0 > 15 > 30 > 90 > 45 ≈ 60 > 75°, and for first VCT process, the order is 15 > 75 > 45 > 30 > 60°. Here, a factor of 5 is multiplied by the actual opacities of all γ except γ ) 15° to magnify the curves. Overall, different processes occur with sufficient probability at different numbers of contributing partial waves. Rotationally Summed State-Selective Differential Cross Section. State-to-state differential cross sections (DCSs) are computed through the relation given below

dσ (jV f V') ) dω

λ

(jV f V') ∑ (2λ + 1)-1 dσ dω

(7)

λ

where

( |

dσλ 1 (jV f V') ) 2 dω 4kjV

) ∑ (2l + 1)P (cos θ l

l

|

2

l cm)Aλ(jV

f V')

(8)

where λ is the index of the Legendre polynomials and in the present case, λ ) 0, 2, 4, ..., N*; N* e 2(N - 1) with N ) 7. Aλl are the unknown coefficients found by solving a set of simultaneous linear equations involving Legendre polynomials. The computation was carried out at collision energies Ecm ) 9.5 and 23 eV as a function of scattering angle (θcm), and we discuss the results of DCSs only at 23 eV because of the availability of experimental results. In addition to state-to-state DCSs, we add them over all vibrational states (up to 20 states) to compute total DCSs, which can be compared against experimental results. The computed DCSs showed considerable undulatory structures as a function of θcm due to various types of quantum interference and fast oscillation effects. As the experimental findings are not able to resolve such effects,5-9 the computed angular distributions were smoothed by folding them with a Gaussian distributions

Nonadiabatic Dynamics and the H+ + O2 System

dσ ¯ (θ) ) dω

[

(θ - θ) +∆θ exp ∫θ¯θ¯-∆θ 2σ2 ¯

θ

2

]

dσ (θ) dθ dω

J. Phys. Chem. A, Vol. 114, No. 38, 2010 10363

(9)

in the same way as it was done in the earlier theoretical j ) is calculations,35 where ∆θ ) 1.0°, σθ ) 0.33°, and (dσ/dω)(θ the calculated values of DCSs. The assumed uncertainty in the related experiments was ∼0.5°, and it corresponds to the value at full width at half-maximum (FWHM) in the energy loss spectra.7 Since the energy resolution in the proton/hydrogen energy loss spectra was not sufficient to resolve the individual rotational transitions in the vibrational excitation spectra, it would be appropriate to compare the smoothed state-to-state DCS data with those of experiments. Hereafter, we report only the results obtained with smoothed DCSs for both IVE and the VCT processes. We present in Figure 9 the results obtained by us and those of Gianturco et al.21 and Sidis et al.20 along with the experimental

Figure 10. The total DCS and state-to-state DCS as a function of θcm for VCT channel at Ecm ) 23 eV are compared with experimental results7 and two other earlier theoretical works.20,21

Figure 9. The total DCS and state-to-state DCS as a function of θcm for IVE channel at Ecm ) 23 eV are compared with experimental results7 and two other earlier theoretical works.20,21

results.7 All of the theoretical values are plotted on the absolute scale, and no normalization was done. In the experiments, the measurements of DCSs could be reported only on the relative scale. Therefore, we normalized the experimental data with the present theoretical data with respect to the theoretical datum value taken at θcm ) 11.72° (experimental rainbow maximum for the total DCS). The present work is relatively in good agreement with the experimental results compared to the previous works. After the rainbow maximum, the computed total DCS matches with the experimental data closely. This work predicts the angular dependency of DCS quite well unlike the case of Gianturco et al., where there is an increase after the initial decrease in the state-to-state DCS. In the previous works, DCSs corresponding to V′ ) 3 and 2 are shifted to a large extent relative to the experimental data, and such a large deviation is not observed in the present computed results. The abovementioned improvements over the previous results demonstrate the use of high-quality PESs in the dynamical computations.

10364

J. Phys. Chem. A, Vol. 114, No. 38, 2010

Xavier

TABLE 3: Measured Relative Transition Probabilities7 for Vibrational Excitation of O2 by Proton Impact at Elab ) 23.7 eV for Various Scattering Angles Compared with Our Computed Values for the Inelastic Channel, up to W′ ) 5 Levels θcm (deg) 3.09 4.12 5.15 6.18 7.21 8.24 9.28 10.31 11.34

expt(IVE)

theory(IVE)

P(0 f V′)

P(0 f V′)

V′ ) 0

1

2

3

4

0.786 0.668 0.598 0.542 0.504 0.489 0.458 0.430 0.402

0.170 0.240 0.268 0.292 0.308 0.279 0.300 0.317 0.313

0.037 0.063 0.095 0.109 0.127 0.138 0.151 0.155 0.169

0.007 0.020 0.035 0.043 0.044 0.058 0.055 0.062 0.060

0.009 0.004 0.014 0.014 0.027 0.025 0.026 0.041

5

V′ ) 0

1

2

3

4

5

0.003 0.010 0.010 0.011 0.015

0.776 0.693 0.643 0.607 0.575 0.548 0.529 0.512 0.490

0.169 0.216 0.223 0.225 0.244 0.261 0.264 0.265 0.280

0.023 0.045 0.075 0.091 0.088 0.088 0.102 0.114 0.112

0.012 0.014 0.016 0.025 0.036 0.039 0.037 0.038 0.046

0.009 0.012 0.012 0.012 0.013 0.017 0.021 0.021 0.019

0.003 0.004 0.006 0.009 0.010 0.010 0.010 0.011 0.012

TABLE 4: Computed Vibrational Transition Probabilities for the Charge-Transfer Process at Ecm ) 23 eV Compared with the Experimental Values7 for the Scattering Angle θcm ) 0-11° θcm (deg) 0.00 1.03 2.06 3.09 4.12 5.15 6.18 7.21 8.24 9.28 10.31 11.34

expt(VCT)

theory(VCT)

P(0 f V′′)

P(0 f V′′)

V′′ ) 0

1

2

3

4

5

6

V′′ ) 0

1

2

3

4

5

6

0.055 0.073 0.144 0.158 0.168 0.151 0.134 0.178 0.168 0.150 0.124 0.121

0.135 0.172 0.278 0.292 0.318 0.343 0.334 0.312 0.294 0.267 0.242 0.231

0.203 0.238 0.261 0.258 0.272 0.295 0.273 0.247 0.240 0.259 0.279 0.274

0.280 0.257 0.190 0.196 0.156 0.149 0.165 0.164 0.174 0.186 0.193 0.207

0.237 0.184 0.096 0.070 0.065 0.047 0.067 0.069 0.080 0.091 0.124 0.121

0.075 0.067 0.023 0.024 0.022 0.014 0.027 0.030 0.045 0.048 0.039 0.046

0.015 0.009 0.010

0.004 0.044 0.067 0.077 0.093 0.110 0.120 0.116 0.110 0.104 0.110 0.120

0.024 0.104 0.120 0.123 0.184 0.233 0.230 0.210 0.201 0.225 0.248 0.250

0.038 0.093 0.099 0.126 0.194 0.234 0.224 0.210 0.220 0.245 0.260 0.251

0.044 0.062 0.075 0.110 0.129 0.136 0.150 0.160 0.170 0.163 0.156 0.155

0.030 0.036 0.052 0.081 0.070 0.055 0.064 0.080 0.078 0.067 0.064 0.070

0.006 0.045 0.110 0.122 0.078 0.045 0.040 0.040 0.034 0.030 0.031 0.036

0.000 0.001 0.001 0.002 0.003 0.004 0.005 0.007 0.007 0.008 0.008 0.008

The results of this work and those of previous works20,21 together with the experimental results for the VCT channel are presented in Figure 10. The numbers in the y-ordinate denote powers of 10 in both cases. The experimental data have been normalized at θcm ) 9.027° (the experimental rainbow maximum in the total DCS). The overall layout of the total DCS and the state-to-state DCS seem to suggest that this work agrees more closely with experiment than the previous works. The ordering of DCS from V′′ ) 0 to 6 matches with the ordering shown by experiment. The presence of a hump at V′′ ) 5 is well-explained but at a slighlty extended region (around θcm ) 12°). All of the state-to-state DCSs are closely spaced due to the narrow energy gap between the corresponding transitions. The presence of the primary rainbow and the nearby data is well-accounted for with

good numerical agreement, but the presence of a secondary rainbow near θ ≈ 0° due to a second possible scattering potential is not explained successfully in this work. Overall, quantitatively, the previous works are not satisfactory, while the present work provides relatively good qualitative and quantitative matching. Transition Probability. The state-to-state transition probability for the IVE process, P0fV′(θcm), and that for VCT process, P0fV′′(θcm), have been computed as a function of scattering angle (θcm) from their respective state-to-state DCS values using the set of formulas given below

P0fV'(θcm) )

TABLE 5: VCC-RIOSA Absolute State-to-State Integral Cross Section (σ) for the IVE Channel and the VCT Channel at Ecm ) 9.5 and 23 eV integral cross section (Å2) V′(V′′)

0 1 2 3 4 5 6 7 8 9 ∑V′(V′′)

9.5 eV IVE

VCT

IVE

VCT

4381.07 49.95 20.99 11.13 6.11 3.18 1.80 1.38 1.11 1.12 4477.84

1.42 1.55 1.34 0.74 0.62 0.54 0.49 0.42

2671.78 59.05 18.52 10.96 6.31 4.32 3.13 2.24 1.73 1.45 2779.511

2.14 3.14 2.75 1.74 1.20 0.96 0.78 0.68

7.12

13.39

dσ (0 f V') | ∑ dω θ

Vmax

V'*0

P0fV''(θcm) )

23 eV

|

dσ (0 f V') dω θcm

|

cm

dσ (0 f V'') dω θcm dσ (0 f V'') | ∑ dω θ

(10)

Vmax

V''*0

cm

The denominator involves the summation over state-to-state DCSs up to 20 vibrational states of O2 and O2+ molecules, and the numerator involves the state-to-state DCS as a function of θcm for a specific vibrational state, and the division between the two gives us the state-to-state transition probability. It is computed at 9.5 and 23 eV, and we present only the data at 23 eV. Table 3 lists the data for the IVE channel for the first six transitions along with experimental data. It accounts for the

Nonadiabatic Dynamics and the H+ + O2 System

J. Phys. Chem. A, Vol. 114, No. 38, 2010 10365

Figure 11. Average vibrational energy transfer ∆Evib(θcm) (in eV) for the IVE (left) and the VCT (right) channels at Ecm ) 9.5 and 23 eV.

angular dependency for all of the transitions, and the numerical matching is found to be fair, though the computed values are on the lower side of the magnitude relative to the experimental ones. The agreement seems to be good for higher-order transitions like 0 f 4, 0 f 5, and so forth. Table 4 shows the data for the VCT channel for the first seven transitions together with experimental values. Like the IVE channel, here too the angular dependency and numerical agreement seem to be fair, with the computed values being less than the experimental ones, and in higher-order transitions, the improvement in agreement is encouraging. In sum, the low magnitudes bring out the fact that the interaction potential between the GS and the first ES is computed with an equal weightage given to all four electronic states, and this will reduce the importance of interaction between the two lowest states. Average Vibrational Energy Transfer. The average (mean) vibrational energy transfer ∆Evib(θcm) (in eV) is computed using the equations given below ∞

∆Evib(θcm) ) ∆Evib(θcm) )

∑ P0fV'(θcm)∆E(0 f V')

V'*0 ∞

∑

V''*0

(11) P0fV''(θcm)∆E(0 f V'')

This is just the sum over the product of the transition probability at a specific vibrational state and the corresponding vibrational energy of either O2(V′) or O2+(V′′), as the case may be, as a function of θcm. The computed average vibrational energy

transfer ∆Evib(θcm) for both types of processes are shown in Figure 11 for both collision energies. The results for the IVE channel are compared with the available data from the experiments.6,7 At Ecm ) 9.5 eV, the agreement between theory and experiment is fairly good, even though there is a large deviation near θcm ≈ 20°. At Ecm ) 23 eV, the experimental trend of variation with scattering angle differs in IVE and VCT entirely. In IVE, the trend is to increase with angle, and the computed curve in the present work reflects the same; numerically, it is very close to the experimental data points, unlike the curve made by Gianturco et al., which is far away from the experimental data even after multiplying the actual data by a factor of 0.66. However, the plotted data in this work are the actual data, free from any multiplication. In VCT, the angular dependence is a sharp decrease followed by gradual increase; the previous work failed to explain this trend, but the computed curve in this work, even though overestimating the experimental data, predicts the observed experimental trend nicely. Therefore, the present calculation supports non-Franck-Condon (FC) nature of the vibronic mechanism because of the ability to predict the observed trend correctly. Like in IVE, the previously computed curve was not the actual data but the data multiplied by a factor of 1.5 for magnification purpose, and the data in the present work is the actual data free from any multiplication. Integral Cross Section. Rotationally summed state-to-state integral cross sections (ICSs) for both the IVE and VCT channels are reported as a function of the vibrational levels of target diatomic molecules O2 and O2+ at Ecm ) 9.5 and 23 in Table 5. There are no previously reported theoretical and experimental results for the ICSs with which to compare. The

10366

J. Phys. Chem. A, Vol. 114, No. 38, 2010

ICS for the IVE channel decreases exponentially, whereas that of the VCT channel passes through a maximum at V′′ ) 1 followed by exponential decay, with the magnitude of the IVE process being higher than that of the VCT process. Also, comparison of the ICSs at Ecm ) 9.5 and 23 eV shows that the ICSs remain almost the same for both of these energies in the IVE channel, whereas they differ widely in the VCT channel. We can conclude that the VCT is more probable at Ecm ) 23 eV than that at Ecm ) 9.5 eV. VI. Summary and Conclusion A new set of global ab initio adiabatic and diabatic PESs for the four lowest electronic states of the H+ + O2 system have been obtained at the MRCI level of accuracy employing the Dunning’s cc-pVTZ basis set in the Jacobi coordinates. The nonadiabatic interactions have been analyzed in terms of the magnitudes of NACMEs and the coupling potential matrix elements. We determined the equilibrium geometrical properties of bound triatomic HO2 ion in the GS and compared them with the earlier theoretical results. Our results are found to be in excellent agreement with those obtained from earlier high-level ab initio calculations which focused on the description of the potential well of the GS PES, thus lending credence to our ab initio calculations. We have taken the 2 × 2 quasi-diabatic potential matrix (describing the GS and the first ES and their coupling) out of the 4 × 4 quasi-diabatic potential matrix. We have modeled the experimentally observed IVE and VCT processes as twostate processes and carried out quantum dynamics. We performed quantum dynamics at two different collision energies of Ecm ) 9.5 and 23 eV using only the coupled GS and first ES PESs along with the associated coupling potentials within the VCC-RIOSA framework to obtain the five different scattering properties, namely, the opacity function, differential cross sections, transition probability, average vibrational energy transfer, and integral cross section. We observed that there is an improvement in the computed collision attributes as compared with those of earlier reported results. The improvement is much better for the IVE channel as compared to the VCT channel, where the predictions for the former are in good agreement with the experimental results. However, some discrepancies still exist between theory and experiment for the VCT channel, which, it is hoped, will be sorted out in future work. Acknowledgment. The author thanks Prof. Sanjay Kumar for encouragement and support and acknowledges the CSIR, New Delhi, for the Senior Research Fellowship. The financial assistance by IIT Madras in procuring MOLPRO through an interdisciplinary research grant in the chemical physics program is also gratefully acknowledged. Finally, the author thanks the reviewers for patiently reading the entire manuscript and providing useful suggestions which helped in improving the quality of the manuscript. References and Notes (1) Herbst, E. Chem. Soc. ReV. 2001, 30, 168–176. (2) Widicus Weaver, S. L.; Woon, D. E.; Ruscic, B.; McCall, B. J. Astrophys. J. 2009, 697, 601. (3) Udseth, H.; Giese, C. F.; Gentry, W. R. J. Chem. Phys. 1974, 60, 3051–3056. (4) Hermann, V.; Schmidt, H.; Linder, F. J. Phys. B 1978, 11, 493– 506. (5) Krutein, J.; Linder, F. J. Chem. Phys. 1979, 71, 599–604.

Xavier (6) Gianturco, F. A.; Gierz, U.; Toennies, J. P. J. Phys. B: At. Mol. Phys. 1981, 14, 667–677. (7) Noll, M.; Toennies, J. P. J. Chem. Phys. 1986, 85, 3313–3325. (8) Niedner, G.; Noll, M.; Toennies, J. P.; Schiler, C. J. Chem. Phys. 1987, 87, 2685–2694. (9) Niedner-Schatteburg, G.; Toennies, J. P. AdV. Chem. Phys. 1992, 82, 553–646. (10) Grimbert, D.; Sidis, V.; Sizun, M. Chem. Phys. Lett. 1994, 230, 47. (11) (a) Bruna, P. J.; Peyerimhoff, S. D.; Buenker, R. J. Chem. Phys. 1975, 10, 323. (b) Gianturco, F. A.; Lamanna, U. T.; Ignazzi, D. Chem. Phys. 1980, 48, 387. (c) Dhilip Kumar, T. J.; Kumar, S. J. Chem. Phys. 2004, 121, 191. (d) Dhilip Kumar, T. J.; Saieswari, A.; Kumar, S. J. Chem. Phys. 2006, 124, 034314(1). (12) Saieswari, A.; Kumar, S. J. Chem. Phys. 2008, 128, 124306(1). (13) (a) Saieswari, A.; Kumar, S. Chem. Phys. Lett. 2007, 449, 358. (b) Saieswari, A.; Kumar, S. J. Chem. Phys. 2007, 127, 214304(1). (c) Saieswari, A.; Kumar, S. J. Chem. Phys. 2008, 128, 064301(1). (14) (a) Gianturco, F. A.; Kumar, S.; Schneider, F. Chem. Phys. 1996, 211, 33. (b) Gianturco, F. A.; Kumar, S.; Ritschel, T.; Vetter, R.; Zu¨licke, L. J. Chem. Phys. 1997, 107, 6634. (15) Staemmler, V.; Gianturco, F. A. Int. J. Quantum Chem. 1985, 28, 553–564. (16) Gianturco, F. A.; Staemmler, V. Intermolecular Forces; Pullman, Reidel: Dordrecht, The Netherlands, 1981. (17) Vazquez, G. J.; Buenker, R. J.; Peyerimhoff, S. D. Mol. Phys. 1986, 59, 291–316. (18) Grimbert, D.; Lassier-Givers, B.; Sidis, V. Chem. Phys. 1988, 124, 187–204. (19) Schneider, F.; Zu¨licke, L.; DiGiacomo, F.; Gianturco, F. A.; Paidarova´, I.; Pola´k, R. Chem. Phys. 1988, 128, 311–320. (20) Sidis, V.; Grimbert, D.; Sizun, M.; Baer, M. Chem. Phys. Lett. 1989, 163, 19–22. (21) Gianturco, F. A.; Palma, A.; Semprini, E.; Stefani, F.; Baer, M. Phys. ReV. A 1990, 42, 3926–3939. (22) Saieswari, A.; Kumar, S. J. Chem. Phys. 2008, 128, 154325(1)– 154325(10). (23) Saieswari, A.; Kumar, S. J. Chem. Sci 2007, 119, 423–431. (24) Saieswari, A.; Kumar, S. J. Chem. Sci. 2009, 121, 797–803. (25) Werner, H. J.; et al. MOLPRO, version 2002.6, a package of ab initio programs; http://www.molpro.net (2002). (26) Dunning, J. T. H. J. Chem. Phys. 1989, 90, 1007–1023. (27) (a) Werner, H.-J.; Knowles, P. J. J. Chem. Phys. 1985, 82, 5053– 5063. (b) Knowles, P. J.; Werner, H.-J. Chem. Phys. Lett. 1985, 115, 259– 267. (28) (a) Werner, H.-J.; Knowles, P. J. J. Chem. Phys. 1988, 89, 5803– 5814. (b) Knowles, P. J.; Werner, H.-J. Theor. Chim. Acta. 1992, 84, 95– 103. (c) Knowles, P. J.; Werner, H.-J. Chem. Phys. Lett. 1988, 145, 514– 522. (29) (a) Buenker, R. J.; Peyerimhoff, S. D. Theor. Chim. Acta 1974, 35, 33. (b) Buenker, R. J. Int. J. Quantum Chem. 1986, 29, 435. (c) Krebs, S.; Buenker, R. J. J. Chem. Phys. 1995, 103, 5613. (30) Hirsch, G.; Bruna, P. J.; Buenker, R. J.; Peyerimhoff, S. D. Chem. Phys. 1980, 45, 335. (31) (a) Smith, F. T. Phys. ReV. 1969, 179, 111–123. (b) Mead, A.; Truhlar, D. J. Chem. Phys. 1982, 77, 6090–6098. (c) Sidis, V. AdV. Chem. Phys. 1992, 82, 73. (d) Pacher, T.; Cederbaum, L. S.; Ko¨ppel, H. AdV. Chem. Phys. 1993, 84, 293. (e) Adhikari, S.; Billing, G. D. AdV. Chem. Phys. 2002, 124, 143. (f) Baer, M. AdV. Chem. Phys. 2002, 124, 39–142. (g) Baer, M. Phys. Rep. 2002, 358, 75–142. (h) Child, M. S. AdV. Chem. Phys. 2002, 124, 1–38. (i) Worth, G. A.; Robb, M. A. AdV. Chem. Phys. 2002, 124, 355–431. (j) Jasper, A. W.; Zhu, C.; Nangia, S.; Truhlar, D. G. Faraday Discuss 2004, 127, 1–22. (k) Ko¨ppel, H. Faraday Discuss. 2004, 127, 35–47. (l) Vertesi, T.; Bene, E.; Vibok, A.; Halasz, G. J.; Baer, M. J. Phys. Chem. A 2005, 109, 3476–3484. (32) Simah, D.; Hartke, B.; Werner, H.-J. J. Chem. Phys. 1999, 111, 4523. (33) Schinke, R.; McGuire, P. Chem. Phys. 1978, 31, 391–412. (34) Le Roy, R. J. Chemical Physics Research Report, Report CP-55R; University of Waterloo: Waterloo, 1996. (35) Baer, M.; Niedner-Schatterburg, G.; Toennies, J. P. J. Chem. Phys. 1989, 91, 4182–4196. (36) van Lenthe, J. H.; Ruttink, P. J. A. Chem. Phys. Lett. 1978, 56, 20–24. (37) Raine, G. P.; Schaefer, H. F., III. J. Chem. Phys. 1984, 80, 319– 324. (38) Robbe, J. M.; Monnerville, M.; Chambaud, G.; Rosmus, P.; Knowles, P. J. Chem. Phys. 2000, 252, 9–16. (39) Huang, X.; Lee, T. J. J. Chem. Phys. 2008, 129, 044312/1–044312/14.

JP107098V