Diffusion-Limited Aggregation: A Model for Pattern Formation
Recent insights from this well-studied model have led to many new applications--from river networks to oil recovery, and from electrodeposition to string theory.

Thomas C. Halsey

Nature confronts us at every turn with patterns--whether the stately spiral shapes of galaxies and hurricanes or the beautiful symmetries of snowflakes and silicon. A host of processes can play a role in forming natural patterns, though they usually involve an interaction between the transport and the thermodynamic properties of the matter and radiation involved.

Typically, convection dominates the transport, in both terrestrial and astrophysical contexts. A classical example is Rayleigh-Bénard convection. The instabilities and patterns generated in a fluid that is convectively transporting heat have implications in contexts as far-flung as laboratory fluid dynamics and solar physics.

In many natural settings, however, convection simply cannot occur. In those cases, diffusion usually dominates the transport. Consider the formation of river networks, frost on glass, or veins of minerals in geologic formations. Similarly, convection plays no role in many patterns in laboratory settings--for example, during ion deposition, electrodeposition, or other solidification processes.

The patterns occurring in this type of system have some general features, which are captured by a number of simple models. The most famous of these models is diffusion-limited aggregation. 1 DLA was originally introduced by Tom Witten and Len Sander as a model for irreversible colloidal aggregation, although they and others quickly realized that the model is very widely applicable. Recent progress in our understanding of DLA has hinged on scaling studies in nonequilibrium statistical physics. Those studies have advanced dramatically in recent years, due in no small part to innovative applications of renormalization group techniques. Yet, many aspects of DLA remain puzzling to specialists.

The basic concept

To understand the basics, consider colloidal particles undergoing Brownian motion in some fluid, and let them adhere irreversibly on contact with one another. Suppose further that the density of the colloidal particles is quite low, so one might imagine that the aggregation process occurs one particle at a time. We are then led to the following model.

Fix a seed particle at the origin of some coordinate system. Now introduce another particle at a large distance from the seed, and let it perform a random walk. Ultimately, that second particle will either escape to infinity or contact the seed, to which it will stick irreversibly. Now introduce a third particle into the system and allow it to walk randomly until it either sticks to the two-particle cluster or escapes to infinity. Clearly, this process can be repeated to an extent limited only by the modeler's patience and ingenuity (the required computational resources grow rapidly with n , the number of particles).

A diffusion aggregation cluster in two dimensions
Figure 1
The clusters generated by this process are 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. An arriving random walker is far more likely to attach to one of the tips of the cluster shown in figure 1a than to penetrate deeply into one of the cluster's "fjords" without first contacting any surface site. Thus the tips tend to screen the fjords, a process that evidently operates on all length scales. Figure 1b shows the "equipotential lines" of walker probability density near the cluster, confirming the unlikelihood of random walkers penetrating the fjords.

The example of Hele-Shaw flow

The preceding model is quite interesting, but its general relevance is not immediately apparent, even for colloidal aggregation at finite concentration. To illustrate the model's generality, let us consider a very different problem: Hele-Shaw fluid flow. 2

In a thin cell, or in a porous medium, a fluid's velocity is proportional to the pressure gradient,


where k is the permeability in a porous medium and m is the viscosity of the fluid. If the fluid is incompressible, then taking the divergence of equation 1 yields the Laplace equation,


Suppose that into such a fluid we inject a second, immiscible fluid of much lower viscosity--the result is Hele-Shaw flow. An example beloved of the oil and gas industry is the injection of water into highly viscous oil in a porous rock (such as sandstone), which is a practical form of secondary oil recovery. Because of its low viscosity, the injected fluid's pressure can be set to a constant. Then the flow of the more viscous fluid is determined by equation 2 with a constant-pressure boundary condition, and its velocity is given by equation 1--which thus also determines the velocity of the interface between the two fluids.

A diffusion aggregation cluster in two dimensions
Figure 2
An experimental realization is displayed in figure 2 . A high-viscosity light-colored hydrophobic fluid (2.5% hexadecyl end-capped polymer) was confined to a space 0.4 mm thick between two glass plates 40 cm across. Water (colored dark) was then injected. The branched structure clearly resembles a smeared-out version of the DLA simulation shown in figure 1. Remarkably, the mathematical descriptions of the two problems are almost identical. For Hele-Shaw flow, the pressure field satisfies the Laplace equation with constant-pressure boundary conditions, and the velocity of the interface between the two liquids is proportional to the gradient of the pressure. For DLA, the probability density of the randomly walking particle satisfies the Laplace equation, with the cluster's surface providing a surface of constant probability density. In this case, the probability of growth (not the growth rate) at the surface is given by the gradient of this probability density. Thus DLA is a stochastic version of the Hele-Shaw problem.

