Vibrational corrections to linear and nonlinear static electric properties of polyatomic molecules at non-optimum reference geometry

Chemical Physics 260 (2000) 1±10

www.elsevier.nl/locate/chemphys

Vibrational corrections to linear and nonlinear static electric properties of polyatomic molecules at non-optimum reference geometry Victoria E. Ingamells a,b, Manthos G. Papadopoulos b,*, Andrzej J. Sadlej c,d,1 b

a Clarendon Laboratory, Department of Physics, University of Oxford, Parks Road, Oxford OX1 3PU, UK Institute of Organic and Pharmaceutical Chemistry, National Hellenic Research Foundation, 48 Vasileos Constantinou, Athens G-11635, Greece c Department of Quantum Chemistry, Nicolaus Copernicus University, 7 Gagarin St., PL-87100 Toru n, Poland d Department of Theoretical Chemistry, University of Lund, Box 124, S-22100 Lund, Sweden

Received 22 May 2000

Abstract We show how vibrational contributions to the static electric properties of polyatomic molecules may be computed with respect to arbitrary reference geometries. The derived geometry corrective term that compensates for the e€ect of using a non-optimum geometry, is a simple function of the energy gradient and ®rst-order nuclear derivative of the electric property, summed over the number of normal vibrational modes of the molecule. We illustrate the possible orders of magnitude of the geometry correction for the major electric properties of the water molecule, computed at di€erent levels of electronic structure theory: SCF, MP2, CCSD and CCSD(T). Ó 2000 Elsevier Science B.V. All rights reserved.

1. Introduction The recognition of the importance of vibrational contributions beyond the usual vibrational averaging of molecular properties [1±4] is relatively recent [5,6]. Over the past decade, the study of what is called the pure vibrational contribution to second- and higher-order electric properties of molecules became, an area of intensive research, thanks mostly to the pioneering investigations of Bishop [7,8]. The original perturbation, theoretical derivation of analytic expressions for vibrational corrections by Bishop et al. [9], assumes that the chosen reference geometry corresponds to a minimum on the Born±Oppenheimer (BO) energy surface. To satisfy this condition, one needs to use optimized molecular geometries. Although for a given computational method, the calculation of the BO surface may not pose any serious problems, a comparison of a series of di€erent methods would require several geometry optimizations. For sizable *

Corresponding author. Tel.: +30-1727-3892; fax: +30-1727-3831. E-mail addresses: [email protected] (V.E. Ingamells), [email protected] (M.G. Papadopoulos), [email protected] (A.J. Sadlej). 1 Also corresponding author. Tel.: +4856-6114760; fax: +4856-6542477. 0301-0104/00/$ - see front matter Ó 2000 Elsevier Science B.V. All rights reserved. PII: S 0 3 0 1 - 0 1 0 4 ( 0 0 ) 0 0 2 5 4 - 8

2

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

