Ind. Eng. Chem. Res. 2002, 41, 1189-1194
1189
MATERIALS AND INTERFACES Influence of Nonuniform Kinetics on Particle Aggregation Bin Lin, Radhakrishna Sureshkumar,* and John L. Kardos Department of Chemical Engineering and Materials Research Laboratory, Washington University, St. Louis, Missouri 63130
Many pathways for the synthesis of novel materials involve aggregation processes that occur in the presence of diffusion barriers and spatially nonuniform kinetics (i.e., the ratio of the timescales of diffusion and reaction depends on the cluster size). In this study, we simulate clusters that grow by diffusion and spatially nonhomogeneous kinetics on two-dimensional square lattices. The influence of nonuniform kinetics on the fractal dimension, particle density distribution, and average diffusion path length is investigated. It is found that if the simulations are performed for a sufficiently long time, the fractal dimension and packing density are not significantly influenced. Specifically, the decrease in the fractal dimension and the increase in the packing density as compared to their values in the case of diffusion-limited aggregation are relatively small. However, the effective diffusion path length prior to attachment is significantly enhanced, especially for smaller cluster sizes. The implications of these observations for particle aggregation processes are discussed. Introduction Diffusion-limited aggregation (DLA) is a fractal growth model exhibiting rich structural complexity.1 DLA was originally introduced by Witten and Sander2,3 as a model for irreversible colloidal aggregation. Since then, the DLA framework has been used in a wide variety of applications, ranging from modeling of river networks, oil recovery, and electrodeposition to string theory.4,5 DLA is only a limiting case of a class of more generic aggregation processes where the balance between diffusion and reaction controls the free particle concentration near a growing surface. DLA occurs when diffusion becomes the rate-limiting step and, consequently, the surface concentration of free particles falls to zero. There has been increasing interest in the growth mechanisms governing the cluster structure in DLA. Several simple generalizations6-8 of the DLA model have been developed to take into account a uniform kinetic constraint by introducing a new parameter, called the sticking or reaction probability Ps. When the diffusing free particle reaches a site that is nearest to the aggregate, it has now only the probability Ps of sticking. For small reaction probability, the aggregate starts to grow as in the Eden model,9 in a quasi-compact fashion, and then, above a characteristic size, the fractal character shows up, and the scaling laws for Ps ) 1 are recovered.10 Meakin11 and Shitamoto et al.12 investigated the deposit profiles in diffusion-controlled deposition on a surface with a homogeneous kinetic constraint. Duffy and Darby13 simulated the deposition structure on a pipe wall from turbulent gas flow. It has been found that decreasing the sticking probability has little effect * Corresponding author. E-mail:
[email protected]. Phone: 1-314-935-6082. Fax: 1-314-935-7211.
on the fractal dimension, although the aggregate structure becomes more compact. In many processes of technological relevance, the timescales of diffusion and reaction are comparable. For instance, Lin et al.14,15 have experimentally investigated the aqueous electropolymerization of pyrrole on carbon fibers. Their experimental data show that the process is initially controlled by the surface reaction kinetics and becomes diffusion limited (i.e., weight gain due to electrocoating increases with the square root of time) only for sufficiently long times. Kane et al.16 developed a reaction-diffusion model for the semiconductor nanocluster growth within polymer films. In their model, nanocluster growth is governed by the ratio of the diffusion timescale τD to the reaction timescale τR. These timescales are dependent on the aggregate size. When the cluster size is small, the reaction timescale is much larger than the diffusion timescale, and the whole process is reaction controlled. When more particles attach to the fractal and the cluster size becomes bigger, τD increases while τR decreases. Consequently, diffusion becomes much slower, and finally, the nanocluster growth stops when τD . τR. These studies indicate the presence of size-dependent kinetics during the cluster formation. Systematic research on the influence of spatially nonhomogeneous kinetics on the cluster structure is lacking. In this work, we have investigated the influence of nonuniform kinetics on the structural characteristics of clusters formed in the presence of diffusion barriers. Stochastic Simulations A Monte Carlo simulation algorithm has been developed on the basis of DLA theory. In the model process, particles moving in space undergo a random walk until
10.1021/ie0104114 CCC: $22.00 © 2002 American Chemical Society Published on Web 02/02/2002
1190
Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002
they reach a live site of the growing fractal. They interact with the live site according to the value of sticking probability Ps. This value is derived from the balance between diffusion and reaction at the point of attachment of the growing fractal. This extended version of DLA with the sticking probability is equivalent to the Laplacian growth model with “the third”,8 or the mixed Robin-type boundary condition, because DLA is isomorphic to the Laplacian growth.17 We consider a first-order surface reaction with a lumped surface kinetic rate ke on the surface of the aggregate corresponding to the mixed boundary condition given by
∂C ∂l
keC ) -D
(1)
where C, D, and l denote the concentration, the diffusion coefficient, and the direction of diffusion, respectively. This equation can be expressed in dimensionless form
∂C* + PsC* ) 0 ∂l*
(1 - Ps)
(2)
where l* and C* denote the dimensionless coordinate l/a and dimensionless concentration C/C0, respectively, where a is a reference length, and C0 is a reference concentration. Here, the sticking probability Ps is interpreted as the probability of reaction once the particle reaches a live site on the aggregate and is defined as
kea Da ) Ps ≡ kea + D 1 + Da
(3)
To reduce statistical errors, a specific number of stochastic simulations should be performed for any given set of parameters using different seeds supplied to the random number generator, and the average values corresponding to these simulations should be used. The stochastic errors can be estimated by means of the χ2 test18 2
(4)
In this study, a more general form will be used to represent the nonhomogeneous Damkohler number. This form is given by
Da ) Krn
on the lattice until it reaches a site adjacent to the seed and becomes part of the growing cluster with a sticking probability Ps. A third particle is then introduced at the outside boundaries and undergoes random walk until it also becomes a part of the growing cluster with a sticking probability Ps. Note that the value of Ps changes with the distance from the central seed particle r according to eqs 3 and 5. Periodic boundary conditions were applied in both lattice directions, and 12 000 particles were used in each simulation. A typical cluster is shown in Figure 1. Results and Discussion
where Da is the Damkohler number defined as Da ≡ kea/D. In many processes, the value of Da is dependent on the distance from the central seed particle r. Specifically, for the growth of metal nanoclusters in polymeric films investigated by Kane et al.,16 this relationship was modeled as
Da ∝ r5
Figure 1. Two-dimensional DLA cluster of 12 000 particles for Ps ) 1.0.
(5)
where K and n are the two controlling parameters in our simulation. Note that when Ps , 1 (Da , 1), the whole process is reaction controlled. In this limiting case, the Brownian free particles visit the aggregate sites many times before sticking on it. On the other hand, as Ps f 1 (Da . 1), the whole process is controlled by the diffusive transport of particles from the bulk to the fractal. This limiting case represents the classical DLA model with the perfectly absorbing boundary. From eqs 3 and 5, it is very clear that, for larger values of n and K, the sticking probability reaches the value of unity (corresponding to DLA) quickly for smaller values of the cluster radius. Simulations were carried out on a two-dimensional square lattice of size Lx × Ly. We start with a single seed particle at the center of the square lattice. A site on the outside boundaries is randomly chosen using a uniform random number generator. A second particle is released from that site and undergoes a random walk
χ (m) )
2 r )k 1 i (ψm-1(ri) - ψm(ri))
∑
kri)1 ψm-1(ri) + ψm(ri)
(6)
where ψm-1(ri) and ψm(ri) represent the average value of the normalized cluster particle distance distribution function ψ(r) at the distance ri measured from the central seed particle corresponding to m - 1 and m stochastic simulations, respectively. Detailed information about the cluster particle distance distribution and density will be discussed later in this paper. From eq 6, it is clear that the larger the value of χ2, the higher the dissimilarity between ψm-1 and ψm. Figure 2 shows the influence of the number simulations m on the values of χ2 with Ps ) 1.0. The statistical errors significantly decrease with increasing m, reaching very small values for m > 10. Therefore, results presented for any given set of parameters in this study will correspond to averages over 15 stochastic simulations performed using different seeds supplied to the random number generator. DLA clusters are among the most widely studied fractal objects. The fractal nature of the DLA clusters has been analyzed by analytical, experimental, and computational methods. Many analytical treatments, based on mean-field theories,19,20 continuum approximations,21 or renormalization-group methods,22-26 have been developed to calculate the fractal dimension.
Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1191
Figure 2. Influence of number of simulations m on the values of χ2 with Ps ) 1.0.
However, almost all of these theoretical attempts do not take into account the influence of sticking probability. By numerical simulations, Ball et al.27 and Meakin6,7 have determined the value of the fractal dimension in up to eight spatial dimensions. The fractal dimension d relates the number of particles N with the size R of the cluster via the relationship N ∼ Rd for large N. Therefore, the values of the fractal dimension can be easily obtained from the double-logarithmic plot of N versus R. In this study, the radius of gyration (Rg) is used to represent the characteristic size of the cluster. Earlier studies6,7,13 on the scaling behavior of DLA clusters have verified that the fractal dimension d is insensitive to the sticking probability over the range 0.1 e Ps e 1.0. However, these studies are based on constant sticking probability throughout the whole process. If we consider nonhomogeneous kinetics during the particle aggregation, the results are different depending on the value of K and n. Table 1 shows the values of the fractal dimension obtained for different values of n and K. The influence of n on fractal dimension is negligible for n > 1 because the whole process becomes diffusion-limited for relatively small values of r. On the other hand, we do observe an influence of K on d when n ) 1. The fractal dimension decreases to 1.50 (≈10% difference as compared to the data for Ps ) 1.0) when K ) 0.01. Therefore, the fractal dimension of DLA clusters is not universal in the sense that it depends on the influence of cluster size on the kinetics. From Figure 1, it is clear that the cluster generated is both highly branched and fractal. The cluster’s fractal structure arises because the faster growing parts of the cluster shield the other parts, which therefore become less accessible to incoming particles.4 In this study, the particle distance distribution function ψ(r), given by the number of aggregate particles with a specific distance r measured from the central seed particle normalized by the total number of simulation particles, is used to characterize the cluster structure. The influence of n and K on ψ is shown in Figures 3 and 4, respectively. We can see (Figure 3) that the fractal structure is insensitive to n for n > 1. However, for n ) 1, there are certain differences between the distributions obtained for different values of K and those obtained for a uniform Ps value of 1.0, especially for small r. The structure of the cluster at the later stages should be
Figure 3. Influence of n on ψ(r): (1) Ps ) 1.0; (2) K ) 0.025, n)2; (3) K ) 0.025, n ) 4.
Figure 4. Influence of K on ψ(r): (1) Ps ) 1.0; (2) K ) 0.025, n ) 1; (3) K ) 0.01, n ) 1.
dependent on that at the earlier stages. Generally, an arriving random walker is far more likely to attach to one of the tips of the cluster than to penetrate deeply into the cluster without first contacting any surface site. With smaller sticking probability, the random walker has more freedom to move around, and the probability of attaching to an inside cluster particle becomes larger. Therefore, smaller sticking probabilities do lead to the formation of a denser structure. When n > 1, the sticking probability reaches the value of 1, corresponding to DLA, very rapidly. Consequently, the cluster structure is not significantly influenced, as shown in Figure 3. Similar conclusions can be drawn from the spatial cluster density profiles shown in Figures 5 and 6. Generally, the cluster density decreases continuously with increasing distance from the central seed particle because of the Brownian motion of the free particles and the fractal characteristics of the cluster. The density function, F(r), is defined as the number of sites occupied by aggregate particles normalized by the total number of lattice sites at a specific r value. Figures 5 and 6 show the influence of n and K on F(r). For n > 1, the value of n has little effect on the fractal structure except for very small values of r (r e 2) (see Figure 5). However, there are certain differences among the density profiles obtained for different values of K when n ) 1 (see Figure 6). To understand the dynamics of cluster formation, we need to examine the time evolution of the aggregation
Table 1. Fractal Dimension d with Nonhomogeneous Kinetic Constraint Ps ) 1.0
K ) 0.025 n)5
K ) 0.025 n)4
K ) 0.025 n)3
K ) 0.025 n)2
K ) 0.2 n)1
K ) 0.1 n)1
K ) 0.05 n)1
K ) 0.025 n)1
K ) 0.01 n)1
1.65 ( 0.05
1.67 ( 0.09
1.64 ( 0.09
1.68 ( 0.09
1.65 ( 0.09
1.62 ( 0.09
1.61 ( 0.09
1.59 ( 0.09
1.56 ( 0.07
1.50 ( 0.05
1192
Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002
Figure 5. Influence of n on F(r): (1) Ps ) 1.0; (2) K ) 0.025, n ) 2; (3) K ) 0.025, n ) 4.
Figure 8. Influence of K on λ(r): (1) Ps ) 1.0; (2) K ) 0.05, n ) 1; (3) K ) 0.025, n ) 1; (4) K ) 0.01, n ) 1.
Figure 6. Influence of K on F(r): (1) Ps ) 1.0; (2) K ) 0.025, n ) 1; (3) K ) 0.01, n ) 1.
Figure 9. Influence of n on L: K ) 0.025.
Figure 7. Influence of n on λ(r): (1) Ps ) 1.0; (2) K ) 0.025, n ) 4; (3) K ) 0.025, n ) 2; (4) K ) 0.025, n ) 1.
process. However, direct estimation of timescales is not feasible from the aforesaid Monte Carlo simulations. In this study, the average diffusion path length λ(r) in lattice units is used to represent the relative “diffusion time” for the incoming particles. This information may be interpreted as a pseudo-time that is indicative of the dynamics of the process. From Figures 7 and 8, we find that the average diffusion path length λ(r) decreases significantly with increasing cluster size. Because of the smaller values of sticking probability and the smaller occupancy fraction in the lattice during the early stages of the cluster formation, the incoming particles need to travel much longer distances before they finally attach to the cluster. As the cluster size increases, the sticking probability becomes larger, and the cluster occupies progressively more lattice sites. Therefore, the average diffusion path length of the incoming particles gradually decreases. Figures 7 and 8 also show the influence of n and K on λ(r). Note that all the simulations are
performed by keeping the lattice size and the final occupancy fraction constant. Hence, the lattice size influences all of the cases uniformly. Consequently, the relative differences seen among the curves in Figures 7 and 8 arise because of the variations in the kinetic parameters (i.e., n and K). It is very clear that the rate of aggregation depends markedly on the values of n and K, and therefore on the sticking probability. Earlier results (Figure 3) have already shown that the cluster structure is insensitive to n for n > 1. However, with decreasing n, the average diffusion path length increases significantly, and the attachment of the incoming particle takes “longer”. This is because, at the same distance from the central seed particle r, the sticking probability is lower with smaller value of n; therefore, it becomes less probable for free particles to attach to the aggregate. During the later stages of the cluster formation, the whole process becomes diffusion-limited, and the influence of n on sticking probability gradually disappears. Consequently, λ(r) obtained for the different n values becomes independent of n. Similar conclusions can be drawn for the influence of K on λ(r), shown in Figure 8. Finally, the area under the λ(r) versus r curve between r ) 2 and r )10, denoted by L, was calculated for different values of K and n. The influence of n and K on L is shown in Figures 9 and 10, respectively. For smaller values of K and n, the sticking probability increases to unity much more slowly with increasing cluster size, leading to larger values of L. This effect becomes less pronounced for larger values of K and n. Conclusions In this work, we have investigated DLA in the presence of spatially nonhomogeneous kinetics. The
Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002 1193
Figure 10. Influence of K on L: n ) 1.
kinetic constraint can be easily and self-consistently incorporated into the classical DLA-based Monte Carlo simulations by using a nonuniform, cluster size-dependent sticking probability model inspired by experimental observations. Earlier studies6,7 have demonstrated that a uniform reaction probability Ps * 1 does not significantly influence the fractal scaling laws of the aggregate. The conclusions reached from this study indicate the following: (a) if sufficiently long time simulations are performed (i.e., if the process timescale is very large), the fractal dimension and density of the cluster are only mildly influenced by nonuniform reaction probability; and (b) the short time dynamics (e.g., path length prior to aggregation) can be significantly influenced by the size-dependence of the reaction rate associated with attachment. Conclusion b further motivates dynamic (e.g., Brownian dynamics) simulations to explore the influence of size-dependent kinetics on aggregate morphology in situations where the process timescale is smaller than the (typically glacial) timescale required to reach kinetics independent states in simulations. Such processes could occur in many situations. Examples include nanoparticle formation in semiconductor films16 and the aggregation and sintering of metal oxide clusters containing a limited number of particles.28 Similarly, surface bound nucleation on irregular surfaces, where the formation of initial nucleation sites is kinetically controlled and further growth occurs via diffusion, presents another interesting scenario. In such cases, the size dispersion of the nucleating regions will depend on the original distribution of the sites (see, for example, ref 29). Similarly, surface-bound electropolymerization discussed by Lin et al.14 could also be subject to dynamic effects arising from size-dependent kinetics. Such effects have not yet been systematically investigated using dynamic simulations. Hence, they present avenues for future research. Acknowledgment B.L. gratefully acknowledges financial support from the Boeing Company. Nomenclature a ) reference length C ) concentration of the free particles C0 ) reference concentration C* ) dimensionless concentration, C/C0 d ) fractal dimension D ) diffusion coefficient of the free particles Da ) Damkohler number, ) kea/D ke ) lumped surface kinetic rate
K ) one of the controlling parameters in the simulation l ) the direction of diffusion l* ) dimensionless coordinate, l/a L ) the area under the λ(r) versus r curve between r ) 2 and r )10 m ) number of simulations n ) one of the controlling parameters in the simulation N ) the number of particles in the cluster Ps ) sticking or reaction probability r ) specific distance measured from the central seed particle R ) characteristic size of the cluster Rg ) radius of gyration ψ(r) ) normalized cluster particle distance distribution function F(r) ) density function λ(r) ) average diffusion path length in lattice units τD ) diffusion timescale τR ) reaction timescale
Literature Cited (1) Vicsek, T. Fractal Growth Phenomena; World Scientific: Singapore, 1992. (2) Witten, T. A.; Sander, L. M. Diffusion-Limited Aggregation, a Kinetic Critical Phenomenon. Phys. Rev. Lett. 1981, 47, 1400. (3) Witten, T. A.; Sander, L. M. Diffusion-Limited Aggregation. Phys. Rev. B: Condens. Matter 1983, 27, 5686. (4) Halsey, T. C. Diffusion-Limited Aggregation: a Model for Pattern Formation. Phys. Today 2000, 53, 36. (5) Sander, L. M. Diffusion-Limited Aggregation: A Kinetic Critical Phenomenon? Contemp. Phys. 2000, 41, 203. (6) Meakin, P. Diffusion-controlled Cluster Formation in Two, Three, and Four Dimensions. Phys. Rev. A: At., Mol., Opt. Phys. 1983, 27, 604. (7) Meakin, P. Diffusion-controlled Cluster Formation in 2-6 Dimensional Space. Phys. Rev. A: At., Mol., Opt. Phys. 1983, 27, 1495. (8) Nagatani, T. Laplacian Growth Phenomena with the Third Boundary Condition: Crossover from Dense Structure to DiffusionLimited Aggregation Fractal. Phys. Rev. A: At., Mol., Opt. Phys. 1989, 40, 7286. (9) Eden, M. Proceedings of the Fourth Berkeley Symposium on Mathematical Statistics and Probability, 1961; Vol. 4, p 223. (10) Jullien, R.; Botet, R. Aggregation and Fractal Aggregates; World Scientific: Singapore, 1987. (11) Meakin, P. Diffusion-Controlled Deposition on Fibers and Surfaces. Phys. Rev. A: At., Mol., Opt. Phys. 1983, 27, 2616. (12) Shitamoto, H.; Okutani, O.; Nagatani, T. Deposit Profiles in Two-Dimensional Electrochemical Deposition. J. Phys. D: Appl. Phys. 1998, 31, 1155. (13) Duffy, G. J.; Darby, M. I. Simulated Structure of Particulate Deposition from a Model of Turbulent Gas Flow in two Dimensions. J. Phys. D: Appl. Phys. 1991, 24, 1665. (14) Lin, B.; Sureshkumar, R.; Kardos, J. L. Electropolymerization of Pyrrole on PAN-based Carbon Fibers: Experimental Observations and a Multiscale Modeling Approach. Chem. Eng. Sci. 2001, 56, 6563. (15) Lin, B.; Sureshkumar, R.; Kardos, J. L. Electropolymerization Coatings on Carbon Fibers: Stochastic Simulations and Experiments. Presented at 2000 AIChE Annual Meeting, Los Angeles, CA, 2000; Paper 216i. (16) Kane, R. S.; Cohen, R. E.; Silbey, R. Semiconductor Nanocluster Growth within Polymer Films. Langmuir 1999, 15, 39. (17) Pietronero, L.; Wiesmann, H. J. Stochastic Model for Dielectric Breakdown. J. Stat. Phys. 1984, 36, 909. (18) Press, W. H.; Teukolsky, S. A.; Vetterling, W. T.; Flannery, B. P. Numerical Recipes in FORTRAN, the Art of Scientific Computing; Cambridge University Press: New York, 1992. (19) Muthukumar, M. Mean-Field Theory for Diffusion-Limited Cluster Formation. Phys. Rev. Lett. 1983, 50, 839. (20) Tokuyama, M.; Kawasaki, K. Fractal Dimensions for Diffusion-Limited Aggregation. Phys. Lett. 1984, 100A, 337.
1194
Ind. Eng. Chem. Res., Vol. 41, No. 5, 2002
(21) Nauenberg, M. Critical Growth Velocity in DiffusionControlled Aggregation. Phys. Rev. B: Condens. Matter 1983, 28, 449. (22) Gould, H.; Family, F.; Stanley, H. E. Kinetics of Formation of Randomly Branched Aggregates: a Renormalization-Group Approach. Phys. Rev. Lett. 1983, 50, 686. (23) Nakanishi, N.; Family, F. Large-Cell Monte Carlo Renormalization of Irreversible Growth Processes. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 32, 3606. (24) Nagatani, T. A Renormalisation Group Approach to the Scaling Structure of Diffusion-Limited Aggregation. J. Phys. A: Gen. Phys. 1987, 20, L381. (25) Pietronero, L.; Erzan, A.; Evertsz, C. Theory of Fractal Growth. Phys. Rev. Lett. 1988, 61, 861. (26) Wang, X. R.; Shapir, Y.; Rubinstein, M. Improved Kinetic Renormalisation Group Approach to Diffusion-Limited Aggregation. J. Phys. A: Gen. Phys. 1989, 22, L507.
(27) Ball, R.; Nauenberg, M.; Witten, T. A. Diffusion-Controlled Aggregation in the Continuum Approximation. Phys. Rev. A: At., Mol., Opt. Phys. 1984, 29, 2017. (28) Yang, G.; Biswas, P. Computer Simulation of the Aggregation and Sintering Restructuring of Fractal-like Clusters Containing Limited Numbers of Primary Particles. J. Colloid Interface Sci. 1999, 211, 142. (29) Fransaer, J. L.; Penner, R. M.; Brownian Dynamics Simulation of the Growth of Metal Nanocrystal Ensembles on Electrode Surfaces from Solution. I. Instantaneous Nucleation and Diffusion-Controlled Growth. J. Phys. Chem. B, 1999, 103, 7643.
Received for review May 9, 2001 Accepted December 12, 2001 IE0104114
