Controlled Self-Assembly of Nanocrystalline Arrays Studied by 3D

ARTICLE pubs.acs.org/JPCC

Controlled Self-Assembly of Nanocrystalline Arrays Studied by 3D Kinetic Monte Carlo Modeling Abuhanif K. Bhuiyan,*,†,‡ Steven K. Dew,*,† and Maria Stepanova*,†,‡ † ‡

Department of Electrical and Computer Engineering, University of Alberta, Edmonton, Alberta, Canada T6G2V4 National Institute for Nanotechnology NRC, Edmonton, Alberta, Canada T6G2M9 ABSTRACT: Fabrication of self-assembled arrays of nanocrystals (NCs) by physical vapor deposition (PVD) is a promising technique rated highly for its potential for various electronic, photonic, and sensing applications. However, the self-assembly process is not straightforward to control and direct in a desired way. A detailed understanding of how to control the size, shape, and density of self-assembled NCs by varying the accessible PVD process conditions, such as deposition rate, duration, or temperature, is critical for the potential of self-assembled nanofabrication to be fully realized. In this paper, we report a systematic kinetic Monte Carlo modeling that explicitly represents PVD synthesis of selfassembled metallic NCs on a crystalline substrate. We investigate how varying the duration of deposition, deposition rate, temperature, and substrate wetting conditions may affect the morphologies of arrays of self-assembled metallic islands and compare our results with previously reported experimentally observed surface morphologies generated by PVD and theoretical studies.

1. INTRODUCTION Fabrication of self-assembled nanocrystals (NCs) by physical vapor deposition (PVD) is a promising technique rated highly for its potential for various electronic, photonic, and sensing applications.16 Employing naturally occurring self-assembly of nanostructures may result in very efficient and inexpensive fabrication techniques. The challenge, however, is that self-assembly processes are not straightforward to control and direct in a desired way. Before subtle interactions at the level of individual atoms could be efficiently controlled by varying the accessible PVD process conditions, such as deposition rate or temperature, a detailed and systematic understanding of the intimate relationship between individual atomic interactions and mesoscopic properties at larger scales should be developed. The major processes of PVD nanocrystalline synthesis, the nucleation and growth of nanocrystalline islands on a substrate, have been extensively addressed in the literature.712 Figure 1 illustrates the growth of Cu islands on a crystalline Cu substrate during a PVD process.13,14 Figure 1a shows an early postnucleation stage when tiny metallic islands have just arisen, and Figure 1b,c demonstrates how the morphology evolves. It can be seen that multilayer structures arise even before the first monolayer is complete and some of the new nucleation events take place on top of previously deposited layers. It is possible to identify a characteristic rectangular shape that is observed both at the stage of island nucleation (Figure 1a) and during the subsequent morphology evolution (Figure 1b,c) The mounded structure and the square or rectangular symmetry of the morphology are apparent also in Figure 1df, illustrating later stages of the deposition. The evolution of the mounded morphology, r 2011 American Chemical Society

comprising the coarsening of the mounds as well as merging of their bases, is evident in the images. The basic physical mechanisms of thin film deposition are known in significant detail.7,10,15 Initially, the incident atoms randomly impinge on the substrate where they are adsorbed and can subsequently re-evaporate or diffuse over the surface.16 The process of nucleation then takes place when isolated adatoms form tiny nanoislands on the substrate. The minimum number of atoms that makes a nucleus stable is known as the critical nucleation number.1721 Only the stable nuclei would steadily increase in size, giving rise to an array of NCs. At room temperature, the size of NCs increases largely due to the arrival of new adatoms at the surface, and also through coalescence of neighboring NCs.8 If the deposition continues, the bases of NCs merge and a multilayer thin film forms with a roughened surface. Mechanical, electrical, and optical properties of surface nanostructures strongly depend on the sizes, shapes, and densities of the fabricated islands,22,23 and thus a detailed understanding of the various deposition conditions (deposition rate, temperature, materials involved, etc.) on the morphology of the growing NCs is important to controlling the functionality and reliability properties of devices to be fabricated. Detailed numerical modeling of self-assembled synthesis of NCs exploring the physics behind the NCs’ fabrication is a subject of constantly growing interest. Kinetic rate equations addressed widely in literature8,9 describe the evolution of growing islands, whereas scaling laws derived from the rate equations relate the average characteristics, Received: June 20, 2011 Revised: August 9, 2011 Published: August 22, 2011 19557

dx.doi.org/10.1021/jp205791t | J. Phys. Chem. C 2011, 115, 19557–19568

The Journal of Physical Chemistry C

ARTICLE

Figure 1. STM images (ac) of Cu on Cu(001) for a size of 100 100 nm at T = 296 K with coverage θ, (a) 0.08 and (b) 1.0 ML at R = 0.16 ML/s and (c) 8 ML at R = 0.2 ML/s, obtained by physical vapor deposition.13 The STM images (df) of Cu on Cu(100) at R = 0.0208 ML/s and T = 299 K produced by physical vapor deposition. Images (d) and (e) show an area of 100 100 nm, whereas (f) is 500 500 nm.14 Reproduced with kind permission from Springer Science+Business Media B.V. and the American Physical Society.

such as the number density of stable nuclei, with the various process parameters.811 Thus, the density of stable islands in a steady-state condition can be expressed by p R Esurf exp N ¼ An0 ð1Þ v kT where N is the density of stable islands, A is a dimensionless constant, n0 is the density of adsorption sites at the surface, R is the deposition rate in ML/s, v is the fundamental frequency, T is temperature, Esurf is the activation energy for surface mobility, and the power p is a parameter depending on the critical nucleation size i. Under steady-state conditions, the generic expression for 2D islands is p = (i/(i + 2)) and for 3D islands it is p = (i/(i + 2.5)).9,10 An estimate for the average distance between the islands (δ) derived from the scaling law (eq 1) indicates that the average island distance and deposition rate are related by δ µ 1/Rp/2. When the island coarsening takes place, the average distance δ increases with time t. In general, the average distance follows a power law with the exponent m depending on the coarsening mechanism.2430 δðtÞ µ t m

ð2Þ