molecules and high-level computational methods, the requirements imposed by the perturbation formulae of Bishop et al. [9] may be dicult to satisfy. One needs to stress that repeated geometry optimization is relatively easy at low levels of theory, e.g. in Hartree±Fock SCF and MP2 methods. However, at the level of coupled cluster schemes the geometry reoptimization becomes the main time-consuming factor. Hence, the use of molecular geometries optimized by using low-level theories is quite common. On the other hand, these methods are usually insucient for adequately accurate calculations of molecular electric properties. Combining low-level geometries with higher level property data leads then to the consideration of the geometry di€erence contribution to the evaluated vibrational corrections. Similar comments are also valid with regard to the use of di€erent basis sets for geometry optimization and for calculations of molecular electric properties. In our recent paper [10], we have presented an extension of the original formalism proposed by Bishop et al. [9] to the case of incompletely optimized, i.e. non-equilibrium molecular geometries. Our treatment assumes that the chosen molecular geometry used for the calculation of vibrational corrections is suciently close to the equilibrium geometry, so that the major di€erence in the Taylor expansion of the BO potential is given by the non-zero force contribution. A comparison with the formulae derived by Bishop et al. [9] exposes a second-order corrective term which accounts for the e€ect of using a non-optimum energy structure. This term ± denoted pgc for property p ± was shown [10] to be a simple function of forces acting on the nuclei at the selected geometry and property derivatives with respect to nuclear coordinates. Our extension of the perturbation approach to vibrational corrections to molecular properties was illustrated by pilot calculations of the electric properties for the BH molecule. It is obvious, however, that the pgc term will be of particular interest for polyatomic molecules. In such cases, a comparative study of vibrationally corrected electric properties with di€erent computational methods leads to a number of reference geometries which, according to the original theory of Bishop et al. [9], would have to be considered separately. The present paper reports on calculations of pgc for several electric properties of the water molecule computed by di€erent routine computational methods of quantum chemistry at optimum and non-optimum reference geometries. Although this molecule supports only three vibrational modes, and its geometry optimization is rather easy for any high-level computational method, the use of a single nonoptimum geometry in all calculations o€ers obvious advantages; this point will be documented in our paper. The data reported will also show the magnitude of error made in calculations at non-optimum geometries under the assumptions of the perturbation approach of Bishop et al. [9]. General formulation of the theory of vibrational corrections to properties at non-equilibrium geometries has been presented in our earlier paper [10], and only its brief summary is given in Section 2. Our calculations of vibrational corrections to the dipole moment, dipole polarizability and hyperpolarizabilities of the water molecule will be presented in Section 3. Their discussion, and the main conclusions of our study of the pgc contribution to the static electric properties of water, are given in Section 4. 2. The geometry change correction: a survey of theory In order to make this report complete, we summarize here some of the results obtained earlier [10]. We consider a polyatomic molecule in its rotationless …J ˆ 0† and lowest vibrational …n ˆ 0† state. Let Q denote the corresponding normal coordinates and Q0 , the selected (arbitrary) geometry of the molecule. The minimum energy geometry obtained either experimentally or by the use of some computational method will be denoted by Qe . We start with the zeroth-order vibrational problem, which is a set of independent harmonic oscillators:  X 1 1 Pr2 ‡ x2r Q2r wn…0† ˆ Wn…0† w…0† …1† n ; 2 2 r

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

3

where Pr is the momentum conjugate to normal coordinate Qr , xr denotes the frequency of the rth mode, and Wn…0† and w…0† n are the eigenenergy and eigenfunction of the nth vibrational state of the system. Eq. (1) is equivalent to the assumption that all forces vanish at Q0 ˆ Qe . However, if the reference geometry does not correspond to the energy minimum …Q0 6ˆ Qe †, then 

oE…Q; 0† oQr

 6ˆ 0;

r ˆ 1; 2; . . . ;

…2†

Q0

where E(Q,0) is the electronic energy calculated at zero external electric ®eld, F ˆ 0. Upon embedding the molecule in a homogeneous static electric ®eld, the electronic energy becomes 1 1 1 E…Q; F† ˆ E…Q; 0† ÿ la …Q†Fa ÿ aab …Q†Fa Fb ÿ babc …Q†Fa Fb Fc ÿ cabcd …Q†Fa Fb Fc Fd ÿ


; 2 3! 4!

…3†

where la ; aab ; babc ; dabcd ; . . . denote cartesian components of the dipole moment and di€erent polarizability tensors [11]. Hence, for F 6ˆ 0, the nuclear potential is modi®ed to account for the ®eld dependence: E…Q; F† ÿ E…Q0 ; F† ˆ

1 1X 2 2 X 1 …k† E …Q; F†; xr Qr ‡ 2 r k! kˆ1

…4†

where E(Q0 ,F) comprises all pure electronic contributions evaluated at Q0 .With respect to vibrational properties, these terms are constant and can be neglected. The relevant electric ®eld dependence of the nuclear potential which enters the ®eld-dependent vibrational equation is represented in terms of di€erent electric properties which themselves depend on nuclear coordinates. The ®eld-dependent term adding to the right-hand side of Eq. (4) is, therefore, a (Q,F)dependent multiple perturbation to the zeroth-order hamiltonian (1): 1 X 1 kˆ1

k!

E…k† …Q; F† ˆ

1 X X 1 …k;m† V …Q ; 0†Qp1 Qp2


Qpk Fa1 Fa2


Fam : k! p1 ;p2 ;...;pk mˆ0 m! p1 ;p2 ;...;pk ;a1 ;a2 ;...;am 0

1 X 1 kˆ1

…5†

The expansion coecients in Eq. (5) are de®ned as 

…Q0 ; 0† Vp…k;m† 1 ;p2 ;...;pk ;a1 ;a2 ;...;am

om‡k E…Q; F† ˆ …1 ÿ dk;2 dm;0 † oFa1 oFa2 . . . oFam oQp1 oQp2 . . . oQpk

 :

…6†

Q0 ;Fˆ0

