A theoretical model for the epitaxial growth of Si(100) surface

A TH EO M ETICAL М О Ш Ѣ ІЮ Ж Т Ш

E P IT A X IA L GRCßWTU O F SI (1 0 0 ) SIT'âFAGS

'Гх^ о V ,i»> s*4 Μ «rtl iш т в "''•’\··>\Γ^ ·ί: %-.У^ГУ. B G x m i < 4 . ^ · ι Λ . * * - ι * 4т ' ¡ W ' ЧѴ< Ч i o ч ' . w w i i Ä V i i . « S . ·é^^

« ¿ «¿^/¿іікмі іііі . А A »» i ^ . iÁi<¿ «іо" '<»»/' ¿L&MJ ’4 оМі ^ vhÁ· JmJtjtfJUİJİ Á!^'á^U.ji ч . А Ό / ¿4^

FQïe THE I9SGEBS OF

■

■

■«(.«■ ,4 Çr'-n'Sî'.'S!* A i? •аі'*'=Г'?5»?СА'®’_{•ΜΐΐΜΑ-ώ} μ··4μ«·' ^ A U . İİİHIV '*» > ОІІ . W '■ »Î «if ί « > ■'■ Хѵ'.Ькні4

<^.5?Д ■^ψΛ .iÉ«r* 4 Ъ n h àM i 'Чл·

;■ á .) ·■ ■ I Ή '■/

A THEORETICAL MODEL FOR THE

EPITAXIAL G R O W TH OF SI(IOO) SURFACE

A THESIS

SUBMITTED TO THE DEPARTMENT OF PHYSICS AND THE INSTITUTE OF ENGINEERING AND SCIENCE

OF BILKENT UNIVERSITY

IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF

MASTER OF SCIENCE

By

Sefa Dag

January 2000

é i i , S

b34

S¡..ooa

I certify that I have read this thesis and that in my opinion it is fully adequate, in scope and in quality, as a dissertation for the degree of Master of Science.

Prof. Salim Çıracı (Supervisor) I certify that I have read this thesis and that in my opinion it is fully adequate, in scope and in quality, as a dissertation for the degree of Master of Science.

Prof. Cemal Yalabık

I certify that I have read this thesis and that in my opinion it is fully adequate, in scope and in quality, as a dissertation for the degree of Master of Science.

— 3 - 'G S S 2 .S -S C

Prof. Şinasi Ellialtioğlu

Approved for the Institute of Engineering and Science:

/ / f-/ 06

Prof. Mehmet Baray,

Abstract

A THEORETICAL MODEL FOR THE EPITAXIAL

GROWTH OF SI(IOO) SURFACE

Sefa Dağ

M. S. in Physics

Supervisor: Prof. Salim Çıracı

January 2000

In this thesis, we investigated atomic scale mechanisms of growth on the Si(lOO) surface. First principles quantum molecular dynamics calculations are performed, where self-consistent field electronic structure calculations within the density functional theory, and ionic relaxations carried out concomitantly. By using this method, we first determined the binding energy and binding sites of ad­ atoms and ad-dimers. In the second stage, we simulated different configuration of dimer construction. In the third stage, we analyzed temperature dependent behaviour of dimers and we investigated translational and rotational motion of dimers. Finally, we presented a new theoretical model for ID epitaxial growth.

K e y w o r d s : Silicon, dimer bonds, temperature, epitaxial growth, ab-initio

özet

SI(IOO) YÜZEYİ ÜZERİNDE EPİK ATM AN BÜYÜME

MODELİ

Sefa Dağ

Fizik Yüksek Lisans

Tez Yöneticisi: Prof. Salim Çıracı

Ocak 2000

Bu çalışma, genel olarak Si(lOO) yüzeyinde atomik skaladaki büyüme

mekanizmalarını ele almaktadır. Yüzeyin yapısal özelliklerinin ve sıcaklığa

göre davranışının irdelenrnesinde “ y46-/nîh'o-Kuvantum Moleküler Dinamik” hesapları kullanılmaktadır. Bu hesaplarda “Yoğunluk Fonksiyonel” yaklaşımına dayalı kendi içinde tutarlı elektronik yapı hesapları ve iyonik durulma

hesapları aynı anda yapılmaktadır. Bu yöntemi kullanarak, ilk etapta OA^

sıcaklığında üzerinde tek atom bulunan ve daha sonra sadece dimer bulunan yüzeyler için yeni atomların en uygun bağ konligürasyonlarmm ve konumlarının bulunması ile ilgili çalışmayı ele aldık, ikinci etapta sonlu sıcaklıktaki yüzeyde ayrı iki atomun yüzey üzerinde değişik konfigürasyonlarda dimer oluşumları, üçüncü etapta ise oluşan dimerlerin sıcaklığa bağlı davranışları ve bunun sonucunda öteleme ve dönü hareketleri, son etapta da tek boyutta “epikatman” büyümenin oluşumu ve rnodellenrnesi incelenmiştir.

Anahtar

sözcükler: Silisyum, dimer bağlar, sıcaklık, epikatman büyüme, ab-initio

Acknowledgement

I would like to present my greatest respect and thanks to Prof. Salim Çıracı for his helps, suggestions, guiding, and advices.

I would like to give my endless thanks to Çetin Kılıç for his helps and beneficial discussions.

And I would like to give my thanks to Altug Ozpineci for his helps in the preparation of my thesis.

Contents

Abstract i

Özet i

Acknowledgement i

Contents i

List of Figures iii

List of Tables vi

1 IN T R O D U C T IO N 1

2 A T O M IC STR U CTU R E OF Si(lOO) SURFACE 9

2.1

Surface S tr u c tu r e ... 9

2.2

Surface Properties: Atomic M o t i o n ... 11

2.3 Thermodynamic Quantities About S u r fa c e ... 13

3 M O D ELLIN G THE IN TE R A C TIO N S 17

3.1 In trod u ction ... 17

3.2

Molecular Dynamics M e th o d ... 18

3.3 First Principles Molecular Dynamics 19

3.3.1 A pproxim ation s... 19

3.3.3 Exchange-Correlation Energy and LDA and GGA For­

malisms ...

22

3.3.4 Bloch’s Theorem and k-point s a m p lin g ... 23

3.3.5 Plane-wave represantation of Kohn-Sham equations . . . . 24

3.3.6 Electron-ion Interactions... 25

3.4 Car-Parrinello L a g ra n g ia n ... 28

3.4.1 Molecular Dynamics M e th o d ... 30

3.5 Relaxation of the ionic s y s t e m ... 32

3.5.1 T h e O th e r M e th o d s : Conjugate Gradient and Steepest Descent 36 4 H O M O E P IT A X IA L G R O W T H ON TH E Si(lOO) SURFACE 39 4.1 Binding Configurations of Si Ad-atom and Ad-dimer on the Si(lOO) S u rface... 39