In several papers, the power m = 1/4 was derived by solving a continuity equation describing a capillary-force-induced island coalescence.29,30 The powers p and m in eqs 1 and 2 may, however, depend on the critical nucleation size and temperature. Such simple scaling laws are used extensively to interpret experimental data. Keeping all other conditions constant, the effect of one process parameter, such as the deposition time (t), deposition rate (R), or temperature (T), can be observed. Such simple dependencies have proven to work well at relatively early

stages of deposition. However, the detailed morphologies that develop at later stages of the process, when the nuclei grow into nanocrystalline islands, cannot be described by the density and/or the average distance of nuclei of alone and demand a thorough investigation of the morphologies as a function of multiple process parameters, such as the deposition time14,31 and deposition rate,25,3234 crystalline orientation,35,36 and material24,32,37 and temperature27,35,36,3840 of the substrate. An understanding of the synergy of these factors and their impact on the morphology is still incomplete and beyond the capabilities of simple kinetic theories and scaling laws. More comprehensive simulation methods are required to realize rational in silico-aided design of nanofabrication processes. Several simulation techniques are employed to model the detailed NC morphologies. The usefulness of a numerical technique depends on the surface area, time scale, and the nature of the properties to be predicted. For example, molecular dynamics (MD) simulations4143 provide detailed information of every atom’s trajectory, kinetic Monte Carlo (KMC) simulations4451 can handle bigger systems while employing a stochastic approach to describe diffusion of deposited atoms, and analytic coarsegrained models5254 can virtually cover unlimited areas at the expense of limitations in the level of detail at the atomic scale . Out of these simulation techniques, KMC provides an advantage of combining power, flexibility, and numerical efficiency. In our previous work,44,45 we reported a KMC simulation of the impact of major PVD process parameters, such as the deposition rate, temperature, and duration on the size, shape, and areal density of two-dimensional projections of growing NCs. At early stages of deposition when the surface coverage is yet incomplete, such a 2D modeling approach is representative of the texture trends 19558

dx.doi.org/10.1021/jp205791t |J. Phys. Chem. C 2011, 115, 19557–19568


ARTICLE

Figure 2. Model for the fcc lattice. The gray rectangles indicate the allowed positions for occupation. (a) Layer 1. (b) Layer 2. (c) Cross-sectional view.

that are unavailable from the scaling laws. However, it does not capture the entire 3D mound morphologies observed experimentally.14,25,3133,55 Multilayer pyramid-like NCs can develop even at relatively low overall surface coverage,12,56 making necessary more detailed 3D modeling studies of the NCs’ self-assembly. In this work, we extend our KMC simulation model introduced previously44,45 to three dimensions and present the results of the simulation, which we compare with pertinent analytic models and experimental reports. In the 3D model, the morphologies can be characterized by the island separation (or density), the lateral size, height, and 3D shape of NCs. One of the objectives is to demonstrate a capability to control the morphology of the growing NCs by tuning the process conditions relevant to PVD processes, such as sputtering or evaporation. In section 2, we describe the KMC model; in section 3, we discuss the results of the simulations, and section 4 summarizes the conclusions.

Deposited atoms can occupy only allowed (gray) sites in the upper layer at the given (x,y) location. When deposited atoms diffuse, they also jump only into an allowed nearest-neighbor (x,y) position, which can, however, be at a different height level h than the initial position. Individual deposition and diffusion events are not correlated, and both occur at randomly selected (x,y) locations. The deposition and diffusion processes are represented by the corresponding numbers of Monte Carlo steps Ndep and Ndiff. For a given period of time t, the number of attempted deposition steps is defined by Ndep = tR, where R is the deposition rate in monolayers (MLs). The number of attempted Monte Carlo diffusion steps is defined as Ndiff = tDadatom/d2, where D is the surface diffusivity of adatoms and d is the average diffusion jump distance adopted equal to the distance between nearest-neighbor atoms in the lattice. The diffusivity can be represented by E ð4Þ Dadatom ¼ D0 exp kT

2. THE KMC MODEL Our 3D simulation of self-assembled NC synthesis employs the solid-on-solid (SOS) model.5759 We represent the substrate with an fcc lattice model structure in a (100) plane. The metal atoms deposited directly onto the substrate (1st layer) are allowed to occupy alternate (x,y) locations, as shown in gray in Figure 2a. In the second layer, the allowed positions change (see Figure 2b,c). The third layer repeats the symmetry of the first layer. In contrast to the previous 2D model,44,45 the present 3D model represents the morphology in a more detailed way by describing the NC’s height, h(x,y). In this 3D fcc lattice, an atom has four nearest-neighbor atoms in the same horizontal layer and four atoms in diagonal locations in the lower and upper layers each. The maximum number of nearest neighbors is, therefore, 12, as is the maximum number of bonds per atom. Binding with atoms located farther than the nearest-neighbor locations is not considered. The simulations were performed for L L = 140 140 sites with periodic boundary conditions. We assume that the distance between nearest-neighbor atoms in the lattice, d, is equal to 0.25 nm, resulting in the size of the simulated region of L = 25 nm. The height of a site is defined by pffiffiffi ð3Þ hðx, yÞ ¼ Nh d= 2

where D0 is the pre-exponential factor, T is temperature, and E is the average activation energy of adatom diffusion for Cu(100), which is adopted equal to 0.46 eV for Cu.60 At each attempted diffusion step, an initial (x,y) position before the jump, “100 , is selected randomly. If the position is occupied, one of its eight (x,y) nearest neighbors is selected randomly and labeled “2”. If the position “2” is available for occupation, the adatoms’ jumping from position “1” to position “2” is simulated employing the Metropolis algorithm with the energy function F + Fs, where F and Fs represent the binding energies of the deposited adatom with other adatoms and the substrate, respectively. The corresponding probability of an adatom’s jump from position “1” to position “2” is given by ðF2 þ Fs2 Þ ðF1 þ Fs1 Þ ð5Þ P12 ¼ min 1, exp kT

where Nh is the number of layers at a given location. In the following discussion, the average height, h, is defined by the sum of maximum heights at all occupied (x,y) sites divided by the total number of the occupied sites. In our SOS model, at each Monte Carlo deposition step, an (x,y) location is selected randomly.