These are recognized as the kth order nuclear coordinate derivatives of the mth order electric property of the system. Note that Eq. (5) contains the ®eld-independent …k ˆ 1; m ˆ 0† contribution due to the non-zero forces acting at Q0 (2). Applying the usual perturbation theory leads to di€erent types of vibrational contributions to the electric properties of the system. They can be derived as derivatives of the pth order perturbation corrections, W0p …F†, to the vibrational energy, W00 …0†, of Eq. (1), with respect to the electric ®eld strength. The ®rst-order corrections W01 …F† have been found to give the usual zero-point vibrational contributions [10]. Of interest here is the second-order …p ˆ 2† contribution, which provides a general form for the second-order vibrational correction to the mth order electric property [10]:

4

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10 …2;m† W0;a1 ;a2 ;...;am

…2†

om W0 …F† ˆÿ oFa1 oFa2 . . . oFam ˆ

1 X 1

! Fˆ0

1 X X 1 X X

1

D E …0† …0† w jQ Q . . . Q jw p1 p2 pk 0 s …0†

k! p1 ;p2 ;...;pk lˆ1 l! q1 ;q2 ;...;ql s6ˆ0 Ws…0† ÿ W0 m D EX 1 …0† Pa1 a2 ...at at‡1 ...am Vp…k;t†  ws…0† jQq1 Qq2 . . . Qql jw0 1 ;p2 ;...;pk ;a1 ;a2 ;...;at t!…m ÿ t†! tˆ0 kˆ1

 …Q0 ; 0†Vq…l;mÿt† …Q0 ; 0†; 1 ;q2 ;...;ql ;at‡1 ;at‡2 ;...;am

…7†

where Pa1 a2 ...at at‡1 ...am produces all possible permutations of the coordinate labels. This formula combines di€erent orders of the coordinate dependence of ®eld-independent energies (either t ˆ 0 or t ˆ m) with di€erent orders of the coordinate dependence of electric properties. In particular, terms with either k ˆ 1 and t ˆ 0 or l ˆ 1 and t ˆ m, i.e. those arising from Eq. (2), will in general give a non-vanishing contribution to electric properties of any order. These terms are responsible for the e€ect which results from the use of the non-optimum geometry. From the general formula (7), one can derive explicit expressions for the major electric properties of any molecule [10], namely the dipole moment, polarizability and ®rst and second hyperpolarizabilities. We will restrict our interests to vibrational corrections to the dipole moment …m ˆ 1†, dipole polarizability …m ˆ 2† and ®rst …m ˆ 3† and second …m ˆ 4† hyperpolarizabilities in the vibrational ground state …n ˆ 0†. For the lowest-order …k ˆ l ˆ 1† term, in the expansion with respect to normal coordinates, Eq. (7) becomes D ED E XXX 1 …2;m† …0† …0† …0† …0† w jQ jw jQ jw w W0;a1 a2 ...am ˆ p q 0 0 s s …0† …0† ÿ W0 p q s6ˆ0 Ws m X 1 …1;t† …l;mÿt† Pa1 a2 ...at at‡1 ...am Vp;a …Q0 ; 0†Vq;a …Q0 ; 0†: …8†  1 a2 ...at t‡1 at‡2 ...am t!…m ÿ t†! tˆ0 The usual selection rules for matrix elements of Qp in the basis of the harmonic oscillator functions combined with the standard formulae for the non-vanishing matrix elements, result in the following (lowest-order) expression for the geometry correction to the vibrational contribution to the dipole moment:     X 1  X h oE ola …2;1† …1;1† Pa Vp…1;0† Vp;a P …9† ˆ lgc ‡ Va…1;1† V …1;0† ˆ ÿ xÿ2 W0;a ˆ a p a : oQ h  x 2x oQ p p p p p p The polarizability tensors a, b and c follow in a similar way with m ˆ 2, 3, 4, respectively.   X1 oE oaab ola olb …2;2† 2 …0;0† xÿ2 P ; ÿ ˆ agc W0;ab ˆ ÿ ab ab ‡ ‰l Š p oQ oQ oQ 2 oQ p p p p p …2;3† W0;abc

…2;4† W0;abcd

ˆÿ

X1 p

ˆÿ

2

 xÿ2 p Pabc

X1 p

2

1 oE obabc ola oaab ‡ 3 oQp oQp oQp oQp

 xÿ2 p Pabcd



ˆ bgc abc ‡ ‰laŠ

…0;0†

;

1 oE ocabcd 1 ola obbcd 1 oaab oacd ‡ ‡ 12 oQp oQp 3 oQp oQp 4 oQp oQp

