- Email: [email protected]
A reduced dimensionality QM study of the BO+H2→HBO+H reaction: tunneling in polyatomic reactions
26 February 2002
Chemical Physics Letters 353 (2002) 446–454 www.elsevier.com/locate/cplett
A reduced dimensionality QM study of the BO þ H2 ! HBO þ H reaction: tunneling in polyatomic reactions H. Szichman 1, M. Gilibert, M. Albertı, X. Gimenez, A. Aguilar
*
Departament de Quımica Fısica and Centre Especial de Recerca en Quımica Te orica, Universitat de Barcelona, Martı i Franqu es 1, 08028 Barcelona, Spain Received 14 November 2001; in final form 12 December 2001
Abstract We report in the present Letter a three-dimensional, quantum mechanical, infinite order sudden approximation (IOSA) study of the H2 þ BO combustion reaction using a recently reported, six-dimensional potential energy surface (PES) for the adiabatic ground electronic state H2 BO. Total reactive probabilities, cross-sections and rate constants are thus computed. These last ones compare relatively well with the experimental data, the agreement being much more improved with regard to previous quasi-classical trajectories (QCT) predictions. Furthermore, it is found that computed reactive cross-sections by both models compare well only beyond a translational energy of 0.8 eV. Ó 2002 Elsevier Science B.V. All rights reserved.
1. Introduction Tunneling in molecular processes is an ubiquitous quantum mechanical (QM) feature. Electron tunneling is technically exploited in instrumentation, such as the Scanning Tunneling Microscope, as well as in several solid-state devices. In chemical reactions, electron as well as heavy particle tunneling spreads over all types of systems, e.g from gas-phase elementary triatomic reactions to complex biomolecular processes [1]. In spite of this, its
*
Corresponding author. Fax: +34-93-402-1231. E-mail address: [email protected] (A. Aguilar). 1 On leave from: Department of Physics and Applied Mathematics, Soreq NRC, Yavne 81800, Israel.
ultimate nature is still a subject of controversy [2], and so is its understanding in terms of either a pure QM context [3] or a semiclassical framework [4], being the latter a truly hot topic in the recent literature. One may argue that such understanding is not necessary at all, for tunneling naturally emerges from the solution of the Schr€ odinger equation, and thus it is fully taken into account, as any other type of QM feature. However, these solutions are only available through numerical computations, being thus severely hampered by the exponential growth of the computational effort with the number of degrees of freedom. Consequently, a physical insight on tunneling in multidimensional, complex systems might only be available through the proper, physically sound extension of low-dimensionality results.
0009-2614/02/$ - see front matter Ó 2002 Elsevier Science B.V. All rights reserved. PII: S 0 0 0 9 - 2 6 1 4 ( 0 2 ) 0 0 0 4 6 - 5
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
The situation becomes even more complicated if tunneling takes place concurrently to classically allowed processes, or when the transition from tunneling to classical regimes does not take place under clear-cut conditions. This is very often the case in chemical reactions, where the passage from below-the-barrier to over-the-barrier collisions may be masked by averaging over angular momentum, as well as over initial states and velocities. A common approach, to single out the tunneling contribution, consists then in comparing QM results to their classical counterpart. A limiting case is obtained if the total energy is below that corresponding to the barrier to reaction, then the full QM contribution may be entirely attributed to tunneling, although this way of proceeding has its subtleties [2]. It is especially desirable to have available results on chemical reactions where the computational study involves more than two mathematical dimensions, since multidimensional tunneling appears to be of major physical relevance, but at the same time it is perhaps the toughest aspect to deal with. In addition, this data may serve for the accurate refining of potential energy surfaces (PES), given the critical sensitivity of the amount of tunneling with the barrier height and width. The present work finds its spot in this latter context. The accurate characterization of the BO þ H2 ! HBO þ H title reaction dynamics looks especially adequate for our purpose, since it appears to be profoundly mediated by H-atom tunneling in the threshold region. This feature was first pointed out by Page [5,6], where an enhancing factor of 22 at 300 K (as well as 1.3 at 1300 K) was found for a transition-state theory (TST) rate constant calculation, after inclusion of unidimensional, Eckart barrier tunneling corrections. Going beyond the useful TST requires solving the complete dynamics on a full PES. This was done classically (i.e., without tunneling) by Sogas et al. [7], where a 6D PES was developed and reaction selectivity, as well as the role of the spectator bond [8], was studied by means of quasi-classical trajectories (QCT). The ideal situation is the use of a full-dimensional QM method, calculate cross-sections and rate constants, check with experimental informa-
447
tion and estimate the tunneling contribution at threshold by comparing to existing QCT results. The exact, comprehensive 6D QM treatment of four-atom systems is currently a most challenging problem of scattering dynamics. Full results for the benchmark H2 þ OH system are already available [9], even though important aspects concerning the PES accuracy are still open [10]. However, as it also happens with polyatomic systems, as soon as one considers heavier masses or the PES gets increasingly involved (the paradigm being the existence of minima in the strong interaction region [11]), the computational overload becomes exceedingly heavy, so that one must resort to approximate methods. In spite of their obvious shortcomings these methods critically check the relevance of the approximations on which they are based and might be fairly extended to more complex systems. Tetraatom reactive systems have been extensively explored, in the last decade, by means of QM reactive scattering approaches [12–22]. The use of negative imaginary potentials (NIPs) [23] to decouple arrangement channels (ACs) [24] has improved the achievements, among which we may note an exact 6D (J ¼ 0) calculation for the H2 þ OH reaction [16], a 6 mathematical dimensional (MD) coupled-states calculation of crosssections for the H2 O þ H reaction [20] and several 5MD coupled-states calculations that yielded state-selected cross-sections and rate constants [14,15,19]. Recently, a number of reactive tetraatom systems within the infinite order sudden approximation (IOSA) have been considered, being mainly 3D calculations with the three Jacobi angles treated as parameters randomly selected. Despite its relative simplicity, this approximation yielded results which are close to those obtained with more accurate methods [20]. This method has been applied to the O þ O3 ! 2O2 reaction [17], for which the calculated rate constants fitted the experimental data very nicely, and to the H2 þ OH $ H2 O þ H reactions [18], whose crosssections were in a reasonable agreement with 5/ 6MD coupled-states calculations [19,20]. Particularly, it should be mentioned that a recent study on the H2 þ CN reaction was treated by this method as well [25]. This study shows that the QM-IOSA
448
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
calculated cross-sections approach the exact 6MD results [26] better than any other QM [27] or QCT predictions [28]. Here, a three-dimensional QM study of the title reaction, within the IOSA is performed. The present study differs from the previous ones, on tetraatom systems, mainly in two aspects. First, the reactant AC is treated by using a polar-averaged PES [29,30] and, second, a discrete grid for the selection of the two remaining reagent Jacobi angles is used [21,30,31], rather than a random rule, so that the statistical uncertainty related to these parameters is expected to be reduced. Among the relevant new features, it will be shown here that a marked anisotropy of the interaction potential, as a function of the orientation angle, might be not as severe as for triatomic reactions, especially at the relatively low translational energies corresponding to the tunneling regime. The remainder of the Letter is organized as follows. Section 2 briefly accounts for the theoretical fundamentals, mainly stressing the original parts with respect to previous treatments. Section 3 shows the main numerical details, while results are shown in Section 4. Finally Section 5 concludes.
2. Dynamical QM method and PES As described elsewhere [13,30], an AC approach, which considers three different, independent regions (reagents, products and the region of close-interaction) is adopted. Thus, in the general case three different Schr€ odinger equations must be solved. However, a further and important simplification can be introduced by considering only the reagents AC, albeit only quantities related to the total flux may be calculated. In this case, one of the possibilities consists in solving at first the Schr€ odinger equation for the Hamiltonian associated to the reagent unperturbed region (Hr ). Its solution, wtk , leads to the total wavefunction wt , solution of total Hamiltonian (H), as a sum of two terms: the unperturbed wtk and the perturbed part, vtk . This unperturbed term follows from the inhomogenous Schr€ odinger equation in the closeinteraction region,
ðE HI Þvtk ¼ Vk wtk ;
ð1Þ
where Vk is the interaction potential and HI ¼ H þ VNIP , a phenomenological procedure to account for flux absorption in the products AC, by means of a NIP VNIP . Equivalently to previous reports on similar diatom–diatom reactions [13,30], a Jacobi coordinate system has been used in the present work, to describe the ACs of the four-atom system. Thus, the diatom–diatom (reagent) channel is described by three radial distances and three Jacobi angles (see Fig. 1 from [30]). However, the calculations reported here differ from the former by the use of the 5D polar angle-averaged PES Z 1 p U ðRr1 r2 jc1 c2 Þ ¼ U ðRr1 r2 c1 c2 /Þ d/: ð2Þ p 0 Since it is assumed that, being the reaction collinearly dominated, the average does not mask any relevant dynamical feature. The calculation of wk0 (and vk0 ) within a 3D IOSA approach, is based on the solution of the following Hamiltonian [21,29] H ¼
h2 o2 h2 o2 h2 o2 r r R 1 2 2l1 r1 or12 2l2 r2 or22 2lR oR2
þ
h2 j1 ðj1 þ 1Þ h2 j2 ðj2 þ 1Þ þ 2l1 2l2 r12 r22
þ
h2 J ðJ þ 1Þ þ U ðr1 r2 Rjc1 c2 Þ; R2 2l
ð3Þ
where li and ji (i ¼ 1; 2) are the reduced masses and the rotational quantum numbers, respectively, related to the two diatomics, J is the total angular momentum quantum number and l is the corresponding reduced mass of the whole diatom + diatom system. Note that /, the out-of-plane angle, is missing. The expression of the (unperturbed) PES to be used (so that Hr ¼ Tr þ W ) is given by W ðr1 r2 Rjc1 c2 Þ ¼ v1 ðr1 Þ þ v2 ðr2 Þ þ wðRjc1 c2 Þ;
ð4Þ
where vi ðri Þ, i ¼ 1; 2 are the asymptotic vibrational potentials of the two diatoms and w is the distortion potential in the diatom–diatom AC. Expressions for vi ðri Þ are directly obtained from the PES at R ! 1. As for w, the following expression is adopted
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
wðRjc1 c2 Þ ¼ U ðr1e r2e Rjc1 c2 Þ;
ð5Þ
where rie , i ¼ 1; 2, are the diatomic equilibrium distances. Note that these quantities are obtained from the equilibrium properties of the H2 and BO diatomic molecules. The function vk0 is derived by solving Eq. (1) in the reagents AC. For that purpose, the range of the reagents vibrational coordinate(s) are enlarged so as to comprise the relevant reactive regions and include the necessary decoupling NIPs. In the title reaction, only one AC is open [5–7], so that the BO bond remains unbroken through the full collision process. This consideration requires adding just two negative imaginary terms to the real Hamiltonian: a vibrational term along the H–H distance r1 , and a translational term along R, the distance joining the diatoms’ center of mass VI ðr1 ; r2 ; RÞ i½vIr1 ðr1 Þ þ vIR ðRÞ :
ð6Þ
The addition of the NIPs to the real averaged potential U converts the scattering problem into a bound system problem, and hence makes vk0 expandable in terms of square integrable L2 functions [32]. These functions are chosen here as localized functions for the translational components and adiabatic basis sets for the vibrational ones. Thus vJk0 ðr1 r2 Rc1 c2 jj1 j2 Þ 1 X J ¼ a gðRjnÞf ðr1 r2 Rc1 c2 jj1 j2 jntÞ; r1 r2 R nt nt
ð7Þ
where gðRjnÞ represents the translational component which is chosen to be a standard Gaussian function. Regarding f ðr1 r2 Rc1 c2 jj1 j2 jntÞ, this is an eigenfunction of a 2D Schr€ odinger equation, which is solved for fixed values of R, c1 and c2 , by defining two vibrational adiabatic basis sets which allow to build-up and diagonalise the implied Hamiltonian. The whole solution process leads to the complete set of non-reactive probabilities (inelastic S-matrix), from which reactive ones are extracted by simply substracting from unity. However, the computation of total reaction probabilities requires a prior averaging of the reaction probabilities so obtained over the c1 and c2 angles
J ðEtr ; k0 Þ ¼ Preact
1 4
Z
449
Z
1
1
d cos c1 1
1
ðEtr ; k0 jc1 c2 Þ:
J d cos c2 Preact
ð8Þ
Then, the QM total reactive cross-sections are calculated by using p rr ðEtr ; k0 Þ ¼ 2 k ðEtr Þð2j1 þ 1Þð2j2 þ 1Þ X J ð2J þ 1ÞPreact ðEtr ; k0 Þ; ð9Þ
J
where kðEtr Þ is the standard wavenumber for the whole diatom + diatom system, which is defined by k 2 ðEtr Þ ¼ ð2lÞ=h2 Etr . As it is known, once the energy dependence of the reactive cross-sections is obtained (Eq. (9)), the corresponding thermal rate constant can be easily calculated. The PES developed for the H2 þ BO ! H þ HBO reaction has already been presented in full detail elsewhere [7], so that here we will limit ourselves to briefly describe its most relevant features. The PES is expanded as an analytical function of the standard Sorbie-Murrell type [33], which includes, as it is well-known, mono-, bi-, tri-, and tetraatomic potential forms, for a tetraatom system. Among the six possible triatomic fragments, only the term corresponding to the HBO molecule was included, whereas a double four-body term was included to characterize the tetraatomic term. The other five possibilities for the triatomic fragments were excluded, since they are not energetically accessible in the present study. The resulting PES shows a collinear TS placed at 8:62 kcal mol1 above reactants, showing a good description of the reaction exothermicity [5]. It also shows the presence of two four-atom minima, one of C2v symmetry and the other showing a trans geometry. As previously discussed, no direct reaction path has been found connecting these minima with the reagents asymptote [7]. 3. Numerical details The quantum dynamical computations of nonreactive probabilities for the H2 þ BO collisional process have been carried out over the range of
450
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
translational energies 0:050 6 Etr =eV 6 0:500. The initial ro-vibrational state of both diatoms was considered to be at the lowest state, for the main purpose is the characterization of the bare tunneling reactivity. Specific calculations would show that at room temperature (T ¼ 300 K), the most populated state is the ortho-H2 species with j1 ¼ 1 ðH2 ðm1 ¼ 0; j1 ¼ 1ÞÞ, whereas for BO is j2 ¼ 7 (BO ðm2 ¼ 0; j2 ¼ 7Þ). According to our previous QCT study, no significant new features are expected from moderate rotational excitations of the reactants, so that the present conclusions might be extended to more physically representative conditions [7]. Adequately solving Eq. (1) (for a given wk0 ) required dividing the R-translational axis (1.70–7.0 ) into up to 50 equidistant sectors. Attached to A each of these sectors one translational gaussian, and a set of twofold adiabatic vibrational basis functions have been used. The number of such functions varies from one sector to another, the maximum value being constrained by a cutoff energy value of 2.2 eV [13,32]. These considerations led to solving about 600 complex equations, in order to obtain the aJnt coefficients of Eq. (7). Two linear-ramp NIPs were used in the calculations: one along the translational coordinate R ; 0:8 eVÞ and with the parameters ðDRI ; vI Þ ð1 A one along the reactive bond r1 with the parameters ; 0:1 eVÞ. Here DRI and Dr1I are ðDr1I ; v1I Þ ð1 A the R and r1I intervals along which the NIPs differs from zero, and vI and v1I are the corresponding NIP heights. The full wavefunction was calculated and along the r1 ðH2 Þ and r2 ðBOÞ ranges 0.4–2.2 A , respectively. 0.8–1.6 A Finally, as already stated, the present calculation is carried out using a fixed rectangular grid for (cos c1 ; cos c2 ). It was established that with this procedure a higher statistical accuracy is obtained [21,30,31]. For this purpose the range of cos c1 ðH2 Þ was chosen as (0:0; 1:0) and the results multiplied by two to take into account the interchangeability of the two H. The range for cos c2 ðBOÞ was found to be shorter, namely (0:625; 1:0). These angle ranges stand for the influence of orientation on the reactive probability, determining the basis for establishing the accuracy of the present IOS methodology.
4. Results 4.1. Distribution functions of reactive probabilities and cross-sections In order to learn about orientation effects on the reactivity of the title reaction, in Fig. 1 are shown the iso-contours of the angular distribution (over (cos c1 ; cos c2 )), of computed cross-sections by the present model for different values of the translational energy Etr . The cross-sections range of the iso-contours runs from a maximum value rrmax to a minimum value rrmin in linear step rrstp . From a first sight it is quite clear that the reaction probabilities are manifestly more sensible to the directions of c2 (the angle between the diatom BO and the translational distance R) than with c1 (the angle that H2 forms with R). Also it is inferred from the same figure that the title reaction has a predominantly collinear geometry, due certainly to the similar structure of the HBO molecule [5]. A further inspection shows that this collinearity is more evident for larger energies than lowers and that the effect of tunneling comes indistinctly from the directions of c1 . A similar scope on the collinear structure of the title reaction can be observed in Fig. 2 where the distribution of averaged reactive cross-sections over the incidence collision angle c1 is given. The distribution is shown as depending on the remaining collision angle c2 , being possible to observe that the departure of collinearity is less than 15°. Fig. 3 shows the calculated total averaged reactive probabilities for the title reaction as a function of the total angular momentum (opacity function). The results are presented for values of Etr ranging from 0.05 to 0.5 eV at intervals of 0.05 eV. In all cases the calculated opacity functions show a smooth behavior with a maximum at J ¼ 0 and decreasing to zero at a cutoff Jc value which increases with increasing energy (for example, as seen in the same figure, Jc ¼ 75 for Etr ¼ 0:5 eV, while Jc ¼ 45 for Etr ¼ 0:15 eV. The present QM calculated cross-sections for the H2 þ BO ! H þ HBO reaction as a function of the translational energy are shown in Fig. 4 compared with similar QCT results obtained using the same PES at identical assigned ro-vibrational
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
451
Fig. 1. Iso-contours of reactive cross-section distributions per square solid angle unit for the title reaction at some characteristic translational energies in Etr in eV. The cross-sections of the iso-contours are given in bohr2 per square solid angle units. The different levels in this figure are indicated by the sequence (see text) Etr , rrmax , rrmin , rrstp ; as follows: (a) 0.30, 1, 0.2, 0.1; (b) 0.35, 10, 1, 1; (c) 0.40, 10, 1, 1; (d) 0.45, 20, 4, 2.
Fig. 2. Distribution of averaged reactive cross-sections (in bohr2 =X) over the incidence collision angle c1 . The distribution is shown as depending on the remaining collision angle c2 . (X indicates here a solid angle unit).
Fig. 3. Opacity functions calculated as a function of the total angular momentum. Curves have been calculated for values of Etr ranging from 0.50 to 0.05 eV at steps of 0.05 eV, which can be identified at J ¼ 0 by descending order from top to bottom.
452
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
Fig. 4. A comparison of the cross-section as a function of the translational energy for the process H2 þ BO ! H þ HBO for the QM and QCT theories mentioned in this Letter. Curves have been computed for values of Etr ranging from 0.05 to 0.80 eV. (a) Full line (QM). (b) Full dots () (QCT).
initial states. As seen from this figure, there is no apparent threshold for this reaction according to the QM values, which seems to be in contradiction with the QCT [7] and ab initio data [5,6]. These data appear to indicate the presence of a classical barrier which is in agreement with the height of the potential energy barrier of the title reaction. We have indeed observed such a barrier in our calculation of the sectorial vibrational energy levels along the evolution on the translational distance, but predominantly along the collinear orientation of the H2 þ BO system. In other directions, this barrier increases, the potential becoming then more repulsive and consequently the reactivity diminishing appreciably. In such conditions the tunneling effects may become of great importance. 4.2. Thermal rate constants for the reaction H2 þ BO ! H þ HBO In Fig. 5 is shown the calculated QM thermal rate coefficient for the title reaction as a function of temperature over the range 625 6 T=K 6 1250. Also included for comparison are the available experimental estimates, the results of QCT calculations and the best results of Page which include tunneling effect from a monodimensional approximation [6]. In order to compare better from a quantitative point of view in Table 1 are listed the same rate constants k for the process H2 ðm1 ¼ 0; j1 ¼ 0Þ þ BOðm2 ¼ 0; j2 ¼ 0Þ as appear-
Fig. 5. Logarithmic plot of the rate constant as a function of the inverse of T.
ing in Fig. 5, for some selected temperatures, together with corresponding QCT and experimental values. Among the aspects to be highlighted, it is observed that: (a) the QCT rate constants are more than one order of magnitude lower than those obtained from QM calculations at low temperatures, (b) the QM/QCT rate constant ratio decreases as the temperature increases, (c) at the lowest temperatures, QM calculated rate constants are close to respective experimental ones, increasing the difference between them when the temperature goes up. The difference in reactivity between QM and QCT calculations found in the range of low
Table 1 Rate constants k in 1013 cm3 s1 for the process H2 ðm1 ¼ 0; j1 ¼ 0Þ þ BOðm2 ¼ 0; j2 ¼ 0Þ for some selected temperatures. Comparison is done of present results with corresponding QCT and experimental values T/K
Expa
QCTb
QMc
QM/QCT
689 775 857 900 1029
0.97 1.60 2.28 3.64 9.00
0.0197 0.0474 0.1040 0.9357 2.5437
0.73 0.96 1.32 1.75 2.06
37 20 13 1.8 0.8
a
Ref. [6]. Ref. [7]. c Present. b
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
energies, is due mainly to tunneling through the potential barrier. This effect is expected to be important in the present H-atom transfer reaction, mainly at low temperatures. A measure of the tunneling effect as a function of temperature can be inferred from the QM/QCT rate constant ratio, which values are reported in the last column of Table 1. As can be seen, this ratio decreases quickly as temperature increases. Thus, from a certain temperature (at the highest value), the QM rate constants (with inclusion of tunneling effect) are lower than the QCT ones. This kind of behaviour has been also observed in the theoretical study of Page [5], whereas it has been pointed out, from a conventional TST, tunneling effect has been included as a simple quantum correction for motion along the reaction coordinate. Page predicts that tunneling may enhance the rates by a factor of 22 at 300 K to a factor of 1.3 at 1000 K. Our 3D QM calculations show higher enhancements of reactivity that, at low temperatures, approaches the rate constant to the experimental value. Returning to Fig. 5, considering our QM results, we may point-out that although they fall nicely close to the experimental data at lower temperatures, we found that the dynamical properties of the inferred Arrhenius curve, are not sufficiently accurate to reproduce rate constant values at higher temperature. Moreover, the extrapolation of both experimental and 3D QM results seems indicate that at lowest temperatures than those used in the experiment, QM rate constant should be slightly higher that the corresponding experimental result. In consequence, at least according with the present results, both the shape and the height of the potential barrier should be modified. Nevertheless, the decrease of the height barrier must be lesser than the expected from QCT calculations.
5. Summary In this Letter a 3D quantum mechanical study of the title system using a six-dimensional tetratomic potential energy surface is presented. The numerical treatment was done within the non-reactive infinite order sudden approximation, where at each calculation the two frozen angles (c1 ; c2 )
453
were discretely selected through a two-dimensional grid and the third, polar angle /, was eliminated by using a /-averaged potential. The final reactive probabilities were obtained by means of a due averaging of the partial results over the quantities odinger equations were (cos c1 ; cos c2 ). The Schr€ solved in the reagents channel and the exchange channel were eliminated by using NIPs. Although considering the reagents channel, the remnant reaction processes were assigned to one only open channel. The calculated reactive cross-sections are quite an improvement with respect to QCT estimates as can be seen from the QM calculations of the rate constant values which are closer to the experimental ones.
Acknowledgements This work was done under the auspices of the ‘‘Profesores Visitantes de IBERDROLA’’ program and the Spanish DGICYT Project (Reference PB97-0919) and the Generalitat de Catalunya project 2000SGR 00016. Thanks are also due to the Centre de Supercomputaci o de Catalunya (CESCA) for the computational facilities used in carrying out the present calculations.
References [1] W.H. Miller, J. Phys. Chem. A 105 (2001) 2492. [2] N.E. Henriksen, F.Y. Hansen, G.D. Billing, Chem. Phys. Lett. 330 (2000) 139. [3] C.L. Lopreore, R.E. Wyatt, Phys. Rev. Lett. 82 (1999) 5190. [4] E.J. Heller, J. Phys. Chem. A 103 (1999) 10433. [5] M. Page, J. Phys. Chem. 93 (1989) 3639. [6] N.L. Garland, C.T. Stanton, H.H. Nelson, M. Page, J. Chem. Phys. 95 (1991) 2511. [7] J. Sogas, M. Albertı, X. Gimenez, R. Say os, A. Aguilar, J. Phys. Chem. 101 (1997) 8877. [8] J. Sogas, M. Albertı, X. Gimenez, A. Aguilar, Chem. Phys. Lett. 347 (2001) 451. [9] U. Manthe, F. Matzkies, J. Chem. Phys. 113 (2000) 5725. [10] A. Rodrıguez, E. Garcıa, J.M. Alvari~ no, A. Lagana, Chem. Phys. Lett. 345 (2001) 219. [11] E.J. Rackham, F. Huarte-Larra~ naga, D.E. Manolopoulos, Chem. Phys. Lett. 343 (2001) 356. [12] N. Balakrishnan, G.D. Billing, J. Chem. Phys. 104 (1996) 4005, and references therein.
454
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
[13] H. Szichman, M. Baer, J. Chem. Phys. 101 (1994) 2081. [14] D.H. Zhang, J.Z.H. Zhang, J. Chem. Phys. 103 (1995) 6512. [15] U. Manthe, T. Seideman, W.H. Miller, J. Chem. Phys. 101 (1994) 4759. [16] D. Neuhauser, J. Chem. Phys. 100 (1994) 9272. [17] H. Szichman, A.J.C. Varandas, M. Baer, J. Chem. Phys. 102 (1995) 3474. [18] H. Szichman, M. Baer, Chem. Phys. Lett. 242 (1995) 285. [19] D.H. Zhang, J.Z.H. Zhang, J. Chem. Phys. 100 (1994) 2697. [20] D.H. Zhang, J.C. Light, J. Chem. Phys. 104 (1996) 4544. [21] H. Szichman, M. Baer, H.-R. Volpp, J. Wolfrum, J. Phys. Chem. A. 102 (1998) 10455. [22] H. Szichman, M. Baer, H.-R. Volpp, J. Wolfrum, J. Chem. Phys. 111 (1999) 567. [23] M.S. Child, Mol. Phys. 72 (1991) 89, and references therein. [24] A. Baram, I. Last, M. Baer, Chem. Phys. Lett. 212 (1993) 649, and references therein.
[25] D. Troya, H. Szichman, A. Aguilar, M. Gonzalez, in preparation. [26] W. Zhu, J.H.Z. Zhang, Y.C. Zhang, L.X. Zhan, S.L. Zhang, D.H. Zhang, J. Chem. Phys. 108 (1998) 3509, and references therein. [27] T. Takayanagi, G.C. Schatz, J. Chem. Phys. 106 (1997) 3227. [28] M.A. ter Horst, G.C. Schatz, L.B.J. Harding, Chem. Phys. 105 (1996) 558. [29] H. Szichman, M. Baer, J. Chem. Phys. 105 (1996) 10380. [30] H. Szichman, M. Baer, H. Nakamura, J. Chem. Phys. 107 (1997) 3521. [31] H. Szichman, M. Gilibert, M. Gonzalez, X. Gimenez, A. Aguilar, J. Chem. Phys. 113 (2000) 176. [32] I. Last, A. Baram, H. Szichman, M. Baer, J. Phys. Chem. 97 (1993) 7040. [33] J.N. Murrell, S. Carter, S.C. Farantos, P. Huxley, A.J.C. Varandas, Molecular Potential Energy Functions, John Wiley & Sons, New York, 1984.
Chemical Physics Letters 353 (2002) 446–454 www.elsevier.com/locate/cplett
A reduced dimensionality QM study of the BO þ H2 ! HBO þ H reaction: tunneling in polyatomic reactions H. Szichman 1, M. Gilibert, M. Albertı, X. Gimenez, A. Aguilar
*
Departament de Quımica Fısica and Centre Especial de Recerca en Quımica Te orica, Universitat de Barcelona, Martı i Franqu es 1, 08028 Barcelona, Spain Received 14 November 2001; in final form 12 December 2001
Abstract We report in the present Letter a three-dimensional, quantum mechanical, infinite order sudden approximation (IOSA) study of the H2 þ BO combustion reaction using a recently reported, six-dimensional potential energy surface (PES) for the adiabatic ground electronic state H2 BO. Total reactive probabilities, cross-sections and rate constants are thus computed. These last ones compare relatively well with the experimental data, the agreement being much more improved with regard to previous quasi-classical trajectories (QCT) predictions. Furthermore, it is found that computed reactive cross-sections by both models compare well only beyond a translational energy of 0.8 eV. Ó 2002 Elsevier Science B.V. All rights reserved.
1. Introduction Tunneling in molecular processes is an ubiquitous quantum mechanical (QM) feature. Electron tunneling is technically exploited in instrumentation, such as the Scanning Tunneling Microscope, as well as in several solid-state devices. In chemical reactions, electron as well as heavy particle tunneling spreads over all types of systems, e.g from gas-phase elementary triatomic reactions to complex biomolecular processes [1]. In spite of this, its
*
Corresponding author. Fax: +34-93-402-1231. E-mail address: [email protected] (A. Aguilar). 1 On leave from: Department of Physics and Applied Mathematics, Soreq NRC, Yavne 81800, Israel.
ultimate nature is still a subject of controversy [2], and so is its understanding in terms of either a pure QM context [3] or a semiclassical framework [4], being the latter a truly hot topic in the recent literature. One may argue that such understanding is not necessary at all, for tunneling naturally emerges from the solution of the Schr€ odinger equation, and thus it is fully taken into account, as any other type of QM feature. However, these solutions are only available through numerical computations, being thus severely hampered by the exponential growth of the computational effort with the number of degrees of freedom. Consequently, a physical insight on tunneling in multidimensional, complex systems might only be available through the proper, physically sound extension of low-dimensionality results.
0009-2614/02/$ - see front matter Ó 2002 Elsevier Science B.V. All rights reserved. PII: S 0 0 0 9 - 2 6 1 4 ( 0 2 ) 0 0 0 4 6 - 5
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
The situation becomes even more complicated if tunneling takes place concurrently to classically allowed processes, or when the transition from tunneling to classical regimes does not take place under clear-cut conditions. This is very often the case in chemical reactions, where the passage from below-the-barrier to over-the-barrier collisions may be masked by averaging over angular momentum, as well as over initial states and velocities. A common approach, to single out the tunneling contribution, consists then in comparing QM results to their classical counterpart. A limiting case is obtained if the total energy is below that corresponding to the barrier to reaction, then the full QM contribution may be entirely attributed to tunneling, although this way of proceeding has its subtleties [2]. It is especially desirable to have available results on chemical reactions where the computational study involves more than two mathematical dimensions, since multidimensional tunneling appears to be of major physical relevance, but at the same time it is perhaps the toughest aspect to deal with. In addition, this data may serve for the accurate refining of potential energy surfaces (PES), given the critical sensitivity of the amount of tunneling with the barrier height and width. The present work finds its spot in this latter context. The accurate characterization of the BO þ H2 ! HBO þ H title reaction dynamics looks especially adequate for our purpose, since it appears to be profoundly mediated by H-atom tunneling in the threshold region. This feature was first pointed out by Page [5,6], where an enhancing factor of 22 at 300 K (as well as 1.3 at 1300 K) was found for a transition-state theory (TST) rate constant calculation, after inclusion of unidimensional, Eckart barrier tunneling corrections. Going beyond the useful TST requires solving the complete dynamics on a full PES. This was done classically (i.e., without tunneling) by Sogas et al. [7], where a 6D PES was developed and reaction selectivity, as well as the role of the spectator bond [8], was studied by means of quasi-classical trajectories (QCT). The ideal situation is the use of a full-dimensional QM method, calculate cross-sections and rate constants, check with experimental informa-
447
tion and estimate the tunneling contribution at threshold by comparing to existing QCT results. The exact, comprehensive 6D QM treatment of four-atom systems is currently a most challenging problem of scattering dynamics. Full results for the benchmark H2 þ OH system are already available [9], even though important aspects concerning the PES accuracy are still open [10]. However, as it also happens with polyatomic systems, as soon as one considers heavier masses or the PES gets increasingly involved (the paradigm being the existence of minima in the strong interaction region [11]), the computational overload becomes exceedingly heavy, so that one must resort to approximate methods. In spite of their obvious shortcomings these methods critically check the relevance of the approximations on which they are based and might be fairly extended to more complex systems. Tetraatom reactive systems have been extensively explored, in the last decade, by means of QM reactive scattering approaches [12–22]. The use of negative imaginary potentials (NIPs) [23] to decouple arrangement channels (ACs) [24] has improved the achievements, among which we may note an exact 6D (J ¼ 0) calculation for the H2 þ OH reaction [16], a 6 mathematical dimensional (MD) coupled-states calculation of crosssections for the H2 O þ H reaction [20] and several 5MD coupled-states calculations that yielded state-selected cross-sections and rate constants [14,15,19]. Recently, a number of reactive tetraatom systems within the infinite order sudden approximation (IOSA) have been considered, being mainly 3D calculations with the three Jacobi angles treated as parameters randomly selected. Despite its relative simplicity, this approximation yielded results which are close to those obtained with more accurate methods [20]. This method has been applied to the O þ O3 ! 2O2 reaction [17], for which the calculated rate constants fitted the experimental data very nicely, and to the H2 þ OH $ H2 O þ H reactions [18], whose crosssections were in a reasonable agreement with 5/ 6MD coupled-states calculations [19,20]. Particularly, it should be mentioned that a recent study on the H2 þ CN reaction was treated by this method as well [25]. This study shows that the QM-IOSA
448
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
calculated cross-sections approach the exact 6MD results [26] better than any other QM [27] or QCT predictions [28]. Here, a three-dimensional QM study of the title reaction, within the IOSA is performed. The present study differs from the previous ones, on tetraatom systems, mainly in two aspects. First, the reactant AC is treated by using a polar-averaged PES [29,30] and, second, a discrete grid for the selection of the two remaining reagent Jacobi angles is used [21,30,31], rather than a random rule, so that the statistical uncertainty related to these parameters is expected to be reduced. Among the relevant new features, it will be shown here that a marked anisotropy of the interaction potential, as a function of the orientation angle, might be not as severe as for triatomic reactions, especially at the relatively low translational energies corresponding to the tunneling regime. The remainder of the Letter is organized as follows. Section 2 briefly accounts for the theoretical fundamentals, mainly stressing the original parts with respect to previous treatments. Section 3 shows the main numerical details, while results are shown in Section 4. Finally Section 5 concludes.
2. Dynamical QM method and PES As described elsewhere [13,30], an AC approach, which considers three different, independent regions (reagents, products and the region of close-interaction) is adopted. Thus, in the general case three different Schr€ odinger equations must be solved. However, a further and important simplification can be introduced by considering only the reagents AC, albeit only quantities related to the total flux may be calculated. In this case, one of the possibilities consists in solving at first the Schr€ odinger equation for the Hamiltonian associated to the reagent unperturbed region (Hr ). Its solution, wtk , leads to the total wavefunction wt , solution of total Hamiltonian (H), as a sum of two terms: the unperturbed wtk and the perturbed part, vtk . This unperturbed term follows from the inhomogenous Schr€ odinger equation in the closeinteraction region,
ðE HI Þvtk ¼ Vk wtk ;
ð1Þ
where Vk is the interaction potential and HI ¼ H þ VNIP , a phenomenological procedure to account for flux absorption in the products AC, by means of a NIP VNIP . Equivalently to previous reports on similar diatom–diatom reactions [13,30], a Jacobi coordinate system has been used in the present work, to describe the ACs of the four-atom system. Thus, the diatom–diatom (reagent) channel is described by three radial distances and three Jacobi angles (see Fig. 1 from [30]). However, the calculations reported here differ from the former by the use of the 5D polar angle-averaged PES Z 1 p U ðRr1 r2 jc1 c2 Þ ¼ U ðRr1 r2 c1 c2 /Þ d/: ð2Þ p 0 Since it is assumed that, being the reaction collinearly dominated, the average does not mask any relevant dynamical feature. The calculation of wk0 (and vk0 ) within a 3D IOSA approach, is based on the solution of the following Hamiltonian [21,29] H ¼
h2 o2 h2 o2 h2 o2 r r R 1 2 2l1 r1 or12 2l2 r2 or22 2lR oR2
þ
h2 j1 ðj1 þ 1Þ h2 j2 ðj2 þ 1Þ þ 2l1 2l2 r12 r22
þ
h2 J ðJ þ 1Þ þ U ðr1 r2 Rjc1 c2 Þ; R2 2l
ð3Þ
where li and ji (i ¼ 1; 2) are the reduced masses and the rotational quantum numbers, respectively, related to the two diatomics, J is the total angular momentum quantum number and l is the corresponding reduced mass of the whole diatom + diatom system. Note that /, the out-of-plane angle, is missing. The expression of the (unperturbed) PES to be used (so that Hr ¼ Tr þ W ) is given by W ðr1 r2 Rjc1 c2 Þ ¼ v1 ðr1 Þ þ v2 ðr2 Þ þ wðRjc1 c2 Þ;
ð4Þ
where vi ðri Þ, i ¼ 1; 2 are the asymptotic vibrational potentials of the two diatoms and w is the distortion potential in the diatom–diatom AC. Expressions for vi ðri Þ are directly obtained from the PES at R ! 1. As for w, the following expression is adopted
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
wðRjc1 c2 Þ ¼ U ðr1e r2e Rjc1 c2 Þ;
ð5Þ
where rie , i ¼ 1; 2, are the diatomic equilibrium distances. Note that these quantities are obtained from the equilibrium properties of the H2 and BO diatomic molecules. The function vk0 is derived by solving Eq. (1) in the reagents AC. For that purpose, the range of the reagents vibrational coordinate(s) are enlarged so as to comprise the relevant reactive regions and include the necessary decoupling NIPs. In the title reaction, only one AC is open [5–7], so that the BO bond remains unbroken through the full collision process. This consideration requires adding just two negative imaginary terms to the real Hamiltonian: a vibrational term along the H–H distance r1 , and a translational term along R, the distance joining the diatoms’ center of mass VI ðr1 ; r2 ; RÞ i½vIr1 ðr1 Þ þ vIR ðRÞ :
ð6Þ
The addition of the NIPs to the real averaged potential U converts the scattering problem into a bound system problem, and hence makes vk0 expandable in terms of square integrable L2 functions [32]. These functions are chosen here as localized functions for the translational components and adiabatic basis sets for the vibrational ones. Thus vJk0 ðr1 r2 Rc1 c2 jj1 j2 Þ 1 X J ¼ a gðRjnÞf ðr1 r2 Rc1 c2 jj1 j2 jntÞ; r1 r2 R nt nt
ð7Þ
where gðRjnÞ represents the translational component which is chosen to be a standard Gaussian function. Regarding f ðr1 r2 Rc1 c2 jj1 j2 jntÞ, this is an eigenfunction of a 2D Schr€ odinger equation, which is solved for fixed values of R, c1 and c2 , by defining two vibrational adiabatic basis sets which allow to build-up and diagonalise the implied Hamiltonian. The whole solution process leads to the complete set of non-reactive probabilities (inelastic S-matrix), from which reactive ones are extracted by simply substracting from unity. However, the computation of total reaction probabilities requires a prior averaging of the reaction probabilities so obtained over the c1 and c2 angles
J ðEtr ; k0 Þ ¼ Preact
1 4
Z
449
Z
1
1
d cos c1 1
1
ðEtr ; k0 jc1 c2 Þ:
J d cos c2 Preact
ð8Þ
Then, the QM total reactive cross-sections are calculated by using p rr ðEtr ; k0 Þ ¼ 2 k ðEtr Þð2j1 þ 1Þð2j2 þ 1Þ X J ð2J þ 1ÞPreact ðEtr ; k0 Þ; ð9Þ
J
where kðEtr Þ is the standard wavenumber for the whole diatom + diatom system, which is defined by k 2 ðEtr Þ ¼ ð2lÞ=h2 Etr . As it is known, once the energy dependence of the reactive cross-sections is obtained (Eq. (9)), the corresponding thermal rate constant can be easily calculated. The PES developed for the H2 þ BO ! H þ HBO reaction has already been presented in full detail elsewhere [7], so that here we will limit ourselves to briefly describe its most relevant features. The PES is expanded as an analytical function of the standard Sorbie-Murrell type [33], which includes, as it is well-known, mono-, bi-, tri-, and tetraatomic potential forms, for a tetraatom system. Among the six possible triatomic fragments, only the term corresponding to the HBO molecule was included, whereas a double four-body term was included to characterize the tetraatomic term. The other five possibilities for the triatomic fragments were excluded, since they are not energetically accessible in the present study. The resulting PES shows a collinear TS placed at 8:62 kcal mol1 above reactants, showing a good description of the reaction exothermicity [5]. It also shows the presence of two four-atom minima, one of C2v symmetry and the other showing a trans geometry. As previously discussed, no direct reaction path has been found connecting these minima with the reagents asymptote [7]. 3. Numerical details The quantum dynamical computations of nonreactive probabilities for the H2 þ BO collisional process have been carried out over the range of
450
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
translational energies 0:050 6 Etr =eV 6 0:500. The initial ro-vibrational state of both diatoms was considered to be at the lowest state, for the main purpose is the characterization of the bare tunneling reactivity. Specific calculations would show that at room temperature (T ¼ 300 K), the most populated state is the ortho-H2 species with j1 ¼ 1 ðH2 ðm1 ¼ 0; j1 ¼ 1ÞÞ, whereas for BO is j2 ¼ 7 (BO ðm2 ¼ 0; j2 ¼ 7Þ). According to our previous QCT study, no significant new features are expected from moderate rotational excitations of the reactants, so that the present conclusions might be extended to more physically representative conditions [7]. Adequately solving Eq. (1) (for a given wk0 ) required dividing the R-translational axis (1.70–7.0 ) into up to 50 equidistant sectors. Attached to A each of these sectors one translational gaussian, and a set of twofold adiabatic vibrational basis functions have been used. The number of such functions varies from one sector to another, the maximum value being constrained by a cutoff energy value of 2.2 eV [13,32]. These considerations led to solving about 600 complex equations, in order to obtain the aJnt coefficients of Eq. (7). Two linear-ramp NIPs were used in the calculations: one along the translational coordinate R ; 0:8 eVÞ and with the parameters ðDRI ; vI Þ ð1 A one along the reactive bond r1 with the parameters ; 0:1 eVÞ. Here DRI and Dr1I are ðDr1I ; v1I Þ ð1 A the R and r1I intervals along which the NIPs differs from zero, and vI and v1I are the corresponding NIP heights. The full wavefunction was calculated and along the r1 ðH2 Þ and r2 ðBOÞ ranges 0.4–2.2 A , respectively. 0.8–1.6 A Finally, as already stated, the present calculation is carried out using a fixed rectangular grid for (cos c1 ; cos c2 ). It was established that with this procedure a higher statistical accuracy is obtained [21,30,31]. For this purpose the range of cos c1 ðH2 Þ was chosen as (0:0; 1:0) and the results multiplied by two to take into account the interchangeability of the two H. The range for cos c2 ðBOÞ was found to be shorter, namely (0:625; 1:0). These angle ranges stand for the influence of orientation on the reactive probability, determining the basis for establishing the accuracy of the present IOS methodology.
4. Results 4.1. Distribution functions of reactive probabilities and cross-sections In order to learn about orientation effects on the reactivity of the title reaction, in Fig. 1 are shown the iso-contours of the angular distribution (over (cos c1 ; cos c2 )), of computed cross-sections by the present model for different values of the translational energy Etr . The cross-sections range of the iso-contours runs from a maximum value rrmax to a minimum value rrmin in linear step rrstp . From a first sight it is quite clear that the reaction probabilities are manifestly more sensible to the directions of c2 (the angle between the diatom BO and the translational distance R) than with c1 (the angle that H2 forms with R). Also it is inferred from the same figure that the title reaction has a predominantly collinear geometry, due certainly to the similar structure of the HBO molecule [5]. A further inspection shows that this collinearity is more evident for larger energies than lowers and that the effect of tunneling comes indistinctly from the directions of c1 . A similar scope on the collinear structure of the title reaction can be observed in Fig. 2 where the distribution of averaged reactive cross-sections over the incidence collision angle c1 is given. The distribution is shown as depending on the remaining collision angle c2 , being possible to observe that the departure of collinearity is less than 15°. Fig. 3 shows the calculated total averaged reactive probabilities for the title reaction as a function of the total angular momentum (opacity function). The results are presented for values of Etr ranging from 0.05 to 0.5 eV at intervals of 0.05 eV. In all cases the calculated opacity functions show a smooth behavior with a maximum at J ¼ 0 and decreasing to zero at a cutoff Jc value which increases with increasing energy (for example, as seen in the same figure, Jc ¼ 75 for Etr ¼ 0:5 eV, while Jc ¼ 45 for Etr ¼ 0:15 eV. The present QM calculated cross-sections for the H2 þ BO ! H þ HBO reaction as a function of the translational energy are shown in Fig. 4 compared with similar QCT results obtained using the same PES at identical assigned ro-vibrational
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
451
Fig. 1. Iso-contours of reactive cross-section distributions per square solid angle unit for the title reaction at some characteristic translational energies in Etr in eV. The cross-sections of the iso-contours are given in bohr2 per square solid angle units. The different levels in this figure are indicated by the sequence (see text) Etr , rrmax , rrmin , rrstp ; as follows: (a) 0.30, 1, 0.2, 0.1; (b) 0.35, 10, 1, 1; (c) 0.40, 10, 1, 1; (d) 0.45, 20, 4, 2.
Fig. 2. Distribution of averaged reactive cross-sections (in bohr2 =X) over the incidence collision angle c1 . The distribution is shown as depending on the remaining collision angle c2 . (X indicates here a solid angle unit).
Fig. 3. Opacity functions calculated as a function of the total angular momentum. Curves have been calculated for values of Etr ranging from 0.50 to 0.05 eV at steps of 0.05 eV, which can be identified at J ¼ 0 by descending order from top to bottom.
452
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
Fig. 4. A comparison of the cross-section as a function of the translational energy for the process H2 þ BO ! H þ HBO for the QM and QCT theories mentioned in this Letter. Curves have been computed for values of Etr ranging from 0.05 to 0.80 eV. (a) Full line (QM). (b) Full dots () (QCT).
initial states. As seen from this figure, there is no apparent threshold for this reaction according to the QM values, which seems to be in contradiction with the QCT [7] and ab initio data [5,6]. These data appear to indicate the presence of a classical barrier which is in agreement with the height of the potential energy barrier of the title reaction. We have indeed observed such a barrier in our calculation of the sectorial vibrational energy levels along the evolution on the translational distance, but predominantly along the collinear orientation of the H2 þ BO system. In other directions, this barrier increases, the potential becoming then more repulsive and consequently the reactivity diminishing appreciably. In such conditions the tunneling effects may become of great importance. 4.2. Thermal rate constants for the reaction H2 þ BO ! H þ HBO In Fig. 5 is shown the calculated QM thermal rate coefficient for the title reaction as a function of temperature over the range 625 6 T=K 6 1250. Also included for comparison are the available experimental estimates, the results of QCT calculations and the best results of Page which include tunneling effect from a monodimensional approximation [6]. In order to compare better from a quantitative point of view in Table 1 are listed the same rate constants k for the process H2 ðm1 ¼ 0; j1 ¼ 0Þ þ BOðm2 ¼ 0; j2 ¼ 0Þ as appear-
Fig. 5. Logarithmic plot of the rate constant as a function of the inverse of T.
ing in Fig. 5, for some selected temperatures, together with corresponding QCT and experimental values. Among the aspects to be highlighted, it is observed that: (a) the QCT rate constants are more than one order of magnitude lower than those obtained from QM calculations at low temperatures, (b) the QM/QCT rate constant ratio decreases as the temperature increases, (c) at the lowest temperatures, QM calculated rate constants are close to respective experimental ones, increasing the difference between them when the temperature goes up. The difference in reactivity between QM and QCT calculations found in the range of low
Table 1 Rate constants k in 1013 cm3 s1 for the process H2 ðm1 ¼ 0; j1 ¼ 0Þ þ BOðm2 ¼ 0; j2 ¼ 0Þ for some selected temperatures. Comparison is done of present results with corresponding QCT and experimental values T/K
Expa
QCTb
QMc
QM/QCT
689 775 857 900 1029
0.97 1.60 2.28 3.64 9.00
0.0197 0.0474 0.1040 0.9357 2.5437
0.73 0.96 1.32 1.75 2.06
37 20 13 1.8 0.8
a
Ref. [6]. Ref. [7]. c Present. b
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
energies, is due mainly to tunneling through the potential barrier. This effect is expected to be important in the present H-atom transfer reaction, mainly at low temperatures. A measure of the tunneling effect as a function of temperature can be inferred from the QM/QCT rate constant ratio, which values are reported in the last column of Table 1. As can be seen, this ratio decreases quickly as temperature increases. Thus, from a certain temperature (at the highest value), the QM rate constants (with inclusion of tunneling effect) are lower than the QCT ones. This kind of behaviour has been also observed in the theoretical study of Page [5], whereas it has been pointed out, from a conventional TST, tunneling effect has been included as a simple quantum correction for motion along the reaction coordinate. Page predicts that tunneling may enhance the rates by a factor of 22 at 300 K to a factor of 1.3 at 1000 K. Our 3D QM calculations show higher enhancements of reactivity that, at low temperatures, approaches the rate constant to the experimental value. Returning to Fig. 5, considering our QM results, we may point-out that although they fall nicely close to the experimental data at lower temperatures, we found that the dynamical properties of the inferred Arrhenius curve, are not sufficiently accurate to reproduce rate constant values at higher temperature. Moreover, the extrapolation of both experimental and 3D QM results seems indicate that at lowest temperatures than those used in the experiment, QM rate constant should be slightly higher that the corresponding experimental result. In consequence, at least according with the present results, both the shape and the height of the potential barrier should be modified. Nevertheless, the decrease of the height barrier must be lesser than the expected from QCT calculations.
5. Summary In this Letter a 3D quantum mechanical study of the title system using a six-dimensional tetratomic potential energy surface is presented. The numerical treatment was done within the non-reactive infinite order sudden approximation, where at each calculation the two frozen angles (c1 ; c2 )
453
were discretely selected through a two-dimensional grid and the third, polar angle /, was eliminated by using a /-averaged potential. The final reactive probabilities were obtained by means of a due averaging of the partial results over the quantities odinger equations were (cos c1 ; cos c2 ). The Schr€ solved in the reagents channel and the exchange channel were eliminated by using NIPs. Although considering the reagents channel, the remnant reaction processes were assigned to one only open channel. The calculated reactive cross-sections are quite an improvement with respect to QCT estimates as can be seen from the QM calculations of the rate constant values which are closer to the experimental ones.
Acknowledgements This work was done under the auspices of the ‘‘Profesores Visitantes de IBERDROLA’’ program and the Spanish DGICYT Project (Reference PB97-0919) and the Generalitat de Catalunya project 2000SGR 00016. Thanks are also due to the Centre de Supercomputaci o de Catalunya (CESCA) for the computational facilities used in carrying out the present calculations.
References [1] W.H. Miller, J. Phys. Chem. A 105 (2001) 2492. [2] N.E. Henriksen, F.Y. Hansen, G.D. Billing, Chem. Phys. Lett. 330 (2000) 139. [3] C.L. Lopreore, R.E. Wyatt, Phys. Rev. Lett. 82 (1999) 5190. [4] E.J. Heller, J. Phys. Chem. A 103 (1999) 10433. [5] M. Page, J. Phys. Chem. 93 (1989) 3639. [6] N.L. Garland, C.T. Stanton, H.H. Nelson, M. Page, J. Chem. Phys. 95 (1991) 2511. [7] J. Sogas, M. Albertı, X. Gimenez, R. Say os, A. Aguilar, J. Phys. Chem. 101 (1997) 8877. [8] J. Sogas, M. Albertı, X. Gimenez, A. Aguilar, Chem. Phys. Lett. 347 (2001) 451. [9] U. Manthe, F. Matzkies, J. Chem. Phys. 113 (2000) 5725. [10] A. Rodrıguez, E. Garcıa, J.M. Alvari~ no, A. Lagana, Chem. Phys. Lett. 345 (2001) 219. [11] E.J. Rackham, F. Huarte-Larra~ naga, D.E. Manolopoulos, Chem. Phys. Lett. 343 (2001) 356. [12] N. Balakrishnan, G.D. Billing, J. Chem. Phys. 104 (1996) 4005, and references therein.
454
H. Szichman et al. / Chemical Physics Letters 353 (2002) 446–454
[13] H. Szichman, M. Baer, J. Chem. Phys. 101 (1994) 2081. [14] D.H. Zhang, J.Z.H. Zhang, J. Chem. Phys. 103 (1995) 6512. [15] U. Manthe, T. Seideman, W.H. Miller, J. Chem. Phys. 101 (1994) 4759. [16] D. Neuhauser, J. Chem. Phys. 100 (1994) 9272. [17] H. Szichman, A.J.C. Varandas, M. Baer, J. Chem. Phys. 102 (1995) 3474. [18] H. Szichman, M. Baer, Chem. Phys. Lett. 242 (1995) 285. [19] D.H. Zhang, J.Z.H. Zhang, J. Chem. Phys. 100 (1994) 2697. [20] D.H. Zhang, J.C. Light, J. Chem. Phys. 104 (1996) 4544. [21] H. Szichman, M. Baer, H.-R. Volpp, J. Wolfrum, J. Phys. Chem. A. 102 (1998) 10455. [22] H. Szichman, M. Baer, H.-R. Volpp, J. Wolfrum, J. Chem. Phys. 111 (1999) 567. [23] M.S. Child, Mol. Phys. 72 (1991) 89, and references therein. [24] A. Baram, I. Last, M. Baer, Chem. Phys. Lett. 212 (1993) 649, and references therein.
[25] D. Troya, H. Szichman, A. Aguilar, M. Gonzalez, in preparation. [26] W. Zhu, J.H.Z. Zhang, Y.C. Zhang, L.X. Zhan, S.L. Zhang, D.H. Zhang, J. Chem. Phys. 108 (1998) 3509, and references therein. [27] T. Takayanagi, G.C. Schatz, J. Chem. Phys. 106 (1997) 3227. [28] M.A. ter Horst, G.C. Schatz, L.B.J. Harding, Chem. Phys. 105 (1996) 558. [29] H. Szichman, M. Baer, J. Chem. Phys. 105 (1996) 10380. [30] H. Szichman, M. Baer, H. Nakamura, J. Chem. Phys. 107 (1997) 3521. [31] H. Szichman, M. Gilibert, M. Gonzalez, X. Gimenez, A. Aguilar, J. Chem. Phys. 113 (2000) 176. [32] I. Last, A. Baram, H. Szichman, M. Baer, J. Phys. Chem. 97 (1993) 7040. [33] J.N. Murrell, S. Carter, S.C. Farantos, P. Huxley, A.J.C. Varandas, Molecular Potential Energy Functions, John Wiley & Sons, New York, 1984.