Structure and Energetics of the Weakly Bound NH3···H2O Complex

8528

J. Phys. Chem. A 1999, 103, 8528-8536

Structure and Energetics of the Weakly Bound NH3‚‚‚H2O Complex J. Sadlej* Department of Chemistry, UniVersity of Warsaw, Pasteura 1, 02-093 Warsaw, Poland, and Drug Institute, Chelmska 30/34, 00-725 Warsaw, Poland

R. Moszynski Department of Chemistry, UniVersity of Warsaw, Pasteura 1, 02-093 Warsaw, Poland, and Laboratoire de Chimie The´ orique, UMR 7551 CNRS/ULP, Institut Le Bel, UniVersite´ Louis Pasteur, 4 rue Blaise Pascal, F-67008 Strasbourg Cedex, France

J. Cz. Dobrowolski Drug Institute, Chelmska 30/34, 00-725 Warsaw, Poland, and Industrial Chemistry Research Institute, Rydygiera 8, 01-793 Warsaw, Poland

A. P. Mazurek Drug Institute, Chelmska 30/34, 00-725 Warsaw, Poland ReceiVed: May 14, 1999; In Final Form: June 29, 1999

Optimal structures, interaction energies, and harmonic vibrational frequencies of the NH3‚‚‚H2O complex have been determined from the second-order Møller-Plesset perturbation theory and DFT/B3PW91 calculations with the aug-cc-pVTZ basis set. Some properties of the complex have been calculated at the MP4 and CCSD(T) levels. The tunneling motion of the water molecule around its c inertial axis was studied, and the barrier to the exchange of the bound and free hydrogen atoms was determined as 1267 cm-1. The nature of the intermolecular interactions in the complex was investigated by symmetry-adapted perturbation theory (SAPT). As revealed by the SAPT analysis the main binding contributions are the electrostatic and induction components. The calculations of the vibrational frequencies and infrared intensities for this complex are presented to facilitate the frequency assignments of the experimental spectra.

I. Introduction The structures of numerous weakly bound complexes have been recently examined spectroscopically using molecular beam electric resonance, Fourier transform microwave spectroscopy, high-resolution vibrational spectroscopy, vibrational predissociation spectra, and matrix isolation infrared spectra.1-3 The challenge is to obtain both the structures and the potential functions. Indeed, a successful determination of the potential energy surfaces for van der Waals and hydrogen-bonded complexes has thus far been possible only for atom-molecule complexes (see, for instance, refs 4-6). Only recently have powerful dynamical techniques been developed for refinement of the potential surfaces of systems as large as the water dimer.7 A complementary approach to obtain structural and potential function information for weakly bound complexes is the ab initio calculation. These methods are becoming accepted tools for studies of the structure and properties of molecular complexes. They are especially useful in those cases in which the structures of different conformers are difficult to determine from the spectroscopic data. Although the accuracy of the ab initio calculations is still below the state-of-the-art accuracy of the spectroscopic data, the calculations provide information on the shape of the potential energy surfaces without any initial assumptions. Moreover, calculations can help to explain why the complex is more stable in a specific configuration.

