J. Phys. Chem. 1983, 87,3865-3872
Because the z axis is directed toward the vacuum from the slab, a decrease in Dpdcorresponds to a decrease in the work function, and vice versa. The dipole moment D,, from mixing s and p orbitals is proportional to the work function change and has been calculated similarly. This calculation for the energy levels, occupation values,
3865
and-dipole moment values is performed at a series of points in k space. For the (100) lattice, the reciprocal lattice consists of a square, and we sampled over l/s of the area. Similarly, for the (111) lattice, of the area of the Brillouin zone was examined. We_averagedall values with appropriate weighting a t these k space points.
A Formalism for Self-consistent Treatment of “Perfect-Pairing’’ Electron Correlation Effects in Crystalline Systems T. H. Upton Corporate Research Science Laboratories, Exxon Research and Engineering Company, Linden, New Jersey 07036 (Received December 21, 1982; I n Final Form: M r c h 22, 1983)
Exact equations for the total energy and eigenvalues are obtained for semiinfinite three-dimensional crystalline wave functions which describe an important class of electron correlations: those between antiparallel spins in closed-shell singlet (bonding) electron pairs. A local (Wannier) orbital representation is used to describe each bond pair and electron correlation effects are introduced via double excitations localized at each cell. The resulting self-consistent wave functions are formally analogous to molecular generalized valence bond (GVB) wave functions. It is shown that, when expressed in a Bloch form, the effect of the double excitations is to couple states of different k within pairs of filled bands. It is also shown that the energy expressions may be decomposed into a translationally invariant portion, and a set of “self“-Coulomband exchange terms confined to an arbitrary origin cell. This decomposition allows a rapidly convergent summation procedure to be used in evaluating the energies which properly accounts for screening of long-range interactions. While only GVB wave functions are explicitly treated, the development may be applied to any MC SCF wave function in which excitations may be confined to individual unit cells.
I. Introduction In recent years, the array of techniques that have been directed toward the treatment of electronic structure problems in quantum chemistry has become enormous. In almost all cases, these methods have been developed in such a way as to maximize their effectiveness for problems with molecular boundary conditions, that is problems where $(r) 0 as r 0. As a result, their application to systems with long-range order and semiinfinite extent is generally less than optimum, and, in the case of true first-principles techniques, not possible without major simplification (such as truncation of the material) or reformulation of the equations. Unfortunately, this prevents the detailed examination of a rather large set of problems of considerable current interest in chemistry, such as those in the areas of surface chemistry and catalysis (e.g., the effects of surface reconstruction or impurities on thermodynamic quantities in chemisorption), molecular crystals (e.g., electronic excited states and long-range interactions), and the chemical physics of anisotropic materials (e.g., polyacetylenes, layered compounds, and thin films). To be sure, there are methods where extensions to problems with periodic boundary conditions is straightforward, and indeed those whose origins lie in solid-state physics. In all of these, however, this flexibility is a result of simplifications that are made in the exact nonrelativistic Hamiltonian to reduce the treatment of nonlocal exchange and electron correlation effects to local problems. This simplification has proven valuable in many applications, but for problems involving bond formation and the interaction of systems with net spin, a more exact treatment is required. In recent years, considerable effort has been
- -
directed toward the development of techniques capable of providing exact solutions for the Hartree-Fock (HF) equations that result from attempting to evaluate the total energy of a single determinant wave function of an infinite periodic system using an “exact” Hamiltonian.’+ The methods used have been varied and have included approaches that involve the Fourier transform of the l / r l z potentia1,l of the exchange interaction,2 and of the wave f ~ n c t i o n .Real ~ space techniques have utilized both the integration of the full lattice Coulomb and exchange fields over a single unit cell,4 and the inverse in which the potentials experienced by a single cell are integrated over the full l a t t i ~ e .Each ~ of these techniques makes extensive use of the translational invariance of the closed-shell (or high-spin) Hartree-Fock Hamiltonian in simplifying the problem. While the solution of the crystalline Hartree-Fock problem represents a significant achievement and is a necessary step in the development of accurate first-principles techniques for such systems, the wave functions and energies that result are of limited value for problems requiring detailed energetic information. If we leave aside concerns about the well-documented fundamental flaws (1) See, for example, F. E. Harris, Theor. Chem.: Adu. Perspect., 1, 147 (1975). (2) N. E. Brener and J. L. Fry, Phys. Reu. B, 17, 506 (1978). (3) A. Mauger and M. Lannoo, Phys. Rev. B, 15, 2324 (1977). (4) (a) J. L. Calais and G. Sperber, Int. J. Quantum Chem., 7 , 501 (1973); (b) G. Sperber and J. L. Calais, ibid.,7,521 (1973); (c) G. Sperber, ibid., 7, 537 (1973). (5) T. H. Upton and W. A. Goddard 111, Phys. Reu. B, 22, 1534 (1980). (6) (a) S. Suhai, J. Chem. Phys., 73, 3843 (1980); (b) C. Pisani and R. Dovesi, Int. J. Quantum Chem., 17, 501 (1980).
0022-365418312087-3865$01.50/0 0 1983 American Chemical Society
3866
The Journal of Physical Chemistry, Vol. 87,
No. 20, 1983
in the Hartree-Fock density of states for extended systems,’ the single-configuration wave function is still only marginally valuable as a quantitative tool for their detailed study. In an earlier study, Upton and Goddard5 contended that treatment of the HF problem using real space basis-set expansion techniques provided a natural foundation from which extensions to more complex wave functions and Hamiltonians might be easily considered. In this paper, a first step in this development is presented. It will be shown that, if one defines a basis of orthonormal functions localized at each site or unit cell (formally Wannier functions), painvise correlation effects may be introduced into the problem through the formation of generalized valence bond (GVB) “perfect-pairing” wave functions in which each GVB pairs is localized a t a particular site in the lattice. The “perfect-pairing’! wave function is sufficient to describe in a self-consistent manner the degree to which the two electrons in each individual singlet pair will seek to minimize their repulsive interactions. The exclusion of such electron correlation effects is known to be a major source of error in simple H F treatments. While formulating the problem in this way leads in principle to a total wave function that is a linear combination of an infinite number of Slater determinants, translational symmetry may be used to reduce the variational problem to one that formally resembles that of a single localized pair.
11. Symmetry Considerations We define first a semiinfinite periodic array of identical unit cells each of which contains one or more atoms. For simplicity, we require that each cell possess an inversion center. We identify an origin within the lattice such that each cell j may be located by a vector R, from the origin. The exact position or number of atoms within the cell at this point is unimportant, it is only necessary to define a basis of m normalized Gaussian or Slater orbitals in each cell with exponents a (denoted pa)that is of sufficient size to describe the orbital amplitudes within the cells. It is convenient to project this atomic orbital basis into lattice symmetry or Bloch functions of the form
where N is the number of cells (and is large), and the transformation variable k is a wave vector ranging from 0 to a value depending on the particular lattice geometry (Le., all values within the first Brillouin zone (BZ)). The transformation is precisely analogous to the process of changing from a local to symmetry orbital basis through the application of Wigner projection operators in molecular systems. When defined in this way, the Bloch functions are related to the original localized atomic orbitals by the inverse transform BZ
@ =
Cexp(-ik.Rj)cpi: k
(2)
and thus we have m Bloch functions associated with each value of the wave vector k. Clearly, both the local and Bloch basis orbitals depend on r; however, for simplicity of notation we will not explicitly indicate this dependence. A GVB wave function for this system may be written as (7) See, for example, J. Callaway, “Quantum Theory of the Solid State”, Academic Press, New York, 1974. (8) F. W. Bobrowitz and W. A. Goddard I11 in ‘Modern Theoretical Chemistry: Methods of Electronic Structure Theory”, H. F. Schaefer 111, Ed., Plenum, New York, 1977, Vol. 3, pp 79-127.
Upton
~a($~(k~)rC/,(k,)...$~(km~)rC/~(k~)~ ..$2(kmax)... [~111~,12(R1) -X121~2l2(Rl)I [~211~112(R2) - X221~212(R2)l... [ X N ~ ~ ’ J ~ I ~-(X RN ~ z) i & i ~ ( RX~ ) l [X112~122(R1) - X1228222(Rl)l... [X.vipoip2(R,v)- X ~ z p ~ z p ~ ( R(3) ~)l~l
~ L A T=
where N is the number of cells (and is large) and A is the singlet spin function appropriate for the wave function. The orbitals +(k) are closed-shell (CS) orbitals. As eigen vectors of a closed-shell (Hartree-Fock) Hamiltonian (which is totally symmetric with respect to lattice spacegroup operations) these orbitals possess the symmetry of the lattice and may be expanded in Bloch basis as m
$,(k) [ = $,(k,r)I = fln-1’2(k)CC;(k)&
(4)
n
with coefficients that depend on k, on the basis function a, and the root n of the m X m Hamiltonian matrix obtained for each value of k. The normalization of these one-electron Bloch functions is m
N
l
The orbitals Olp(Rj)are natural orbitals and are most easily expanded in the terms of the local A 0 basis as m N
olp(Rj) [= olp(Rj,r)l = Jn,p-”2(Rj)CCC$(Ra+ Rj)&+, a a
(6)
with coefficients which depend on the A 0 function CY and its site a relative to site j , as well as the natural orbital ( I ) and GVB pair ( p ) indices. These natural orbitals are related to the usual nonorthogonal GVB orbitals x,(R,).* (xjlp201p2(Rj)- k,vpe,,p2(Rj))(aP- 01. = (xp(Rj)xp’(Rj)+ Xp’(R,)Xp(Rj))(aP- Pa) ( 7 ) The orbitals Blp(Rj) for given values of 1 and p are identical by symmetry and orthonormal for all sites j . As a result, the coefficients are independent of j and, in subsequent discussion the ceh designation, j will be dropped from these terms. Similarly, the expansion coefficients C;(Rj + R,) will be the same at each site (to within a sign) and the normalization m h’
.mlp(Rj) = CCCP,(R,)C%(R,)(~~+jI~+j) =Alp
(8)
a,Oa,b
is also independent of j . Again, as indicated parenthetically in eq 4 and 6, the orbitals $,(k) and B,,(R,) depend explicitly on r; however, except where necessary, we will identify each orbital exclusively by its location vector Rj (and pair indices 1 and p ) , or its wave vector k (and band index n). In eq 3 there are p pairs of electrons in each cell that have been expanded in this manner, and we may retain the cell and pair designations on each of the nonorthogonal orbitals Xp(Rj) as well. A simple example of a lattice defined in this manner would be a hypothetical array of ethylene molecules with one molecule per cell. Including correlation effects via the perfect-pairing approach only in the s-orbital of each molecule would yield p = 1 and n = 7 (a Is orbital on each carbon, four CH bonds, and the CC u-bond define the $,(k)) and the GVB orbitals in eq 3 would appear explicitly as [X181r2(Rj)- hz,e,,2(Rj)](aP - Pa) = [xr(Rj)x,’(Rj) + x,’(Rj)x,(Rj)l(aP - Pa) (9) where x,(R,) and x,’(R,) denote nonorthogonal rr-orbitals
“Perfect-Pairing” Electron Correlation Effects in Crystalline Systems
that become localized predominantly on the individual carbon atoms in order to reduce Coulombic repulsions and minimize ionic contributions to the wave function. As is suggested by the above discussion, the problem of optimizing both CS and GVB orbitals for crystalline systems requires that one determine those areas where a formulation in terms of Bloch orbitals or in terms local orbitals would be most advantageous. In the following sections, we will attempt to demonstrate that a formulation based on the use of local interactions throughout possesses certain practical advantages. 111. Total Energies and Eigenvalues in a Mixed Bloch a n d Local Orbital Representation The total electronic energy for the general lattice GVB wave function (eq 3) may be written cs cs ELAT = C(2hnn(k,k) + C [2Jnn’(k,k’) - Knn!(k,k’)]l kn
k‘n‘
The Journal of Physical Chemistry, Vol. 87, No. 20, 1983 3067
Olp(Rj,r) = Olp(R:,r
+ Rj’ - RJ)
(12)
and must form complete orthonormal sets satisfying (~lp(R;)IOpp4Rj’))= Jnip6i~!6ppJ;y
As such, they may (in crystalline problems) be identified as Wannier orbitals, with uniform partial occupation Alp2 for all members of each set of orbitals OLp(RJ).We may expand the members of a set in Bloch functions (using discrete normalization) as BZ
m
Oip(RJ)= N-1~2Cexp(ik.RJ),Nip-1~2(k) CC$,(k)& k
(13)
a
where the sum ranges over all values of k within the first Brillouin zone, and the expansion coefficients are related to the local orbital coefficients by N
Nip-1/2C$,(k)= N-1/2JnLp-112CC&(R,) exp(ik.R,)
(14a)
a
(loa)
BZ
JnLp-’~2C$,(R,)= N-’/2C,NLp-1/2(k)C;(k)exp(-ik.R,) k
(14b) N pairs 2
[2J~~i~~dR;,Rj’) - K L ~ ~ ~ ~ R11; ,( R 1 0~~’ )) CS N GVB
+ 2 CC C iX~p~[2Jnip(k,Rj)- Knip(k,Rj)ll (10d) kn I
pl
where the expression has been partitioned into k dependent terms for the CS Bloch orbitals (loa), purely local terms for the origin GVB orbitals (lob), interactions between GVB orbitals at the origin and other cells (lOc), and a set of field terms requiring a mixture of both representations (10d). Quantities appearing in this expression have the usual meanings, with one-electron matrix elements given by hnn((k,k’) = ($n(k)[hl$nt(k’)) hiprp4R;,Rj’) = (OLp(R;)lhlOitp4Rj’))
(11)
where
A n expansion of this type in terms of k states suggests that we might view the wave function as if it was a Bloch wave function consisting of complete bands of orbitals in which all states were partially occupied. The differences here are significant, however: for any value of p, all orbitals in bands identified by different values of 1 will be coupled irrespective of k. This may be seen by substituting the expression 14b for each orbital Olp(R,)in eq 3 and collecting the resultant determinants in orders of k. The Hamiltonians appropriate for such an expansion will not be totally symmetric, rather they will depend on k and thus link orbitals of different k to produce totally symmetric integrands. This is in marked contrast to the single-determinant Bloch wave function (where the Hamiltonian is totally symmetric and only links orbitals with the same value of k), and provides insight into the way in which the inclusion of electron correlation leads to deviations from the Bloch picture. Expanding the lattice electron density due to closed-shell orbitals in (3) we find
N
h = -Y2V2 - CZ,(lR, - rl)-’ (2
S:f =
and two-electron matrix elements Jnn,(k,k’) = ($n(k)lJnl(k’)l$n(k)) = ($n(k)$n(k) l$n((k’)$n((k’)) Knn&,k’) = (+n(k)lKnT&’)l$n(k)) = ($n(k)$n