…10†

…11† 

…0;0† 2 …0;0† ‰a Š : ˆ cgc abcd ‡ ‰lbŠ

…12†

In all these formulae, the permutation operators Pab . . . produce a sum of all terms resulting from permutations of the cartesian coordinate labels. For each property, the terms following the geometry cor-

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

5

rection may be identi®ed as the harmonic contributions which are already given by the Bishop±Kirtman perturbation equations [9]. A more detailed derivation of the formulae presented in this section can be found in our earlier paper [10]. 3. Computational details and results We have selected the water molecule to illustrate the magnitude and importance of the geometry correction to static electric properties of polyatomic molecules. The present study covers its dipole moment, dipole polarizability, and ®rst and second hyperpolarizabilities which are computed with the most widely used methods of computational quantum chemistry, i.e. the SCF, MP2, CCSD, and CCSD(T) methods [12]. It should be stressed that we do not attempt to obtain the most accurate data for the properties studied. The present goal is to study the geometry change contribution which accompanies the use of nonoptimum reference geometries. All calculations of electric properties reported in this paper have been carried out with the so-called polarized (Pol [13]) basis sets. The medium-sized Pol sets ([10.6.4/5.3.2] for oxygen and [6.4/3.2] for hydrogen) are suciently ¯exible for dipole moment and dipole polarizability calculations. Also, in the case of hyperpolarizabilities, their performance is acceptably good [14]. The optimized molecular geometries have also been calculated with the Pol basis sets. The required nuclear derivatives of the electric properties were obtained numerically [15] by taking di€erences of dipole moment, polarizability and hyperpolarizabilities with respect to nuclear displacements about the given reference geometry. In turn, the values of l, a, b and c were determined by taking ®nite di€erences of the electric-®eld-perturbed electronic energy. Such calculations force an extremely tight convergence on the electronic energy, particularly in the estimation of accurate ®rst and second derivatives of b and c. At all four levels of theory, the electronic energies were converged to a criterion of 12 decimal places, which was found sucient to determine reliable electric properties via the ®nite ®eld approach. A suitable ®eld strength was selected following several preliminary calculations [16], where numerical results were compared with analytically determined data at the SCF level. A value of 0.003 atomic units was ®nally chosen as our ®nite-®eld value. For the non-optimum geometries and certain directions of the external electric ®eld, the CCSD and CCSD(T) energies were found to be insuciently accurate for the numerical evaluation of the required derivatives. For this reason, three components of czpva …xxxx; xxyy and xxzz† are left undetermined. It should be stressed, however, that the absence of these data in no way a€ects the main issues discussed in this paper. Our results are shown in Tables 1±4 for the SCF, MP2, CCSD and CCSD(T) methods, respectively, where the geometry correction is listed in the penultimate column of each table. For each computational method considered in this study, we determine the optimum energy structure of water and present the corresponding electronic pel , zero-point vibrational average (ZPVA) pZPVA and pure vibrational pPV contributions. The geometry correction formulae are derived under the assumption of the harmonic approximation. In view of this, the pure vibrational contributions presented in the tables include only the harmonic terms arising from the Bishop±Kirtman perturbation theory; we refer to those terms given in Eqs. (10)±(12) for the polarizabilities, and ®rst and second hyperpolarizabilities, respectively. Similarly, the ZVPA corrections only include the contribution due to a harmonic potential approximation. We note that the vibrational corrections, including pure and ZPVA contributions, to the electric properties of water have been computed previously by Bishop et al. [17]. They also apply the Bishop±Kirtman formulae but include mechanical and electrical anharmonic perturbation terms. Consequently, a direct comparison of our present results with those of Ref. [17] only re¯ects the di€erence in approximation. However, we observe that the two sets of results show similar orders of magnitude. This is expected since the harmonic perturbation terms tend to provide the largest contribution to the vibrational correction, particularly for l and a. Therefore, we may be assured that our results are correct within the harmonic approximation.

6

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

Table 1 SCF results: the electronic and (harmonic) ZVPA and pure vibrational contributions to the electric properties l, a, b and c of H2 Oa   Optimum: r…OH† ˆ 0:9450 A, Non-optimum: r…OH† ˆ 0:9861 A, …pel …Q0 † ÿ pel …Qe †† \HOH ˆ 105:8° \HOH ˆ 103:5° pel

pZPVA

pPV

pel

pZPVA

pPV

pgc

lz

0.769

ÿ0.002

0.0

0.791

ÿ0.002