Recently, the NH3‚‚‚H2O complex was examined by microwave and tunable far-infrared laser spectroscopy.8,9 Also, the spectroscopic constants were published, obtained from microwave and radio frequency spectra, observed by molecular beam electric resonance spectroscopy.10 The detailed microwave study by Stockman et al.8 provides the following structural features: the equilibrium structure of the complex has the heavy atoms collinear with the hydrogen of water bound to the nitrogen of the NH3 molecule. However, the hydrogen bond in N‚‚‚HO-H is nonlinear (the C3 axis of the ammonia molecule forms an angle of 23° with the N‚‚‚O complex axis). The bond length was determined to be 2.983 Å. Both the nearly free internal rotation of NH3 and the tunneling of the H2O protons were observed. These authors also estimated the barrier to NH3 rotation about the C3 axis as 10.5 ( 5.0 cm-1 and the water tunneling splitting as 700 cm-1 (in a one-dimensional model). The latter is a measure of the barrier to the H2O proton exchange. If one assumes that only the water molecule is moving in the tunneling process, the barrier height was estimated as 840 cm-1. Finally, the ground-state dipole moment obtained from the Stark effect measurements and the nitrogen quadrupole coupling constants were also published.10 Another group of papers is devoted to the matrix isolation studies of NH3‚‚‚H2O in the Ar, Ne, and Kr matrices. Nelander et al.l1,l2 and Yeo and Ford13 registered vibrational frequencies

10.1021/jp991581l CCC: $18.00 © 1999 American Chemical Society Published on Web 09/30/1999

Structure and Energetics of NH3‚‚‚H2O of the ammonia and water intramolecular modes and several intermolecular modes of the complex. Previous ab initio potential energy surfaces of NH3‚‚‚H2O were reported in a few papers.14-19 Del Bene found (HF/STO basis set) that the ammonia molecule is free to rotate about its own C3 axis in the complex.14 Yeo and Ford (MP2/6-31G**)16 and Latajka and Scheiner (MP2/6-31G**)18 characterized the trans structure of the complex as a minimum, whereas cis and bifurcated structures are transition states with one imaginary frequency.16 By contrast, Skurski and Gutowski found the cis form for the minimum and the trans form as a transition state.19 Accurate prediction of the dissociation energy remains an open problem, because the experimental dissociation energy is not known. Inclusion of the correction for the basis-set superposition error (BSSE) together with the zero-point vibrational energy calculated in a good basis set is one reason for the present paper investigations. The second controversy in the literature concerns the infrared spectra. Nelander et al.l1,12 measured a band at 19.5 cm-1 and interpreted it as a torsion of the linearly (N‚‚‚H-O) hydrogen bonded water molecule about the a-axis of the dimer. However, Stockman et al.,8 on the basis of the far-infrared measurements, assigned this band to “a superposition of some or all the observed rotation-tunneling bands, which fall between 20 and 22 cm-1”. The authors of ref 8 also stressed the important discrepancies between the ab initio calculations18 and the matrix isolation work.11 For instance, the out-of-plane ammonia bend is assigned at roughly twice the calculated ab initio value.13,18 The study of molecular complexes including NH3 is still a challenging problem for theory. The goal of this paper is fourfold: (i) to compute the characteristic points on the potential energy surface of the NH3‚‚‚H2O complex and to study the effect of the rotation of the water molecule on the shape of this surface; (ii) to achieve a better understanding of the nature of the interaction in the complex by analyzing the components such as the electrostatic, exchange, induction, and dispersion computed by symmetry-adapted perturbation theory (SAPT)30,31 [This method provides information about the interplay of the four major physical contributions governing the bonding in van der Waals and hydrogen-bonded complexes. Our analysis will address the question of the origin of the differences in the relative stabilities of the three isomers: trans, cis and bifurcated.]; (iii) to check the applicability of the density functional theory (DFT) to predict structure and energetics of various conformers of the same complex [Indeed, it is known that DFT often offers an electron correlation comparable to MP2. However, when calculations concern different structures of a hydrogen-bonded complex, the accuracy of the DFT results may not be the same for different forms of the complex. Thus, a comparison of the DFT and SAPT results will address the question of the reliability of the DFT method for structures governed by different components of the interaction energy.]; (iv) to resolve, despite quite a few computational16,18,19 and experimental investigations,9,11,12 some important discrepances between theory and experiment concerning the vibrational frequencies. This paper is organized as follows. In Section II we briefly describe the methods of calculations used in this work. In Section III we present and discuss the results of our calculations. Finally, Section IV concludes our paper. II. Methods of Calculations A. Supermolecular Calculations. The supermolecular calculations have been done with the Gaussian 9420 program. Full

J. Phys. Chem. A, Vol. 103, No. 42, 1999 8529 geometry optimizations have been carried out using both the second-order Møller-Plesset perturbation theory (MP2) and the density functional approach (DFT/B3PW91). B3PW91 denotes the three-parameter functional of Becke21 with the semilocal correlation potential provided by the Perdew-Wang expression.22 This method gives accurate harmonic frequencies for molecules as well as for hydrogen-bonding systems. The geometry optimizations have been performed for the isolated subunits and for the entire complex at the MP2 and DFT/ B3PW91 levels of theory without taking into account the basis set superposition error in the optimization procedure. The vibrational frequencies have been computed within the harmonic approximation at the MP2 level. It is well-known that the geometries of the monomers in the complex differ slightly from their optimal geometries when isolated. Therefore, ab initio calculations of the interaction energies, binding energies, and dissociation energies require some care. The interaction energy of the complex, ESM int , has been obtained by subtraction of the energies of complex and of the monomers (the superscript “SM” denotes the method used in the actual calculation). Here, the geometries of the monomers correspond to those in the complex and, thus, are deformed from their equilibrium values. In the calculations of the interaction energies we employed the full basis of the dimer; that is, we corrected the computed interaction energies for the basis-set superposition error using the prescription of Boys and Bernardi.23 This is a rather commonly accepted procedure to obtain reliable interaction energies,24,25 despite some controversial arguments reported in the literature.26 To compute the binding energies DSM e , the interaction energies have to be corrected for the so-called deformation energies Edef,27-29 that is, for the energy needed to deform the monomers from their optimal equilibrium geometries to their geometries in the complex. The latter quantity was computed on the basis of the monomer. (The meaning of the “SM” superscript is the same as for the were interaction energy.) Finally, dissociation energies DSM 0 obtained from the binding energies by adding the correction for the zero-point vibrational motion. In our work, the interaction CCSD(T) , and the respective energies of the complex, EMP4 int and Eint MP4 values of the binding energies, De and DCCSD(T) , and the e MP4 CCSD(T) dissociation energies, D0 and D0 , were calculated at MP2 optimized geometries of the complexes. The zero-point correction (ZPE) was calculated at the MP2 level. The dipole moment and nuclear quadrupole coupling constants at the HF and MP2 levels were calculated as analytical derivatives of the corresponding energies. The MP4 and CCSD(T) dipole moments were obtained by using the finite field technique (the field strength was set equal to 0.001 au). B. Perturbation Theory Calculations. To get more insight into the physical nature of the NH3‚‚‚H2O interaction, we performed SAPT calculations for the three conformers: trans, cis, and bifurcated. The SAPT approach employed in the present paper is the same as in previous papers (see refs 30-33 and, for the most recent applications, refs 34 and 35). The intermowas represented by the sum of lecular interaction energy ESAPT int the first- and second-order polarization and exchange contributions36 (1) (2) (2) (2) ESAPT ) E(1) int pol + Eexch + Eind + Edisp + Eexch-ind + (2) (2) + Eexch-disp (1) Eexch-def

where E(1) pol is the electrostatic energy, fully accounting for the overlap (penetration) of the monomer charge distributions, E(2) ind

8530 J. Phys. Chem. A, Vol. 103, No. 42, 1999 and E(2) disp denote the nonexpanded induction (with response) (2) (2) (2) and dispersion energies, and Eexch-ind , Eexch-def , and Eexch-disp are the second-order exchange-induction (with response), exchange-deformation, and exchange-dispersion energies. The contributions appearing on the right-hand side of eq 1 were evaluated using the many-body techniques developed in refs 37-44 (see also ref 36 for a review). The exchangedeformation term is defined as in refs 45 and 46, was computed from the supermolecule Hartree-Fock interaction energy, and may be viewed as part of the exchange energy due to the exchange terms in the perturbation theory equations for the wave functions.46 The computational scheme adopted in this work was the same as in refs 32, 34, and 35. C. Computational Details. There are number of considerations in the selection of a basis set suitable to describe intermolecular interactions. First, it is necessary to reproduce various electric properties of the monomers as accurately as possible. This will ensure a correct description of the long-range interactions in the complex. Second, the sensitivity of these electric properties to the approach of the partner’s orbitals should be minimized to diminish the basis set extension effects. We decided to use the aug-cc-pVTZ atomic basis.47 As shown by Xantheas et al.,28,29,48 this basis set accurately reproduces the geometries, frequencies, and electric properties of isolated molecules and their complexes. For computational reasons the SAPT calculations have been performed on a smaller aug-ccpVDZ basis, supplemented with the 1f function on the nitrogen and oxygen (with exponents 0.3640 and 0.5, respectively) and the 1d function on the hydrogen (with exponent 0.2470). For the trans conformer the MP4 interaction energy in this basis amounts to -6.131 kcal/mol, and it differs by only 3% from the value computed in the full aug-cc-pVTZ basis. Because in the present work the SAPT method was not used to determine the potential energy surface and the characteristic points on this surface, but just to provide a qualitative view of the nature of the interaction, we decided that the use of the smaller basis is legitimate, the 3% disagreement with the large basis set results being of no importance for our qualitative discussion. We considered the three most important conformers of the complex NH3‚‚‚H2O, in which the ammonia molecule is the electron donor: trans, cis, and bifurcated (see Figure 1a-c). Additionally, a bifurcated-type structure, denoted bifurcatedperpendicular (see Figure 1d), was investigated as a transition configuration for the water molecule rotation. For the H2NH‚‚‚OH2 complex, it was not possible to locate a minimum for a structure with the NH3 molecule serving as the hydrogen bond donor and water as the hydrogen bond acceptor. This failure is not surprising given the fact that the ammonia is a stronger hydrogen bond acceptor but a weaker donor, whereas water is a stronger donor but a weaker acceptor.16 III. Results and Discussion A. Geometry and Energetics of the Complex. Ab initio calculations at the MP2/aug-cc-VTZ and B3PW91/aug-cc-VTZ levels predicted a global minimum corresponding to the trans conformer (Figure 1a). The cis conformer (Figure 1b) is a saddle point (one imaginary frequency). This result was obtained when the optimization procedure was run until the maximum force and the maximum displacement were converged to 10-6 and 10-5, respectively. However, when the optimization process was stopped in a step with the maximum force and displacement equal to 10-3 and 1.6 × 10-3, respectively, all frequencies were positive. One large negative eigenvalue calculated for the

