5.2 Atomistic models and rate equations	149

	60		±1
				3D: = 0
)	40
j
G(
energychange,				2D:	9= 0
	20
		1
Free
	0
	±20	2
	0	20	40	60	80
			Cluster size, j

Figure 5.4. Free energy of nucleation DG(j) for 3D (full line) and 2D (dashed line) clusters. These curves are to scale for the surface free energy term X54, and for Dm or Dm9521, 0, 1 and 2. DG(j), X and Dm are all in the same (arbitrary) units. However, if the unit is taken as kT, then the scale of free energy is the same as used by Weeks & Gilmer (1979) as shown in ®gure 1.12 and discussed in section 1.3 (replotted after Venables et al. 1984).

energy in an atomistic calculation is very similar to ®gure 5.4, but has a discrete character which can show secondary minima at particularly stable sizes, which are sometimes referred to as magic clusters. An early example of such a calculation is given by Frankl & Venables (1970). However, an atomistic model should be consistent with the macroscopic thermodynamic viewpoint in the large-i limit. To ensure this is not trivial, and most models don't even try; if I rather emphasize this, it is because I am attempting to sort this out in my research papers (Venables 1987). In other words, there are (at least) two traditions in the literature; it would be nice to unify them convincingly.

5.2Atomistic models and rate equations

5.2.1Rate equations, controlling energies, and simulations

We have considered simple rate equations for adatom concentrations in section 1.3, and in problem 1.2, adding a diVusion gradient in problem 1.3. Now we need to add nonlinear terms to describe clustering and nucleation of 2D or 3D islands. These equations are governed primarily by energies, which appear in exponentials, and also by fre-

quency and entropic preexponential factors.

150 5 Surface processes in epitaxial growth

The most important energies are indicated in ®gure 5.3(c): Ea and Ed control desorption and diVusion of adatoms. These processes are linear, and have been discussed in chapter 1; Ej and Ei are binding energies which control clustering ± these processes are non-linear. In the simplest three-parameter model, we can build the cluster energies out of pair bonds of strength Eb. Without this simplifying assumption, we can't make explicit predictions; but with it, we can develop models which describe nucleation and growth process over a large range of time and length scales. This is the main advantage of such `mean ®eld' models (Venables 1973, 1987). They are known not to describe ¯uctuations very well, so various quantities, such as size distributions of clusters, are not described accurately. In current research, using fast computational techniques such as `Kinetic Monte Carlo' (KMC), the early stages can be simulated on moderate size lattices. These KMC simulation results using the same assumptions can then be used to check whether mean ®eld treatments work for a particular quantity.

The emergence of computer simulation as a third way between experiment and theory is clearly a growth area of our time. To make progress in this area, one has to start with the simplest models, and stick with them until they are really understood. You need to beware generating more heat than light, and in particular of generating special cases which may or may not be of real interest. As we will see later in this chapter, the number of important parameters can become alarmingly large. Simulations can however be very illuminating, and may suggest inputs for simple models that one hadn't thought of. Animations are immediately appealing, and if Spielberg can do it, why shouldn't we? The problem lies only in the subsequent claims for correspondence with reality; then a measure of self-discipline is needed, both from the author and the reader. An extensive list of methods and suitable warnings are provided by Stoltze (1997).

5.2.2Elements of rate equation models

To develop an atomistic model, we consider rate equations for the various sized clusters and then try to simplify them. If only isolated adatoms are mobile on the surface, we have

dn1/dt5R (or F) ± n1/ta 2 2U1 2 oUj,	(5.5)
and for larger clusters
dnj /dt5Uj21 2 Uj ( j$ 2),	(5.6)

where Uj is the net rate of capture of adatoms by j-clusters. This is not very useful yet, since we need expressions for the Uj, and we need the simpli®cation introduced by the idea of a critical nucleus size. In its simplest form, this means that (a) we can consider all clusters of size .i to be `stable', in that another adatom usually arrives before the clusters (on average) decay; the reverse is true for clusters of below critical size; and (b) these subcritical clusters are in local equilibrium with the adatom population.

The ®rst consideration leads to de®ning the stable cluster density, or nucleation density nx, via the nucleation rate

5.2 Atomistic models and rate equations	151
dnx/dt5oj $i (Uj 2 Uj 11)5Ui 2 . . .,	(5.7)

since all the other terms cancel out in pairs. The . . . means that we can add other terms, such as the loss of clusters due to coalescence, in various approximations. The second consideration leads to arguments about detailed balance, and the Walton relation, named after a paper where local equilibrium was ®rst discussed in this context (Walton 1962). These detailed balance arguments lead to all the Uj, for j,i, being zero separately, and hence dnj /dt50. Note that this is not the same as a steady state argument where dnj /dt50 because Uj215Uj; it is more stringent.

A typical form of the Uj contains both growth and decay, and in local equilibrium these two terms are numerically equal; the growth term due to adding single adatoms by diVusion to j2 1 clusters is of the form sj Dn1nj21, where s is known as a capture number. The decay term has the form 2ndnj exp{2(Ed1DE)/kT}, where DE is the binding energy diference between j and j2 1 clusters. The key point is that if there is local equilibrium, then the ratio

nj /nj215n1Cexp(DE/kT),

(5.8)

where C is a statistical weight, which is a constant for a particular size (and con®guration) cluster. Note that this equilbrium must not depend on D, which is only concerned with kinetics. This argument can then be cascaded down through the subcritical clusters, yielding the Walton relation

n	5(n	)j o	C	(m)exp (E	(m)/kT),	(5.9)
j	1	m	j	j

where (m) denotes the mth con®guration of the j-sized cluster. This formula gives essentially the equilibrium constant, in the physical chemistry sense, of the polymerization reaction j adatoms 1 j-mer. It can thus be derived using classical statistical mechanics on a lattice, with N0 sites.

In the above formulae, ML units have been used for nj for simplicity, but sometimes the N0 is put in explicitly. In that case the nj are areal densities, and we need nj /N0 and n1/N0 in the above equation. You may note that at low temperature, we would only need to consider the most strongly bound con®guration, because of the dominant role of the exponential in (5.9). However, the critical cluster size is large typically at high temperature, so we need to be on our guard. If i51, at low temperature, the above discussion is not required anyway, since pairs of adatoms already form a stable cluster, and so are part of nx.

At this point, we do have something useful, because we can simplify the rate equa-

tions down to just two coupled equations, namely
dn1/dt5R ± n1/ta 2 (2U11oj ,iUj )2 sxDn1nx,	(5.10)

where the term in brackets is almost always numerically unimportant, and the last term describes the capture of adatoms by stable clusters, and can be written as n1/tc, and

dnx /dt5siDn1ni 2 Ucl.

(5.11)

In (5.11) the assumption of local equilibrium for ni, which is only a ®rst approximation, will make the ®rst term explicit, and proportional to the (i11)th power of the

152 5 Surface processes in epitaxial growth

adatom concentration: this is very non-linear if the critical nucleus size is large! The second term in (5.11) is typically due to coalescence of islands; this Ucl is proportional to nxdZ/dt, where Z is the coverage of the substrate by the stable islands. Thus dZ/dt is related to the (2D or 3D) shape of the islands, and how they grow (Venables 1973, 1987, Venables & Price 1975, Venables et al. 1984); all these details are not repeated here, but in the simplest 2D island case the last two loss terms in (5.10) equal NadZ/dt, where Na is the 2D density of atoms in the deposit. Note that Z is measured in ML, and is therefore dimensionless.

The capture numbers are related to the size, stability and spatial distribution of islands. The simplest mean ®eld model, which was worked on long ago (Venables 1973, Lewis & Anderson 1978, Stoyanov & Kaschiev 1981), and which others have worked on more recently (Bales & Chrzan 1994, Brune et al. 1999), is referred to as the uniform depletion approximation; it considers a typical cluster of size k immersed in the average density of islands of all sizes. Then one can set up an ancilliary diVusion equation for the adatom concentration in the vicinity of the k-cluster (size speci®c), or x-cluster (the average size cluster), which has a Bessel Function solution. This model gives exactly in the incomplete condensation (re-evaporation dominant) limit

sk52pXkK1(Xk)/K0(Xk) and	(5.12a)
sx52pXK1(X)/K0(X),	(5.12b)

where Xk25rk2/Dta and X25rx2/Dta, rk and rx being the corresponding island radii, and K0 and K1 the Bessel functions. For complete condensation the mean ®eld expressions are the same, but the arguments of the Bessel functions contain tc instead of ta; in general we should use t as de®ned in the next section. In complete condensation, these capture numbers are just functions of the coverage of the substrate by islands, Z.

The details of these capture numbers are the subject of problem 5.2, but for the moment we need to remember that they are simply numbers, with si in the range 2±4 and sx often in the range 5±10. Using these expressions one can compute the evolution of the nucleation density with time, or more readily with Z, as the independent variable, as ®rst proposed by Stowell et al. in the early 1970s. There are also other approximations for the various ss, such as the lattice approximation. However, the appropriateness of any of these depends on the spatial correlations between islands which develop as nucleation proceeds. The key point of Bales & Chrzan's (1994) paper was to show, for the particular case of i51 in complete condensation, that (5.12) with tc plus the other small terms as the argument, is the correct expression in the absence of spatial correlations, in agreement with their KMC simulations.

5.2.3Regimes of condensation

The above reasoning needs a bit of explaining to make it clear; if you are interested in the details it may be worth doing problems 5.1 and 5.2 at some point. Let us start by focusing our attention on the rate equation for the adatoms, where we can write as

dn	/dt5R2 n	/t , where t215t	211t	211t 21	1. . .	(5.13a)
1	1	a	n	c

	5.2 Atomistic models and rate equations		153
	Arrival (R)	Evaporation (τa )
	Nucleation	nl	Capture
			(τc )
	(τn )

nx

n i

Figure 5.5. Schematic illustration of the interaction between nucleation and growth stages. The adatom density n1 determines the critical cluster density ni; however, n1 is itself determined by the arrival rate R in conjunction with the various loss processes which have characteristic times as described in the text (Venables 1987).

We can de®ne the various time constants by comparison with (5.10), to obtain

t215(2U		1o	j,i	U	)/n	and t215s Dn		.	(5.13b)
n	1			j	1	c	x x

It is reasonable that the nucleation term (in brackets) is almost always numerically unimportant, as we have already convinced ourselves that Uj is close to zero for subcritical clusters. The ratio r5ta/tc5sxnxDta then determines whether we are in the complete (..1) or incomplete (,,1) condensation regime. So at high temperatures, t ta and at low temperatures (and/or high deposition rates), t tc. This is set out pictorially in ®gure 5.5, where the diVerent reaction channels are illustrated. It is useful to think of competitive capture; equations (5.13a,b) describe processes in which all the adatoms end up somewhere, and the diVerent competing processes (or channels) add as resistances in parallel.

Often, the condensation starts out incomplete, and becomes complete by the end of the deposition. This is the initially incomplete regime. If the diVusion distance on the surface is so short that only atoms which impinge directly on the islands condense, then we have the extreme incomplete regime. In each of these limiting cases, the two coupled equations ((5.10) for n1 and and (5.11) for nx ) can be evaluated explicitly, and give the nucleation density of the form nx ,Rp exp(E/kT), with p and E dependent on the regime considered. The formulae are given in table 5.1, and can be explored further via problem 5.1. Perhaps the most important regime for what follows is complete condensation. Re-evaporation is absent in this regime; one can notice from table 5.1 that the adsorption energy Ea does not appear in the expression for the cluster density.

In the complete condensation regime the (non-linear) interplay between nucleation and growth is most marked. In general, it is clear that if one includes diVerent processes, then one would expect to get diVerent power laws and energies. For example, Markov (1996) and Kandel (1997) explore diVerent models to study the role of surfactants in promoting layer growth during complete condensation, and obtain diVerent power laws which could be tested by experiments. There have been many discussions