J . Phys. Chem. 1984,88, 4880-4886
4880
97 kcal/mol relative to hydroxide binding. As discussed above, the weakening of the metal-oxygen bond will drop the proton affinity of the water, allowing the proton on the water to drift toward the Glu-Thr subunit. This will have a synergistic effect on the bicarbonate binding. In the limit of complete proton transfer, the binding energy of the bicarbonate drops to 40.3 kcal/mol. In effect, the proton shift regenerates a hydroxyl on the metal to compensate for the loss of bicarbonate; the local electrostatic potential around the zinc remains unchanged throughout. The internal transfer of the bicarbonate proton localizes negative charge on the oxygen nearest the zinc. This increases the binding energy of bicarbonate to 220.8 kcal/mol in the L3MOH2complex and to 65.5 kcal/mol in the limit of proton transfer to the Glu-Thr subunit. Dissociation of the bicarbonate would be accompanied by solvation and we can suggest tbe relative hydrogen bond energies involved: A single water to cafbon dioxide hydrogen bond is 3.05 kcal/mol while a hydrogen bdnd to bicarbonate is 18 kcal/mol. The large number of hydrogen bonds possible thus corresponds to a large energy of solvation. Conclusions These calculations illustrate a number of salient points about the mechanism of carbon dioxide hydration. The proton affinity of metal bound water is consistent with a zinc hydroxide mechanism. The proton on this water can exchange with solvent, possibly assisted by proton hopping through the His 64-water chain. The proton relay (Glu 106-Thr 199) proposed by Kannan8 is not needed to deprotonate the water during the approach of substrate. The metal can accomodate a fifth ligand, water, and the dissociation energies are consistent with observed rates of ligand exchange. The fifth water is weakly bound to the zinc hydroxide complex; it is more strongly bound in the diaquo complex. Metal to ligand distances correlate inversely with the binding energies reflecting the strong Coulombic nature of the interacfipn. The
relaxation of bond distances upon protonation of the hydroxide (or after nucleophilic attack) is an important part of the reaction mechanism. The addition of substrate (carbon dioxide) results in the formation of an adduct shown in Figure 3. The metal-oxygen bonds relax, shrinking the metal-water bond and lengthening the metal hydroxide as negative charge moves out onto the bicarbonate. Direct formation of carbonic acid via a proton transfer from the fifth coordinate water ligand is not observed. However, a water-mediated proton transfer of the bicarbonate as proposed by Pocker (Figure 4) is favorable. There is a strongly cooperative coupling between the loss of bicarbonate product and a proton transfer from the fifth coordinate water to the Glu 106-Thr 199 subunit. As the bicarbonate leaves, the proton affinity of the remaining water drops, allowing the proton to drift to the threonine. This shift regenerates the hydroxide, further weaking the metal-bicarbonate bond. In effect, one unit of negative charge has moved from the hydroxide to the bicarbonate while a second has moved from the glutamate to the water. The average electrostatic potential around the zinc remains essentially constant. (Note that the shrinkage of the metal-water bond prevents the bicarbonate from forming a bidentate ligand with a correspondingly high binding energy.) The catalytic cycle is completed with the addition of water and the release of a proton. This step involves multiple proton transfers and a significant barrie2 is expected. His 64 can aid this step by providing multiple pathways for proton release.
Acknowledgment. The authors thank the NIH, G M 26462, for generous financial assistance. Registry No. OH-, 14280-30-9; H 2 0 , 7732-18-5; HCOO-, 71-47-6; Im, 288-32-4; HOCOO-, 71-52-3; [L3M-OH]+, 91389-06-9; [L,M(OH),], 91389-07-0; [L3M(OH)(OHJ]+, 91409-22-2; [L,M(OH)(qCOOH)], 9 1409-23-3; [L3M(OH2)(OCOOH)],9 1389-08-1; HC03-, 71-52-3; carbonic anhydrase, 9001-03-0.
Model Potential Method in Molecular Calculations Sigeru Huzinaga,* Mariusz Klobukowski, Department of Chemistry, University of Alberta, Edmonton, Alberta, Canada T6G 2G2
and Yoshiko Sakai College of General Education, Kyushu University, Ropponmatsu, Fukuoka, Japan 810 (Received: February 13, 1984)
The theoretical background of the model potential method has been reviewed with emphasis on practical aspects of preparation of atomic model potentials. The newly constructed basis sets and model parameters were used in molecular calculations (F2, IF, Iz, KrF,, XeF,, XeF,, CaO, ScO, CrO, NiO, NiCO, and Ni(CO),). Consistently satisfactory results of these calculations render strong support to the view that the model potential method is a reliable and economical tool for molecular structure calculations.
Introduction In an attempt to widen the scope of reliable molecular calculations, a number of approximate methods have been proposed. The methods, known as pseudopotential, effective potential, or model potential, are aimed at simplifying the treatment of chemically inert core electrons. (A recent review of the methods covered thp period up to 1977.' Contributions from various laboratories report later d e v e l ~ p m e n t . ~ ~Such ~ ' ) a method has (1) R. N. Qixon and I. L. Robertson, "Theoretical Chemistry", Vol. 3, R. N. Dixon and I. L. Robertson, Eds., The Chemical Society, London, 1978, p 100.
0022-3654/84/2088-4880$01.50/0
been developed during the past decade in our l a b ~ r a t o r y . ~ - ~ Recently, our model potential (MP) approach has been modified so as to allow for approximate treatment of relativistic effect^.^ (2) P. A. Christiansen, K. Balasubramanian, and K. S. Pitzer, J . Chem. Phys. 76, 5087 (1982). (3) I. Ortega-Blake, J. C. Barthelat, E. Costes-Puech, E. Oliveros, and J. P. Daudey, J . Chem. Phys., 76, 4130 (1982). (4) L. S. Bartell, M. J. Rothman, and A. Gavezzotti, J. Chem. Phys., 76, 4136 (1982). ( 5 ) V. Bonifacic and S. Huzinaga, J. Chem. Phys., 60, 2779 (1974). (6) Y.Sakai, J . Chem. Phys., 75, 1303 (1981). (7) Y. Sakai and S. Huzinaga, J. Chem. Phys., 76, 2537 (1982). (8) Y. Sakai and S . Huzinaga, J . Chem. Phys., 76, 2552 (1982).
0 1984 American Chemical Society
Model Potential Method in Molecular Calculations The aim of the present paper is to review the M P method as used in atoms and molecules. During the past year we have prepared the M P parameters together with small economical basis sets for nine atoms and we have tested their performance in molecular calculations. The obtained results are reported here. The process of preparation of the model potentials is similar to the process of preparation of small, economical Gaussian basis sets for molecular calculations. In both cases the process is performed via atomic calculations. However, the quality of a given model potential-or of a small basis set-cannot be judged by its performance in the atomic calculations alone. The simplicity of the prepared computational tool requires that its quality and usefulness have to be judged in test molecular calculations. The test calculations may force us to take some corrective actions: splitting the valence part of a minimal basis set, introducing the scaling factors for the exponential parameters, or, in the case of the model potentials, modifying values of some parameters. (Preparation of small basis sets and of model potentials is similar also in that not too many research groups are willing to perform this tedious task.) The present review of the model potential approach and its applicntions is organized as follows. First, the atomic model potential formalism is presented, followed by a discussion of some important issues (choice of the core functions, choice of the valence functions, and introduction of some relativistic effects). Next, the molecular formulation of the method is reviewed, and certain aspects of the molecular calculations are discussed (choice of diffuse and/or polarization functions and core-valence interactions in the molecular environment). Finally, the results of the molecular calculations are presented for a variety of molecular systems. Very recently, our method has been applied in studies of correlation energy in diatomic halides and halogen hydrides'O and in studies of interaction energies in weakly bonded aggregates of the noble-gas atoms." These results have been presented at this symposium and will be published elsewhere. Atomic Model Potential Formalism A detailed discussion of the theoretical background of the M P method has been recently given by Sakai and Huzinaga.7 Here we present the results of the discussion. (The atomic units will be used throughout, unless specified otherwise). The atomic M P Hamiltonian for the N , valence electrons of an atom a may be written as a sum of the one-electron and the two-electron parts: Np
Np
The two-electron term describes the Coulomb repulsion between the valence electrons. The one-electron part contains, apart from the kinetic energy and the nuclear attraction energy terms, the characteristic model potential terms which approximate the valence-core interaction and which prevent the valence orbitals from collapsing onto the core space occupied by the NF core electrons:
Cf"
where Z" represents the atomic number of atom a , AIa,a p , and B; are the M P parameters, and vla is allowed to take values 0 and 1.7,12The coordinates Fai of electron i in the coordinate frame (9) M. Klobukowski, J . Comput. Chem., 4, 350 (1983). A number of papers on the relativistic effects may be found in the recently published volume "Relativistic Effects in Atoms, Molecules, and Solids", G . L. Malli, Ed., Plenum, New York, 1983. (10) S.Huzinaga, M. Klobukowski, and Y. Sakai, paper presented at the 8th Canadian Symposium on Theoretical Chemistry, Halifax, Canada, Aug 1983. (1 1) J. Andzelm, S. Huzinaga, M. Klobukowski, and
E. Radzio-Andzelm, Mol. Phys., in press. (12) G. Hojer and J. Chung, I n t . J . Quantum Chem. 14, 623 (1978).
The Journal of Physical Chemistry, Vol. 88, No. 21, 1984 4881 centered on the atom a are abbreviated as i. In the analytical (matrix) formulation of the M P method both the core and the valence orbitals are expanded in terms of the basis (Gaussian) functions:
(3) (4)
Thus, the M P terms of the Hamiltonian ha(i) (eq 2) give rise to merely two new types of integrals:
It is worth mentioning here that in the determinantal formulation of the M P method the last two terms of the one-electron Hamiltonian (eq 2) have to perform collectively triple duty: they have to represent the local Coulomb interactions, to approximate the nonlocal exchange terms, and to act as a genuine projector. In the M P approach an atom is described by a set of core functions, {$/), a set of valence orbitals, (+pa), and a set of the M P parameters, Ala, are, and B / . We will discuss the choice of the orbitals first. The core atomic orbitals could have been taken as any available Gaussian representation. However, it is only recently that such functions have been systematically prepared for all atoms up to Rn.13 Furthermore, for atoms beyond Rn such functions are not available. Thus, having in mind also the quasi-relativistic extension of the method, we decided to use the Gaussian expansions obtained by the least-squares fitting to the numerical radial functions. The latter may be obtained at the present time without any difficulties for all the atoms by using both the nonrelativi~tic'~ and the quasi-relativi~tic'~ formulations. The details of our fitting technique have been described e l ~ e w h e r e . ~ The valence atomic orbitals (Le. 5;" and C,,, of eq 4) are also determined by the least-squares fits to the numerical radial functions. This approach allows for great flexibility in representing the M P valence orbitals. First of all, we are able to retain the nodal structure of the orbitals (it should be noted that often a single Gaussian is sufficient for proper location of inner nodes). Next, we can increase the accuracy of the valence part to meet requirements of a particular system under study. For example, we can use a sufficiently large number of Gaussians to accurately describe the outer parts of the valence orbitals, which is necessary in studies of weak interactions in the van der Waals systems." Furthermore, the same procedure may be used for both nonrelativistic and quasi-relativistic formulations of the M P method. Once the analytical representations of both core and valence orbitals have been chosen, the M P Fock matrix becomes a functional of the M P parameters Ala, ala, and B/: F(A*,a*,B")C, = SC,e,
(6)
where C, denotes the set of the valence expansion coefficients for (eq 4), S is the matrix of the valence basis function orbital overlaps, and e, denotes an eigenvalue. The values of the Bea parameters are expressed in terms of the respective reference orbital energies, Fca, as
+,
B," = -Fcqca (7) At the initial stage of the M P parameters optimization the values of F,"'s are usually taken to be equal to L7 The remaining parameters, A and a, are determined from minimization of the differences between the reference valence orbital shapes and energies and the ones generated from the M P eigenvalue problem (eq 6). The minimized expression has the form (13) S . Huzinaga, J. Andzelm, M. Klobukowski, E. Radzio-Andzelm, Y. Sakai, and H. Tatewaki, "Gaussian Basis Sets for Molecular Calculations", Elsevier, Amsterdam, 1984. (14) Ch. Froese-Fischer, "The Hartree-Fock Method for Atoms", Wiley, New York, 1977. (15) R. D. Cowan and D. C. Griffin, J . Opt. Soc. Am., 66, 1010 (1976).
4882 The Journal of Physical Chemistry, Vol. 88, No. 21, 1984
TABLE I: SCF Results for the Spectroscopic Constants of the F2 Molecule'
P
where w,,and W, denote appropriately chosen weight factors, 6, is the orbital energy difference
6 , = It, - Z,,I and A,, is the orbital function difference
(9)
k
where P denotes the radial part of the valence orbital and the summation runs over points generated by the numerical atomic HF program.16 It should be noted here that using the relative orbital function difference A i
A i = CI(Pp&(r)- Pfik(r.)j/P#&(r.)12 k
( 1 Ob)
allows for securing excellent accuracy of the valence orbitals in their outermost region, thus providing good quality functions for studies of the weak interactions." It should be noted here that the exponents of the Gaussian functions, 3;", are not reoptimized in the atomic M P calculations. To sum up, the parameters of the M P method for an atom CY are determined as follows: (a) the core orbitals 6"and Cpu from LS fitting; (b) the valence orbitals 3;" from LS fitting and C:, from the eigenvalue problem (eq 6); (c) the M P parameters AIa and utu from minimizing the deviation u (eq 8); (d) the M P parameters B," (or, alternatively, F;) from the LS fitted valence basis functions. (The dependence of u on F," is usually small. Therefore, they are readjusted in molecular S C F calculations on simple diatomic molecules.) The relativistic extension of the M P method9 is based on the reference radial functions and orbital energies computed via the relativistic Hartree-Fock (RHF) method of Cowan and Griffin.15 In the R H F method the mass-velocity and the Darwin terms of the Pauli Hamiltonian are added to the Hartree-Fock differential equations, which are subsequently solved in a self-consistent fashion. This technique allows for incorporation of the major direct and indirect relativistic effectsl7-l9 into the orbital functions while allowing for retention of the commonly used nonrelativistic formulation of the energy level calculation. The spin-orbit effects, which are omitted in this approach, could be treated as a perturbation, possibly simultaneously with the electron correlation effects. The numerical radial functions calculated by the R H F method are used together with the corresponding orbital energies in determination of the quasi-relativistic M P (QRMP) parameters in exactly the same way as for the nonrelativistic M P ( N R M P or MP) case.
Molecular Model Potential Formalism The molecular M P Hamiltonian consists of the sum of the atomic Hamiltonians (eq 1) and the internuclear repulsion term
where N = XaNa is the total number of the valence electrons in a molecule (the sums Earun over all the atomic centers in the molecule), hu(i) is given by eq 2, and the internuclear repulsion term, P@(R,,),has been approximated aszo V'fl(R,B) = (2"- N,")(Z@- N,B)/R,p
Huzinaga et al.
(12)
where R,, is the distance between the nuclei CY and p. It is worth mentioning here that implementation of the M P formalism may be very easily accommodated within any existing all-electron (AE) molecular calculation program: the post-HF parts of the program and the two-electron packages need not be (16) Ch. Froese-Fischer, Compur. Phys. Commun., 4, 107 (1972). (17) J.-P. Desclaux, Int. J. Quantum Chem., 6, 25 (1972). (18)P. Pyykko, Adu. Quantum Chem., 11, 353 (1978). (19) K. S.Pitzer, Acc. Chem. Res., 12, 271 (1979). (20)V. Bonifacic and S. Huzinaga, J. Chem. Phys., 62, 1509 (1975).
R,/A
w,/cm-' DJeV
tw
AE
HF
EVW
TMB
exptl
1.330 1273 -1.41
1.337 1259 -1.12
1.323 1257 -1.37
1.32 1638
1.418
1.42 916.6 1.68
'Abbreviations: tw, this work using F-1 model potential for fluorine; AE, all-electron S C F calculations from ref 24; HF, Hartree-Fock results from ref 25; EVW, effective core potential of Ewig and Van Wazer;26TMB, pseudopotential of Teichteil, Malrieu, and Barthelat;27 exptl, experimental values from ref 28. TABLE II: SCF Results for the Spectroscopic Constants of HF"sb
ReIA w,/cm-' DeleV
tw
AE
exptl
0.903 4503 4.19
0.898 4494 4.26
0.915 4138 6.12
'Basis set for hydrogen was ( 4 ~ )contracted ~ ~ to (31).23 *Abbreviations: tw, this work using F-1 model potential for fluorine; A€?, all-electron SCF calculations from ref 24; exptl, experimental values from ref 28. changed; the new one-electron integrals are the three-center generalizations of the integrals given by eq 5a,b; the derivatives with respect to the nuclear coordinates (for the analytic gradient method of the geometry optimization) are calculated easily.z1 The QRMP molecular calculations are done using the quasirelativistic model potential parameters and valence basis functions determined in the procedure described in the previous section. In order to ensure good quality of the molecular (MP or AE) calculations, we need polarization functions and diffuse functions. The exponents of the polarization functions, awl,used in the present study were obtained by maximizing the value of the overlap integral Such polarization functions usually provide sufficient quality of the S C F results for strong chemical interactions. A better choice is to use two polarization functions determined from the overlap criterion. However, if we want to perform calculations beyond the HF level, then some optimization of the polarization function(s) in the molecular environment may be needed. In the present study we employed only a single polarization function. The diffuse functions are usually needed in studies on molecular anions, on Rydberg states, on correlation energy, and on weak interactions. In the M P studies on CaO reported here we needed such functions, obtained via the fitting technique where the excited-state numerical radial functions were used as the reference functions. The least-squares fitting based on A+' (eq lo), with a proper choice of the radial points k , has been successfully used to generate some needed diffusion function exp0nents.l' In our M P MCSCF calculatiodOthe additional diffuse functions were optimized in the atomic M P calculations on the respective anions.
Results of Molecular Calculations The model potential method has been recently applied to studies on molecules belonging to three very different groups: the halogen and interhalogen molecules, compounds containing noble-gas atoms and fluorine, and the transition-metal compounds. For the systems containing heavy atoms the relativistic effects were investigated. The present section shows the results of these calculations. The model potential parameters and basis sets are not presented here, but they are available from the authors upon request. A standard basis set for h y d r ~ g e n ~has * ? been ~ ~ used in calculations on HF molecule. Only the S C F results are shown, and they are compared with other S C F results. (21) K.Kitaura, private communication. (22)S.Huzinaga, J . Chem. Phys., 42, 1293 (1965). (23)T.H. Dunning and P. J. Hay in "Modern Theoretical Chemistry", Vol. 3, H. F. Schaeffer 111, Ed., Plenum Press, New York, 1977,p 1 .
Model Potential Method in Molecular Calculations
The Journal of Physical Chemistry, Vol. 88, No. 21, 1984 4883
TABLE 111: SCF Results for the Spectroscopic Constants of IFa AE SG tw EVW exDtl 1.99 1.910 1.936 1.86 1.918 R,/A 756 610 698 637 w,/cm-' 711 Abbreviations: AE, results of all-electron S C F calculations from ref 24; SG, results of all-electron S C F calculations using small, balanced Gaussian basis sets (ref 31); tw, this work using F-l model potential for fluorine and 1-2 model potential for iodine; EVW, effective core potential of Ewig and Van Wazer;26 exptl, experimental values from ref 28.
F,, HF, IF, I,, I,-, and I,+. The results for the diatomic fluorine and for hydrogen fluoride are shown in order to assess the quality of the M P for the atom which will be used in further studies. For the light fluorine atom the M P approximation works very well. The results for F, and HF are given in Tables I and 11, respectively, and compared with results of AE24v25and other effective core potential calculation^^^^^^ and with the experimental data.28 The present results for both Fz and H F are quite close to the ones obtained in the AE calculations. This finding has been known since the first papers on the M P method. However, when the M P method is applied to such small systems, the gain in computer time is minimal on the S C F level and is even smaller on the post-HF level calculations. The main benefit of using the model potentials for the first-row atoms lies in avoiding the basis set superposition e r r 0 r s ~ ~(BSSE) 3 ~ ~ which would be present in the AE calculations if only a small number of Gaussians was used for description of the 1s orbital. The BSSE originate from improving description of atomic shells (mainly 1s shell) by basis functions centered on other atomic sites of a molecule. Thus, the composite system is described better than its components, which leads to errors in bonding energies and charge di~tribution.~' Replacement of the 1s shell, which is the major contributor to BSSE, by a model potential eliminates a substantial part of this error. A side benefit is that, due to removal of the 1s shell, more Gaussian functions may be used for description of the valence shells. The results of SCF calculations on I F (Table 111) show that the two model potentials used in the calculation provide a balanced description of the valence shell bonding effects in the interhalogen compound. The calculated bond length of 1.94 A is only 1%larger than the value obtained in the all-electron calculations using a comparable valence basis set (SG). The model potential used for iodine, 1-2, has been carefully chosen in studies described below. The I, molecule was a difficult test for the model potential method. The first M P calculation on the molecule gave a bond length of 2.806 A,32i.e. nearly 0.1 A larger than the H F value of 2.720 A.33 This discrepancy has been attributed to difficulties in modeling the large interactions between the valence shell and , ~ ~ difficulties the outermost core shell on the M P l e ~ e l . ~These may be solved in three ways. Firstly, the outermost core electrons may be included e ~ p l i c i t l y .Secondly, ~ the outermost core may be treated in the frozen orbital (FO) a p p r o ~ i m a t i o n .Thirdly, ~~ the value of the M P parameter F (eq 7) may be adjusted so as to correct for the M P deficiencies. The three approaches give the (24) P. A. Straub and A. D. McLean, Theor. Chim. Acta, 32,227 (1974). (25) G. Das and A. C. Wahl, J . Chem. Phys., 44, 87 (1966). (26) C. S. Ewig and J. R. Van Wazer, J . Chem. Phys., 63,4035 (1973). (27) Ch. Teichteil, J. P. Malrieu, and J. C. Barthelat, Mol. Phys., 33, 181 (1977). (28) K. P. Huber and G. Herzberg, "Molecular Spectra and Molecular Structure", Vol. IV, Van Nostrand Reinhold, New York, 1979. (29) S. F. Boys and F. Bernardi, Mol. Phys., 19, 553 (1970). (30) J. J. Kaufman, private communication. (31) J. Andzelm, M. Klobukowski, and E. Radzio-Andzelm, J . Comput. Chem., in press. (32) 0.Gropen, S. Huzinaga, and A. D. McLean, J. Chem. Phys., 73,402 (1980). (33) A. D. McLean, 0. Gropen and S. Huzinaga, J . Chem. Phys., 73, 396 (1980). (34) A. Stramberg, 0. Gropen, and U. Wahlgren, to be submitted for publication. (35) A Stremberg and 0. Gropen, to be submitted for publication.
TABLE IV: Properties of I1 in the Ground Electronic State 'Eg+' HF SG 1-4 FO 1-2 1-30 exDtl R,/A 2.678 2.704 2.732 2.727 2.757 2.736 2.667 w,/cm-' 236.2 236.2 227.1 222.0 226.4 224.8 214.5 1OZb0,/ 6.24 5.25 esu-cm2 "Abbreviations: HF, near-HF results (ref 33) obtained with f-type polarization function (the bond length increases to 2.720 A when the f function is removed from the basis set"); SG, S C F results obtained with a new Gaussian basis set (ref 31); 1-4, M P results from present study using 1-4 model potential; FO, frozen-orbital results from ref 35; 1-2, M P results from present study using 1-2 model potential; I-3Q, Q R M P results from present study; exptl, experimental values from ref 28.
TABLE V: Electron Affinities and Ionization Potentials of L o s b HF 1-2 I-3Q TMB exptl EA, I2('C,) 1.04 1.28 1.20 1.7 A 0.07 12-(22")
-
The vertical values in electronvolts. bAbbreviations: H F , near-HF results from ref 33; 1-2, M P results from present study using 1-2 model potential; I-3Q, Q R M P results from present study; TMB, pseudopotential results from ref 27; exptl, experimental values from ref 33. 'Value taken from ref 27.
bond lengths of 2.732, 2.727, and 2.757 A, respectively (cf. Table IV). It seems that the three methods give results of comparable quality; however, it only the third one which does not increase the computational effort in the molecular calculations. For the iodine molecule we have also prepared a quasi-relativistic model potential. The results of the Q R M P calculations (column I-3Q, Table IV) indicate that the relativistic effects are responsible for a slight bond shortening of 0.02 A. The remaining difference between the QRMP value of Re and the experimental value, which is equal 0.07 A, may be attributed to the neglect of the correlation effects and omission of the f-type polarization functions. Tables V and VI illustrate results of M P calculations on ions of 12. The calculations were done for IF in the state and for the 1,' ion in both ,rI6 and 226states. The M P parameters and basis sets were the same as employed in the calculations on the neutral iodine molecule (Table IV). The electron affinities and the vertical ionization potentials reported in Table V were all calculated by using the ASCF method. The M P values are quite close to the HF ones. The differences between the computed values and the experimental data may be due to both basis set effects and correlation energy corrections. The equilibrium properties of the ions are shown in Table VI. The overall agreement between the results of the AEHF calcul a t i o n ~and ~ ~the results of the present M P calculations is very good considering that the M P parameters were the same as in the calculations on the ground state of I,. KrF,, XeF,, and XeF4. The electronic structure of krypton difluoride has been studied on the all-electron level by Bagus et who used a large STO basis set, and by Collins et al.,37who employed a Gaussian basis of STO-3G quality. In the model potential calculations we first treated the 3d electrons explicitly as the valence electrons (Kr-1 model potential), and then we buried them among the core electrons, reducing at the same time the value of the parameter F (eq 7) (Kr-2 model potential). The Kr-1 model potential calculation compares very well with the all-electron calculation: the Kr-F bond length is virtually the (36) P. S. Bagus, B. Liu, and H. F. Schaeffer 111, J. Am. Chem. Soc., 94, 6635 (1972). (37) G. A. D. Collins, D. W. J. Cruickshank, and A. Breeze, J . Chem. Soc., Faraday Trans. 2, 70, 393 (1974).
4884
The Journal of Physical Chemistry, Vol. 88, No. 21, I984
Huzinaga et al.
TABLE VI: Properties of Ions of ItUvb
I2-(*zJ method HF
Re 3.313 3.392 3.350 3.281
1-2 I-3Q
exptl
I2+(*II,)
we
6,
106.8 107.7 110.1 109.0
-14.65 -13.90
Re 2.593 2.646 2.625
12+(zz,)
we
8,
227.0 265.3 261.3
15.75 14.60
Re 3.053 3.196 3.180
we
Be
120.7 115.5 109.7
25.89 27.00
in A, w, in cm-', and Be in esu.cm2. bAbbreviations: HF, near-HF values from ref 33; 1-2, MP results from present work using 1-2 model potential; L3Q, QRMP results from present work using I-3Q model potential; exptl, experimental values from ref 33.
- 0.4
TABLE VII: SCF Results for KrK ( D , P EXT Kr- 1 Kr-2
AE
- 0.5
(A) R K ~ - F I A 1.811 1.822
1.896
- 0.6
-0.468' -0.531 -0.599 -0.676 -0.731 -1.216
- 0.7
1.813
5a2u=--===z
(B) Valence Orbital Energies/auc -0.541 -0.567 -0.648 -0.711 -0.766 -1.241 -1.560 -1.579 -3.984 -4.003 -4.012
-0.549 -0.569 -0.656 -0.714 -0.769 -1.246 -1.573 -1.594 -3.998 -4.015 -4.024
-0.541 -0.570 -0.652 -0.713 -0.767 -1.243 -1.568 -1.588
Abbreviations: EXT, AE calculations using extended STO basis set (ref 36); Kr-1, MP results from present work using Kr-1 model potential; Kr-2, MP results from present work using Kr-2 model potential; AE, AE calculations using STO-3G basis set (ref 37). bResults obtained at 1.9 A. cCalculated at R K ~=F1.852 A (3.5 bohrs). (I
same in both cases, and the largest difference in the orbital energies is 0.015 au (Table VII). The more economical Kr-2 model potential gives almost the same results as the ones obtained with Kr- 1. The geometry predicted by Bagus et al.36and by our model potential calculations differs by more than 0.05 A from the experimental data (1~ 7 and5 1A89 ~ A39), ~ The Kr-F bond length predicted by Collins et aL3' is closer to the experimental one; however, this seems to be the result of a poorer basis set used in those calculations. The bonding in the xenon difluoride molecule has been studied theoretically for the first time by Basch et who used a small Gaussian basis set in performing SCF calculations at a Xe-F distance of 2.00 A. A more detailed study was done by Bagus et al.,41 who determined geometry and binding in calculations which evaluated correlation energy using an extended STO basis set. Recently, a pseudopotential study has also been published.42 We have performed four sets of S C F calculations on XeF2. Firstly, the all-electron calculations have been done in order to establish a good reference point for the MP calculations (Bagus et al.41did not report results of the one-configuration SCF calculations). Secondly, the nonrelatilvistic MP calculations were done in which the 4d electrons were treated explicitly (Xe-1). Next, similar quasi-relativistic calculations were performed (Xe-3Q). Finally, the nonrelativistic MP calculations were repeated, this time with the 4d electrons included into the core and the value of F (eq.7) reduced to 1.0 (Xe-2). The results, collected in Table VIII, again show the ability of the MP method to simulate quite accurately the results of the AE calculations (columns 2-5). The shorter bond length predicted by Bagus et aL4' may be due to both a large basis set (f polarization functions included) and (38) W. Harshbarger, R. K. Bohn, and S. H. Bauer, J . Am. Chem. SOC., 89, 6466 (1967). (39) C . Murchinson, S. Reichman, D. Anderson, J. Overend, and F. Schreiner, J . Am. Chem. SOC.,90, 5690 (1968). (40) H. Basch, J. W. Moskowitz, C. Hollister, and D. Hankin, J. Chem. Phys., 55, 1922 (1971). (41) P. S.Baaus. B. Liu, D. H. Liskow, and H. F. Schaeffer 111. J . Am. Chem.' Soc., 97,7216 (1975). (42) L. S. Bartell, M. J. Rothman, C. S. Ewig, and J.-R. Van Wazer, J . Chem. Phys., 73, 367 (1980).
- 0.8 - 0.5 - 1.0 - 1.i - 1.2 -1.3
- 1,4 - 1.5 -1.6
- 1.7 E&U.
3.4
3.6
3.8
4.0
A
Figure 1. Valence orbital energies of XeF,.
to performing the calculation on the two-configuration level. The relativistic bond contraction is quite small in XeF, (0.003 A). The agreement between the quasi-relativistic bond length (1.978 A) and the experimental one (1.977 f 0.002 A in gas phase43) is fortuituous, as in the calculations the correlation effects were not included and the f polarization functions were missing in the basis set. The bond length predicted by the EP method of Ewig4, is too long compared with both the AE results and the experimental value. Almost all the methods predicted the same electronic configuration of XeF,: it is only our AE calculations (column 2, Table VIII) which puts the 5?r, level above the loa, one (this happen to be the same ordering as the one obtained in the extended calculations by Bagus et a1.:' however, those calculations are the two-configuration SCF). For the more complex XeF4 molecule no Bagus-quality calculations were performed. The only AE SCF calculations reported by Basch et aL40 employed a relatively small Gaussian basis set in a single-point calculation at RXWF= 1.95 A. Bartell et al.,42 using the effective potential method of Ewig, obtained the Xe-F bond length equal to 2.037 in D4,, symmetry. In our MP calculations we have used Xe-2 and F-1 model potentials, arriving at Xe-F bond length of 1.957 A. Again, the agreement with the experimental value of 1.95 is probably due to simultaneous neglection of the correlation, relativistic, and basis set extension .&42343
(43) S. Reichman and F. Schreiner, J . Chem. Phys., 51, 2355 (1969). (44) M. Yoshimine, J . Phys. SOC.Jpn., 25, 1100 (1968). (45) K. D. Carlson, K. Kaiser, C . Moser, and A. C. Wahl, J. Chem. Phys., 52, 4678 (1970).
The Journal of Physical Chemistry, Vol. 88, No. 21, 1984 4885
Model Potential Method in Molecular Calculations TABLE VIII: SCF Results for XeF, (Dah)' EXT AE1 1.989
1.965
Xe- 1 (A) 1.98 1
Xe-3Q
Xe-2
EP
1.993
2.068
-0.537 -0.499 -0.632 -0.657 -0.71 1 -1.037 -1.547 -1.552
-0.350 -0.466 -0.5 17 -0.601 -0.632 -1.010 -1.413 -1.417
AE2
RXeF/8,
1.978
(B) Valence Orbital Energies/aub loa, 5H" 357s 4H" 6 UlJ 9% 5 uu 8%
-0.576 -0.500 -0.640 -0.663 -1.028 -1.020 -1.544 -1.527
-0.509 -0.514 -0.617 -0.655 -0.701 -1.042 -1.535 -1.539
-0.537 -0.520 -0.637 -0.665 -0.716 -1.046 -1.549 -1.554
-0.534 -0.529 -0.633 -0.665 -0.715 -1.112 -1.544 -1.551
-0.470 -0.499 -0.587 -0.636 -0.675 -1.008 -1.482 -1.487
"Abbreviations: EXT, results of two-configuration S C F using extended STO basis (ref 41); A E l , results of AE calculations from present work (The basis set used was (22s, 17p, 10d/7s, 4p) contracted to [6s, 5p, 3d/3s, 2pl. The minimum total energy was -7430.4652 au.); Xe-1, M P results from present study using Xe-1 model potential; Xe-3Q, M P results from present study using Xe-3Q quasi-relativistic counterparts of Xe-1 model potential; Xe-2, M P results from present study using Xe-2 model potential; EP, results of the effective potential calculations from ref 42 (The E P orbital energies were calculated a t RXep = 2.00 k, (3.78 bohrs).); AE2, results of AE calculations from ref 40. (Only a single S C F calculation was done at RXeF = 2.00 8, (3.78 bohrs).) bCalculated at RXeF = 2.01 k, (3.8 bohrs). TABLE IX: SCF Results for the ]EtState of CaO" HF-1 HF-2 Ca-1 ReIA 1.840 1.855 1.842 w,/cm-' 903.7 803.0 847.4
MP 0.735
"Abbreviations: HF-1, near-HF results from ref 44; HF-2, near-HF results from ref 45; Ca-1, M P results from present work using Ca-1 model potential for C a and the oxygen model potential from ref 7; MP, M P results obtained by 0. Gropen et al.46 TABLE X SCF Results for the *E+State of ScO" STO GTO sc-1 R,IA w,/cm-'
1.62 1373
1.69 974
1.63 1047
sc-2 1.52 750
"Abbreviations: STO, AE results from S T O basis set (ref 48); GTO, AE results from G T O basis set (RVV) (ref 47); Sc-1, M P results from present study using Sc-1 model potential for scandium and 0 - 1 model potential for oxygen (3p electrons treated in the valence electron set); Sc-2, M P results from present study using 0 - 1 model potential for oxygen and Sc-2 model potential for scandium. (In the latter the 3p electrons were treated as the core electrons.).
effects. Variation of the orbital energies with the Xe-F distance D4,, symmetry is shown in Figure 1. Metal Oxides: CaO, ScO, CrO, and NiO. In the case of the metal oxides careful preparation of the model potential parameters and basis sets turned out to be crucial in securing good quality results. In calculations on the ' Z + ( u 2 )state of CaO the model potential Ca-1, which leads to excellent agreement with the AE results (Table IX), required an extended basis set. The diffuse two p functions and three d functions were obtained by fitting to the numerical 4p and 3d radial functions calculated from 4 ~ 4 p ( ~ and P) 4~3d(~D configurations ) of calcium. The oxygen model potential was taken from ref 7. The careful choice of the basis set functions substantially improved the M P results obtained by Gropen et al.46 The transition metal oxides pose new difficulties for the model potential calculations: even though the 3p shell electrons are inert chemically and are well described by the frozen core approximation:' the model potential method cannot describe adequately the 3 p 3 d nonlocal exchange interactions. The early M P calculations of Gropen et al.46indicate a very short bond distance for the scandium oxide molecule. In our own M P calculations on the *E+state of this molecule (Table X), the "collapse" does not occur and a reasonable value of the bond length (1.52 A) is obtained. The M P calculations in which the recalcitrant 3p electrons are treated explicitly give again excellent agreement with the AE values. It is worth noting here that the AE calculations predict (46) 0. Gropen, U.Wahlgren, and L. Pettersson, Chem. Phys., 66,459 (1982). (47) L. Pettersson and U. Wahlgren, Chem. Phys., 69, 185 (1982).
TABLE XI: SCF Results for the IEt State of CrO" AE Cr-1 ReIA 1.77 1.74 w,/cm-' 598
P2 1.60
Abbreviations: AE, AE results from ref 46; Cr- 1, M P results from present study using Cr-1 model potential for chromium; P2, P2-type model potential results from ref 46. TABLE XII: SCF Results for the ]E+ State of NiO" DrODertV AE MP 1u/au -1.120 -1.122 rd/au -0.658 -0.674 -0.594 6/au -0.574 2u/au -0.396 -0.396 -0.347 HP/au -0.343 1.64 1.62 418, we/cm-' 892 717 a Abbreviations: AE, all-electron results from present work; MP, M P results from present work using Ni-1 and 0 - 4 model potentials for nickel and oxygen, respectively.
TABLE XIII: SCF Results for the 'E State of NiCO",b Property 1u/au 2a/au 3u/au *P/au 6/au rd/au 4a/au
AE -1.495 -0.763 -0.621 -0.605 -0.395 -0.367 -0.316
-1.526 -9.780 -0.628 -0.614 -0.392 -0.374 -0.314
MP
R,(Ni-CO) 18,
1.713
1.689
a R c q = 1.15 k, (fixed). bAbbreviations: AE, the AE results from present stqdy; MP, the M P results using Ni-1, 0 - 1 , and C-1 model potentials for nickel, oxygen, and carbon, respectively.
different values of the bond length, which bracket the experimental value of 1.67 In the case of chromium oxide (Table XI) our careful choice of the valence basis functions, accompanied by the reduction of the value of F from 2 to 1.2, results in satisfactory +greement of the M P bond length with the AE value. The 3p electrons were included into the core. Our result is a significant improvement over the value reported by Gropen et al.46 The model potential parameter values and basis set functions for the nickel atom were determined in atomic M P calculations on the 4 ~ ' 3 d ~ ( ~ofDthe ) atom. Our molecular AE reference calculations were done using the (43321/431*1*/311) basis set for Ni and the (421/31) basis set for oxygen. The reduced, soft-core value of F (= 1.0) brings a b o u t very good agreement between the results of AE and MP calculations (Table XII), both
Huzinaga et al.
The Journal of Physical Chemistry, Vol. 88, No. 21, 1984
4886
I
I
I
,
I
I
. --, 1.74
&JW.
1.7'4
184
1.94
2.04
A
i.84
1.94
Figure 3. Variation of atomic Mulliken charges in Ni(CO),,.
Figure 2. Valence orbital energies of Ni(CO)A.
for the bond length and orbital energies. Nickel Carbonyls: NiCO and Ni(CO)4. The studies on NiCO in the lX+ state were undertaken in order to provide a good quality model potential for nickel which could be subsequently used in calculations on Ni(C0)4. The AE calculations have been performed first, giving us the molecular reference results. In the calculations the C-0 distance was kept fixed at 1.15 A, the experimental C-0 bond length in Ni(CO)+ The minimum of the AE total energy was found at RNiqO = 1.713 A (Table XIII). In the model potential calculations 0 - 1 and C-1 model potentials have been used for oxygen and carbon atoms, respectively. In the first M P calculation, which employed the Ni-3 model potential, the Ni-CO distance turned out to be 0.16 A too long (RNIqO = 1.875 A). Reduction of the F value to 1.2 (Ni-2) resulted in shortening the Ni-GO bond length to 1.813 A, still larger than the AE value. The final model potential for the nickel atom, Ni-1, brings about good agreement with the AE bond length and orbital energies (Table XIII). The geometry of NiCO has been recently investigated@by using Goddard's effective core potential method.49 The computed Ni-CO bond length equals 1.75 A, in reasonable agreement with our AE and MP results. The LCGTO-Xa method gives 1.65 and 1.15 A for the equilibrium internuclear distances Ni-CO and C-O, respecti~ely.~~ In the past decade a number of calculations were performed on the nickel tetracarbonyl m o l e c ~ l e . ~In~ -all~ ~those calculations, which aimed at determining the electronic structure of the molecule, the geometrical structure was fixed at the one found o f 0.03 A = 1.82 Rc4 = 1.15 experimentally,60i.e. R N ~ =~ 1.82 f 0.02 A in Td symmetry. (48) A. B. Rives and R. F. Fenske, J . Chem. Phys., 75, 1293 (1981). (49) C. F. Melius, B. D. Olafson, and W. A. Goddard 111, Chem. Phys. Lett., 28, 457 (1974). (50) B. I. Dunlap, H. L. Yu, and P. R. Antoniewicz, Phys. Reu. A, 25, 7 ( 1982). (51) I. H. Hillier and V. R. Saunders, Mol. Phys., 22, 1025 (1971). (52) K. H. Johnson and U. Wahlgren, Int. J . Quantum Chem., 6, 243 (1972). (53) J. Demuynck and A. Veillard, Theor. Chim. Acta, 28, 241 (1973). (54) I. H. Hillier and M. F. Guest, Mol. Phys., 27, 215 (1974). (55) E. J. Baerends and P. Ros, Mol. Phys., 30, 1735 (1975). (56) R. Osman, C. S. Ewig, and J. R. Van Wazer, Chem. Phys. Lett., 39, 27 (1976); 54, 392 (1978). (57) B.-I.Kim, A. Adachi, and S. Jmoto, J . Electron Spectrosc. Relaf. Phenom., 11, 349 (1977). (58) J. Demuynck, Chem. Phys. Lett., 45, 74 (1977). (59) Y.Sakai and S.Huzinaga, J . Chem. Phys., 76, 2552 (1982). (60) L. E. Sutton, Ed., Spec. Pub1.-Chem. Soc., No. 11, 1958.
Our previous single-point model potential studySgresulted in the electronic structure of Ni(C0)4 identical with the one obtained in AE calculation^.^^*^^ However, the nickel model potential used in that studyS9failed to give a reasonable description of the geometry of tbe molecule. In the present investigation of the equilibrium geometry of Ni(CO), we have used the best nickel model potential (Ni-1) prepared in the studies on NiCO. The carbon-oxygen distance was fixed at 1.15 A,6oand the respective model potentials used for these atoms were C-1 and 0 - 1 . The equilibrium Ni-CO distance was found to be 1.978 A. Although this value is 7% larger than the experimental one, it is not unreasonable in view of the limitations of the basis set used (particularly no f functions on nickel) and neglect of the electron correlation effects. Another nickel model potential, Ni-4, in which only four d Gaussians were used gives R e N , q O = 1.965 A. Figures 2 and 3 illustrate variation of the orbital energies and the Mulliken charges on atoms with Ni-CO distance. In the region investigated there is no crossing of the orbital energy curves; thus, assignment of the photoelectron spectra may be done at any Ni-CO distance.
Conclusions It has been shown that the model potential method is capable of producing results which are very close to the results of the all-electron calculations (whenever such are available) for a variety of molecules. The method also allows one to estimate the relativistic effects in systems where the spin-orbit coupling is small. For molecules containing a number of small atoms the use of the model potential to replace the inner core electrons may not result in a substantial saving in computational time, but it provides a good way of eliminating the basis set superposition error originating from the core regions. For molecules containing heavy elements, the saving is substantial and the method offers an excellent possibility to explore a vast region of chemistry untouched so far computationally. In particular, studies on geometrical configurations of such molecules, as well as investigation of reaction paths for heavy-atom systems, become feasible. Acknowledgment. Part of this work has been supported by a grant from the Natural Sciences and Engineering Research Council of Canada. We are grateful to the Computing Services and M.K.) and at the Kyushu at the University of Alberta (S.H. University and Institute for Molecular Science, Okazaki, Japan (Y.S.), for grants for computer time which made this work possible. Registry No. NiO, 1313-99-1; NiCO, 33637-76-2.