Sadlej et al.

Figure 1. Geometries of the H3N‚‚‚HOH complex corresponding to the (a) trans, (b) cis, (c) bifurcated, and (d) bifurcated-perpendicular conformers.

bifurcated structure (Figure 1c) at the MP2 level indicated that this structure also corresponds to a saddle point. In Table 1 the fully optimized geometries of the H3N‚‚‚HOH complex at the MP2 and DFT levels are reported for the trans, cis, and bifurcated conformers. For the cis conformer we present the results obtained from the converged optimization (i.e., for which one frequency is imaginary). The microwave experiment predicts a structure with an O‚‚‚N distance equal to 2.983 Å and with the ammonia C3 symmetry axis forming an angle of 23.1° with the O‚‚‚N axis. The theoretical MP2 values for the distance and angle quoted above are 2.922 Å and 22°, respectively, in good agreement with the experimental results. A small disagreement of 0.06 Å in the O‚‚‚N distance may be due to the ground-state vibrational averaging. The trans and cis conformers differ only slightly by the values of the bond lengths and the bond angles. The two angles that define the orientations of the ammonia molecule with respect to the N‚‚‚H hydrogen bond are quite different, though. The geometries of the water and ammonia subunits almost do not change upon complexation. Both sets of geometries (MP2 and DFT/B3PW91) for trans and cis structures differ only slightly.

Structure and Energetics of NH3‚‚‚H2O

J. Phys. Chem. A, Vol. 103, No. 42, 1999 8531

TABLE 1: Geometrical Parameters (in Angstroms and Degrees) of the H3N‚‚‚HOH Complex Calculated at the MP2 and B3PW91 Levels for the Trans, Cis, Bifurcated, and Bifurcated-Perpendicular Conformersa trans MP2 r (N1‚‚‚O5) r (N1H2) r (N1H3) r (O5H7) r (O5H6) R (H7O5H6) R (H2N1‚‚‚O5) R (H3N1‚‚‚O6) R (H3N1H4) R (H2N1H3) R (N1‚‚‚O5H6) R (N1‚‚‚O5H7) δ (H2N1‚‚‚O5H6) δ (H2N1‚‚‚O5H7) δ (H3N1‚‚‚O5H6) δ (H3N1‚‚‚O5H7) a

cis B3PW91

2.9219 1.0130 1.0126 0.9746 0.9605

2.9181 1.0131 1.0127 0.977 0.959

MP2 2.9230 1.0123 1.0127 0.9747 0.9603

104.8 99.4 117.7 107.0 106.9 111.5 6.8

105.4 101.2 116.5 107.3 107.2 111.4 6.7

104.7 124.3 105.4 106.8 106.9 6.7 111.4

-180.0 -180.0 65.2 65.2

-179.9 -179.9 64.1 64.1

0.0 0.0 -123.6 -123.6

bifurcated B3PW91 2.920 1.013 1.013 0.976 0.959

MP2 3.2411 1.0178 1.0129 0.9624 0.9624

bifurcated-perpendicular

B3PW91 3.3120 1.0131 1.0132 0.9609 0.9609

MP2 3.2430 1.0134 1.0128 0.9630 0.9623

B3PW91 3.3119 1.0129 1.0130 0.9619 0.9610

105.3 121.9 106.2 107.2 107.3 6.5 111.3

100.6 110.5 113.3 106.3 106.5 50.4 50.3

101.7 110.5 112.8 105.7 106.8 50.9 50.8

100.7 113.1 112.0 105.2 106.3 50.3 50.3

101.2 112.9 111.9 105.1 106.2 50.6 50.6

-0.03 -0.03 -123.1 -123.1

-89.0 91.0 151.6 -28.4

-89.0 90.1 151.5 -28.5

0.0 180.0 120.1 -59.8

0.0 180.0 120.1 -59.9

Global energy minimum corresponds to the trans conformation, whereas the cis and bifurcated structures correspond to saddle points.

TABLE 2: Interaction Energies, Binding Energies, and Dissociation Energies (in Kilocalories per Mole) for the Trans, Cis, and Bifurcated Conformers of the NH3‚‚‚H2O Complex Calculated by Various Methods trans

cis

bifurcated

EHF int EB3PW91 int EMP2 int EMP4 int ECCSD(T) int

-4.468 -5.839 -6.365 -6.323 -6.285

-4.472 -5.847 -6.361 -6.320 -6.282

-2.410 -2.276 -3.095 -3.179 -3.219

DHF e DB3PW91 e DMP2 e DMP4 e DCCSD(T) e

-4.037 -5.664 -6.240 -6.219 -6.168

-4.091 -5.677 -6.238 -6.216 -6.165

-2.353 -1.767 -3.017 -3.101 -3.161

DHF 0 DB3PW91 0 DMP2 0 DMP4 0 DCCSD(T) 0

1.787 3.390 3.990 3.969 3.918