4

.

1.1

Single atom adsorption on the Si(lOO) surface: Adatom . . 39

4

.

1.2

Addimers on the Si(lOO) su rfa ce... 42

4.2 Dimer F o r m a tio n ... 45 4.3 Surface K in e tics... 47 4.3.1 Dimer R o t a t io n ... 52 4.4 Epitaxial G r o w t h ... 54 5 CO N CLU SIO N 60

11

List of Figures

1.1 Positions of adatoms on the Si(001)-(2 x

1

) surface are designed

by M,H,B)D·· according to the labeling of Ref. 5. Adatom-adatom

combinations forming A , B, C , D type dimers are indicated by

arrows. Large, me dium and small filled, empty circles denote, respectively, the adatom (or ad-dimer), raised and lowered surface

dimer atoms and subsurface atoms...

8

2.1 Si(lOO) surface structures with symmetric dimers[p(l x

2

)] and

asymmetric (buckled) dimers[p(2 x 2) and c(2 x 4 ) ] ... 10

2.2 Potential energy of an adatom in a physisorbed state on a planar

surface as a function of its distance z from the surface...

12

2.3 Time-line plots. At the top are shown two plots of the ad-dimer

position along the row as a function of time acquired at 65 and

128 C° [Ref.

20

]... 14

2.4 Arrhenius plot of the logarithm of the diffusion rate versus 1/kT

from which the activation energy barrier is extracted [Ref. 20]. 15

2.5 Atom-track data taken at 128 C°. Site occupation probabilities

and the extracted relative free energies of the eleven sites are measured from six data sets, total time is 196 s. Error bars reflect

the statistics related to the number of hops at each site [ref.

20

]. 16

3.2

3.3

3.4

3.1 An illustration of the full all-electronic (AE) wavefunction and

electronic potential (solid lines) plotted against distance, r, from the atomic nucleus. The corresponding pseudowavefunction and pseudopotential are plotted(dashed line). Outside a given radius,

Гс, the all electron and pseudo electron values are identical. . . . 27

schematic representation of the evolution of coefficients { c } , Kohn- Sham Hamiltonian Я, and the total energy, in the final two time

steps of the molecular-dynamics method... 32

Flow chart describing the computational procedure for the

calculation of the total energy of a solid with molecular dynamics 33

A schematic view of the potential surface of E[4>i,Ri\. A typical

ionic relaxation trajectory is shown by the solid line which must lie close to the Born-Oppenheimer surface (dashed line) for the Hellmann-Feynman theorem to be implemented without incurring

large errors... 35

Top view of clean surface of Si: top layer is Si dimers layer and

bulk is saturated with hydrogens (bottom atoms)... 40

Optimized structure of Si(lOO) surface with an ad-atorn at (a)

M-site and (b)H -site... 41

The top views of the Si ad-dimers on the Si(

100

) surface. The

atomic positions are obtained from the T

=0

K, ab-initio molecular

dynamics calculations: (a) A configuration , (b ) В configuration

, (c) C configuration , (d) D configuration... 43

4.4 Total energy change of A- and B- with respect to different cutoff

energies within L D A ... 45

4.5 Total energy difference between A- and B- with respect to different

cutoff energies within LDA ... 46

4.6 Total energy change of A- and B- with respect to different cutoff

energies within P B E ... 47

4.7 Total energy difference between A- and B- with respect to different

cutoff energies within PBE 48

4.1

4.2

4.3

4.8 Velocity change of Si(113) and Si(114) with respect to time 49

4.9 Energy change of system along construction of D dimer... 49

4.10 Schematic diagram for translational motion of A dimer ( top of

the dimer row) at 450 K ... 51

4.11 A 65A X 65A STM image at 360 K showing a dimer adsorbed on

top of a dimer row. The rows are diagonal in the image. Since

dimers readily rotate between configurations A and B at this

temperature, the resulting image is a mixture of the two. Borowsky

et al. [Ref. 19]... 52 4.12 Change of center of mass position of ad-dimer at T=300 K and

T=600 K ... 54

4.13 Change of geometric structure of surface in rotational motion at

300 K. First, ad-dimer is B configuration and after that ad-dimer

rotate to A configuration: (a)0 ps,(b)0.06 ps, (c)0.2 ps, (d)0.25

ps, (e)0.34 ps, (f)0.39 ps... 55

4.14 Schematic representation of atoms and electronic charge density

during the rotational motion of ad-dirner atoms at T —600 K. . . 58

4.15 (a ,b ,c ) Nucleation of epitaxial dimei's perpendicular to the

existing rows when A and D are aligned at T = 0 K. (d ) Schematic

top-view representation of the reconstructed Si(lOO) surface with

new A D A D A D .. dimer row. The top three layers and new row

are shown. The vertical seperation between the up and down atoms,

List of Tables

4.1 Total energies in Hartree and relative energies in m e V o f A , B , C and D type dimers calculated within LDA and G G A by using cutoff energy 12 R y... 42

4.2 Total energies (Hartree) of A- and B- site dimers calculated within

LDA and GGA by using different cutoff energies... 44

4.3 Total Energy of System for a different ad-atom and ad-dimer configurations with |k + Gc|^ = 15Ry. in G G A formalism... 57

Chapter 1

IN TR O D U C TIO N

Since the beginning of the history of mankind, human has always been a part

of the nature. Every time, human who is living with the nature was very

curious about the events his vicinity, and he has always tried to understand these events. Now in all parts of our life we can see the result of this, from space shuttle to microchips; the technology. Improvement of technology can only go further when human continues to explore the nature. Nature comprise the fundamental knowledge so that technology is constructed on this basis like a building. Therefore, fundamental elements of knowledge will be important here. In our world there are a lot of interesting materials which have different and amazing features. With the reason of these properties, these materials can be used in the different fields of technology. They can even open new technological areas, for example, semiconducting elements which are bulding stones of microelectronic technology. They can be activated with electrical and optical effects. With these properties they can be used in various applications in electronics. In nature, there are different kinds of semiconductor materials; one of them, silicon, is important for all technological and scientific research. First we can ask ourselves: Why silicon or where would we be without silicon? This abundant element comprises 27 percent of the earth’s crust and provides the starting ingredients for nearly all the technologies that have transformed our life. In computing and communicating, continuing innovation in smaller and faster circuitry, with

CHAPTER 1. INTRODUCTION

Si as a base material, is taking us to the different though into the information age. Silicon remains the choise for many aplications, despite other materials that promise faster electronics, because it is durable and reliable. Whether it is to solve important problems in science or for the rapidly evolving new world of communications and entertainments, we need faster circuitry, smaller and more perfect chips. With the development of the world, technology must overcome the necessities. Actually, surface investigations starts from this problem. Because, the above properties we talked about, will be related to the surface behaviours of a solid. So that the properties of silicon surface will be fundamental knowledge for the other surface studies. Along with our theoretical studies, we will investigate the events that can occur on surface and homoepitaxial growth of silicon. Even now, the reason of epitaxial growth is a question mark in brains.

Experimentally, by using Molecular Beam Epitaxy^ (M BE) technique the Si(100) surface can be grown on a silicon or on another material as substrate. But althougth this technique is a most powerful technique, during the process there can be imperfections on the grown surface. This will be important for device quality. Reasons of these phehomena relate with the surface structure. Studies in atomic dimensions will reflect these reasons. Growth on Si(lOO) surface occurs in one dimension. And these rows destroy the binding conflguration of surface dimers. This kind of growth occur only on the Si(lOO) surface. Scanning tunneling microscopy experiments will be helpful to reveal a novel form of diffusion by a

simple molecule, the Si dimer, on the Si(

100

) surface. The Si(

100

) surface is the

basis for most lithographic and nanofabrication methods, and has been subject of intensive experimental and theoretical investigations. Generally, this surface

is in the dimerized structures. But there has been considerable controversy

concerning the question of whether the dimers are symmetric or buckled (i.e., the two atoms forming a dimer have different height). The picture of assymmetric dimers (buckled) originally proposed by Chadi.^

In all surface studies thermodynamic equilibrium conditions must be considered because in order to control surface morphologies it is important to achieve an understanding of the thermodynamics and various kinetic processes.

CHAPTER 1. INTRODUCTION

The equilibrium surface morphologies are controlled by thermodynamics, that is, the binding or configurational free energies, while the rates at which dynamic events occur on the surface are controlled by details of the activation barriers, kinetics. With this we can say that all surface events are stochastic. On the surface one event occurs only at a certain probability. On the surface all simple molecules want to choose most stable positions. At these choices the total energy for the system will be minimum. Theoretically, these points could be found by using first principle total energy calculations (total energy in these calculations have been introduced by Car and Parrinello^). Generally all programs achieve to find minimum energy of system in any configurations. To begin, it is necessary to deal with the addition (adsorbsion) of one atom case on the Si(lOO) surface. Because, on the surface, evaluation of simple molecules like dimers depend on the atomic behaviour at any enviromental condition. First studies with one atom are to find stable positions, so that the calculated energy defines the binding geometry of the ad-atom to the substrate. Brocks and Kelly^ have obtained energy surface for the ad-atorn; they have described several saddle points on this energy surface which determine the activation barriers for the diffusion of the ad-atom along

the Si (

100

) surface. In total energy calculations there are several approximation

that can be used. One of the important approximation is the “local density

approximation” (LDA)^ in density functional formalism. Kelly et all’ have used

LDA in their total energy calculations. Hovewer, using the “gradient dependent functional formalism” (GGA)® instead of LDA was more reliable for most results. They have found two adsorption positions on the surface. One of them is the H-point which is between two dimers on the row. The other point is the M- point(absolute minimum) on the trough. These points were stable points for one atom. Namely, when the atom is adsorbed at one of these positions, total energy becomes minimum relative to small displacements of the adsorbate. Kelly

et al. have found that the total energy for the system when the atom is at the H-point was 0.25eV higher than the absolute minimum so that equilibrium

population of this site will not be significant. Yamasaki et alJ have also treated

CHAPTER 1. INTRODUCTION

points to be

0

.

2

eV. Their results are consistent with the result of Kelly’s. After

the binding configuration of Si adatom on Si(lOO) surface is determined the next step is generally to study the diffusion of the ad-atom.

Theoretically Kelly et al} and experimentally Mo et al.^ found that the