Here, F = yn, where y is the energy of a bond between two adatoms and n is the number of adatomadatom bonds in a given position. We set F = 0 for locations where the number of bonds n is less than a critical number nc, n e nc, where nc is adopted to be 2. In the absence of defects at the surface of the substrate, the adatomsubstrate binding energy Fs can be represented as 4ys for the first layer and is equal to 0 for upper layers. Here, ys is a nominal energy of an adatomsubstrate bond, which may differ from y. The number of bonds for an adatom with its neighbors n depends on the location,6163 and so does the energy function F. The substrate binding energy, Fs, accounts for the possible difference of the binding at the early stages submonolayer of deposition and/or for adatoms located at 19559



ARTICLE

Figure 3. Simulated surface morphologies demonstrating Cu islands’ growth for deposition times of 5, 30, and 60 s at T = 300 K and R = 0.086 ML/s or 0.2 nm/s. (a) The 3D view of the substrate. (b) The 2D view of the same images. (c) The 2D heightheight correlation function and (d) the 1D heightheight correlation function of the above substrates.

the first surface layer. As a result, in the present model, the jumping probability P12 for adatoms is not constant throughout the simulation, but depends on the location. A discussion of a variability of the surface environment of an adatom and of the corresponding changes in the energetic can be found elsewhere.6163 The energy of bonds between adatoms, y, was derived from the cohesive energy, which provides a reasonable approximation for calculating the bond-breaking energy.6466 For Cu, the cohesive energy is equal to 3.52 eV.67 Assuming that the cohesive energy represents the binding of a surface atom with eight nearest neighbors, our estimate of the energy per bond based on the measured cohesive energy for Cu provides approximately 0.44 eV. The labels “1” and “2” in eq 5 denote the energy function in the corresponding locations before and after the attempted jump, respectively. Because the binding energies F and Fs, which describe attractive interatomic binding, are negative values, the probability P given by eq 5 is equal to 1 for jumps resulting in an increasing number of bonds, and is less than 1 otherwise. The entire process of simultaneous depositiondiffusion is repeated for a number of time cycles, typically totaling 1 s. Time in our KMC simulation is scaled proportionally to the number of deposition steps, taking into account the selected deposition rate. For a given deposition rate, the number of deposition steps Ndep is determined by the deposition time, whereas the number of attempted diffusion steps Ndiff depends on both time and

temperature. The described realistic definitions of the diffusion and deposition processes allows us to study the impact of the deposition rate, time, and temperature, as independent factors determining 3D surface morphology during PVD deposition. Varying the substrateadatom binding is also important to understanding how the change in the substrate binding energy affects the morphology of deposited nanostructures.

3. RESULTS AND DISCUSSION 3.1. Evolution of Surface Morphology during Deposition of Nanocrystals. Our first example presents a simulation of

deposition of Cu atoms on a Cu(100) substrate at the rate of R = 0.086 ML/s, which is equivalent to 0.2 nm/s, at room temperature. Figure 3 shows snapshots of the surface morphology after 5, 30, and 60 s of deposition. Row (a) presents 3D images of the morphology, and row (b) shows the corresponding 2D projections. As one could expect,47 the islands nucleate spontaneously at the surface and their size increases with time. At some point, the islands start merging while retaining the pyramidal surface morphology, similarly to the experimentally observed behavior;13,14,26 see also Figure 1. At room temperature, the growth occurs mainly because of arrival of deposited atoms, as well as by the merging of neighbor islands (coalescence). Another mechanism that may contribute to the increase of an island’s size, 19560



ARTICLE

Figure 4. Average island size (S) and interisland distance (δ) as functions of time (t) at the deposition rate of R = 0.086 ML/s and T = 300 K.

primarily at elevated temperatures, is due to the growth of larger islands at the expense of smaller ones, also known as Ostwald ripening.6870 To characterize the surface morphologies, we employed the islands’ density (N) and lateral size (S). We define the latter by the island’s average surface diameter,44 calculated as 4 times the total surface area of islands divided by their total surface perimeter. Furthermore, we have estimated the interisland distance from the heightheight correlation function Cðr, sÞ ¼

∑ðhðx, yÞ h̅ Þðhðx þ r, y þ sÞ h̅ Þ

ð6Þ

Here, h(x,y) is the height of a location of reference defined by eq 3, h(x + r, y + s) is the height at a location {x + r, y + s}, and h is the average height value. Row (c) in Figure 3 shows the examples of C(r,s) functions representing the symmetry of the patterns, and row (d) presents the normalized C(r,0) and C(0,s) dependencies. The first large minima in these dependencies are representative of half of the average distance (δ/2) between islands. In the correlation function shown in Figure 3di, oscillations are visible resulting from the periodic fcc structure of the uncovered substrate. At later stages of deposition when the coverage increases (Figure 3d-ii,d-iii), this oscillatory component vanishes. In Figure 4, the average size and distance of Cu clusters are shown as functions of deposition time. It can be seen that both the island size and the interisland distance increase with time. At early stages of deposition when stable islands form, their lateral size increases rapidly with the arrival of new adatoms. If the deposition continues, stable islands keeps growing both laterally as well as vertically (see Figure 3a), and the evolution in the lateral dimensions somewhat slows down after the initial stage, which can be seen in Figure 4. One could expect that the regimes of slower growth (after approximately 15 s in Figure 4) corresponds the best to the conditions of applicability of the analytical scaling law for the interisland distance (eq 2). This power law dependence with m = 0.28 is also shown in Figure 4. The power law function follows closely the interisland distance δ(t) that we derived from the minima of the correlation functions. However, the average size dependence S(t) is reasonably close to the theoretical prediction only for deposition times less than approximately 60 s, after which the dependence S(t) deviates from both δ(t) and the power law. The reason is that, with time, the islands’ bases start to overlap, which affects the estimates of the average island’s size S, but not the intermound distance δ. Thus, the described dependencies S(t) and δ(t) provide complementary descriptions

