Ab initio and kinetic calculation for the abstraction reaction of atomic O(3P) with Si2H6

12 March 2002

Chemical Physics Letters 354 (2002) 291–297 www.elsevier.com/locate/cplett

Ab initio and kinetic calculation for the abstraction reaction of atomic O ð3PÞ with Si2H6 Qingzhu Zhang, Shaokun Wang, Jianhua Zhou, Yueshu Gu

*

School of Chemistry and Chemical Engineering, Shandong University, Jinan, Shandong 250100, PR China Received 9 October 2001; in final form 3 January 2002

Abstract The hydrogen abstraction reaction of O ð3 PÞ with Si2 H6 has been studied theoretically for the first time. Two transition states of 3 A00 and 3 A0 symmetries have been located for this abstraction reaction. Geometries have been optimized at the MP2 level with 6-311G+(d) basis set. A modified G3MP2 method has been used for the final singlepoint energy calculation. Based on the ab initio data, the rate constants have been calculated over a wide temperature range of 200–3000 K using canonical variational transition state theory (CVT) with small curvature tunneling effect (SCT). The calculated CVT/SCT rate constants match well with the experimental value. Ó 2002 Elsevier Science B.V. All rights reserved.

1. Introduction Disilane is an important material in plasma Chemical Vapor Deposition (CVD) and in semiconductor device process [1]. The reaction of disilane with atomic oxygen Si2 H6 þ O ð3 PÞ ! products is a likely initial step in systems where atomic oxygen is present, e.g. plasma enhanced CVD or photochemical vapor deposition (CVD) from Si2 H6 = N2 O3 mixtures [2]. However, despite its importance, kinetics studies of this reaction have been very limited; only one experimental study is on record. In 1993, Taylor and Marshall [3] measured the rate

*

Corresponding author. Fax: +86-531-8564464. E-mail address: [email protected] (Y. Gu).

constant of this reaction, and they obtained a value is ð6:0  1:0Þ  1012 cm3 molecule1 s1 at 295 K. To our best knowledge, little theoretical attention has paid to the reaction of O ð3 PÞ with Si2 H6 . We have initiated a theoretical study of the application of ab initio electronic calculations combined with the variational transition state theory for the hydrogen abstraction reaction of atomic O ð3 PÞ with Si2 H6 . The reasons for initiating such a work are twofold. First, the reaction mechanism and kinetic nature of this reaction are essential input data for computer-modeling studies directed towards obtaining an understanding of the factors controlling CVD processes. Second, the theoretical study for this hydrogen abstraction reaction from Si2 H6 by O ð3 PÞ is particularly challenging because the approach of the O ð3 PÞ to Si2 H6 with Cs symmetry proceeds over two potential energy surfaces (PESs), 3 A0 þ 3 A00 ,

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 1 1 8 - 5

292

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297

generated by the pseudo-Jahn–Teller effect. Indeed, for Cs symmetry the irreducible representa3 0 tion is A þ 3 A00 for reactants and 0 00 3 A þ 3 A þ 1 A0 þ 1 A00 for products, and therefore, both asymptotes adiabatically correlate through the PESs 3 A0 and 3 A00 in Cs . In this Letter, we focus on the H-atom abstraction reaction of O ð3 PÞ with Si2 H6 . We have first calculated geometrical parameters, frequencies, and energies of the reactant, transition states, and products for this hydrogen abstraction reaction. In a second step, we have carried out the kinetics calculations using the canonical variational transition state theory (CVT) [4–7] with the small curvature tunneling (SCT) approximation [8].

of the reaction of O with Si2 H6 . We have made two modifications in our G3MP2 calculations: (1) the geometries are obtained at the MP2/6-311G+(d) level instead of MP2(FULL)/6-31G(d) and (2) the ZPE and vibrational frequencies are obtained at the MP2/6-311+G(d) level. All relative energies quoted and discussed in the present Letter include ZPE corrections with unscaled vibrational frequencies, except in Table 2 and Figs. 2 and 3, where we also show and discuss some values without ZPE. The initial information obtained from our ab initio calculations allowed us to calculate the variational rate constants including the tunneling effects. All the kinetical calculations have been carried out using the PO L Y R A T E 7.8 program [11].

2. Computation methods

3. Result and discussion

