Chemical Physics Chemical Physics 182 ( 1994) 39-5 1
Ab initio studies on the mixed heterodimers of ammonia and hydrogen cyanide Soumitra Chattopadhyay
a,b,Patricia L.M. Plummer b,c
LIScience Department, Columbia College, Columbia, MO 65216, USA b Department of Physics and Astronomy, University of Missouri-Columbia, Columbia, MO 65211, USA ’ Department of Chemistry, University of Missouri-Columbia, Columbia, MO 65211, USA
Received 21 June 1993; in final fom 7 December 1993
Abstract Ab initio molecular structure calculations have been used to examine the association of NH3 with HCN to form a heterodimer. Because of the potential importance of electron correlation in weakly associated complexes, the structure and energy of the
heterodimers and the associated monomers are reported at both Hartree-Fock and Mgller-Plesset levels of calculation. Other calculated properties including harmonic vibrational frequencies and rotational constants are also reported for each system. Several stationary points on the potential surface were located and characterized. In the two most stable configurations, the monomer units are associated through a hydrogen bond with HCN acting as the hydrogen donor. The structure predicted to be the most stable has a linear hydrogen bond and C3” symmetry in accordance with the experimentally observed symmetric top spectrum. The stability of this structure is 5.28 kcal/mol (with MP2FULL/6-31+ G( 3df, 2p) and 5.01 kcal/mol with counterpoise correction). The second structure having only C, symmetry is identified as a transition state. A third stable hydrogen bonded structure has HCN acting as the hydrogen acceptor and is the most weakly bound of the three.
1. Introduction
Hydrogen bonding has been a subject of interest to scientists for most of this century, but the weakness of the bonding, together with the general flexibility of the hydrogen bond has made direct study difficult. It has been especially difficult to detect a C-H.. .X bond as a C-H group is a weak proton donor. Only since the late 70’s have the experimental and theoretical methods become sufficiently sensitive to provide quantitative as well as qualitative characterization of hydrogen bond properties in different systems. Specifically, the development of molecular-beam techniques has made high-resolution rotational spectroscopic methods available to probe the structure of 0301-0104/94/$07.00 0 1994 Elsevier Science B.V. All rights reserved SSDIO301-0104 (93)E0445-2
weakly bound molecular complexes #I. These studies have provided the geometry of many such systems and, as a result, a much enhanced insight into the nature of the bonding in these systems. However many questions remain unanswered, especially as to the role of the “chemically useful” practice of partitioning of the interaction energy into electrostatic, charge transfer, polarization, dispersion, and exchange. Prior to the current generation of computers, there was considerable interest in decomposing the Hartree-Fock SCF interaction energy in this fashion. The motivation for such partitioning was two-fold: first, to increase the understanding of the association mechanism in the weakly ‘I For a discussion on microwave and infrared spectroscopy see refs. [ 1,2], respectively.
40
S Chattopadhyay,
P.L.M. Plummer/Chemical
bound molecular complexes, and second, to design empirical models which could predict structures without the need to carry out the ab initio calculations. Additionally, the limited basis sets which could be used at the time gave rise to large basis set superposition errors and large uncertainties in both the structure and the stability of weakly bound systems. It was found by Morokuma et al. [ 31 and Kollman et al. [43 that for many weakly bound dimers, the electrostatic energy dominates the binding. In many hydrogen-bonded dimers, this dominance could be explained by the fact that the exchange repulsion energy was nearly balanced by the attractive polarization and charge-transfer terms. Legon and Millen [ 51 in reporting their experimental results, proposed an empirical electrostatic model to predict the shape of H-bonded dimers. Later Buckingham and Fowler [6] used more accurate multipole expansions to predict the structure of a variety of weakly bound molecular dimers, including both Hbonded and non-H-bonded complexes. The models based on electrostatic interactions have had considerable success in predicting the shapes of many molecular dimers. However these models are less reliable in predicting which bonding pattern is preferred when more than one isomer is possible. The electrostatic models have also been less successful in predicting the structure when one or more of the monomer units contained multiple bonds. In addition, most of the electrostatic models require keeping the intramolecular distances between the heavy atoms fixed at their experimental monomeric geometries. This restriction limits the usefulness of the models when experimental data is lacking or when cluster formation distorts the monomers. The general reliability of the early electrostatic models was criticized [ 71 because of the large errors in many experimentally determined low-order moments and in the inaccuracies arising from using truncated distributed multipole expansions. There were also sizable errors in the models based on multipole expansions from small basis set SCF calculations. The advances of the decade in both the computational hardware and software available to investigate small van der Waals’ and hydrogenbonded clusters, makes the time right to reinvestigate the bonding in some weakly bound clusters using ab initio quantum-chemical techniques including postHartree-Fock analysis. This study continues our examination of the structure and stability of van der Waals’ clusters. Specifically we examine the association of
Physics 182 (1994) 39-51
HCN with NH,. Both molecules are capable of serving as donor or acceptor molecules in the hydrogen bond and thus this system provides a variety of bonding possibilities. In an early experimental study, Legon and Millen [ 51 measured the rotational gas-phase spectra of several hydrogen-bonded heterodimers (dimers formed between the monomers of two different molecules) of HCN. Later Buckingham and Fowler [ 61 applied their electrostatic model for the prediction of structures of hydrogen-bonded dimers, to the dimer between NH, and HCN. According to Buckingham and Fowler [ 61, some systems may exhibit a stable anti-hydrogen bonded dimer instead of, or in addition to, a hydrogen bonded structure. In a previous study, we [S] have shown that such a dimer exists between HCN and SOZ. But for ammonia and hydrogen cyanide Buckingham and Fowler considered only the dimer having HCN as the donor molecule. Four experimental studies [9-121 have examined the HCN...NH, system. The complex H,N.HCN, which is a solid at room temperature, belongs to a special class of weakly bound binary complexes whose subunits form crystalline solids if equal molar quantities of them are taken, In the condensed phase ammonia and hydrogen cyanide exist as an ion pair but no experimental evidence has been found for the ion pair in the gas phase. The first experimental study [ 91 of the system was an examination of the vapor over NH&N by low resolution infrared spectroscopy. Three fundamental frequencies of the NH,--HCN complex were identified. The structure of the complex was assumed to be C,, with the HCN unit strictly linear. In the second study [ lo], Fraser et al. examined the microwave rotational spectrum using the molecular beam electric technique. They observed a symmetric top spectrum for the complex which they identified as a hydrogen bonded dimer with HCN serving as the H donor. In the most recent study [ I 11, Bohn and Andrews measured the infrared spectrum of the complex in an argon matrix. They also found a C,, symmetry for the dimercomplex. Smith and Yarwood [ 121 used Fourier transform IR to reexamine the heterodimer in the gas phase. They state that the dissociation energy of NH3 . HCN is unknown. Their attempt to dissociate the dimer using an IRmicrowave double resonance experiment produced no evidence of dissociation. Using experimental results for
S. Chattopadhyay, P.L.M. Plwnmer/ Chemical Physics 182 (1994) 3%51
HCN .HF and (HCN),, they suggest a dissociation energy of = 3.4 kcal/mol. For the system involving ammonia and hydrogen cyanide, previous theoretical studies of both homodimers and heterodimers have been reported [ 13-191. The majority of these theoretical studies were carried out in the seventies and early eighties using either the minimal basis sets or not including polarization functions. The previous ab initio studies [ 15-181 for the NH,-HCN heterodimer were limited to the experimentally identified hydrogen bonded structure. In addition, none of these earlier studies included electron correlation effects which have been shown to be important in the description of weakly bound systems. In a recent review, Hobza and Zahradnik [ 191 state “. . .the role of correlation contribution to the interaction energy is very important. There are no types of van der Waals’ molecules for which the correlation contribution could be neglected.” These authors also suggest that many-body perturbation treatment such as the Moller-Plesset method [20] “is the most promising for application to extended van der Waals’ molecules”. The present study was undertaken to confirm the most stable structure of the heterodimer and to examine the potential surface for other stable or metastable states. Of particular interest was whether a stable anti-hydrogen bonded structure for the dimer between HCN and NH3 could be identified. In addition, we wanted to examine, for the gas-phase dimer, the importance of the role played by the electrostatic contribution in the long-range interactions and in hydrogen bond formation. We also wanted to include electron correlation effects to increase the accuracy of the calculated vibrational frequencies for the complexes as well as improve the description of the bonding region.
41
grams employed for these studies were GAUSSIAN 86-92 [ 23,421. In addition, we have carried out the MP4 calculations, including contributions from single, double, triple and quadruple excitations, at the optimized geometry obtained from the MP2 calculations. Energies at the configuration interaction with single and double excitations including the inner shells (CISD) were also obtained. To estimate the stability of the dimers, the energy was compared with the sum of the energies of the optimized monomers. Because basis set superposition errors can be large compared to the interaction energy in van der Waals’ dimers, counterpoise calculations [24] at the HF and MP2 levels for the optimized dimer structure were carried out for the fragment monomers. The sum of these energies, when compared to the dimer energy, serves to put a lower limit on the dimer stability. For the frequency calculations at both the HF and MP2 levels, the analytical calculations of the second derivatives were used. Because of the possibility of the existence of multiple minima and transition states for weakly bound complexes, the Hessian matrix was examined for all stationary point geometries located during the optimization. The calculations were carried out on IBM-438 1 with FPS- 164 attached processor and IBM3090 machines at the University of Missouri Columbia and IBM-3090 at IBM-Los Angeles and IBM-Cambridge Scientific Centers. For displaying the molecules, we have used the software package PLUTO [ 251. The structures for the monomers and the most stable CsV dimer structure were reoptimized at the HF and MP2 (both at the frozen core and the full) levels using the larger 6-3 1 lG( d, p) basis set [ 261. This basis contains additional polarization functions on the hydrogens. Thus it provides additional flexibility for the valence electrons and in the hydrogen bonding region.
2. Method of calculation 3. Results and discussions We have employed the gradient optimization techniques to fully optimize the geometries of the heterodimers at the Hartree-Fock (HF) level. The monomers were also fully optimized at this level. The potential surface searches were carried out using the 6-31G(d) basis set [ 21,221. The structures corresponding to stationary points were reoptimized in the Moller-Plesset [ 201 (MPn) calculations at the second order including all orbitals in the correlation (MP2FULL). The pro-
Since our ultimate goal is to extend these calculations to the study of larger clusters (containing up to ten molecules), the size of the basis set is a critical component of the investigation. It must be large enough to provide a realistic description of the structure and relative stability of the cluster while being small enough to permit the larger clusters to be investigated. However it is well documented that the stability is over-estimated
42
S. Chattopadhyay. P.L.M. Plummer/Chemical
if the basis is too limited and the predicted structure may be incorrect as well. As an example we note that, in the case of the dimer H,N.. .HCl, the use of the 43 1G basis set [ 271 predicted an incorrect ionic structure, whereas the use of 6-31G(d) basis set gave the correct hydrogen-bonded structure [ 281. For the description of molecular clusters containing lone pair electrons, multiple bonds or hypervalent atoms, d functions appear essential [ 291. Frisch et al.‘s [ 17 ] calculations on the homodimer of NH,, and Latajka and Scheiner’s [ 281 calculation on H,N.. .HCl found that the 6-3 1G(d) basis set provided a good description of the interaction of dimers containing ammonia. The studies of George et al. [ 301, and Chattopadhyay and Plummer [ 81 find the 6-3 lG(d) basis set yields geometrical parameters for HCN in reasonable agreement with the experiment and MP2 harmonic frequencies within 1% to 4% (about 5 to 80 cm-‘) of the experimental harmonic frequencies. Even so, the dipole moments calculated for the heterodimers of HCN are usually slightly higher than their experimental values (but within 10%) andMP2/6-31G(d) lengthsformultiple bonds such as CN are generally overestimated in both the monomers and the dimers [ 30.3 I]. Tables 1 and 2 give the geometries and energies calculated with the6-3lG(d) andthe6-3llG(d,p) basesforthemonomers and the most stable heterodimer. Examination of the results obtained with these two basis sets found no qualitative and little quantitative differences of the geometries or stabilities at either the I-IF or the MP2 levels of calculation, The MP2 harmonic frequencies calculated with the 6-3lG(d) basis set agree with experiment somewhat better than those obtained with the larger basis set. Additionally the 6-3 lG( d) basis is one of the smallest basis which can be used with confidence in post-Hartree-Fock calculations to include electron correlation. For weakly bound systems, the first polarization functions play the dominant role in the correlation energy. Thus, it is not surprising that the additional polarization functions do not significantly alter the results. For the van der Waals’ complexes associated by hydrogen bonds the major part of the stabilization originates from the Hartree-Fock interaction energy. The 6-31G basis with heavy atom polarization functions provides the required accuracy for the molecular geometries and for the monomer dipole moments especially if MP2 electron density is used for dipole calculation. Thus we believe that the results
Physrcs 182 (1994) 39-51
obtained with the 6-31G(d) basis can be used with confidence for the dimers and this basis has the advantage of being extendable to the larger clusters as well. The 6-31G(d) basis was used to search the HF potential surface. Both the hydrogen- and the antihydrogen-bonded structures were investigated. Jones et al. [ 93 had suggested that this latter structure would be more rigid than the hydrogen-bonded structure but found no conclusive evidence of its existence in their experiment. Our potential surface search also failed to find a stable anti-hydrogen-bonded structure for the dimer. On the other hand, three stable (relative to dissociation to monomers) hydrogen-bonded isomers were identified. One of the structures had the ammonia serving as the hydrogen donor while in the other two it was the acceptor. Of the two structures having ammonia as the acceptor m the hydrogen bond, one was identified as a transition state after examining the frequencies associated with it (see Fig. 1 and Table 6). The H bond stabilization energy for the transition state structure at the HF/6-3lG(d) level was found to be 5.79 kcal/mol with respect to the optimized structure for the fragments and is reduced by = 1 kcal/mol to 4.82 kcal/mol by counterpoise calculations. The other structure (Fig. 2)) which was a minimum on the potential surface, had a binding energy at the HF/6-3 lG(d) level of 7.93 kcal/mol (counterpoise corrected to 6.97 kcal/mol) and at the MP2FULL/6-3lG(d) level of 9.29 kcal/mol. (When corrected for zero-point energies the stability of the dimer is 7.2 kcal/mol at the MP2 level.) We note here that the hydrogen bond energy predicted by the electrostatic model of Buckingham and Fowler [ 61 for this structure of the dimer was 9.35 kcal/mol. In two previous ab initio studies using different basis sets, the corresponding hydrogenbond energy was found to be 7.9 kcal/mol (by Vishveshwara [ 16 J using STO-3G basis set) and 9.0 kcal/ mol respectively (by Kollman et al. [ 151 using 4-3 IG) Somasundaram et al. [ 181 estimate a hydrogen-bond energy of 7 kcal/mol at the double zeta level but do not report a stabilization energy for the dimer using their larger basis set. The third structure (having HCN as the acceptor molecule) has an association energy of 2.05 kcal/mol at the HF level when compared to the optimized structures for the monomers and 1.1 kcal/mol when counterpoise corrections were taken into account. The stability at the MP2FULL level is 2.71 kcal/mol (1.81 kcal/mol when zero-point ener-
S. C~tto~byay,
P.L.M. Plummeri Ckmical Physics 182 (1994) 39-B
Table 1 Comparison of bond lengths and bond angles between this study and previous theoretical (All bond lengths are in .&and bond angles in deg) Bond iengthsiaogles
Reference level
Ammonia
Hydrogen
cyanide
and experimental
43
results at HF and MPLFULL levels.
NHIImHCN (TS3
NH,.
NH,. NCH
a
C-N
HF/6-3lG(d) MP2/6-3X(d) MP2/6-31iG(d, p) MP216-31 +G(3df, 2p) previous talc. ’ experimental
1.1325 1.1761 b i.1697 1.1677 1.1242 1.156d. 1.153”
1.1336
1.134 1.1769 1.1707 1.1686 1.1256
1.132 1.1753
C-H
HF/6-3lG(d) MP2/6-31G(d) MP2/6-311G(d, p) MP2/6-31 +G{3df, 2p) previous talc. c experimental ’
1.0590 1.0695 b 1.0671 1.0652 I .Q587 1.0640
1.065
1.0721 1 .Og72 1.0854 1 .oa04
1.0593 1.0693
N2-H
HF/6-3lG(d) MP2/6-31G(d) MP2/6-31lG(d,p) MP2/6-31 +G(ldf, previous cak. ’ experimental d
2~)
Nl-H4
HF/6-3lG(d) MP2/6-31G(d) MP2/6-31lG(d, p) MP2/6-31+ G( 3df, 2p) previous talc. ‘rf experimental s
LHlNlH2
HF/6-3lG(d) MP2/6-31G(d) MP2/6-311G(d, p) MP216-31 +GfJdf, 2~) previous talc. ’ experimental d
LHINlH4
HF/6-31G(d) MP2/6-31C(d) MP2/6-311G(d, p) MP2/6-31 +G(3df,2p)
1.0026 1.0167 1.0131 1.0109 1.Oao 1.012
107.2 106.4 106.1 106.2 108.2 106.7
a For labeling of atoms, refer to Figs. l-3. hRef. [S]i. c Ref. [32]. d Ref. [33]. e Ref~ [34]. f There are some other calculated values available for this hydrogen-bond the best one is mpotied here. See refs. [ 1%19]. a Ref. [8]
gies are taken into account). The final structure for this conformer is shown in Fig. 3, The electrostatic model of Kollman et al. [ 1.51 predicts the stability of this isomer to be 2.3 kcal/mol (1.2 kcallmol when scaling factor is applied).
1.0702 1.0041
1.0033 1.0177 1.0141 1.0179 1.0516
1.0029 1.0178
2.35 I
2.146 2.0536 2.0553 2.0583 2.220 2.16
2.6115 2.4233
106.36
107.0 106.4 106.4 106.4 107.5
106.9 106.1
125.2
111.9 112.2 112.3 112.2
length, obtained from ab initio and other theoretical
calculations,
only
The detailed bond length and bond angles for the stable structures are reported in Table 1. The table shows the optimized geometry predicted with the 631G(d) and 6-311G(d, p) basis sets at the HF and MP2 levels of calculation. The values are compared with experimental results and earlier calculations,
S. Chattopadhyay. P.L.M. Plummer/Chemica~ Phystcs 182 (1994) 39-J-1
44
Table 2 Energies obtained from Hartree-Fock in hartree) a
and Moller-Plesset
NH3 HF/6-31G(d) CP/HF/6-31G
calculations
on the systems NHa, HCN, NH,. HCN and NH,. NCH. (All energies are
HCN
NH3.HCN
(TS)
- 56.18385 -56.18531
- 92.87520 - 92.87526
- 56 20999
- 92.895 16
MP2/FC/6-31G (d) CP/MP2FC/6-31G (d)
- 56.3542 I -56 35655
-93.15893 -93.15916
- 149.52342
- 149.52774
MP3/FC/6-31G (d) MP4D/FC/6-31G (d) MP4DQ/FC/6-31G(d) MP4SDQ/FC/6-31G (d) MP4SDTQ/FC/6_31G(d) CISD, 4/K/6-31G (d)
- 56.36595 - 56.36987 - 56.36823 - 56.36894 - 56.37126 -56.36194
-93.15772 -93.16753 -93 16152 -93.16630 - 93.18069 - 93.13996
-
-
MP2/FULL/6-31Gld) MP3/FULL/6-31G(d) MP4D/FULL/6-31G(d)
- 56.35738 - 56 36894 - 56.37286
-93.16680 -93.16544 -93.17535
- 149.53467 - 149.54459 - 149.55817
- 149.53912 - 149.54883 - 149.56237
MP4DQ/FULL/6-31G (d) MP4SDQ/FULL/6_31G (d) MP4SDTQ/FULL/6-31G(d) CISD, 4/FULL/6-31G (d)
-5637119 - 56.37190 - 56.37424 - 56.36448
-93.16911 - 93.17392 - 93.18842 -93.14691
-
-
MP2/FULL/6-311G
-5642788
- 93.23736
- 56.43945 - 56.4406 1
- 93.24288 -93.24319 _
HF/6-311G
(d) (d, p)
(d. p)
MP2/FULL/6-31+G (3df. 2p) CP/MP~FULL/~-~I+G(~~~,~P) MP4/FULL/631+G(3df,2p)
- 149.06826
NH, ’ HCN - 14907168
NH,. NCH - 149.06232
- 149.11693
149.53379 149.54734 149.53961 149.54516 149.56206 149.47812
149 55027 149.55586 149.57294 149.48699
149.53789 149 55140 149.54362 149.54917 149.56619 149.48197 - 149.52850 _
149.55442 149.56001 149.57721 149 49089
- 149 67920 _
- 149.69345 - 149.73965
a Stability (m kcal /mol ) of the dimer at 6-3 1G (d ) level with respect to: ( i ) tbe optimized monomer geometry: (a) for the TS HF structure: 5.79, (b) for the optnnized HF structure: 7.93, (c) for the optimized MP2 structure: 9.29, (d) corrected for zero-point energy: 7.34, (ii) the counterpoise calculations: (a) for the TS structure: 4.83, (b) for the optimized HF structure: 6.97, (c) for the optimized MP2 structure: 7.55, (d) corrected for zero-point energy 5.60. Stability (in kcal/ mol ) of the dimer at 6-311G (d, p) level with respect to: (i ) the optimized monomer geometry: (a) for the optimized HF structure: 7.40. (b) for the optimized MP2 structure. 8.76. Stability (m kcal/mol) of the dimer at 631 + G (3df, 2~) level with respect to: (i) the optimized monomer geometry. (a) for the optimized MP2 structure: 6.97, (b 1 counterpoise MP2. 6.06, (c) counterpoise and corrected for zero-point energy: 4.37
where available. The optimized geometry (for the most stable structure) obtained in this study, both at HF and MP2 levels, has a CaVsymmetry within the accuracy of theoretical calculations. In the third hydrogen-bonded structure, the ammonia acts as a hydrogen donor. The nitrogen of the C=N was the hydrogen acceptor. This configuration was the least stable of the three and is characterized by a very long hydrogen bond. The monomer fragments are only minimally distorted from their isolated geometries. This lack of stability and a long hydrogen bond are not surprising given the reluctance of nitrogen to serve as a proton donor. (Fraser et al. [ 35 ] state: “It appears that NH has little propensity for proton donation.. .” ) What is, perhaps, surprising is the existence of this structure for the heterodimer as
a local minimum on the potential surface. A second transition state having this bonding pattern was identified but was found to have an energy higher than the isolated monomer units. Detailed results for this second TS are not reported in this paper. The energies of the monomers and the dimers are reported in Table 2. Post-Hartree-Fock calculations were considered necessary for the dimers to account for the electron correlations involved. The geometry was first optimized at the MP2 level of calculation for all three structures with all electrons used in the correlation. With that geometry, single-point MP3 and MP4 calculations (with SDTQ excitations) were carried out. The Moller-Plesset calculations at all levels were done both with frozen core (FC) and full corre-
S. Chattopadhyay, :K
tat12 TRmSl
P.L.M. Plummer / Chemical Physics 182 (1994) 39-51
7, ON STATE
N2
Cl
8
"4
a
n2
WI
“1
Fig. I. A transition state for the dimer between ammonia and hydrogen cyanide.
determined at the HF level. As in our earlier report, we have scaled the frequencies obtained at the HF level by a factor of 0.89 as suggested by Nobes et al. [ 371. The values for the rotational constants from both HF and MP2 calculations are reported in table 3. Except in the case of B0 for ammonia, where the difference was found to be = 8% at the HF level, all other calculated values for the constants differ by less than 5% from their experimental counterparts. The vibrational harmonic frequencies for the monomers are reported in table 4. As expected the values obtained for these quantities at the MP2 level are closer to the experimental values than those obtained from the HF calculations. In table 5, we have reported the values of the harmonic vibrational constants obtained from the HF calculations for the stable dimer, their scaled values, their MP2 counterparts, the values obtained from previous calculations and also the experimental values wherever they are available. We also tried to identify the individual monomer contributions to the dimer spectrum. In Table 6, the frequencies associated with the stable transition state of dimer and the other minimum energy structure are reported. As no experimental values for these quantities are available for these structures, only the unscaled numbers are given. NAL
NHS+CH
OPTI U 7X0
Fig. 2. The most stable structure for the dimer between ammoniaand hydrogen cyanide.
lations (FULL). The results of the energies obtained at the HF and MP2, MP3 and MP4 levels are all included in Table 2 for comparison. Harmonic vibrational constants were calculated for the monomers and for the structures of the dimer at both HF and MP2 level. It has been suggested by Hehre et al. [36] and substantiated in a previous study by Chattopadhyay and Plummer [ 81 that MP2 frequencies are in better agreement with experiment than those
45
STRUCTURE
N2
Cl
H
N4
Fig. 3. A third hydrogen-bonded structure for the dimer between ammonia and hydrogen cyanide with one of the hydrogens of ammonia bonded to the nitrogen of the hydrogen cyanide subunit.
S. Chattopadhyuy, P. L.M. Plummer / Chemical Physrcs I82 (I 994) 39-51
46 Table 3 Comparison
of rotational constants.
(All frequencies
are m GHz)
System ammonia
hydrogen
HF/6-31G(d) expenmental d MP2/6-3lG(d) MP2/6-31 lG(d, p) MP2/6-3 1 + G( 3df, 2p)
305.7 _ 295.1 291.2 298.8
305.7 283.3 295.7 297.2 298.8
HF/6-31G(d)
cyanide
45 584 44.341 43 004 43.400 43 586
43 004 43.400 43.586
116.1
3.15
3 I1
192.7 _
3 04 3.02 f .00024
3 04 _
182.0 192.1 198.8 190.9 191.5
3 07 2.95 3.10 3.10 3.06
2.95 3.10 3 10 3 06
expertmental * MP2/6-31G(d) MP2/6-311G(d, p) MP2/6-3 1 + G( 3df, 2p) NH,.HCN(TS)
HF/6-31G(d)
NH,. HCN
HF/6-31G(d) experimental
h
experimental ‘ prevtous talc. ’ MP2/6-31G(d) MP2/6-311G(d,p) MP2/6-3 I+ G( 3df. 2p)
192.6 1859 189.3 191.4 191.8 45 584 _
“Ref. [33]. ‘Ref. [IO]. c Ref. [ 121 (note that these values are reported for the V, excited state of the dimer). “Ref [32].
In Table 7 we report the dipole moments of the monomers as well as those for the dimers obtained from the ab initio calculations to show the change accompanying dimer formation. Both the experimental and calculated Table 4 Comparison In cm-‘)
of calculated
harmonic vibrattonal
frequencies
dipole moments for the dimer exceed the sum of the monomer moments. We see from the Table that the ab initio calculations overestimate the dipole moment for both the monomers and the dimer, when compared with
for ammonia and hydrogen cyanide with experimental
values. (All frequenctes
are
System
Reference d
01
eJ2
wr
04
Ea (kcal/mol)
ammonia
HF/6-31G(d) (scaled MP2/6-31G(d) MP2/6-31 lG(d, p) MP2/6-3 1 + G( 3df. 2p) experimental h
3690 3284 3505 3532 3518 3337
1209 1076 1160 112s 1080 968
3823 3402 3661 3676 3670 3444
1850 1647) 1755 1678 1684 1627
23 22
HF/6-31G(d) (scaled MP2/6-31G(d) MP2/6-31 lG(d, p) MP2/6-3 1 + G( 3df, 2p) expertmental b experimental ’ (harmonic)
2445 2176 2046 2027 2039 2097
892 794 734 767 736 714
3676 3272) 3509 3495 3485 3312
2128
727
3440
hydrogen
cyanide
“Scalefactoris0.89forHF/6-31G(d).
bRef. [33].
‘Ref.
[38].
22 16 21.97 21.88
11.29 10.04 10.09 10.00
S. Chattopadhyay, P.L.M. Plummer/Chemical
Physics 182 (1994) 39-51
47
Table 5 Calculated harmonic vibrational frequencies at HF and MP2 levels for NHs. HCN (unscaled) and corresponding monomer frequencies calculated at MP2 level and their comparison with previously calculated and experimental values. (All frequencies are in cm-‘) a Harmonic frequency No.
oat HF level b No. 1
oat MP2 level b No. 1
Previous theoretical calculations”
Experimental ref. [ 1 l] d
Corresponding frequency ref. [9]
monomer
Remarks
ref. [ 121
No. 3
at MP2 level b
experimental
e
No. 1 No. 2 No. 3 192 3 4.5 6,7 8 9, 10 11 12 13 14, 15
126(2) 168 312(2)
123(2) 190 308
132(2) 158 283(2)
1062(2) 1275 1846(2) 2412 3689 3491 3818(2)
1008(Z) 1215 1749(2) 2038 3506 3251 3656(2)
953(2) 1133 1680(2) 2024 3507 3251 3656(2)
113 146 302 lOOO(2) 1203 799(2) 2382 3682 3439 3799(2)
‘bend” “stretch” “shear”
140 f 25 279 878 * 1055 2077 3305 3022 *
1040 2086 3110
3110.5*0.2
734 1160 1755 2046 3505 3509 3661
767 1125 1678 2027 3522 3495 3676
736 1080 1684 2039 3518 3485 3670
727(721) 968(974) 1627 2128(2093) 3337 3440(3305) 3444
0, of HCN 0, of NH, o, of NH, 0, of HCN o, of NH, 0, of HCN o, of NH,
’ Note that the degenerate frequencies are grouped together, &(MP2/6-31G(d)) = 34.30 kcal/mol. b Basis No. 1 is 6-31G(d), basis No. 2 is 6-31 lG(d, p) and basis No. 3 is 6-31 +G(3df, 2~). ’ Previous calculations were done using the TZ + 2P basis set, ref. [ 183. dThe asterisk refers to values obtained from the argon spectrum. ’ Experimental values for ammonia are obtained from ref. [ 331 and those for HCN are from ref. [ 381. Values in parentheses (table 1) , in argon matrix.
the experimental values. Prior to GAUSSIAN 92 [ 421 only SCF electron densities were used in the calculation of the dipole moments, even for MP4 calculations. In GAUSSIAN 92, use of MP2 electron density became an option. Use of MP2 density improved the comparison with experiment for all basis sets employed. The closest value for the dipole moment for the dimer from any previous ab initio calculation was found by Vishveshwara [ 161 using the STO-3G basis set, although the other parameters obtained from that study, such as the hydrogen-bond length or the charge densities were far from being realistic. The charge distribution in the dimers and monomers for the different levels of calculation are shown in table 8. For the Csv dimer, at the MP2 level, the hydrogens in the ammonia subunit have equal charges on them. This suggests that this level of theory provides a realistic description of the system. It can be seen that in the dimer, the NH3 subunit has a charge of about + 0.02 e at both HF and MP2 levels. Thus, the complex formation is accompanied by a minimal charge transfer
are from ref. [ 1 l]
between the monomer units and is not a charge transfer complex. However, as expected, both subunits have been polarized on complex formation with a greater change in the HCN subunit. The presence of the ammonia thus increase the basicity of the CN group by about 20% relative to the isolated molecule. The polarization is less in the transition state and symmetry of charges on the hydrogens (of ammonia) is reduced. In the third structure, the complex formation is accompanied by even less (about 0.007 e and 0.009 e at HF and MP2FULL levels, respectively) charge transfer between the subunits. Thus these calculations show ammonia to be a more effective Lewis base than a Lewis acid when interacting with HCN. Examination of the orbital energies for the most stable Csv structure reveals a shift of = 6% relative to the monomers, decreasing for NH3 and increasing for HCN. For NH3 functioning as a Lewis acid (H donor), the energies shift considerably less in the opposite direction. In the condensed phase, ammonia and hydrogen cyanide exist as an ion pair: NH: CN- . The experimental
S. Chattopudhyay,
48
P.L.M. Plummer/Chemical
Table 6 Calculated harmonic vibrational frequencies (unscaled) with 63 1G(d) basis set for NH,. HCN (transition state at HF level) and NH,. NCH (at HF and MP2 levels) dimers. (All frequencies are in cm-‘) Mode No.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
Frequency for NH,. HCN ’ at HF level
Frequency
for NH,. NCH b
at HF level
at MP2 level
- 54.7 9.4 139.7 252.4 287.8 903.5 988.3 1292.8 1850.1 1851 2 2477.9 3608.8 3679.8 3802.9 3806 8
28 38 77 144 187 893 894 1245 1852 1863 2441 3678 3688 3813 3823
30 40 104 187 227 729 729 1199 1758 1778 205 1 3503 3512 3654 3656
d &(TS) =35.599 kcal/mol. hEo(HF)=35.115kcal/molandE,(MP2)=33.105kcal/mol.
[9-l 21 of this system found no evidence of the ion pair in the gas phase or in the argon matrix. In the course of this study, we made a limited investigation of the energy associated with a proton transfer reaction. Holding the NH, and CN groups fixed at their optimized positions in the stable C,, dimer, the energy was studies
Physics 182 (1994) 39-51
calculated for the movement of hydrogen along the bond from CN to within 0.8 ,& of the nitrogen. A maximum in energy was found for a N-H distance of 1.3 A. For that geometry, the energy was optimized relaxing all of the constraints with the exception of the NC...N distance. Examination of the Hessian matrix for the structure thus obtained showed it to be a transition state having an energy of = 44.5 kcal/mol above that of the most stable geometry. Thus, it appears from this limited study that the energy barrier for the hydrogen transfer for the gas-phase dimer is substantially above the energy of the separated monomers. Examination of the electron distribution for this metastable state showed that this hydrogen was no more positive than the hydrogens from the ammonia. Thus this reaction path is better described as hydrogen transfer, rather than a proton transfer reaction. For NH, reacting with hydrogen halides, large cluster complex formation mvolves proton transfer [43]. In the gas phase, estimates for the enthalpy for NH, + HX -+ NH: X- range from 8 + 6.5 to 17 + 6.5 kcal/mol as X changes from I to Cl. Since many properties of HCN are very similar to HX with certain cyano compounds, they are often referred to as pseudo halogens. Thus we expect that proton transfer may also occur in larger clusters in the NH,-HCN system. The cluster size at which the proton transfer process becomes energetically favorable is of special interest. We [ 44.451 are currently in the process of extending these calculations to larger clusters to investigate the process.
Table 7 Calculated dipole moments of the monomers and the dimers obtained from ab initio calculations (in units of D )
and comparison
with the experimental
values
System
From HF
From MP2 d
From previous calculation
Experimental
HCN/6-31G (d) /6-3llG(d,p) /6-31+G(3df,2p) NH,/6-31G(d) /6-311G(d,p) /6-31+G(3df,2p) NH3.HCN (TS) /6-31G (d) NH,. HCN /6-31G (d) /6-311G (d, p) /6-31+G(3df, 2p)
3.209 3.205 (3.032) 1.922 1.728 (1.593)
3.209 (2.950) 3.244
3.25 ’
2.985 ’
1.969(1 967) 1.802
2.30 d
1.471 e
4.18 ‘, 5.03 g, 5.37’
5.2608 i. 0.005 ’
5.529 5.780 5.821 (5.541)
5.778 6.192 (5.948) 6.059
a Values in parentheses are computed with MP2 electron density - an option not available until G92. ‘Ref [40]. dRef. [15] ‘Ref [41]. ‘Ref [13]. gRef. [141. hRef. I161. ‘Ref [39]
“Ref.
[lo]
S. Chattopadhyay,
P.L.M. Plummer/
Chemical Physics 182 (1994) 3%51
49
Table 8 Charge distribution on the atoms (in terms of electronic charge) a System
Reference
HCN
HF/6-31G(d) MP2/6-31G(d) MP2/6-31 lG(d, p)
NH3
NH?. HCN (TS) NH3 HCN (Cw) NH,. NCH
onN1
on Hl
on H2
on H3
on H4
on C
on N2
0.313 0.320 0.205
0.066 0.062 0.084
- 0.379 -0.381 - 0.289
HF/6-31G(d) MP2/6-31G(d) MP2/6-31 lG(d, p)
- 0.996 - 0.989 - 0.556
0.332 0.330 0.185
0.332 0.330 0.185
0.332 0.330 0.185
HF/6-31G(d) MP2/6-31G(d)
- 1.022 - 1.022
0.352 0.354
0.345 0.344
0.343 0.343
0.350 0.364
0.054 0.048
-0.421 -0.431
HF/6-31G(d) MP2/6-31G(d) MP2/6-311G(d,
- 1.025 - 1.039 - 0.582
0.348 0.354 0.211
0.352 0.354 0.211
0.344 0.354 0.211
0.362 0.415 0.261
0.046 0.020 0.011
- 0.427 -0.458 - 0.325
- 1.021 - 1.022
0.324 0.320
0.324 0.320
0.366 0.373
0.319 0.328
0.081 0.079
- 0.393 - 0.398
HF/6-31G(d) MP2/6-31G(d)
p)
’ Refer to figs. l-3 for the labeling of the atoms.
In addition to the calculations described above, another set was performed to find out if there exists a stable T-shaped structure with the HCN subunit perpendicular to the C3” axis of NH3 (acting as a mirror plane). Such a structure is reported [ 461 to be preferred by the acetylene dimer ( ( HC=CH)2) whose association is partly stabilized by hydrogen bonding. A Tshaped structure is the only stable minimum reported [ 471 for the (H,) 2 dimer. An additional motivation in our search for such a structure was to see if we could find another local potential minimum separated from the most stable structure by a saddlepoint. However, despite an extensive search of the potential surface, we were not able to locate such a state. We conclude that the transition state located in this study represents a barrier to dissociation by bending of the hydrogen bond. A generally flexible hydrogen bond as suggested by Fraser et al. [ lo] and Larson et al. [48] is consistent ‘with our characterization of the potential surface. We find it relatively flat for small variations of the planar angle between the two monomers; this plateau is surrounded by several metastable transition states, one of which we have fully characterized in this report.
4. Conclusions The structure of the most stable heterodimer of NH3 with HCN is predicted to have a hydrogen bond with
C,, symmetry and HCN as the hydrogen donor. The calculated length of the hydrogen bond and the harmonic frequencies reported are in reasonable agreement with experimental observations. At the HF/ 6-31G(d) level, the stabilization energy for this conformer is predicted to be = 8 kcal/mol compared with 5 kcal/mol [ 181 and 6.29 kcal/mol[32] for (HCN),. As with many hydrogen bonded complexes, the stability calculated at the HF level is comparable to that predicted from the electrostatic interaction of the monomers. Inclusion of electron correlation increases the predicted stabilization energy by = 15%, while inclusion of zero-point vibrational energy reduces it by =20%. In addition, two other structures have been characterized and predicted to be stable. No stable anti-hydrogen bonded structure was found. The relative stabilities of the heterodimers at the HF level compared to possible dissociation limits [ 491 are shown schematically in Fig. 4. The most stable heterodimer [50] of NH, with HNC and the associated asymptote are included for comparison. Thus our calculations at the HF level with 6-31G(d) basis set predict all three hydrogen bonded states to be stable, not only with respect to the monomer fragments, but also with respect to the counterpoise calculations, hydrogen transfer, or proton transfer and electron transfer processes. The ground state of the heterodimer exhibited minimal charge transfer. When N rather than C is the hydro-
50
S. Chattopadhyay, P.L.M. Plummer/Chemical NW-+cN*
NH4+ +
(174.1-l)
cm-(16n.lkalhml)
Physics 182 (1994) 39-51
Acknowledgement The authors wish to thank the IBM, Inc. for a substantial grant of computer time.
NIB + HNC (19.6 kdhnol)
Note added MtmNC(c3v) (92kL-alhd) NH3 +HCN NH3NCH
(7.9.kc&mI)
(5.9kdhd)
NHWCN(-l-S)W kdhol) NHtHCN(ar) (O-refemce) DIMERS
ASYMPTOTES
Fig. 4. A comparison of binding energies for different systems involving ammonia. hydrogen cyamde and hydrogen tsocyamde.
gen donor, as in the isocyanide-ammonia dimer, the charge transfer is about twice as great and the donated hydrogen is 30% more positive [ 501. Correspondingly the basicity of the CN group in isocyanide is 30% more than the monomer compared with a 20% increase for HCN. The investigation of the structure and stability of the heterodimers formed in this system continues our quest for a better understanding of bonding in weakly bonded molecular complexes. The eventual goal of these studies is the development of intermolecular potential functions which can describe the energy of interaction in terms of the relative separation and orientation of the monomer fragments. A more immediate goal is the characterization of stable stationary points on the potential surface. Calculations such as these are necessary complements to experimental studies. Nesbitt [ 5 1 ] has discussed the presence of multiple minima in the potential surface. In order to decide whether isomers are kinetically stable and distinct or whether they are better described as a vibrationally excited state having large amplitude motions along the interconversion coordinates depends on the height of the barriers separating the minima. Since spectral analysis cannot provide a unique determination of the intermolecular potential energy surface for weakly bound complexes, calculations such as we have presented can provide insights not presently available from experimental data.
Geometries obtained from the 6-3 1G(d) basis at the MP2 level using all electrons (MP2FULL) are used in Gl [ 521 and G2 [ 531 to predict the molecular energies and properties. The agreement between the calculated equilibrium geometries and experiment is excellent [ 541. Concern expressed by a referee about this basis lead us to reoptimize the geometries of the monomers and the most stable dimer at the MP’LFULL level with a larger basis set. This basis set (6-3 It G( 3df, 2p) ) adds diffuse functions and additional higher angular momentum terms to the heavy atoms and a second set of p functions on the hydrogens. As shown in Table 1, the equilibrium geometries are only slightly changed from those obtained with the 6-3 1G(d) basis set. Inclusion of zero point vibrational energies and counterpoise corrections predict a stability of the dimer of between 5.01 kcal/mol for the largest basis set to 5.60 kcal/mol for the smallest one used, the 6-3 lG(d) Thus, we again reiterate our belief that the use of 6-31G(d) with correlation is justified and can provide valuable insights when used in the study of larger cluster systems. References [ 11J.S. Muenter, in: Structure and dynamics of weakly bound
[2] [3]
[4]
[5] [6] [7] [8] [9]
molecular systems. ed. A. Weber, NATO AS1 Series, Vol. C3 ( 1987). D.J. Nesbitt, Chem. Rev. 88 (1988) 843. K. Morokuma, J. Chem Phys. 55 (1971) 1236, H. Umeyamaand K. Morokuma. J. Am. Chem Sot 98 ( 1976) 4400.99 (1977) 1316. P. Kollman, J. Am. Chem. Sot. 99 ( 1977) 4875; P. Kollman and L.C. Allen, J. Am. Chem. Sot. 93 ( 1971) 4991. A C Legon and D.J. Millen, Chem. Sot 73 (1982) 71, Chem. Rev. 86 (1986) 635 A.D Buckingham and P W. Fowler, Can. 1. Chem. 63 ( 1985) 2018 F.A. Bamcchi, W. Reiher and W. Klemperer, J. Chem. Phys. 79 (1983) 6428. S. Chattopadhyay and P L.M. Plummer, J. Chem. Phys. 93 (1990) 4187. W.J. Jones, R.M. See1 and N. Sheppard, Spectrochim. Acta 25A (1969) 385
S. Chattopadhyay,
P. L.M. Plummer / Chemical Physics I82 (1994)
[ 101 G.T. Fraser, K.R. Leopold, D.D. Nelson Jr., A. Tung and W. Klemperer, J. Chem. Phys. 80 (1984) 3073. [ 111 R.B. Bohn and L. Andrews, J. Phys. Chem. 93 ( 1989) 3974. [ 121 I.W. Smith and G. Yarwood, Mol. Phys. 64 (1988) 627. [ 131 A. Goel and C.N.R. Rao, Trans. Faraday Sot. 67 (1971) 2828. [ 141 D. Bonchev and P. Cremaschi, Theoret. Chim. Acta 35 ( 1974) 69. [ 151 P. Kollman, J. McKelvey, A. Johansson and S. Rothenberg, J. Am. Chem. Sot. 97 (1975) 955. [ 161 S. Vishveshwara, Chem. Phys. Letters 59 (1978) 26. [ 171 M.J. Frisch, J.A. Pople and J.E. Del Bene, J. Phys. Chem. 89 (1985) 3664. [ 181 K. Somasundaram, R.D. Amos and N.C. Handy, Theoret. Chim. Acta ( 1986) 491. [ 191 P. Hobza and R. Zahradnik, Chem. Rev. 88 ( 1988) 871. [20] C. Mller and M.S. Plesset, Phys. Rev. 46 (1934) 618; J.A. Pople, J.S. Binkley and R. Seeger, J. Quantum Chem. Symp. 10 (1976) I; J.A. Pople, R. Krishnan, J.S. Binkley and R. Seeger, Intern. J. Quantum Chem. 14 (1978) 545; R. Krishnan, M.J. Frisch and J.A. Pople, J. Chem. Phys. 72 ( 1980) 4244. [21] W.J. Hehre, R. Ditchfield and J.A. Pople, J. Chem. Phys. 56 (1972) 2257. [22] P.C. Hariharan and J.A. Pople, Theoret. Chim. Acta 28 ( 1973) 213. [23] M.J. Frisch, J.S.Binkley, H.B. Schlegel, K. Raghavachari, C.F. Mel@, R.L. Martin, J.J.P. Stewart, F.W. Bobrowicz, C.M. Rohlfing, L.R. Kahn, D.J. DeFrees, R. Seeger, R.A. Whiteside, D.J. Fox, E.M. Fluder and J.A. Pople, GAUSSIAN 86 (Carnegie-Mellon Quantum Chemistry Publishing Unit, Pittsburgh, PA, 1986). [24] M.J. Frisch, J.E. Del Bene, J.S. Binkley and H.F. Schaefer Ill, J. Chem. Phys. 84 (1986) 2279; D.W. Schwenke and D.G. Truhlar, J. Chem. Phys. 82 (1985) 2418. [25] W.D.S. Motherwell, PLUTO, Program for plotting molecular and crystal structures (Cambridge Univ. Press, Cambridge, 1976). [26] R. Krishnan, J.S. Binkley, R. Seeger and J.A. Pople, J. Chem. Phys. 72 (1980) 650. [27] R. Ditchfield, W.J. Hehre and J.A. Pople, J. Chem. Phys. 54 (1971) 7241. [28] Z. Latajkaand S. Scheiner, J. Chem. Phys. 82 (1985) 4131. -[29] S. Rothenberg and H.F. Schaefer Ill, J. Chem. Phys. 53 (1970) 3014. [30] P. George, M. Trachtman and C.W. Bock, J. Comput. Chem. 2 (1981) 334. [31] L.L. Lohr, J. Phys. Chem. 89 (1985) 3465.
39-51
51
[32] S.-K. Rhee and A. Karpfen, Chem. Phys. 120 (1988) 199. [ 331 G. Herzberg, Electronic spectra of polyatomic molecules (Van Nostrand Reinhold, New York, 1966). [ 341 G. Winnewisser, A.G. Maki and D.R. Johnson, J. Mol. Spectry. 39 (1971) 149. [35] G.T. Fraser, F.J. Lovas, R.D. Suenram, D.D. Nelson and W. Klemperer, J. Chem. Phys. 84 (1986) 5983. [36] W.J. Hehre, L. Radom, P. von R. Schleyer and J.A. Pople, Ab initio molecular orbital theory (Wiley, New York, 1986). [37] R.N. Nobes, J.A. Pople, L. Radom, N.C. Handy and P.J. Knowles, Chem. Phys. Letters 138 (1987) 481. [38] W. Quapp, J. Mol. Spectry 125 (1987) 122. [39] P.K. Pearson, G. Blackmun, H.F. Schaefer Ill, U. Wahlgren and B. Ross, Astrophys. J. 184 (1973) L 19. 1401 A.G. Maki, J. Phys. Chem. Ref. Data 3 (1974) 221. [41] M.D. Marshall and J.S. Muenter, J. Mol. Spectry 85 (1981) 322. [42] M. Frisch, J. Foresman and A. Frisch, GAUSSIAN 92 (Gaussian Inc., Pittsburgh, PA, 1992). [43] J.T. Cheung. D.A. Dixon and D.R. Herschback, J. Phys. Chem. 92 (1988) 2536. [ 441 P.L.M. Plummer and S. Chattopadhyay, Z. Physik D 26 ( 1993) S329. [45] P.L.M. Plummer, S. Chattopadhyay, A.M. Turner and B.M. Simpson, Z. Physik D 26 (1993) S 326; S. Chattopadhyay and P.L.M. Plummer, Z. Physik D 26 ( 1993) S 332. [46] D.G. Prichard, R.N. Nandi and J.S. Muenter, J. Chem. Phys. 89 ( 1988) 115. [ 471 B. Schneider, P. Hobza and R. Zahradnik, Theoret. Chim. Acta 73 (1988) 201. [48] J.W. Larson, J.E. Szulejko and T.B. McMahon, J. Am. Chem. Sot. 110 ( 1988) 7604. [49] R.A. Whiteside, M.J. Frischand J.A. Pople,eds.,The CarnegieMellon quantum chemistry archive, 3rd Ed. (Department of Chemistry, Carnegie-Mellon University, Pittsburgh, PA, 1983). [50] S. Chattopadhyay and P.L.M. Plummer, in: Physics and chemistry of finite systems, Vol. 1, eds. P. Jena, S.N. Khanna and B.K. Rao, NATO AS1 Series (Kluwer, Dordrecht, 1992) p. 211. [Sl] D.J. Nesbitt, in: The structure of small molecules and ions, eds. R. Naaman and Z. Vager (Plenum Press, New York, 1988). [52] J.A. Pople, M. Head-Gordon, D.J. Fox, K. Raghavachari and L.A. Curtiss, J. Chem. Phys. 90 (1989) 5622. [53] L.A. Curtiss, K. Raghavachari, G.W. Trucks and J.A. Pople, J. Chem. Phys. 94 ( 1991) 7221. [54] D.J. DeFrees, B.A. Levi, S.K. Pollack, W.J.Hehre, J.S. Binkley and J.A. Pople, J. Am. Chem. Sot. 101 (1979) 4085.