of the surface morphology, and only δ(t), but not S(t), can be represented by scaling law (eq 2) with m = 0.28 at advanced stages of the deposition. At room temperature, the morphology determined experimentally for Cu islands grown on Cu(100) can be described by m = 0.2514 and that for Fe/Fe(100) has been matched by m = 0.16.25 At increased temperatures, an Fe/Mg(100) system fits the value of m = 0.25 at 400 K,26 the Ge/Ge(100) morphology can be represented by m = 0.424at 500 K, and Rh/mica corresponds to m = 0.3327 at 700 K. Theoretically, a universal value of the coarsening exponent m = 0.25 has been derived by Mullins30 by analyzing a kinetic equation describing adatom diffusion on a mounded surface. More recently, the value of power exponent m = 0.25 was confirmed by Hunt29 considering a more detailed kinetic model of surface mobility on a mounded landscape with an accounting for the Schwoebel barriers. However, subsequently, Siegert28 has demonstrated that, although in most cases a 0.25 value of the exponent could be expected, it may adopt values between 0 and 0.33 depending on the material and islands’ geometry. Our estimate of m = 0.28 falls into the interval between the generic value of m = 0.25 and the upper theoretic limit of 0.33 and is in a reasonable agreement with the reported experiments. Overall, the predicted morphology of growing NCs is compatible with both the theory as well as with experimentally reported morphologies, both at an early postnucleation stage (see e.g., Figure 1) as well as at later stages when the deposited structures contain multiple layers of atoms. 3.2. Control of NC Synthesis by Varying the Deposition Rate. As the next step, we explore the capability of our model to demonstrate control over the morphology of self-assembled arrays of NCs during PVD synthesis by varying the conditions of deposition. One of the major process parameters in PVD of nanocrystals is the deposition rate.14,25,3234,71,72 Figure 5a shows the 3D images representing the evolution of the morphology with time at various deposition rates, and Figure 5b shows 2D images for the same examples. From Figure 5, it is evident that there are gradual changes in size and density and interisland distance due to the change in deposition time t and rate R. By varying these two parameters, a significant control over the NC fabrication process can be reached. Below, we investigate the effect of deposition rate on the interisland distance and density for two different process conditions. One deposition process occurs at a fixed time of deposition, and the other one employs a fixed amount of deposited material θ (coverage). This approach is in the spirit of our earlier 2D simulation reported in ref 44; however, in contrast to the previous 2D study where the coverage θ was kept relatively low to ensure submonolayer deposition regimes, here, we do not have this constraint and consider larger amounts of deposited material that result in a pronounced mounded morphology. Figure 6 shows the effect of deposition rate R on N and δ for a constant time of t = 30 s and constant coverage for θ = 2.6 ML, by varying deposition rate in the same limits from R = 0.0051 ML/s to R = 0.102 ML/s. In the figure, circles denote the constant time process, triangles denote the constant coverage process, and the dotted lines show power law approximations. Figure 6a shows the dependence of the island density (N) with deposition rate (R). The density of islands is determined by counting their number on a 25 nm 25 nm surface area accounting for the periodic boundary conditions. At constant deposition time and low deposition rate R, the island density N increases with the increase of R, reaching a maximum between 19561



ARTICLE

Figure 5. (a) Effect of deposition rate and time on the 3D surface morphology. Gradual changes of NCs’ size and shape observed due to change in the parameter R (ML/s) and t (s) at T = 300 K, y s = y = 0.39 eV (3D view). (b) Effect of deposition rate and time on the surface morphology. Gradual changes of NCs’ size and shape observed due to change in the parameter R (ML/s) and t (s) at T = 300 K, y s = y = 0.39 eV (2D view).

19562



Figure 6. (a) The average island density and (b) average interisland distance, as functions of the deposition rate R, for a fixed time of t = 30 s and for a fixed coverage of θ = 2.6 ML, at T = 300 K.

0.01 and 0.02 ML/s, then decreasing smoothly with the further increase of R. For R > 0.08 ML/s, the decrease in N becomes slightly sharper. The reason for this behavior might originate from the synergy of several elementary processes; the arrival of new adatoms, resulting in the island nucleation and growth; and further coalescence of growing islands. As a result, the growth of existing stable islands and nucleation of new islands occur simultaneously, and the morphology is determined by these processes. As a result, at low deposition rate R with fixed t, a small number of stable islands arise, leaving a vast surface area islandfree. At higher R, the island density N keeps growing as long as the amount of deposited material θ increases.44 The maximum number of islands is observed at coverage θ ∼ 0.50.7 ML. A further increase of R results mostly in an increase of the islands’ height h, whereas the number of islands N rather starts decreasing by merging of neighboring islands. When deposition occurs at constant coverage θ, the density N increases steadily with increasing R because a facilitated nucleation at higher R provides a larger number of stable nuclei per unit area. In Figure 6a, we also show the dependence of N(R) from the scaling law described by eq 1. The corresponding exponent value p was estimated to be approximately 0.60, which is close to the asymptotic value of 0.55 based on the critical number of bonds adopted in our model (i = 2) with the correction for 3D geometry.9,10 Figure 6b presents the dependence of the average interisland distance δ on the deposition rate R. For a constant deposition time t, the distance δ increases with R, the increase being more pronounced at low deposition rates. This trend is consistent with the corresponding decreasing part of the dependence of the island density N(R) in Figure 6a and is explained by the prevailing coalescence of islands, as discussed earlier. In the case of a constant amount of deposited material θ, with low R, the interisland distance δ decreases with R because higher deposition

ARTICLE