diffusion of the Si adatoms is very anisotropic; they claim that it should be directly observable by STM experiments at temperatures somewhat below room temperature. It is of the order of 1000 times faster along the surface dimer rows than across them. In first-principle calculations of the energy barriers that a

Si adatom adsorbed on Si(

100

) at various symmetry sites, has to overcome is

found to be

0.6

eV for the migration along the dimer rows and 1.0 eV across the

rows. The surface diffusion coefficients at temperatures between 348 and 573 K can be measured by STM, which yields the activation energy for diffusion to be 0.67 ± O.OSeV.

When individual atoms are deposited onto the surface at room temperature they diffuse rapidly and readily to find another atom with which to form a stable cluster,i.e., an ad-dirner.^ Dimers are simple molecules on the surface. Ab initio

calculations by both Kelly et and Yamasaki et al7 have emphasized the

significance of dimers, and they claim that the dimers are stable, i.e., a dimer has a lower energy than two seperate monomers. Experimental and theoretical

calculations showed that there are four possible dimer configurations on Si(

100

).

Two of them, A and B , are on the surface dimer rows, the others, C and D , are on the trough (see fig. 1.1) First studies on dimers aimed at finding absolute minimum energy configuration of the dimer on the surface. Two contradictory

results were obtained: Ab initio calculations done by both Kelly and Yamasaki

showed that dimer configurations on the row was more stable then the other

configurations, whereas Zhang et al}^ have emphasized that dimer on the trough

was more stable relative to the dimer on the row. In the latter work, the

interactions between the Si atoms are described by Stillinger-Weber^^ potential.

Experimental studies by both Bedrossian^“* and Mo et al}^ supported their

results; the STM images showed that the residence time for a dimer on the trough was higher than the residence time for dimer on the row. But this observation

CHAPTER 1. INTRODUCTION

does not necessarily imply anything about the absolute minimum, absolute minimum does not necessarily mean a higher activation energy. Although both

Yamasaki et al. and Kelly et al. use ab initio calculations, their results are

different from each other due to the differences in the assumed surface structure.^

Recent ab-initio calculations by Yamasaki et al.^ indicate that configuration A

is the lowest in energy, followed by B , C, and D , with relative energies 0, 80,

180, and 760 meV, respectively. But Kelly et al.^^ proposed that relative energies

which they denoted as B , A , C , D , were 0.0, 10, 310, and 1110 meV, respectively. Experimentally, configuration A is an order of magnitude more populated than

that of B . Experimentally, Swartzentruber et al.^^ estimated the energy difference

to be 59 ± 9 meV.

Diffusion of Si addimers on the Si(lOO) surface have attracted numerous experimental and theoretical investigations because it plays an essential role in the homoépitaxial growth of silicon films. All calculations have suggested that diffusion of ad-dimers on Si(lOO) surface is anisotropic: They prefer diffusion along the top of the dimer row.^® Theoretically diffusion of dimers have been

investigated by using two step optimization in Ref. 7. They calculated the

energy barriers of two diffusion processes along and across the substrate dimer rows, and found that along the row the energy barrier is 1.44 eV, but by using atom exchange mechanism they observed the activation barrier decreases to the 1.28 eV for the dimer diffusion along the row. However, this value is still

