A modified Stillinger-Weber potential for modelling silicon surfaces

surface s c i e n c e ELSEVIER

Surface Science 366 (1996) 177-184

A modified Stillinger-Weber potential for modelling silicon surfaces P . C . L , S t e p h e n s o n 1, M . W . R a d n y , P . V . S m i t h * Physics Department, University of Newcastle, Callaghan NSW2308, Australia

Received 8 May 1995; accepted for publication 23 April 1996

Abstract

The widely used Stillinger-Weber potential for silicon interactions has been modified to provide an accurate description of the Si(lll)q7 x 7) surface including the highly reactive adatom and rest-atom sites. This modified potential also provides a good representation of bulk silicon, and the Si(001)-(1 x 1), Si(001)-(2 x 1), Si( 111)-( 1 x 1t and Si( 111)-(2 x 1) surfaces. Above the melting temperature of 1683 K, however, the original Stillinger-Weber potential is probably superior. Keywords: Adatoms; Computer simulations; Low index single crystal surfaces; Molecular dynamics; Semi-empirical models and model calculations; Silicon; Surface relaxation and reconstruction; Surface structure, morphology, roughness, and topography

1. Introduction

Many molecular dynamics (MD) simulations use analytic potentials to model the interactions between atoms. A commonly employed potential for describing silicon interactions is the Stillinger-Weber (SW) potential [1]. This potential has been used to simulate bulk silicon [ ! ], the (001) and (111) surfaces of silicon [ 2 - 5 ] , and the etching of the Si(001)-(2 x 1) surface by fluorine [ 5 - 8 ] . While providing a good representation of the (1 x 1) and ( 2 x 1) (001) and (111) surfaces of silicon, this potential fails to provide an accurate description of the adatoms of the S i ( l l l ) - ( 7 x 7 ) surface [9]. The latter is important because initial

* Corresponding author. Fax: + 61 49 216907; e-mail: [email protected] i Present address: Physics Laboratory, University of Kent at Canterbury, CT2 7NR, UK.

chemisorption on the (7 x 7) reconstructed Si(111) surface occurs predominantly at the dangling-bond adatom and rest-atom sites. Recent ab initio calculations, for example, have shown that fluorine and chlorine atoms are readily chemisorbed near the adatom and rest-atom sites of the Si(111)-(7 x 7) surface and form the various SiFx and SiClx species which are the precursors to silicon surface etching [ 10,11]. Thus, in order to simulate processes such as the etching of the S i ( l l l ) - ( 7 x 7 ) surface by atomic fluorine and chlorine, it is essential that the adatom and rest-atom sites, and their immediate environment, be modelled accurately. The aim of the present work, therefore, is to modify the original SW potential to provide a good description of the highly complex Si( 111 )-(7 x 7) surface, in addition to the bulk solid and silicon surfaces already well modelled by the SW potential. This will provide the basis for reliable M D simulations of the S i ( l l l ) - ( 7 x 7 ) surface and the various reac-

0039-6028/96/$15.00 Copyright © 1996 Elsevier Science B.V. All rights reserved PH S0039-6028 (96) 00801-1

P. C.L. Stephenson et al./Surface Science 366 (1996) 177-184

178

tions that occur on this very important surface. In particular, it will enable simulations of the interaction of fluorine with the Si(lll)-(7 ×7) surface analogous to the original SW potential calculations for the (001) surface of silicon [ 5-8]. Such calculations constitute a natural extension of the abovementioned total-energy calculations and should permit these halogen chemisorption and etching processes, which are very important technologically, to be studied dynamically at a variety of different temperatures and exposures.

2. The Stiilinger-Weber potential The SW potential for interactions between silicon atoms can be written as the sum of two-body and three-body interaction terms, denoted f2 and f3, respectively [1]. In terms of the energy and length units e=2.1672 eV and a=2.0951 A, these can be written as

{A(Br~4-1)exp[~(rij-a)-t],rij<_a f2(rij) =

diamond structure as the most stable periodic arrangement of atoms at low pressure, and that the melting point and fluid structure, as calculated from molecular dynamics simulations, are in agreement with experimental results. No surface properties were considered in this fitting procedure, although some surfaces are reasonably well described by this SW potential. Calculations with this potential provide a good description of the rest-atom sites on the Si(lll)-(7 x 7) surface but predict the adatoms to be approximately 0.6 ,~ higher than the empirically determined LEED positions (see Table 1). This results in significantly larger bond lengths and bond angles for the adatom sites than the values of around 2.35 A and 60 ° obtained from fitting the LEED data [12]. The adatom binding energies yielded by this potential are also much smaller than the results obtained from first-principles calculations (see Table 2). It thus follows that some modification of the SW potential is required to provide an adequate description of both the rest-atom and adatom sites of the Si(lll)-(7 x 7) surface.

otherwise

(1) and

3. The modified potential

f 3(ri,rj,rk) = h(rij,rik,Ojik) + h(rij,rkj,Oijk) + h(rik,rkj,Oikj),

(2) where

h(r,s,O) 2 [b(cos0 + k) z - c] exp{~ [(r - a)- 1 + ( s - a)- 13}, r and s_
(3) In these expressions, rij is the distance (in reduced units) between atoms i and j, centred at vectors ri and rj, and O~jk is the angle at rj subtended by the lines to the other two points, i and k. For the original Stillinger-Weber potential the parameters have the values A=7.049 556 277, B=0.602 224 558 4, ~ = 1.0, a = 1.80, 2 = 21.0, b = 1.0, k = 1/3, c = 0 and ~ = 1.20. The criteria for selecting the values of the above parameters were that they give the

In order to determine this modified Stillinger-Weber potential, a number of different bulk and surface properties of silicon have been considered. These include some of the bulk lattice dynamical properties, the equilibrium geometries of the (111)-(2 x 1), (111)-(7 x 7) and (001)-(2 x 1) surfaces, and the adatom and rest-atom binding energies on the S i ( l l l ) - ( 7 x 7 ) surface. Attention was also paid to ensuring that the diamond structure was the most stable phase. To obtain a good description of all of these properties it was necessary to vary both the two-body and three-body terms (Eqs.(1)-(3)). The resulting parameter values are A=6.8932, B=0.3535, (=0.7872, a = 1.80, 2=21.0, b=0.0980, k=0.0500, c = - 0 . 0 3 2 4 and ~ = 0.5778. The results predicted by this modified SW potential for the different properties of silicon are presented below, and compared with the values obtained from the original SW potential, other theoretical calculations, and experiment.

P. C L, Stephenson et al./Surface Science 366 (1996) 177 184

179

Table 1

Geometrical properties (in .~) characterising the adatom, rest-atom, and dimer configurations of the Si(111)-(7 x 7) reconstructed surface (see Fig. 1) as determined by LEED [12], ab initio pseudopotential calculations [13,14], and the original and modified Stillinger-Weber empirical potentials (present work); the adatom and rest-atom inter-atomic distances li and interlayer spacings d~ have been determined by averaging over the four adatoms and two rest-atoms which lie along the long diagonal of the Si( 111)-{7× 7) unit cell (see Fig. 1 of Ref, [ 13]); the dimer 11 and 12 values have been averaged over all nine dimers in the Si(lll )-(7 x 7) unit cell

Adatom

Rest-atom

Dimer

Parameter

Original SW potential

Modified SW potential

Ab initio pseudopotential

LEED

11 12 13 dl d2 d3 ll lz di d2 I~ Iz

2.713 3.638 2.446 1.749 1.277 2.288 2.387 3.963 0.679 2.403 2.426 2,455

2.483 3.730 2.43 l 1,280 1,060 2,385 2.365 3.805 0.894 2.382 2.398 2.455

2,490 3.630, 3.618 2.382 1.344 1.130 2.270 2.379 3.774 0.955 2.413 2.442 2.431

2.349 3.448 2.376 1.235 1.293 2.155 2.364 3.909 0.704 2.365 2.450 2.490

Table 2 Binding energies for the adatoms and rest-atoms of the Si(111)-(7 × 7) surface for both the original and modified potentials; also presented are the corresponding values obtained from experiment [ 15], and from semi-empirical AM1 [ 18] and ab initio H F-DFT calculations

Adatom A Adatom B Rest-atom

Original SW potential

Modified SW potential

AM1

3.36 3.19 6.39

5.86 5.56 7.05

6.76 6.12 7.85

Experiment

HF-DFT

7.80

6.34" 6.34 a 7.60 b

" From Ref. El6]. b From Ref. [17].

3.1, Silicon surfaces 3.1.1. The Si(111)-(7 x 7) surface To determine the equilibrium geometries of the Si( 111 )-(7 x 7) reconstructed surface resulting from the original and modified Stillinger-Weber potentials, molecular dynamics calculations employing the Verlet algorithm were performed. These calculations were carried out using an eight-layer slab and applying periodic boundary conditions to the (7 x 7) rhombic surface unit cell. The L E E D structure of Tong et al. [12] was employed as the starting geometry, and the positions of the atoms of the six topmost layers was allowed to vary while keeping the atoms of the remaining two b o t t o m layers fixed at their ideal (1 x 1) positions.

The geometrical parameters characterising the resulting minimum-energy configurations are presented in Table 1 and illustrated in Fig. 1. Also tabulated are the corresponding L E E D [12] and ab initio pseudopotential results 1-13,14]. It can be seen from Table 1 that the adatom geometry resulting from the modified SW potential is in very good agreement with the ab initio pseudopotential results. The original SW potential, on the other hand, produces a substantial outward movement of both the adatom and its first layer atoms, and hence yields values for 11, dl and d2 which are considerably larger than the first-principles results. The average height of an adatom above the ideal tetrahedral position is predicted by the pseudopotential calculations to be 0.48 [13] and 0.51

180

P. C.L. Stephenson et al./Surface Science 366 (1996) 177-184

adatom site (side view)

rest-atom site (side view)

surface dimer (top view)

Fig. 1. The rest-atom, adatom and dimer structures of the Si(lll)-(7x 7) surface, showing the geometrical parameters presented in Table 1.

[14]. These values are to be compared with 1.09 for the original SW potential, 0.57A for the modified potential, and the LEED value of 0.31 ,~. The rest-atom geometry predicted by the modified SW potential is also in very good agreement with the ab initio pseudopotential results, whereas the original SW potential is more consistent with the LEED data in predicting a substantially smaller value of the interlayer spacing dl, and a larger value for /2. The average height of a rest-atom above the ideal tetrahedral position has been determined by the pseudo.potential calculations to be 0.27 [13] and 0.20A [14]. The corresponding empirical potential and LEED values are 0.21 (modified SW potential), 0.07 (original SW potential), and 0.01 ,~ (LEED [12]). The original and modified potentials both yield Si( 111)(7 x 7) dimer geometries, in good agreement with the LEED and ab initio pseudopotential values. To determine the reconstruction energy for the S i ( l l l ) - ( 7 x 7 ) surface, the difference in energy between the optimally relaxed Si(111 )-(1 x 1 ) surface and the geometry-optimised Si(111)-(7×7) surface was calculated for both the original and modified SW potentials. The resulting energy values were 0.47 eV per (1 × 1 ) surface unit cell for the original SW potential, and -0.05 eV per (1 x 1) surface unit cell for the modified potential. The new potential thus predicts the (7 x 7) reconstruction of the Si(lll) surface to be more stable than the (1 x 1) geometry, in agreement with experiment, whilst the original SW potential incorrectly favours the (1 × 1) topology. The binding energies appropriate to an adatom

A (with one six-member and two adjacent fivemember rings), an adatom B (with one five-member and two six-member rings), and a rest-atom of the Si(lll)-(7 x 7) surface (see Fig. 1 of Ref. [13]) were also determined for both potentials. This was done by first minimising the energy of the system, and then repeating the same calculation with the adatom (or rest-atom) removed. The difference in these energies gives the required binding energy. The resulting values are presented in Table 2, from which it can be seen that the rest-atom binding energy predicted by the modified potential is 0.75eV lower than the 7.8 eV binding energy obtained from experiments on the Si(lll)-(1 x 11 surface [15], and 0.55 and 0.80eV below the theoretically determined values [17,18]. The adatom binding energies predicted by the modified potential are also only 0.48 and 0.78 eV below the ab initio value of 6.34eV [16], whereas those yielded by the original SW potential are much lower, at 3.0 eV below the ab initio value. To determine the accuracy with which the twc empirical SW potentials actually model the interactions at the adatom and rest-atom sites of the Si (111 )-(7 × 7) surface, calculations of the variation of the energy E with the height h of the adaton~ or rest-atom above the surrounding surface were also performed. Due to the 3.77 ,~ cut-off charac. terising both potentials, only the Si12 and Si19 clusters shown in Figs. 2a and 3a need to be considered in order to include all the adatomdependent and rest-atom-dependent two-body and three-body interactions, respectively. The resulting

P, C.L. Stephenson et al./Surface Science 366 (1996) 177 184

181

adjoin

041 .... • 0.3

\

:J/

0.4

'i modified SW j ". ofiginalSW ....... /' / '.. HF o /

0.3

i

i

i

i

i

!

modified SW original SWHF o

o

°

// o ,,,

b13

0.2

&

0.2 ',

0.1

o // /

', <>

0.1 0 -0.3 -0.2 -0.1 (b)

/

o

,,,'

~ ,

°'/// ", o 'X?

0

0.1

0.2

height (/~)

Fig. 2. (a) An adatom and its 11 nearest-neighbour atoms on the Si( 111 )-(7 x 7) surface. The geometry is that obtained from 6-31G*(3df,2p) HF calculations using the SisHv cluster formed from the shaded silicon atoms and allowing only the height of the adatom to vary from the LEED values. (b) The energy of an adatom as a function of its height with respect to its equilibrium position for the original and modified SW potentials. Also plotted are the corresponding SisHv:6-31G*(3df,2p) HF values. Each curve has been shifted so that its minimum lies at the origin.

<>"7 o/,// ',,o

0.3 0

I

I

o// ^~

-0.3 -0.2 -0.1 (b)

/

0

0.1

I

I

0.2

0.3

height (A)

Fig. 3. (a) A rest-atom and its 18 nearest-neighbour atoms on the Si(111)-(7 x 7) surface. The geometry is that obtained from 6-31 G* (3df,2p) HF calculations using the Si? H t 5 cluster formed from the shaded silicon atoms and allowing only the height of the rest-atom to vary from the LEED values. (b) The energy of a rest-atom as a function of its height with respect to its equilibrium position for the original and modified SW potentials. Also plotted are the corresponding SivH15:6-31G* HF values. Each curve has been shifted so that its minimum lies at the origin.

E(h)

curves for both potentials are shown in Figs. 2b and 3b. In these calculations, all of the atoms were kept fixed at their LEED-determined positions [ 12] except for the adatom or rest-atom, which was moved in the direction normal to the surface. Also plotted in Figs. 2b and 3b are the results of corresponding Hartree-Fock calculations. These have been performed using an SisHv (adatom) and SiTH:5 (rest-atom) cluster and

employing the GAUSSIAN92 program [-19] with the 6-31G*(3df,2p) basis set. These clusters have been formed by taking the shaded silicon atoms in Figs. 2a and 3a and adding hydrogen atoms along the directions to the next-nearest neighbour silicon atoms, at a distance of 1.5 A. Calculations with the smaller 3-21G* and 6-31G* basis sets and a larger Si12H15 adatom cluster, were also performed to

182

P. C.L. Stephensonet al./Surface Science366 (1996) 17~184

confirm the convergence of these results with respect to both the size of the cluster and the range of the basis functions. For each of the curves in Figs. 2b and 3b the energy is relative to its own minimum, and the height of the adatom or rest-atom is with respect to the corresponding equilibrium position. It can be seen that both the original and modified SW potentials provide a fairly good description of the energy of an Si(lll)-(7 x 7) adatom or rest-atom as a function of its height above the surrounding surface. For the adatom, the predictions of the modified potential are in almost perfect agreement with the ab initio HF results, whereas those of the original SW potential are somewhat too steep. For the rest-atom site, the original SW potential yields an E(h) curve which gives a reasonably good description of the HF results in the vicinity of the minimum but is too flat for large negative values of the height h, while the modified SW potential energy curve follows the overall shape of the HF curve more closely, but is too shallow for all h. 3.1.2. The Si(001)-(1 x I ) surface and the dimerised Si(001)-(2 x I) surface Both the original and modified SW potentials provide a good description of the Si(001)-(1 x 1) and dimerised Si(001)-(2 x 1) surfaces. The dimer bond lengths of the Si(001)-(2x 1) surface were 2.404 and 2.379 A for the original and modified SW potentials, respectively. Both of these values are close to the accepted value of 2.40 ~, [20]. The corresponding dimerisation energies are - 1.68 and -2.18 eV. Both potentials thus correctly predict that the dimerised surface is energetically more favourable than the Si(001)-(1 x 1) surface. The reconstruction energy predicted by the modified SW potential, however, is in slightly better agreement with the value of approximately - 2 . 0 eV determined by state-of-the-art ab initio pseudopotential calculations [21]. 3.1.3. The S i ( l l l ) - ( l x I ) surface and the n-bonded Si(111)-(2 x 1) surface Both the original and modified SW potentials produce re-bonded chain reconstructions for the S i ( l l l ) - ( 2 x l ) surface, analogous to those predicted by ab initio calculations [22]. The respective reconstruction energies are found to be 0.35 and

0.36 eV per surface atom. Both these values are positive, suggesting that the (2 x 1) surface is energetically less favourable than the (1 x 1) surface. While this is contrary to what is observed experimentally, such predictions are characteristic of most of the analytic potentials for silicon [23,24]. The Si(lll)-(1 x 1) surface is reproduced well by both potentials, although neither predicts any significant relaxation of the surface. 3.2. Bulk silicon

In order to determine the relative stability of the different possible bulk phases of silicon, the totalenergy versus density curves have been calculated for the diamond, simple-cubic (sc),//-Sn [25], bodycentred cubic (bcc) and face-centred cubic (fcc) structures. These E(p) curves are plotted in Fig. 4 for both the original and modified StiUinger-Weber potentials. The structure with the lowest energy, for both the original and modified SW potentials, was found to be the diamond structure with a density of p=0.460o--3. This corresponds to a nearestneighbour distance of 2.351 A, in agreement with the accepted value [26]. The corresponding minimum energy values are -4.33 and -5.27 eV per atom for the original and modified SW potentials, respectively. Both of these values are in reasonable agreement with the experimental value for the cohesive energy of silicon of 4.63 eV per atom. For the original SW potential, the diamond structure was found to be the lowest-energy configuration up to a density of approximately 0.6a -3, but becomes essentially degenerate with the sc structure at higher densities. For the modified potential, the diamond structure was always the lowest energy phase. Neither potential yields the experimentally observed phase transition to the fl-Sn structure at higher densities, although the original SW potential does predict almost identical energies for the/~-Sn and diamond structures near p=0.56a -3. For densities below this value, the diamond structure becomes increasingly more stable relative to the other structures, for both potentials. The diamond structure is thus the clearly preferred bulk phase at high temperatures. Phonon-mode frequencies were calculated from the difference in energy between the relaxed bulk-silicon system, and the same system with

P. C.L. Stephenson et al./Surface Science 366 (1996) 177 184

Table 3 Some Si phonon frequencies (in THz), elastic constants cn and c12 (in 10l° N m-2), and bulk modulus B1. (in 10 m N m 2), calculated using the original and modified Stillinger-Weber potential. Also shown are the corresponding experimental results at 296 K 1-27,28]

0 b.c.c. \ \

-1 fl-Sn e~ ©

r~

f.e.e.

-2

\~i<

-3

s.e

-4

•

,

..

-5 0.3

i

i

i

0.4

0.5

0.6

(al

0.7

density (or-3)

oI -2

,/: ....

',,

f.c.c.

/ iiiI

,'

-3 ~""

"""

/ / / b.c.c.

0.3

0.4

0.5

183

.."" , , y

(DLA (L) O)LO (L) O~TA(L) com (L) OJLAO ( X ) OOTA(X) O93o (X) ~OL,rO{F) c~l q2 Br

Experiment

Original SW

Modified SW

11.35 12.60 3.43 14.68 12.32 4.49 13.90 15.53 16,58 6.39 9.88

13.66 11.95 4.79 17.10 13.24 6.77 15.94 18.73 14.99 7.73 10.15

11.61 9.15 3.34 14.54 11.23 4.04 13.55 15.73 11.03 6.55 8.04

the original SW potential. The ctl elastic constant and the bulk modulus, however, are determined more accurately by the original potential. In order to determine an estimate for the melting temperature predicted by the modified potential, the procedure of Stillinger and Weber [1] was repeated to calculate the mean potential energy per atom as the temperature is systematically increased from 0 K. The melting temperature was found to be approximately 3800 K. This is to be compared with the experimental value of 1683 K [26] and the approximate value of 2370 K predicted by the original potential [1]. 3.3. Silicon clusters

(b)

0.6

0.7

density (a -a)

Fig. 4. The total potential energy (in eV per atom) for (a) the original Stillinger-Weber potential, and (b) the modified potential, as a function of density in reduced units (er-3), for five possible structural phases of Si.

displacements to the silicon atoms corresponding to each phonon mode. Also calculated were the elastic constants c n and q2, and the bulk modulus. All of these results are shown in Table 3. The values of the phonon-mode frequencies and the c12 elastic constant calculated using the modified potential are clearly better than those derived using

Ab initio calculations [29] show that the minimum energy configuration for a n Si 4 molecule is a rhombus. For the Si 6 duster, the lowest energy structure is predicted to be an edge-capped trigonal bipyramid structure with C2v symmetry. The equilibrium topologies of an Si4 molecule using the modified potential were found to include the rhombus as a stable structure, while those of the original potential did not. The Si 6 edge-capped trigonal bipyramid structure is predicted to be a local minimum by the original SW potential, but transforms into the closely related regular octahedron to yield the lowest-energy Si 6 structure for the modified potential.

184

P. C.L. Stephenson et al./Surface Science 366 (1996) 177-184

4. Conclusion In the present work the widely used Stillinger-Weber potential for silicon has been modified to provide a good representation of the Si(lll)-(7 x 7) surface. This has been achieved by retaining the original SW functional form and simply varying the parameters in the two- and three-body potentials. The resulting model potential yields results for the geometry and energetics of the Si(lll)-(7 × 7) surface which are in much better agreement with both experiment and ab initio calculations than the original SW potential. The solid-state silicon systems modelled accurately by the SW potential - bulk silicon and the Si(001) and Si(lll)-(1 × 1) and (2 x 1) surfaces - are also well described by the new potential. This modified potential, however, seems less satisfactory for modelling silicon at high temperatures (> 1600 K), and hence simulations of liquid silicon would probably best be performed using the original SW potential. This is not surprising, as the original SW potential was formulated with just such simulations in mind. The main objective of the present work, however, has been to develop an empirical potential which can be used to study the processes of silicon surface reconstruction, chemisorption and desorption, all of which occur predominantly at temperatures below 1000 K. For such studies, the modified potential appears to be superior.

Acknowledgements The authors acknowledge the Australian Research Council and the University of Newcastle for financial support during the course of this work. Time on the VP2200 supercomputer at the Supercomputer Facility of the Australian National University is also gratefully acknowledged.

References [1] F.H. Stillinger and T.A. Weber, Phys. Rev. B 31 (1985) 5262.

[2] F.F. Abraham and I.P. Batra, Surf. Sci. 163 (1985) L752. [3] K. Khor and S. Das Sharma, Phys. Rev. B 36 (1987) 7733. [4] P.C. Weakliem and E.A. Carter, J. Chem. Phys. 98 (1993) 737. [5] P.C. Weakliem, Z. Zhang and H. Metiu, Surf. Sci. 336 (1995) 303. [6] F.H. Stillinger and T.A. Weber, Phys. Rev. Lett. 62 (1989) 2144. [7] T.A. Schoolcraft and B.J. Garrison, J. Vac. Sci. Technol. A 8 (1990) 3496. [8] L.E. Carter, S. Khodabandeh, P.C. Weakliem and E.A. Carter, J. Chem. Phys. 100 (1994) 2277. [9] X-P. Li, G. Chen, P.B. Allen and J.Q. Broughton, Phys. Rev. B 38 (1988) 3331. [10] M.W. Radny, P.V. Smith and Pei-Lin Cao, Surf. Sci, to be published. [ 11 ] M.W. Radny, P.V. Smith and Pei-Lin Cao, Surf. Sci, to be published. [ 12] S.Y. Tong, H. Huang, C.M. Wei, W.E. Packard, F.K. Men, G. Glander and M.B. Webb, J. Vac. Sci. Technol. A 6 (1988) 615. [13] K.D. Brommer, M. Needels, B.E. Larson and J.D. Joannopoulos, Phys. Rev. Lett. 68 (1992) 1355. [ 14] I. Stich, M.C. Payne, R.D. King-Smith and J.-S. Lin, Phys. Rev. Lett. 68 (1992) 1351. [15] J. Dieleman, Le Vide Les Couches Minces, Suppl. 218 (1983) 3. [16] S. Wang, M.W. Radny and P.V. Smith, in preparation. [17] P.J. Van den Hock, W. Ravenek and E.J. Baerends, Phys. Rev. B 38 (1988) 12508. [18] P.V. Smith and Pei-Lin Cao, J. Phys.: Condens. Matter 7 (1995) 7125. [19] M.J. Frisch, G.W. Trucks, M. Head-Gordon, P.M.W. Gill, M.W. Wong, J.B. Foresman, B.G. Johnson, H.B. Schlegel, M.A. Robb, E.S. Replogle, R. Gomperts, J.L. Andres, K. Raghavachari, J.S. Binkley, C. Gonzalez, R.L. Martin, D.J. Fox, D.J. Defrees, J. Baker, J.J.P. Stewart and J.A. Pople, GAUSSIAN92 (Revision E.3) (Gaussian Inc., Pittsburgh, PA, 1992). [20] Z. Jing and J.L. Whitten, Surf. Sci. 274 (1992) 106. [21] H. Balamane, T. Halicioglu and W.A. Tiller, Phys. Rev. B 46 (1992) 2250. [22] J.E. Northrup and M.L. Cohen, J. Vac. Sci. Technol. 21 (1982) 333; Physica B 117/118 (1983) 774. [23] J. Tersoff, Phys. Rev. B 37 (1988) 6991. [24] J. Tersoff, Phys. Rev. B 38 (1988) 9902. [25] J.C. Jamieson, Science 139 (1963) 762. [26] O. Madelung, Landolt-B6rnstein: Numerical Data and Functional Relationships in Science and Technology, New Series III/22a Eds. O. Madelung and M. Schulz (Berlin, Springer, 1987) p. 14. [27] G. Dolling, in: Inelastic Scattering of Neutrons in Solids and Liquids II (IAEA, Vienna, 1963) p. 37. 1-28] H.J. McSkimin and P. Andreatch, Jr., J. Appl. Phys. 35 (1964) 2161. [29] K. Raghavachari, J. Chem. Phys. 84 (1986) 5672.