As one can see from Table 2, which presents the energetics of the complexes, the trans conformer is the global minimum. The fact that the cis structure differs only marginally from the trans supports the conclusion of the earlier theoretical and experimental findings that the ammonia molecule in the complex is almost free to rotate about its C3 axis. Inclusion of the correlation effects at the MP2 level is important and increases the interaction energy from --4.468 to --6.323 kcal/mol for trans form, that is, by ∼40%. The MP4 and CCSD(T) calculations change slightly this result. For all structures considered in Table 2 the convergence of the Møller-Plesset expansion for the interaction energy appears to be good, the MP2 approximation providing over 98, 98, and 96% of the CCSD(T) result for the trans, cis, and bifurcated geometries, respectively. One may note that the performance of the B3PW91 functional is rather good for the trans and cis conformers and very bad for the bifurcated structure. Indeed, the interaction energies for the trans and cis geometries computed by B3PW91/ DFT differ by ≈7% from the CCSD(T) results. By contrast, for the bifurcated structure this error is as large as 30%. See the discussion of the SAPT results for an explanation of this strange performance of the B3PW91/DFT method. The defor-

mation energies of the monomers do not introduce substantial changes in the binding energies. This is consistent with the marginal geometry changes of the monomers upon complexation. Trends similar to those observed for the interaction energy remain for the dissociation energy. However, it is worth noting that, similar to the Hartree-Fock method, the B3PW91 functional predicts the cis conformer to be more stable than trans. This is in contrast to the MP2, MP4, and CCSD(T) results. The experimental dissociation energy of this complex is not known. Extensive calculations of the binding energies of van der Waals complexes indicate that our value in the aug-cc-VTZ basis set is underestimated, mainly due to an underestimation of the dispersion component of the interaction energy, by ∼1015%. To provide the best estimate of the binding and dissociation energies, we recomputed the interaction energy at the MP2 level in a larger aug-cc-pVQZ basis and corrected it for higher order correlation effects, as well as deformation energies and zero-point corrections, computed in the smaller aug-cc-pVTZ. This gave -6.399 kcal/mol for De and 4.149 kcal/mol for D0. We believe that the latter values represent the best estimates of the binding and dissociation energies of the NH3‚‚‚H2O complex. B. Tunneling Motion of the Molecules in the Complex. Let us discuss now the internal motion of the molecules in the complex. First we discuss the tunneling motion of the water molecule. The H2O tunneling motion is strongly coupled to the ammonia internal rotation. Experimentally, the height of the barrier (V3) was obtained by fitting the parameter V3 in the potential term of the nuclear motion Hamiltonian to the observed microwave spectra.8 The tunneling pathway corresponds to a partial rotation of the H2O unit about an axis normal to the H2O molecular plane. However, in contrast to the nearly free rotation of the ammonia molecule about its C3 axis, the barrier for the water tunneling motion about its c-axis in the a-b plane of the complex, leading to the exchange of the hydrogen atoms, is estimated as 704 cm-1 (assuming a simple one-dimensional model).8 When the motion of the water subunit only was considered in the tunneling process, the barrier height was estimated as 840 cm-1.8 Following this analysis we examined the tunneling motion of water around the c-inertial axis (this motion corresponds to the change of the N1‚‚‚O5H6 angle in Figure 1a). The results are presented in Figure 2. Taking into account the minimum energy path, which involves the optimization of all the geo-

8532 J. Phys. Chem. A, Vol. 103, No. 42, 1999

Figure 2. Energy profile for the tunneling process of the water molecule calculated at the MP2 level. The cis conformer is on the left, the trans on the right, and the bifurcated-perpendicular in the middle of the x-axis.