rates produce more stable nuclei per unit area, resulting in a decrease in the interisland distance. The power law describing the average distance δ at fixed θ, as derived from eq 1, reads δ µ (1/(Rp/2)). The dependence in Figure 6b is reasonably described with p equal to approximately 0.54, which matches very closely the asymptotic estimated value of 0.55. Thus, for constant coverage θ, we find that the average density N and distance δ of islands can be described by the scaling law (eq 1). In contrast, the scaling law expressed by eq 1 and its derivatives are inapplicable to the constant t dependencies in Figure 6a,b. From the results presented, it is evident that the variation of deposition rate is a powerful control factor, and a process engineer can determine density, growth, and separation by setting deposition rate conditions during NC fabrication. Although the scaling laws can provide a reasonable approximation for constant-coverage deposition processes, in practice, the fabrication of nanocrystals may depend on a multitude of various factors not all of which can be captured by simple analytical dependencies. An efficient, yet more flexible, KMC simulation technique can help to comprehend the synergism of the various contributing factors. Furthermore, 3D KMC simulations reported here provide important information on the detailed morphology of the surface for better understanding of the trends of deposition, as well as for efficient in silico-aided PVD process design. This process, however, depends on other conditions, such as the temperature and substrate type, whose impact should be accounted for properly. 3.3. Effect of Temperature and Surface Binding Energy. Temperature of the substrate is another major control factor for NC fabrication.18,25,27,3234,36,39,73 From PVD practice, it is known that the islands’ size and density are extremely sensitive to temperature. Despite extensive literature on this subject, few numerical studies address simulation of the NC synthesis in practical temperature regimes. Here, we report 3D KMC simulations of the temperature dependence for the range from 225 to 375 K, which is immediately relevant to PVD applications. In our model, temperature has a direct effect on the surface diffusivity of adatoms and also influences the reversibility of the model; see eqs 4 and 5, respectively. In eq 5, representing the probability of reversible atomic jumps, the temperature dependence comes in proportion to the energy function, which, in turn, is determined by the nominal adatom’s binding energies with other adatoms, y, and with the substrate, ys. The ratio ys/y represents the affinity of NCs to the substrate, a large ys/y corresponding to strong affinity. In this section, we employ our 3D KMC model to investigate the temperature dependence of the surface morphologies for various wetting conditions represented by various ys/y ratios. Figure 7a,b shows the surface morphologies after 30 s of deposition at temperatures varying from 225 to 350 K, for a set of representative substrate binding energy ratios from 0 to 2. It can be seen that the number of islands decreases and the size increases with temperature, because the increased surface mobility at higher temperatures facilitates coarsening of islands both by coalescense and by Ostwald ripening. It is also observed that the change in binding energy has an effect on surface coverage, height, and orientation of the NCs formed. By varying the ratio ys/y, we also investigate the influence of the wetting conditions on the surface morphology at various temperatures. In the limiting case of ys = 0, pillar-like islands are formed with a relatively high aspect ratio. With an increase of the ratio ys/y from 0 to 2, the islands acquire a more regular pyramid-like shape and their bases are more pronouncedly aligned along the [100] 19563



ARTICLE

Figure 7. (a) Effect of temperature and the adatomsubstrate binding on the surface morphology. Gradual changes of NCs' size and orientation observed due to change in temperature T and parameter ys/y at t = 30 s, R = 0.086 ML/s, y = 0.39 eV (3D view). (b) Effect of temperature and the adatomsubstrate binding on the surface morphology. Gradual changes of NCs' size and orientation observed due to change in temperature T and parameter ys/y at t = 30 s, R = 0.086 ML/s, y = 0.39 eV (2D view).

direction. We have shown simulations with the ratio ys/y up to 2; however, for ys/y higher than 3, the influence of the wetting conditions on the size, shape, and orientation of the islands becomes less pronounced. As an example, Figure 8 shows the interisland distance δ and the average height h as functions of the ratio ys/y at T = 300 K. When ys = 0, the interisland distance δ is relatively low because a larger number of higher pillars arise at

this condition. With an increae of ys, the interisland distance increases slightly, which is accompanied by a strong decrease in the islands’ height. The latter, however, stabilizes when the ratio ys/y reaches 1, and the reason might be the merging of atoms with the increase of wetting condition. The temperature dependencies of the islands’ size, height, and interisland distance, for the example, of ys/y = 1, are shown in 19564



Figure 8. Average interisland distance (a) and average island height (b) as functions of the binding energy ratio (ys/y) at T = 300 K, for t = 30 s, R = 0.086 ML/s, y = 0.39 eV.

Figure 9. All of the dependencies are increasing functions of temperature, reflecting the trend of the increase of the NCs' size, which is also seen in Figure 7. Such a temperature dependence is in agreement with reported experiments18,25,27,3234,36,39,73 and results from an increased surface mobility that facilitates the NCs' coarsening at higher temperatures. The trend seen in Figure 9a of the distances δ increasing faster than the estimated lateral size S at high temperatures is also an indirect consequence of the island’s coarsening. To maintain an equilibrium shape, large islands formed at higher temperatures tend to grow in height (Figure 9b), resulting in relatively distant nanocrystals separated by an uncovered surface; see also the high-temperature morphologies in Figure 7a,b. Figure 10 presents the density of islands as a function of temperature for three different values of ys. As one could expect, all the dependencies show that a lesser number of islands per unit area are observed at higher temperatures. The graph also shows that a crossover of two Arrhenius-like regimes with different activation energies exists similar to the predictions of our 2D model47 and experimental results.25,32,74 We attribute this crossover to the increasing impact of reversible jumping at increased temperatures. At lower temperatures, the nucleation and growth of NCs is an irreversible process, whereas at higher temperatures, thermally activated bond-breaking events occur more often and the transition from largely irreversible to reversible (e.g., somewhat slower) island growth takes place. Figure 10 shows the average density of islands as a function of temperature (1000/T) for t = 30 s, R = 0.86 ML/s, and three different values of ys. Symbols represent the numerical results, and lines show the Arrhenius approximations.

ARTICLE

Figure 9. Average island size and interisland distance (a) and average island height (b) as functions of temperature (T), at t = 30 s, R = 0.86 ML/s, ys = y = 0.39 eV.

Figure 10. Average density of islands as a function of temperature (1000/T) for t = 30 s, R = 0.86 ML/s, and three different values of ys. Symbols represent the numerical results, and lines show the approximate slopes.

Figure 10 also illustrates the effect of binding energy on the crossover temperature Tx. We observe that, with ys = 0, the crossover occurs at Tx = 255 K; at ys = 0.5y, Tx = 272 K (not shown); at ys = y, Tx = 278 K; and at ys = 2y, Tx = 285 K. Thus, increasing the ratio ys/y from 0 to 2 shifts the crossover temperature from 255 to 285 K. The observed increase of Tx with ys can be explained by higher substrate binding energies impeding the reversible bond-breaking process, thereby shifting the stronger coarsening regimes toward higher temperatures.

4. SUMMARY (i) We have reported a systematic solid-on-solid KMC simulation of PVD synthesis of self-assembled metallic nanocrystals 19565