higher than the experimental value. According to Borowsky et al.^^ diffusion on

the row is activated just slightly above the room temperature. Experimentally Swartzentruber have found that the diffusion barrier for silicon ad-dimers on the row is 0.94 ± 0.09 eV at 400K.^° This activation barrier for the diffusion along the top of the substrate rows is higher than the measured barrier for the Si monomer diffusion, 0.7 eV. On the other hand, diffusion barrier across the

dimer row is 1.8 eV.^ Ho et al.^^ studied the diffusion from the trough to the

top of the dimer row using the tight-binding molecular dynamics and ab initio

total energy calculations. The relative energies of the four principle ad-dimer

CHAPTER 1. INTRODUCTION

meV, respectively, for the A , B, C , D addimer configurations. The energy

barrier for the diffusion of the ad-dimer from C to B along this path is 1.72 eV.

Borowsky et al}^ observed that at 450 K dimers are able to diffuse across the

dimer rows of the substrate. At this temperature, average hopping rate for the

dimer was 0.4 ±

0.2

Hz. Experimentally, using the average residence time in the

trough, the activation energy is found to be 1.36 ± 0 .0 6 eV for leaving the trough at 450 K. Third possibility for dimer diffusion path is diffusion along trough. Goringe and Bowler^^ propose that, as for row diffusion, trough diffusion requires the dimer to move one atom at a time. This mode of diffusion gives a barrier of 1.28 eV. Atom tracking studies*^ of trough diffusion confirm the presence of long jumps in this diffusion channels. The dimer along the trough at a much slower rate of 6.5 x 10“ ^ Hz at 433 K. Experimentally, the energy barrier for trough diffusion was found to be 1.27 ± 0.08 eV.

The other experimental analysis was rotation of ad-dimers on the row between

A and B configurations. The possibility of transition either A —> B or B —> A

through a rotation was pointed out based on the arguments of energetics and

energy barriers between these two configurations.^^ Experim entally, it was

shown that the rotational transition occurs on the time scale of seconds at

room temperature. Swartzentruber et al}^ obtained qualitative data on the

dynamics of rotation by performing atom tracking with STM. Borowsky e t

observed dimer diffusion on the row with rotation of dimer between A and B

configurations at 360 K. Yamasaki had found energy difference between A and

B configuration was very low. So that dimer can make flip flop m otion between

A and B configurations. The measured rotational barrier for the A —> B (B —+

A) configuration is 0.7 ± 0.08eY(0.64 ± O.OOeY) with a frequency prefactor of

1010.8± 1.4 sec- 1

(10

10.7± 1.3sec ^).

Epitaxial growth of Si(lOO) surface has been challenging for most theoretical and experimental studies. According to the experimental results, epitaxial growth

could not be observed below 500K. Pearson et al?"^ observed island growing

(ID growth) with an STM on the Si(lOO) surface at 536 K. They observed ID row formed perpendicularly to surface dimer row. Bedrossian^'* investigated

CHAPTER 1. INTRODUCTION

the nucléation mechanisms of silicon ad-dimers at 400 K. Actually dimers can nucleate at this temperature. But this row will not be epitaxial. This study showed metastable lines of dimers arranged end to end. And in this 1 ine dimers sitted on the trough as a D configuration. As we stated before, the Si(lOO) surface reconstructs to form rows of dimerized atoms as a ID. And these rows will be perpendicular to the substrate dimer rows.

All these Si(lOO) surface studies explained the diffusion, adsorbsion, exchange and growth mechanisms. By considering all these facts in this thesis we will try to explain these events within time scale of picosecond. Rotation of dimer on the row will b e investigated. Diffusion event for a dimer on the row will be explained. Last thing to be delth with will be the epitaxial growth. How dimers

can nucleate as a ID row, and what is the mechanism for periodical A and

D site rows? Is there any interaction between A and D dimers? How does

this interaction effect the growth? We will try to answer these questions. In

our studies we used ab-initio molecular dynamics technique to be able to find

energy and binding configurations of atoms or dimers on the Si (

100

) surface by

CHAPTER 1. INTRODUCTION

Figure 1.1: Positions of adatoms on the Si(001)-(2 x

1

) surface are designed by

according to the labeling of Ref. 5. Adatom-adatom combinations forming A , B , C , D type dimers are indicated by arrows. Large, me diurn and small filled, empty circles denote, respectively, the adatom (or ad-dimer), raised and lowered surface dimer atoms and subsurface atoms.

Chapter 2

A TO M IC STRUCTURE OF

Si(lOO) SURFACE

2.1

Surface Structure

The Si(100) surface has been the subject of a number of experimental and theoretical studies. The microscopic structure of the Si(lOO) surface was first investigated 30 years ago in a low-energy electron diffraction study (LEED) by

Schlier and F a rn sw orth .E v id en ce of a 2 x

1

structure (Fig. 2.1) was observed

and a reconstruction based on the formation of dimers was proposed. Subsequent investigations showed higher order diffraction spots, intensity and sharpness that was highly dependent on sample treatment.

Theoretical modeling^*"’^® together with experimental observations eventually

established the surface dimer for an important feature of the Si(

100

) surface,

but the origin of higher order diffraction spots remained controversial. It was recognized that buckled dimers could account for higher-order diffraction. In empiracal tight binding calculations, Chadi^ represented that the buckled dimer

structure and the most favorable geometry for the Si(lOO) surface is c(4 x

2

) with

respect to the symmetric (

2

x

1

) surface structure. It was shown that the buckling

involves charge transfer between two atoms forming a dimer. Althought every silicon atom makes tetrahedral bond, one of the buckled dimer atoms makes

CHAPTER 2. ATOMIC STRUCTURE OE SI(IOO) SUREACE

10

hybridization and the other makes sp^ hybridization. Later, all STM studies and

pseudopotential calculations showed that dimers can be found either

(2

x

1

), or

c(4 X 2) structures depending on the temperature. Wolkow^^ has shown, however,

that the number of buckled dimers increases as the temperature is reduced and that for T'=120 K the c(2 x 4) reconstruction occurs on most of the surface.

« ® m d--- € • • s /1 p --- ¥ • • m 7... ^ { } K) • m • 7 ^ • i 7 \ ) • • y.7 K © --- ( 7 K ^--- i J K ) --- i 7 X/ )--- © IK 1x2) p(2x2> i 1 7 ---• a---• • \ / 7 V7---^ \ / 1 i 7 t • • 7 \ • • 2 7 t )---1 7---V Ik---6 c(2x4)

Figure 2.1: Si(lOO) surface structures with symmetric dimers[p(l x

2

)] and

asymmetric (buckled) dimers[p(2 x 2) and c(2 x 4)]

Pseudopotential density-functional calculations have examined the energetics of dimer buckling[27-30]. In general, these calculations have found the buckled and symmetric dimer surfaces to be close in energy: the difference in energy is