metrical parameters for all the points except for the N1‚‚‚O5H6 angle, we found a barrier of 1267 cm-1. This height of the barrier is considerably higher than the experimental estimates.8 One should keep in mind, however, that these estimates neglect the coupling between various intermolecular degrees of freedom in the complex. Given the fact that the recorded far-infrared spectra reveal several lowenergy vibrational modes corresponding to large-amplitude motions in the complex, we believe that such a decoupling is not fully justified, and our ab initio result may be more realistic. It is possible that further lowering of the rotation barrier can be obtained by introducing more polarization functions to the basis set. We also investigated the energetics of the NH3 rotation about the C3 axis of ammonia. This rotation is described by the dihedral angle d ) (H6O5‚‚‚N1H2) (cf. Figure 1). This angle is equal to -180° and 0° for the trans and cis forms, respectively. The MP2 energy difference trans - cis is equal to 0.41 cm-1 only. Taking into account the minimum energy path connecting the cis and trans forms (which involves fully optimized structures except for the dihedral angle d), we found a barrier of 2.2 cm-1 for d ) 24°. The experimental estimate of this barrier height is 10.5 ( 5.0 cm-1.8 Again, the agreement between theory and experiment is not perfect. However, this barrier is very difficult to determine, both theoretically and experimentally, due to the extreme flatness of the potential energy surface in this region. C. SAPT Calculations: The Nature of the Interaction. The SAPT calculations provide us with a decomposition of the interaction energy into physically different contributions (see Table 3). The total SAPT interaction energy amounts to -6.222 kcal/mol for the trans form, which compares well with the value of -6.131 kcal/mol computed at the MP4 level. As one can see, the main binding contribution comes from the electrostatic, induction, and dispersion energies, with the electrostatic component being the most important. These attractive energies have to be compensated mainly by the first-order exchange interaction (higher than the electrostatic component) and also by the exchange-induction interaction. The exchange-dispersion effect is much smaller. The decomposition of the interaction energy for the cis conformer is very similar to the decomposition for the trans form. To obtain an estimate of which components contribute to the water tunneling barrier, some additional SAPT calcula-

Sadlej et al. tions have been performed for the bifurcated-perpendicular saddle point. As one can see from Table 3 (last row), the total interaction energy for this structure is smaller by almost twice as compared to the trans and cis structures. The main contribution to the barrier comes from the electrostatic effect, followed by the exchange, dispersion, and induction effects. The exchange repulsion decreases considerably (∼3 times) in the bifurcatedperpendicular form as compared to the trans and cis ones. However, the bonding in the bifurcated form is largely due to the dispersion energy and much less to the induction and electrostatic. Summarizing, symmetry-adapted perturbation theory indicates that the hydrogen bonding in the trans complex is mainly caused by strong electrostatic and induction interactions quenched by the first-order exchange-repulsion interactions and the secondorder exchange-induction effects. The stability of the bifurcated form is largely due to the dispersion attraction and a much smaller exchange-repulsion effect. It is interesting to note that the B3PW91 functional works rather well for those configurations which are mainly governed by the electrostatic and induction forces. However, when the dispersion interaction starts to play an important role, the performance of this functional deteriorates drastically. These results clearly show that the DFT methods can be used with trust only for those configurations of the complexes for which the classical electrostatic and induction effects are predominant. These findings are in agreement with the results of a theoretical analysis reported in ref 49. D. Frequency Analysis. The internal coordinates used in the normal-mode analysis for potential energy distribution (PED) are presented in Table 4. The vibrational harmonic frequencies and IR intensities of the three conformers of H3N‚‚‚HOH were calculated at the MP2 and B3PW91 levels for the respective optimized geometries (given in Table 1). The results are presented in Table 5. The PED analysis is presented in Table 5 for the trans form only. A single imaginary frequency is predicted for the bifurcated structure (Figure 1c) at the MP2 level. Single and double imaginary frequencies are calculated for the cis and bifurcated structures at the DFT/B3PW91 level, respectively. Apart from one imaginary frequency, the calculated frequencies and intensities are very similar for both the trans and cis structures. Table 6 lists changes induced in the subunit parameters by the formation of the H-bonded complex. In this table are reported the harmonic frequencies for the trans conformer and the ratio Acom/Amon, that is, the effect of the complexation on the intensities. The computed vibrational frequencies are in most cases closer to the experimental values than any of the three other calculations known in the literature.16,18,19 The main reason for the improvement is the better quality basis set we used. Let us first discuss the intramolecular modes of the H3N‚‚‚ HOH complex (Tables 5 and 6). The vibrational modes of the monomers remain localized in ammonia and water molecules. The monomer H2O frequencies in the aug-cc-pVTZ basis set are 1629, 3821, and 3947 cm-1. Formation of the complex with NH3 destroyed the C2V symmetry of HOH. The two stretching modes located at 3947 and 3821 cm-1 for monomer HOH occur at 3905 (free OH) and 3586 (H-bonded OH) cm-1; that is, they are shifted by -42 and -235 cm-1. The calculated large red shift of -235 cm-1 is in agreement with the shift of -223 cm-1 observed in the Ar matrix.11,12 The strong 100-fold intensity enhancement agrees with the results for similar systems and agrees with the earlier calculations for this system.16,17 The HOH bending frequency of 1678 cm-1 is shifted toward higher

Structure and Energetics of NH3‚‚‚H2O

J. Phys. Chem. A, Vol. 103, No. 42, 1999 8533

TABLE 3: Components of the Interaction Energy (in Kilocalories per Mole) of the H3N‚‚‚HOH Complex Computed by SAPT for the Trans, Cis, and Two Bifurcated Structures trans cis bifurcated bifurcated-perpendicular

(1) Eelst

E(1) exch

E(2) ind

(2) Eexch-ind

E(2) disp

(2) Eexch-disp

(2) Eexch-def

ESAPT int

-11.698 -11.677 -4.504 -4.494

13.036 12.994 3.501 3.481

-5.931 -5.909 -1.543 -1.535

3.381 3.367 0.993 0.986

-4.127 -4.117 -1.887 -1.880

0.758 0.756 0.265 0.263

-1.643 -1.634 -0.132 -0.131

-6.222 -6.219 -3.307 -3.309

TABLE 4: Internal Coordinates Used in the Normal Modes Analysis (Atom Numbering as in Figure 1) for the Trans Configuration of the H3N‚‚‚HOH Complex coordinate

definition

explanation

notation

S1 S2 S3 S4 S5 S6 S7 S8 S9 S10 S11 S12 S13 S14 S15

r5,6 r1,4 + r1,2 - r1,3 r1,2 - r1,3 - r1,4 r5,7 r1,2 + r1,3 + r1,4 β7,5,6 - β3,1,4 + β2,1,4 + β2,1,3 β2,1,4 - β3,1,4 - β2,1,3 β7,5,6 + β3,1,4 - β2,1,4 - β2,1,3 τ4,1,2,3 τ2,1,7,5 β2,1,7 - β1,7,5 r1,7 τ1,5,2,4 β2,1,7 + β1,7,5 τ1,7,5,6 + τ3,1,7,5

ν (O5H6) ν (N1H4) + ν (N1H2) - ν (N1H3) ν (N1H2) - ν (N1H3) - ν (N1H4) ν (O5H7) ν (N1H2) + ν (N1H3) + ν (N1H4) β (H7O5H6) -β (H3N1H4) + β (H2N1H4) + β (H2N1H3) β (H2N1H4) - β (H3N1H4) - β (H2N1H3) β (H7O5H6) + β (H3N1H4) - β (H2N1H4) - β (H2N1H3) τ (H4N1H2H3) τ (H2N1H7O5) β (H2N1H7) - β (N1H7O5) ν (N1H7) τ (N1O5H2H4) β (H2N1H7) + β (N1H7O5) τ (N1H7O5H6) + τ (H3N1H7O5)

νfree (OH) νas,1 (NH3) νas,2 (NH3) νHB (OH) νs (NH3) β1 (H2O-NH2) βa (NH3) β2 (H2O-NH2) τ (Ninv) τ1 (HB) β1 (HB) ν (HB) τ2 (HB) β2 (HB) τ3 (HB)

TABLE 5: Calculated Frequencies (cm-1), IR Intensities (km/mol), and ZPE (kcal/mol) for the H3N‚‚‚HOH Complex at the MP2 and B3PW91 Levels for the Trans, Cis, and Bifurcated Conformers trans MP2 no. 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 ZPE

ν

cis B3BPW91

I

15 66.1 169 29.7 179 43.2 201 26.7 460 82.2 716 86.4 1094 146.8 1653 55.5 1655 16.6 1678 1.7 3499 2.1 3588 570.3 3642 17.2 3646 15.8 3907 85.7 37.341

ν

I

51 68.7 182 10.5 194 38.0 205 40.4 473 80.6 730 87.0 1083 144.5 1649 54.3 1660 18.3 1673 7.4 3489 4.4 3536 659.4 3611 12.1 3614 11.0 3901 57.8 37.245

MP2 ν

I

-8 168 30.2 179 42.0 201 25.1 456 82.1 718 88.4 1093 147.3 1660 55.6 1667 17.5 1673 0.5 3499 1.4 3588 568.0 3643 17.3 3646 15.8 3907 85.7 37.315

frequency from the monomer value of 1629 cm-1 without the intensity change. In contrast to the water molecule, the complex modes associated with the intramolecular modes of ammonia are changed very little from their monomer values (Tables 5 and 6). The insensitivity of the NH3 frequnecies to the complexation is entirely consistent with the mentioned lack of appreciable changes in its geometry. The monomer values are 1035, 1668 (degen), 3505, and 3652 (degen) cm-1. Complexation of NH3 with HOH breaks the C3V symmetry, causing a splitting of the degenerate vibrational modes of the NH3 asymmetric stretching and bending. Asymmetric stretching is affected very little by the formation of the complex. One of the stretching modes is lowered by 4 cm-1 while the other is lowered by 9 cm-1, although their intensities are doubled. The degenerate bending mode 1668 cm-1 in the monomer is split into two compontents with frequencies of 1669 and 1653 cm-1, respectively. The lower component exhibits substantial intensity increase (4-fold) relative to the corresponding monomer mode. The lowest monomer mode, 1035 cm-1, is shifted to higher frequency by

bifurcated B3BPW91 ν

I

-16 179 16.5 191 37.3 201 33.1 467 80.6 730 89.1 1081 145.2 1655 51.4 1657 18.9 1669 9.3 3489 3.9 3538 656.7 3612 12.4 3614 10.8 3899 57.9 37.315

MP2 ν

B3BPW91 I

-313 4 0.0 97 3.0 112 0.8 174 22.2 438 221.3 1076 141.3 1631 151.6 1669 16.0 1670 15.8 3497 1.3 3639 11.9 3641 12.0 3831 18.1 3924 56.5 36.318

ν

I

-29 -49 94 0.6 105 2.8 179 22.9 418 228.3 1064 138.4 1636 159.6 1661 17.3 1602 17.3 3486 1.7 3607 7.2 3608 7.4 3835 16.8 3911 46.4 36.123

58 cm-1 in the complex and is described by us as a mainly torsional ammonia mode. This is in good agreement with the matrix isolation result (76 cm-1)13 and the tunable microwave laser studies (71 cm-1).9 Because of the difficulty in the detection of bands in the farIR region, the calculations can play a particularly important role in predictions and interpretations of this part of spectrum. In Tables 5 and 6 we report the six intramolecular frequencies. In some cases, the highest frequency vibration is predicted to be the in-plane wagging of the proton donor. Our PED analysis classified the 716 cm-1 band as out-of-plane wagging of the proton donor. The second frequency, 462 cm-1, is characterized as in-plane wagging of HOH. The calculated values of 716 and 462 cm-1, both almost of the same intensity, are close to the experimental values in the Ar matrix: 600 and 440 cm-1, respectively. The intermolecular stretching band is generally of quite low frequency and intensity, making its observation rather difficult. However, for this complex particularly good agreement is noted for the stretching N‚‚‚H-O mode: 201 cm-1 (calcd) and 202 cm-1 (exptl). The two next modes with relatively low

8534 J. Phys. Chem. A, Vol. 103, No. 42, 1999

Sadlej et al.

TABLE 6: Comparison of the Infrared Frequencies (cm-1) and Intensities (km/mol) for the H3N‚‚‚HOH Complex Calculated at the MP2 Level for the Trans Conformer with the Theoretical and Experimental Results of Yeo and Ford Intramolecular this work ν

I

3905 3648 3643 3586 3500 1678 1669 1653 1093

86 16 17 570 2 2 17 55 147

this work PED% free

96ν (OH) 99νas,1 (NH3) 98νas,2 (NH3) 94νHB (OH) 99νs (NH3) 80β1 (H2O - NH2) 98βa (NH3) 83β2 (H2O + NH2) 86τ (Ninv)

∆νa -42 -4 -9 -235 -5 49 1 -15 58

Yeo-Ford

exp

I(com)/I(mon)

∆νa

I(com)/I(mon)

νa

0.6 2.0 2.0 104.0 0.7 0.03 1.2 4.0 1.0

-45 0 2 -167 3 60 -4 -13 36

1.5 45.0 17.1 87.6 17.4 37.0 1.3 3.7 0.9

-31 -223 -1 29 76

Intermolecular this work

a∆ν

Yeo-Ford

exp

ν

I

PED%

ν

I

∆νa

716 462 201 178 169 54

86 83 26 47 30 62

91τ 1 (HB) 87β 1 (HB) 75ν (HB), 21β 2 (HB) 83τ 2 (HB) 77β 2 (HB), 22ν (HB) 92τ 3 (HB)

712 460 202 172 161 30

184 110 19 44 40 77

600 440 202 180 10

) νcomplex - νmonomer.

