- Email: [email protected]
Polarizable continuum model calculations including electron correlation in the ab initio wavefunction
Journal of Molecular Structure (Theochem) , 279 (1993) 223-23 I Elsevier Science Publishers B.V., Amsterdam
223
Polarizable continuum model ca .culations including electron correlation in the ab initio wave ‘unction F.J. Olivares de1 Valle, M.A. Aguilar and S. Tolosa Departamento de Quimica Fisica, Universidad de Extremadura, 06071 Badajoz (Spain) (Received 11 May 1992)
Abstract Some calculations are presented using a recently developed version of the ab initio polarizable continuum model that includes electron correlation in the representation of the solute. Calculations were carried out to estimate the self-consistent field electronic solvation energies and one-electron properties of several molecules. Comparison is made of results at several levels of approximation.
Introduction
In the current methodology of quantum mechanics, one topic of particular interest because of possible applications in the fields of chemistry and biology is the study of molecular systems affected by the presence of a solvent. It is common to use a continuous polarizable solvent, ignoring its microscopic structure and possible short-range influences (e.g. charge transfer, charge exchange) on the solute [l-3]. This is an extremely simple model, and although not giving a satisfactory response for some thermodynamic quantities (e.g. ionic hydration energies), it has proved useful in explaining a large number of phenomena associated with solvation in both equilibrium and non-equilibrium situations. Its simplicity in describing the solvent allows one to focus on studying the solute, employing techniques (Moller-Plesset (MP), configuration interaction (CI), multiconfiguration selfconsistent field (MCSCF)) that would be prohibi-
Correspondence to: F.J. Olivares de1 Valle, Departamento de Quimica Fisica, Universidad de Extremadura, 06071 Badajoz, Spain.
tively difficult if other models in which the solvent is represented in greater detail were used. One fairly extensively used model is the ab initio polarizable continuum model (PCM), developed initially by Tomasi and co-workers [4]. This model (from now on denoted MST) considers the total solute-created electrostatic potential in an exact form, i.e. with no necessity to truncate. It also allows one to work with solutes immersed in arbitrarily shaped cavities, and takes into account the effect of the presence of the solvent on the solute wavefunction. The model has been modified over the last few years with respect to the construction of the cavity [5-81, the nature of the interaction that one wishes to represent [9-131, and the description of the solute itself [ 14,151. In this last extension, our group presented a perturbational version of the model which takes into account the effects of electron correlation in the solvation process at two levels: (a) at the energetic level, computing the electron correlation energy with a Hartree-Fock solvent reaction potential (PTE approximation); (b) at the density matrix level, including the electron correlation modification of the first-order solute density matrix and calculating the energy up to first order (PTD approximation) or higher
224
F.J. Olivares de1 Valle et al./J. Mol. Strut.
orders (PTDE approximation). We have developed the PTDE model in the framework of the perturbation-independent response density matrix, which includes (via the coupled-perturbed HartreeFock equations) all effects of orbital relaxation due to electron correlation. We also use the response density matrix to compute some induced one-electron properties in solution. In this paper, we present some comparative results relative to the effect that considering the correlation correction in the wavefunction has on the calculated properties at the different levels of approximation. The paper is organized as follows. First, details of the method of calculating the components of the solvation energy are given, including electron correlation in the density matrix. Then the results for the energies and one-electron properties are presented, with mention of the most relevant conclusions. Methods and approximations The model adopted in the present work to represent the solute lies within the context of the models known as second generation or quantum mechanical versions of the PCM [16]. This method makes explicit use of the Schrddinger equation
where l?” is the Hamiltonian for the non-perturbed solute and I& represents the effective solutesolvent interaction potential. Details of the method can be found in the original work [4]. We have recently developed [14] a method of treating the effects of electron correlation in solvated molecules: the MST method for a cavity of general shape was reformulated by applying standard Rayleigh-Schriidinger perturbation theory using the Fock operator as the zeroth-order Hamiltonian and the Slater determinant of SCF spin orbitals as the zeroth-order wavefunction. The incorporation of the post-SCF procedure in the continuum model could introduce significant changes in the results: the basic units of the independent electron model theory are modified in a solution by
(Theochem)
279 (1993) 223-231
the solvent reaction field, and the electric field of the solute depends on the charge distribution, which is in turn influenced by electron correlation. In order to break down the mutual effects of correlation and polarization into component parts, we have developed the theory at three separate levels. At the first level designated as PTE, the internal energy of the solute is computed by the expansion (1) in vacua, and by Esor(PTE) = f Ek{@ k=O
(2)
in the presence of the dielectric medium. Here 9 and D denote the first-order density matrix in vacua and in a solution respectively, and (521 (or {D}) the dependence of each internal energy component on this matrix. The suffix k refers to the order of the perturbation and m to the maximum order of the expansion (in our case m < 3). The Hartree-Fock energy is obtained when m= 1: ~vAc(SCF) = &‘o{% + &‘,(%
(3)
Eso,(SCF) = Eo {D) + E, {D1
(4)
In vacua, the zeroth-order component may be decomposed into two contributions: the zerothorder zero-body contribution PO associated with the nuclear repulsion energy d,,,, and the zerothorder one-body contribution, associated with the sum of the one-electron orbital energies &‘A{&@} = c. ,n,q .( 9 } (n, IS ’ t h e occupation number for the ith molecular orbital). We include (91 in the notation of this last component because the density matrix is used in its evaluation. The first-order two-body contribution is 8, {9} = -+&9ijij (where the $yij = (zj ]]ij) represent the antisymmetrized twoelectron molecular integrals for all the ij occupied spin orbital pairs). The SCF internal energy in eqn. (4) may be decomposed in the same form as &vVAC(SCF), except that in this case the different components are modified by the terms derived from the perturbing
F.J. Olivares de1 Valle et a/./J. Mol. Struct.
(Theochem)
potential V;,, . Thus, Ei = 6: + Qqn, where Qqn represents the electrostatic interaction between the solute nuclei and the virtual charge distribution (or image) on the cavity’s inner surface. This term may be thought of as the nuclear contribution to the solvation energy. Two terms may be associated with the modification of the zeroth-order one-body contribution by the effect of the solvent. One term is not present in the vacuum and is directly related to the solvent reaction field and the perturbed density matrix (see ref. 4): (DFNT) (where V,, is the matrix with elements (x,1 I&IX”), x being a one-electron basis function). The other term derives from using the perturbed density matrix D in the calculation of si. The averaged value (DP&,) may be thought of as the electron-static contribution to the zeroth-order one-body component. The modification of the first-order twobody contribution by the effect of the reaction field is only of the electron-polarization type and depends exclusively on the change in the density matrix. The higher order components (m = 2,3) in eqns. (1) and (2) are associated with the electron correlation energy. In vacua, the Msller-Plesset secondorder contribution may be written as (5) where %ij&= &&#= (ij Ilab)(&i + Ej - E, - &b)-‘, and the occupied spin orbitals are labelled with the hole indices i, j, k,. . . and unoccupied spin orbitals with the particle indices a, b, c, . . . (later, p, q, r, . . . are used as general indices). The third-order contribution to eqn. (1) may be decomposed into three parts: g3{9}
=
SC
+
1 ij abed
9&zbxubcd1cdij
6 11 %job ijkl ab
$k,q
-
225
279 (1993) 223-231
cc 2ijabfkbic90ckj ijk abc
up to first-order terms in the expansions: (7) G~L(PTD) = JWV
(8)
Here the notation {DC} emphasizes the dependence of the solute internal energy components E, and E, on the correlated density matrix to order u (2 2) over the molecular basis via the changes in the fi,, potential: the inclusion of electron correlation effects in the wavefunction produces changes in the electric field created by the solute in its environment which can induce non-negligible changes in I&,. The correlated density matrix to a predetermined order may be obtained (in the molecular orbital basis) by adding to the SCF density matrix, @cr = 6, (Vi,j), the total ath density matrix corrections: gTo= E cgn 0
= gscr + @oRR
(9)
To second order, gcoRR E G@? and the occupiedoccupied blocks are given [17,18] by {gZ)ij = 41
(10)
sjkab9ikob
kab
Because .Ypqrsare integrals over the spin molecular orbital basis, in practical calculations it is preferable to make use of transformed spin-free expressions over the spatial molecular orbitals basis. Then eqn. (10) can be written as G% >I./= - 1 WK4B(2%K4s - %KBA) KAB +
-zJKBA(2%KBA
-
=%KAB)l
(11)
where the integrals xlJAB and sIJAs are defined in a similar fashion to $&$, and L!&&. The elements of the virtual-virtual block may likewise be written as (92)~~
=
1 [~IJACWJBC CIJ +
=%bkl
At the second level of approximation, designated as perturbation theory on the density distribution alone (PTD), we compute the solute-solvent interaction potential J&r, taking into account correlation effects in the wavefunction and including only
+ E, {D,)
-
%JCA(2=%CB
-
~IJCB>
-%JBC)I
(12)
The corresponding elements of the occupiedvirtual block are obtained from {g2jIJ and {&}AB by solving the perturbation-independent linear equation (see for example refs. 17 and 18) @@B.&J
-
EB)ABJ,A,
=
xAI
(13:
226
F.J. Olivares del Valle et ral./J. Mol. Struct.
where ABJ,AIis given by A AI,BJ
=
1 -
c2$AIBJ
-
9ABJI
-
xAJBI)(EI
-
&A)-’
(14) The {XA,}are intermediate quantities made up of molecular integrals and orbital energies: xA,
=
[ll
+
121+ L31+ L41
(1%
where KJ
PI E C(~*}B~~B~CA
-
(17)
XBlac
BC
L31= - c
[y,BKJ(29KJAB
-
dKJBA)
JKB
-
L41= 1
9,BJK(9KJAB
+
[9CBAJ(2%JCB
%UBA)l
-
(18)
=%JBC!)
JBC
-
yCBJA(9,CB
+
=%JBC)l
(19)
Once gz (or D, in the solution case) has been constructed, the matrix gToT (or DToT) is readily transformed to the atomic orbital basis (in our case by obtaining the natural orbitals and the new occupation numbers) and can be used for the calculation of I&T, thus introducing the electron correlation correction into the primary aspects of the PCM. Thus, with respect to the original version (eqn. (4)), the PTD approximation (eqn. (8)) supposes a modification of the results which has its origin in the different electron distribution (represented by its density matrix) used to calculate the virtual electron distribution CJ(S)generated on the cavity surface. Finally, at the third level of approximation (PTDE) the internal energy of the molecule in solution is
are the free energy terms Cp,2Jw4~ relating to the variation of the correlation energy (to order p) including the new potential I$.,, calculated from the correlated density matrix DToT . In
where
(Theochem)
279 (1993) 223-231
sum, the correlated density matrix refers to the situation in solution, taking into account the mutual effects of the solvent reaction potential and the fluctuation potential. In both the PTD and PTDE approximations, the matrix DToT is calculated iteratively: starting from the density matrix in vacua, one calculates the solute-solvent interaction V,NT; this is then used in resolving the Schrodinger equation, when one obtains a new density matrix (in the MP context) that includes the modifications due to correlation and to the effect of the dielectric, and allows the cycle to be re-initiated. Four or five cycles are usually enough to achieve convergence in the energy. At this point it must be made clear that, although I&T includes all the components of the solutesolvent interaction (electrostatic, dispersion, exchange etc.), only the electrostatic component is considered in the present work. This approximation is made with a view to facilitating the interpretation of the results and is based on the results of earlier work which analysed the influence of the dispersion component [19]. This component has a considerable influence on the energy aspect. However, it has only a very small effect on the density matrix of the solute, and hence on its oneelectron properties. Once a model has been adopted to represent the solute and the solvent, one has to specify the interactions included in the interaction potential &,T. In general, for a non-microscopic description of the solute, we use a multipole expansion of the potential. This is the case for the historical model. The expansion must be truncated at some fixed order, accepting a compromise between the quality of the results and the complexity of the method. Reference 16 gives a brief selection of the methods relating to the classical continuum electrostatic approach which use a multipole expansion for version. &NT ; ref. 20 gives a more sophisticated In our case, the potential c,, is not obtained from an analytical parametrix expression, but rather by a numerical procedure which is iterative and convergent (for more details, see ref. 4). By making use of eqns. (1) and (2) or eqn. (20)
F.J. Olivares de1 Valle et a/./J. Mol. Struct.
(Theochem)
279 (1993) 223-231
227
and its vacuum equivalent, one may obtain the variation AE,, of the solute internal energy at any of the three levels of approximation, together with its corresponding components (AE,, is a measure of the electrostatic contribution to the solvation energy). However, this quantity is not employed directly [21] because it is not an energy of primary interest. The magnitude (which is well defined thermodynamically) that interests us is
We compute AG,, using eqn. (24). It may be expressed [ 13,141 as
AG,, = AE,, - N&n +)
(21)
which is the free energy of formation of the set of charges off in the cavity of the dielectric. Briefly, the standard free energy of solvation for M, AGL, may be estimated by means of the expression AG; = G;’ - GM = A(PV) + AGF= + Aq
(25) with AG:’ = AGO”i- AG;
(26)
AGO”= +Q4,,
(27)
AG,’= xAei - $(DP&T)
(28)
AGF’= - +cA-“iiu G
(29)
Finally, the contributions AGT’and AG$’ are given by (see eqns. (5) and (6)) (30)
(22)
where the solute-solvent interaction energy [ 1l] AGET is (omitting the contribution from dispersion
AGP’ = AE3”+ AE3”+ AE,4
[191) AGENT = Gw + AG,,
Results and discussion
(23)
The term AG,, is the change in the internal electronic energy of M as it passes from the gaseous phase into solution (we consider that at the temperature T of our conditions, all molecules are in the ground state) AG,, = G;;’ - E,i””
(31)
(24)
Included in Giy’ are the proper solute electronic energy itself, the energy coming from the electrostatic solute-solvent interaction (sum of two contributions, one nuclear, the other electronic), and the polarization energy of the solvent because of the solute. In eqn. (22), the term Aq contains the entropy or thermal contributions associated with the change in the rotational, translational and vibrational energies under the effect of the solvent [13], and AGZ in addition to depending on temperature, is also affected by the molecular geometry in so far as this has an influence on AGENT.In this paper we assume that there is no geometrical distortion on passing from the gaseous to the condensed phase [15,22].
The solutes that we used to analyse the changes produced in the free energy on taking into consideration correlation in the wavefunction are H,CO, HFCO, F&O, F&S, H,CS, HClCO and FClCO. Their geometries are given in Table 1. The basis set used in al the calculations was 6-3 1G and the dielectric constant was E = 78.3 (value equivalent to that of water at 25°C). . The cavities were constructed using a set of interpenetrating spheres of variable radius according to the procedure proposed previously by our group [5-81, which takes into account the population supported by each atom and the basis set that is being employed. The radii used in each case are given in Table 2. The use of this procedure permits the generation of cavities adapted to the local characteristics of each atom in each molecule. This can be seen, for example, in the case of the carbon atom for the series H,CO, HFCO and F&O. This atom undergoes a reduction of the electron population as the hydrogen atoms are progressively substituted by their corresponding fluorine atoms. Thus, in this series, the value of the net charge on the carbon
228
F.J. Olivares de1 Valle et al./J. Mol. Struct.
(Theochem)
279 (1993) 223-231
TABLE 1 Geometrical
parameters”
Parameter
H,COb
HFCO’
F, COd
HClCO”
FClCO’
1.182 1.090
1.173
F,CSp
H,CSh
Bond length 1.210 1.116
R co R HC R FC
1.181 1.095 1.338
_ _
R CC, R CS
Bond angle L HCH L HCF L FCF L HCCl L FCC1 L HCO L HCS L FCO L FCS L ClCO
_ _
116.5 _ _
_ 109.9 _
_ 121.7
1.315 _ _
_ _
_
_
107.6 _
1.087 1.334 1.765 _
_ _
110.5
_ 108.8
127.3
_ _ _
“Bond lengths in kngstrsms;
1.170
122.8 _
_ _
126.5 _
126.2 _ _
123.1
123.7 _ 127.5
1.315 1.725 1.589 _ _ 107.1 _ _ _ _ 126.5
_ 1.611 116.5 _ _ _ _ 121.8 _ _
bond angles in degrees. bRef. 23. ‘Ref. 24. dRef. 25. ‘Ref. 26. ‘Ref. 27. gRef. 28. hRef. 29.
atom changes from + 0.14e in H&O to 1.14e in F$O. This fact is clearly reflected in the value that is used for R,, which goes from 1.9 A in the first case to 1.5 I$ in the second (Table 2). The use of
fixed radii (equal to the van der Waals radii, RVDW, or scaled up by a certain quantity k x RVDW)takes no account of this fact, each carbon atom in any of the solutes under study contributing with an
TABLE 2 Charges and radii”,b Parameter
H&O
HFCO
F,CO
HClCO
FClCO
0.14 -0.58 0.21
0.63 - 0.53 0.32
1.14 - 0.45 _ _
0.57 - 0.42 _ _
- 0.42
-0.34 _
0.120 - 0.46 0.33 -i _
F&S
H,CS
Charges :: Q” z:
Q,,
0.00
-0.35 0.20
0.53 _
- 0.43 _
_ 0.13 - 0.33 _
0.26 - 0.03 _
1.71 _
2.13 _
Radius Rc & & & R, Rc, “In gngstrsms. bObtained with the PTE model.
1.70 1.90 0.90 _
1.50 1.89 _ _
1.88 1.88 0.92 _
1.70 1.87 _ _
1.80 _
1.75 _
2.77
1.74 2.70
2.79 1.74 _
0.95 2.88
3
- 10.89
- 7.33
AGsc,
AG,, 3.40
-7.11
- 10.51
- 0.59
3.99
- 144.10
133.60
H,COC
3.55
- 10.26
- 13.81
- 0.79
4.33
-98.22
84.42
HFCOb
3.28
- 9.77
- 13.05
- 0.73
4.02
-94.60
81.50
HFCO”
1.21
- 4.83
- 6.04
- 0.77
1.98
54.29
- 60.33
F,COb
“In kcalories per mol. bObtained using the PTE model. ‘Obtained
3.56
- 0.82
AG,
AAGd
4.17
AG,
139.56
- 150.45
AG,
AGo
H,COb
Variations in free energiesa
TABLE
0.66
- 0.77
- 1.43
0.36
0.30
27.10
- 28.54
HCICOb
1.07
- 0.56
- 1.63
0.53
0.54
75.84
- 77.47
HClCO’
1.91
- 2.27
-4.18
-0.12
2.03
- 78.65
74.47
FCICOb
using the PTD model. dAAG = AG,,,
1.09
-4.39
- 5.48
-0.77
1.86
47.38
- 52.57
F&O’
- AG,.
1.93
- 2.28
- 4.21
-0.13
2.05
- 78.38
74.18
FClCO’
2.16
- 5.94
-8.10
- 0.65
2.81
- 72.89
64.79
F,CSb
1.99
- 5.82
-7.81
-0.51
2.50
- 72.30
64.49
F,CS
0.72
- 1.97
- 2.69
- 0.30
1.02
27.78
- 30.47
H,CSb
0.55
- 1.67
2.22
- 0.22
0.77
27.43
- 29.65
H,CS
230 TABLE
F.J. Olivares de1 Valle et al./J. Mol. Strut.
(Theochem)
279 (1993) 223-231
4
Variation in molecular properties”,b AP”
H,COd
AQc AQo AQH AQs AQF AQ,,
- 0.009
0.045
0.023
0.120
0.070
-0.134 0.067
- 0.225 0.085
-0.106 0.117
-0.147 0.104
- 0.053
-0.145
- 0.033
- 0.078
- 0.001
- 0.060
H,CO”
HFCOd
HFCO”
F,COd
F,CO” 0.266
MC
0.004
0.009
0.011
0.200
0.025
0.056
A&
0.106
- 0.037
- 0.029
- 0.074
- 0.009
- 0.069
A&
0.030
0.023
0.061
0.055
AA A&
- 0.006
0.200
0.004
0.024
F,CSd 0.018
F&Se 0.154
- 0.019
- 0.014
0.000
- 0.065
0.007
- 0.002 0.002
0.035
Ap
0.442
0.679
0.354
0.426
0.089
0.149
0.038
0.793
2.086
0.085
0.920
- 0.364
0.079
- 0.210
“In atomic units. bAll values represent the electronic contribution variation. “AP = P,, - P,,.
HClCO
0.011
0.078
-0.071 0.093
-0.162 0.142
- 0.043
- 0.058
0.008
0.001
- 0.019
- 0.098
0.047
0.047
- 0.005
0.508
0.287
0.531
- 0.079
1.071
- 0.235 - 0.059
A&,
A
HCICOd
-0.189 0.936
except for the dipole moment which represents the total
dObtained using the PTE model. “Obtained using the PTD model.
identical sphere to the resultant cavity. The above can be generalized to each of the atoms involved. A more details discussion of these questions may be found in refs. 5-8. Table 3 lists the components of the free energy of solvation that we obtained using the PTE and PTD approximations, calculated according to eqn. (26) and (29)-(31). Without going into detail (which we postpone until a future article dealing with a greater number of diverse solutes) we can draw the following general conclusions: (1). In no case are the third-order contributions to the free energy of solvation especially significant, and they hardly vary from one approximation to another. One must not ignore, however, the fact that AG, is a compensated sum of quite diverse contributions (from two, three and four bodies), and some of them individually can undergo major changes. An analysis of these variations can be found in ref. 8. It is therefore enough to consider up to the contribution AG2 in order to extract a value near to the convergence of AGioL. As a general rule, AG3 is a stabilizing contribution tending to raise the absolute value of AGioL. (2). The second-order contribution to the free
energy of solvation, AG2, is quite significant in the solutes we studied, reaching a value near 45% in the case of H,CS. In every case, the contribution tends to reduce the absolute value of AGzoL. It is therefore advisable to consider this contribution in those calculations of solvation which handle free energies. Its importance is possibly critical in the study of reactive processes in which the saturated nature of the initial and final products is different, as the value of AGz (and, in some cases, AG,) is usually influenced by the degree of saturation of the molecule under study. We have also found a certain correlation between the variations of AG2 and the polarizability of the solute [14]. (3). The third conclusion worthy of note is, in our judgement, the modest significance of including electron correlation in the wavefunction. The variations that appear in the value of AGioL are less than 0.5 kcal mall’ , smaller than the margin of accuracy for this type of model. We found, in a more complete sample of solutes, that this behaviour is general and always an order of magnitude below that found using the PTE approximation. Finally, in Table 4 we list the results for the variation of a selection of one-electron properties.
F.J. Olivares de1 Valle et al/J.
Mol. Struct.
(Theochem)
279 (1993) 223-231
In general, one can say that the inclusion of correlation in the density matrix leads to variations which are more significant in this class of properties than in the energies. Thus, one appreciates an increase in the net charge on the atoms and in the size of the electron distribution, and, lastly, an accentuated polarity of the molecule which can reach major variations, as occurs in the solute HClCO, whose dipole moment varies by 0.7 D in the PTE approximation and 1.33 D in the PTD. ,
10 11 12 13 14 15 16
Acknowledgement
17
This research was partially sponsored by the Direction General de Investigation Cientifica y Tknica, Projects 1017/87 and 930/90.
18 19 20
References M. Born, Z. Physik., 1 (1920) 45. L. Onsager, J. Am. Chem. Sot., 58 (1936) 1486. J.G. Kirkwood, J. Chem. Phys., 2 (1934) 351. S. MiertuS, E. Scrocco and J. Tomasi, Chem. Phys., 55 (1981) 117. M.A. Aguilar, M.A. Martin, S. Tolosa and F.J. Olivares de1 Valle, J. Mol. Struct. (Theochem), 166 (1988) 313. M.A. Aguilar and F.J. Olivares de1 Valle, Chem. Phys., 129 (1989) 439. F.J. Olivares de1 Valle, M.A. Aguilar, S. Tolosa, J.C. Contador and J. Tomasi, Chem. Phys., 143 (1990) 371. F.J. Olivares de1 Valle and M.A. Aguilar, J. Comput. Chem., 13 (1992) 115. D. Rinaldi, B.J.C. Cabral and J.L. Rivail, Chem. Phys. Lett., 125 (1986) 495.
21
22 23 24 25 26 27 28 29
231
F. Floris and J. Tomasi, J. Comput. Chem., 10 (1989) 616. J. Langlet, P. Claverie, J. Caillet and A. Pullman, J. Phys. Chem., 92 (1988) 1617. M.A. Aguilar and F.J. Olivares de1 Valle, Chem. Phys., 138 (1989) 327. F.J. Olivares de1 Valle and M.A. Aguilar, J. Mol. Struct. (Theochem), 280 (1993) 25. F.J. Olivares de1 Valle and J. Tomasi. Chem. Phys., 150 (1991) 139. M.A. Aguilar, F.J. Olivares de1 Valle and J. Tomasi, Chem. Phys., 150 (1991) 151. J. Tomasi, R. Bonaccorsi, R. Cammi and F.J. Olivares de1 Valle, J. Mol. Struct. (Theochem), 234 (1991) 401. G.W. Truck, E.A. Salter, C. Sosa and R.J. Bartlett, Chem. Phys. Lett., 147 (1988) 359. G.W. Truck, E.A. Salter, J. Noga and R.J. Bartlett, Chem. Phys. Lett., 150 (1988) 37. F.J. Olivares de1 Valle and M.A. Aguilar, in preparation. M.E. Davis and J.A. McCammon, Chem. Rev., 90 (1990) 509. R. Bonaccorsi, C. Ghio and J. Tomasi, in R. Carbo (Ed.), Current Aspects of Quantum Chemistry, 1981, Elsevier, 1982, pp. 407-426. R. Bonaccorsi, R. Cammi and J. Tomasi, J. Comput. Chem., 12 (1991) 301. K. Tagaki and T. Oka, J. Chem. Sot. Jpn., 18 (1974) 1174. R.F. Miller and R.F. Curl, J. Chem. Phys., 34 (1961) 1847. V.W. Laurie, D.T. Pence and R.H. Jackson, J. Chem. Phys., 37 (1962) 2995. H. Takeo and C. Matsumura, J. Chem. Phys., 64 (1976) 4536. H. Oberhammer, J. Chem. Phys., 73 (1980) 4310. A.J. Careless, H.W. Kroto and B.M. Landsberg, Chem. Phys., 1 (1973) 371. P.H. Turner, L. Halonen and I.M. Mills, J. Mol. Spectrosc., 88 (1981) 402.
223
Polarizable continuum model ca .culations including electron correlation in the ab initio wave ‘unction F.J. Olivares de1 Valle, M.A. Aguilar and S. Tolosa Departamento de Quimica Fisica, Universidad de Extremadura, 06071 Badajoz (Spain) (Received 11 May 1992)
Abstract Some calculations are presented using a recently developed version of the ab initio polarizable continuum model that includes electron correlation in the representation of the solute. Calculations were carried out to estimate the self-consistent field electronic solvation energies and one-electron properties of several molecules. Comparison is made of results at several levels of approximation.
Introduction
In the current methodology of quantum mechanics, one topic of particular interest because of possible applications in the fields of chemistry and biology is the study of molecular systems affected by the presence of a solvent. It is common to use a continuous polarizable solvent, ignoring its microscopic structure and possible short-range influences (e.g. charge transfer, charge exchange) on the solute [l-3]. This is an extremely simple model, and although not giving a satisfactory response for some thermodynamic quantities (e.g. ionic hydration energies), it has proved useful in explaining a large number of phenomena associated with solvation in both equilibrium and non-equilibrium situations. Its simplicity in describing the solvent allows one to focus on studying the solute, employing techniques (Moller-Plesset (MP), configuration interaction (CI), multiconfiguration selfconsistent field (MCSCF)) that would be prohibi-
Correspondence to: F.J. Olivares de1 Valle, Departamento de Quimica Fisica, Universidad de Extremadura, 06071 Badajoz, Spain.
tively difficult if other models in which the solvent is represented in greater detail were used. One fairly extensively used model is the ab initio polarizable continuum model (PCM), developed initially by Tomasi and co-workers [4]. This model (from now on denoted MST) considers the total solute-created electrostatic potential in an exact form, i.e. with no necessity to truncate. It also allows one to work with solutes immersed in arbitrarily shaped cavities, and takes into account the effect of the presence of the solvent on the solute wavefunction. The model has been modified over the last few years with respect to the construction of the cavity [5-81, the nature of the interaction that one wishes to represent [9-131, and the description of the solute itself [ 14,151. In this last extension, our group presented a perturbational version of the model which takes into account the effects of electron correlation in the solvation process at two levels: (a) at the energetic level, computing the electron correlation energy with a Hartree-Fock solvent reaction potential (PTE approximation); (b) at the density matrix level, including the electron correlation modification of the first-order solute density matrix and calculating the energy up to first order (PTD approximation) or higher
224
F.J. Olivares de1 Valle et al./J. Mol. Strut.
orders (PTDE approximation). We have developed the PTDE model in the framework of the perturbation-independent response density matrix, which includes (via the coupled-perturbed HartreeFock equations) all effects of orbital relaxation due to electron correlation. We also use the response density matrix to compute some induced one-electron properties in solution. In this paper, we present some comparative results relative to the effect that considering the correlation correction in the wavefunction has on the calculated properties at the different levels of approximation. The paper is organized as follows. First, details of the method of calculating the components of the solvation energy are given, including electron correlation in the density matrix. Then the results for the energies and one-electron properties are presented, with mention of the most relevant conclusions. Methods and approximations The model adopted in the present work to represent the solute lies within the context of the models known as second generation or quantum mechanical versions of the PCM [16]. This method makes explicit use of the Schrddinger equation
where l?” is the Hamiltonian for the non-perturbed solute and I& represents the effective solutesolvent interaction potential. Details of the method can be found in the original work [4]. We have recently developed [14] a method of treating the effects of electron correlation in solvated molecules: the MST method for a cavity of general shape was reformulated by applying standard Rayleigh-Schriidinger perturbation theory using the Fock operator as the zeroth-order Hamiltonian and the Slater determinant of SCF spin orbitals as the zeroth-order wavefunction. The incorporation of the post-SCF procedure in the continuum model could introduce significant changes in the results: the basic units of the independent electron model theory are modified in a solution by
(Theochem)
279 (1993) 223-231
the solvent reaction field, and the electric field of the solute depends on the charge distribution, which is in turn influenced by electron correlation. In order to break down the mutual effects of correlation and polarization into component parts, we have developed the theory at three separate levels. At the first level designated as PTE, the internal energy of the solute is computed by the expansion (1) in vacua, and by Esor(PTE) = f Ek{@ k=O
(2)
in the presence of the dielectric medium. Here 9 and D denote the first-order density matrix in vacua and in a solution respectively, and (521 (or {D}) the dependence of each internal energy component on this matrix. The suffix k refers to the order of the perturbation and m to the maximum order of the expansion (in our case m < 3). The Hartree-Fock energy is obtained when m= 1: ~vAc(SCF) = &‘o{% + &‘,(%
(3)
Eso,(SCF) = Eo {D) + E, {D1
(4)
In vacua, the zeroth-order component may be decomposed into two contributions: the zerothorder zero-body contribution PO associated with the nuclear repulsion energy d,,,, and the zerothorder one-body contribution, associated with the sum of the one-electron orbital energies &‘A{&@} = c. ,n,q .( 9 } (n, IS ’ t h e occupation number for the ith molecular orbital). We include (91 in the notation of this last component because the density matrix is used in its evaluation. The first-order two-body contribution is 8, {9} = -+&9ijij (where the $yij = (zj ]]ij) represent the antisymmetrized twoelectron molecular integrals for all the ij occupied spin orbital pairs). The SCF internal energy in eqn. (4) may be decomposed in the same form as &vVAC(SCF), except that in this case the different components are modified by the terms derived from the perturbing
F.J. Olivares de1 Valle et a/./J. Mol. Struct.
(Theochem)
potential V;,, . Thus, Ei = 6: + Qqn, where Qqn represents the electrostatic interaction between the solute nuclei and the virtual charge distribution (or image) on the cavity’s inner surface. This term may be thought of as the nuclear contribution to the solvation energy. Two terms may be associated with the modification of the zeroth-order one-body contribution by the effect of the solvent. One term is not present in the vacuum and is directly related to the solvent reaction field and the perturbed density matrix (see ref. 4): (DFNT) (where V,, is the matrix with elements (x,1 I&IX”), x being a one-electron basis function). The other term derives from using the perturbed density matrix D in the calculation of si. The averaged value (DP&,) may be thought of as the electron-static contribution to the zeroth-order one-body component. The modification of the first-order twobody contribution by the effect of the reaction field is only of the electron-polarization type and depends exclusively on the change in the density matrix. The higher order components (m = 2,3) in eqns. (1) and (2) are associated with the electron correlation energy. In vacua, the Msller-Plesset secondorder contribution may be written as (5) where %ij&= &&#= (ij Ilab)(&i + Ej - E, - &b)-‘, and the occupied spin orbitals are labelled with the hole indices i, j, k,. . . and unoccupied spin orbitals with the particle indices a, b, c, . . . (later, p, q, r, . . . are used as general indices). The third-order contribution to eqn. (1) may be decomposed into three parts: g3{9}
=
SC
+
1 ij abed
9&zbxubcd1cdij
6 11 %job ijkl ab
$k,q
-
225
279 (1993) 223-231
cc 2ijabfkbic90ckj ijk abc
up to first-order terms in the expansions: (7) G~L(PTD) = JWV
(8)
Here the notation {DC} emphasizes the dependence of the solute internal energy components E, and E, on the correlated density matrix to order u (2 2) over the molecular basis via the changes in the fi,, potential: the inclusion of electron correlation effects in the wavefunction produces changes in the electric field created by the solute in its environment which can induce non-negligible changes in I&,. The correlated density matrix to a predetermined order may be obtained (in the molecular orbital basis) by adding to the SCF density matrix, @cr = 6, (Vi,j), the total ath density matrix corrections: gTo= E cgn 0
= gscr + @oRR
(9)
To second order, gcoRR E G@? and the occupiedoccupied blocks are given [17,18] by {gZ)ij = 41
(10)
sjkab9ikob
kab
Because .Ypqrsare integrals over the spin molecular orbital basis, in practical calculations it is preferable to make use of transformed spin-free expressions over the spatial molecular orbitals basis. Then eqn. (10) can be written as G% >I./= - 1 WK4B(2%K4s - %KBA) KAB +
-zJKBA(2%KBA
-
=%KAB)l
(11)
where the integrals xlJAB and sIJAs are defined in a similar fashion to $&$, and L!&&. The elements of the virtual-virtual block may likewise be written as (92)~~
=
1 [~IJACWJBC CIJ +
=%bkl
At the second level of approximation, designated as perturbation theory on the density distribution alone (PTD), we compute the solute-solvent interaction potential J&r, taking into account correlation effects in the wavefunction and including only
+ E, {D,)
-
%JCA(2=%CB
-
~IJCB>
-%JBC)I
(12)
The corresponding elements of the occupiedvirtual block are obtained from {g2jIJ and {&}AB by solving the perturbation-independent linear equation (see for example refs. 17 and 18) @@B.&J
-
EB)ABJ,A,
=
xAI
(13:
226
F.J. Olivares del Valle et ral./J. Mol. Struct.
where ABJ,AIis given by A AI,BJ
=
1 -
c2$AIBJ
-
9ABJI
-
xAJBI)(EI
-
&A)-’
(14) The {XA,}are intermediate quantities made up of molecular integrals and orbital energies: xA,
=
[ll
+
121+ L31+ L41
(1%
where KJ
PI E C(~*}B~~B~CA
-
(17)
XBlac
BC
L31= - c
[y,BKJ(29KJAB
-
dKJBA)
JKB
-
L41= 1
9,BJK(9KJAB
+
[9CBAJ(2%JCB
%UBA)l
-
(18)
=%JBC!)
JBC
-
yCBJA(9,CB
+
=%JBC)l
(19)
Once gz (or D, in the solution case) has been constructed, the matrix gToT (or DToT) is readily transformed to the atomic orbital basis (in our case by obtaining the natural orbitals and the new occupation numbers) and can be used for the calculation of I&T, thus introducing the electron correlation correction into the primary aspects of the PCM. Thus, with respect to the original version (eqn. (4)), the PTD approximation (eqn. (8)) supposes a modification of the results which has its origin in the different electron distribution (represented by its density matrix) used to calculate the virtual electron distribution CJ(S)generated on the cavity surface. Finally, at the third level of approximation (PTDE) the internal energy of the molecule in solution is
are the free energy terms Cp,2Jw4~ relating to the variation of the correlation energy (to order p) including the new potential I$.,, calculated from the correlated density matrix DToT . In
where
(Theochem)
279 (1993) 223-231
sum, the correlated density matrix refers to the situation in solution, taking into account the mutual effects of the solvent reaction potential and the fluctuation potential. In both the PTD and PTDE approximations, the matrix DToT is calculated iteratively: starting from the density matrix in vacua, one calculates the solute-solvent interaction V,NT; this is then used in resolving the Schrodinger equation, when one obtains a new density matrix (in the MP context) that includes the modifications due to correlation and to the effect of the dielectric, and allows the cycle to be re-initiated. Four or five cycles are usually enough to achieve convergence in the energy. At this point it must be made clear that, although I&T includes all the components of the solutesolvent interaction (electrostatic, dispersion, exchange etc.), only the electrostatic component is considered in the present work. This approximation is made with a view to facilitating the interpretation of the results and is based on the results of earlier work which analysed the influence of the dispersion component [19]. This component has a considerable influence on the energy aspect. However, it has only a very small effect on the density matrix of the solute, and hence on its oneelectron properties. Once a model has been adopted to represent the solute and the solvent, one has to specify the interactions included in the interaction potential &,T. In general, for a non-microscopic description of the solute, we use a multipole expansion of the potential. This is the case for the historical model. The expansion must be truncated at some fixed order, accepting a compromise between the quality of the results and the complexity of the method. Reference 16 gives a brief selection of the methods relating to the classical continuum electrostatic approach which use a multipole expansion for version. &NT ; ref. 20 gives a more sophisticated In our case, the potential c,, is not obtained from an analytical parametrix expression, but rather by a numerical procedure which is iterative and convergent (for more details, see ref. 4). By making use of eqns. (1) and (2) or eqn. (20)
F.J. Olivares de1 Valle et a/./J. Mol. Struct.
(Theochem)
279 (1993) 223-231
227
and its vacuum equivalent, one may obtain the variation AE,, of the solute internal energy at any of the three levels of approximation, together with its corresponding components (AE,, is a measure of the electrostatic contribution to the solvation energy). However, this quantity is not employed directly [21] because it is not an energy of primary interest. The magnitude (which is well defined thermodynamically) that interests us is
We compute AG,, using eqn. (24). It may be expressed [ 13,141 as
AG,, = AE,, - N&n +
(21)
which is the free energy of formation of the set of charges off in the cavity of the dielectric. Briefly, the standard free energy of solvation for M, AGL, may be estimated by means of the expression AG; = G;’ - GM = A(PV) + AGF= + Aq
(25) with AG:’ = AGO”i- AG;
(26)
AGO”= +Q4,,
(27)
AG,’= xAei - $(DP&T)
(28)
AGF’= - +cA-“iiu G
(29)
Finally, the contributions AGT’and AG$’ are given by (see eqns. (5) and (6)) (30)
(22)
where the solute-solvent interaction energy [ 1l] AGET is (omitting the contribution from dispersion
AGP’ = AE3”+ AE3”+ AE,4
[191) AGENT = Gw + AG,,
Results and discussion
(23)
The term AG,, is the change in the internal electronic energy of M as it passes from the gaseous phase into solution (we consider that at the temperature T of our conditions, all molecules are in the ground state) AG,, = G;;’ - E,i””
(31)
(24)
Included in Giy’ are the proper solute electronic energy itself, the energy coming from the electrostatic solute-solvent interaction (sum of two contributions, one nuclear, the other electronic), and the polarization energy of the solvent because of the solute. In eqn. (22), the term Aq contains the entropy or thermal contributions associated with the change in the rotational, translational and vibrational energies under the effect of the solvent [13], and AGZ in addition to depending on temperature, is also affected by the molecular geometry in so far as this has an influence on AGENT.In this paper we assume that there is no geometrical distortion on passing from the gaseous to the condensed phase [15,22].
The solutes that we used to analyse the changes produced in the free energy on taking into consideration correlation in the wavefunction are H,CO, HFCO, F&O, F&S, H,CS, HClCO and FClCO. Their geometries are given in Table 1. The basis set used in al the calculations was 6-3 1G and the dielectric constant was E = 78.3 (value equivalent to that of water at 25°C). . The cavities were constructed using a set of interpenetrating spheres of variable radius according to the procedure proposed previously by our group [5-81, which takes into account the population supported by each atom and the basis set that is being employed. The radii used in each case are given in Table 2. The use of this procedure permits the generation of cavities adapted to the local characteristics of each atom in each molecule. This can be seen, for example, in the case of the carbon atom for the series H,CO, HFCO and F&O. This atom undergoes a reduction of the electron population as the hydrogen atoms are progressively substituted by their corresponding fluorine atoms. Thus, in this series, the value of the net charge on the carbon
228
F.J. Olivares de1 Valle et al./J. Mol. Struct.
(Theochem)
279 (1993) 223-231
TABLE 1 Geometrical
parameters”
Parameter
H,COb
HFCO’
F, COd
HClCO”
FClCO’
1.182 1.090
1.173
F,CSp
H,CSh
Bond length 1.210 1.116
R co R HC R FC
1.181 1.095 1.338
_ _
R CC, R CS
Bond angle L HCH L HCF L FCF L HCCl L FCC1 L HCO L HCS L FCO L FCS L ClCO
_ _
116.5 _ _
_ 109.9 _
_ 121.7
1.315 _ _
_ _
_
_
107.6 _
1.087 1.334 1.765 _
_ _
110.5
_ 108.8
127.3
_ _ _
“Bond lengths in kngstrsms;
1.170
122.8 _
_ _
126.5 _
126.2 _ _
123.1
123.7 _ 127.5
1.315 1.725 1.589 _ _ 107.1 _ _ _ _ 126.5
_ 1.611 116.5 _ _ _ _ 121.8 _ _
bond angles in degrees. bRef. 23. ‘Ref. 24. dRef. 25. ‘Ref. 26. ‘Ref. 27. gRef. 28. hRef. 29.
atom changes from + 0.14e in H&O to 1.14e in F$O. This fact is clearly reflected in the value that is used for R,, which goes from 1.9 A in the first case to 1.5 I$ in the second (Table 2). The use of
fixed radii (equal to the van der Waals radii, RVDW, or scaled up by a certain quantity k x RVDW)takes no account of this fact, each carbon atom in any of the solutes under study contributing with an
TABLE 2 Charges and radii”,b Parameter
H&O
HFCO
F,CO
HClCO
FClCO
0.14 -0.58 0.21
0.63 - 0.53 0.32
1.14 - 0.45 _ _
0.57 - 0.42 _ _
- 0.42
-0.34 _
0.120 - 0.46 0.33 -i _
F&S
H,CS
Charges :: Q” z:
Q,,
0.00
-0.35 0.20
0.53 _
- 0.43 _
_ 0.13 - 0.33 _
0.26 - 0.03 _
1.71 _
2.13 _
Radius Rc & & & R, Rc, “In gngstrsms. bObtained with the PTE model.
1.70 1.90 0.90 _
1.50 1.89 _ _
1.88 1.88 0.92 _
1.70 1.87 _ _
1.80 _
1.75 _
2.77
1.74 2.70
2.79 1.74 _
0.95 2.88
3
- 10.89
- 7.33
AGsc,
AG,, 3.40
-7.11
- 10.51
- 0.59
3.99
- 144.10
133.60
H,COC
3.55
- 10.26
- 13.81
- 0.79
4.33
-98.22
84.42
HFCOb
3.28
- 9.77
- 13.05
- 0.73
4.02
-94.60
81.50
HFCO”
1.21
- 4.83
- 6.04
- 0.77
1.98
54.29
- 60.33
F,COb
“In kcalories per mol. bObtained using the PTE model. ‘Obtained
3.56
- 0.82
AG,
AAGd
4.17
AG,
139.56
- 150.45
AG,
AGo
H,COb
Variations in free energiesa
TABLE
0.66
- 0.77
- 1.43
0.36
0.30
27.10
- 28.54
HCICOb
1.07
- 0.56
- 1.63
0.53
0.54
75.84
- 77.47
HClCO’
1.91
- 2.27
-4.18
-0.12
2.03
- 78.65
74.47
FCICOb
using the PTD model. dAAG = AG,,,
1.09
-4.39
- 5.48
-0.77
1.86
47.38
- 52.57
F&O’
- AG,.
1.93
- 2.28
- 4.21
-0.13
2.05
- 78.38
74.18
FClCO’
2.16
- 5.94
-8.10
- 0.65
2.81
- 72.89
64.79
F,CSb
1.99
- 5.82
-7.81
-0.51
2.50
- 72.30
64.49
F,CS
0.72
- 1.97
- 2.69
- 0.30
1.02
27.78
- 30.47
H,CSb
0.55
- 1.67
2.22
- 0.22
0.77
27.43
- 29.65
H,CS
230 TABLE
F.J. Olivares de1 Valle et al./J. Mol. Strut.
(Theochem)
279 (1993) 223-231
4
Variation in molecular properties”,b AP”
H,COd
AQc AQo AQH AQs AQF AQ,,
- 0.009
0.045
0.023
0.120
0.070
-0.134 0.067
- 0.225 0.085
-0.106 0.117
-0.147 0.104
- 0.053
-0.145
- 0.033
- 0.078
- 0.001
- 0.060
H,CO”
HFCOd
HFCO”
F,COd
F,CO” 0.266
MC
0.004
0.009
0.011
0.200
0.025
0.056
A&
0.106
- 0.037
- 0.029
- 0.074
- 0.009
- 0.069
A&
0.030
0.023
0.061
0.055
AA A&
- 0.006
0.200
0.004
0.024
F,CSd 0.018
F&Se 0.154
- 0.019
- 0.014
0.000
- 0.065
0.007
- 0.002 0.002
0.035
Ap
0.442
0.679
0.354
0.426
0.089
0.149
0.038
0.793
2.086
0.085
0.920
- 0.364
0.079
- 0.210
“In atomic units. bAll values represent the electronic contribution variation. “AP = P,, - P,,.
HClCO
0.011
0.078
-0.071 0.093
-0.162 0.142
- 0.043
- 0.058
0.008
0.001
- 0.019
- 0.098
0.047
0.047
- 0.005
0.508
0.287
0.531
- 0.079
1.071
- 0.235 - 0.059
A&,
A
HCICOd
-0.189 0.936
except for the dipole moment which represents the total
dObtained using the PTE model. “Obtained using the PTD model.
identical sphere to the resultant cavity. The above can be generalized to each of the atoms involved. A more details discussion of these questions may be found in refs. 5-8. Table 3 lists the components of the free energy of solvation that we obtained using the PTE and PTD approximations, calculated according to eqn. (26) and (29)-(31). Without going into detail (which we postpone until a future article dealing with a greater number of diverse solutes) we can draw the following general conclusions: (1). In no case are the third-order contributions to the free energy of solvation especially significant, and they hardly vary from one approximation to another. One must not ignore, however, the fact that AG, is a compensated sum of quite diverse contributions (from two, three and four bodies), and some of them individually can undergo major changes. An analysis of these variations can be found in ref. 8. It is therefore enough to consider up to the contribution AG2 in order to extract a value near to the convergence of AGioL. As a general rule, AG3 is a stabilizing contribution tending to raise the absolute value of AGioL. (2). The second-order contribution to the free
energy of solvation, AG2, is quite significant in the solutes we studied, reaching a value near 45% in the case of H,CS. In every case, the contribution tends to reduce the absolute value of AGzoL. It is therefore advisable to consider this contribution in those calculations of solvation which handle free energies. Its importance is possibly critical in the study of reactive processes in which the saturated nature of the initial and final products is different, as the value of AGz (and, in some cases, AG,) is usually influenced by the degree of saturation of the molecule under study. We have also found a certain correlation between the variations of AG2 and the polarizability of the solute [14]. (3). The third conclusion worthy of note is, in our judgement, the modest significance of including electron correlation in the wavefunction. The variations that appear in the value of AGioL are less than 0.5 kcal mall’ , smaller than the margin of accuracy for this type of model. We found, in a more complete sample of solutes, that this behaviour is general and always an order of magnitude below that found using the PTE approximation. Finally, in Table 4 we list the results for the variation of a selection of one-electron properties.
F.J. Olivares de1 Valle et al/J.
Mol. Struct.
(Theochem)
279 (1993) 223-231
In general, one can say that the inclusion of correlation in the density matrix leads to variations which are more significant in this class of properties than in the energies. Thus, one appreciates an increase in the net charge on the atoms and in the size of the electron distribution, and, lastly, an accentuated polarity of the molecule which can reach major variations, as occurs in the solute HClCO, whose dipole moment varies by 0.7 D in the PTE approximation and 1.33 D in the PTD. ,
10 11 12 13 14 15 16
Acknowledgement
17
This research was partially sponsored by the Direction General de Investigation Cientifica y Tknica, Projects 1017/87 and 930/90.
18 19 20
References M. Born, Z. Physik., 1 (1920) 45. L. Onsager, J. Am. Chem. Sot., 58 (1936) 1486. J.G. Kirkwood, J. Chem. Phys., 2 (1934) 351. S. MiertuS, E. Scrocco and J. Tomasi, Chem. Phys., 55 (1981) 117. M.A. Aguilar, M.A. Martin, S. Tolosa and F.J. Olivares de1 Valle, J. Mol. Struct. (Theochem), 166 (1988) 313. M.A. Aguilar and F.J. Olivares de1 Valle, Chem. Phys., 129 (1989) 439. F.J. Olivares de1 Valle, M.A. Aguilar, S. Tolosa, J.C. Contador and J. Tomasi, Chem. Phys., 143 (1990) 371. F.J. Olivares de1 Valle and M.A. Aguilar, J. Comput. Chem., 13 (1992) 115. D. Rinaldi, B.J.C. Cabral and J.L. Rivail, Chem. Phys. Lett., 125 (1986) 495.
21
22 23 24 25 26 27 28 29
231
F. Floris and J. Tomasi, J. Comput. Chem., 10 (1989) 616. J. Langlet, P. Claverie, J. Caillet and A. Pullman, J. Phys. Chem., 92 (1988) 1617. M.A. Aguilar and F.J. Olivares de1 Valle, Chem. Phys., 138 (1989) 327. F.J. Olivares de1 Valle and M.A. Aguilar, J. Mol. Struct. (Theochem), 280 (1993) 25. F.J. Olivares de1 Valle and J. Tomasi. Chem. Phys., 150 (1991) 139. M.A. Aguilar, F.J. Olivares de1 Valle and J. Tomasi, Chem. Phys., 150 (1991) 151. J. Tomasi, R. Bonaccorsi, R. Cammi and F.J. Olivares de1 Valle, J. Mol. Struct. (Theochem), 234 (1991) 401. G.W. Truck, E.A. Salter, C. Sosa and R.J. Bartlett, Chem. Phys. Lett., 147 (1988) 359. G.W. Truck, E.A. Salter, J. Noga and R.J. Bartlett, Chem. Phys. Lett., 150 (1988) 37. F.J. Olivares de1 Valle and M.A. Aguilar, in preparation. M.E. Davis and J.A. McCammon, Chem. Rev., 90 (1990) 509. R. Bonaccorsi, C. Ghio and J. Tomasi, in R. Carbo (Ed.), Current Aspects of Quantum Chemistry, 1981, Elsevier, 1982, pp. 407-426. R. Bonaccorsi, R. Cammi and J. Tomasi, J. Comput. Chem., 12 (1991) 301. K. Tagaki and T. Oka, J. Chem. Sot. Jpn., 18 (1974) 1174. R.F. Miller and R.F. Curl, J. Chem. Phys., 34 (1961) 1847. V.W. Laurie, D.T. Pence and R.H. Jackson, J. Chem. Phys., 37 (1962) 2995. H. Takeo and C. Matsumura, J. Chem. Phys., 64 (1976) 4536. H. Oberhammer, J. Chem. Phys., 73 (1980) 4310. A.J. Careless, H.W. Kroto and B.M. Landsberg, Chem. Phys., 1 (1973) 371. P.H. Turner, L. Halonen and I.M. Mills, J. Mol. Spectrosc., 88 (1981) 402.