typically H.leV! dimer or less. Dabrowski and Scheffler,^^ using well-converged

plane wave basis set, found the

2

x

1

buckled dimer surface to be lower in

energy then the symmetric dimer surface by O.leV/dimer. Roberts and Needs^°

found that the p{2 x

2

) buckled dimer surface is lower in energy than the

symmetric dimer surface by D Æ eVjdimer. These results suggested that within

CHAPTER 2. ATOMIC STRUCTURE

OF

SI(IOO) SUREACE

10

hybridization and the other makes sp^ hybridization. Later, all STM studies and

pseudopotential calculations showed that dimers can be found either

(2

x

1

), or

c(4 X 2) structures depending on the temperature. Wolkow^^ has shown, however,

that the number of buckled dimers increases as the temperature is reduced and that for T=120 K the c(2 x 4) reconstruction occurs on most of the surface.

w

a----i• ) ( p---- 9 • • / / / / /a f k9 • • • • \ e 9 V • m 7 \9 m • ©----i 9 ^ )----ii----( 7 \ 9 i----© ,HU2) p(2x2) M i1 ?---^ • • 'I { d----o • • ^ i \ y/ / O '/Y i¥ 7 V7 \9 ¥ 1 • • s / 7 t / • • M f t ^----i 7 t ) J----© c(2x4)

Figure

2

.

1

: Si(lOO) surface structures with symmetric dirners[p(l x

2

)] and

asymmetric (buckled) dimers[p(2 x 2) and c(2 x 4)]

Pseudopotential density-functional calculations have examined the energetics of dimer buckling[27-30]. In general, these calculations have found the buckled and symmetric dimer surfaces to be close in energy: the difference in energy is

typically DAeV!dimer or less. Dabrowski and Scheffler,^^ using well-converged

plane wave basis set, found the

2

x

1

buckled dimer surface to be lower in

energy then the symmetric dimer surface by O.leVfdimer. Roberts and Needs^*^

found that the p(2 x

2

) buckled dimer surface is lower in energy than the

symmetric dimer surface by 0.09eVfdimer. These results suggested that within

CHAPTER 2. ATOMIC STRUCTURE OE SI(IOO) SUREACE

11

favorable. The total-energy calculations for the c(2 x

4

) and p(2 x

2

) surfaces

supported this a s s e rtio n .A cco rd in g to these calculations, the computed ground

state of the flat Si(OOl) surface was the c

(2

x 4) reconstruction, which induces

an energy gain of 0.2AeV per surface dimer, and a buckling of 0.72Â. The

p(2 X 2) reconstruction was nearly degenerate energetically with the c(2 x 4),

so that its energy was higher by only SrneV per surface d i m e r . I n this case

the buckling is 0.7lA. Northrup^^ found that also excited-surface energy for the

c(2 X 4) structure is significantly lower than that of the symmetric dimer surface at

room temperature. Surface state dispersions and band gaps are calculated. The results support the existence of correlated dimer buckling at room temperature. The equilibrium c(2 x 4) surface obtained from total-energy calculations is 0.14 eV /dim er lower in energy than the 2 x 1 symmetric dimer surface, exhibits a

dimer buckling of 0.69A, and has a surface energy of

1

.

39

e F /( l x

1

).

2.2

Surface Properties: Atomic Motion

A long-range van der Waals type interaction is observed to grow as the molecule approaches the surface but the attractive interaction gives way to a repulsive one cis the distance between the two decreases further. However, because of some chemical interactions, such as the filling of the molecule’s antibonding orbitals by the electrons at the surface or the dumping of its electrons into the empty surface bonding levels (chemical adsorbtion), turns on at shorter distances, the repulsion does not continue to increase. The interaction first becomes less repulsive and then turns strongly attractive, giving rise to a potential maximum. An adsorbed molecule (physically adsorbed) is bound to the surface via a rather weak van der Waals type bond. This bond involves no charge transfer from the substrate to the adatom or vice versa. Rather, the attractive force is provided by the instantaneous dipole moments of the adatom and its nearest neighbour surface atoms. The interaction can be described by the potential energy diagram shown

in fig. 2.2 An incoming molecule with kinetic energy Ek has to lose at least this

CHAPTER 2. ATOMIC STRUCTURE OE SI(IOO) SUREACE

12

lattice phonons in the substrate, and the molecule is then equilibrated in a state of oscillation in the potential well of depth equal to the binding or adsorption

energy Ea. V r ^Repulsive t:. Atiraciivc V

Figure 2.2: Potential energy of an adatom in a physisorbed state on a planar surface as a function of its distance z from the surface.

When atoms or dimers come onto Si(100) surface, their motion can reveal interesting results. Actually atomic movement strongly depends on temperature and surface potential. On the surface all ad-atoms or dimers continuously migrate until they constitute a stable structure on the surface. Here the important thing is the behaviour of atoms until stable configuration. Direction of movement for atoms or dimers is arbitrary (random) under the effect of temperature. But, generally atoms or dinlers choose minimum energy barrier to pass. Passage of atoms throught the barrier is dictated by the diffusion rate. Diffusion rate can

CHAPTER 2. ATOMIC STRUCTURE OF SI(IOO) SURFACE

13

be written as,

1

, Ea.

r = - = ^ „ x e x p ( - — ).

(2.1)

In this equation Uq is frequency prefactor which only depends on surface structure.

And Ea is activation energy. Diffusion rate for atoms or simple molecules is

inversely proportional to mean residence time r. By experimental determination

of residence time for atom, activation energy can be found. For example,

according to Borowsky’s^^ results mean residence time for silicon ad-dimer which is on the row is 675 ±151 sec. The mean residence time on the trough was 116 ±2 0 sec. These results show us that the dimer spends more of its time on the row.

We can see that diffusion rate for a silicon atom will increase with te m p e ra tu re ,b e ca u se the probability for diffusion over a barrier is related to temperature, fig. 2.3. So that according to this figure the mean residence time decreases from 13.1 to 0.11 s as temperature is increased.

And also one could find activation energy and frequency prefactor by using

the data in fig. 2.3. The activation barrier Ea is measured as the slope of the

Arrhenius plot of the logarithm of the transition rate versus 1/kT (Fig. 2.4).

2.3

Thermodynamic Quantities About Surface

In general, adsorption includes the binding to the surface of an atom or molecule arriving either from the vapour outside the solid or from interior by diffusion process. There are three macroscopic properties of a surface of a solid containing

more than one chemical element. These are the surface tension, the surface

stress, and the specific surface free energy. The surface tension

7

is defined as the reversible work done in creating a unit area of new surface at constant

temperature, volume, and total number of molecules. Under conditions of