0.0

0.022

0.022

axx ayy azz

7.792 8.955 8.337

0.021 0.143 0.099

0.0 0.140 0.859

7.891 9.366 8.652

0.021 0.153 0.104

0.0 0.143 0.756

0.111 0.462 0.348

0.099 0.411 0.315

bzxx bzyy bzzz

0.515 ÿ9.009 ÿ4.094

ÿ0.031 ÿ0.603 ÿ0.475

ÿ0.436 3.408 8.754

0.626 ÿ10.407 ÿ5.120

ÿ0.042 ÿ0.676 ÿ0.547

ÿ0.176 4.519 10.353

0.121 ÿ1.595 ÿ1.135

0.111 ÿ1.398 ÿ1.026

cxxxx cxxyy cxxzz cyyyy cyyzz czzzz

1216 291 326 447 273 760

ÿ3.7 2.5 2.8 32.2 8.7 21.5

20.2 28.3 15.1 220.7 73.3 58.5

1226 297 331 499 291 787

ÿ4.1 3.4 3.2 35.9 9.6 22.2

23.3 33.5 18.9 287.0 98.6 92.4

12 7 5 60 21 33

10 6 5 52 18 27

a

Values are reported at the optimized SCF/Pol geometry …Qe †, and at an arbitrary (non-optimum) geometry …Q0 † that corresponds to the MP2/Pol optimized structure. The last two columns compare the geometry correction pgc for each property with the di€erence in the electronic contributions at the optimum and non-optimum geometries.

Table 2 MP2 results: the electronic and (harmonic) ZVPA and pure vibrational contributions to the electric properties l, a, b and c of H2 Oa   Optimum: r…OH† ˆ 0:9861 A, Non-optimum: r…OH† ˆ 0:9450 A, …pel …Q0 † ÿ pel …Qe †† \HOH ˆ 103:5° \HOH ˆ 105:8° pel

pZPVA

pPV

pel

pZPVA

pPV

pgc

lz

0.738

ÿ0.002

0.0

0.720

ÿ0.003

0.0

ÿ0.019

ÿ0.018

axx ayy azz

9.642 10.246 9.913

0.029 0.126 0.096

0.0 0.117 0.700

9.483 9.846 9.578

0.031 0.126 0.099

0.0 0.0 0.843

ÿ0.143 ÿ0.356 ÿ0.303

ÿ0.159 0.400 ÿ0.335

bzxx bzyy bzzz

ÿ2.294 ÿ11.047 ÿ9.138

0.072 ÿ0.227 ÿ0.094

ÿ0.322 2.407 7.896

ÿ2.428 ÿ9.995 ÿ8.100

0.055 ÿ0.285 ÿ0.151

ÿ0.670 1.521 6.456

ÿ0.005 0.951 0.979

ÿ0.046 1.052 1.038

cxxxx cxxyy cxxzz cyyyy cyyzz czzzz

2187 559 586 657 427 1237

7.5 6.6 6.8 31.7 6.5 21.4

52.6 51.5 24.6 315.3 121.3 128.0

2135 544 574 605 411 1206

8.7 7.8 5.6 27.6 5.6 23.5

46.6 43.8 17.1 253.9 97.8 85.5

ÿ47 ÿ13 ÿ10 ÿ46 ÿ13 ÿ22

ÿ52 ÿ15 ÿ12 ÿ52 ÿ16 ÿ31

a

Values are reported at the optimized MP2/Pol bond length (Qe ), and at an arbitrary (non-optimum) bond length (Q0 ) that is the optimized SCF/Pol geometry. The last two columns compare the geometry correction pgc for each property with the di€erence in the electronic contributions at the optimum and non-optimum geometries.

A similar set of electronic and vibrational contributions are reported at the selected non-optimum geometry, together with the non-zero geometry corrective term, pgc . At the SCF level, the non-optimum geometry refers to the MP2 optimized structure, whereas the correlated methods use the SCF optimized geometry. This choice corresponds to a relatively common strategy in which the geometry optimization is

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

7

Table 3 CCSD results: the electronic and (harmonic) ZVPA and pure vibrational contributions to the electric properties l, a, b and c of H2 Oa   Optimum: r…OH† ˆ 0:9676 A, Non-optimum: r…OH† ˆ 0:9450 A, …pel …Q0 † ÿ pel …Qe †† \HOH ˆ 103:5°

\HOH ˆ 105:8°

pel

pZPVA

pPV

pel

pZPVA

pPV

pgc

lz

0.736

