A fundamental evolution equation is developed to describe the distribution of areas of capture zones (CZs) associated with islands formed by homogeneous nucleation and growth during submonolayer deposition on perfect flat surfaces. This equation involves various quantities which characterize subtle spatial aspects of the nucleation process. These quantities in turn depend on the complex stochastic geometry of the CZ tessellation of the surface, and their detailed form determines the CZ area distribution (CZD) including its asymptotic features. For small CZ areas, behavior of the CZD reflects the critical island size, i. For large CZ areas, it may reflect the probability for nucleation near such large CZs. Predictions are compared with kinetic Monte Carlo simulation data for models with two-dimensional compact islands with i = 1 (irreversible island formation by diffusing adatom pairs) and i = 0 (adatoms spontaneously convert to stable nuclei, e.g., by exchange with the substrate).

## I. INTRODUCTION

Homogeneous nucleation and growth of islands during submonolayer deposition on perfect low-index extended surfaces have been explored for decades as a route for self-assembly of arrays of supported nanoclusters or “islands.”^{1–3} This process involves a competition between simultaneous nucleation and growth of stable islands mediated by terrace diffusion of deposited atoms. Nucleation is often characterized by specifying a critical size i (measured in atoms) above which islands are stable and effectively immobile on the timescale of island formation. Thus, i = 1 corresponds to irreversible island formation where dimers and larger clusters are stable. For i > 1, island formation is reversible as sub-stable clusters of i or less atoms can dissociate. For i = 0, one regards isolated diffusing adatoms as being able to spontaneously convert to stable nuclei (e.g., by exchange with the substrate) about which stable islands then grow.^{4}

Theoretical and experimental analyses have focused on behavior of the mean island density, *N*_{isl}, and of the island size distribution (ISD).^{1–3} However, there has also been a long-standing appreciation of the significance of capture zones (CZs) associated with individual islands.^{5,6} The key idea underlying the CZ concept is that the surface can be tessellated into a distribution of CZs, one per island, where the area of the CZ for each island determines its growth rate.^{3,5–11} Recently, attention has turned to analysis of the CZ area distribution (CZD).^{12–21} A heuristic formulation by Pimpinelli and Einstein suggested that the CZD has a Gaussian tail for large areas, and that the small area behavior is related to the critical size, i, characterizing nucleation.^{14}

Our strategy for analysis of these surface deposition systems is to develop exact beyond-mean-field evolution equations for the basic quantities of interest: the ISD, the CZD, or related quantities such as the joint probability distribution (JPD) for island sizes and CZ areas.^{9,12,13,16} While exact analysis is not possible, these equations generally provide a solid foundation for reasonable approximate treatments. Also from this approach, one immediately recognizes the importance of characterizing spatial aspects of the nucleation process and of the resulting island distribution. This includes description of correlations between island size and island separation (or, equivalently, the dependence of CZ area on island size). Critically, here for analysis of the CZD, it is essential to provide a detailed description of various subtle spatial aspects of nucleation. Indeed, the central importance of the latter has been noted previously both in analytic theories based on the JPD,^{12,13} and in the development of an efficient “geometry based simulation” (GBS) procedure to describe island formation.^{22–24} Of particular significance is the feature that nucleation is enhanced by a higher adatom density. It therefore occurs preferentially far away from stable islands, the edges of which are sinks for diffusing adatoms, and thus nearby CZ boundaries.^{3,22–24}

In Sec. II, we describe our model for island formation, develop the fundamental equation for evolution of the CZD, and indicate expected scaling behavior of the solution. In Sec. III, we describe in detail key quantities characterizing the spatial aspects of nucleation. The scaling form of the fundamental evolution equation for CZD is described in Sec. IV, as well as the anticipated asymptotic form of its solutions for small and large CZ area. Scenarios for creating CZs with small areas are described in detail in Sec. V. Simulation results for i = 1 and i = 0 are presented and analyzed in Sec. VI. Further discussion and conclusions are presented in Sec. VII. A brief sketch of the basic formulation presented here and limited results were provided in a previous Comment.^{16}

## II. EVOLUTION OF CZ AREA DISTRIBUTION (CZD)

### A. Model prescription and basic behavior

Our atomistic model for island formation during submonolayer deposition includes the following basic ingredients: (i) Random deposition at rate *F* per site on a square lattice of adsorption sites. (ii) Hopping of isolated adatoms to empty nearest-neighbor (NN) sites with rate *h* per direction. (iii) Aggregation of adatoms into two-dimensional (2D) islands, which are stable above a critical size, i (measured in atoms). Nucleation requires aggregation of i + 1 atoms for i ≥ 1. Nucleation is mediated by exchange of single adatoms with the substrate for i = 0. Sub-stable islands can be regarded as in quasi-equilibrium with isolated adatoms. (iv) Periphery diffusion is facile so that islands relax to achieve compact square shapes during growth. (v) Atoms deposited on top of islands are instantaneously transported to the island edge and hop down. See Appendix A for further details.

We will sometimes comment on behavior for analogous “point-island” models where all islands occupy just a single lattice site but carry a label to indicate their size.^{25} These models avoid some complications seen in the above more realistic models for which the “footprint” of the 2D compact islands expands during growth. Basic features of the ISD, CZD, and JPD are similar for point and compact islands at lower coverage, although we find that some quantities characterizing the details of nucleation differ significantly.

We define the coverage, θ, from the deposition time, *t*, via θ = *Ft* in monolayers (ML). Then, the theory below is applicable to the regime of lower coverage, θ ≤ 0.1 ML, where the typical linear dimension of islands is well below the mean island separation. Let *N*_{s} denote the density of islands of s adatoms, where densities are measured per adsorption site. Thus, *N*_{1}, *N*_{i}, and *N*_{isl} = ∑_{s>i} *N*_{s} denote the densities of adatoms, critical clusters, and stable islands with size s > i, respectively. The Walton relation for quasi-equilibrated critical clusters implies that *N*_{i} ∝ (*N*_{1})^{i}.^{1,26} For i = 0, one naturally sets *N*_{0} ≈ 1. For i ≥ 1, in terms of the nucleation rate, *K*_{nuc} ∝ *hN*_{1}*N*_{i}, and the rate of aggregation, *K*_{agg} ∝ *hN*_{1}*N*_{isl}, of adatoms with stable islands, one has that^{1–3}