The relation between Hele-Shaw and DLA is even more subtle than this, however. In 1984, Boris Shraiman and David Bensimon analyzed the growth of the surface in the Hele-Shaw problem in two dimensions, and reached the surprising conclusion that the problem is, in a mathematical sense, ill-posed. 3 An arbitrary initial surface will generate singular cusps within a finite time after the initiation of growth, a mathematical reflection of the so-called Mullins-Sekerka instability in solidification. Thus, one must add some other physical effect, such as surface tension, to our model of the Hele-Shaw problem to hold these mathematical singularities at bay. In DLA, by contrast, the finite particle size prevents the appearance of any such singularities.

A diffusion aggregation cluster in two dimensions
Figure 3
In colloidal aggregation, the particles diffuse, while in Hele-Shaw flow, the fluid's pressure diffuses. In each case, the growth of the interface is sufficiently slow that we can use the Laplace equation rather than the diffusion equation to model the diffusing field. This suggests that the Laplacian model might be useful for general pattern formation problems in which diffusive transport controls the growth of a structure. This is indeed the case: DLA, or some variant of DLA, has been used to model phenomena as diverse as electrodeposition, surface poisoning in ion-beam microscopy, and dielectric breakdown. 4 Figure 3 shows a mineralogical example, in which a deposition process on a rock surface has led to beautiful dendritic patterns.

DLA, fractals, and multifractals

DLA clusters are among the most widely known and studied fractal objects. The fractal dimension D connects the number of particles n with the size r of the cluster: n = r D . In two dimensions, one finds D » 1.71, and in three dimensions, D » 2.5. Numerical simulations have determined D in up to eight spatial dimensions, with the result 5 that in high numbers of spatial dimensions d , the cluster fractal dimension D ® d - 1.

However, in two dimensions, where DLA has been most completely studied, its fractal nature is curiously fragile. For example, the fractal dimension is sensitive to the lattice structure of the problem. Thus, if one performs the succession of random walks, and grows the cluster without an underlying lattice, one obtains the aforementioned D = 1.71. However, if one studies precisely the same problem on a square lattice, one finds, 6 for large clusters, that D crosses over to a value of 3/2. One of the few rigorous results on the fractal properties of DLA is the bound D ³ 3/2 in two dimensions, proved by Harry Kesten. 7

In addition, the fractal dimension of DLA appears to depend weakly on the geometry of the simulation. The result D = 1.71 is obtained for radial growth from a seed. However, for growth from a surface, or in a channel, one obtains a result closer to D = 1.67, a small but robust difference from the radial growth case that seems to persist to the asymptotic growth limit. 8 DLA clusters also exhibit "multifractality," a property of the growth probabilities on the surface of the cluster. 9 Consider a cluster of n particles. The i th particle has a probability p i that the next particle to arrive at the cluster will attach to it. The probability measure defined on the surface of a 2D cluster by the { p i } is termed the "harmonic measure," due to its relationship with the theory of analytic functions. The probabilities p i are distributed over a wide range, being relatively large at a cluster's outer tips and quite small deep within the fjords. It is thus natural to examine the scaling of the moments of this probability distribution. We can define a scaling function s ( q ) as an exponent,


The existence of a nontrivial function s ( q ) implies multifractality, which can be interpreted as each particular range of growth probability d p being associated with a different fractal dimensionality. Nevertheless, the cluster as a whole still has a unique fractal dimension--the maximum over the fractal dimensions of all the possible ranges of growth probability.

For deterministic problems, the multifractal scaling function s ( q ) often exists even for negative values of q . In those cases, the sum over probabilities in equation 3 is dominated by the very small values of p i . For a stochastic problem such as DLA, one might be skeptical about the existence of such negative- q scaling behavior, which can be easily disrupted by fluctuations. Several researchers have explored the breakdown of scaling for negative values of q ; in general, the precise manner of the breakdown depends on the details of averaging the summation over the stochastic ensemble of DLA clusters. 10

The multifractal exponents corresponding to the harmonic measure have recently been computed exactly by Bertrand Duplantier, using quantum gravity techniques, for a variety of "equilibrium" fractals in two dimensions, such as percolation or Ising clusters, and Brownian walks. 11 The results agree qualitatively with the scenario envisioned for DLA clusters, including the breakdown of the formalism at sufficiently negative values of q . Alas, there is no indication as yet that these techniques can be extended to nonequilibrium problems.

Scaling laws for DLA