intensities, at 178 and 169 cm -1, one of which is found in the experiment at 180 cm-1, are associated mainly with the wagging of the ammonia molecule (in- and out-of-plane). The predicted intensities of these two modes are nearly equal to one another. The smallest frequency of 54 cm-1, which we assigned as the torsional motion about the hydrogen bond, is a relatively intense band. Summarizing, the agreement between our calculated harmonic frequencies and the experimental measurements is indeed encouraging, despite the flatness of the intermolecular potential and the anharmonicity of the low-frequency modes. The calculated shifts are general in agreement with those calculated by Yeo and Ford16 with some exceptions, though. For instance, the blue shift of 58 cm-1 for the inversion mode of ammonia (no. 7 in Table 5) is higher than that of ref 16 by 36 cm-1. This point was one of the disagreements between the IR measurements11,12 and the calculations mentioned by Stockman et al.8 By contrast, the calculated shifts for the water molecule are like those calculated by Yeo and Ford16 except for the intensity change of the HOH bending mode, which is significantly lower in the complex than in the monomer. However, our work does not provide any direct confirmation of the conjecture8 that the band at 20 cm-1 measured by Engdahl and Nelander11,12 can be ascribed as a superposition of some tunneling bands.8 E. Molecular Properties. In Table 7 we report the values of the total dipole moment for the trans, cis, and bifurcated geometries of the complex, computed at various levels of the theory. An inspection of this table shows that the correlation effects are not very important. The Hartree-Fock value overestimates the CCSD(T) result by 4%, whereas the error of the MP2 method is only 0.6%. The convergence of the MøllerPlesset expansion appears to be satisfactory. The B3PW91/DFT result is very close to the Hartree-Fock value. Similar trends are observed for the component of the dipole moment along the a-axis of the dimer. Our best value of the total dipole moment agrees fairly well with the experimental results of refs 8 and 10, the error of the CCSD(T) result with respect to the experimental data of refs 8 and 10 being 8 and 14%, respec-