Ab initio calculations have been carried out using GA U S S I A N 94 programs [9]. The geometries of the reactant, transition states and products have been optimized at the MP2/6-311+G(d) level. Throughout the Letter, MP2 and QCISD(T) denote the unrestricted versions, UMP2 and UQCISD(T). The vibrational frequencies have been calculated at the same level of theory in order to determine the nature of the stationary points, the zero-point energy (ZPE), and the thermal contributions to the free energy of activation. An intrinsic reaction coordinate (IRC) calculation confirmed that the transition state connects the designated reactants and products. At the MP2/6311+G(d) level, the minimum energy path (MEP) has been obtained with a gradient step size of 0.02 amu1=2 bohr in mass-weighted Cartesian coordinates. The force constant matrices of the stationary points and selected non-stationary points near the transition state along the MEP have been also calculated. Although the geometrical parameters and the frequencies of various species can be determined satisfactorily at the MP2/6-311+G(d) level, the energies obtained at this level may not be accurate enough for the subsequent kinetic calculations. Therefore, a reliable and inexpensive theory level, G3MP2 [10], was used for the energy calculation

The optimized geometries of the reactant, transition states, and products are shown in Fig. 1. The vibrational frequencies and the ZPE of the reactant, products, and transition states are listed in Table 1. The zero-point-inclusive and zeropoint-unclusive barrier heights (DH0TS and V TS ,

Fig. 1. MP2/6-311+G(d) optimized geometries for stationary points. Distances are in angstrom, and angles are in degree. The values in parenthesis are experimental data [25]. The values in italics are calculated at MP2/6-311G(2df, p) level.

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297 Table 1 Vibrational frequencies (in cm1 ) and zero-point energies (ZPE in kcal/mol) for the reactant, products, and transition states involved in the reaction of O ð3 PÞ with Si2 H6 Species

Frequencies

ZPE

Si2 H6

2272 2272 2267 2263 2263 2257 982 982 965 965 965 887 662 662 448 391 391 149 2179 2179 2155 2155 2152 2145 940 940 929 929 909 843 625 625 434 379 379 125

31.51

Si2 H5

2278 2270 2269 2256 2249 977 976 963 910 666 630 452 435 408 145

25.57

OH

3776 3

00

5.40

TSð A Þ

2283 2275 2270 2268 2263 975 972 959 914 862 821 655 639 479 414 402 388 124 65 42a 2065i

28.69

TSð3 A0 Þ

2395 2284 2275 2271 2264 1124 976 961 951 921 859 769 668 604 475 437 412 389 70 64a 2061i

30.26

The values in italics are experimental data [13]. The lowest-frequency vibrations of the transition states are considered as internal rotations. a

Table 2 The calculated potential barriers DH0TS and V TS and the reaction enthalpy DH0 at various theory levels including the MP2 ZPE correction at 0 K (the values are in kcal/mol) Levels

DH0TS ðV TS Þ 3

00

A

DH0 3

0

A

MP2/6-31G(d) 10.60 (13.42) 12.03 (13.27) )8.62 QCISD(T)/6-31G(d) 8.04 (10.87) 9.45(10.71) )6.34 G3MP2 3.03 4.57 )12.78 DH0TS is the zero-point-inclusive barrier height evaluated at the conventional transition state, V TS is the classical (i.e., zeropoint-unclusive) barrier height at the conventional transition state, and DH0 is the enthalpy of reaction at 0 K.

respectively), and the reaction enthalpy DH0 calculated are summarized in Table 2. 3.1. Reliability of the calculation Since unrestricted Hartree–Fock (UHF) reference wavefunctions are not spin eigenfunctions for

293

