J Phys. Chem. 1988, 92, 300-306
300
Vibrational-Rotational Analysis of ab Initio Potential Energy Surfaces for Symmetric-Top Molecules: Application to Ammonia I sotopomers Leonard 0. Hargiss Lever Research, Inc.. Edgewater, New Jersey 07020
and Walter C. Ermler* Department of Chemistry and Chemical Engineering, Stevens Institute of Technology, Hoboken, New Jersey 07030 (Received: February 27, 1987; In Final Form: June 15, 1987)
An ab initio potential energy surface for ammonia has been used to calculate vibrational-rotational spectroscopic constants. A complete quartic Taylor series fit in both internal and normal coordinates was derived. Spectroscopicconstants were calculated for deuteriated isotopomers by using normal coordinate force constants, equilibrium geometries, and normal mode eigenvectors within the framework of second-order perturbation theory. Anharmonic spectroscopic constants for ammonia compare favorably with other ab initio calculations. Normal mode force constants are in reasonable agreement with empirically obtained values.
Introduction
With the increasing refinement of a b initio techniques for calculating potential energy surfaces of polyatomic molecules, there have been many studies concerning the problem of reconciling ab initio energies with observed (mainly spectroscopic) data. Early work was performed by Eckart’ and Wilson and Howard,2 who developed the vibration-rotation Hamiltonian, and Nielsen,3 who used second-order perturbation theory to calculate the relationships between internal coordinate force constants, normal coordinate force constants, and spectroscopic term formulas for a number of molecular configurations. With the more recent a b initio methods for calculating molecular energies in terms of internal coordinates, it is natural to use these relationships for the comparison of calculated energies with observed quantities. In a series of papers, a procedure to calculate normal coordinate energy and property expansions from internal coordinate potential energy and property surfaces has been developed on the basis of perturbation theory.”7 The resulting computer programs accept as input energy and property surfaces in internal coordinates. The program used here, SURVIB,develops an expansion in symmetry-adapted internal coordinates about a calculated equilibrium point, using least-squares polynomial regression technique^.^ The normal mode analysis follows the procedure of Gwinn,s in which numerical second derivatives define the force constant matrix relative to mass-weighted Cartesian displacement coordinates. Eigenvectors and eigenvalues of the force constant matrix then define the Cartesian to normal coordinate transformation and the harmonic frequencies. A potential energy surface relative to normal coordinates is then calculated, based directly on energy values derived from the expansion initially calculated in terms of internal coordinates. In this way, the program accomplishes a proper nonlinear transformation between internal coordinate and normal coordinate descriptions of the energies or properties. The normal coordinate expansion constants and eigenvectors are then available for calculating spectroscopic constants, in the case of energy calculations, and vibrational expectation values, including anharmonicity effects, in the case of other properties. (1) Eckart, C. Phys. Reu. 1935, 47, 552. ( 2 ) Wilson, E. B., Jr.; Howard, J. B. J. Chem. Phys. 1936, 4 , 260. (3) Nielsen, H. H. Reu. Mod. Phys. 1951, 23, 90. (4) Krohn, B. J.; Ermler, W. C.; Kern, C. W. J. Chem. Phys. 1974, 60, 22. ( 5 ) Ermler, W. C.; Krohn, B. J. J. Chem. Phys. 1977, 67, 1360. ( 6 ) Ermler, W. C.; Rosenberg, B. J.; Shavitt, I. In Comparison of ab Initio Quantum Chemistry with Experiment f o r Small Molecules; Bartlett, R. J., Ed.: Reidel: New York. 1985. (7) Harding, L. B.; Ermler, W. C. J. Comput. Chem. 1985, 6 , 13. (8) Gwinn, W. D. J. Chem. Phys. 1971, 55, 411.
0022-365418812092-0300$01.50/0
The earliest extensive application of a precursor of SURVIB was to the water m o l e c ~ l e . ~ ”Expressions were derived by using second-order perturbation theory for energy and property anharmonicity corrections to the purely vibrating three-atom molecule and applied to H 2 0 and several isotopic variants. Later papers explored refinements of the procedure, including the introduction of the nonlinear transformation s ~ h e m e Applications .~ to additional properties included electric field gradients, dipole moments, and second moments, and higher quality a b initio calculations led to corresponding improvements in agreement with experiment. These papers presented evidence that disagreement with experimental values for the harmonic frequencies was principally due to error in the ab initio surface itself, since energies calculated from the internal coordinate expansion were demonstrably very close to the input energy surface. In this work, the procedure has been revised to accommodate symmetric-top molecules, again with the approximation of second-order perturbation theory, following the development of Mills.9 The present study reports internal and normal coordinate energy expansion coefficients, as well as spectroscopic constants, for ammonia and its deuteriated analogues. The ammonia molecule has its shortcomings as a candidate for this procedure, given its low “umbrella” inversion barrier of about 2000 cm-’. For this reason, some treatments have separated the vibration along the inversion coordinate ( u 2 ) from the other vibrations, treating it separately by a nonrigid invertor Hamiltonian.1° As a subject ammonia is compelling, however, since its normal coordinate force constants and spectroscopic constants are the subject of a broader literature than most symmetric tops.”,‘2 The noninversion motions, fortunately, show relatively small splitting from the inversion potential (1 cm-I or less for the lowest excited states)1° and thus may be expected to exhibit more normal behavior. Calculations
Ab Initio Potential Energy Surfaces. A 198-point potential energy surface for ammonia was calculated by using GAUSSIAN 82, a program for a b initio electronic structure calculation^.^^ (9) Mills, I. M. In Molecular Spectroscopy: Modern Research; Rao, K. N., Mathews, C. W., Eds.; Academic: New York, 1972. (10) Bunker, P. R.; Kraemer, W. P.; Spirko, V. Can. J . Phys. 1984, 62, 1801. (11) Gaw, J. F.; Handy, N. C. Chem. Phys. Lett. 1985, 121, 321. (12) Hoy, A. R.; Mills, I. M.; Strey, G. Mol. Phys. 1972, 24, 1265. (13) Binkley, J. S.; Frish, M.; Raghavachan, K.; Defrees, D.; Schlegel, H. B.; Whiteside, R.; Fluder, E.;Seeger, R.; Pople, J. A. GAUSSIAN-82, Release A ; Carnegie-Mellon University: Pittsburgh, PA, 1982.
0 1988 American Chemical Society
Ammonia Isotopomers Each point took about 25 min of CPU time on a VAX 11/780 computer. The calculations for NH3 used a 6-31G** basis set and the Hartree-Fock approximation plus fourth-order MerllerPlesset perturbation theory to estimate electron correlation. These features will be discussed briefly. The asterisks refer to (1) a set of six-component uncontracted (that is, single-exponent) 3d Cartesian Gaussian polarization functions for the heavy atoms (Li through F) and (2) a set of three-component uncontracted p-type polarization functions for each hydrogen (pX,pY,PA. The correlation energy was calculated according to MdlerPlesset perturbation theory.I4 The formulas employ the Hartree-Fock operator as the zero-order Hamiltonian and apply perturbation formulas to calculate contributions to the total energy based on the ground-state S C F wave functions. The GAUSSIAN 82 program calculates Merlle~Plessetcorrections to fourth order. This involves sums over single, double, and quadruple excitations, including second-, third-, and fourth-power products of two-electron antisymmetrized integrals in the numerator and first, second, and third powers of eigenvalues differences in the denominator. Analysis of Potential Energy Surfaces. The program used to perform the normal mode analysis is a modification of SURVIB.I5 This program has four main sections. Firstly, given input in the form of internal coordinates and energies, SURVIB calculates the equilibrium nuclear configuration and potential energy expansion in internal coordinates to a desired degree. Secondly, a conventional normal mode analysis is performed by (a) calculating second-degree force constants relative to mass-weighted Cartesian displacement coordinates by finite differences, (b) solving for the eigenvectors in mass-weighted Cartesian form, and (c) sorting the eigenvectors into vibrational and nonvibrational modes by eigenvalue. Thirdly, SURVIB creates an appropriate grid of points in normal coordinate space in increments of about 0.0005 au in the vicinity of the chosen stationary point and obtains the energy for each point by back-transforming into internal coordinates and using the internal coordinate energy expansion. Fourthly, the energy expansion in normal coordinates is obtained by performing a polynomial regression on the energy values from the previous section, obtaining normal mode force constants to a desired degree. Finally, spectroscopic constants and property expectation values are calculated by using formulas of Krohn et aL4 and M i l k 9 In the first section, SURVIB improves the economy of the calculation by symmetry adapting the internal coordinate expansion. Each term in the expansion associates an energy value with a product of powers of displacements from equilibrium in the internal coordinates. We can take as an example the term in p A k p d . SURVIB begins with a nuclear configuration of the same symmetry as the stationary configuration of the molecule studied, displaces coordinate pA by k units and coordinate pB by 1 units, and then calculates the corresponding three principal moments of inertia. The program does the same for every subsequent term in the expansion. If for any other term the three moments of inertia match the first, then this other term pCmpo" is treated identically with the first term in the polynomial regression, leading to a common expansion constant. In this way, the quartic expansion for ammonia is reduced from 210 to 51 terms. An advantage of the program is that the internal coordinate energy and property expansion constants, being isotopically invariant, need be calculated only once. Calculations for isotopic variants (isotopomers) of a molecule are then done with the same internal coordinate expansions. In the first section of the program, terms in bond length are usually expanded in the Simons-Parr-Finlan form, in powers of ( R - R,)/R.I6 In the third section, the energies are calculated (14) See, for example, Szabo, A.; Ostlund, N. S. Modern Quantum Chemistry: Introduction to Advanced Electronic Structure Theory; Macmillan: New York, 1982. (15) Harding, L. B.; Ermler, W. C . SURVIB: Potential Energy and Property Surface Fitting and Vibrational Analysis, Including Anharmonicity,f o r Polyatomic Systems, Program 505; Quantum Chemistry Program Exchange,
1986.
The Journal of Physical Chemistry, Vol. 92, No. 2, 1988 301 in internal coordinate space, so the energy is exact on the basis of the expansion in internal coordinates. The normal mode expansion constants are thus more accurate than would be possible in a linear coordinate transformation. Basis for the Perturbation Treatment. The starting point for the perturbation treatment of the complete vibration-rotation Hamiltonian is the harmonic oscillator approximation. The 3 N displacement coordinates pp of an N-atom molecule are defined in a Cartesian frame, such that pf = 0 at the equilibrium position of each nucleus in the internal potential field: pp
=rp-r&,
i = 1 , ..., N , a = x , y , z
(1)
The energy is obtained as a MacLaurin series expansion in pf:
where a,@, y = x, y , z ; i, j , k = 1, ..., N. The first-degree terms (dE/dpy) are zero by the condition that pf = 0 at the equilibrium positions of all nuclei. Through diagonalization of the matrix of mass-weighted second-degree Cartesian force constants, the orthonormal set of mass-weightd normal coordinate eigenvectors I;u may be obtained, such that (3) and
PP = C [ ( 1 / m , ' i 2 ) (rC l L Q , ) l
(4)
I
where i = 1, ..., N CY = x , y , z,; and r = 1, ..., (3N - 6).* It is also useful to define a dimensionless normal mode displacement coordinate q, = Q , ( ~ $ C W , / ~ )In~ these / ~ . ~coordinates *~ the potential energy of the molecular configuration may be expressed as
v = (1/2!)Cw,q,2 + ( 1 / 3 ! ) c c C 4 , s , q A A f + (1 /4!)CCCC+rsruq~sqrqu + ... ( 5 ) where w, are the harmonic frequencies and &,and &fu are the anharmonic force constants, all in cm-I. Force constants are commonly given in dynecm-' or N-m-'. If these are divided by nuclear masses expressed in g or kg, respectively, they have units of s-~. This is the basis for the convenience of mass-weighted coordinates. If the mass-weighted force constant kyp = ( m i m j ) - 1 / 2 ( d 2 E / d p y disp ~in) dyn.cm-'.g-' or N. m-l-kg-', then w will be in rad/s. If k is in atomic units of hartree.bohr-2, then w will be in rad/atu, where atu is the atomic time unit, 2.41888 X 1 O - l ' ~ . The common convention is that w is given in units of cm-I, obtained from rad/s using (27rc)-' = 5.30884 X It is noted that 1 rad/atu = 219480.4 cm-]." The energy levels of the vibrating-rotating symmetric-top molecule are given by the term formula
E = Go + Cw,(u,
+ y2) + ~ w , ( u r+ 1) +
C C d u , + 1/2)(us + Y2) + CCx,,(u, + f / Z ) ( U , + 1) + CCXrt,(u,+ 1)(u,l + 1) + CCgrJJp + Erot ( 6 ) where u, are positive integers, r and s are nondegenerate mode indices, t are (two-fold) degenerate mode indices, I, (the vibrational angular momentum quantum number) takes values from u, to -ur in increments of 2, and E,,, is the pure rotational energy (the formula for asymmetric tops simply omits all terms with t indi~es).~ The rotational energy formula for the symmetric top is covered in detail by Herzberg and need not be discussed here.I8 We give (16) (a) Simons, G.; Parr, R. G.; Finlan, J. M. J . Chem. Phys. 1973, 59, 3229. (b) Simons, G. J . Chem. Phys. 1974, 61, 369. (17) Values for constants from Cohen, E. R.; Taylor, B. N. J . Phys. Chem. Re$ Data 1973, 4, 663.
Hargiss and Ermler
302 The Journal of Physical Chemistry, Vol. 92, No. 2, 1988
criteria are that a: = 0 for all t2, and that all E,!,, = 0. This is equivalent to requiring mode t2 to be antisymmetric with respect to the xz plane. With these conditions satisfied, all t l and t2 eigenvectors are specified to within a sign. Note that with these restrictions all &2,2,2 = 0, since ql: is r 1 1 antisymmetric with respect to a symmetry operation of the (A similar formula holds for A , in the case of prolate symmetric molecule. The rule can be stated: All normal coordinate force tops.) The formula for B, is similar, omitting the term in qJ,. constants containing an odd power of any t2 are zero, except those The asymmetric-top formulas also omit this last term, being that contain odd powers in both t2 and t12. The sign of t 2 is thus identical in form for all axes. The a; and g, constants, together truly arbitrary. In general, t , does have a third-degree normal with the rotational constants, determine the rotational energy of coordinate force constant; therefore the choice of sign of the t l the molecule for a limited number of (usually low-lying) energy eigenvector will affect the sign of all force constants with odd states and have been published for a number of molecules. powers of t , . The perturbation analysis is simplified by defining certain Two general effects that figure strongly in vibration-rotation convenient constants, based on the normal mode eigenanalysis. spectra of low-lying states are Coriolis interactions and Fermi The inertial tensor expansion coefficients a:@appear in the perresonances. We have incorporated means for identifying and turbation formulas for vibration-rotation coupling constants. treating these effects into SURVIB. These are considered in turn. These are calculated from the atomic masses and mass-weighted Coriolis interactions in molecular spectra are of two kinds: Cartesian to normal coordinate eigenvectors: between nondegenerate vibrations and between the two vibrational coordinates making up a degenerate pair. The interaction between aLp = -2Xm,'I21;.Ji = -2Cmil/21;pi, (aP = x y , y z , z x ) displacement along the two coordinates defining a degenerate mode (84 and a rotation about the axis that couples them is particularly strong. In effect, instead of two independent modes in which atoms and move in approximately straight lines, the modes are free to combine. The Coriolis force deflects the vibration into a cyclical a:, = & ~ , ' / ~ ( l ; ~-y I;&, ~ (CY& = xyz, y z x , zxy) (8b) trajectory, such that the vibrational angular momentum about where cyi, pi,and y i are the equilibrium distances of particle i from the frame approaches in value the molecule-fixed angular mothe center of mass along the respective inertial axes. We also mentum K . Thus instead of defining separate vibrational quantum which may be considered define the Coriolis coupling constants numbers for the two degenerate vibrations, there is one vibrational as the projection of the cross-product of axis a with eigenvector quantum number ut and an angular momentum quantum number r onto eigenvector s: I,, which ranges from -ut to +u, in increments of 2. It is this 1, that is responsible for splitting the degeneracy in K values by f2(A{,),IlK in the rotational energy formula for the symmetric top.I8 In general, the rigid (prolate) symmetric top rotational energy is given by BJ(J + 1) + (A, - Bo)@. In the nonrigid case, The purpose of the perturbation analysis is to connect the however, because of this strong Coriolis interaction, there is the formulas for the internal coordinate energy expansion (eq 5) with term linear in K given above. Mills identifies {, with the Coriolis constants in the empirical term formula (eq 6). This calculation coupling constant S,l,r2;Hougen gives a more complicated forbegins with the complete nonrelativistic vibration-rotation Hamm ~ l a .Since ~ ~ we are concentrating on vibrational spectra, the iltonian, first derived by Wilson and Howard2 and Darling and details need not concern us here. Dennison.I9 Its first extensive application to molecules in the The K splitting is thus a first-order effect in the perturbation perturbation form was given by Nielsen in 1951.3 Watson pubexpansion, and the I , quantum numbers have their own anharlished a concise derivation in 1968.20 monicity constants, q,,,, in the energy term formula (see eq 6). The perturbation analysis was carried out, in the contact transformation form, by Nielsen and co-workers in the 1 9 4 0 ~ , ~ The I, also slightly affect the rotational constants about the symmetry axis; this is quantified by the I-doubling constant q, (eq using matrix elements evaluated by Shaffer*' and commutation 7). (Note that this I-type doubling constant q, is unrelated to the relations established by Herman and Shaffer.22 These make it dimensionless normal coordinates q, and q1 in eq 5, above. The possible to diagonalize the Hamiltonian, so far as is possible, by similarity in notation is unfortunate but these symbols have become using unitary transformations, then evaluating the off-diagonal established.) matrix elements remaining in terms of the quantum numbers, The calculated values for the q, of the symmetric tops are and inertial tensor expansion coefficients Coriolis coefficients presented together with the a:, to which they bear a formal a&. The spectroscopic constants in the empirical term formulas resemblance, in Table 111. are found by examining the coefficients of the quantum numbers Vibrational angular momentum is also present whenever a in the sum of the perturbation energies.23 The formulas for molecule rotates about an axis for which two modes have a nonzero symmetric tops were first made explicit by Mills.9 Coriolis coupling constant. Small interactions lead to a negative Two sets of perturbation formulas are given by M i l k g One term in the vibration-rotation coupling constant a:; large interset gives the anharmonicity constants x,, x,,, g,,, and gtl, in terms actions may lead to "Coriolis resonance", in which the Coriolis of the energy expansion constants in eq 5 and the Coriolis coupling term must be eliminated because of a near-zero denominator. constants. The second set of formulas gives the vibration-rotation Fermi resonances are another important effect in vibrational coupling constants a:, a;, and qr in terms of the equilibrium spectra. In this paper, the term "Fermi resonance" has a slightly rotational constants, Coriolis coupling constants, third- and restricted meaning. In the perturbation formulas for xrSand g,, fourth-degree constants from eq 5 and the inertial tensor expansion there are a number of terms with denominators like (4w; - w:) coefficients aaB. or [(w: + w,4 + w:) - 2(0,2w: w:q2 w120;)]. These Mills' formulas for the symmetric top require that the degenw, or w, w, + cor, denominators approach zero when 2w, erate vibrational modes assume a particular ~ r i e n t a t i o n . The ~ respectively. We note that, depending on the closeness of the harmonic frequencies to the observed frequencies, the denomi(18) Herzberg, G. Molecular Spectra and Molecular Structure I I Innators may or may not reflect an actual Fermi resonance. In any frared and Raman Spectra of Polyatomic Molecules; Van Nostrand Reinhold: case, however, the offending terms must be removed from the New York, 1945. calculation and the constant involved must be designated as having (19) Darling, B. T.; Dennison, D. M. Phys. Rev. 1940, 57, 128. (20) Watson, J. K. G. Mol. Phys. 1968, IS, 419. been modified. only the rotational constant C, in terms of the vibration-rotation coupling constants CY:, the I-doubling constant q,, and the vibrational angular momentum quantum number I,, discussed below:
e,,
cs,
+
(21) Shaffer, W. H. Rev. Mod. Phys. 1944, 16, 245. (22) Herman, R.; Shaffer, W. H. J . Chem. Phys. 1948, 16, 453. (23) Califano, S. Vibrational Stares; Wiley: London. 1976.
-
+
(24) Hougen, J. T. J . Chem. Phys. 1962, 47, 5410
-
Ammonia Isotopomers
The Journal of Physical Chemistry, Vol. 92, No. 2, 1988 303
TABLE I: Ammonia Polynomial Expansion in Terms of Internal Coordinateso term RNHI RNH2 RNH3 LHlNH2 LH 1NH3 0 0 1 0 0 0 0 0 0 0 2 1 1 0 0 0 0 3 0 0 4 2 0 0 0 1 0 1 0 5 0 1 0 0 1 6 0 0 7 1 0 0 0 0 0 2 8 0 1 1 9 0 0 0 0 0 0 0 10 3 0 1 0 0 11 2 0 2 0 0 1 12 0 0 0 2 0 13 0 1 0 1 1 14 0 0 1 15 1 1 1 1 1 0 0 16 0 2 0 1 0 17 1 0 1 1 0 18 1 0 1 0 19 0 0 0 0 1 0 20 0 3 0 0 0 21 2 1 0 0 22 0 0 1 1 0 0 23 4 0 0 0 24 0 0 0 0 3 1 25 0 1 0 3 0 26 0 0 0 27 3 0 0 0 0 2 2 28 1 1 0 0 29 2 1 0 1 0 30 2 0 0 1 2 1 31 0 0 0 2 1 32 0 0 2 2 0 33 1 1 0 0 2 34 1 0 0 2 0 35 0 0 0 0 36 2 1 0 1 1 37 1 2 0 1 1 0 38 1 1 1 0 1 39 2 0 1 1 0 40 1 1 0 0 1 41 3 0 0 0 42 1 2 1 0 0 43 1 2 0 0 1 0 44 1 1 0 0 45 1 1 0 1 0 0 46 0 0 0 47 1 0 4 0 0 0 0 48 3 1 0 0 49 0 2 2 0 0 0 50 1 0 0 2 51 0
LH2NH3 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 0 1 2 0 0 1 0 0 0 1 0 0 0 0 1 0 0 1 2 0 0 0 0 1 0 0 1 1 2 3 0 0 0 1
coeff -56.398 100 0.000 000 0.000 000 0.868 300 -0.012500 0.057 090 0.01 1 030 0.077 790 -0.01 1 370 -0.172 200 -0.01 5 070 0.064 740 0.009 068 -0.014 040 -0.087 440 0.023 910 -0.020 490 0.012 7 10 0.023 190 0.007 525 -0.025 280 -0.003 552 0.062 440 -0.337 100 -0.045 430 0.031 060 0.023 570 0.001 452 -0.020 880 -0.066 620 -0.027 210 0.068 790 -0.022 740 0.069 860 0.054 690 0.003 754 -0.025410 0.083 470 0.023 040 0.016 170 0.002 827 0.016220 0.019860 0.001 728 0.014280 0.005 887 -0.002 204 -0.010300 -0.002 800 0.021 280 -0.01 5 290
std dev 0.000 000 0.000 002 0.000 001 0.000075 0.000 094 0.000 022 0.000056 0.000006 0.000 008 0.000 963 0.001 301 0.000612 0.001 492 0.004 5 19 0.000 9 16 0.001 146 0.000 124 0.000 3 13 0,000 268 0.000 249 0.000 020 0.000 043 0.000 128 0.009 210 0.025 720 0.007 491 0.035 870 0.043 580 0.048 030 0.009 808 0.050 310 0.044 6 10 0.001 821 0.004 243 0.008 342 0.006 056 0.041 880 0.003 999 0.01 1090 0.007 987 0.022030 0.000 601 0.000 758 0.003 814 0.003 190 0.003 459 0.001 117 0.000059 0.000 165 0.000 276 0.000 307
'Ab initio 6-31G** MP4 quartic force field. The integers refer to exponents of internal coordinate displacements for each term in the energy expansion (see eq 2). The units are as follows: energy, hartrees; distance, bohrs; angle, radians. The units of bond length are expressed as powers of (R - R,)/R. Geometry at equilibrium is RNH = 1.035 A; LHNH = 105.9'.
It is relatively straightforward to treat Fermi resonances to first order. A secular determinant is constructed, with the zero-order energies along the diagonal and matrix elements of the cubic anharmonicity constants of eq 5 with respect to the zero-order wave functions of the two states constituting the off-diagonal elements. The solutions are the energies of the perturbed levels.ls It is possible to superimpose this calculation on the perturbation analysis used here, provided that the anharmonic perturbations taken into account here are elsewhere removed from the spectroscopic constants. Such a calculation requires substitution terms to be inserted into the spectroscopic constant formulas of eq 6 for each set of resonating frequencies. Such terms have been published for the asymmetric-top case by Califano for the anharmonicity constants xij" and by Papousek and Aliev for Coriolis resonances in the a;.25To our knowledge, no substitution terms have been rigorously derived for the symmetric top. Currently the recommended procedure is simply to eliminate terms with resonance denominators from the CY, and x, formulas and to designate them as modified constants with a symbol such as
In general, perturbation theory is complicated when near degeneracies arise in certain combinations of zero-order energy levels. In cases where large interactions occur between vibrational states, a better approach may be to use variational methods, in which the wave functions are explicitly expressed in normal coordinate space as linear combinations of rigid-rotator and harmonic-oscillator wave functions and the solutions found by eigenanalysis.26 Such methods should be able to cope with Fermi resonances but are currently applied predominantly to three-atom molecules because of unfavorable scaling with molecule size.6
Results and Discussion Table I lists the 5 1 symmetry-unique internal coordinate force constants of ammonia calculated from the 198-point 6-31G** MP4 (25) Papousek, D.; Aliev, M. R. Molecular Vibrational Rotational Specfra;Elsevier: Amsterdam, 1982. ( 2 6 ) Carney, G . D.; Sprandel, L. L.; Kern, C . W. Adu. Chem. Phys. 1978, 37, 305.
304
Hargiss and Ermler
The Journal of Physical Chemistry, Vol. 92, No. 2, 1988
TABLE I 1 Reduced Normal Coordinate Force Constants of Ammonia (cm-') NH, ND3 6-31G** emoirical" 6-31G** 3546.17 3503.0 WI 1140.82 1030.0 866.71 w2 3687.25 3591.6 2711.76 "3 1733.19 1689.9 1257.58 w4 813.37 -1333.25 -1285.4 4111 -60.42 147.94 239.0 4112 488.86 329.83 445.7 4122 -624.43 402.1 1 -516.4 $222 -1352.47 -1341.6 -845.95 6l3a3a = 6 l 3 b 3 b -13.72 -29.8 20.08 613a4a = d'13Mb 156.67 243.9 96.48 414a4a = d'14Mb 253.81 330.9 -168.39 b23a3a = 623b3b -80.76 -105.1 -42.60 d'23a4a = '#'23b4b -31.78 232.6 12.72 424a4a = 624b4b -978.20 -974.9 -629.69 b3a3a3a = -63a3b3b -122.85 -73.1 -76.93 43a3a4a = -63a3b4b = -43b3a4b -1 3 1.47 -141.7 -88.94 43a4a4a = -63a4b4b = -@)Ma46 -105.75 -96.2 -58.8 44a4a4a = -44a4b4b 415.52 418.8 225.73 4llll -69.22 18.70 41112 -267.26 -290.6 -149.97 41122 189.50 -99.26 dl222 59.13 47.08 -36.7 42222 436.76 235.56 438.3 4113a3a = 4113b3b 17.91 -3.44 4 l l 3 a 4 a = d'113b4b -80.74 -157.25 -174.7 6114a4a = 6114b4b -89.41 45.27 6123a3a = d'123b3b 60.06 -29.89 6123a4a = "123b4b 14.88 33.48 4124a4a = 4124b4b 318.24 175.10 613a3a3a = -613a3b3b 24.95 40.68 413a3a4a = -4l3b3b4a = -'#1' 3a3b4b 53.52 95.15 613a4a4a = -613a4Mb = -613Ma4b 34.09 13.93 614a4a4a = -414a4Mb -317.39 -336.5 -174.68 4223a3a = 4223b3b -6 1.94 -26.78 4223a4a = 6223Mb -22.35 51.6 -7.57 4224a4a = 4224b4b -76.42 45.55 $23a3a3a = -423a3b3b -26.16 55.60 423a3a4a = -423b3b4a = -423a3b4b -9.05 7.24 423a4a4a = -423b4a4b = -423a4b4b -7.19 4.01 624a4a4a = -624a4b4b 668.16 365.74 688.7 43a3a3a3a = 363s3a3b3b = @3b3b3b3b -9. I 7 16.27 #3a3a3a4a = 343a3a3b4b = 3'#3' a3b3b4a
Wavenumber 4000
I
~
~~
63b3b3Mb = 43b3Mb4b
jd'3a4a4b4b
-1 18.57 83.35 -284.94 -17.11
64b4b4b4b
-35.77
63a3a4a4a d'3a3Ma4b
43b3b4a4a 43a4a4a4a
=
= =
43a3a4Mb 3d'3b4a4a4b
=
63b4b4b4b
64a4a4a4a
= 344a4a4b4b =
~~
-109.0 -266.8 30.0
-59.66 46.68 -152.99 -4.32
-1 1.05
'Reference 12. potential energy surface. Normal coordinate force constants are given in Table I1 and N H , and ND, spectroscopic constants in Table 111. In these and in subsequent tables the indices have the following meanings: for C3, isotopomers, mode 1 is the symmetric stretch (A,), mode 2 is the umbrella bend (A,), mode 3 is the H N H symmetric stretch (E), and mode 4 is the H N H bend (E). For C, isotopomers, mode 1 (A') has N-H symmetric stretch character, mode 2 (A') has N-X stretch character, where X is the heavier hydrogen isotope, modes 3 and 4 (A') have symmetric scissors bend character, and modes 5 and 6 (A") have antisymmetric stretch and bend character, respectively. The calculation of the normal energy expansion is not symmetry adapted. Where certain relationships must be satisfied by the force constants, for example 4 r ~ 4 t ~ 4 r = 1 434t14t14t24~2 tl = 4r2@r24,Z4t2 for the degenerate modes of a symmetric top, small numerical discrepancies have been averaged out in the reported values. There is an inherent difficulty in comparing ab initio spectroscopic constants with experiment. It is possible to insert empirical harmonic frequencies into the analysis to check the accuracy
3500 3000
2500
1
Hydrogen Mass (a.u.)
Figure 1. Fundamental frequencies of NX3 vs X mass. Wavenumber 4000
I
1000 5001.008
1.259
1.511
1.762
2.01
Hydrogen Mass (a.u.) Figure 2. Fundamental frequencies of NH2X vs X mass.
of the higher force constants (these frequencies are most practically introduced at the beginning of the normal mode analysis, as the wi affect many details of computation). We have chosen to present purely ab initio constants to be consistent with previous work. The introduction of empirical parameters into SURVIB presents no particular difficulty, and has in fact been done at one stage in the analysis of formaldehyde.' Ammonia fared more poorly than molecules studied earlier in terms of the quality of the internal coordinate fit. The rootmean-square residual errors for the internal coordinate potential expansion were about 1 X loe6hartree for formaldehyde and 4 X lo6 for ammonia, for similar values of the geometric distortions. We attribute this to the strong anharmonicity of the inversion coordinate for ammonia; the inversion barrier is quite low (about 2000 cm-I) and makes it impossible to approximate the potential energy with a Taylor series for large distortions. Obviously for higher vibrational levels the approximation must break down completely, and the vibrational energies will no longer follow the perturbation formulas. The table of vibrational frequencies (Tables IV) includes estimated error values. These are based on the residual error in fitting energies to the internal coordinates and indicate approximately the maximum error that can be attributed to the Taylor series fit. Additional disagreement with experimental values must be attributed either to inaccuracy in the potential energy surface itself or breakdown of the perturbation theory approach. Table V lists spectroscopic constants for the asymmetric tops N H 2 D and NHD2. These isotopomers are interesting subjects for study, as they suffer from no Fermi interactions among their normal modes at the MP4 level of approximation. They also exhibit no Coriolis resonance problem since they have no neardegenerate frequencies of like symmetry. Consequently, our spectroscopic constants should be more accurate for these two molecules than for the symmetric-top isotopomers. With SURVIB,we can solve for the spectroscopic energies for a purely theoretical molecule. We have plotted fundamental frequencies for some pseudoisotopic variants of ammonia (Figures 1-3). When the C3, symmetry of a symmetric top is perturbed, say by formally increasing the mass of an atom in the x z plane,
The Journal of Physical Chemistry, Vol. 92, No. 2, 1988 305
Ammonia Isotopomers TABLE I11 Spectroscopic Constants of NH3 and ND3 (cm-I) NH3 LHNH, deg H-N, 8, w1
w2 w3 w4
XI I XI2 XI3 x14
x22 x23 x24
x33
7.34
x44 g33 g34 g44
4C
: ; 4 ag
:4; 93
94
6-31G** 105.9 1.035 3546.1 1140.8 3687.2 1733.1 -28.61 16.99 -107.04 14.59 -47.52 18.68 -7.15 -48.36 -15.26 -16.43 15.36 -3.71 9.79 0.0449 0.0509 0.0697 -0.0086 0.1557 0.1113 0.1624 -0.2383 0.0355 0.2255
empirical'
DZ + db
exptcSd
107.0 1.025 3503 1030 3591.6 1689.9 -29.81 8.86 -1 14.96 24.48 -41.77 17.06 35.66 -52.25 -21.44 -26.16 14.5 -6.9 13.2 0.023 0.037 -0.005 -0.003 0.201 0.089 0.218 -0.137 -0.026 0.255
107.6 1.003 3706.3 1187.8 3851.2 1841.3 -23.2 11.3 -92.3 118.1 -54.3 15.8 -11.2 -41.4 -14.2 -42.3 14.6 -8.2 34.3 0.044 0.062 0.015 0.064 0.127 0.110 0.160 -0.239 -0.1 190 0.3363
106.7 1.012 3506 1022 3577 1691 20.6 -92?? -6.7 32.4 -10.7 -18.5 -17.3 -8.8 2.7 0.078 0.098 -0.009 0.066 0.135 0.015 0.176 -0.230
6-31G**
ND3 D Z db
2534.2 866.7 2711.7 1257.5 -13.60 9.59 -54.27 48.50 -26.90 12.37 -3.48 -27.44 -8.77 -2.65' 9.41 -1.99 0.33c 0.0199 0.0180 0.0248 -0.0037 0.0510 0.0591 0.0584 -0.0892 0.0190 0.0795
2648.0 902.6 2836.3 1336.8 -11.5 5.6 -47.3 -62.6 -3 1.3 10.2 -5.2 -24.0 -9.0 9.3 9.1 -3.9 -12.8 0.020 0.022 0.007 0.021 0.043 0.057 0.058 -0.088 -0.037 0.1219
+
exptcsd
5.9 -45.8 17 -5?? -13.6 -18??
0.052 0.029 -0.0002 0.025 0.032 0.054 0.065 -0.080
'Reference 12. bReference 11. cBenedict, W. S.;Plyler, E. K. Can. J . Phys. 1957, 35, 1235; J . Chenz. Phys. 1960, 32, 32. dReference 27. eThese constants have been modified to omit terms containing Fermi resonance effects (see Results and Discussion). TABLE IV: Fundamental Frequencies of Ammonia as a Function of Hydrogen Mass (cm-I)" A. N H 3 through ND3
modeb w3 w1 w4 w2
N1H7 3498.1 (1.6) 3405.0 (6.0) 1681.3 (0.5) 1064.2 (1.2)
N1.2SH 3177.3 (1.1) 3075.8 (3.3) 1520.5 (1.2) 978.9 (1.7)
Nl.5H3 2939.9 (1.1) 2833.7 (4.0) 1399.2 (0.5) 916.0 (0.8)
N' 15H7 2761 .O (0.9) 2650.5 (3.4) 1310.0 (0.4) 869.6 (0.9)
ND, 2609.1 (0.5) 2506.1 (2.2) 1263.4c (0.5) 826.6 (0.8)
B. N H 3 through NHzD modeb
N1H7
w5
w4
3498.1 (1.6) 3498.1 (1.6) 3405.0 (6.0) 1681.3 (0.5) 1681.3 (0.5) 1064.2 (1.2)
modeb
N1H3
WI
w2 w3 w6
WI
w5 w2 *6
w3 w4
3498.1 (1.6) 3498.1 (1.6) 3405.0 (6.0) 1681.3 (0.5) 1681.3 (0.5) 1064.2 (1.2)
NH21,25H
NH,1,5H
NH,1.75H
NH,D
3498.6 (1.8) 3434.9 (4.0) 2876.6 (3.6) 1659.9 (1.7) 1525.2 (1.5) 1019.2 (2.0)
3502.6 (1.1) 3432.7 (3.1) 2699.3 (2.5) 1655.4 (0.6) 1480.0 (0.7) 1006.1 (0.9)
3502.0 (0.9) 3432.5 (3.4) 2552.0 (2.1) 1653.4 (0.6) 1445.2 (0.3) 994.6 (0.9)
C. NH3 through NHD2 NH',25H2 NH' 5H2 3469.9 (1.3) 3471.0 (1.1) 2906.7e (0.6) 3183.0 (0.8) 3097.8 (3.6) 2889.5 (3.4) 1617.8 (0.6) 1608.9' (0.4) 1433.OC(1.4) 1536.8 (0.3) 1008.4 (0.8) 970.3 (0.7)
NH1.l5H2 3460.4 (1 3 ) 2748.4 (0.7) 2643.2f (2.7) 1534.7 (0.5) 1345.6 (0.5) 930.5 (0.7)
NHD2 3462.7 (2.4) 2604.0 (1.2) 2481.1 (3.2) 1510.3 (0.9) 1283.7 (0.7) 904.6 (1.2)
3502.3 (1.4) 3441.7 (3.6) 3113.2d (3.8) 1671.9 (1.3) 1591.2d (1.2) 1040.3 (1.5)
'Exact masses used for calculation are l,2sH = 1.259394 amu, lsoH = 1.510964 amu, and I.l5H = 1.762533 amu. Estimated error based on residual error in internal coordinate fit is given in parentheses (see Analysis of Potential Energy Surfaces). bSpectroscopic numbering in parts B and C refer to asymmetric-top conventions. 'Spectroscopic constants used in this calculation have been modified by removing terms corresponding to a Fermi resonance between modes 1 and 4 (wl N 2w4). "Spectroscopic constants used in this calculation have been modified by removing terms corresponding to a Fermi resonance between modes 2 and 6 (w2 2 4 . 'Spectroscopic constants used in this calculation have been modified by removing terms corresponding to a Fermi resonance between modes 3, 5, and 6 (w5 N w 3 + w6). /Spectroscopic constants used in this calculation have been modified by removing terms corresponding to a Fermi resonance between modes 2 and 3 (u2 N 2w3). the particular pair of degenerate modes with motions symmetric and antisymmetric to this plane are forced apart in energy. This pair corresponds to the linear combinations denoted tl and t2. The symmetric or t l eigenvector, since it lies on the xz plane, exhibits a change in frequency with a change in mass of this atom, whereas
the t 2 eigenvector remains substantially unchanged (see Figures 2 and 3). The most appropriate comparison for our calculations is a recent report by Gaw and Handy," who used analytical formulas to obtain third-degree potential constants from an SCF a b initio
Hargiss and Ermler
306 The Journal of Physical Chemistry, Vol. 92, No. 2, 1988 Wavenumber
I
4000
2500 -
w1
1500
+
2oooQ 1000 Wl
5001,006
1.259 1.511 1.762 Hydrogen Mass (8.u.)
2.014
Figure 3. Fundamental frequencies of NHX2vs X mass.
TABLE V: Spectroscopic Constants of NHID and NHD2 (cm-’) A. NHZD wl w2 xII
xi2 x,3 ~ 1 4 XIS
XI6
~ 2 2
“I
“2 “3 “4
“5 “6
3598 2649
w3
-40.3 -14.5 -3.7 14.8 -156.8 -9.3 -41.7
~ 2 3
A axis 0.2237 0.0757 -0.0308 0.0443 0.1875 -0.3550
w4
~ 2 4
~ 2 5 ~ 2 6 x33
~ 3 4
x3*
1699 1052
~5
-5.0 17.8 -10.4 -14.7 -11.1 -8.7 -22.2
xjg ~ 4 4
~ 4 5 ~ 4 6
~ 5 5
x56 x66
B axis 0.0724 0.1503 -0.0483 -0.1009 0.0892 -0.0666
al ( ~ 2 “3
a4
“5 “6
3687 1478
w6
-7.7 -36.3 17.5 -10.5 -46.7 -11.2 -3.3
C axis 0.0146 0.0359 -0.0 14 1 0.0356 0.0010 0.0809
a1 “2 “3 “4
“5 “6
B. NHD2 Wl 02
XI I XI2 x13 x14 x16 x77
“1 “2
a3 CY^ “5 016
constants are more accurate (see Table 111). Note that while our harmonic frequencies lie between those of Gaw and Handy and the experimental values, the experimental geometry lies almost exactly between their S C F minimum and our MP4 minimum. Their bonds are too short and the pyramid too high, whereas our bonds are too long and the pyramid too flat. Gaw and Handy compared a simple double { basis set to a double {set with a d function on the nitrogen. The first led to a nearly planar equilibrium geometry, as would be expected due to the lack of polarization f ~ n c t i o n s . ’Results ~ from the latter ( D Z d) are similar to ours, both because they too use Mills’ perturbation formulas and because of similarities in the higher force constants. Their equilibrium geometry is comparable to ours: 1.003 A vs our 1.035 8, and 107.6’ vs our 105.9’. Their values and g, from this surface of spectroscopic constants wr, xrs, g,,, a,”, are included in Table 111 for comparison. The most serious discrepancy with Gaw and Handy in the table, for both N H 3 and ND, (the latter for which we have found no experimental values), is in the xi4 constants. Gaw and Handy found a strong Fermi resonance between modes 1 and 4 of NH3. Both studies find that 2w4 - w1 is smaller for ND, than for NH,. This term comprises a denominator in the formulas for x44and g4., thus making the perturbation analysis invalid. In accordance with accepted practice, these terms have been deleted in our calculations of our ND3 spectroscopicconstants. The contributions of these terms were +15.7 cm-’ to the x44and -1.1 cm-l to the g., This does not explain the anomalous xi4 value, however, since the resonance denominator does not appear in its calculation. The other constants for both molecules seem quite comparable. Hoy, Mills, and Strey12used the reverse procedure, calculating the values of or,x,, and a: from an empirical force field derived from algebraic transformation of internal to normal coordinates in an assumed potential form (Table 111). This study employed the internal coordinate force field of Morino, Kuchitsu, and Y a m a m ~ t owhich , ~ ~ included third- and fourth-degree terms in N-H bond stretching. The anharmonicity constants in internal coordinate space were obtained from data from diatomic systems. Due to the well-known problem of indeterminacy in XY3 molecules, it is not possible to determine all the higher order internal coordinate potential constants from spectroscopic data (there are fewer values of x, and a,“ than there are constants in the expansion). The authors solved this problem by assuming all higher order internal coordinate constants to be zero except those associated with the N-H stretching motion. Despite these simplifications, given that they are applying the same formulas to experimentally obtained potential constants, their results are a useful cross check. Some of the agreement, particularly in the diagonal constants, is quite good (see Table 11).
1318 952
3644 2590
~3
-80.7 -10.6 -8.9 16.0 -11.9 -24.3 -21.0
~ 2 3
A axis 0.1344 0.0615 -0.2307 0.0214 0.1371 -0.0274
~4
~ 2 4
~ 2 5 X26 ~ 3 3 ~ 3 4 ~
1
5
w5 w6
-42.4 8.5 -82.1 -7.4 4.0 -19.3 -12.8
x36 x44
x45 x46
x55 x56 X66
B axis “1
“2
a3 014
CY^ CY6
0.0651 0.0773 -0.0529 -0.0589 0.0836 -0.0402
CY] “2 “3 “4
“5 016
2712 1550 -2.1 -28.2 13.0 -0.8 -27.5 -11.9 -8.2
C axis -0.0010 0.0304 0.0648 0.0263 0.0117 -0.0196
calculation and calculated fourth-degree constants from finite differences. These calculations are thus far possible only at the S C F level, where there is a sizable positive error in the harmonic frequencies. This is well-known from past analyses of potential energy surfaces. They found, however, that the higher degree force
Acknowledgment. This research was partially supported by the National Science Foundation under grant no. CHE-8214689 and by the New Jersey Commission on Science and Technology. Registry No. NH,, 7664-41-7; ND,, 13550-49-7 (27) Morino, Y . ;Kuchitsu, K.; Yamamoto,S.Spectrochim. Acta, Part A 1968, 24A, 3 3 5 .