ÿ0.003

0.0

0.719

ÿ0.003

0.0

ÿ0.017

ÿ0.017

axx ayy azz

9.322 10.084 9.642

0.021 0.137 0.097

0.0 0.087 0.670

9.188 9.692 9.326

0.024 0.133 0.097

0.0 0.091 0.805

ÿ0.122 ÿ0.352 ÿ0.282

ÿ0.134 ÿ0.392 ÿ0.316

bzxx bzyy bzzz

ÿ1.989 ÿ11.336 ÿ7.730

0.054 ÿ0.376 ÿ0.190

ÿ0.407 2.016 7.627

ÿ2.046 ÿ10.161 ÿ6.831

0.046 ÿ0.403 ÿ0.213

ÿ0.683 1.328 6.454

ÿ0.038 1.064 0.828

ÿ0.057 1.175 0.899

cxxxx cxxyy cxxzz cyyyy cyyzz czzzz

2062 543 554 613 403 1125

8.7 9.2 9.0 38.6 7.7 29.2

39.7 44.6 19.8 313.5 119.7 119.3

2024 528 545 563 388 1105

±b ±b ±b 29.3 6.4 20.8

35.8 38.4 13.9 250.1 94.7 80.8

ÿ37 ÿ13 ÿ9 ÿ44 ÿ13 ÿ17

ÿ38 ÿ15 ÿ9 ÿ50 ÿ15 ÿ20

a

Values are reported at the optimized CCSD/Pol bond length (Qe ), and at an arbitrary (non-optimum) bond length (Q0 ) that is the optimized SCF/Pol geometry. The last two columns compare the geometry correction pgc for each property with the di€erence in the electronic contributions at the optimum and non-optimum geometries. b Value undetermined due to insucient accuracy in the CCSD energies leading to the second-order derivative of c. Table 4 CCSD(T) results: the electronic and (harmonic) ZVPA and pure vibrational contributions to the electric properties l, a, b and c of H2 Oa   Optimum: r…OH† ˆ 0:9694 A, Non-optimum: r…OH† ˆ 0:9450 A, …pel …Q0 † ÿ pel …Qe †† \HOH ˆ 103:5° \HOH ˆ 105:8° pel

pZPVA

pPV

pel

pZPVA

pPV

pgc

lz

0.728

ÿ0.003

0.0

0.719

ÿ0.003

0.0

0.018

0.016

axx ayy azz

9.541 10.269 9.847

0.023 0.139 0.098

0.0 0.082 0.655

9.387 9.835 9.500

0.026 0.135 0.099

0.0 0.086 0.797

ÿ0.138 ÿ0.384 ÿ0.311

ÿ0.154 ÿ0.434 ÿ0.346

bzxx bzyy bzzz

ÿ2.686 ÿ11.986 ÿ9.027

0.068 ÿ0.353 ÿ0.144

ÿ0.456 1.791 7.351

ÿ2.719 ÿ10.679 ÿ8.000

0.052 ÿ0.401 ÿ0.192

ÿ0.748 1.095 6.152

0.002 1.167 0.961

0.033 1.307 1.027

cxxxx cxxyy cxxzz cyyyy cyyzz czzzz

2237 593 605 665 437 1227

3.4 4.3 5.9 32.0 16.1 27.0

44.5 48.0 21.0 332.9 128.5 130.4

2186 575 592 608 420 1201

±b ±b ±b 29.8 6.4 20.6

39.6 40.6 13.9 260.7 100.4 86.9

ÿ46 ÿ16 ÿ11 ÿ50 ÿ15 ÿ19

ÿ51 ÿ18 ÿ13 ÿ57 ÿ18 ÿ26

a

Values are reported at the optimized CCSD(T)/Pol bond length (Qe ), and at an arbitrary (non-optimum) bond length (Q0 ) that is the optimized SCF/Pol geometry. The last two columns compare the geometry correction pgc for each property with the di€erence in the electronic contributions at the optimum and non-optimum geometries. b Value undetermined due to insucient accuracy in the CCSD(T) energies leading to the second-order derivative of c.

carried out at some low-level approximation and the selected molecular features are computed at that geometry by using high-level theories.

8

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

All SCF and MP2 electronic structure calculations were performed with the C A D P A C program [18]. The program [12] was used to perform the coupled-cluster CCSD and CCSD(T) calculations. The potential energy and electric property derivatives were determined using a numerical derivative algorithm developed by one of us [19]. GAUSSIAN