open-shell species [12], we monitored the expectation values of hS 2 i in the MP2/6-311+G(d) optimization. Values of hS 2 i are always in the range 0.750–0.756 for doublets and in the range 2.000– 2.046 for triplets. After spin annihilation, the values of hS 2 i are 0.750 for doublets and 2.000 for triplets, which are the exact values for a pure doublet and for a pure triplet. Thus, spin contamination is not severe. This suggests that a single determinant reference wavefunction is suitable for this system for the level of theory used in the optimization. To clarify the general reliability of the theoretical calculations, it is useful to compare the predicated chemical properties of the present particular systems of interest with experimental data. As shown in Fig. 1, the geometric parameters of Si2 H6 and OH optimized at MP2/6-311+G(d) are in good agreement with the available experimental values. From this result, it might be inferred that the same accuracy could be expected for the calculated transition state geometries, but such an inference would be unjustified because transition states are much harder to calculate. In order to check the dependence of the ab initio results on the basis set, we performed the MP2 calculation with the flexible 6-311G(2df, p) (computationally more expensive) for the reactant, transition states and products. The optimized geometrical parameters are also shown in Fig. 1. Compared with the results at MP2/6-311+G(d), it can be seen that the extension of basis set, 6-311G(2df, p), does not cause observable change. As can be seen from Table 1, the vibrational frequencies of Si2 H6 agree well with the experimentally observed fundamentals, and the maximum absolute error is 93 cm1 . These good agreements give us confidence that the MP2/6-311+G(d) theory level is adequate to optimize the geometries and to calculate the frequencies. In order to choose a reliable theory level to calculate the energy, we have calculated the dissociation energy of the Si–H bond for the Si2 H6 . The experiments give the ground-state dissociation energy D0 , and so the calculated equilibrium dissociation energies have been corrected for the ZPE. The values at MP2/6-31G(d), QCISD(T)/631G(d), G3MP2 levels are 75.52, 77.42 and 87.67

294

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297

kcal/mol. The experimental result is 88.9 kcal/mol [13]. It can be seen that the calculated value at the G3MP2 level is in good agreement with the experimental value. Both the MP2 and the QCISD(T) levels underestimate the Si–H bond dissociation energy. Therefore, the G3MP2 theory is a good choice to calculate accurate energies for the title system. 3.2. Reaction mechanism 3.2.1. Geometry The reactant disilane has two stationary geometries: the staggered and the eclipsed. The energy of the eclipsed form is 3.75 kJ/mol higher than that of the staggered form at the G3MP2//MP2/6311+G(d) level. So in this Letter the staggered form is assumed to be the stationary geometry of disilane. As mentioned above, the hydrogen abstraction reaction of O ð3 PÞ with Si2 H6 proceeds over two PESs, of symmetries 3 A0 and 3 A00 , generated by a pseudo-Jahn–Teller effect. This kind of situation has been studied previously for O ð3 PÞ þ CH4 , where, due to the higher symmetry, it is a true Jahn–Teller effect [14]. At MP2/6-311+G(d) level, two transition states with 3 A00 and 3 A0 symmetries were located, whose geometrical structures are shown in Fig. 1. The transition states of 3 A00 and 3 0 A symmetries are denoted as TS1 and TS2 , respectively. Population analysis shows that the halffilled p-orbital of the 3 A00 symmetry is in HSiHO plane and the half-filled p-orbital of the 3 A0 symmetry is perpendicular to this plane. It can be seen from Fig. 1 that TS1 and TS2 have nearly identical geometries except for the orientation angle of Si– H–O. The O ð3 PÞ atom abstracts one of the H atoms of Si2 H6 either from the trans-position of the other H atom via TS1 on the 3 A00 surface or from cis-position via TS2 on the 3 A0 surface. The breaking Si–H bonds are elongated by 9.8%, while the forming O–H bonds are longer than the equi by 36.8% for the transition librium value of 0.97 A 3 00 3 0 states of A and A symmetries. Therefore, these two transition states are reactant-like, and the H abstraction reactions from Si2 H6 by O atom via early transition states. This rather early character in these transition states is in accordance with the

low reaction barriers and the high exothermicities of this reaction, in keeping with Hammond’s postulate [15]. To get a clearer view of the 3 A00 and 3 A0 surfaces, we have carried out electron calculations at the MP2/6-311+G(d) level, varying the angle Si– H–O in relation from its respective optimized transition state geometry. Fig. 2 shows the calculations of the PESs in Cs symmetry, where all values (angles and energies) are referred to the most stable 3 A00 structure, and the internuclear distances dðSi–HÞ and dðH–OÞ are fixed. The calculated PESs are very flat. In order to confirm further that these two transition states of 3 A00 and 3 A0 symmetries connect the designated reactants and products, we traced the reaction paths as follows. With a step size of 0.05 amu1=2 bohr, the IRC has been calculated at the MP2/6-311+G(d) level from the transition states to the reactants and the products, respectively. For these two transition states of 3 A00 and 3 A0 symmetries, the profiles of the changing of the bond lengths along the reaction path are almost overlapping. The breaking bond Si–H is almost unchanged from s ¼ 1 to 0:5 amu1=2 bohr where it is close to the value of the reactant, and it stretches after s ¼ 0:5 amu1=2 bohr. The forming H–O bond shortens rapidly from reac-