Multifractality is an interesting formal property in its own right, but its special interest for DLA lies in the existence of scaling laws connecting the multifractal properties of the probabilities to the fractal dimension of the cluster. 12 The first, and best established, of these laws was found by Nikolai Makarov. He showed that for any continuous curve in two dimensions, the harmonic measure has an "information dimension" of one. Translated into our notation, this implies that


which is in good agreement with numerical results.

A second scaling relation was proposed by Leonid Turkevich and Harvey Scher. Consider the particle of the cluster that is farthest from the center. One might expect that the cluster radius will grow only if a particle attaches to this "tip" particle; a process for which the next arriving particle will have a probability p tip . Since in this event the maximum radius r max will grow by roughly the particle size a , it follows that d r max /d n ~ p tip a . Given p tip as a function of either r or n , and the supplementary assumption that all radii--including the maximum radius--of the cluster scale in the same way with n , this equation can be integrated to give the dependence of r on n .

Let us suppose that p tip is the maximum over the set of all growth probabilities of the particles in the cluster. Then its scaling can be extracted from the multifractal behavior of the growth probability distribution, connecting the asymptotic behavior of this distribution with the fractal dimension. The result is the Turkevich-Scher scaling relation,


Relaxing our assumptions leads to an inequality, in which the dimension is greater than or equal to the right-hand side of equation 5. But therein lies a puzzle: It is the inequality, not the equality, that is satisfied by numerical results.

An additional scaling relation is the "electrostatic" scaling relation that I proposed. It originates in a formula for the change in the capacitance C of a surface with small changes in the surface geometry. Since the growth of a cluster by the addition of particles results in a succession of relatively small changes in the surface, one can convert this formula into a form relevant for DLA:


In two dimensions, this yields S i p 3 i µ 1/ n . Equivalently, comparing with equation 3, s (3) = 1. In higher dimensions, this argument yields a modified scaling relation connecting D , d , and s (3). That relation agrees with numerical results.

Theoretical approaches to DLA

Naturally, the richness of DLA has attracted a number of theoretical attempts at a comprehensive analysis. Challenges and puzzles, however, abound. One difficulty facing all such attempts has been the absence of an easily identifiable small parameter that would allow a perturbation analysis. DLA seems to yield fractal structures in which fluctuations are important up to arbitrarily high spatial dimensions; there is no upper critical dimension, above which mean-field theory would be valid. In fact, mean-field theory for DLA predicts D = d - 1, which only appears to be true in the limit of infinite spatial dimensionality. Also, as a nonequilibrium problem, DLA has no obvious relationship to the class of problems--mostly related to equilibrium statistical mechanics--that can be solved in two dimensions by conformal field-theory techniques.

The self-similarity of DLA clusters suggests that their structure might be determined by a renormalization group approach. Several proposals for applying real-space renormalization methods to DLA have indeed been made. Probably the most sophisticated and successful has been the "fixed scale transformation" of Luciano Pietronero and his coworkers. 13 Although not a real-space renormalization group in the classical sense, it is based on the enumeration of real-space configurations, and uses a transformation between scales. It gives good results for the fractal and multifractal properties of DLA (and a number of other statistical physics problems) in two dimensions. It also shares with real-space renormalization groups the lack of a small perturbative parameter.

My coworkers and I have taken an entirely different approach. 14 A noticeable feature of DLA is the way that branches screen one another simultaneously on a variety of length scales. In two dimensions without a lattice, DLA typically has four or five large branches, which are more or less stable. At smaller length scales, however, branches compete in a never-ending vicious cycle of precarious survival. In fact, for any two neighbors among these smaller branches, at most one will survive as the cluster grows. The death of branches as they are screened by their neighbors is balanced by the creation of new branches via microscopic tip-splitting processes.

This picture of DLA growth led to the "branched growth model," in which the competition--on all length scales--between branches is represented as a dynamical system. The overall cluster dynamics is then represented as a large family of coupled dynamical systems running simultaneously.

This approach allows approximate but quite detailed solutions for the cluster dynamics and fractal properties in all dimensions. Results for D from the branched growth model are in excellent agreement with numerical results, especially in high dimensions. This approach also allows one to compute multifractal properties; those results also agree with simulations. Finally, this approach is especially well suited for computing the topological self-similarity of the clusters (see the box on page 39).

The Hastings-Levitov approach

Recent work on DLA has been dominated by a new formulation of the problem in two dimensions, due to Matthew Hastings and Leonid Levitov. 15 Although it has always been known that the Hele-Shaw problem in two dimensions has a natural conformal representation, Hastings and Levitov were the first to generate an elegant representation of DLA growth as a problem in iterated conformal maps. Their formulation has revived interest in the DLA problem, by making available the powerful tools of analytic function theory.

