2995
erne for the Rapid Computation of Molecular Electronic
:F-LCAO-( ST0)-MO Formalism, and Its Application Superoxide Molecule
e
Billingsley, II,I and Carl Trindle*z
Department of Chemistry, University of Virginia, Charlottesville, Virginia 28901
(Received M a y 86, 1872)
aublicotim costs assisted by the National Science Foundation
A STO- LCAO-SCF computational scheme is described, in which one-center integrals and most two-center inlegrals are evaluated by fast analytic methods, while less tractable integrals are approximated by Gaussian (hxparrsion of the XTO’s. The method is distinguished from the method of Cook and Palmieri by the use of +diegeneralized G-function techniques of Dyatkina and Klimenko so that a greater fraction of all integrals may be evaluated exactly. Illustrative results are reported for formaldehyde, ethylene, and allyl cation, mth comparison to exact STO-SCF computations and the LEDO approximate method. Since our method of Computation is designed for the rapid exploration of potention surfaces for reactive species, we examine ,ortioris of the Li-02 surface. Two stable forms are found, an isosceles triangular species observed by Andrews n matrix-isolation studies and a bent species analogous to HOZ.
Introduction The choice of a hasis set in which to express molecular wave function is a compromise between the accuracy of a Slates-type, basis3 arid the convenience of a GausOne particularly appealing device is sian-type basis t o employ Slatcr-type orbitals (STO’s) to evaluate those integrals which are rapidly done by exact methods,5 and to approximate the more difficult intcgrals. Cook and PnlmierP h v e developed a scheme in which Gaussian expansionb of various accuracy represent STO’s, and multicenter integrals over STO’s are approximated by sums of multicenter integrals over Gaussian funclions; this is called a “mixed basis” computation since the basis functions entering onecenter integrals and Iwo-center coulomb integrals are Gaussian-series approximations to STO’s. The errors attending this type of computation have been carefiiliy di~cussed;~ in brief, the variation principle does not apply t,o result 5, of mixed basis computations. Since Gaussian funci ions decay rapidly in comparison vith Slater functions and multicenter integrals contain major contributions from the long-distance region of the basis functions, such integrals are consistently underestima teci by 1he Gaussian expansion. Charge has a tendenc) l o migrate close to nuclei, into regions of 103~-iiuclear--s-lectronicpotential, since the opposing electron-electron rtlpulsion is not fully recognized; the total energy is lowered by this computational artifice. One of us (F. P.Pi.) has developed a method for approxirnating the drficult multicenter integrals, called the Limited I3xpansion of Diatomic Overlap technique A two-centcr charge distribution is approximated by a series of one-center charge distributions,
with the expansion coeficients chosen so to reproduce the values of exactly computable multicenter integrals. This demand leads to an often ill-conditioned matrix equation which must be solved by special numerical techniques. Errors in this type of calculation may arise from the bad manners of the central matrix equation, or the poor representation of two-center charge distributions by superpositions of one-center charge distributions; the latter error is marked when very tightly bound orbitals occur as in fluorine and oxygen systems. However, the most significant disadvantage of the LEDO method is that the single KO usually associated with a hydrogen atom must be augmented with 2s and 2p functions in order that multicenter integrals involving hydrogen be accurately estimated. I n systems with many hydrogen atoms, the effective basis set rapidly becomes unmanageably large. We wish t o describe a variant of these two approximation methods which avoids some of the disadvantages of each. Our proposal employs Gaussiam representation of STQ’s in order Lo approximate certain (1) Present address, Chemistry Division, U. S. Department of Commerce, National Bureau of Standards, Washington, D. C. (2) To whom correspondence should be directed. A. P. S!oan Fellow 1971-1973. (3) J. C. Slater, Phvs. Rev., 36, 57 (1930). (4) S. F. Boys, Proc. Rou. Soc., Ser. A , 201, 125 (1950). (5) K. Ruedenherg, C. C. J. Roothaan, and W. Jaunzemis, J . Chem. Phys., 24, 201 (1956). (6) D. B. Cook and P. Palmieri, M o l . Phya., 17, 271 (1969). (7) D. B . Coolr, P. D . nacre, J. L. Dodds. and &I. Elder, Theor. Chim. Acta, 22, 167 (1971). (8) F . P. Billingsley, 11, and J. E. Bloor, Chem. Phys. Lett., 4, 48 (1969) ; F. P . Billingsley, II, Dissertation, University of Virginia, 1970.
T h e Journal
of
Physical Chemistry, Val. 76, N o . 8 i p197B
F. P.RILLINGSLEY, 11, AND CARLTRINDLE
2996 multicenter integrals, but by adopting the generalized C-function techniques of Klimenko and Dyatkina9 we are able t o evaluate exactly a larger proportion of integrals than are done in previously reported mixed basis computations. At the same time we avoid the artificial inflation of the basis set size necessary in the IAEDQtreatment of systems with many hydrogen atoms. ‘The accuracy of our proposed variant is tested by comparison with an exact STO minimal basis computation and LED0 results for HzCO with a closely related STO cornputation on GH4, and with a mixed basis study of allyl cation.] Since the intended use of this method i s the study of molecular potential functions, we present pieliminary results of geometry computations on the IIJiOz system.’O
Details of the
Table I: Comparison of Exact, LEDQ8and Mixed Basic Computations in a h!Iinimal Basis of Rlater Orbitals on Formaldehyde
_----__ 2
Comparison of LE 0, Mixed Basis, and Exact Minimal Basis S Computations on Formaldehyde Tablc 1 contains total energies, occupied-140 energies, and atomic charges for I-IzCO. Aung, Pitzer, and ChanI6 report exacl computations in a minimal basis of Slater-type orbitah with exponents prescribed by Slator’s rules. Our mixed basis results arc slightly but consistently superior t o the LEDO results, perhaps because t,he tightly drawn in orbitals of oxygen provide an insufficient representation of niulticcnter charge distributions in th-, LEDO computation. The tenThe Journal of Physical C’memistrg, Vol. 76, N o . 21, 1972
0.0 0.0 0.91655 -0.91655
0.0 0.0 0.0 0.0
Mired basis
Exact18
- 113 .449
- 20.589 - 11.357 -1.369 -0.837 -0.675 -0.671 -0,470 -0.385 0.247 0.615 0.746
Basis Computation
The familiar formulation of the molecular orbital equations in the LCAO approximation due to Roothaanll is tlre basis for all computations. Different Fock matrices for different spins are used for opensheil systems, foilowing Pople and Nesbet.Iz The exact methods of Ruedenberg, Roothaan, and Jaunzemisj allow svalua ,ion of all one-center one- and twoelectron integrals, and with the generalized C functions described by Rlimernko arid D y a t l ~ i n a ,all ~ two-center integrals except Gxchange integrals. Three-center nuClem attraction integrals and all remaining t%-o-electron integrals are computed by replacing STO’s by the Gaussian series of S t ~ r n a r t , ’and ~ evaluating the resulting integrals over Gaussian orbitals by standard methods, l 4 These expansions are not obtained variationa1ly1&but match the STQ’s in a least-squares may. They should represent the long-distance region of the A 0 better than energy-optimized expansions do, and provide better estimates of the multicenter integrals. he naclear attraction integrals is so fast Ey ernploy the most accurate six-Gausof thc STO, while the three-Gaussian expansion is used for the more numerous and more time-consuming multicenter two-electron integrals. Per cent errors IIZ the approximations to very small integrals secin large when the Lhree-Gaussian expansion is used, but wra error in the integral evaluation exceeds
Geometry--------Y
0.821
113 4482 1.0031 - 20 588 ~
I
-11.351 -1.369 -0.829 -0.674 -0.574 -0.468 -0.385 0.247 0,602 0.754 0.818
z
0.0 -1.2709 8 0.52917 0.52917 LED08
- 113.471 - 20.558 -11.359 -1.356 -0.857 -0.673 -0.566 -0.477 -0.398 0.244 0.600
0.761 0 . a57
Charges
-0.162 -0.080 0.121
-0.179 -0.079 0.129
-0.148 -0.082 0.116
dency of computations approximating repulsion integrals to provide anomalously low energies is indicated by the -0.022 au error in the LEDO energy.
Comparison of Mixed Basis and Exact Computations on Ethylene Pallte and Lipscombl’ report an ab initio study of ethylene in which the minimal basis of Slater orbitals with exponents chosen by Slater’s rules for carbon and an STO with an exponent of 1.2 for hydrogen. We do not present a LEDO computation for this system since the number of hydrogens is sufficiently large that LEDO computations are much slower than mixed basis computations. Exact and mixed basis results for ethylene are set forth in Table 11. Some errors associated with (9) N. M. Klimenko and M . E. Dyatkina, Zh. S t ~ u k fKhim., . 6 , 604, 755 (1965); also N. M. Klimenko, E. N . Krylova, P;.X. Mikhaleva, Y. I. Churikov, and SI.E. Dyatkina, ibid., 6 , 407 (1965). (10) L. Andrews, J . Cliem. Phys., 50, 4288 (1969). (11) C. C. J. Roothaan, Rea. Mod. Phys., 23, 69 (1951). (12) J. A. Fople and R . K. Nesbet, J . Ciiem. Phys., 22, 571 (1954). (13) R. F. Stewart, ibid., 5 2 , 431 (1970). (14) I. Shavitt, Methods Comput. Phys., 2, 45 (1963). (15) K. 0-ohata, H. Taketa, and 9. Huainaga, J . Phys. Soc. Jap., 21, 2306 (1966) (16) S.Aung, R. M. Pitser, and 5 . Chan, J . Ciiem. Ph,ys., 45, 3457 (1966). (17) W. Pallre and W. Lipscomb, J . Amel. Chrm. Soc., 88, 2384 I
(1966).
COMPVTATTON G F MOLECULAR ELECTRONIC ENERGY mixed basis computations become apparent. The kinetic energy i13 anomalously high, but the electronic potential compcnsates in the opposite direction, so that the total electronic energy is in error by only -0.0145 au. A 8 a result of t,he contraction about the nuclei, the virial ratio -2T/V i s higher for the mixed basis compul ation than for the exact computation.
z
2/
z
13 rk 0.92664
0 0
10.675 Iir1.210 A(PL -
STO/G
STO/G)
-11.2875 - 1k ,2868 - 9.0144 -0.7u23 - 0.6438 --0.56 16 -0.5061 - 0.3709 0.2426 0.5888 0,6206 0,6395 0.8453 0.8917
-11.2807 -11.2795 -1.0166 0.7711 - 0.6400 -0.5713 - 0.5146 - 0.3688 0,2443 0.5899 0,6092 0.6283 0.8396 0.9126
- 0.0068 - 0.0073
77.52’75
77.6633
-. 111.0830
- 7’7 8362 0.99302 -0.280 I
- 111.2312 - 77.8486 0.99881 - 0.2946
Table 111: Allyl Cationa
Etotal
Ethylene Coordinates, d
-
orbitals. We used single STO’s with Slater exponents, while Cook, et al., used SCF-atomic ork)itads. We chose only the three-Gaussian approximation t o AO’s entering multicenter integrals.
Gzometry CCC = 118‘ Rc-c = 1.39 A Rc-H = 1.09 A HCH = 120’
Table 11 : A S g s t m with Numerous Hydrogen Atoms: Ethylene Molecule
PL,U
2997
$0.0022 - 0.0112 - 0.0038 $0.0097 $0.0085 - 0.0021 -0.0017 - 0.0031 - 0,0114 0.0112 0.0057 - 0.0209 -0.1358
- 0.1482 0.0124 - 0.00079 $0. 015
In spite of the errors in the kinetic-potential partitioning in the mixed basis computation, the small errors in total energy (-0.014 au), orbital energies (0.01 au), and atomic charges (0.01 electrons) suggest that the mixed basis method may be very useful in the rapid and appr oximate study of potential surfaces.
Comparison of the Proposed Mixed Basis Method with a Previous I\/iixed Basis Computation We have suggested that our variant of the mixed basis method nnay be more reliable than the method discussed in rciil 7 since a greater fraction of the large integrals can be done exactly by the techniques of Klimenlco and D ~ a t k i n a . ~Table I11 contains results of our mixed basis computation on allyl cation, where the errors of tne mixed basis scheme become serious. A significant difference between our computation and that of ref 7 is thtvt oar programs are not equipped t o incorporate sums of STO’s to represent the atomic
la1 1bz 2a1 3aL 2bz 4al 5al 3bz 4bz 681
1biw 1atr 2bi,
STQ/G
3 G7
4G7
-116.7146 - 11.704 -11.704 - 11.618 -1.421 -1.186 -1.045 -0.991 -0.917 -0.834 -0.774 -0.707 -0.252 -0.002
-115.9221 -11.770 -11.770 - 11.709 -1.432 -1.263 -1.095 -0.983 -0.958 -0.926 .- 0.856 -0.762 -0.326 -0.089
-115.8031 - 11.780 -11.780 - I1 “711 - 1.426 -1.257 -1.098 -0.995 -0.950 -0.884 -0.852 -0.750 -0.324 -0.098
HB
Charges 6.15 6.04 0.79 0.80 0.69
- 2T/V
Virial Ratio 0.9957
Interior C Terminal C
H, Hz
6.31 6.26 0.66 0.65 0.61
6.29 6.24 0.67 0.65 0.63
a Column heading 3G indicates that all multicenter integrals are approximated by representing STO’s by a series of three Gaussians; 4G means a series of four Gaussiaiis was used. STO/G refers to the present work, in which only a portion of multicenter integrals are approximated, by use of a three Gaussian series for each STO.
While the 1\10energies in the mixed basis results of Cook, et al., are consistently depressed relative to our M O energies, presumeably due to their superior carbon AO’s, the ordering scheme agrees, As we expect, passage from a three-Gaussian to a €our-Gaussian representation generally raises the total energy in Cook’s computation, in this case by 0.1 stu. The supposition that Cook’s four-Gaussian energy i s within 0.1 au of the exact value is of interest since oiir three-Gaussian cnergy lies 0.1 au above Cook’s energy. However, it seems more reasonable to suppose that our threeGaussian result is closely comparable l o the previous scheme’s four-Gaussian result.
Mixed Basis Study of the Li-02 System Andrews has reported the infrared spectrum of LiOz isolated in an inert matrix.1° The spectrum is consistent with a Cz, arrangement, while the small value of the 0-0 stretching force constant suggests an ionic The Journal of Physical Chemistru, Vol. 7 0 , iVo. 21, 1572
F. P. BILLINGSLEY, 11: AND CAIKTRINDLE
2998
-1.0
0
1.0
2.0
30
50
40
Figure 2. Contours of the t o t d density for the C, L i 0 ~system show that the Li 2s A 0 provides an extension of one lobe of the 0 2 ?r* MO, reducing the 0-0 a antibonding and establishing an 0-Li bond.
-30
-20
-1.0
0
IO
2.0
30
Figure 1. (a) A contour plot of the total density for the Cz, Li02 system shows that the system is basically a Li’ core which strongly polarizes the charge distribution of an 0 2 species, in agreement with the bonding picture proposed by Andrews. Inner contour(s) represents a density of 0.5 electron/ (bohrp; the contour8 form a geometric series, ( 1 / 2 ) n + L electron/(bohr)a. Black dots represent nuclei; scale is in bohrs. (b) Contours of the density associated with highest occupied MO in Cz,tion show that polarization is accomplished by transfer from an 0 2 - a* 3x0 to the Li p orbital of like symmetry.
interaction in which a Li electron is transferred to an O2 antibonding a orbital. To illuEAxate the intended application of our fast mixed basis method of computation, we have examined a portion of the LiO, potential surface. The known Cz, species was given close attention. The computed equilibrium geometry, force constants, charges, virial ratio, and dissociation energy are collected in Table IV, with the experimentally obtained counterparts where available. The geometry is well reproduced, though the force constants are substantially in error. The atomic (;\iuJikenls) charges are consistent with the qualitative view that bonding is accompanied by charge transfer from Li to 02. Since the 14ulliken populatioris can be misleading in a system with orbitals of such different exponents, we provide a contour plot of the charge density in the molecular plane (Figure 1). The low-force constant for motion of the Li atom parallel with the 0-0 axis (the “b motion”) suggests The Journal ofPkyeicn1 Chemistry, Vol. 76, N o . $1, 197s
the possibility of an interconversion to B system in which Li is bound to only a single oxygen atom. We adopted the computed 0-0 and Li-0 distances from the Czuspecies and traced the change in energy as the 0-0-Li angle was altered. The least energy a.ngle was 135”; the bonding in the open species is similar t o bonding in HQ219(Figure 2 ) .
Table IV: Geometry, Force Constants, and Computed Quantities for Li0f Theoretical results RO-Li
Ro-o F(Li-02 a-stretch), mdyn/A F ( 0 - 0 stretch) F(Li-02 b stretch), Edissaoiation
1.72 8 (0.05j 1.30 A (0.04) 5.2
19.2 1.3 0.103 ail 64 kcal/mol - 1.56.547 BLI 1.00570 0.197 electrons - 0 098 electrons
Experiments1 estimates10
I . 77 I 33 ~
A (0.07) A (0.06)
I.lS(0.01) 5,59(0.06) ca. 50 kcal/mol
I
a CzVgeometry species; note very small computed force constant opposing b symmetric motion of Li. This is the motion which leads to either dissociation or a stable bent G , LiOt.
The energy of the open (C,) LiQz system is comparable to the energy of the known Cposystem; since we did not explore variations in 0-0 and Li-0 bond lengths, it is very likely that wc liavve not found the precise minimum in energy of the open system. Thus (18) R. S.Mulliken, J . Chem. P h p . , 23, 1833, 1841, 2338 (1955). (19) D. H. Liskow and H. F. Schaeffer, 111, J . Amw. Chem. Soe., 93, 6734 (1971).
CQMPDTATIO N IJF MOLECULAR ELECTRONIC ENERGY
I
I
3 1.3 ,
1,5
1.4
"_
I
I
1.6
1.7
I
1.8
I
1.9
Figure 3. The dissociation of L i 0 ~such that CzVsymmetry is maintained requires the population of an a A state ~ lying above the 2Ezbound state, The motion of I,i parallel to the 0-0 axis is a permitted dissociation route, as well as a route to the bound Lion system of C, symmetry.
we expect that the open species will prove to be more stable than the isosceles triangular molecule. The iact that the CvZvrnolecule was discovered first may be due to the means of preparation, the greater intensity of infrared transitions in a system with largely ionic bonding relathe to a system with predominately covalent bonding, or the kinetic stability of the CZvsystem. This kinetic stability arises from the fact that dissociation or isomerization from the C2, molecule requires a change in electronic state from the ground 2Bzstate to dissociative 2A2state (see Figure 3). The easiest way t o isomerize or dissociate the Cz,, Li02 is to move the Li atom so that the symmetry is broken; the analogy betw-etm this motion and the motion of carbenes adding to double bonds is notable. The experimental 6tudy of photoisomelrizntion of CzvLiOzto the C, species is contemplated. 2o
~
~
~
i
~
g
Multicenter integral evaluation is the most timeconsuming part of m y ab initio computation. It is relevant t o compare computation times in the LEDO approximation and also with a pure Gaussian basis to ab i n i t i o times. Timing comparisons are not often easy clue t o differences in machine cycle times, differing usage of peripheral storage (which affects times dras-
2999 tically), and differences in details of programming which can change run times considerably. On our Burroughs BFi500, which has a cycle time comparable with an IBM 360/SO, the rigorously determined one- and twocenter integrals require about 0.040 seelintegral. The LED0 estimate of multicenter iniegrals* requires 0.009 sec/integral regardless of the size of the system. Billingsley* estimates this t o be roughly t w o orders of magnitude faster than the rigorous method of Wahl,21 and roughly comparable with the ox^^^^^ P method of Newton.22 We ha8vestressed Billingsley's timing results for his LED0 approximation since all our experience shows that the STO/G method reported here competes in speed for small systems and is superior in speed for systems with many hydrogen atoms, vrrhere many auxiliary p functions must be artificially introduced in the LEDQ algorithm. The question remains whether a pure-Gaussian cal. culation competes with the STO/G ~ ~ t h o dShillady has developcd an ab initio Camsinn lobe program for the Burroughs B5500. If a sufficiently small Gaussian basis i s chosen so that times of the S'TO/G and Gaussian cornputations match, the energy evaluated by the Ctiussian lobe method is far inferior tJothe STO/C energy. Typically 95% of the EIarLree-Fock limit i~ produced by Gaussian lobe calculation in the time required by ~~~/~ to produce 99% oF the Eari,rw-Fock energy.
Summary An approximate, near-ab initio computational system is described which is suitable for rapid exploration of molecular potential functions. Results of electronic structure computations on the LiOz system are presented for the lcZv species, which is experimentally known, and for a C, form predicted t o be comparable in energy with the known form.
Acknowledgments. This work was supported by the lJnited States National Science Foundation through Grant No. GP-3081'7 t o C. T. Postdoctoral support for F. P. B. was provided by the University of Virginia. ~ o ~ ~ u tcosts a t were ~ o defrayed ~ by a grant from the Frederick Gardner Cottrell Fund of the Research Corporation, and Grant No. 1830-Ca of the Petroleum Research Fund of the American Chemical Society. (20) L. bndrews, private oornninnication. (21) A. C . Wahl and R. H. Land, Int. J. Quanfum Chem., 1, 375 (1967); J . Claem. Phgs., 50, 4725 (1969). (22) M. D. Newton, W. A. Lathsn, W. J. Helire, and J. A. Pople, ihid., 5 2 , 4064 (1070).
The Joiirnal of Phvsical Chemistry, Vol. '76, N o . $1, 107.9