The Journal of Physical Chemistry C on a crystalline substrate. The model of adatom mobility comprises reversible jumps between allowed fcc sites on a surface of a three-dimensional landscape. The surface wetting conditions can be varied by changing the energy of binding between adatoms and the substrate. The surface diffusion and deposition, although occurring simultaneously, are represented by independent elementary events. The time of deposition is scaled proportionally to the number of deposition events, taking into account the selected deposition rate. The model has been validated by comparison with benchmark experiments reported in the literature, the analytical scaling laws, and published simulations, including our earlier 2D modeling of nanocrystalline synthesis. (ii) Employing the model, we investigated the capabilities to control the self-assembling surface morphology by varying the major PVD process conditions, such as the duration of deposition, deposition rate, and temperature, in the regimes that are directly applicable to the realistic deposition experiments. Distinct from scaling laws and 2D simulations, which are efficient mostly for early deposition stages, the present 3D KMC simulations allow handling of advanced deposition processes when most of the substrate represents a mounded three-dimensional landscape. This allowed, for example, better understanding of the various ways to control surface morphology by varying the deposition rate, without constraining the analysis by low surface coverage. The reported trends of the temperature dependence of the three-dimensional morphologies in the practically relevant temperature regimes of 225375 K would also be useful for researchers to efficiently optimize PVD processes. (iii) Although the present model considers a single (001) crystalline orientation of an fcc substrate, the possibility to vary the wetting conditions allows predicting changes in the morphology when the binding to the substrate is varied. In this model, columnar growth of round-shaped islands is facilitated by poor wetting conditions, whereas a stronger affinity to substrate promotes more regular pyramidal shapes. Surprisingly, the strength of binding to the substrate also seems to have an influence on the hightemperature coarsening, which is facilitated under poor wetting conditions. (iv) The lattice KMC approach reported in this paper has demonstrated a high level of flexibility to capture a complex interplay of the multitude of various factors on which the fabrication of nanocrystals may depend. We expect that the model and its results reported here will facilitate the PVD process design by replacing trial-anderror experiments with numerical prediction. (v) One of the limitations of the present model, however, is the absence of defects—interstitials, grain boundaries, mismatched substrate lattices, and others—which are encountered in reality. A detailed accounting for lattice defects is beyond the framework of the SOS model adopted in the present study. Further development of the lattice-based PVD modeling on which we have been working includes substrates with point and extended defects as well as an extension to account for the Schwoebel barrier at the step edges. A more comprehensive solution would be provided by an off-lattice approach, where adatoms are not restricted to occupy sites in a specific lattice, but migrate to locations of minimized

ARTICLE

energy. Although the present model has proven to be reasonably accurate to explore the PVD process control parameters, there is ample room for the model refinement, which is becoming increasingly important as processes of nanocrystalline self-assembly progress toward deep nanoscale.

’ AUTHOR INFORMATION Corresponding Author

*E-mail: [email protected] (A.K.B.), [email protected] (S.K.D.), [email protected] (M.S.).

’ REFERENCES (1) Dinh, L. N.; Chase, L. L.; Balooch, M.; Siekhaus, W. J.; Wooten, F. Optical properties of passivated Si nanocrystals and SiOx nanostructures. Phys. Rev. B 1996, 54, 5029–5037. (2) Robinson, I. K.; Pfeiffer, F.; Vartanyants, I. A.; Sun, Y.; Xia, Y. Enhancement of coherent X-ray diffraction from nanocrystals by introduction of X-ray optics. Opt. Express 2003, 11, 2329–2334. (3) Chen, J. H.; Wang, Y. Q.; Yoo, W. J.; Yeo, Y. C.; Samudra, G.; Chan, D. S.; Du, A. Y.; Kwong, D. L. Nonvolatile flash memory device using Ge nanocrystals embedded in HfAlO high-ktunneling and control oxides: Device fabrication and electrical performance. IEEE Trans. Electron Devices 2004, 51, 1840–1848. (4) Lee, J. J.; Kwong, D. L. Metal nanocrystal memory with high-k tunneling barrier for improved data retention. IEEE Trans. Electron Devices 2005, 52, 507–511. (5) Xu, N. S.; Huq, S. E. Novel cold cathode materials and applications. Mater. Sci. Eng., R 2005, 48, 47–189. (6) Tjong, S. C. Nanocrystalline Materials: Their Synthesis-StructureProperty Relationships and Applications; Elsevier Science Ltd.: Amsterdam, 2006. (7) Venables, J. A. Nucleation and growth of thin films: recent progress. Vacuum 1983, 33, 701–705. (8) Venables, J. A.; Spiller, G. D. T.; Hanbucken, M. Nucleation and growth of thin films. Rep. Prog. Phys. 1984, 47, 399. (9) Venables, J. A. Nucleation and growth processes in thin film formation. J. Vac. Sci. Technol., B: Microelectron. Nanometer Struct.— Process., Meas., Phenom. 1986, 4, 870. (10) Venables, J. A. Atomic processes in crystal growth. Surf. Sci. 1994, 299300, 798–817. (11) Ratsch, C.; Venables, J. A. Nucleation theory and the early stages of thin film growth. J. Vac. Sci. Technol., A 2003, 21, S96. (12) Evans, J. W.; Thiel, P. A.; Bartelt, M. C. Morphological evolution during epitaxial thin film growth: Formation of 2D islands and 3D mounds. Surf. Sci. Rep. 2006, 61, 1–128. (13) Biham, O.; Furman, I.; Mehl, H.; Wendelken, J. Diffusion, Nucleation and Growth on Metal Surfaces. Quantum Dots: Fundamentals, Applications, and Frontiers; Springer: Dordrecht, The Netherlands, 2005; pp 5570. (14) Zuo, J. K.; Wendelken, J. F. Evolution of mound morphology in reversible homoepitaxy on Cu(100). Phys. Rev. Lett. 1997, 78, 2791–2794. (15) Venables, J. A. Nucleation growth and pattern formation in heteroepitaxy. Physica A 1997, 239, 35–46. (16) Reichelt, K. Nucleation and growth of thin films. Vacuum 1988, 38, 1083–1099. (17) Stoyanov, S.; Kashchiev, D. Thin film nucleation and growth theories: A confrontation with experiment. Curr. Top. Mater. Sci. 1981, 7, 69141. (18) Ratsch, C.; Zangwill, A.; Smilauer, P.; Vvedensky, D. D. Saturation and scaling of epitaxial island densities. Phys. Rev. Lett. 1994, 72, 3194–3197. (19) Amar, J. G.; Family, F. Critical cluster size: Island morphology and size distribution in submonolayer epitaxial growth. Phys. Rev. Lett. 1995, 74, 2066–2069. 19566