TABLE 7: Dipole Moment µtot, Its Component along the a Axis µa (in Debye), and the Nuclear Quadrupole Coupling Constant at the 15N Nucleus χaa (15N) (in Megahertz) for the Trans, Cis, and Bifurcated Conformers of the H3N‚‚‚HOH Complex Calculated at Various Levels of the Theory µHF tot µB3PW91 tot µMP2 tot µMP4 tot µCCSD(T) tot µexp tot , ref µexp tot , ref

8 10

trans

cis

bifurcated

3.600 3.581 3.478 3.444 3.455 3.022 3.202

3.605 3.580 3.483

4.047 3.834 3.862

µHF a µB3PW91 a µMP2 a µMP4 a µCCSD(T) a µexp a , ref 8 µexp a , ref 10

-3.367 -3.317 -3.270 -3.242 -3.249 2.9763 2.9766

-3.605 -3.316 -3.483

-4.047 -3.834 -3.862

15 χHF aa ( N) 15 χMP2 aa ( N) exp 15 χaa ( N)

-4.477 -3.726 -3.143

-4.472 -3.723

-4.581 -3.831

tively. Both papers8,10 also report the component of the dipole moment along the a-axis. Here, the experimental results agree one with the other to four significant figures. The level of agreement between theory and the two experiments in this case is about the same as for the total dipole moment, the theoretical result being 9% off the experimental value. The error in the computed dipole moment may be due to the basis set unsaturation or the neglect of the vibrational averaging effects. To check the first conjecture, we recomputed the total dipole moment and its component along the a-axis in the aug-cc-pVQZ basis. We assumed that the convergence of the Møller-Plesset expansion is basis-set independent, so that the MP2 result should be close to the converged CCSD(T) result in this basis. Our calculation in the aug-cc-pVQZ basis gave a total dipole moment of 3.486 D and its a component of 3.278 D. Comparison with the values computed in the aug-cc-pVTZ basis shows that both

Structure and Energetics of NH3‚‚‚H2O quantities are rather basis-set independent, so the theoretical results should be realistic. Thus, we tend to believe that the disagreement between theory and experiment is due to the vibrational average included in the experimental data. Herbine and Dyke10 also measured the nuclear quadrupole coupling constants for the 15N nuclei. Table 7 displays this parameter calculated at the Hartree-Fock and MP2 levels. Similar to the case of the dipole moment, the agreement between theory and experiment is not perfect. Morevoer, the theoretical value in the aug-cc-pVQZ basis (3.835 MHz) is rather close to the result computed in the aug-cc-pVTZ basis. Again, the disagreement between theory and experiment is probably caused by the vibrational average. One could argue that higher correlation effects could be important. However, similarly as the dipole moment, the nuclear quadrupole coupling constants are first-order properties and are essentially defined by the electron density. The good convergence of the Møller-Plesset expansion for the dipole moment suggests that the electron density of NH3‚‚‚H2O is not very sensitive to the correlation beyond the MP2 level, so the MP2 values of the nuclear quadrupole coupling constants should be rather accurate.

J. Phys. Chem. A, Vol. 103, No. 42, 1999 8535 largely due to the dispersion attraction and a much smaller exchange-repulsion effect. 6. Density functional methods give a rather accurate description of the intermolecular interactions for those structures of the complex that are dominated by the classical electrostatic and induction effects (trans and cis conformers). For the bifurcated structures, for which the dispersion contribution is large, the DFT/B3PW91 calculations are not appropriate. Thus, care should be taken when the DFT methods are applied in the calculations of the potential energy surfaces of hydrogen-bonded complexes. 7. The calculated dipole moment and nuclear quadrupole coupling constant agree fairly well with the experimental data.8,10 The observed disagreement is probably due to neglect of the vibrational average in our calculations. Acknowledgment. This work was supported by a KBN grant through the University of Warsaw (Grant BW-1418/10/98). J.S. thanks the Interdisciplinary Center for Computational Modeling, University of Warsaw, for providing computer time. References and Notes

IV. Conclusions In the present study we investigated the structural, energetic, and spectroscopic properties of the H3N‚‚‚H2O complex. Our results can be summarized as follows: 1. The equilibrium geometry of the complex corresponds to the trans structure with a nonlinear hydrogen bond formed by the nitrogen atom of ammonia and a hydrogen atom of water. The computed N‚‚‚O distance and the angle between the C3 axis of ammonia and the N‚‚‚O axis are in good agreement with the experimental data from the microwave measurements.8 Our best estimates for the binding and dissociation energies are -6.399 and 4.149 kcal/mol, respectively. 2. Investigation of the internal rotation dynamics in the complex revealed an almost free rotation of the ammonia molecule about its own C3 axis, in agreement with the experimental findings.8 Our best estimate for the barrier separating the two equivalent minima is 2.2 cm-1. 3. Similar study of the tunneling motion of the water molecule in the a-b plane of the dimer allowed us to estimate the barrier to the water proton exchange as 1267 cm-1, to be compared with the experimental barrier of 700 cm-1 derived from a onedimensional mode.8 The calculated height of the barrier is considerably higher than the experimental estimates. These estimates neglect the coupling between various intermolecular degrees of freedom in the complex. Given the fact that the recorded far-infrared spectra reveal several low-energy vibrational modes corresponding to large-amplitude motions in the complex, we believe that such a decoupling is not fully justified, and our ab initio result may be more realistic. 4. Despite a fairly weak character of the hydrogen bond in the H3N‚‚‚HOH complex, a substantial red shift (of 235 cm-1) and intensification (a factor of 104) was found for the OH stretch of water, in agreement with experiment (225 cm-1). The spectrum of NH3 is almost unaffected. In the present study we reproduced with good accuracy all essential features of the recorded infrared spectra in matrices. 5. SAPT calculations indicate that the bonding in the complex is mainly caused by strong electrostatic and induction interactions, quenched by the first-order exchange repulsion and second-order exchange-induction effects. The dispersion contribution is relatively less important and represents ∼30% of the induction energy. The stability of the bifurcated form is