Fig. 2. Energy profile for the 3 A0 and 3 A00 potential energy curves with respect to their minima at the MP2/6-311+G(d) level. For these calculations, the other internal coordinates were fixed at their respective values in the minimum. The zero value on the x-axis corresponds to the Si–H–O angle of the 3 A00 transition state.

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297

tants and almost reaches the equilibrium bond length of OH at s ¼ 0:7 amu1=2 bohr. Other bond lengths are almost unchanged during the reaction process. Therefore, the transition states of 3 A00 and 3 0 A symmetries both connect the reactants (Si2 H6 and O) with the products (Si2 H5 and OH). The geometric changes mainly takes place in the region from s ¼ 0:5 to 0:7 amu1=2 bohr.

295

dissociation energy at these levels, we think the G3MP2 level is the most reliable level, and the MP2/6-31G(d) and QCISD(T)/6-31G(d) levels yield strong overestimates for the potential barrier. Therefore, in this Letter, we calculate the potential barrier and the reaction enthalpy at the G3MP2 level for the kinetics calculation. 3.3. Kinetics calculation

3.2.2. Frequency As shown in Table 1, at the MP2/6-311+G(d) level, both transition states have one and only one imaginary frequency; the values of the imaginary frequencies are large, which implies that the quantum tunneling effect may be significant and may play an important role in the calculation of the rate constants. Direct inspection of the transition state low-frequency mode indicates that, in both cases, the mode of the lowest frequency is a hindered internal rotation instead of a small-amplitude vibration. These modes were removed from the vibrational partition function for the transition state and the corresponding hindered rotor partition function QHR ðT Þ, calculated by the method devised by Truhlar [16], was included in the expression for the rate constant. 3.2.3. Energy Table 2 summarizes zero-point-inclusive and zero-point-unclusive barrier heights (DH0TS and V TS , respectively), and the reaction enthalpy DH0 relative to the reactants at the G3MP2 level. Results obtained at MP2 and QCISD(T) levels are also listed for comparison purposes. At MP2/6-31G(d) and QCISD(T)/6-31G(d) levels, the potential barriers without ZPE contribution are almost identical. After inclusion of the MP2 ZPE, the potential barrier of the 3 A00 surface is lower 1.43 kcal/mol at MP2/6-31G(d) and 1.41 kcal/mol at QCISD(T)/6-31G(d) than that of the 3 0 A surface. At G3MP2 level, the potential barrier of the 3 A00 surface is lower 1.54 kcal/mol than that of the 3 A0 surface. This means that the 3 A00 transition state is more stable than the 3 A0 one. As we raise the calculation level ðMP2 ! QCISDðTÞ ! G3MP2Þ, the barrier height is lowered drastically for both surfaces. Taking into account the calculated results of the Si–H bond

The reaction dynamics for the 3 A00 and 3 A0 surfaces have been calculated independently. The two surfaces show similar features and further discussion in this Letter will refer only to 3 A0 surface. The MEP is calculated at the MP2/6311+G(d) level, and the energies along the MEP are refined by the G3MP2//MP2 method. The classical potential energy VMEP and the groundstate vibrational adiabatic potential energy VaG , which are functions of the IRC s, are important dynamical parameters. The changes of VMEP and VaG on the 3 A00 potential surface with the reaction coordinate are s shown in Fig. 3. From s ¼ 1 to 1:0 amu1=2 bohr, the potential curves change very slowly. From s ¼ 1:0 amu1=2 bohr, the curves increase gradually and decrease quickly on passing the transition state. It can be seen that the VMEP and VaG curves are similar in shape, and their maximum positions are almost the same at the G3MP2//MP2 level. So the ZPE, which is the difference of VaG and VMEP , is practically constant as s varies. In order to analyze this behavior in greater