The Journal of Physical Chemistry C (20) Bartelt, M. C.; Evans, J. W. Scaling analysis of diffusionmediated island growth in surface adsorption processes. Phys. Rev. B 1992, 46, 12675–12687. (21) Villain, J.; Pimpinelli, A.; Tang, L.; Wolf, D. Terrace sizes in molecular beam epitaxy. J. Phys. I 1992, 2, 2107–2121. (22) Rao, C. N. R.; Kulkarni, G. U.; Thomas, P. J.; Edwards, P. P. Size-dependent chemistry: Properties of nanocrystals. Chem.—Eur. J. 2002, 8, 28–35. (23) Law, M.; Luther, J. M.; Song, Q.; Hughes, B. K.; Perkins, C. L.; Nozik, A. J. Structural, optical, and electrical properties of PbSe nanocrystal solids treated thermally or with simple amines. J. Am. Chem. Soc. 2008, 130, 5974–5985. (24) Van Nostrand, J. E.; Chey, S. J.; Hasan, M. A.; Cahill, D. G.; Greene, J. Surface morphology during multilayer epitaxial growth of Ge(001). Phys. Rev. Lett. 1995, 74, 1127–1130. (25) Stroscio, J. A.; Pierce, D. T.; Stiles, M. D.; Zangwill, A.; Sander, L. M. Coarsening of unstable surface features during Fe(001) homoepitaxy. Phys. Rev. Lett. 1995, 75, 4246–4249. (26) Th€urmer, K.; Koch, R.; Weber, M.; Rieder, K. H. Dynamic evolution of pyramid structures during growth of epitaxial Fe (001) films. Phys. Rev. Lett. 1995, 75, 1767–1770. (27) Tsui, F.; Wellman, J.; Uher, C.; Clarke, R. Morphology transition and layer-by-layer growth of Rh(111). Phys. Rev. Lett. 1996, 76, 3164–3167. (28) Siegert, M. Coarsening dynamics of crystalline thin films. Phys. Rev. Lett. 1998, 81, 5481–5484. (29) Hunt, A. W.; Orme, C.; Williams, D. R. M.; Orr, B. G.; Sander, L. M. Instabilities in MBE growth. Europhys. Lett. 1994, 27, 611. (30) Mullins, W. W. J. Appl. Phys. 1957, 28, 333. J. Appl. Phys. 1959, 30, 77. (31) Biham, O.; Furman, I.; Karimi, M.; Vidali, G.; Kennett, R.; Zeng, H. Models for diffusion and island growth in metal monolayers. Arxiv preprint cond-mat/9611095; 1996. (32) Zhang, C.-M.; Bartelt, M. C.; Wen, J.-M.; Jenks, C. J.; Evans, J. W.; Thiel, P. A. Submonolayer island formation and the onset of multilayer growth during Ag/Ag(100) homoepitaxy. Surf. Sci. 1998, 406, 178–193. (33) G€unther, S.; Kopatzki, E.; Bartelt, M. C.; Evans, J. W.; Behm, R. J. Anisotropy in nucleation and growth of two-dimensional islands during homoepitaxy on” hex” reconstructed Au(100). Phys. Rev. Lett. 1994, 73, 553–556. (34) Ernst, H. J.; Fabre, F.; Lapujoulade, J. Nucleation and diffusion of Cu adatoms on Cu(100): A helium-atom-beam scattering study. Phys. Rev. B 1992, 46, 1929–1932. (35) Elliott, W. C.; Miceli, P. F.; Tse, T.; Stephens, P. W. Temperature and orientation dependence of kinetic roughening during homoepitaxy: A quantitative X-ray-scattering study of Ag. Phys. Rev. B 1996, 54, 17938–17942. (36) Mo, Y. W.; Kleiner, J.; Webb, M. B.; Lagally, M. G. Activation energy for surface diffusion of Si on Si(001): A scanning-tunnelingmicroscopy study. Phys. Rev. Lett. 1991, 66, 1998–2001. (37) Kopatzki, E.; G€unther, S.; Nichtl-Pecher, W.; Behm, R. Homoepitaxial growth on Ni(100) and its modification by a preadsorbed oxygen adlayer. Surf. Sci. 1993, 284, 154–166. (38) Stroscio, J. A.; Pierce, D. T.; Dragoset, R. A. Homoepitaxial growth of iron and a real space view of reflection-high-energy-electron diffraction. Phys. Rev. Lett. 1993, 70, 3615–3618. (39) Bardotti, L.; Stoldt, C. R.; Jenks, C. J.; Bartelt, M. C.; Evans, J. W.; Thiel, P. A. High-resolution LEED profile analysis and diffusion barrier estimation for submonolayer homoepitaxy of Ag/Ag(100). Phys. Rev. B 1998, 57, 12544–12549. (40) Costantini, G.; Buatier de Mongeot, F.; Boragno, C.; Valbusa, U. Temperature dependent reentrant smooth growth in Ag(001) homoepitaxy. Surf. Sci. 2000, 459, L487–L492. (41) Wolf, D.; Yamakov, V.; Phillpot, S. R.; Mukherjee, A.; Gleiter, H. Deformation of nanocrystalline materials by molecular-dynamics simulation: Relationship to experiments? Acta Mater. 2005, 53, 1–40. (42) Erkoc-, S-. Stability of gold clusters: Molecular-dynamics simulations. Physica E 2000, 8, 210–218.

ARTICLE