constant temperature and volume the work done in creating new unit area of surface is the H e lm h o ltz free e n e rg y for that area and this is defined as being the specific surface free energy. The surface stress, on the other hand, is a more complex quantity because it involves the work done in deforming a surface.

CHAPTER 2. ATOMIC STRUCTURE OF SI(IOO) SURFACE

14

0^

0

-5

0 -15

CL

< 2 0 ^ ';r '

10

-O 0 -

1

-

10

.

£ -

20-1

J j —

1

^.--- - p

|1

I

1'30

J * j 'l i 3 » C | ¡

1

;) / ' H < \ t

200

'ı;*^ı

30D

4CiO

rl^,

r 5fli0

1

. ,‘4|i7·'

w

' '*1

0

uc .■ ^

0

, , ' 5

tome (sec)

20

Figure 2.3: Time-line plots. At the top are shown two plots of the ad-dirner

position along the row as a function of time acquired at 65 and 128 C° [Ref. 20].

On a surface, free energy can change from one site to another. For this reason occupation probability will differ. According to fig. 2.5, atoms generally choose sites which has a low free e n e r g y . T h i s is quite ordinary result. Because every canonical system wants to make its free energy minimum. So the total energy is minimum. Result of this is that the atoms or dimers want to settle for the most stable configuration of surface.

Temperature is another powerful factor controlling kinetic processes and growth on a surface. Actually not every surface show ID growth like Si(lOO) surface. They can display multilayer growth depending on temperature. For a

growth, critical island size Rc is defined at which a second layer nucletites on top

CHAPTER 2. ATOMIC STRUCTURE OE SI(IOO) SUREACE

15

Temperature (K)

400 380

360

340

320

300

Figure 2.4; Arrhenius plot of the logarithm of the diffusion rate versus

1

/k T from

which the activation energy barrier is extracted [Ref. 20].

will nucleate a second layer before coalescence, giving rise to multilayer growth. But with an increasing temperature, this problem can be surmounted. Because with an increasing temperature, atoms can diffuse and hop from the diffusion barrier at the island edge without the other atoms nucleating on top of the islands. As a result of this, layer-by-layer growth can be realized. Generally, with temperature effect surfactants are used to make layer-by-layer growth. Crucial role of a surfactant is to reduce the spacing between islands, so that coalescence occurs before the islands ever reach the critical radius.

All the above parameters will effect the growth process. But now we want to discuss the initial stages of ID growth on the Si(lOO) surface. Actually, there are a lot of ideas about this subject. The explanation of this curious shape anisotropy

has been controversial. On the S'i(lOO) — (

2

x

1

) surface, the long axis of growth

is perpendicular to the substrate dimer rows. It is puzzling that the islands grow perpendicular to the direction of fast diffusion. Because the mobility is high along the rows, one might expect that islands would capture atoms from a greater

CHAPTER 2. ATOMIC STRUCTURE OF SI(100) S UREACE

16

S '

0.14 i

rzz i

I

_{C 0.10}

^

3.08 H

3.C6

4

5 5

7

H S 10 11

SlteNuniier

QOO IS jX

Figure

2

.

5

: Atom-track data taken at 128 C°. Site occupation probabilities and

the extracted relative free energies of the eleven sites are measured from six data sets, total time is 196 s. Error bars reflect the statistics related to the number of

hops at each site [ref.

20

].

distance in this direction. This would lead to a growth predominantly along the

rows.^^ Mo^® et al. proposed a kinetic model based on sticking anisotropy to

explain the island shape anisotropy. If island ends are stickier than island sides, arriving ad-atorns stick preferentially to the ends, resulting in the growth of many

ID islands. Metiu et al}^ proposed an alternative model based on an exchange

mechanism. Metiu suggested that an adatom arriving on the side of a Si island may displace an existing island atom to the top of the island. The displaced atom diffuses rapidly along the top of the island (the direction of fast diffusion) until reaching an end, where it can fall over the edge cind stick.

Chapter 3

MODELLING THE

IN TER ACTIO N S

3.1

Introduction

In order to carry out a theoretical study of complex materials, it is important to model the interactions between the ions in these structures as accurately as possible, but still keeping the calculations computationally feasible. Actually,

there are several methods to model these interactions. One of them is the

molecular dynamics technique. With the more powerful computers, the molecular dynamics technique has become an effective tool in the physics of condensed matter systems. In the molecular dynamics method, the forces acting on particles in a cell are found and the classical Newtonian equations of motion are solved numerically. The largest part of a molecular dynamics simulation is the evaluation of the forces which are required to be known in order to find the relaxed ionic positions.

To model interactions in which no a priori bonding information is known,

an approach at a more fundamental level is required. One must turn to the formidable task of solving the Schrôdinger equation for the electrons (in fact it is the Kohn-Sham equations that are solved where the many electron interactions are approximated by a local potential). This method will be described in later

CHAPTER 3. MODELLING THE INTERACTIONS

18

sections where the bonding topologies of the structures under consideration are unknown. There are obviously great advantages in using a molecular dynamics method where the only specification of the atomic numbers of the ions present

(so called ab initio methods) are required. Before an evaluation of the forces

on the ions can be performed, a massive minimization calculation is required in the extremely large phase space of the basis set of the electronic wavefunctions. So that speed of calculations strongly depends on the number of atom in the system under study. Actually, for large number of atoms, for example, at the order 10®, this method will be useless. Instead of this, generally, an empirical potential method is used. The empirical potentials allow this, because the forces are calculated by the evaluation of a simple function.

3.2

Molecular Dynamics Method

Whether the interatomic forces are calculated by an empirical potential or ab-

initio techniques, the methodology for moving the ions in a molecular dynamics simulation is the same. Each method produce the forces, F¿ exerted on particle i by the rest of the system. Throughout this work, periodic boundary conditions

are used on a cell containing N ions. Some initial starting configuration R¿ is

assumed and F¿ is calculated for each i using ab-initio technique. The initial

positions R and velocities R are given and equations of motion

Fv = niiRi (3.1)

are then integrated numerically (for example by using Verlet algorithm) over a

small time step At. The size of A t for atomic systems is of the order of sec =

1

fs . This allows the positions to be updated and the process is repeated. One of

the more intuitive method of finding the minimum energy structure is to extract the kinetic energy (thus eventually finding a OK configuration) by quenching the system. That is, removing kinetic energy by stopping the ions when the force exerted on them is in the opposite direction to their velocities. This method will continue until the forces are sufficiently small, where the ions are in a local

CHAPTER 3. MODELLING THE INTERACTIONS

19

minimum of energy. Introducing temperature into the system is also desirable. An intuitive method for stabilizing the temperature is to scale the velocities of the particles at every timestep so that total kinetic energy of the system is given by f h e T per degree of freedom. Although this ensures that the temperature is constant at every time step by adding or removing energy from the system, it is not physically correct. Temperature is a statistical quantity here and therefore an average over many time steps is required to define the temperature. A method