Fig. 3. The potential energy ðVMEP Þ and vibrationally adiabatic potential energy curves ðVaG Þ as functions of s at G3MP2 level for the 3 A00 potential energy surface.

296

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297

Fig. 4. Changes of the generalized normal-mode vibrational frequencies as functions of s at the MP2/6-311+G(d) for the 3 A00 potential energy surface.

detail, we have shown the variation of the generalized normal mode vibrational frequencies of the hydrogen abstraction reaction in Fig. 4. In the negative limit of s, the frequencies are associated with the reactants ðSi2 H6 þ OÞ, while in the positive limit of s, the frequencies are associated with the products ðSi2 H5 þ OHÞ. For the sake of clarity, the vibrational frequencies can be divided into three types: spectator modes, transitional modes and reactive modes. The spectator modes are those that undergo little change and sometimes remain basically unchanged in going from reactants to the transition state. The transitional modes appear along the reaction path as a consequence of the transformation from free rotation or free translations within the reactant or the product limit into real vibrational motions in the global system. Their frequencies tend to zero at the reactant and the product limit and reach their maximum in the saddle point zone. The reactive modes are those that undergo the largest change in the saddle point zone, and, therefore, they must be related to the breaking/forming bonds. For the hydrogen abstraction reaction by O ð3 PÞ on the 3 A00 potential surface, mode 6 that connects the frequency of Si–H stretching vibration of Si2 H6 with the frequency of the O–H stretching vibration of OH is the reactive mode. Modes 19 and 20 are transitional modes, and other modes are spectator modes. From s ¼ 1:0 to 1:0 amu1=2 bohr, the reactive modes drop dramatically, this behavior is similar to that found in other hydrogen abstrac-

tion reactions [17–20]. If changes in other frequencies were small, this drop could cause a considerable fall in the ZPE near the transition state, but the drop of the reactive mode is compensated partially with the transitional modes. As a result, the ZPE shows very little change with s, and the classical potential energy VMEP and the ground-state vibrational adiabatic potential energy VaG curves are similar in shape. The CVT with a SCT, which has been successfully applied for several analogous reactions [8,14,21–24], is an effective method to calculate the rate constants. In this Letter, we have used this method to study the title reaction over a wide temperature range of 200–3000 K. In order to calculate the rate constants, 30 points near the transition state region along MEP were selected – 15 points on the reactant side and 15 points on the product side. Fig. 5 shows the total rate constant, which is obtained as the sum of the calculated CVT/SCT rate constants for the 3 A00 and 3 0 A surfaces, and the experimental values against the reciprocal of the temperature for the reaction of O ð3 PÞ with Si2 H6 . The calculated CVT/SCT rate constants exhibit typical non-Arrhenius behavior. Therefore, a three-parameter formula is fitted in units of cm3 molecule1 s1 as follows: kðT Þ ¼ 1:53  1017 T 2:25 expð834=T Þ: Comparing the calculated value with the limited experimental result for the title reaction, it can be

Fig. 5. Rate constants as function of the reciprocal of the temperature ðT Þ in the temperature range of 200–3000 K for the whole reaction of O ð3 PÞ with Si2 H6 . j is the experimental value.

Q. Zhang et al. / Chemical Physics Letters 354 (2002) 291–297

seen that the CVT/SCT rate constant is in good agreement with the experimental value at 295 K. This good agreement testifies that the hydrogen abstraction mechanism is reasonable for this reaction. No other addition product such as 3 OSi2 H6 structure, or insertion products such as 3 H3 SiOSiH3 or 3 Si2 H5 OH, was located, so that at present ab initio theory supports abstraction as the main reaction mechanism in accord with the recommendation of Taylor and Marshall [3].

4. Conclusion In this Letter, we have studied the reaction of O ð3 PÞ with Si2 H6 using parameterized ab initio theory at the G3MP2 level and the canonical variational transition state theory (CVT) with the small curvature tunneling approximation. Rate constants were reported over the temperature range of 200–3000 K. Several major conclusions can be drawn from this calculation: 1. This reaction proceeds via a direct hydrogen abstraction mechanism. There is no theoretical evidence for the insertion or addition mechanism. 2. For this H atom abstraction reaction, the approach of the O ð3 PÞ to Si2 H6 with Cs symmetry proceeds over two PESs, 3 A0 and 3 A00 , generated by a pseudo-Jahn–Teller effect. 3. The calculated rate constants exhibit typical non-Arrhenius behavior. 4. The calculated CVT/SCT rate constant is in good agreement with the limited experimental data. Acknowledgements The authors thank Professor Donald G. Truhlar for providing the PO L Y R A T E 7.8 program. This work is supported by the Research Fund for the Doctoral Program of Higher Education of China.