(43) Heino, P.; H€akkinen, H.; Kaski, K. Molecular-dynamics study of copper with defects under strain. Phys. Rev. B 1998, 58, 641–652. (44) Bhuiyan, A. K.; Dew, S. K.; Stepanova, M. Kinetic Monte Carlo simulation of metallic nanoislands grown by physical vapor deposition. CiCP 2011, 9, 49. (45) Bhuiyan, A. K.; Stepanova, M.; Dew, S. K. Trends of nanoclusters’ growth by physical vapor deposition studied by atomistic simulation. MRS 2009, 1177, Z09-05. (46) Davis, G. D.; Venables, J. D. In Surface Treatments of Metal Adherends; Chaudhury, M., Pocius, A.V., Eds.; Surfaces, Chemistry and Applications; Elsevier Science B.V.: Amsterdam, 2002; pp 9471008. (47) Gilmer, G. H.; Huang, H.; de la Rubia, T. D.; Dalla Torre, J.; Baumann, F. Lattice Monte Carlo models of thin film deposition. Thin Solid Films 2000, 365, 189–200. (48) Nurminen, L.; Kuronen, A.; Kaski, K. Kinetic Monte Carlo simulation of nucleation on patterned substrates. Phys. Rev. B 2000, 63, 35407. (49) Park, S.; Kim, Y.; Won, T. Kinetic Monte Carlo study on boron diffusion posterior to pre-amorphization implant process. Microelectron. Eng. 2009, 86, 430–433. (50) Wang, L.; Clancy, P. Kinetic Monte Carlo simulation of the growth of polycrystalline Cu films. Surf. Sci. 2001, 473, 25–38. (51) Zhang, P.; Zheng, X.; Wu, S.; Liu, J.; He, D. Kinetic Monte Carlo simulation of Cu thin film growth. Vacuum 2004, 72, 405–410. (52) Friedrich, L. J.; Dew, S. K.; Brett, M.; Smy, T. Thin film microstructure modelling through line-segment simulation. Thin Solid Films 1995, 266, 83–88. (53) Hildebrand, M.; Mikhailov, A. S. Mesoscopic modeling in the kinetic theory of adsorbates. J. Phys. Chem. 1996, 100, 19089–19101. (54) Schulze, T. P.; Smereka, P. Coupling kinetic Monte-Carlo and continuum models with application to epitaxial growth. J. Comput. Phys. 2003, 189, 197–211. (55) Van der Vegt, H. A.; Huisman, W. J.; Howes, P. B.; Vlieg, E. The effect of Sb on the nucleation and growth of Ag on Ag(100). Surf. Sci. 1995, 330, 101–112. (56) Liu, S. J.; Wang, E. G.; Woo, C. H.; Huang, H. Threedimensional SchwoebelEhrlich barrier. J. Comput.-Aided Mater. Des. 2000, 7, 195–201. (57) Oleksy, C. Solid-on-solid model of overlayer-induced faceting. Surf. Sci. 2004, 549, 246–254. (58) Pan, E.; Sun, M.; Chung, P. W.; Zhu, R. Three-dimensional kinetic Monte Carlo simulation of prepatterned quantum-dot island growth. Appl. Phys. Lett. 2007, 91, 193110–193110. (59) Searson, P. C.; Li, R.; Sieradzki, K. Surface diffusion in the solid-on-solid model. Phys. Rev. Lett. 1995, 74, 1395–1398. (60) Pomeroy, J. M.; Brock, J. D. Critical nucleus phase diagram for the Cu(100) surface. Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 73, 245405.1245405.7. (61) Karimi, M.; Tomkowski, T.; Vidali, G.; Biham, O. Diffusion of Cu on Cu surfaces. Phys. Rev. B 1995, 52, 5364. (62) Breeman, M.; Barkema, G. T.; Boerma, D. O. Binding energies and stability of Cu-adatom clusters on Cu(100) and Cu(111). Surf. Sci. 1995, 323, 71–80. (63) Antczak, G.; Ehrlich, G. Surface Diffusion: Metals, Metal Atoms, and Clusters; Cambridge University Press: New York, 2010. (64) Uppenbrink, J.; Wales, D. J. Structure and energetics of model metal clusters. J. Chem. Phys. 1992, 96, 8520. (65) Eckstein, W. Computer Simulation of Ion-Solid Interactions; Springer-Verlag: Berlin, 1991; p 296. (66) The description of the SRIM code can be found at http://www. srim.org/. (67) Kambe, K. Cohesive energy of noble metals. Phys. Rev. 1955, 99, 419–422. (68) Penn, R. L.; Banfield, J. F. Morphology development and crystal growth in nanocrystalline aggregates under hydrothermal conditions: Insights from titania. Geochim. Cosmochim. Acta 1999, 63, 1549–1557. (69) Yang, H. G.; Zeng, H. C. Preparation of hollow anatase TiO2 nanospheres via Ostwald ripening. J. Phys. Chem. B 2004, 108, 3492–3495. 19567



ARTICLE

(70) Leite, E. R.; Giraldi, T. R.; Pontes, F. M.; Longo, E.; Beltran, A.; Andres, J. Crystal growth in colloidal tin oxide nanocrystals induced by coalescence at room temperature. Appl. Phys. Lett. 2003, 83, 1566. (71) Linderoth, T. R.; Mortensen, J. J.; Jacobsen, K. W.; Lægsgaard, E.; Stensgaard, I.; Besenbacher, F. Homoepitaxial growth of Pt on Pt(100)-hex: Effects of strongly anisotropic diffusion and finite island sizes. Phys. Rev. Lett. 1996, 77, 87–90. (72) Swan, A. K.; Shi, Z. P.; Wendelken, J. F.; Zhang, Z. Fluxdependent scaling behavior in Cu(100) submonolayer homoepitaxy. Surf. Sci. 1997, 391, L1205–L1211. (73) Bott, M.; Hohage, M.; Morgenstern, M.; Michely, T.; Comsa, G. New approach for determination of diffusion parameters of adatoms. Phys. Rev. Lett. 1996, 76, 1304–1307. (74) Brune, H.; Bales, G. S.; Jacobsen, J.; Boragno, C.; Kern, K. Measuring surface diffusion from nucleation island densities. Phys. Rev. B 1999, 60, 5991–6006.

19568


Controlled Self-Assembly of Nanocrystalline Arrays Studied by 3D

Recommend Documents