4. Discussion and conclusions An inspection of our results presented in Tables 1±4 clearly indicates the importance of the pgc term for all studied electric properties. In this respect, the geometry change contribution to the dipole moment makes an exception. The magnitude of lgc is quite small, being of the order of 10ÿ2 at each level of theory. However, one should note that this small contribution is larger than that which follows from vibrational averaging. At the level of the SCF approximation, the pgc contribution to the dipole polarizability evaluated at the MP2 optimum geometry is found to be commensurable with the pure vibrational correction and is de®nitely larger than the pZPVA term. For b and c tensors, the relative importance of the pgc contribution appears to diminish in proportion to the increase of the magnitude of the pPV correction. A similar pattern is observed for MP2, CCSD, and CCSD(T) data calculated at the optimized SCF geometry. The corrective term increases from around a tenth to one atomic unit for both the polarizability and ®rst hyperpolarizability. The largest value for the geometry correction is found for the second hyperpolarizability. At the MP2 and coupled-cluster levels, the magnitude of cgc is between 2 and 50 atomic units. This signi®es the di€erence in the SCF and correlated optimum geometries. The last column of each table lists the di€erence in the electronic property at the optimum and nonoptimum geometries. This value should be compared with the geometry correction. In our previous paper, we illustrated how properties determined at an arbitrary geometry may be used to obtain properties at the optimum geometry. Computing the total electric property with respect to a given arbitrary reference geometry Q0 gives ptotal …Q0 † ˆ pel …Q0 † ‡ pZPVA …Q0 † ‡ pPV …Q0 †:

…13†

However, we recall from perturbation theory that at a non-optimum geometry the second-order vibrational correction included a term due to the non-zero forces acting on the molecule. Thus, in terms of the perturbation expansion we have the following approximate relation: ptotal …Q0 †  pel …Qe † ‡ pZPVA …Q0 † ‡ ‰pgc …Q0 † ‡ pPV …Q0 †Š:

…14†

Equating (13) and (14) gives pel …Q0 †  pel …Qe † ‡ pgc …Q0 †:

…15†

Hence, it turns out that the geometry correction is approximately equal to the di€erence between the computed electronic contributions at the non-optimum and optimum geometries. From the data presented in Tables 1±4, we ®nd that the agreement between the two sets of numbers is good at all levels of theory. The discrepancy tends to widen for some components of the c tensor. This may be partly due to some numerical inaccuracies in the property derivatives calculated using the ®nite di€erence method. Another explanation may follow from pragmatic assumptions of the present theory: the

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

9

present formulae for the geometry correction are obtained under the assumption of the harmonic approximation. An extension to include anharmonicities in the potential energy may provide better agreement. The result expressed by Eq. (15) is not surprising if one recalls that the exactly calculated total property value in the given vibrational state must be the same for any choice of the reference geometry. Thus, its partition into electronic and vibrational contributions is completely arbitrary. However, with the vibrational contribution calculated by means of the ®nite order perturbation theory this invariance is lost. The given order in the perturbation approach and the level of anharmonicity which is taken into account do not represent exactly the same level of approximation at di€erent reference geometries. The di€erence in the electronic contribution is, within the assumed approximation for the calculation of the property surface, computed exactly. The geometry correction term follows from the lowest-order approximation in vibrational theory. One should bear in mind that the magnitude of the corrective term is primarily an indication of the validity of a chosen non-optimum geometry. Thus, the use of more elaborate formulae for pgc is not imperative. For values of pgc , which are small compared to the total vibrational correction, the present perturbation expressions may be employed to estimate the total vibrational contribution which accounts for the use of the non-optimum molecular geometry. A relatively small value of the pgc term means that the current reference geometry and the true equilibrium structure predicted by the given computational method are quite close to each other. This is usually the case for MP2 geometries used in higher-level calculations of vibrational corrections to molecular electric properties. The present study supplements our earlier exploratory calculations and demonstrates how vibrational corrections to electric properties of any molecule may be computed at non-optimum geometries at relatively low cost. In the low-order approximation recommended in this paper, the corrective term is easy to calculate and provides useful information concerning the appropriateness of the assumed reference geometry of the molecule. The present method avoids the use of a multitude of di€erent reference geometries and permits the identi®cation of the meaning of the vibrational correction at the selected geometry. From this point of view, our extension of the Bishop±Kirtman perturbation scheme becomes a useful tool in computation of vibrationally corrected molecular electric properties. If one takes into account that the geometry optimization at coupled-cluster levels of theory is quite time consuming and costly, then the use of a single reference geometry becomes quite an advantage. The same conclusion applies to calculations in which the geometry data is obtained using a standard basis set, e.g. from the G A U S S I A N library, and the property data is evaluated with a larger, speci®cally tailored basis set. Finally, let us add that the theory presented and exempli®ed in this paper deals with static electric properties. How much of the vibrational correction survives at non-zero frequencies depends on both the frequency range and the character of the studied frequency-dependent property. For predominantly used frequency ranges which go far beyond the range of vibrational frequencies, the pure vibrational contributions tend to become quite small and negligible in comparison with the pure electronic contribution [15,20]. Then, the choice of reference geometry and its contribution to the vibrational correction become of secondary importance.