Consider a cluster of n particles. The Riemann mapping theorem assures us that there exists a conformal map, w = F n ( z ), that maps a unit circle in the complex z -plane onto the surface of the cluster in the physical w -plane. If the exterior of the unit circle is mapped onto the cluster exterior, then it follows that F n is an analytic function in the exterior of the unit circle. Conformal mapping then tells us that the angular distance between two points on the circle circumference in z -space is proportional to the total growth probability along the arc connected by the images of those two points in the physical space.

A diffusion aggregation cluster in two dimensions
Figure 4
Hastings and Levitov gave a simple algorithm for the construction of the function F n ( z ) corresponding to a given cluster. Suppose that a function, f l , q ( z ), giving a "bump" on the unit circle, corresponds to the attachment of one particle of size l at an angular position q in the z -plane (see figure 4 ). Then if F n is the map for an n -particle cluster, the map for an ( n + 1)-particle cluster, where the last particle is added at the image of the angular position q , is given by F n ( f l , q ( z )). Iteration then allows the determination of the cluster map from the "one-particle" maps f . This is a concrete realization of what mathematicians refer to as a "stochastic Loewner process." Curiously, such processes have recently been used in the rigorous proof of some of Duplantier's results for multifractal scaling.

Hastings used the conformal representation of DLA growth to perform a momentum-space renormalization group calculation for the DLA dimension in two dimensions. Although the result, D = 1.7, was highly accurate, the calculation suffered from the ever-present defect of perturbative approaches to DLA: It was based on a small parameter that wasn't small.

The Hastings-Levitov algorithm lets us reproduce not only DLA, but also more general models. If the growth positions, q , in z -space (the pre-image of the physical space) are chosen randomly, we get DLA. But other choices are possible. An interesting choice is q n = 2 pW n for the angle q n of the n th particle, with W a constant. Although perhaps unphysical, this fully deterministic choice allows one to explore the significance of randomness in the DLA model. Benny Davidovitch and his colleagues 16 have shown that for irrational values of the parameter W , the Hastings-Levitov model leads to branched structures qualitatively similar to DLA, but with significantly higher values of the fractal dimension D , as shown in figure 5 . Scaling functions for the Hastings-Levitov model were computed, both for the stochastic (DLA) case and for the quasiperiodic ( W irrational) case; the functions gave numerically accurate results for the dimensions of both types of cluster in two dimensions.

DLA, string theory, and beyond

A diffusion aggregation cluster in two dimensions
Figure 5
The most surprising recent development suggests a possible relationship between Hele-Shaw growth and string theory. 17 The starting point for this development is the remarkable fact that Hele-Shaw growth conserves the harmonic moments of the exterior domains. Those moments are defined by


where the integral is over the exterior of the growing structure, z = x + iy is the ordinary complex variable, and divergences in the integral are suitably regularized. Of course, C 0 varies as the Hele-Shaw pattern grows, but all of the other C k are fixed during the growth. Thus, the problem of determining the patterns created by Hele-Shaw growth is equivalent to determining the families of curves with different values of C 0 but fixed values of the other C k .

This problem, in turn, can be related to a set of equations known as the "integrable Toda hierarchy," which also appear in 2D quantum gravity, and hence in string theory. In this relation, the parameters C k become the degrees of freedom of this integrable hierarchy. Furthermore, it is known that a particular solution of the Toda hierarchy is related to the statistical mechanics of Hermitian N ´ N matrices (which, in the large- N limit, is also believed to reproduce the scaling behavior of 2D quantum gravity). It is precisely in the N ® ¥ limit that the Toda hierarchy maps exactly onto the pure Hele-Shaw problem; this suggests a strong, yet still obscure mathematical relationship between the latter problem and string theory.

Next year marks the 20th anniversary of the Witten-Sander model, which opened the door to the wonderful physics of diffusion-limited aggregation, and revived interest in the classical problem of Hele-Shaw growth. Those beautiful structures seemed, at the outset, likely to be understood by then-conventional techniques. Instead, there remains a certain amount of mystery. Our understanding of the phenomenology of DLA has certainly become quite sophisticated, and the new techniques of Pietronero, Hastings and Levitov, and others have afforded new insights. In addition, it appears that there are deep connections between Hele-Shaw growth--and thus DLA--and 2D quantum gravity. Such newly discovered connections to other problems of theoretical physics suggests that the next 20 years are liable to be full of surprises.

1 . T. A. Witten Jr, L. M. Sander, Phys. Rev. Lett. 47 , 1400 (1981). P. Meakin, Phys. Rev. A 27 , 1495 (1983). T. Vicsek, Fractal Growth Phenomena , World Scientific, Singapore (1989).

