Highly correlated systems: Structure, binding energy and harmonic vibrational frequencies of ozone

Volume 158, number

3.4

CHEMICAL

PHYSICS LETTERS

HIGHLY CORRELATED SYSTEMS: STRUCTURE, AND HARMONIC VIBRATIONAL FREQUENCIES Krishnan

RAGHAVACHARI,

AT&T Bell Laboratones,

9 June 1989

BINDING ENERGY OF OZONE

Gary W. TRUCKS

Murray Hill, NJ 07974, USA

John A. POPLE and Eric REPLOGLE Drpartmmt of Chemistry, Carnegie-Melion Received 14 February

University. Pittsburgh, PA 152 13, USA

1989: in final form 22 March 1989

The structure, binding energy and harmonic vibrational frequencies of 0, are calculated using the quadratic configuration interaction (QCI) method with a variety ofbasis sets. This method reproduces the geometry and vibrational frequencies ofozone fairly accurateiy using both restricted (RHF) and unrestricted Hartree-Fock (UHF) starting points. The asymmetric stretching frequency of ozone, which is extremely scnsitivc to the level of electron correlation, is calculated with a deviation of about 10% from expcnment. The symmetric frequencies are calculated with errors of l-596. The vibrational frequencies do not depend greatly on the level of correlation if a UHF starting point is used. The atomization energy is also accurately calculated, using a composite method where the effects of larger basis sets are evaluated using fourth-order perturbation theory.

1. Introduction It is well known that traditional single configuration Hartree-Fock theory #’ provides an inadequate description of the ground electronic structure of ozone [ 2-7 1. This is well recognized as due to the mixing of the ground state configuration ...( lb1)2(4b2)2(6a,)2( 1a2)2 with a significant amount of the low-lying doubly excited configuration ...(lb.)2(4b2)2(6a,)2(2b,)2. Recently, Stanton, Lipscomb, Magers and Bartlett [ 7 ] have calculated the geometry and harmonic vibrational frequencies of 0, with coupled-cluster- (CC) based methods [ 8lo] which include a high degree of electron correlation. They have pointed out that the asymmetric (b,) stretching frequency o3 of 0, is particularly sensitive to the level of electron correlation and that even the CCSDT- 1 method [ lo], which includes the effects of single, double and triple substitutions, obtains a 38% error for w3, yielding a value of only 680 cm - I, significantly smaller than the experimental *’ For a general introduction seeref.

to Hartree-Fock-based

methods,

[l].

0 009-2614/89/$ ( North-Holland

03.50 0 Elsevier Science Publishers Physics Publishing Division )