Acknowledgements The authors acknowledge the ®nancial support obtained via the TMR Network, contract no. FMRXCT96-0047. One of us (A.J.S) acknowledges a partial ®nancial support from the Swedish Science Research Council (NFR, Sweden, contract no. K-AA/KU 01624-311) and the Scienti®c Research Committee (KBN, Poland, contract no. 3 T09A 108 15).

10

V.E. Ingamells et al. / Chemical Physics 260 (2000) 1±10

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12]

[13] [14] [15] [16] [17] [18] [19] [20]

A.D. Buckingham, J. Chem. Phys. 36 (1962) 3096. C.J. Jameson, Bull. Magn. Res. 3 (1980) 3. P.W. Fowler, G. Riley, W.T. Raynes, Mol. Phys. 42 (1981) 1463. P. Lazzaretti, R. Zanasi, A.J. Sadlej, W.T. Raynes, Mol. Phys. 62 (1987) 605. P.K.K. Pandey, D.P. Santry, J. Chem. Phys. 73 (1980) 2899. D.S. Elliott, J.F. Ward, Mol. Phys. 51 (1984) 45. D.M. Bishop, Rev. Mod. Phys. 62 (1990) 343. D.M. Bishop, Adv. Chem. Phys. 104 (1998) 1. D.M. Bishop, J.M. Luis, B. Kirtman, J. Chem. Phys. 108 (1998) 10013. V.E. Ingamells, M.G. Papadopoulos, A.J. Sadlej, J. Chem. Phys. 112 (2000) 1645. A.D. Buckingham, Adv. Chem. Phys. 12 (1967) 107. M.J. Frisch, G.W. Trucks, H.B. Schlegel, P.M.W. Gill, B.G. Johnson, M.A. Robb, J.R. Cheeseman, T.A. Keith, G.A. Petersson, J.A. Montgomery, K. Raghavachari, M.A. Al-Laham, V.G. Zakrewski, J.V. Ortiz, J.B. Foresman, J. Cioslowski, B.B. Stefanov, A. Nanayakkara, M. Challacombe, C.Y. Peng, P.Y. Ayala, W. Chen, M.W. Wong, J.L. Andres, E.S. Replogle, R. Gomperts, R.L. Martin, D.J. Fox, J.S. Binkley, D.J. Defrees, J. Baker, J.J.P. Stewart, M. Head-Gordon, C. Gonzales, J.A. Pople, G A U S S I A N 9 4 (Rev. B2 and C3), Gaussian, Pittsburgh, PA, 1995. A.J. Sadlej, Collect. Czech. Chem. Commun. 53 (1988) 1995. T. Pluta, A.J. Sadlej, Chem. Phys. Lett. 297 (1998) 391. V.E. Ingamells, S.G. Raptis, M.G. Papadopoulos, Chem. Phys. Lett. 307 (1999) 484. H. Reis, S.G. Raptis, M.G. Papadopoulos, submitted for publication. D.M. Bishop, B. Kirtman, H.A. Kurtz, J.E. Rice, J. Chem. Phys. 98 (1993) 8024. R.D. Amos, I.L. Alberts, J.S. Andrews, S.M. Colwell, N.C. Handy, D. Jayatilaka, P.J. Knowles, R. Kobayashi, G.L. Laming, A.M. Lee, P.E. Maslen, C.W. Murray, P. Palmieri, J.E. Rice, E.D. Simandiras, A.J. Stone, M.D. Su, D.J. Tozer, C A D P A C 6 . 0 , The Cambridge Analytic Derivatives Package, Cambridge, UK, 1995. V.E. Ingamells, N U M D E R , a numerical derivative algorithm, National Hellenic Research Foundation, Athens, Greece, 2000. D.M. Bishop, E.K. Dalskov, J. Chem. Phys. 104 (1996) 1004.