297

References [1] C.J. Giunta, J.D. Chapple, R.G. Gordon, J. Electrochem. Soc. 137 (1990) 3237. [2] Y.K. Bhatnagar, W.I. Milne, Thin Solid Films 168 (1989) 345. [3] C.A. Taylor, P. Marshall, Chem. Phys. Lett. 205 (1993) 493. [4] K.K. Baldrige, M.S. Gordor, R. Steckler, et al., J. Phys. Chem. 93 (1989) 5107. [5] A. Gonzalez-Lafont, T.N. Truong, D.G. Truhlar, J. Chem. Phys. 95 (1991) 8875. [6] B.C. Garrett, D.G. Truhlar, J. Phys. Chem. 83 (1979) 1052. [7] B.C. Garrett, D.G. Truhlar, R.S. Grev, A.W. Magnuson, J. Phys. Chem. 84 (1980) 1730. [8] Y.-P. Liu, G.C. Lynch, T.N. Truong, D.-H. Lu, D.G. Truhlar, B.C. Garrett, J. Am. Chem. Soc. 115 (1993) 2408. [9] M.J. Frisch et al., GA U S S I A N 94, Revision E.1, Gaussian, Pittsburgh, PA, 1995. [10] L.A. Curtiss, P.C. Redfern, K. Raghavachari, V. Rassolov, J.A. Pople, J. Chem. Phys. 110 (1999) 4703. [11] R. Steckler, Y.Y. Chuang, P.L. Fast, J.C. Corchade, E.L. Coitino, W.P. Hu, G.C. Lynch, K. Nguyen, C.F. Jackells, M.Z. Gu, I. Rossi, S. Clayton, V. Melissas, B.C. Garrett, A.D. Isaacson, D.G. Truhlar, PO L Y R A T E Version 7.8, University of Minnesota, Minneapolis, 1997. [12] R. Liu, J.S. Francisco, J. Phys. Chem. A 102 (1998) 9869. [13] C.B. Musgrave, S. Dasgupta, W.A. Goddard, J. Phys. Chem. 99 (1995) 13321. [14] J.C. Corchado, J. Espinosa-Garcia, O. Roberto-Neto, Y.Y. Chuang, D.G. Truhlar, J. Phys. Chem. A 102 (1998) 4899. [15] G.S. Hammmond, J. Am. Chem. Soc. 77 (1955) 334. [16] D.G. Truhlar, J. Comput. Chem. 12 (1991) 266. [17] D.G. Truhlar, A.D. Isaacson, J. Chem. Phys. 77 (1982) 3516. [18] J. Espinosa-Garcia, J.C. Corchado, J. Phys. Chem. 100 (1996) 16561. [19] J.C. Corchado, J. Espinosa-Garcia, J. Chem. Phys. 106 (1997) 4013. [20] J. Espinosa-Garcia, J.C. Corchado, J. Phys. Chem. 101 (1997) 7336. [21] V.S. Melissas, D.G. Truhlar, J. Phys. Chem. 98 (1994) 875. [22] H.M. Yin, B.H. Yang, K.L. Han, G.Z. He, J.Z. Guo, C.P. Liu, Y.S. Gu, Phys. Chem. Chem. Phys. 2 (2000) 5093. [23] X. Yu, S.M. Li, Z.F. Xu, Z.S. Li, C.C. Sun, Chem. Phys. Lett. 320 (2000) 123. [24] Y.X. Yu, S.M. Li, Z.F. Xu, Z.S. Li, C.C. Sun, Chem. Phys. Lett. 302 (1999) 281. [25] J. Urban, P.R. Schreiner, G. Vacek, P.R. Schleyer, J.Q. Huang, J. Leszczynski, Chem. Phys. Lett. 264 (1997) 441.