2 . P. G. Saffman, G. I. Taylor, Proc. R. Soc. London, Ser. A 245 , 2312 (1958).

3 . B. Shraiman, D. Bensimon, Phys. Rev. A 30 , 2840 (1984).

4 . R. Brady, R. C. Ball, Nature 309 , 225 (1984). L. Niemeyer, L. Pietronero, H. J. Wiesmann, Phys. Rev. Lett. 52 , 1033 (1984). J. Nittmann, G. Daccord, H. E. Stanley, Nature 314 , 141 (1985). See also J. M. Garcia-Ruiz et al ., eds., Growth Patterns in Physics and Biology , Plenum, New York (1993).

5 . R. Ball, M. Nauenberg, T. A. Witten, Phys. Rev. A 29 , 2017 (1984).

6 . R. C. Ball, R. M. Brady, G. Rossi, B. R. Thompson, Phys. Rev. Lett. 55 , 1406 (1985). P. Meakin, R. C. Ball, P. Ramanlal, L. M. Sander, Phys. Rev. A 35 , 5233 (1987).

7 . H. Kesten, Stoch. Proc. Appl. 25 , 165 (1987).

8 . B. B. Mandelbrot, A. Vespignani, H. Kaufman, Europhys. Lett. 32 , 199 (1995).

9 . T. C. Halsey, P. Meakin, I. Procaccia, Phys. Rev. Lett. 56 , 854 (1986). C. Amitrano, A. Coniglio, F. di Liberto, Phys. Rev. Lett. 57 , 1016 (1986).

10 . M. E. Cates, T. A. Witten, Phys. Rev. A 35 , 1809 (1987). T. C. Halsey, in Fractals: Physical Origin and Properties , L. Pietronero, ed., Plenum Press, London (1989). B. B. Mandelbrot, C. J. G. Evertsz, Nature 348 , 143 (1990). S. Schwarzer, J. Lee, S. Havlin, H. E. Stanley, P. Meakin, Phys. Rev. A 43 , 1134 (1991). R. C. Ball, R. Blumenfeld, Phys. Rev. A 44 , 828 (1991).

11 . B. Duplantier, Phys. Rev. Lett. 82 , 880 and 3940 (1999); 84 , 1363 (2000).

12 . N. G. Makarov, Proc. London Math. Soc. 51 , 369 (1985). L. Turkevich, H. Scher, Phys. Rev. Lett. 55 , 1026 (1985). T. C. Halsey, Phys. Rev. Lett. 59 , 2067 (1987).

13 . L. Pietronero, A. Erzan, C. J. G. Evertsz, Phys. Rev. Lett. 61 , 861 (1988). A. Erzan, L. Pietronero, A. Vespignani, Rev. Mod. Phys. 67 , 545 (1995).

14 . T. C. Halsey, M. Leibig, Phys. Rev. A 46 , 7793 (1992). T. C. Halsey, Phys. Rev. Lett. 72 , 1228 (1994). T. C. Halsey, B. Duplantier, K. Honda, Phys. Rev. Lett. 78 , 1719 (1997). T. C. Halsey, Europhys. Lett. 39 , 43 (1997).

15 . M. B. Hastings, Phys. Rev. E 55 , 135 (1997). M. B. Hastings, L. S. Levitov, Physica D 115 , 244 (1998). B. Davidovitch, H. G. E. Hentschel, Z. Olami, I. Procaccia, L. M. Sander, E. Somfai, Phys. Rev. E 59 , 1368 (1999).

16 . B. Davidovitch, M. J. Feigenbaum, H. G. E. Hentschel, I. Procaccia, Phys. Rev. E 62 , 1706 (2000). B. Davidovitch, I. Procaccia, Phys. Rev. Lett. (in press). B. Davidovitch, A. Levermann, I. Procaccia, Phys. Rev. E (in press).

17 . M. Mineev-Weinstein, P. B. Wiegmann, A. Zabrodin, Phys. Rev. Lett. 84 , 5106 (2000).

18 . R. E. Horton, Bull. Geol. Soc. Am. 56 , 275 (1945). A. N. Strahler, Trans. Am. Geophys. Union 38 , 913 (1957). E. L. Hinrichsen, K. J. Måløy, J. Feder, T. Jøssang, J. Phys. A 22 , L271 (1989). J. Vannimenus, X. G. Viennot, J. Stat. Phys. 54 , 1529 (1989).

Thomas Halsey is a senior staff physicist at ExxonMobil Research and Engineering Co in Annandale, New Jersey.