Integration of (1) for i ≥ 1 reveals a “transient regime” of rapidly growing *N*_{1} for θ below a crossover value^{3} θ^{∗} ∼ (*h*/*F*)^{−2/(i+3)} and a subsequent steady-state regime where adatom gain due to deposition roughly balances the loss due to aggregation so that *F* ≈ *K*_{agg}. The island density scales like *N*_{isl}^{∗} ∼ (*h*/*F*)^{−χ*} with χ^{∗} = (i + 1)/(i + 3) at crossover and like *N*_{isl} ∼ (*h*/*F*)^{−χ} with χ = i/(i + 2) in the steady-state regime.^{1–3} Since χ^{∗} − χ = 2(i + 2)^{−1} (i + 3)^{−1} > 0, nucleation in the steady-state regime dominates that in the transient regime for large enough *h*/*F* when i ≥ 1.^{3} Indeed, this dominance underlies selection of scaling forms for the ISD and CZD.

For i = 0, one must also specify an exchange probability, σ, per hop, so that the nucleation rate becomes *K*_{nuc} ∝ σ*hN*_{1}.^{4} In this case, there are also transient and steady-state regimes, although the characterization of behavior differs.^{3,27} Our version of this model is tailored to represent a “pure” i = 0 mechanism with no nucleation by aggregation of multiple diffusing adatoms. See Appendix A.

### B. CZ area distribution (CZD)

Let *N*_{A} denote the density of stable islands with CZ area *A*, where areas will be measured in units of the surface unit cell area (and distances will be measured in units of surface lattice constant). Also let *A*_{av} denote the mean CZ area. Then, one has that

*r*

_{isl}= (

*A*

_{av})

^{1/2}= (

*N*

_{isl})

^{−1/2}. To assess the evolution of

*N*

_{A}, we introduce the quantity

*P*

_{A}denoting the probability that the CZ of a just-nucleated island overlaps that of a pre-existing island with CZ area

*A*(and thus destroys a pre-existing CZ of area

*A*). Similarly, we introduce

*P*

^{+}

_{A}denoting the probability that formation of the CZ of a just-nucleated island reduces to

*A*the area of the CZ of a pre-existing island which was larger than

*A*(and thus creates a new CZ of area

*A*). We also let

*P*

^{∗}

_{A}denote the probability that the CZ of a just-nucleated island has area

*A*. Then, one obtains

^{16}the exact

*fundamental evolution equation*

Fig. 1 illustrates schematically the terms in this equation. We shall see in Sec. III A that *P*_{A} and *P*^{+}_{A} are related, and one has that ∑_{A}*P*_{A} = ∑_{A}*P*^{+}_{A} = *M*, where *M* denotes the mean number CZs of pre-existing islands overlapped by the CZ of a just-nucleated island. This equality, together with the normalization condition, ∑_{A}*P*^{∗}_{A} = 1, on *P*^{∗}_{A}, ensures consistency of (3) after applying ∑_{A} to both sides, and noting the identity *N*_{isl} = ∑_{A}*N*_{A}. We remark that there are some subtleties in the definition of *P*_{A}, etc., for smaller *A*_{av} where *A* are strongly discrete, and thus there is a chance that the CZ of a just-nucleated island can overlap two or more existing CZs of the same size. These issues are briefly discussed in Appendix B.

Note that the mean number of CZs overlapped by that of a just-nucleated island is *M* ≈ 5.5 for point-island models with i = 1,^{13} and we will find a similar value for compact islands. We emphasize that this result demonstrates the failure of a simple fragmentation picture where each nucleation event is regarded as simply fragmenting an existing CZ into two parts (corresponding to *M* = 1).^{28} The reason for the high *M*-values is that, as noted above, nucleation tends to occur near existing CZ boundaries where overlap with multiple CZs is expected.^{22–24}

Our focus will be on the asymptotic scaling regime of large *A*_{av} or small *N*_{isl} (corresponding to *h*/*F* ≫ 1 for i ≥ 1) where *N*_{A} is proposed to have the form^{12–21}

Unless stated otherwise, ∫dα ranges over 0 ≤ α ≤ ∞. Further comment on two aspects of this scaling are appropriate. *First*, we should note that deviations from the limiting scaling form can be significant for *h/F* ∼ 10^{6}-10^{7} considered here.^{3,29,30} The true limiting or asymptotic shape, *g*, of the distribution may only be achieved for much higher *h/F*. *Second*, the formulation (4) leaves implicit any dependence of g on coverage, θ, given that in this study we just focus on a single θ = 0.1 ML. However, there certainly is a dependence on θ, although we do not expect this dependence to be strong for models incorporating compact islands in the regime for small θ (specifically in the regime θ ≤ 0.1 which we consider). See Ref. 30 for a detailed discussion and analysis of θ-dependence in the context of the ISD. Regarding the CZD, a weak θ-dependence of *g* will require that key quantities characterizing nucleation also have weak dependence on θ. In Secs. III and V, we will offer some comments related to the origin of the θ-dependence for compact islands. From these observations, it is possible that basic behavior of the CZD for small CZ areas could depend on θ, but this is not expected of behavior for large areas. As an aside, we remark that θ-dependence should strictly only be absent in asymptotic scaling functions for point-island models.^{13,25,29,30}

Our main goal is to provide insight into the form of the scaling function *g*(α) based on simulation studies, and also on analysis of (3). The complexity of (3) undoubtedly *precludes* simple forms for this solution, providing a significant caveat for most existing analyses of simulation data for *g*(α). However, just as for the ISD, there is a natural appeal in developing simple fitting-forms which can be conveniently used for experimental data analysis. In this context, generalizing a strategy used in the field of stochastic geometry for analysis of tessellations,^{31–33} it is natural to consider fits to the generalized gamma (GG) distribution of the form^{16}

where

The value of “b” is determined by the unit mean constraint, and that of “a” from normalization. The special case where ν = 2 corresponds to the so-called Generalized Wigner Surmise (GWS) choice in Ref. 14. Another significant component of Ref. 14 was to suggest that β is determined by the critical size i for nucleation, and specifically that β_{GWS} = i + 1^{14} (later changed to β_{GWS} = i + 2^{17}) for island formation in a 2D surface deposition system. It should also be noted that traditional applications to analysis of tessellations in stochastic geometry usually set ν ≈ 1.^{31–33}

## III. DETAILED CHARACTERIZATION OF THE NUCLEATION PROCESS

### A. Characterization of key quantities describing nucleation

Application of (3) to determine the CZ area distribution, *N*_{A}, or the associated scaling function, *g*(α), requires characterization of the quantities *P*^{+}_{A}, *P*_{A}, and *P*^{∗}_{A} describing nucleation. For analysis of *P*^{+}_{A} and *P*_{A} and the relationship between them, we introduce a more detailed quantity, *P*_{A}(*A*_{sn}), giving the probability that the CZ of a just-nucleated island overlaps that of a pre-existing island with CZ area *A* by an amount *A*_{sn}, where *A*_{sn} < *A*.^{34} See Fig. 1. Then, it follows that

Summing both quantities over A immediately yields the identity $ \u2211 A P + A = \u2211 A P A $ (= *M*) mentioned in Sec. II B. See Appendix B.

Next, we note that since *P*_{A} should be proportional to *N*_{A}, it is natural to write

It is also convenient to write *P*_{A}(*A*_{sn}) in the factorized form

We regard *Q*_{A} as an “intrinsic” overlap probability, and *E*_{A}(*A*_{sn}) as characterizing the distribution of overlap areas for CZs of area *A*. In addition, it will be appropriate to introduce the quantity *A*_{subnuc}(*A*) giving the mean area of the portion of the CZ of a just-nucleated island which overlaps that of a pre-existing island with CZ area *A* (or equivalently, the mean area of the portion of the existing CZ of area *A* which is overlapped by the new CZ). It then follows that

Averaging *A*_{subnuc}(*A*) over all CZ areas, *A*, yields the quantity

defining μ_{av} as the average fractional overlap with each existing CZ. In fact, previous analyses for point islands with i = 1 found that μ_{av} ≈ 0.17 and also indicated that *A*_{subnuc}(*A*) ≈ μ_{av}*A* for all *A*.^{12,13} We will find a significantly lower value of μ_{av} for square islands at ∼0.1 ML.

Finally, it is also appropriate to introduce a quantity *A*_{avnuc} = ∑_{A}*AP*^{∗}_{A} corresponding to the average area for the CZs of just-nucleated islands. Note that α_{avnuc} = *A*_{avnuc}/*A*_{av} ≈ 0.97 for point islands with i = 1,^{13} but we will find a significantly smaller value for compact islands at ∼0.1 ML. It is tempting to equate *A*_{avnuc} with the product *MA*_{avsubnuc} = *M*μ_{av}*A*_{av},^{12,13} i.e., to identify α_{avnuc} with *M*μ_{av} ≈ 0.94 for point islands. However, this is only a (reasonable) approximation as it neglects the correlation between the number of CZs overlapped by a new CZ and the area of the portion of the new CZ overlapping each individual existing CZ.

### B. Scaling forms for key quantities describing nucleation

As indicated in Sec. II B, we will focus on scaling behavior, *N*_{A} ≈ (*N*_{isl}/*A*_{av}) *g*(*A*/*A*_{av}), of the CZD for large *A*_{av}. Correspondingly, for the key quantities described above, we introduce scaling functions

and

and

The simplest anticipated form of the scaling function in (12) for the intrinsic overlap probability, *q*(α), might be *q*(α) ≈ c α^{n}. The increase in *q*(α) with α reflects enhanced nucleation near larger CZs where the adatom density is generally higher. For i = 1, a simplistic analytic estimate might suggest n ≈ 3 (cf. Sec. V), but previous limited simulation data for point islands indicated a value of roughly n ≈ 1.5,^{13} and more complex forms with effective n varying with α were also suggested. In (13), *e*(μ) represents the scaling function for the distribution of fractional overlaps μ = *A*_{sn}/*A*.

Applying the above results to determine the key quantities in the fundamental evolution Equation (3), one obtains

In Appendix B, we confirm that these expressions preserve the equality ∑_{A}*P*^{+}_{A} = ∑_{A}*P*_{A}. Since *A*_{av} is large, one can replace the sum ∑_{Asn} in (15) with an integral ∫d*A*_{sn}. After making the change of variable to the fractional overlap μ_{sn} = *A*_{sn}/(*A* + *A*_{sn}), so that (*A* + *A*_{sn})/*A* dμ_{sn} = d*A*_{sn}/(*A* + *A*_{sn}), one obtains

where 〈…〉_{sn} = ∫μ_{sn}*e*(μ_{sn})... denotes a weighted average over fractional overlaps integrating over 0 ≤ μ_{sn} ≤ 1.

Given the somewhat obscure nature of the average in (16), it is instructive to consider a simple, although very crude approximation which greatly simplifies its form. If one assumes that the overlap area *A*_{sn} with a pre-existing CZ of area *A* is always close to *A*_{subnuc}(*A*) ≈ μ_{av}*A*, then the distribution *e*(μ) is sharply peaked around μ = μ_{av}. In the simplest approximation, one could set *e*(μ) ≈ δ(μ − μ_{av}), where δ is the Dirac delta function. Then, (16) reduces to

which imparts a simpler but non-local nature to the fundamental evolution equation analogous to that seen in evolution equations for fragmentation and coalescence processes. It is also readily confirmed that (17) preserves the equality ∑_{A}*P*^{+}_{A} = ∑_{A}*P*_{A}. However, again we emphasize the crude nature of the approximation yielding (17).

## IV. SCALING FORM OF THE EVOLUTION EQUATION FOR CZ AREAS

### A. Scaled evolution equation

Given the forms in Sec. III of various quantities describing nucleation, we recast (3) as an equation for the scaling function, *g*(α), of the CZ area distribution. From (4) and using that *N*_{isl} = 1/*A*_{av}, one obtains

Then substituting (18) and the scaling relations from Sec. III B for *P*_{A}, *P*^{+}_{A}, and *P*^{∗}_{A} into (3), we obtain the *scaled evolution equation*

Again 〈…〉_{sn} = ∫dμ_{sn}*e*(μ_{sn})... denotes a weighted average over fractional overlaps μ_{sn} = α_{sn}/(α + α_{sn}) which can be simplified using (17), but this crude approximation is not pursued here. As a consistency check, integration ∫dα of the left-hand side (LHS) of (19) using integration by parts yields ∫dα*g*(α) = 1. See Appendix B. Integrating the right-hand side (RHS), the first and second terms cancel exactly, and the last term yields ∫dα*p*^{∗}(α) = 1 matching the LHS.

Of course, the solution *g*(α) to (19) depends on the detailed form of the scaling functions for the intrinsic overlap probability, *q*(α), the CZ areas of just-nucleated islands, *p*^{∗}(α), and the distribution of fractional overlaps, *e*(μ). The complex form of (19) *precludes* a simple form such as GWS or GG for *g*(α). This is analogous to the situation for the scaling function for the ISD where the complex form of its exact evolution equation also strictly precludes commonly adopted simple forms.^{3}

### B. Asymptotic forms for solution

Despite the difficulty of exact analysis for (19), it is instructive to explore some approximations which might provide insight into asymptotic behavior of the solution *g*(α). First, *for large* *α*, one might anticipate that the major impact on large CZs is reduction in their size due to overlap by CZs of just-nucleated islands. (Such large CZs are perhaps created early at the end of the transient regime, and surviving into the steady-rate regime.) In this case, it is plausible that the first loss term on the RHS of (19) dominates the gain terms on the RHS, so that

where the integrals are from 0 to α and we utilize the simplified relation *q*(α) ∼ α^{n} introduced in Sec. III B and neglect any α-dependent prefactor. We emphasize the various assumptions or approximations leading to (21): a dominant loss term on the RHS of (19); the simplified form for *q*(α) and the approximation of the integral in the exponential. Nonetheless, this analysis suggests that ν ≈ n in (5) relating behavior of the large area tail of *g*(α) to details of the nucleation process.

Given uncertainties and approximations in the above analysis, it is appropriate to comment on and compare with the quite different heuristic analysis of Pimpinelli and Einstein.^{14} This approach anticipates more universal asymptotic behavior for the stochastic process of island and CZ formation, independent of the details of nucleation. Specifically, the CZ areas are assumed to satisfy a strongly damped Langevin equation subject to an effective force field. For large CZ areas, there is an effective “restoring force” which is roughly linear in area and has the effect of reducing area. This force is regarded as reflecting an effective pressure from surrounding CZs. From the corresponding Fokker-Planck equation, one concludes that the CZD has a Gaussian tail for large areas so that ν = 2.^{35} Since neither (20)-(21) or this heuristic Langevin formulation is rigorous, there remains some uncertainty.

Second, *for small* *α*, given that *q*(α) → 0 as α → 0, the first term on the RHS should be negligible compared to other terms. For simplicity, let us initially also assume that the third *p**-term dominates the RHS, i.e., that “direct creation” of new CZs with small area dominates (but see below). Then, terms on the LHS of (19) must balance the *p**-term yielding

If *p*^{∗}(α) ∼ α^{β}, it follows that also *g*(α) ∼ α^{β} as α → 0, and specifically that^{16}

This analysis is consistent with assumption in previous analysis^{14} that *g*(α) for small α is determined by the form of the distribution of new CZ areas, *p*^{∗}(α), for small α.

However, our picture for the formation of small CZs is quite different from Ref. 14 (see Sec. V), and furthermore is actually more complicated than suggested above. New CZs with small area can also be created “indirectly” by the CZs of just-nucleated islands mostly overlapping those of existing islands with typical areas. This corresponds to the second term on the RHS of (18) with small α and α_{sn} = O(1), and thus does not involve *p*^{∗}(α) for small α. More generally, suppose that the overall rate of creation of new CZs with small area α, combining contributions from the second and third terms on the RHS of (19), has the scaling form *p _{new}*(α) ∼ α

^{β}. Then, it still follows from (19)) that

*g*(α) ∼ α

^{β}, as α → 0.

## V. SCENARIOS FOR CREATING SMALL CZs

The scenario considered by Pimpinelli and Einstein for creating small CZs involves nucleation of new islands at a location on the surface surrounded by multiple pre-existing CZs having small size.^{14,17} However, such regions are extremely rare, their population scaling like a higher power of the population of small CZs, a feature not analyzed in these previous studies.^{36} Instead, we claim the dominant pathway for creating small CZs does not involve pre-existing small CZs, i.e., that new small CZs are primarily created in the *absence* of pre-existing small CZs.^{16} Below, we present and analyze various key scenarios, although a comprehensive analysis of all possibilities is difficult.

### A. Direct creation of new CZs with small area

One scenario is where a small CZ is created by nucleation of an island in between the sides of two (m = 2) closely spaced near-square islands. See Fig. 2(a). We note that the likelihood of this scenario is dependent on the structure of the islands. Here, the finite extent of islands is key in facilitating creation of a small CZ between a nearby existing pair. Also higher coverage, θ, will facilitate this scenario. Another scenario is nucleation of an island in the middle of a triangle of three (m = 3) closely spaced islands. See Fig. 2(b). In both these cases, the small CZ is associated with the just-nucleated island, and thus with *P*^{∗}_{A} or *p*^{∗}(α) described in Sec. III. Of course, one could also consider nucleation in the middle of larger collections of nearby islands (m > 3). However, analysis below indicates that the likelihood of such events is negligible. A semi-quantitative analysis of the *p*^{∗}(α) associated with these scenarios has *two components*. *First*, we must determine the small population, *P*_{m}, of the above “rare” CZ configurations. *Second*, we must determine the probability, *P*_{nuc}, that nucleation occurs in the appropriate restricted sub-region as needed to create a desired small CZ.

For the *first* part of the analysis, consider the density, *N*_{isl}(**r**), of islands with their centers separated from the center of a specified island by a vector separation **r**. One has that *N*_{isl}(**r**) ≈ *C*_{isl}(**r**) *N*_{isl}, where the island-island correlation function satisfies *C*_{isl}(**r**) ∼ [*N*_{1}(**r**)/*N*_{1}]^{i+1}, and *N*_{1}(**r**) represents the adatom density at position *r* from the center of an island.^{1,3} *N*_{1}(**r**) decreases linearly to zero at the island edge,^{1,3} so we expect that *N*_{1}(**r**)/*N*_{1} ∼ r/r_{isl} for closely spaced islands, where r = |**r**| and again r_{isl} is the mean island separation. We conclude that the relative probability to find a pair (m = 2) or trio (m = 3) of closely spaced islands with separations r ≪r_{isl} should scale like

The factor of (*C*_{isl})^{m−1} gives the probability for a specific orientation of the m-tuple of islands, and the additional factor of (r/r_{isl}) accounts for the volume of phase-space for arbitrary rotations.

The *second* factor is the relative probability for nucleation in an appropriate restricted subregion, as indicated above, separated by island edges by a distance of order r. First, we note that the total rate of nucleation, *K*_{CZ}, for a typical CZ scales like

which follows from the steady-state results from Sec. II B. This type of scaling of nucleation rates with linear size was a key component of previous geometry-based simulations of island nucleation and growth.^{22–24} Next, we consider the total nucleation rate for the above-mentioned restricted sub-region of area ∼r^{2} denoted by *K*_{SUB}. Let *N*_{1}(loc) denote the typical value of the adatom density in this subregion. Then, it follows that *K*_{SUB} scales like

Thus, one has that

Combining the results (24) and (27), and using that α ∼ (r/r_{isl})^{2} ≪ 1, the likelihood to create a new small CZ by the above mechanisms scales like

Thus, one has β_{2} = 3i/2 + 3 and β_{3} = 2i + 7/2. For compact islands with finite extent, one expects that contributions from m = 2 (i.e., the smaller β_{m}) will dominate.

### B. Major reduction in size of existing CZs

Another scenario can generate island and CZ configurations similar to those described in Sec. V A, but with a different order of island nucleation. This scenario will be most viable for m = 3, i.e., a triangular triple of nearby islands surrounding a fourth, rather than for m = 2. The idea for m = 3 is that islands could first nucleate creating two corners of the triangle, then a third island nucleates which will end up in the middle, and finally a fourth nucleates constituting third corner. See Fig. 3. In this case, the CZ of the last just-nucleated island is not generally small. However, the CZ of the third nucleated island, which was generally large upon formation, is greatly reduced (i.e., mostly covered) upon nucleation of the fourth island. In this case, the formation of a small CZ is associated with the term *P*^{+}_{A} in d*N*_{A}/d*t*. Specifically, it reflects the small but non-zero weight in the upper end of the fractional overlap distribution, *e*(μ), with μ close to unity.

For a semi-quantitative analysis of the scenario in Fig. 3, the probability to nucleate the first three islands in the above configurations should scale as *P*′_{3} ∼ (r/r_{isl})^{2i+3} just as does P_{3}. However, the last “outer” island will nucleate in a region where *N*_{1}(*r*)/*N*_{1} ∼ r/r_{isl} noting the linear increase of adatom density moving away from the first three islands. This suggests a relative nucleation probability of *P*′_{nuc} ∼ (r/r_{isl})^{2} [*N*_{1}(*r*)/*N*_{1}]^{i+1} ∼ (r/r_{isl})^{i+3}. Using that α ∼ (r/r_{isl})^{2}, the overall rate for formation of the small CZ scaling like

(versus β_{3}(i) = 2i + 7/2 in Sec. V A). Since β′_{3} has the same value as β_{2} given in Sec. V B, this process should also make a significant contribution to the measured β. In fact for point islands, where the m = 2 pathway of Sec. V A is not effective, one expects this alternative pathway for m = 3 to dominate.

## VI. KINETIC MONTE CARLO (KMC) SIMULATION RESULTS FOR SQUARE ISLANDS

Efficient KMC simulation for large values of *h*/*F* utilizes a standard Bortz-algorithm to enable collection of extensive statistics on island configurations. To assess quantities such as *g*(α), CZ tessellations are constructed which here are taken as edge cells, so that CZ boundaries are equidistant from island edges.^{10,11,22–24} However, to collect adequate statistics at a selected coverage for quantities associated with new islands such as *p*^{∗}(α), the following refined approach is implemented:^{9,12,13} (i) we perform a standard simulation to generate an island distribution at the desired coverage; (ii) we continue simulation until one new island is nucleated and measure desired quantities associated with this new island; (iii) we remove the new island and continue simulation until another new island is nucleated and perform desired measurements; (iv) we repeat this process for ∼100 nucleations per configuration generated in (i). This approach is far more efficient than collecting information on just one new island from each simulated configuration.

### A. Behavior for i = 1

We focus below on analysis of data from simulations with i = 1 and *h*/*F* = 10^{6} at θ = 0.1 ML. The total number of islands and thus CZs used to assess *g*(α) is ∼268 500, with average area *A*_{av} = 372.4. With regard to detailed analysis of quantities such as *p*(α), *p*^{∗}(α), *e*(μ), etc., we consider a total of 50 000 new CZs which have an average area of *A*_{avnuc} = 201.6. Key parameters extracted from this analysis include α_{avnuc} = *A*_{avnuc}/*A*_{av} = 0.54, μ_{av} = 0.099, and *M* = 4.65. Note that the α_{avnuc} and μ_{av} values are significantly lower than for point islands. This feature is a natural consequence of the finite extent of the square islands at 0.1 ML reducing the relevant full or partial CZ areas. Also, as for point islands, α_{avnuc} is somewhat above *M*μ_{av} ≈ 0.46.

We start by performing least-square fits of GG forms to *g*(α).^{37} See Fig. 4. First, this is done by selecting the reasonable choices of n-value and fitting β. For *h*/*F* = 10^{6}, this yields β = 3.05 for n = 2, β = 3.48 for n = 1.75, and β = 4.07 for n = 1.5. One can also simultaneously adjust both β and n. For *h*/*F* = 10^{6}, this yields β = 3.84 and n = 1.59. As an aside, from more limited simulation data, we obtain β = 3.55 and n = 1.65 for *h*/*F* = 10^{7}, and β = 3.36 and n = 1.70 for *h*/*F* = 10^{8}.^{38}

The advantage of above fitting approach, which assumes that *g*(α) is described by a GG form, is that one utilizes simulation data for all α. These data are extensive for CZ areas around the average value (α ∼ 1), thereby compensating for limited data in the upper and lower tails of the distribution. However, there are also significant disadvantages. *First*, the fit is expected to be sensitive to deviations for finite *h*/*F* from asymptotic behavior as *h*/*F* → ∞. Variation with *h*/*F* in the peak if the CZD has been assessed previously and significantly impacts estimates of β.^{15} *Second*, a more fundamental concern is that *g*(α) will not be exactly described by a GG form (as indicated above). This feature can introduce some bias and correlation in the resulting values of β and n (e.g., higher n produces smaller β).

The above observations motivate direct unbiased fitting of the lower tail of *g*(α) to assess β, and thus the proposed relationship to i. Here, one must choose the range of α over which to fit. Lower (and narrower) ranges better describe asymptotic behavior as α → 0, but these suffer from limited data. Thus, we give β-values from fitting over different ranges with results shown in Table I. The decrease in estimated β with broader range of α just reflects non-asymptotic behavior as log g(α) versus log α is convex. Thus, extrapolating values for broader ranges to narrower ranges should give the most reliable estimates. Results from analysis of extensive data for *h*/*F* = 10^{6}, and more limited data for higher *h*/*F*, consistently suggest that β is around 3.5.

α (h/F = 10^{6}) | 0.04-0.3 | 0.04-0.4 | 0.04-0.5 | 0.04-0.6 |

β (h/F = 10^{6}) | 3.37 | 3.34 | 3.28 | 3.20 |

α (h/F = 10^{7}) | 0.07-0.3 | 0.07-0.4 | 0.07-0.5 | 0.07-0.6 |

β (h/F = 10^{7}) | 3.30 | 3.28 | 3.19 | 3.09 |

α (h/F = 10^{8}) | 0.1-0.3 | 0.1-0.4 | 0.1-0.5 | 0.1-0.6 |

β (h/F = 10^{8}) | 3.49 | 3.29 | 3.14 | 2.98 |

α (h/F = 10^{6}) | 0.04-0.3 | 0.04-0.4 | 0.04-0.5 | 0.04-0.6 |

β (h/F = 10^{6}) | 3.37 | 3.34 | 3.28 | 3.20 |

α (h/F = 10^{7}) | 0.07-0.3 | 0.07-0.4 | 0.07-0.5 | 0.07-0.6 |

β (h/F = 10^{7}) | 3.30 | 3.28 | 3.19 | 3.09 |

α (h/F = 10^{8}) | 0.1-0.3 | 0.1-0.4 | 0.1-0.5 | 0.1-0.6 |

β (h/F = 10^{8}) | 3.49 | 3.29 | 3.14 | 2.98 |

As indicated in the theoretical development of Secs. II–V *g*(α) behavior is controlled by subtle spatial aspects of the nucleation process. Thus, it is appropriate to also examine simulation results for the key scaling functions *p*^{∗}(α), *p*(α), *q*(α), and *e*(μ). See Fig. 5. Consistent with the (over)simplified analysis of (22) and (23), the distribution, *p*^{∗}(α), of areas of new CZs appears to scale in the same way as *g*(α) for small α. Also, *g*(α) is naturally smaller than *p*^{∗}(α) in the small-α range. However, direct creation of new small CZs associated with *p*^{∗}(α), as described in Sec. V A and as assumed in (23), is not necessarily the dominant pathway. Rather, examination of simulation movies reveals a significant fraction of cases where such small CZs are formed indirectly from a substantial reduction in area of existing CZs, as described in Sec. V B. See Fig. 6. While the associated values of the fractional overlap scaling function, *e*(μ), for larger μ are small, so are values of *p*^{∗}(α) for correspondingly small α ∼ 1 − μ. (We should also note that *e*(μ) is far from the delta function form corresponding to a narrow distribution of fractional overlaps where typically μ ≈ μ_{av}. Instead, small fractional overlaps are common.)

With regard to behavior of *g*(α) for large α, as anticipated, the intrinsic overlap probability, *q*(α), increases non-linearly with α. The precise form is not clear from the limited data, but a simple fit to *q*(α) ∼ α^{n} yields n ≈ 1.5 with considerable uncertainty.

### B. Behavior for i = 0

For the i = 0 model, the island density scales like *N*_{isl} ∼ σ^{1/2} independent of *h*/*F*.^{3,4} Most studies considered the regime of “rapid conversion”^{4,19} with σ ≫ (*h*/*F*)^{−1/2} yielding persistent nucleation in the steady state and a monotonically decreasing ISD.^{3,4} A crossover occurs reducing σ to σ = (*h*/*F*)^{−1/2} and below to growth behavior avoiding persistent nucleation. Thus, in this work, we select σ = (*h*/*F*)^{−1/2} which leads to relatively low *N*_{isl} ∼ (*h*/*F*)^{−1/4}. Our construction and analysis of CZs includes single converted adatoms as well as larger clusters, but not isolated diffusing adatoms. Even at a low θ ∼ 0.1 ML, there is significantly more island coalescence for i = 0 than for i = 1 consistent with the larger island-island correlation function, *C*_{isl}(**r**) ∼ (r/r_{isl})^{i+1}, for small r. Since our theory is geared to the pre-coalescence regime, we also consider behavior for lower θ = 0.04 ML to limit coalescence.

We present a more limited analysis of behavior for i = 0 than for i = 1 focusing on the case *h*/*F* = 10^{6} and σ = 10^{−3} for θ = 0.1 ML. The total number of islands and thus CZs used to assess *g*(α) is ∼346 500 with an average CZ area *A*_{av} = 192.71. Least-squares fitting of a GG form adjusting both β and n yields β = 2.80 and n = 1.42 for θ = 0.1 ML. One obtains β = 3.46 and n = 1.42 for θ = 0.04 ML. Given the limitations or bias in least square fitting to GG forms as discussed for i = 1, we do not give a more extensive analysis here. The reader is referred to Fig. 7 where GG forms are compared with simulation data for 0.1 ML for various choices of β and n. Choosing n = 2 requires β values above 2.5 for a reasonable fit. It appears that the data are better fit by lower n values and higher β closer to 3 which is consistent with β ≈ β_{2}(0) = 3 for i = 0 in (27). As for i = 1, direct fitting of the lower tail of g(α) to assess β avoids bias of an assumed GG form, and is thus most reliable provided there is sufficient simulation data. Results are shown in Table II which suggest again that β is just below 3.0.

α (θ = 0.1 ML) | 0.057-0.3 | 0.057-0.4 | 0.057-0.5 | 0.057-0.6 |

β (θ = 0.1 ML) | 2.74 | 2.64 | 2.53 | 2.42 |

α (θ = 0.04 ML) | 0.056-0.3 | 0.056-0.4 | 0.056-0.5 | 0.056-0.6 |

β (θ = 0.04 ML) | 2.75 | 2.76 | 2.71 | 2.64 |

α (θ = 0.1 ML) | 0.057-0.3 | 0.057-0.4 | 0.057-0.5 | 0.057-0.6 |

β (θ = 0.1 ML) | 2.74 | 2.64 | 2.53 | 2.42 |

α (θ = 0.04 ML) | 0.056-0.3 | 0.056-0.4 | 0.056-0.5 | 0.056-0.6 |

β (θ = 0.04 ML) | 2.75 | 2.76 | 2.71 | 2.64 |

Finally, we remark that analysis of simulation movies suggests that as for i = 1, a significant fraction of small CZs is created indirectly rather than directly. See Fig. 8.

## VII. DISCUSSION AND CONCLUSIONS

Recent extensive interest in the capture zone area distribution (CZD) has been prompted by the proposal of Pimpinelli and Einstein (PE) that the CZD has a GWS form, and furthermore that the exponent, β, describing behavior for small areas is directly related to the critical size.^{14,17,21} Considering *2D systems for* *i* ≥ 1, Shi, Shim, and Amar (SSA) utilized point-island models to critique the original proposal of PE that β_{WS} = i + 1 = 2 for i = 1.^{15} Their simulation data for typical *h*/*F* = 10^{7} fitted to a GWS form instead suggested that β ≈ 3. Furthermore, SSA demonstrated that key features of the scaled CZD, particularly the peak height, had a significant dependence on *h*/*F*. Also, slightly different versions of the point-island model had significantly different peak heights for smaller *h*/*F*. However, extrapolation to *h*/*F* → ∞ removed discrepancies between different versions of the model, as might be expected. This analysis again indicated that β ≈ 3 *assuming* that the CZD, and particularly its peak, is described by a GWS form. This analysis appears consistent with the modified proposal of PE that β_{WS} = i + 2 (= 3 for i = 1).^{17} However, as discussed above, there are limitations to an analysis based on fitting to a GWS (or GG) form which does not correspond to the exact form of the CZD.

Additional studies by Oliveira and Aarão Reis^{18–20} considered the CZD for point, square, and fractal island models utilizing a modified scaling formulation to reduce dependence of the scaled CZD upon h/F. However, this modified formulation was not designed to highlight or analyze behavior for small CZ areas. Point-island behavior for i = 1 is fit reasonably but not perfectly by β ≈ 3 and n = 2, or by β ≈ 4 and n = 1.5. Square island behavior is also fit reasonably by β ≈ 3 and n = 2, but with significant discrepancies in the tails of the CZD. A key focus of these studies was on the form of the right tail of the CZD for large areas. For i ≥ 1, simulation data for point islands and also for square islands at pre-coalescence coverages appear consistent with a Gaussian tail, i.e., with n = 2, or at least with n close to but probably not above 2. For square islands at higher coverages where coalescence is significant, there was a clear deviation from Gaussian towards exponential tails for large CZ areas.

Our own analysis for *2D systems with* *i* ≥ *1* is motivated by our development of an exact evolution equations for the CZD. This equation suggests (but does not rigorously prove) the possibility of a non-Gaussian tail for large areas. It also provides and incorporates a different picture from PE of behavior for small areas. Specifically for the latter, appropriately identifying and analyzing the dominant pathways for creation of small CZ areas (which in our picture do not require pre-existing small CZs), we predict larger β-values than the refined PE prediction of β_{WS} = i + 2. Furthermore, we argue that the most reliable assessment of β comes from direct fitting of CZD data for small α. This contrasts the strategy of previous studies of fitting to a GWS or GG form (which can bias the outcome since the CZD is not exactly described by these forms). Such direct fitting does indicate somewhat larger values of β than β_{WS} = i + 2 = 3 for i = 1, although the values are below those for the pathways which we identified for creation of new small CZs. We note one limitation of our analysis is that data are obtained for relatively small *h*/*F* (as in analyses of PE). Data for higher *h*/*F* are preferred to avoid possible discrete lattice effects for small α, but its acquisition is computationally more demanding. With regard to behavior for large areas, our simulation data are not sufficiently extensive to make a definitive assessment. Behavior is consistent with Gaussian tails (n = 2) but could also be described well by somewhat smaller n for i = 1. Also taking into consideration the analysis of Refs. 18 and 19, one concludes that n is close to 2 (and that it does not exceed 2).

Considering *2D systems for* *i* = 0, our analysis indicates β values significantly above the refined PE prediction of β_{WS} = i + 2 = 2 for i = 0. It also seems likely that the large area tail of the CZD is described by n < 2 rather than by a Gaussian tail. Previous studies for i = 0 with larger exchange rates, σ > > (*h*/*F*)^{−1/2}, and thus persistent nucleation, indicated exponential tails for large areas.^{19} This is certainly reasonable since persistent random nucleation produces quasi-random island distribution for which one expects a CZD closer to that of a Poisson point process. For the latter, a Gamma function form with exponential tail is expected.^{29–31} However, for our choice of lower exchange rate σ = (*h*/*F*)^{−1/2} which eliminates such persistent nucleation,^{3} the form of the large area tail is less clear *a priori*.

In summary, our exact evolution equation for the CZD provides a valuable additional formulation and perspective beyond a heuristic treatment for analysis and understanding of this distribution. In particular, our formulation highlights the importance of subtle spatial details of the nucleation process, and specifically of related quantities such as the intrinsic overlap probabilities and fractional overlap distribution. These quantities have received limited attention and are not considered in heuristic treatments, but their detailed form impacts that of the CZD. Given the complexity of the process of homogeneous island nucleation and growth of islands during deposition, there does not yet exist a definitive characterization and understanding of either the ISD or the CZD. However, this goal is facilitated with the aid of detailed analytic theory such as that developed here, together with more extensive simulation targeting key quantities characterizing the spatial aspects of nucleation. Another challenge is that extensive simulation data are most readily obtained for smaller *h*/*F* = 10^{6}-10^{7}, but one is primarily interested in the scaling limit *h*/*F* → ∞, so reliable extrapolation is also required.^{29}

## Acknowledgments

We thank Tiago Oliveira for valuable comments on the manuscript. Y.H. and J.W.E. were supported for this work by NSF Grant Nos. CHE-1111500 and CHE-1507223. Computations utilized NSF-supported XSEDE resources. M.L. was supported by NSF of China under Grant No. 51271197. The work was performed at Ames Laboratory which is operated for the USDOE by Iowa State University under Contract No. DE-AC02-07CH11358.

### APPENDIX A: ATOMISTIC MODELS FOR ISLAND FORMATION

Our tailored models for island formation with i = 1 involve: (i) Random deposition at rate *F* on a square lattice of adsorption sites with periodic boundary conditions. (ii) Isolated adatoms hop to a NN site with hop rate *h* per direction. (iii) Typically, singly coordinated “edge atoms” can hop at rate *h* to either NN or next-NN edge sites provided they are still coordinated (i.e., they cannot detach from islands). (iv) Multiply coordinated atoms are immobile. (v) Atoms deposited on top of islands immediately undergo a random walk until they reach the edge of the island and hop down to the substrate. Also, atoms landing on isolated adatoms immediately hop down.

Incorporation of facile diffusion of atoms along the island periphery both along straight edges and around corners or kinks suffices to achieve compact islands.^{3} Implementation of rapid transport of atoms landing on top of islands to the first layer avoids second layer population.

One extra constraint not described above is that singly coordinated adatoms in clusters of size s = 2 or 3 are not allowed to undergo hops to next-NN edge sites. This constraint excludes long-range diffusion of these small clusters, which if present would have the effect of changing the scaling behavior of *N*_{isl}.^{3} For clusters of 4 or more atoms, restructuring will produce a rectangular core of doubly coordinated adatoms which is permanently frozen, thus excluding long-range diffusion.

Our version of a “pure” i = 0 model excludes nucleation by aggregation of multiple diffusing adatoms. This is achieved by preventing isolated diffusing adatoms from reaching NN sites, and also exclude deposition on top of or next to isolated adatoms. Island nucleation occurs exclusively by isolated adatoms converting to permanently immobile adatoms at rate σ*h*. Diffusing adatoms aggregate irreversibly with converted immobile adatoms and with clusters of adatoms which form around these, i.e., detachment from such islands is prohibited. Edge diffusion is included as for i = 1 ensuring the formation of compact islands. Atoms deposited on top of islands are also treated as for i = 1. In our construction and analysis of CZs for i = 0, we will count as islands single converted adatoms in addition to larger clusters, but not isolated diffusing adatoms.

We will not present results for i > 1, however we mention that the model for i = 1 can be naturally extended to i = 2 or 3. For example, for i = 2, singly coordinated edge atoms can only detach at rate *h* if their unique neighbor is also singly coordinated (i.e., if they are part of an ad-dimer).^{22} Otherwise detachment is excluded. Doubly coordinated adatoms are again frozen. In these models with i ≥ 1, edge diffusion, as described above, ensures the formation of compact islands. For KMC images of compact islands in such models with i = 1-3, see Refs. 23 and 24.

### APPENDIX B: OVERLAP PROBABILITIES AND CONSISTENCY RELATIONS

First, we comment on the regime of smaller *A*_{av} so that CZ areas *A* are discrete. Then, there is a significant probability that the CZ of a just-nucleated island can overlap k = 2 (or more) pre-existing CZs of the same size, *A*. Then, in determining *P*_{A}, one considers a large number, j = 1 − N, of nucleation events and counts the number k_{j} of overlaps in each case, setting *P*_{A} ≈ ∑_{j} k_{j}/∑_{j} 1 = ∑_{j} k_{j}/N. A similar treatment is implemented for related quantities *P*^{+}_{A} and *P*_{A}(*A*_{sn}).

Second, we show that $ \u2211 A P A + = \u2211 A P A $ exploiting a change of variable *B* = *A* + *A*_{sn},

The same change of variable *B* = *A* + *A*_{sn} applies to simplify

Third, we comment on the consistency of the evolution Equation (19) for the scaling function which follows from integrating over α. After integration, the second term on the RHS of (19) yields

which cancels the first term given the normalization of *e*(μ_{sn}). Finally, we remark that for the most simple δ-function choice of e(μ_{sn}) ≈ δ(μ_{sn} − μ_{av}) (which again has limitations), (19) reduces to the non-local equation

## REFERENCES

*Current Topics in Materials Science*, Volume 3, edited by E. Kaldis (North Holland, Amsterdam, 1979), Chap. 4, pp. 421-462.

A_{sn} corresponds to A_{subnuc} in previous papers.^{10}

If d/dt A = F(A) + η(t) with F(A) ∼−kA for large A, and noise η satisfying 〈η(t)〉 = 0 and 〈η(t)η(t′)〉 = Λδ (t – t′), then one has that P_{A} ∼ exp[2Λ^{−1}∫dAF(A)] ∼ exp[ − kΛ^{−1} A^{2}].

Ref. 14 used the local nucleation rate k(loc) ∼ A^{i+1} to estimate the rate of formation of CZs with small area A in regions of the surface with other CZs of that size. Ref. 17 instead used the total nucleation rate, K(loc) ∼ A k(loc) ∼ A^{i+2} integrated over an area A. The latter can be obtained by integrating the local adatom density ∫_{A} d**r** N_{1}(**r**)^{i+1}, but such analysis is not necessary as N_{1}(**r**) adopts similar values over most of the area. These analyses did not consider the probability of finding the configurations of existing CZs in which they regarded formation of new small CZs as occurring.

Fitting of a GG form to p*(α) yields β = 3.52 selecting n = 2 and β = 4.01 selecting n = 1.75, and β = 4.69 selecting n = 1.5. Adjusting both β and n to fit p* yields β = 3.70 and n = 1.90.

For h/F = 10^{7} we also obtain β = 2.94 for n = 2, and β = 3.92 for n = 1.5. For h/F = 10^{8}, we obtain β = 2.87 for n = 2 and β = 3.83 for n = 1.5.