A density-functional based tight-binding approach to 111-V semiconductor clusters Joachim Elsner Michael Haugk Gerd Jungnickel and. Thomas Frauenheim Technzsche Unzversztat lnstitut fur Physzk,. Theoretzsche Physzk Ill D-09107 Chemnztz Germany We determine here the Hamiltonian and overlap matrix elements necessary in a non-orthogonal two-centre tight-binding (TB) scheme for gallium aluminium arsenic gallium arsenide and aluminium arsenide using a local density-functional method. The repulsive energy required by such methods is described by radial short-range repulsive pair contributions. These are determined with respect to self-consistent cohesive energy calculations performed on diatomic molecules and bulk crystalline forms Without inclusion of any parameters we show that the method can be applied successfully to simulations of small homo- and hetero- nuclear clusters giving an accuracy comparable to ub initzo calculations In addition we prove the transferability of the scheme to the crystalline modifications of As Al Ga GaAs and AlAs by determining the equilibrium structural vibrational and electronic band-structure properties There is great technological interest in 111-V compound semi- conductors in particular in the study of epitaxial growth mechanisms and interface formation which has highlighted a need for an understanding of the properties of such materials at a fundamental level.Theoretical models have been developed to characterize such diverse problems as small homo- and hetero-nuclear clusters' and crystalline surfaces.These methods range from empirically constructed potentials derived from fitting parameters to equilibrium structures to fully self- consistent ub inrtzo schemes ' Although the former are very efficient and thus capable of dealing with extended (amorphous or Crystalline) systems they suffer from a transferability prob- lem. They generally only work well within the regime in which they were fitted and thus are not predictive in structural simulations In contrast the problem of transferability is solved in general by using fully self-consistent ab initio schemes on the basis of density-functional (DF) theory '. These methods can charac- terize very accurately the equilibrium structures and properties Although a lot of progress has been made in applying these methods to ever-larger systems they are still far too slow to examine many interesting problems for example unconstrained surface reconstructions epitaxial growth and nucleation phenomena Therefore more approximate schemes combining the advan- tage of the efficiency of the empirical potentials with the transferability and accuracy of the self-consistent field (SCF) methods are highly desirable In this context tight-binding (TB) models have become very popular recently,' lo providing one of the most accurate alternatives in the determination of the total energy and equilibrium geometry of various systems In particular two-centre-oriented schemes considering only two-centre integrals in the Hamiltonian give results that deviate only slightly from those of more sophisticated methods However the usual procedure of fitting the matrix elements necessary to calculate the band-structure energy is to an arbitrary set of input data It is thus rather complicated and does not guarantee a general transferability to all scale systems Nevertheless good results have been achieved by using empiri- cal TB models for the description of GaAs surface reconstruc- tions" and a-GaAsI2 However it is not clear whether the same TB formulations apply to small GaAs clusters too Our method based on density-functional theory using a minimal basis linear combination of atomic orbitals ( LCAO) framew~rk,'~avoids difficulties arising from an empirical parametrisation Instead we calculate the two-centre Hamiltonian and overlap matrix elements from atom-centred valence-electron orbitals and atomic potentials derived from SCF single-atom calculations within the local-density approxi- mation (LDA). This yields exactly the same energy expression as common parametrised TB schemes but with the important difference of having a well defined procedure for determining the desired matrix elements Recently this method has been applied successfully to the description of the equilibrium structure and properties of various scale C Si and BN struc-tures including interactions with hydrogen,14 l6 as well In addressing possible applications to other 111-V compound semiconductors we now apply the same DF-TB method to Ga As Al GaAs and AlAs clusters and verify the transfer- ability to crystalline forms This paper is organized as follows First we outline briefly the theoretical background of the method followed by a description of how the two-centre Hamiltonian and overlap matrix elements and the short-range repulsive potential are determined In the next section we apply the method to small Ga As Al GaAs and AlAs clusters and present results of the ground-state configurations To verify the transferability of the DF-TB scheme to extended crystalline forms we then calculate the equilibrium bulk properties (cohesive energy bulk modulus lattice constant) and compare our derived band structure with the results of more sophisticated methods Finally we summar- ize our results and give suggestions for further research Method Our method based on the work of Seifert Eschrig and Bieger,13l7 applies the formalism of optimized linear combi- nation of atomic orbitals (0-LCAO) as introduced by Eschrig and Bergert for band-structure calculations '' In this approxi- mation the Kohn-Sham orbitals $ of the system containing 1 atoms are expanded in terms of atom-centred localized basis functions q5v $49= 1Cl"Z4JY-4) (1) lv where I= 1 2 etc is the number of atoms and z= 1 2 etc is the number of Kohn-Sham orbitals By restricting to a minimal basis the q5y represent only valence-electron orbitals which are expanded in terms of Slater-type functions and spherical harmonics As many tests have shown,lg five different values of a and n = 0 1 2 3 form a sufficiently accurate basis set for all elements up to the third row in the periodic table J Muter Chem 1996,6(lo) 1649-1656 1649 Using eqn (2) we perform a self-consistent solution of the modified all-electron single-atom Kohn-Sham equations [T+ vpsat(r)]4v(r)=E,psat(bv(r) (3) vpsat(r)vnucleus(r) + JGartreeCn(r)I + vxcLDACn(r)I += tJ (4) and determine the single-atom electron wavefunctions and potentials. The exchange correlation potential Vx is expressed in terms of the local density approximation as parameterized by Perdew and Zunger 2o.The additional term (r/r,)Nappearing in V(r)in eqn (4) was first introduced by Eschrig et a1 l9 in order to improve the band-structure calculations performed within LCAO It forces the wavefunctions to avoid areas far from the nucleus thus resulting in an electron density that is compressed in comparison to the free atom.The radius r may be optimized to yield the best results however we have found that ro z 2r,, where r, is the covalent radius of the element and N = 2 is an overall good choice Using eqn (1) the Kohn-Sham equations of the many-atom structure are solved within a non-self-consistent treatment in the next step fi$,(V) = &,$,W fi= T+ Kff(4 (5) As an approximation we express the effective one-electron potential V,,Ar) of the many-atom structure as a sum of spherical atomic contributions where V is the Kohn-Sham potential of a neutral pseudo- atom due to its compressed electron density but no longer containing the additional term (~/r~)~(For more details of eqn (6) the reader is referred to the paper of Porezag et a! 14) As a result the Kohn-Sham eqn (5)are transformed into a set of algebraic equations 1CI vr(Hflv'l-E,SpvI11 =0 vl I( = 1,2 etc (number of valence orbitals) (7) 1 = 1,2 etc (number of atoms) where HflVl1= (4jh -R,)IfiI4,(r -RID = (4p(r -R')I4,'k -R' 1) (8) Eqn (7) may also be written as a generalized eigenvalue problem HC = &,SC (9) with the hermitian Hamilton matrix H and the positive definite overlap matrix S are defined by eqn (8) The overlap matrix consists only of two-centre elements which can be calculated in a straightforward manner Consistent with eqn (6) one can neglect several contributions to the Hamiltonian matrix elements HpV,l7 yielding ifp=v and I=I v,' + v,' ) If 1 # 1' otherwise (10) This actually reduces the problem to a two-centre approach As has been shown already by several a~thors,~l-~~ the total energy of the system within this approach can be written with the usual tight-binding equation as a sum over the 'band- structure' energy (sum of the occupied Kohn-Sham orbital energies) and a short-range repulsive two-particle potential 1650 J Mater Chem 1996 6(lo) 1649-1656 Etot({Rk >)=EBS({Rk 1) + Erep({ IRk -R,1 1) = 2 niEt({Rk))+ 2 kp(lR1-R/rl) (11) I k <I where n is the occupation number of Kohn-Sham orbital The short-range repulsive contributions V,,,(R) can be deter- mined easily as the difference of the cohesive energy resulting from self-consistent total-energy calculations on molecular and crystalline reference systems and the related band-structure energy E, for different values of interatomic distances R %p(R) = ELDASCF(R) -E13S(R) (12) From eqn (11) we can now derive the interatomic forces FtJ= -2=1,2 N J=l,2,3w acting on the nuclei As one can show these derivatives can be expressed in the form Here the c are the eigenvectors of the single particle states with energies E Eqn (13) simply involves derivatives of the Hamiltonian and overlap matrices and the repulsive potential Hence for MD simulations one has to diagonalise eqn (9) to calculate the forces at each time step Determination of the tight-binding matrix elements As an example we discuss the construction of interatomic Hamiltonian and overlap matrix elements for GaAs In determining the atomic wavefunctions and potentials as described in eqn (3) and (6) we choose ro = 4 7a0 for gallium and ro=45ao for arsenic as the confinement radius of the additional potential term All matrix elements of eqn (10) to be substituted into the general eigenvalue problem [eqn (9)] are evaluated for differ- ent internuclear distances r ranging from 2 0 to 10 Oao using a distance step of 0 02a0 and tabulated as Slater-Koster (SK) integrals.The corresponding curves for GaAs As and Ga are shown in Fig 1-3 t From these SK integrals which are determined only once the eigenvalue problem (9) for any geometry can be constructed by a uniform coordinate trans- formation 25 Hence we can calculate the required integrals (and their derivatives) from their pre-determined atomic orbitals for every geometry rather than fitting a simpler functional form Finally the repulsive potential is obtained by eqn (12) including SCF cohesive energy data of the diatomic molecule and experimental values for the bulk as refeiences [Fig (4)] Clusters Since we have included the SCF dimer data into the fit of our repulsive energy the diatomic properties of the SCF calculation are reproduced Aluminium and gallium clusters A13 Ga,.Configuration interaction (CT) calculations for the aluminium trimer find that the ground state is an isosceles very close to an equilateral triangle with a bond angle x of ?The related SK tables are available upon request for the authors Fig. 1 (u)Hamiltonian and (b)overlap matrix elements vs interatomic separation H, S for gallium arsenide 562' (ref 26) and 605" (ref 27) and bond lengths r=5u Jones' performed SCF LDA calculations and obtained an equilateral triangle with a bond length of 4 65a We also find that an isosceles triangle (a= 62") is the most stable structure but we obtain a slightly enlarged bond length (r = 5 394,).The linear chain (r = 5 16a,) is higher in energy by 0 3 eV Here again Jones' found a shorter bond length of r =4 86~0 For Ga various authors report many low-lying states with very different geometries Balasubramanian and Feng" find an isosceles triangle (a= 61 2" r = 4 88u,) for the ground state which has been confirmed by Meier et a12' but with slightly larger bond lengths In contrast Jones' predicts an equilateral triangle having considerably shorter bond lengths (r = 4 39u,) to be the ground state Further he reports that the most stable linear form (r = 4 56a0) is about 0 9 eV higher in energy In our calculations we determine a variety of geometries of nearly degenerdte energies also.There is a flat isosceles triangle (3 = 129" r = 9 la,) an equilateral triangle (r = 5 31Uo) the linear chain and a second isosceles triangle (a= 50 5" r = 4 42u0) all of them separated by only 0 04 eV Since the very small energy difference is surely within the error margin of the present method we are therefore not able to predict any configuration to be the most stable one A14 Ga,.For Al we have performed SCF LDA calculations indicating a very soft energy surface with respect to geometrical changes As a consequence any rhombus with angles between 72 and 90 IS stable to reasonable accuracy While Jones indeed determines the rhombus (D,,,bond angle a= 56 5") to be the most stable configuration,' we obtain a square for the ground state (D,,,5 14Uo) Pacchioni and K~utecky,~'and Meier et a129 also predict a rhombic ground state but the square Fig. 2 (a)Hamiltonian and (b)overlap matrix elements us interatomic separation H,, S for arsenic (E = -0 5193Eh,E* = -0 1954Eh) (D4h)is very close in energy (0005 and 0 04 eV higher respect-ively) In accord with Jones,' there are also a lot of higher lying local minima in the energy surface The most stable Ga cluster in our calculations is again a square (r = 5 06a,).This is in between the isoenergetically stable square (r =4 632a0) and rhombus (a= 71 6" r = 4 49Uo) found by Jones' and a square with significantly longer bond lengths (r= 5 31a0)obtained by CI 29 A slightly smaller second square is obtained to be metastable by the DF-TB method with a side length of 4 744 and an energy decrease of 0 19 eV In addition a rhombus (a= 50" AE = 1 eV) forms another metastable structure A] Ga,.. There have been several studies of the aluminium pentamer Whereas CI calculation^^^ and semi-empirical LCAO calculation^^^ found that pyramidal structures [Jahn-Teller distorted CZL:(ref 27) and C4c(ref 31)] were the most stable configurations Petterson et al 32 report a planar (C ) form (with bond lengths constrained to be equal) to be 0 2 eV more stable than the pyramid In contrast Jones' calculations' yield two low-lying structures with almost identical energies a substantially deformed pyramid and a planar structure (C,t) This result is also confirmed by the present calculations where we find that a slightly deformed Al pyramid is the most favourable structure but the planar structure [C Plate 1(a)] is higher in energy by only 0 25 eV In the case of Ga a square capped by an atom on one side [see Plate l(b)] is obtained as the most stable geometry in agreement with SCF-LDA results A planar structure [Plate l(a)] and a pyramid are metastable at slightly higher energies of 0 27 and 0 14 eV respectively Fig. 4 Short-range repulsive-pair potentials for Ga (rcUtoff= 5 2a,) As (rCutoff= 5 14ao)and GaAs (rcutoff= 5 1 8ao) Al Ga,.In our calculations for Al the regular octahedron Jones,' we find that the most stable planar structure consisting 0 (r = 5 07no)forms the most stable structure in accordance of two aligned rhombi with bond angle a M 60" In Plate 1(d) with the ah znztzo predictions of Petterson3' for symmetric is significantly higher (AE = 0 9 eV) in energy structures. Jones' and Jug et found that the octahedron. The situation In Ga IS quite similar to that in A16.We again undergoes a Jahn-Teller distortion to a more stable trigonal identify the regular octahedron 0 (r= 5 13ao) as the most antiprism stable structure almost degenerate with a prism structure. The A metastable prism structure [Plate 1 (c)] in the present most stable planar structure formed by two aligned squares is calculations is obtained at AE = 0.26 eV. In agreement with clearly higher in energy (0.54 eV) 1652 J. Muter Chem. 1996 6( lo) 1649-1656 Al Ga,. As in the calculations of Jug et Raghavachari and Jones,' we confirm the Al ground state to be a C3c structure which can be viewed as a 'capped' antiprism shown in Plate l(e). In the energetic order it is followed by a planar structure [see Plate l(f)] at AE = 1.4 eV.For Ga the ground state is the same as for Al and the planar structure lies 1.33 eV above the ground state. Al Gag Al Ga,. Results on larger clusters Al Al and Ga Gag have been reported by SCF-LDA' and semi-empiri- cal LCA031 calculations. In accordance with Jones we find an Al ground state shown in Plate l(g).. This has a shape similar to a face-centred cubic (fcc) lattice primitive cell but is quite distorted bond angles and bond lengths vary from 50 to SO" and 4.91 to 5.05~~~respectively whereas they are 60" and 5.41~1,in the fcc lattice primitive cell. For Ga our calculations yield a very similar structure but with slightly different bond lengths. Finally the lowest lying Al structure shown in Plate l(h) is in good agreement with the structure predicted by Jones.' Again a similar geometry is stable for Ga confirming the result of Jones.' In Fig.5 we have plotted the dependence of the cohesive energies per atom on the cluster size. As can be seen the DF-TB scheme yields energies which are in a good agreement with the SCF-LDA calculations of Jones.'. The energy devi- ations for some of the medium-sized clusters are probably due to the fact that these energy surfaces are very complex thus permitting a variety of local minima. Arsenic clusters As,.. The arsenic trimer has been identified in gas-phase charge-transfer reaction^.,^. Theoretical studies at the SCF-LDA and HF-CI levels have been performed by Jones,2 Igel-Mann et a/.,' and Balasubramanian et uL3 In agreement with Jones,2 we determine an isosceles triangle (C21; CI z 58" r = 4.40~~)and a nearly degenerate D structure (r= 4.49~~) to be the most stable geometries..The linear chain (r = 4.29~~) lies at significantly (2.1 eV) higher energy. As,. Since As4 is the most prominent component of arsenic vapour between 400 and 850 K,37*38the tetramer has received considerable attention. Gas-phase diffraction measurements of .~~Morino et ~ 1find a tetrahedron with an interatomic separa- tion of r = 4.602 & 0.008~~as the most stable configuration. This is supported by local spin density (LSD) calculations which give 4.613 (ref.40) and 4.56~~~~while CI calculations lead to slightly enlarged bond lengths 4.73 (ref. 41) and 4.67~1,~~~respectively.. The present calculations also yield a tetrahedron with r = 4.652 followed by a 'roof' structure Plate 2(u) at AE = 1.5 eV and a square (r= 4.49~~)which however lies 2.0 eV above the ground state.As,. In accordance with Jones,2 we find that the As ground state is a C2" structure as shown in Plate 2(b).. The slightly distorted pentagon a planar CZ0 structure with bond lengths between 4.40 and 4.49~1,~and the regular pentagon (D,,,r= 4.55~~)are almost isoenergetic metastable configurations at 0.8 eV above the ground state. As,.. The present calculations give several local minima in the energy surface of As,.. The lowest lying structure is a D configuration (trigonal prism).. The regular planar hexagon (D6h) lies 2.74 eV higher.Other local minima are described by a buckled structure shown in Plate 2(c) and a distorted octahedron (bond lengths 4.52-5.12~~)found at AE = 2.1 and 4.5 eV respectively. A%. In agreement with Jones,2 our calculations yield a C structure [Plate 2(d)] corresponding to a local minimum in the energy surface. However a capped prism [C2,; Plate 2(e)] is found to be 0.2 eV more stable. As,. For As we report a cubic oh structure (bond length 4.58~~)to be the most stable geometry. Jones2 predicts that although the cubic structure corresponds to a local minimum a structure of the form Plate 2(f) is clearly more stable. In our calculations this geometry is at 0.40eV above the cubic ground state. 3 4 5 6.7 8 9 cluster size Fig.5Energies of small clusters for Ga (a) and A1 (b). 0,from ref. 1 - this work. Plate 2 (u)-(g) Equilibrium structures for arsenic clusters J. Muter. Chem. 1996 6(lo) 1649-1656 1653 As,. As usual there is a variety of local minima on the energy surface of clusters of this size.. The geometry correspond- ing to the lowest energy determined for As is shown in Plate 2(g). It can be derived from the cubic structure of As by adding a bridging atom with two-fold coordination. Some examples of Ga,As and Al,As clusters There are of course a large number of possible isomers for the binary clusters. We see however that chemical bonding constraints reduce significantly the number of favourable geo- metries. Our calculations were performed with starting geo- metries suggested by Andreoni and some highly symmetric structures.Ga2As2 Al,As,. For the smallest reported heteronuclear GaAs cluster we confirm the ab initio result of Andreoni in obtaining a planar rhombus as the most stable isomer. In this cluster the As-As (4.44~~)bond is clearly favoured and the As-Ga-As bond angle is ca. 51" (ca. 52O in ref. 3). Another planar structure [Plate 3(a)],is a local minimum for Ga,As but the lowest lying geometry of A1,As2 lies 0.21 eV above. The energetically lowest three-dimensional cluster for Ga As is a roof structure at 1.4 eV above the ground state. Our calculations confirm that the same structures also correspond to local minima on the Al,As energy surface. However the energetic ordering is changed..The minimum energy configuration has the geometry shown in Plate 3(a). The Al-1 Al-As and As-As bond lengths are 5.04 5.24 Plate 3 (a)-(c) Stable structures for GaAs blue spheres represent gallium atoms the red spheres arsenic 1654 J. Muter. Chem. 1996 6( lo) 1649-1656 and 4.40a0 respectively.. The planar rhombus the most stable structure for Ga,As is now found to be 0.2eV higher in energy. Ga,As Al,As,. In accordance with Andreoni3 we report the hexamer to have two low lying states (AEz 0.35 eV).. The lower one can be viewed as an edge-capped trigonal bipyramid [Plate 3(h)].. The Ga-Ga bond (I' = 5.83~~)is much longer than the As-As and the Ga-As bonds which are 4.5 and 5.1a0 respectively.. The second structure is shown in Plate 3(c).Here the same large differences between the bond lengths occur. In both structures the Ga-As arrangement is not alternate one As atom having two As nearest neighbours which is consistent with indications from experiments on the singly ionized clusters. For the AI,As clusters we find that the same geometries are stable and the energetic ordering of the isomers is identical. Ga,As Al,As,.. The present calculations predict a slightly distorted octahedron for the ground-state geometry.. The bond lengths vary from 4.7 to 5.1U0 the shorter one corresponding to the As-As bond. A regular octahedron with bond length 5.13~~is 1.47 eV higher in energy. Here as in the previous structures we observe that the Ga-Ga and the Ga-As bonds are of the same length.In our most stable planar structure a slightly distorted hexagon however the Ga-Ga and Ga- As bond lengths differ considerably Ga -As being 0.6~~~ shorter. As in the case of other stoichiometric relations the same isomers in the same energetic ordering were obtained for Al,As also. Ga4As Al,As,. As for Ga,As we find a slightly distorted octahedron to be the ground state of Ga,As,.. The As-As bonds are again shorter than the Ga-As bonds by 0.6~'. Again a variety of local minima exist on the energy surface. The most stable linear structure a chain is determined at 8 eV above the ground state. The lowest lying geometry for Al,As is almost identical to that of Ga,As2. Most of the isomers tested for Ga,As were also found to have the same geometry for Al,As,.. There is however a slight change in their energetic ordering.Bulk properties and band-structure calculations In order to verify the transferability of the density-functional tight-binding approach to bulk structures and properties we have calculated the cohesive energy per GaAs and AlAs dimer in the zinc blende and rocksalt structures as a function of the interatomic distance.. The equilibria for the zinc blende struc- ture determined at an interatomic distance I' = 4.6~~for both compounds confirm very well the experimental results.42. The bulk moduli are calculated to be 7.87 x 10" Pa (GaAs) and 7.37 x 10" Pa (AlAs).. They are in good agreement with the experimental values of 7.56 x 10"Pa (GaAs) and 7.80 x 10'' Pa (A~As).~~ For the rocksalt structure we find the equilibrium interatomic distance at 5.12~~ (AIAs) and 5.13~~ (GaAs) with corresponding decreases of the cohesive energy of 0.72 and 0.66 eV per dimer relative to the zinc blende structure.While the cohesive energy differences between the two com- pounds as well as the energetic difference between the zinc blende and the rocksalt structure are described qualitatively well the absolute energies show the usual over-binding which is known to occur in LDA calculations.. The crystalline phases of aluminium gallium and arsenic are described at a compar- able level of accuracy. Furthermore we make use of the solutions of eqn. (2) as atomic basis functions in a usual tight-binding equation to calculate the band structures of Al As GaAs and AIAs..The resulting valence bands are in good agreement with those calculated with more sophisticated methods. However as can be expected from any sp3 tight-binding scheme our method produces poor results concerning the conduction bands As an example we show the derived band structures of AlAs in Fig 6 and As in Fig 7 In addition we have calculated the phonon densities of states (DOS) for GaAs and AlAs. The results are shown in Fig 8. The eigenvalues of both spectra were broadened with Lorentzians of 15 cm-' width A rough analysis reveals the following properties GaAs.. The vibrational density of states of zinc blende GaAs can be split into two main bands.The high wavenumber region between 200 and 300cm-' is clearly occupied by Ga-As stretching vibrations. The dominating peak is centred at 258 cm-' Bending vibrations can be assigned to the DOS at wavenumbers between 75 and 200cm-' Here we obtain two maxima at 116 and 181 cm-' At these low wavenumbers there are apart from bending vibrations also strong translational components of atomic groups occupying the modes Fig. 8 Vibrational density of states for AlAs (a) and GaAs (h) AIAs. For AlAs we obtain qualitatively the same DOS as -03 for GaAs Here the high wavenumber region arising from Al-As stretching vibrations varies from 250 to 350 cm-' with one principal peak lying at 329cm-' and a smaller peak at -04 289 cm-'. The wavenumbers due to bending vibrations range from 100 to 250cm-' with two maxima determined at 121 and 204 cm-' Again at these low frequencies we observe -05 translational modes of atomic groups 4-06 L G X U.K G Summary We have presented here a density-functional based non-ortho- gonal tight-binding scheme for Al Ga As AlAs and GaAs Fig.6 Calculated band structure for AlAs (E,= -0 13Eh) Fig. 7 Calculated band structure for As (E,= -0 07Eh) Though only two-centre integrals are considered the derived interatomic potential is highly transferable without additional changes to the total energy expression We have applied the method to small clusters of gallium aluminium and arsenic and their binary compounds Apart from a few cluster con- figurations we confirm the most stable geometries and the energetic order of metastable clusters as predicted by self- consistent LDA and ab rnztro quantum chemical calculations In order to verify the transferability to bulk systems we have calculated bulk crystalline properties such as the binding energies bulk moduli and phonon densities of states for GaAs and AlAs.These results are also in a good agreement with experimental data We have therefore demonstrated that despite its simplicity above all in its complete neglect of three- centre integrals the method is highly transferable giving reliable results for geometries cohesive energies phonon fre- quencies and band structures for clusters and solids and hence may be used in predictive MD simulations We are actually applying the method to perform a detailed study of the various surface reconstructions of differently oriented GaAs and AlAs In near future this will serve as the basis for molecular growth J Muter Chem 1996 6(10) 1649-1656 1655 simulations and investigations of interface formation in coupling with other semiconducting materials 19 20 H Eschrig Optimized LCAO Method and the Electronic Structure of Extended Systems Akademie-Verlag Berlin 1988 J P Perdew and A Zunger Phys Rev B 1981,23,5048 We would like to thank all of our group for useful discussions and in particular Dirk Porezag who supported our investi- gations on Al and Ga using an SCF full potential LDA code 21 22 23 24 G Seifert and R 0 Jones Z Phys D 1991,20,77 P Blaudeck T Frauenheim D Porezag G Seifert and E Fromm J Phys Condens Matter 1992,4,6389 W M C Foulkes and R Haydock Phys Rev B 1989,39,12521 D Tomanek and M A Schluter Phys Rev B 1987,36 1208.25 J C Slater and G F Koster Phys Rev 1954,94 1498 References 26 27 H Basch Chem Phys Lett 1987,136,289 T H Upton J Chem Phys 1987,86,7054 1 2 3 4 5 6 7 8 9 10 11 12 13 14 R 0 Jones J Chem Phys 1993,99,1194 P Ballone and R 0 Jones J Chem Phys ,1994,100,4941 W Andreoni Phys Rev B 1992,45,4203 G-X Qian R M Martin and D J Chadi Phys Rev B 1988 38,7649 T Ohno Phys Rev Lett 1992,70,631 E Kaxiras Y Bar-Yam J D Joannopoulos and K C Paudey Phys Rev B 1987,359625 R Car and M Parinello Phys Rev Lett ,1985,55,2471 K Laasonen and R M Nieminen J Phys Condens Matter 1990 2,1509 C H Xu C Z Wang C T Chan and K M Ho J Phys Condens Matter 1992,4,6047 M Menon and K R Subbaswamy Phys Rev B 1993,47,12754 C Mailhiot C B Duke and D J Chadi Surf Sci ,1985,149,366 C Molteni L Colombo and L Miglio Phys Rev B 1994 50 4371 G Seifert H Eschrig and W Bieger Z Phys Chem (Leipzig) 1986,267,529 D Porezag,.Th Frauenheim,. Th Kohler R Kaschner and 28 29 30 31 32 33 34 35 36 37 38 39 40 K Balasubramanian and P Y Feng Chem Phys Lett 1988 146,155 U Meier S D Peyerimhoff and F Grein Z Phys D 1990,17,209 G Pacchioni and J Koutecky Ber Bunsenges Phys Chem 1984 88,242 K Jug H P Schluff. H Kupka and R Iffert J Comput Chem 1988,9,803 L G M Petterson C W Bauschlicher Jr and T Halicioglu J Chem Phys 1987,87,2205 K Raghavachari Bull Am Phys Soc 1990,35,606 G Bachelet D R Hamann and M Schluter Phjs Rev B 1982 26,4 199 G Igel-Mann H Stoll and H Preuss Mol Phys 1993,80,325 K Balasubramanian K Sumathi and D Dai J Chem PIzys 1991 95,3494 J S Kane and J H Reynolds J Chem Phys 1956,25,342 J M Dyke,S Elbe1,A Morrisand J C H Stevens J Chem SOC Faradny Trans 2,1986,82,637 Y Morino T Ukaji and T Ito Bull Chem SOCJpn 1966,39,64 J Andzelm N Russo and D R Salahub Chem Phys Lett 1987 142,169 G Seifert Phys Rev B 1995,51 12947 41 U Meier S D Peyerimhoff and F Grein Chem Phys 1991 15.Th Frauenheim F Welch,. Th Kohler D Porezag G Seifert and 150,331. 16 S Uhlmann Phys Rev B 1995,52,11492 J Widany,. Th Frauenheim,. Th Kohler M Sternberg D Porezag and G Jungnickel Phys Rev B 1996,52,4443 42 43 J Ihm and J D Joannopoulos Phys Rev B 1981,24,4191 Landolt-Bornstein Elastic Piezoelectric and Related Constants vol III/l 1 of New Series Springer-Verlag Berlin 1979 17 18 G Seifert and H Eschrig Phys Status Solidi B 1985,127 573 H Eschrig and I Bergert Phys Status Solidi B 1978,90,621 Paper 6/00703A Received 30th January 1996 1656 J Muter Chem 1996 6(lo) 1649-1656