(1) Atomic and Molecular Clusters; Bernstein, E. R. Ed.; Elsevier: Amsterdam, The Netherlands, 1990. (2) Dynamics of Polyatomic Van der Waals Complexes; NATO ASI Series B; Halberstadt, N., Janda, K., Eds.; Plenum: New York, 1990; Vol. 227/ (3) Nesbitt, D. J. Annu. ReV. Phys. Chem. 1994, 45, 367. (4) Hutson, J. M. Annu. ReV. Phys. Chem. 1990, 41, 123. (5) Cohen, R. C.; Saykally, R. J. J. Phys. Chem. 1992, 96, 1024. (6) Schmuttenmaer, C. A.; Cohen, R. C.; Saykally, R. J. J. Chem. Phys. 1994, 101, 146. (7) Leforestier, C.; Braly, L. B.; Liu, K.; Elrod, M. J.; Saykally, R. J. J. Chem. Phys. 1997, 106, 8527. (8) Stockman, P. A.; Bumgarner, R. E.; Suzuki, S.; Blake, G. A. J. Chem. Phys. 1992, 96, 2496. (9) Fraser, G. T.; Suenram, R. D. J. Chem. Phys. 1992, 96, 7287. (10) Herbine, P.; Dyke, T. R. J. Chem. Phys. 1985, 83, 3768. (11) Nelander, B.; Nord, L. J. Phys. Chem. 1982, 86, 4375. (12) Engdahl, A.; Nelander, B. J. Chem. Phys. 1989, 91, 6604. (13) Yeo, G. A.; Ford, T. A. Spectrochim. Acta A 1991, 47, 485. (14) del Bene, J. E. J. Am. Chem. Soc. 1973, 95, 5460; J. Phys. Chem. 1988, 92, 2874. (15) Diercksen, G. H. F.; Kraemer, W. P.; von Niessen, W. Theor. Chim. Acta 1972, 28, 67. (16) Yeo, G. A.; Ford, T. A. Can. J. Chem. 1991, 69, 632. (17) Ford, T. A. In Molecular InteractionssFrom Van der Waals to Strongly Bound Complexes; Scheiner, S., Ed.; Wiley: New York, 1997; p 181. (18) Latajka, Z.; Scheiner, S. J. Chem. Phys. 1990, 94, 217. (19) Skurski, P.; Gutowski, M. J. Chem. Phys. 1998, 108, 6303. (20) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Robb, M. A.; Cheeseman, J. R.; Keith, K.; Petersson, G. A.; Montgomery, J. A.; Raghavachari, K.; Al-Laham, M. A.; Zakrzewski, V. G.; Ortiz, J. V.; Foresman, J. B.; Cioslowski, J.; Stefanov, B. B.; Nanayakkara, A.; Challacombe, M.; Peng, C. Y.; Ayala, P. Y.; Chen, W.; Wong, M. W.; Andres, J. L.; Replogle, E. S.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Binkley, J. S.; Defrees, D. J.; Baker, J.; Stewart, J. P.; HeadGordon, M.; Gonzalez, C.; Pople, J. A. Gaussian, Inc., Pittsburgh, PA, 1995. (21) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (22) Perdew, J. P.; Wang, Y. Phys. ReV. B 1992, 45, 13244. (23) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553. (24) van Duijneveldt, F. B.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H. Chem. ReV. 1994, 94, 1873. (25) Duijneveldt, F. B. In Molecular InteractionssFrom Van der Waals to Strongly Bound Complexes; Scheiner, S., Ed.; Wiley: New York, 1997; p 81. (26) Turi, L.; Dannenberg, J. J. J. Phys. Chem. 1993, 97, 2488. (27) Zirz, G.; Ahlrichs, R. J. Chem. Phys. 1981, 85, 4989. (28) Xantheas, S. S. J. Chem. Phys. 1994, 100, 7523. (29) Xantheas, S. S.; Dang, L. X. J. Phys. Chem. 1996, 100, 3989. (30) Jeziorski, B.; Moszynski, R.; Szalewicz, K. Chem. ReV. 1994, 94, 1887.

8536 J. Phys. Chem. A, Vol. 103, No. 42, 1999 (31) Jeziorski, B.; Moszynski, R.; Ratkiewicz, A.; Rybak, S.; Szalewicz, K.; Williams, H. L. In Methods and Techniques in Computational Chemistry, METECC-94, Vol. B, Medium Size Systems; Clementi, E., Ed.; STEF: Caglian, 1993. (32) Moszynski, R.; Wormer, P. E. S.; van der Avoird, A. J. Chem. Phys. 1995, 102, 8385. (33) Moszynski, R.; Korona, T.; Wormer, P. E. S.; van der Avoird, A. J. Chem. Phys. 1995, 103, 321. (34) Moszynski, R.; Korona, T.; Wormer, P. E. S.; van der Avoird, A. J. Phys. Chem. 1997, 101, 4690. (35) Heijmen, T. G. A.; Korona, T.; Moszynski, R.; Wormer, P. E. S.; van der Avoird, A. J. Chem. Phys. 1997, 107, 902. (36) Jeziorski, B.; Moszynski, R.; Szalewicz, K. Chem. ReV. 1994, 94, 1887. (37) Szalewicz, K.; Jeziorski, B. Mol. Phys. 1979, 38, 191. (38) Jeziorski, B.; Moszynski, R.; Rybak, S.; Szalewicz, K. In ManyBody Methods in Quantum Chemistry; Lecture Notes in Chemistry 52; Kaldor, U. Ed.; Springer: New York, 1989; p 65. (39) Rybak, S.; Jeziorski, B.; Szalewicz, K. J. Chem. Phys. 1991, 95, 6576.

Sadlej et al. (40) Moszynski, R.; Jeziorski, B.; Ratkiewicz, A.; Rybak, S. J. Chem. Phys. 1993, 99, 8856. (41) Moszynski, R.; Jeziorski, B.; Szalewicz, K. J. Chem. Phys. 1994, 100, 1312. (42) Moszynski, R.; Jeziorski, B.; Rybak, S.; Szalewicz, K.; Williams, H. L. J. Chem. Phys. 1994, 100, 5080. (43) Moszynski, R.; Cybulski, S. M.; Chalasinski, G. J. Chem. Phys. 1994, 100, 4998. (44) Moszynski, R.; Jeziorski, B.; Szalewicz, K. Int. J. Quantum Chem. 1993, 45, 409. (45) Jeziorska, M.; Jeziorski, B.; Cizek, J. Int. J. Quantum Chem. 1987, 32, 149. (46) Moszynski, R.; Heijmen, T. G. A.; Jeziorski, B. Mol. Phys. 1996, 88, 741. (47) Dunning, Jr., T. H. J. Chem. Phys. 1989, 90, 1007. (48) Xantheas, S. S.; Dunning, T. H. J. Chem. Phys. 1993, 99, 8774. (49) Milet, A.; Korona, T.; Moszynski, R.; Kochanski, E. J. Chem. Phys. (submitted for publication).