which allows the instantaneous temperature, T, to fluctuate, but, where < T >

remains constant was suggested by Nose.^^ This couples the temperature to an external heat reservoir, allowing the temperature to become a dynamical variable of the system. It is analogous to attaching a ’spring’ to the temperature so that

T can fluctuate around a mean value. To use Nose thermostat, the Newtonian

equations of motion have to be integrated at each timestep requiring knowledge of the instantaneous force on each particle.

3.3

First Principles Molecular Dynamics

3.3.1 Approximations

In solid state physics to be able to carry out a theoretical study of complex materials, it is necessary to model the interactions between all ions and electrons.

But this is not easy. In order to model interactions in which no a priori bonding

information is known, an approach at a more fundamental level is required. Simplifications are introduced when we deal with electron-electron and electron-

ion interactions. So that density-functional theory (Hohenberg and Kohn, 1964;

Kohn and Sham, 1965) is used to model the electron-electron interactions,

pseudopotential theory (Philips, 1958; Heine and Cohen, 1970) to model the

electron-ion interactions, supercells to model systems with a periodic geometry,

and iterative minimization techniques (Car and Parrinello, 1985; Payne et ai,

1986; Williams and Soler, 1987; Gilan, 1989; Stich et ai, 1989; Teter et al, 1989)

CHAPTER 3. MODELLING THE INTERACTIONS

20

3.3.2 Density Functional Theory

With the exception of the free or nearly free electron gas, it is standart practise in solid state physics to reformulate the many-electron problem in terms of one- electron Hamiltonians. Practically in most cases, there are terms which depend on the ground state one-electron density matrix, so that the problem must be solved by a self-consistent procedure. The Kohn and Sham one-electron Hamiltonian is most widely adopted: it includes, along with the kinetic and Hartree terms,

an effective exchange and correlation potential, ExcV^{r)\, which is in turn a

functional of the electronic density n(r) which is given as;

n (-■) =

2