value [ 111 of 1089 cm-‘. Stanton et al. [7] have performed all their calculations with a restricted Hartree-Fock (RHF) starting point. It is also known [ 2,6] that oLone possesses a significant amount of biradical character with singly occupied x orbitals on the terminal oxygen atoms singlet coupled to each other. One possible description of such a state is based on spin-unrestricted Hartree-Fock theory (UHF) as reported previously by Rama Krishna and Jordan [ 5 1. In this paper, we explore the structure and vibrational frequencies of ozone using both RHF- and UHF-based methods.

2. Computational methods In this work, electron correlation effects have been included by means of the recently developed theory of quadratic configuration interaction (QCI) [ 121. This is a size-consistent single configuration HFbased method which is known to reproduce full configuration interaction results quite well, particularly near equilibrium geometries [ 121. Recently, it has also been shown to be reliable for accurate calculaB.V.

207

Volume 158, number 3,4

CHEMICALPHYSICSLETTERS

tions of bond energies of small molecules including spin-contaminated systems such as CN and HCC radicals [ 131, excitation energies of highly correlated systems such as CZ [ 141, and for the d”s’+ d”+‘s’ excitations of the first-row transition metal atoms [ 151. The QCI method [ 121 has been formulated as an iterative technique in the space of all single and double substitutions (QCISD). In addition, the effects of triple substitutions have been included in a noniterative perturbative manner, leading to a method termed QCISD (T). These methods have been compared in detail with various approximate treatments of coupled-cluster theory [ 141 and this work provides further comparisons between the two related techniques. Several basis sets have been used in this work. The standard 6-3 lG* basis set [ 161 (valence double-zeta plus a set of d-type polarization functions on each 0) has been used to calculate the geometries and vibrational frequencies at a variety of theoretical lcvels. These include HF, MP2 [ 17,181 (second-order Moller-Plesset perturbation theory), QCISD, and QCISD(T) using both RHF and UHF starting points. In addition, effects of larger basis sets were investigated using the 6-3 I G( 2d) basis set [ I] (valence double-zeta plus two sets of d-type functions) and the 6-31 lG* basis set [ 191 (valence triple-zeta plus d) at the RQCISD(T) (restricted QCI) level of theory. Finally, the polarized double-zeta basis set used by Stanton et al. [7] was also used to investigate if there was any special problem associated with that basis set. All the properties at the QCI level were evaluated completely by numerical techniques. In particular, since the frequencies involve double numerical differentiation, it is quite important [ 71 to evaluate the geometrical parameters precisely. The bond lengths at all levels were calculated correctly to 0.002 A and the bond angles to 0.02”. The force constants at the QCI level were computed in terms of internal coordinates and the frequencies then evaluated using the standard FG matrix formalism of Wilson, Decius and Cross [20]. For the RHF/6-3 lG* method, the frequencies determined by this numerical technique using a step size of 0.01 au agreed within 1 cm-’ with those computed analytically 12 1 ] ; therefore, we expect that all reported frequen208

9 June 1989

ties in this paper are also correct to within a few cm-‘. The MP2 frequencies reported in this paper were evaluated using analytical gradient techniques [ 2 I] where only one differentiation was performed numerically. All electrons were correlated in the MP2 calculations whereas only the valence electrons were correlated for the QCI studies. However, as shown previously by Stanton et al. [7], the core electron contributions to the frequencies are very small.

3. Results and discussion We have Iisted the RHF- and UHF-based calculations separately in tables 1 and 2. This facilitates comparisons of our results with those of previous calculations which are all RHF based. 3.1. RHF-based results As noted previously [ 2-7 1, the RHF model is quite inadequate to provide a proper description of the bonding in 0,. The bond length at the RHF/6-3 1G* level is too short by 0.07 A and, correspondingly, the stretching frequencies are too high by ~35% compared to the experimentally derived harmonic frequencies [ 1I ] which are also listed in the tables. As noted by Stanton et al. [ 71, second-order perturbation results (RMP2/6-3 I G*) are quite good for geometries and symmetric (a, ) frcqucncics but totally inadequate for the asymmetric (b,) stretching frequency. The mj(b2) at this level is 2373 cm-‘, compared to the experimental value [ 111 of 1089 cm-‘. Our RMP2/6-31G* value is very similar to that calculated by Stanton et al. [7] using a comparable (4s2pld) basis set. We have not calculated the higher-order perturbation results for the frequencies since they have already been reported by Stanton et al. [ 71. When higher-order electron correlation effects were evaluated by the QCI-based methods, the structures and frequencies are all in very good agreement with experiment. In particular, the c+(bL) frequency is calculated to be 921 and 950 cm- ’ at the RQCISD and RQCISD(T) /6-3 1G* levels, respectively, values which are only 15 and 13% lower than experiment [ 11I. The symmetric stretching frequencies are

Volume 158, number Table

CHEMICAL

PHYSICS LETTERS

Y June 1989

I

Prupertier

calculated

Basis 6-3

3,4

I c-t*

(4s2pld) 6-31G(2d) 6-31 lG*

with

RHF-based

methods a’

Method

Energy b,

&

0,

rfll (a,)

uJz(a,)

w(bz)

RHF RMP2 RQCISD RQCISD(T) RQclSD(T) RQCISD(T) RQCISD(T)

-224.26 144 -224.87615 -224.84761 -224.87905 -224.94285 -224.95396 -224.99053

1.204 1.299 1.216 1.298 1.289 1.284 1.279

119.0 116.3 117.4 116.7 1 17.0 116.6 117.0

1538 1174 1221 1116 II28 1129 1143

849 728 717 680 697 707 712

1454 2373 921 950 934 964 971

1.272

116.8

1135

716

1089

experiment

a) Bond lengths in A, bond angles m deg and frequencies in cm-‘. h1 Optimized total energies in hartree. Only valence electron correlations (see text).

Table 2 Properties

calculated

with UHF-based

have been included

in all correlated

techniques

except MP2

methods a)

Basis

Method

Energy h’

R,

0,

w,(a,)

w(al

6.31G*

UHF UMP2 UQCISD UQCISD(T)

-224.33151 - 224.8 1667 - 224.85472 - 224.86979

1.295 1.306 1.312 1.298

111.6 113.3 114.6 116.3

1137 1060 966 980

718 679 663 678

1163 1078 946 1037

1.272

116.8

1135

716

1089

expcrimcnt

‘i Bond lengths in A, bond angles in de-g and frequencies in cm-‘. bi Optimized total energies in hartree. Only valence electron correlations (see text).

described even better with a deviation of less than 5% at RQCISD(T)/6-31G* level. The surprisingly large difference between our results and those of Stanton et al. [ 7 ] prompted us to investigate whether part of this difference arises from basis set effects. Using the same (4s2pld) basis set employed by Stanton et al. [ 71, the calculated RQCISD (T) frequcncics (table I ) were all very close to the 6-3 lG* results. Thus we conclude that the differences between the two works are all due to the different treatment of electron correlation. Use of larger basis sets improves the results slightly. At the RQCISD(T)/6-31G(Zd) level, the two symmetric frequencies are within 1-2O/6 of experiment and the o)~ is only 11% lower than experiment. At the RQCISD (T) /6-3 11G* level, the frequencies are similar with deviations of less than 1% for wI and o2 and 11% for w3. The larger deviation for w3 shows the importance of correlation effects but all three values should be considered quite satisfactory. The

have been included

in all correlated

1

techniques

w(b)

except MP2

calculated geometrical parameters at the RQCISD(T)/6-3lG* level, 1.279 8, and 117.0” are both in excellent agreement with the corresponding experimental values [22] of 1.272 8, and 116.8”. In addition, the dipole moment calculated using finite field techniques at this level (0.48 D) is also in excellent agreement with the experimental value [23] of 0.53 D. 3.2. UHF-based results The UHF-based results are quite different from the RHF-based values. In particular, even at the UHF/ 6-31G” level, the calculated bond length (1.295 A) is much closer to the experiment than the corresponding RHF value. In addition, the calculated frequencies are also in excellent agreement with experiment with the maximum deviation being only 7O/o for the o3 frequency. Another interesting observation is that inclusion of the effects of electron cor209

Volume 158,

number 3,4

CHEMICAL

PIIYSICS LETTERS

relation has only a relatively small effect on the calculated frequencies. At the UHF, UMP2, UQClSD and UQCISD(T)/6-31G* levels of theory, the agreement is moderately good in all cases, the maximum deviation of any frequency being only about 15%. However, detailed comparison at higher correlated levels reveals that the RHF-based frequencies are better than the UHF-based ones. In particular, the relative ordering of the symmetric and asymmetric stretching frequencies is obtained incorrectly at the UHF level as well at the highest UQCISD (T) levels. However, it is important to note that the simple UHF theory provides a qualitatively correct description of such biradicals and this may be quite useful for a larger molecule where higher-order correlated treatments are not feasible. It is also interesting to note that at the highest level studied, UQCISD (T) /6-3 1G*, the calculated geometry is virtually identical to the RQCISD(T) /6-3 lG* geometry indicating the excellent performance of QCI theory.

4. Binding energy of ozone We have used a recently developed composite model to consider the total binding energy of ozone. In this method, termed Gl theory [ 131, we perform QCISD (T) calculations with a polarized 6-3 11G** basis set [ 19 1, using MP2/6-3 1G* geometries. Effects of diffuse functions and higher polarization functions are computed using fourth-order perturbation theory [ 24,25 1. Zero-point energy corrections as well as a higher-level correction (molecule independent but dependent on the unpaired electron count) are then computed. In the Gl theoretical model, the individual correction terms are assumed to be independent and are added together to yield a total energy which is then used to compute binding energies. The Gl theory was able to calculate the total atomization energies of 31 representative molecules with a mean deviation from experiment of only 1.5 kcal/mol. We have applied the G 1 theory to calculate the total atomization energy of ozone using both RHF and UHF starting points. The RHF-based method is more well defined due to the inclusion of the higher-level correction in the Gl theory which depends on the 210

9 June 1989

unpaired electron count. The calculated value of the total atomization energy using the RHF starting point is 143.7 kcal/mol, which is in excellent agreement with the experimental value of 142.2 kcal/mol. Use of the UHF starting point, with the same unpaired electron count, gives an atomization energy which is too low by 6.7 kcal/mol. Previously, Rama Krishna and Jordan [ 51 have performed perturbation theoretical calculations from a UHF starting point and negative dissociation energies for obtained O,+O,+O. This was attributed to the deficiencies arising from spin contamination in the UHF wavefunction. However, the fact that our RHF- and UHFbased QCI results are quite close illustrates that the UQCISD (T) scheme is able to take into account the principal problems resulting from the spin contamination of the UHF wavefunction.

5. Comparison with previous calculations In table 3, we compare our results obtained at the RQCISD(T) level with those from the literature. In this table, we have listed the results of Lee, Allen and Schaefer [ 61 using a two-configuration CISD calculation, CASSCF results of Alder-Golden et al. [ 41 and the CCSD+T(CCSD) and CCSDT-I results (augmented coupled-cluster methods) of Stanton et al. [ 71. The largest basis set used by Lee et al. [ 61 was of polarized triple-zeta (TZP) quality whereas the other two studies used polarized double-zeta (DZP) basis sets. In order to make reasonable intercomparisons possible, we have included our results with both 6-3 lG* and 6-3 11G* basis sets in table 3. It is clear from table 3 that all the methods describe the symmetric frequencies adequately with errors ranging from 1 to I O”h. However, the am frequency is described well only by the CASSCF and the QCISD(T) methods with errors around 10%. It should be noted that the QCI is a single-configuration-based method whereas the CASSCF is a multicontiguration treatment. The ZR-CISD has an error of about 24%. Surprisingly, the largest errors occur for the augmented coupled-cluster schemes. In particular, the CCSD + T (CCSD ) scheme yields an imaginary w3 frequency (yielding an asymmetric minimum) with the (4s2pld) basis set and a very low

Volume 158, number 3,4 Table 3 Comparison

CHEMICAL

of our results with previous publications Basis set

R,

2R-CISD b, CASSCF ‘) CCSD+T(CCSD) CCSDT-1 d, RQCISD(T) RQCISD(T)

TZP DZP DZP DZP 6-31G’ 6-31 IG’

1.252 1.296 1.293 1.295 1.298 I .279 1.272

experiment

9 June 1989

a’

Method

d,

PHYSICS LETTERS

w(a,)

w(a,)

w(b)

116.7 116.5 117.0 116.6 116.7 117.0

1242 1098 1097 1077 1116 1143

770 689 685 612 680 712

1353 989 I28i, 327 e1 680 950 971

116.8

1135

716

1089

a) Bond lengths in A, bond angles in degand frequencies in cm ‘. h1 Ref. [ 61. cl Ref. [4]. e1 The two values were obtained with a (4sZpld) and a (Ss3p2d) basis, respectively.

frequency (327 cm-‘) with a somewhat larger basis set. Even the computationally demanding CCSDT-1 scheme yields a value of only 680 cm-‘, too low by 38%. An interesting observation is that the effect of triple substitutions appears to be more dominant for CC schemes versus QCI-based methods. The QCISD and CCSD methods (without the triples correction ) yield values of 921 and 1240 cm- ’ for CC)~, respectively, (using DZP basis sets) thus having similar errors though the signs are different. The value of the amplitude for the doubly excited configuration 1 ( a2)2+2(b, )’ as evaluated by the QCISD method (0.22) is also very similar to that given by the CCSD method. However, it appears clear that the effects of triple substitutions are significantly overestimated using the CCSD+T(CCSD) or CCSDT-1 schemes. Thus, while the change in o3 is small on going from QCISD to QCISD( T), the corresponding change in going from CCSD to the augmented schemes is much larger. In this context, it should be noted that both QCISD and CCSD methods already include some contributions of triple substitutions due to the inclusion of product terms such as T,Tz (disconnected triples) [lo] and thus the comparison is not straightforward. We have discussed in another publication [ 141 that computationally the procedures most amenable for large-scale applications are the QCISD(T) and CCSD+T(CCSD) methods. Both methods have an iterative N6 procedure followed by a non-iterative N’ step to evaluate the effects of triples. However, the above results for ozone appear to suggest that the QCISD (T) scheme may have a more balanced treat-

d, Ref.

[ 71

ment of the effects of different substitutions. This appears to be related to the fact that the triples correction formula for the QCISD (T) scheme contains contributions from both singles and doubles whereas the T (CCSD) scheme contains only fourth-order-like contributions from double substitutions. Further investigations of the differences between Ihe QCI and CC schemes for such highly correlated systems are clearly desirable.

Acknowledgement We are indebted to J.F. Stanton, D.H. Magers and R.J. Bartlett for communicating their 0, results prior to publication.

References [ 1] W.J. Hehre, L. Radom, P. von R. Schleyer and J.A. Pople, Ab initio molecular orbital theory (Wiley, New York, 1986). [2] P.J. Hay, T.H. Dunning and W.A. Goddard, J. Chem. Phys. 62 (1975)

3912.

[ 31 W.D. Laid&and H.F. Schaefer III, J. Chem. Phys. 74 ( 1981) 3411. [4] S.M. Alder-Golden, S.R. Langhoff, C.W. Bauschlicher Jr. and G.D. Carney, J. Chem. Phys. 83 (1985) 255. [5 ] M.V. Rama Krishna and K.D. Jordan. Chem. Phys. 1I (1987) 423. [6] T.J. Lee, W.D. Allen and H.F. Schaefer 111, J. Chem. Phys. 87 (1987) 7062. [7] J.F. Stanton, W.N. Lipscomb, D.H. MagersandR.J Bartlett, J. Chem. Phys. 90 (1989) 1077. [8] J. Ciiek, J. Chem. Phys. 45 (1966) 4256; J. Paldus, J. Ciiek and 1. Shavitt, Phys. Rev. A 5 ( 1972) 50.

231

Volume 158, number

3,4

CHEMICAL

[Y] J.A. Pople. K. Krishnan, H.B. Schlegel and J.S. Btnkley, Intern. J. Quantum Chem. 14 (1978) 545. [ 101 R.J. Bartlett and G.D. Purvis III, Intern. J. Quantum Chem. 14 (1978) 561; G.D. Purvis 111and R.J. Bartlett, J. Chem. Phys. 76 (1982) 1910; Y S. Lee, S.A. Kurrhanki and R.J. Bartlett, J. Chem. Phys. 81 (1984) 5906; M. Urban, J. Noga, S.J. Cole and R.J. Bartlett, J. Chem. Phys. X3 (1085) 4041; J. Noga, R.J. Bartlett and M. Urban, Chem. Phys. Letters 134 (1987) 126.

[ 111 A. Barbe, C. Secroun and P. Jouve, J. Mol. Spectry. 49 (1974)

171.

[ 121J.A. Pople, M. Head-GordonandK.

Raghavachan, J. Chem. Phys. 87 (1987) 5968. [ 131J.A. Pnple, M. Head-Gordon, D.J. Fox, K. Raghavachari and LA Curtiss, J. Chem. Phys., to be published. [ 141 K. Raghavachari, J.A. Pople and M. Head-Gordon, in: Many-body methods in quantum chemistry (Springer, Berlin, 1989), to be published; J. Paldus, .I. &ek and B. Jeziorski, J. Chem. Phys., to be published. [ 151 K. Raghavacharr and ti.W. Trucks, J. Chem. Phys., to be published. [ 161P.C. Hariharan and J.A. Pople, Chem. Phys. Letters 66 (1972) 217.

212

9 June 1989

PHYSICS LETTERS

[ I7 ] C. Moller and M.S. Plesset, Phys. Rev. 46 ( 1934 ) 6 18. [ 181 J.A. Pople, J.S. Binklcy and R. Scegcr, Intern. J. Quantum Chem. Symp. IO (1976)

1.

[ 191 R. Krishnan, R. Seeger and J.A. Pople, J. Chem. Phys. 72 (1980) 6.50. [20] E.B. Wilson, J.C. Decius and P.C. Cross, Molecular vibrations (Dover, New York, 198 1) [21 ] J.A. Pople, R. Krishnan, H.B. Schlegel and J.S. Binkley, Intern. J. Quantum Chem. Symp. 13 ( 1979) 255.

[ 221 T. Tanaka and Y. Mot-ino, J. Mol. Spectty. 33 ( 1970) 538. [23] M. Lichtenstein, J.J. Gallagher and S.A. Clough, J. Mol. Spectty. 40 ( 197 1) 10. [24] R. Krishnan and J.A. Pople, Intern. J. Quantum Chem. 14 (1978) 91; R. Krishnan, M.J. Frisch and J.A. Pople, J. Chem. Phys. 72 (1980) 4244: M.J. Frisch, R. Krishnan and J.A. Pople, Chem. Phys. Letters 75 (1980) 66. [ 251 R.J. Bartlett and I. Shavitt, Chem. Phys. Letters 50 (1977) 190; R.J. Bartlett andG.D. Purvis III, J. Chem. Phys. 68 (1978) 2114; R.J. Bartlett, H. Sekino and G.D. Purvis III, Chem. Phys. Letters 98 (1983) 66.