E I * ( (3.2)

The total energy can be written as a sum of several terms;

E[n{r)] = j Vion{r)n{r)dr + ^ j + G[n{r)] (3.3)

is a unit functional of electron density, for a fixed set of atomic nuclei at Ri.

First two terms are the classical Coulomb interaction between electron-ion and electron-electron. G[n(r)] include the terms of kinetic energy of free electrons and exchange and correlation energy T’a,-c[?^(r)].

where

G[n{r)\ = T [n(r)] -h E^c[n{r)]

T[n(r)] = V’* ( ^ ) ( ^ ) V V i ( r · )

(3.4)

(3.5) According to the Hohenberg-Kohn theorem, the total energy function given

by eq.

3.3

is stationary with respect to variations in the ground state charge

density, that is, it is subject to condition

^¿T[n{r)] , , f

=

0

(3.6)

where iXidr) is the functional derivative of the exchange-correlation energy with

CHAPTER 3. MODELLING THE INTERACTIONS

21

variation in the charge density leaves the particle number

N =

J

n{r)dr (3.7)

unchanged. This can be ensured by the condition

J

6n{r)dr =

0

. (

3

.

8

)

Applying the condition of constant particle number to equation gives the result

SE{r) _ 8T[n{r)\ , , [ n (r') + V^°^{r)+ f - I ^ d r +

J

r

—

r

(3.9)

Sn{r) Sn(r)

is a constant, ¡y, which is the Lagrange multiplier associated with the requirement

of constant particle number. Compairing this to the corresponding equation for a

system with an effective potential V(r) but without electron-electron interaction

results in

^F!(r\ f)T\n(r\]

(3.10)

, = ^ + K (r).

6n{r) Sn(r)

it can be seen that the mathematical representations are similar provided that

l

(

r

')

V{r) = V U r ) + f T ^ ^ . d r + t i W ) ·

J

r

—

r\ (3.11)

The effect of this is to allow an indirect variation in n(r) through variation in

the Khom-Sham single particle orbitals, It then follows that the solution can

be found by solving the Schrodinger equation for noninteracting particles moving under the influence of an effective potential T(r)

+

r

(

r

) } * (

r

) =

c

,.

V

-

i

(

r

)

(3.12)

which gives the charge density n(r). The solution of the system of eq. 3.12 leads to the energy and electronic charge density of ground state and all quantities which can be derived from them. The minimum of the Khom-Sham energy functional leads to the ground state charge density of the electronic system where ions are located at the fixed positions {i? ,}. It is only this minimum which has any physical meaning.

CHAPTER 3. MODELLING THE INTERACTIONS

22

3.3.3

Exchange-Correlation Energy and LDA and GGA

Formalisms

Actually, the Coulombic interaction between electrons is dealt with Hartree

approximation. Because of electron-electron interactions, wavefunctions for

electron system will change to But we assume that ^ and free electron

wavefuntion $ are not too dissimilar, so we can calculate approximate interaction

energy, £/i/, by forming < >

Ehf = < ^\Ho + Vh\^ >

where is free Hamiltonian and Vh is interaction potential. We know that

total wavefunction of electron system is antisymmetric, so antisymmetricity will be the reason of spatial seperation between two electrons with the same spin. When there is exchange between these electrons, the Coulombic energy of the electronic system will reduce. The reduction in the energy of the electronic system due to the antisymmetry of the wavefunction is called the exchange energy. For electrons with different spins, spin part of the wavefunctions of these electrons will be antisymmetric but spatial component will be symmetric. According to this result electrons can be found at the same space. But with an effect of repulsion between elctrons, this will be impossible. So that for different spin electrons there

must be correlation between them. Difference between exact energy and Shf is

the correlation energy. The important point that makes this system easier to solve than, for example the Hartree-Fock equations, is that the effective potential is local. Therefore there is no more complexity added in solving Eq. 3.12 than there is in the Hartree approximation. Of course, this is only true if the exchange- correlation energy can be described as a function of the local charge density. A

method of doing so is known as the local density approximation{LDA) defined by

Ceperley and Alder.® In LDA, the exchange-correlation energy per electron at a

point r in the electron gas, exc{f) is equal to the exchange-correlation energy per

CHAPTER 3. MODELLING THE INTERACTIONS

23

point r. It follows that

so that with Exc[n{r)] = j exc[n{r)]n{r)dr _ ^{»(r)ga;c[n(r)]} 8n{r) exc(n(r)l = 4 r [ > ‘ ('·)] (3.13) (3.14) (3.15) where Eq. 3.15 is the assumption that the exchange-correlation energy is purely

local. The other approximation for exchange-correlation is the generalized

gradient approximations (G G A ’s). We know that only the exchange-correlation energy as a functional of the electron spin densities n|(r) and n|(r) must be approximated. The exchange energy for a system of electrons with density n(r)

IS

(3.16) where n(r) is electron density at point r and n i(r, r - f i?) is the density at position r + i? of the exchange hole about an electron at position r. The basis for this approximation is the gradient expansion of the exchange hole, with real-space cutoffs chosen to guarantee that the hole is negative everywhere and represents a deficit of one electron. The GGA clearly takes into account the fact, neglected in the LDA, that the exchange hole in the inhomogeneous system is “off center” with respect to the electron it surrounds. Now above equations become

where s = \An\/2kfTi, kj = and

= (1 -1- 1.296

s

2 -h Ids'* +

3.3.4 Bloch’s Theorem and k-point sampling

(3.17)

(3.18)

It is difficult to handle the infinite number of interacting electrons moving in the static field of an infinite number of ions. Wavefunction for each of the infinite

CHAPTER 3. MODELLING THE INTERACTIONS

24

number of electrons which will extend over the whole space of the solid and the basis set in which the wavefunction wilt be expressed will be infinite. The ions of the perfect crystal have translational symmetry at OK. So, the external potential

is expressed as K (r) = V{r + 1), with I being the translational vector. This is

the requirement needed for the use of Bloch’s theorem. By use of this theorem, it is possible to express the wavefunctions at reciprocal space vectors of a Bravais lattice. With Bloch’s theorem, the wavefunctions can be written as the product of a cell-periodic part and a phase part;

iii(F) = eip(!k.r)./j(r)

(3,19)

Cell periodic part of the wavefunction/·(?) can be expressed by expanding it into a finite number of plane waves whose wavevectors are reciprocal lattice vectors of the crystal.

/.(? ) = (3-20)

G

where G are reciprocal lattice vectors which are defined by G · 1 = 2Trm. The first Brillouin zone can be mapped out by a continious set of points, throughout

that region of reciprocal space (k — space). The occupied states at each k — point

contribute to the electronic potential of the bulk solid. Since the set [k] is dense,

there are an infinite number of k — points in the Brillouin zone at which the

wavefunctions must be calculated. Therefore, if a continium of plane wave basis sets are required, the basis set for any calculation would still be infinite, no matter how small the plane wave energy cut-off is chosen. For this reason electronic states

are only calculated at a set of k —points determined accordind to the shape of the

Brillouin zone. It is therefore possible to represent the electronic wavefunctions

over a region of reciprocal space at a single k — point.

3.3.5 Plane-wave represantation of Kohn-Sham equations

According to the Bloch theorem plane waves are used as a discrete basis set for the electronic wavefunctions. In principle there must be infinite Fourier series. The coefficients of the plane waves, Cj-,k+G) correspond to a kinetic energy of

CHAPTER 3. MODELLING THE INTERACTIONS

25

/2m)\L· + G p . Plane waves with smaller kinetic energy typically have more important role than those with very high kinetic energy. The introduction of a plane wave energy cut-off reduces the basis set to a finite size. This kinetic energy cut-off will lead to an error in the total energy of the system. Nevertheless, it is possible to make this error arbitrarily small by increasing the size of the basis set by allowing a larger energy cut-off.

Another advantage of expanding the electronic wavefunctions in terms of a basis set of plane waves is that Kohn-Sham equations take a particularly simple form. Substitution of Bloch equation into the Kohn-Sham equations gives

y ^ { - ---- | k -| -G | Î G G '+ T i o n ( G — G ' ) - f T e ;e c (G — G ' ) - f T j .c ( G — G ') } c i , k + G ' = £ iC i,k + G >

^ Am

(3.21) It can be seen that in this form the reciprocal space represantation of the kinetic energy is diagonal and the various potentials can be described in terms of their Fourier components. Usual methods of solving the plane wave expansion of the Kohn-Sham equations is the diagonalization of the Hamiltonian matrix. It follows that the size of the Hamiltonian matrix is determined by the energy cut-off

(3.22)

where |Gc| is defined as a magnitude of the maximum value of reciprocal lattice vector. This matrix can be diagonalized by conventional matrix diagonalization techniques, but a more computationally efficient method exists where the plane wave coefficients are treated as dynamical variables.

3.3.6 Electron-ion Interactions

A-The Pseudopotential Approximation

We mentioned earlier that the electronic wave functions can be expanded in terms of plane waves. For a core electrons wavefunction must be expanded interms of very large number of plane waves. For this reason an extremely large plane- wave basis set would be required to perform an all-electron calculation, and a

CHAPTER 3. MODELLING THE INTERACTIONS

26

vast amount of computational time would be required to calculate the electronic wavefunctions. The pseudopotential allows the electronic wavefunctions to be expanded using a much smaller number of plane wave basis states. The most physical properties of solids depend on the valence electrons to a much greater extent than on the core electrons. The pseudopotential approximation exploits this by removing the core electrons and by replacing them and the strong ionic potential by a weaker pseudopotential that acts on a set of pseudo wave functions rather than the true valence wavefunctions

The pseudopotential is also constructed such that the scattering properties of the pseudowavefunctions are identical to the scattering properties of the ion and core electrons. In general, this will be different for each angular momentum component of the valence wavefunction, therefore the pseudopotential will be angular momentum dependent. Pseudopotentials with an angular momentum dependence are called non-local pseudopotentials. If they are independent from angular momentum it is called as a local pseudopotential. A local pseudopotential is a function only of the distance from the nucleus.

B- Kleinman-Bylander Pseudopotential

The most general form of a non-local pseudopotential which is improved by Kleinman and Bylander^® is,

Kv _{E l i " ' ”}

Im

> V i < Y iIm (3.23)

where Yim are spherical harmonics and Vi is the Ith angular momentum

component of the pseudopotential acting on the wavefunction. If there are Npw

plane waves in the expantion of the wavefunction at each Â;-point and there are Nk

À;-points then the evaluation of Vion will require NpiuNk{Npw -

1

- l

)/2

projectors

of the above form to be calculated for each angular momentum component /.

This can be seen as follows: the crystal potential Vcr{r) is obtained by placing a

CHAPTER 3. MODELLING THE INTERACTIONS

27

^scudü \ V"' t s 7"' '1^ V ' _>cuib f / \ y 7Jt

Figure 3.1: An illustration of the full all-electronic (AE) wavefunction and electronic potential (solid lines) plotted against distance, r, from the atomic

nucleus. The corresponding pseudowavefunction and pseudopotential are

plotted(dashed line). Outside a given radius, Tc, the all electron and pseudo

electron values are identical.

incorporates the crystal symmetry, hence

V „{G - G') = J2Ss{G - G') V,,(G - G') (3.24)

where the summation index is over ionic species and the structure factor for each

species is Ss{G — G') = Y^iexp[i{G — G ').7i), where T is lattice translational

vector. The total ion-electron energy is then

= Y ,< i> \ Y h .> V „{G - f f ) < > (3.25)