Fractals And Their Applications _ Liu1986.pdf

SOLID STATE PHYSICS, VOLUME

39

Fractals and Their Applications in Condensed Matter Physics S. H. LIU Solid State Division Oak Ridge National Laboratory Oak Ridge, Tennessee 37831

..................................................

11.

Ill.

IV.

V.

V1.

ground ................................. 2. Fractals in Condensed Matter Physics.. ................................ Measurements of Fractal Dimension. ....................................... 3. Real-Space Measurements .................. .. 4. Small-Angle Light or Neutron Scattering ............................... Growth of Fractal Aggregates ............................................. 5. Computer Simulation of Growth . ............................. 6. Kinetic Theory of Aggregation ........................................ 7. Kinetic Theory of Dendritic Growth ..... .................... 8. Experimental Studies. ................................................ Diffusion on Fractal Networks. .......................... ... 9. Anomalous Diffusion. ........... ............................ 10. Diffusion on the Sierpinski Gasket ..................................... 11. The ac Response of Rough Interfac .... Elastic Properties of Fractal Networks ...................................... 12. Critical Elastic Threshold of Percolation Networks. ...................... 13. Spectral Dimension of Elastic Vibrations . ................ Magnetic Ordering on Fractal Networks .................................... 14. Magnetic Ordering of Ising Systems. ............. ..........

207 208 212 219 219 222 227 228 235 239 24 1 244 245 249 253 251 251 263 261 261

1. Introduction

It is often said that a review article in modern physics is out of date as soon as it is written. In a new and fast growing field such as fractals, one runs the risk of writing a review which is obsolete even before it is written. There are a multitude of reasons for this phenomenon. (1) The idea of fractal touches upon too many areas of natural science for one reviewer to be aware of all the new developments. (2) It takes six months to a year before a new research result appears in print, and before it does the reviewer often does not know of its 207 Copyright @ 1986 by Academic Press. Inc. All rights of reproduction in any form reserved.

208

S. H. LIU

existence. (3) Even a knowledgeable reviewer needs time to digest a piece of work and to put it into perspective. Consequently, the review may not do justice to some important works. (4)New advances will take place after the completion of the review and before its appearance in print. Faced with these severe constraints, the best this reviewer can hope to accomplish is to limit the scope of the article to condensed matter physics and to concentrate on the physical principles and some mathematical techniques. The present article grows out of a series of lectures the author gave at the Oak Ridge National Laboratory to a group of scientists from many disciplines who are interested in applying the idea of fractals to their own fields. The purpose of the course was to bridge the gap between Mandelbrot’s books’.’ and research papers. Experts in the field may find this review lacking in sophistication and depth, but I hope that beginners will find it helpful and not intimidating. Since there is no uniformity of notation in the literature, I have chosen to simplify the mathematical symbols by denoting all dimensions related to fractals by d with appropriate subscripts. The spatial dimension is designated by D. The reader is forewarned about a possible confusion, that in the percolation theory one often uses d for spatial dimension and D for fractal dimension. 1 . MATHEMATICAL BACKGROUND A century ago pure mathematics went through a period of crisis. The issues concern the properties of functions of real variables in general, and in particular the question whether all continuous functions are differentiable. Efforts by mathematicians to construct counter examples, i.e., functions that are continuous everywhere but not differentiable anywhere, resulted in a plethora of “monster curves.” These events hailed the birth of the field of fractals. An example of the monster curve is the Koch snowflake, shown in Fig. 1. One begins with an equilateral triangle. Next, one smaller equilateral triangle is added to the middle of each side to obtain the sixfold symmetrical star of the second stage. In the third stage 12 still smaller triangles are erected, one on every side of the star. The Koch curve is the boundary of the complicated shape when the procedure is repeated indefinitely. The curve is continuous but not differentiable because it makes an infinite number of zigzags between any two points. The length between any two points on the curve is infinity, yet the area bounded by the curve is finite. Other examples of monster curves are found in Mandelbrot’s books.

’ B. B. Mandelbrot, “Fractals:Forms, Chance and Dimension.” Freeman, San Francisco, 1977. B. B. Mandelbrot, “The Fractal Geometry of Nature.” Freeman, San Francisco, 1983.

FRACTALS A N D CONDENSED MATTER PHYSICS

209

FIG. 1. The first four steps in the construction of a Koch snowflake. The boundary of the complicated shape, the Koch curve, has the fractal dimension of In4/ln 3 = 1.26.

Hausdorff in 1919 suggested a way to generalize the notion of dimension, thereby putting the monster curves in a class of their own. His idea is based on scaling, which means measuring the same object with different units of measurement. Suppose we measure a line segment of length L with a measuring stick of length 1. We do this by finding out how many times the measuring stick fits into the segment. The result of the measurement is the ratio Y(1)= L/l.

(1.1)

Then we measure the same segment with another measuring stick of length 1/N. The result of the new measurement is Y ( I / N )= L / ( l / N )= N ' Y ( 1 ) .

(1.2) We can do the same thing with area measurements. If a square of L x L is measured with a 1 x 1 measuring unit, the result is the ratio

&(I)

= L2/12.

(1.3)

If a ( l / N ) x ( l / N )measuring unit is used, the result is

d ( l / N ) = N%d(l). (1.4) The reader can generalize this consideration to the measurement of volumes

210

S. H. LIU

and obtain V ( l / N )= N 3 V ( l ) . (1.5) In every case the dimension appears in the ratio of the results as the exponent of the ratio of the length scales. Applying the same consideration to the Koch snowflake, Hausdorff found it necessary to add the condition that one must not count any detail smaller than the unit of measurement. Physical scientists may understand the unit as the limit of resolution, namely our ability to discern fine details. Suppose at the beginning stage the Koch curve in Fig. 1 is 1 cm on a side. With a resolution of 1 cm we see the curve as a triangle with three equal sides. If the resolution is improved to f an,we begin to see 12 segments as in the second stage. Every time the unit of measurement is reduced by a factor of 3, the number of visible segments increases4 times. Following the rule of the last paragraph, we equate 4 to 3d to obtain the dimension of the Koch snowflake d = In 4/ln 3 = 1.26. Therefore, the strange properties of the snowflake come from the fact that it is not a one-dimensional object. The word fractal was coined by Mandelbrot as a generic name for objects such as the Koch snowflake which possess fractional Hausdorff dimensions. Fractals need not be monster curves. A fractal network called the Sierpinski gasket is constructed according to the rules illustrated in Fig. 2. At every step three equilateral triangles are stacked together to form a larger equilateral triangle. The linear size of the object increases by a factor of 2. The reader may imagine the sequence of steps as the results of observing the same object under increasing magnification, by a factor of 2 at a time, so that more and more internal structure becomes visible. One can calculate the dimension of the object by counting the total length of bonds between vertices. For instance, with the resolution equal to the length of each side of the outer periphery, the total length of the bonds is 3. Improvement of the resolution twofold enables one to discern a total bond length of 9. Therefore, the dimension of the object is d = In 3/ln 2 = 1.58. The Sierpinski gasket is a favorite object of study for theoretical physicists because it is a nontrivial model for a noncrystalline solid, with an atom at

0

1

2

FIG.2. The first three steps in the construction of a Sierpinski gasket. This geometrical object has the fractal dimension of In 3/ln 2 = 1.58.

FRACTALS AND CONDENSED MATTER PHYSICS

21 1

FIG.3. The three-dimensional extension of the Sierpinski gasket, the skewed web, is constructed by the iterative process in this illustration. The object has the fractal dimension of In 4/In 2 = 2.

every vertex and a bond between every pair of nearest neighbors. The dimension of the network can also be deduced from the number of atoms in the sequence of steps in Fig. 2. The numbers are 3, 6, 15, 42,. . . ., with the general formula N. = (3" + 3)/2. For sufficiently large n the number of atoms increases by a factor of 3 every time the linear size doubles. Hence, the dimension is In3/ln2 as obtained earlier. The Hausdorff dimension is a geometrical property of the object, and should be independent of the way it is calculated. The rules of constructing the Sierpinski gasket may be generalized to three dimensions to obtain the skewed web shown in Fig. 3. Four identical tetrahedra are joined together to form a larger tetrahedron at every step. The dimension of this object is d = In 4/ln 2 = 2. This example illustrates that the dimension of a fractal is not necessarily nonintegral. It is usually different from the topological dimension of the object. Fractals may also have dimensions less than 1. An example, called the Cantor bar, is shown in Fig. 4. In the beginning there is a solid bar. In the next step the middle portion of the bar is removed to obtain two shorter bars.

~~

FIG.4. The first five steps in the construction of a Cantor bar. This is an example of a fractal object with a dimension less than one.

212

S. H. LIU

Subsequently the shorter bars are subjected to the same subdivision and the procedure is carried ad injnitum. If at every step the shorter bar is l / a times the length of the original bar, with a > 2, the dimension of the infinitely fragmented object is d = In2/lna < 1. The number of subdivisions may be generalized to N , with N < a. The dimension of this general Cantor bar is d = InN/lna 1. Fractals have the important invariance property called self-similarity under scale transformations, or self-similarity for short. Take the Koch snowflake for example. The last entry in Fig. 1 consists of six prickly protrusions. If a protrusion is magnified threefold, the result is identical to a large portion of the third entry. The same property is seen in the Sierpinski gasket and the Cantor bar. Of course, self-similarity is the necessary and sufficient condition that the Hausdorff dimension can be uniquely defined. In later Parts we will show how the self-similar property can be used to advantage in discussing the physical properties of fractal systems.

-=

2. FRACTALS IN CONDENSED MATTER PHYSICS This subsection is a compendium of those fractals which are most frequently discussed in the literature of condensed matter physics. The more familiar examples of fractals, such as the coastlines, the terrain of mountain ranges, the pock-marked surface of the moon, the distribution of matter in the universe, etc., will not be discussed. a. Path of Random Walk The random walk follows a zigzag path which is much longer than the distance covered. We will show that the path is a fractal and determine its fractal dimension. Consider a one-dimensional random walk. Starting from the origin, the walker may take a step either to the right or to the left by one unit of length. After N steps there are 2N possible final positions. The probability that the walker ends at m steps to the right or left of the origin is P(m,N) = 2-"!/[+(N

- m)]![+(N

+ m)]!.

(2.1)

For sufficiently large N and m, P(m, N) approaches the Gaussian distribution P(m, N) = ( 2 ~ )'/'- exp( - m2/2N). If we increase N by a factor of 4, the length scale of the Gaussian distribution increases by a factor of 2. Since N is the length of the total path, we find that the path length is proportional to the square of the length scale. In other words, the path has the Hausdorff dimension d = 2. The same result holds for

213

FRACTALS AND CONDENSED MATTER PHYSICS

random walks in any spatial dimension. The fact that d = 2 is an alternative way of stating the well known, that the net distance traveled by the random walker is proportional to the square root of the number of steps taken, irrespective of the dimension of the space in which the walk takes place.

b. Self-Avoiding Random Walk If in the random walk process the walker avoids any site he has visited previously, the path of his journey will be different from that of a pure random walk. Clearly the one-dimensional self-avoiding random-walk (SAW) problem is trivial because the walker can only proceed in one direction. In two and three dimensions the problem is not simple, and an enormous amount of literature has been devoted to it.’ The problem is important because it concerns the structure geometry of long-chain molecules, such as polymers. When these macromolecules are dissolved in solvents, the molecules curl up into nearly spherical shapes. The viscosity of the solution, which is determined by the average radius of the curled up molecules, is found to be proportional to a nonintegral power of the molecular weight, which is proportional to the total length? The exponent lies between 0.6 and 0.8 for a large number of systems, compared with the theoretical number 0.5 if the molecular shape is assumed to be that of the path of random walk. This indicates that the molecules have a more open structure. Physically two atoms on a chain can not occupy the same space, and for this reason the molecular shape is expected to resemble the path of SAW. Similar to the random-walk problem the analysis becomes manageable if the walk is restricted on a grid. For small numbers of steps N the number of SAW paths and the mean-square distance traveled ( I ; ) have been enumerated.’~~ Fisher and Hiley studied the ratio’ For large N the ratio is expected to approach the limit 1. The power-law relation between the molecular radius and the molecular weight can be written as ( r 2 ) ’ / 2 oc N”.

(2.4)

It follows that 2v = lim N ( p N - 1). N+ m

’M. E. Fisher and B. J. Hiley, J. Chem. Phys. 34, 1253 (1961), and references cited therein. ’

P. J. Flory, “Principlesof Polymer Chemistry,” Ch. 14. Cornell Univ. Press, Ithaca, New York, 1969. M. E. Fisher and M. F. Sykes, Phys. Rev. 114,45 (1959).

214

S. H. LIU

The limit is determined numerically by plotting the sequence 2vN = N(pN- 1) versus 1/N and extrapolating to 1/N 0. This gives v = 0.6 k 0.01 for a simple cubic grid. The SAW problem can be mapped onto a phase-transition problem and solved by the renormalization-group techniq~e.~.’ De Gennes showed by using the scaling argument that the exponent v in Eq. (2.4)is the same as that which appears in the temperature dependence of the correlation length.6 The value of v calculated by using the E = 4 - d expansion method of Wilson and Fisher’ is in agreement with that calculated by the numerical extrapolation method. The SAW path in two-dimensional space is also of interest because it mimics the molecular configuration of polymer molecules adsorbed on a surface. Domb extended the Fisher and Hiley calculation to two dimensions and found v = 0.75.’ The Hausdorff dimension of the SAW path is given by d = l/v = 1.67 in three dimensions and 1.33 in two dimensions. Stapleton and his co-workers have measured the fractal dimension of the backbone of a large number of proteins from their structure data and found that d lies between 1.2 and 1.8.”*” c. Percolation Clusters

The percolation problem arises from the observation that a mixture of metallic particles and insolating powder is nonconducting until about 30% of the volume is metallic, and at this point a rapid rise in conductivity is obtained. Recent interest in the process is rooted in the suggestion that percolation is implicated in the metal-semiconductor transition seen in some disordered solids.” The mathematical model for percolation is also based on a lattice. A set of bonds between nearest-neighbor sites are taken out at random, and a fluid is prevented from flowing from one site to another if the bond between them is missing. The percolation threshold is reached when the fluid can no longer flow through a large grid. This model is called the bond-percolation model. In the site-percolation model a random set of sites are taken out of the lattice P. G. De Gennes, Phys. Lett. 38A, 339 (1972).

’J. C. Le Guillou and J. Zinn-Justin,Phys. Rev. Lett. 39,95 (1977). * K. G. Wilson and M. E. Fisher, Phys. Rev. Lett. 28,240 (1972).

C. Domb, J. Chem. Phys. 38,2957 (1963). H. J. Stapleton,J. P. Allen, C. P. Flynn, D. G. Stinson, and S. R. Kurtz, Phys. Rev. Lett. 45,1456 (1980). J. T. Colvin, G. C. Wagner, J. P. Allen, and H. J. Stapleton,pers. comm. (1985). ” S. Kirkpatrick, in “Ill-Condensed Matter” (R. Balian, R. Maynard, and G. Toulouse, eds.), p. 324. North-Holland Publ., Amsterdam, 1978. lo

FRACTALS AND CONDENSED MATTER PHYSICS

215

such that the fluid cannot flow through the bonds connected to a missing site. Extensive literature exists on the percolation problem, and the reviews by Kirkpatrick l 2 and Stauffer13 should be consulted both as comprehensi /e summaries of the field and as exhaustive sources of references up to about 1979. Consider the site-percolation problem. Near the percolation threshold the occupied sites form a large connected network and some smaller isolated clusters. The large network, called the percolation network, has been found to be a self-similar object. There has been some controversy with regard to the fractal dimension of this network. The source of the controversy is best described by Sta~ffer.'~ Let us simulate the site-percolation problem on a large number of large chessboards for a fixed occupation probability p. From the resulting configurations we can calculate the probability ns(p)for a cluster to have s occupied sites. The average size of the clusters can be defined as the ratio

Near the percolation threshold the diameter of the cluster is given by the correlation length (. Stanley defines the fractal dimension d by the relati~n'~'' (s) = p'.

(2.7)

Near the threshold (s) a (p - p c ) - y and ( cc (p - pJV, where the critical exponents are defined in Stauffer's article. This gives d' = y/v. Stauffer argues that near the threshold the size of the large clusters is better calculated from

where k is an integer equal to or greater than 2. This quantity has a different singular behavior s y a ( p- P , ) - ( ~ + ~where ), /3 is a third critical exponent. If sg is used instead of (s) in the definition of the fractal dimension, we obtain a larger value d = (/3 + y)/v. The second point of view is strengthened by the following scaling argument. Let N(r) be the number of sites on the percolation network within a radius r. Since the correlation length ( is the only length scale in the system, we must have = rd'Y(r/0,

(2.9)

where f ( x ) is the scaling function which is a constant for x << 1. For x >> 1 we l3

'3n

D. Stauffer, Phys. Rep. 54, 1 (1979). H.E. Stanley, J . Phys. A 10, L211 (1977); also in "Kinetics of Aggregation and Gelation" (F. Family and D. P. Landau, eds.), p. 1. North-Holland Publ., Amsterdam, 1984.

216

S. H. LIU

can calculate N ( r ) by (2.10)

N(r) = P P , ,

where P, is the probability that a site belongs to the infinite cluster and D is the dimension of the lattice. The critical behavior of P, is given by P, cc ( p - P , ) ~cc t-p/v; therefore f ( x ) = x p l Vand (2.11)

N(r) = rD-fl/V(r/t)fl/".

The last result gives d" = D - B/v. Note that d" is the fractal dimension when r / t << 1, i.e., the interior of the cluster, and d" = d on account of the identity y 28 = Dv. Another identity, Dv = P(v l), allows one to write d = D/(1 1/6). The exponents involved in these expressions are tabulated by Stauffer,' and one can readily calculate d and d' and compare them in Table I. The reader will find both sets of values quoted in literature, but the set d is favored by more authors in recent publications because it is believed to describe better the intrinsic property of the percolation network. For instance, the density profile of the percolation network is determined by d.I3 The backbone of the percolation network is what remains in the cluster when all dangling chains are eliminated. Imagine that a current flows from one end of the network to the other end, the backbone is the part of the structure that carries the current. In Fig. 5 we show the backbone of a 2D bondpercolation network. We expect the backbone to have a lower fractal dimension because it is less compact. The numbers reported by Kirkpatrick are d = 1.6 for 2D and 2.0 for 3D.l' Of particular interest in the last few years is the critical conductivity problem near the threshold. Consider the site-percolation model. Let p denote the fraction of occupied sites, and pc the threshold fraction. It is found numerically that near but above pc the conductivity of the cluster behaves like

+ +

+

(2.12)

G(P)cc ( P - PAP,

where p 2 1.2 in 2D and 1.7 in 3D. Efforts to relate the critical conductivity to

TABLE 1. FRACTAL DIMENSION OF PERCOLATION NETWORKS I N D-DIMENSIONAL SPACE

d" d"

1 1

1.9 1.8

2.5 2.1

3.3 2.1 ~~

Values of d determined by Stauffer (Ref. 13). Values of d' determined by Stanley (Ref. 13a).

3.8 2.0

4.0 2.0

FRACTALS AND CONDENSED MATTER PHYSICS

217

FIG.5. The backbone of a percolation network as simulated on a 2D square grid (Kirkpatrick ”).

the random-walk problem on the percolation network will be reviewed in Section VI. In a related problem one deals with the critical mechanical rigidity of a percolation cluster with elastic bonds. Physically the model is realized in the process of gelation, i.e., the condensation of small molecules into a macroscopic polymer network of complex branching structures. d. Random Aggregates

Examples of random aggregates are seen in fine particles, snowflakes, dendrites, high polymers, and many living species. Beautiful pictures of these objects are often found in science magazines as visual examples of fractal

218

S. H. LIU

systems. Kaye pioneered the fractal description of fine particles by observing that the profiles (boundaries of the projected images onto a plane) of a large collection of fine particles from metals to shattered Styrofoam have fractal The deviation of the dimensions between 1 and 2 over a range of ~ca1es.l~ dimension from 1 is a measure of the ruggedness of the profile. Carbonblack or soot particles are formed by random collision of small carbon spheres in the turbulent region of a flame. The particle grows in size by adding new carbon spheres to the existing agglomerate, very much like the addition of new equilateral triangles to the sides of a Koch snowflake of a finite stage. When the length scale is as small as the size of the carbon sphere, the profile ceases to have the fractal property, and its dimension becomes 1. It is also possible for an aggregate to have different dimensions in different ranges of length scale. When a metal particle is attacked by an acid, its profile becomes more rugged with a corresponding increase in dimension. In recent years a considerable amount of research effort has been devoted to the understanding of the kinetics of the growth process in relation to the shape and the dimension of these systems. e. Porous Materials

A porous material consists of solid matter with multiply connected voids interspersed within. In constrast to aggregates, which have very open structures with fractal dimension of the order 2 or less, porous materials often have dimensions between 2 and 3. Examples are rocks15 and coal.I6 The pore structure determines the conductivity of the rock when the pores are filled with water. The conductivity measurement in turn is an important tool in prospecting for water, oil, or other mineral resources. Recent progress in the study of porous materials is mainly in the determination of their fractal structure and dimension.

f. Rough Surfaces Rough surfaces are of great importance in chemical technology. Catalysis takes place in real systems on roughened surfaces in order to increase the effective area. Surface-enhanced Raman scattering is observed on the surfaces of noble-metal electrodes which are electrochemically roughened. Surface roughness also affects the low-frequency impedance across an electrodeelectrolyte interface. An effective way to create a rough surface is to make a l4 l5 l6

B. H. Kaye, in “ParticleCharacterization in Technology”(J. K. Beddow, ed.),Vol. 1, p. 81. CRC Press, Boca Raton, Florida, 1984. A. J. Katz and A. H. Thompson, Phys. Rev. Lett. 54, 1325 (1985). H. E. Bale and P. W. Schmidt, Phys. Rev. Lett. 53, 596 (1984).

FRACTALS A N D CONDENSED MATTER PHYSICS

219

material porous, such as activated charcoal. Pfeiffer and co-workers 17-19 have worked out an ingenious way to measure the fractal dimension of rough surfaces by adsorption, and have reported that the absorbing ability of activated charcoal is correlated with the quantity (3 - d), where d is the dimension of the surface. Typically fractals in nature lack the regularity of the Koch snowflake, the Sierpinski gasket, and the Cantor bar. Nevertheless, natural fractals are selfsimilar in the statistical sense, namely that with a large enough collection, one can magnify a small portion of one member and match it with some other member of the collection.

II. Measurements of Fractal Dimension Many methods have been developed to measure the dimension of fractal systems. The choice of a suitable method for a particular system depends on the nature of the system and the range of length scale. We will divide the method into two classes, the measurement in real space, which is based on the use of yardsticks of various lengths, and the measurement in reciprocal space, which is based on light, x-ray, or neutron scattering.

3. REAL-SPACE MEASUREMENTS Real-space methods for measuring the dimension of fractal objects all make use of the Hausdorff definition of dimension and its consequences. The classic example is Richardson’s measurements of continental coastlines.20Totally unaware of Hausdorff’s work, Richardson discovered that the length depends on the scale of the measurement: the finer the scale, the longer the coastline. When a log-log plot of the length versus the scale was made, a straight line with a slope around -0.1 to -0.2 was obtained. Richardson did not understand the significance of his discovery. The slope of his plot is (1 - d), where d is the fractal dimension of the coastline. The same method was used by Kaye to determine the dimension of the profile of fine particle^.'^ An example of this work is shown in Fig. 6. A variation of this method was used by Pfeiffer et nl. to determine the dimension of rough surfaces of porous adsorbent Molecules of different sizes are used to coat the surface of the material until a monolayer P. Pfeiffer, D. Avnir, and D. Farin, SurJ Sci. 126,569 (1983). P. Pfeiffer, and D. Avnir, J . Chem. Phys. 79, 3558 (1963). D.Avnir, D.Farin, and P. Pfeiffer, J . Chem. Phys. 79,3566 (1983). 2o See the article by B. B. Mandelbrot, Science 156,636 (1967).

220

S. H. LIU

' 5 . 2 0.005 0.01

0.02

0.05

0.1

0.2

A FIG.6. The profile of a carbonblack (soot) particle (KayeI4).The rugged curve has the fractal dimension of 1.18 over one decade of length scale I , which is normalized with respect to the longest projected length of the profile.

is reached. The number of adsorbed molecules N is related to the average radius I by the power law N cc rFd, where d is the dimension of the surface. Weitz and Oliveria" studied the aggregates of fine gold particles. The gold particles are of uniform size, but they form aggregates of up to very large sizes. The transmission electron micrograph of one such cluster is shown in Fig. 7. The number of gold particles in an individual cluster was counted from the TEM image of the cluster, and this number was found to be related to the linear size of the cluster by a power law, as shown in Fig. 8. The power, found to be approximately 1.75, is the fractal dimension of the cluster. Stapleton and co-workers"." determined the dimension of the backbone of protein molecules by drawing concentric spheres of radii R = rnl, where 1 is a length scale and rn is an integer, and counting the number of carbon atoms on the backbone within each sphere. One difficulty with this task is that the fractal dimension determined this way depends somewhat on the choice of the scale. The authors varied the scale until the dimension is minimized, and associated the minimum value with the dimension of the protein backbone. An automated procedure using scanning electron microscopy has been developed by Krohn and Thompson for the fractured surfaces of porous materials such as s a n d ~ t 0 n e .The l ~ image of the surface is digitized, and the number of geometrical features of various sizes are counted automatically. The principle of the measurement is the same as the coastline example, i.e., that the finer the feature the more numerous they are expected to be. They have D. A. Weitz and M. Oliveria, Phys. Rev. Lett. 52, 1433 (1984).

FRACTALS AND CONDENSED MATTER PHYSICS

22 1

FIG. 7. The transmisssion electron microscope image of a gold colloidal aggregate (Weitz and Oliveria'').

. .

FIG.8. The number of gold particles N plotted versus the length scale L of the colloidal aggregate in Fig. I. The slope of the log-log plot gives the fractal dimension d = 1.75. The unit of L is 14.5 nm (Weitz and OIiveriaz1).

222

S. H. LIU

reported results for the fractal dimension of sandstone samples ranging from 2.5 to 2.9. 4. SMALL-ANGLE LIGHTOR NEUTRON SCATTERING The principle of the small-angle-scatteringtechnique is that the differential scattering cross section is related to the Fourier transform of the charge distribution (light or x-ray scattering) or the mass distribution (neutron scattering) of the object under study. In case the object is a fractal, its nonintegral dimension will manifest itself in the differential cross section. Since the method is an important tool for studying the structure of solid materials in the micrometer scale, we will take some care to develop the theory, especially to point out the distinction between volume and surface scattering. The latter issue is particularly relevant for fractal systems. The experimental setup is schematically represented in Fig. 9. A beam of photons or neutrons is scattered by the experimental object, and the intensity of the scattered beam depends on the scattering angles 8 and 9. The dependence on the azimuthal angle 4 is averaged out either by the symmetry or by the random orientation of the object. The scattering angle is related to the momentum transfer q by q = 2ksin(8/2),

(4.1) where k is the wave vector or momentum of the incident particles. The ratio of the intensity of the scattered beam to the flux of the incident beam is the differentialcross section a(q), which is the quantity to be measured. We need to develop the theory to relate a(q) to the charge or mass distribution of the scattering object. It is helpful to illustrate some of the basic concepts by considering the scattering problem for a uniform sphere. The scattering potential is V(r) = V,, = 0,

r < a, r > a,

where a is the radius of the sphere. In the Born approximation the scattering amplitude f ( q ) is the Fourier transform of V(r) f ( q ) cc

J V(r)eiq*'d3r

= 4nVO[sin(qa)- qacos(qa)]/q3.

(4.3)

For qa << 1, we can expand the right-hand side of Eq. (4.3) to obtain f ( 4 ) = f(O)(l - R 2 q 2 / 6 ) ,

(4.4)

223

FRACTALS AND CONDENSED MATTER PHYSICS

l NCl DENT

BEAM

-k \ FIG.9. The schematic diagram of a small-angle x-ray or neutron scattering experiment.

where the forward scattering amplitude f(0) cc V0(2xa3/3)is a measure of the total amount of matter, and R Z = 3a2/5is the radius of gyration of the sphere. For qa >> 1, the differential cross section is given by d q ) = lf(4)l2

= 47W-o cos(qa)12(S/q4)9

(4.5)

where S = 4na2 is the area of the sphere. The differential cross section oscillates around a q-4 tail whose amplitude is a measure of the area of the sphere. The cross-over point from small-q to large-q behavior occurs at q = l/a. Therefore, for a sphere of micrometer size or larger and beam particles with wavelengths in the angstrom range, the momentum transfer q is very small compared with the particle momentum k, and the scattering angle of importance is very small. Consider a general distribution of matter with density p(r).22*23We can write I-

In the near-forward direction we expand the exponential factor in powers of q. The leading term is I-

which is the total amount of matter. The linear term in q makes no contribution if the distribution has inversion symmetry or if an average over 22

23

A. Guinier and G. Fournet, “Small-AngleScattering of X-Rays,” p. 5. Wiley, New York, 1955. G . Porod, in “Small Angle X-Ray Scattering” (0.Glatter and 0.Kratky, eds.), p. 18. Academic Press, New York, 1982.

224

S. H. LIU

random orientation is made. The quadratic term is

after angular averaging. The last integral is the second moment of the distribution. Therefore,

k

f ( 4 ) = f(O)( 1 - R 24

19

(4.9)

where R is the radius of gyration of the object. These results hold for q << 1/R regardless of the nature of the distribution, including fractals. For q > 1/Rit is more convenient to study the quantity o(q)= If(q)12 cc

where p2(r)=

s

s

p2(r)eiq”d3r,

(4.10)

+ r’)d3r’

(4.1 1)

p(r’)p(r

is the density-density correlation function. The meaning of p2(r) is as follows. Consider the distribution as a collection of particles, then p2(r) is the probability of finding two particles separated in space by r. For a distribution of uniform density, we find p2(r) = 1 except for surface effects. For a fractal network, such as the path of a self-avoiding random walk or an aggregate in Fig. 7, we expect p2(r) to decrease with increasing r. We choose an arbitrary point on the network and construct a spherical shell of radius r and thickness dr. The amount of matter enclosed in the shell is proportional to 4nrd- dr. The density-density correlation function is this quantity divided by the volume of the shell, 4nr2 dr. Thus p2(r) cc Mrd-3,

(4.12)

where M is the total mass or charge of the object resulting from the integration. Consequently, o(q) K MT(d - 1) sin[n(d - 1)/2]/qd.

I‘

(4.13)

Schaefer et al. have used light and x-ray scattering to study the colloidal aggregates of silica particles.24 Their data, shown in Fig. 10, indicate d = 2.1. Sinha et al. have reported a neutron-scattering experiment on a different kind of silica aggregates with d = 2.5.25 24 25

D. W. Schaefer, J. E. Martin, P. Wiltzius, and D. S. Cannell, Phys. Reo. Lett. 52, 2371 (1984). S. K. Sinha, T. Freltoft, and J. Kjems, in “Kinetics of Aggregation and Gelation”(F. Family and D. P. Landau, eds.), p. 87. North-Holland Publ., Amsterdam, 1984.

FRACTALS AND CONDENSED MATTER PHYSICS

225

t v) 2 W

c

I

0.000 1

0.0 1

0.00 1

0.1

K (l/R)

FIG.10. Combined small-angle x-ray-scattering (SAXS) and light-scattering results on a colloidal aggregate of silica particles. The slope of the log-log plot gives the fractal dimension d = 2.1 over most of range of momentum transfer. Surface scattering from individual particles gives rise to the slope - 4 in the large-momentum-transfer end (Schaefer et a1.24).

The result in Eq. (4.13) applies as long as d c 3. For d = 3 the sine factor vanishes so that there is no q 3 contribution. We have learned from the uniform sphere example that the leading term is actually 4-4 which comes from surface scattering. In the following we will derive the theory for surface scattering for a solid of arbitrary shape. Consider a uniform solid of volume V enclosed by a smooth surface of area S. It is clear that pz(0) = p z V . For finite r the value of p2(r) is represented by the overlap in volume between the object and its identical copy displaced by r, shown in Fig. 11. The volume of the overlap

FIG.1 1 . The calculation of the surface-scattering contribution from a solid body of arbitrary shape.

226

S. H. LIU

region is smaller than the total volume by an amount measured by the area of the boundary, i.e., p2(r) a V -

s-

Ir nldS = V - aSr,

(4.14)

where n is the unit normal of the surface and a is a geometric factor of the order 1. Furthermore r is bounded by D, where D is the caliper of the boundary. Then Eq. (4.10)gives o(q) a p 2

joD

( V - aSr)eiq"d3r.

The integrals can all be carried out and the result contains many oscillatory terms except one, which comes from the term containing S. The oscillatory terms are not observed if we study a collection of objects of various sizes. The nonoscillatory term has the q-4 dependence a(q) = 8np2aS/q4.

(4.15)

If the solid is bounded by a fractal surface, Bale and Schmidt argued that Eq. (4.14)should be replaced by16 p2(r) = V - Nr3-d,

(4.16)

where d is the Hausdorff dimension of the surface. Recall that the quantity subtracted from V is the volume of the nonoverlapping part in Fig. 11. We estimate this volume by covering it with cubes of volume r3. We find that the smaller the cube the more fine detail on the boundary can be resolved so that the number of cubes needed to cover the nonoverlapping part goes like Nr-d, where N is a proportionality constant characteristic of the fractal surface. Using this result, we find the surface contribution to the differential cross section to be a(q) a nNT(5 - d)sin[n(d - 1)/2]/q6-d.

(4.17)

This formula resembles Eq. (4.13) for the volume scattering case, but has totally different physical meaning. Since d < 3 in both cases, the volumescattering part gives a q-" dependence with a < 3, but the surface-scattering part gives a > 3. Thus we can deduce the nature of the fractal object from whether the exponent of the q dependence is greater than or less than 3. Bale and Schmidt reported that lignite coal from North Dakota exhibits a -3.5 power-law q dependence in its x-ray-scattering cross section (see Fig. 12), and concluded that the material is a porous solid whose pores have rough surfaces with d = 2.5.16

227

FRACTALS AND CONDENSED MATTER PHYSICS

10'0 :

10s

7

t v)

El08 b

z 0

1

2 107 Ic

a

5: W

L 10' +

9

W

K

105

104

104

103

102

10"

SCATTERING ANGLE IRADIANSI

100

FIG.12. The small-angle x-ray-scattering intensity for a lignite coal. The slope of the log-log plot has the slope of - 3.5, which indicates that the surface of the pores has the fractal dimension d = 2.5 (Bale and SchmidtI6).

111. Growth of Fractal Aggregates In this Part we will summarize our current understanding of how fractals are formed in nature. This is by far the most active area of fractal research. An international conference has been devoted to this topic, and review talks have been featured in many other recent conferences on disordered systems. The

228

S. H. LIU

problem is being studied in three ways: computer simulation, kinetic equations, and experiments, in descending order of activity. As a result, there are more theoretical predictions than experimental results to verify them. It seems feasible at this time to review only the basic ideas and leave out most of the theoretical details. Interested readers are urged to consult the reviews by Stanley26 and S a n d e ~ - . ~ More ’ * ~ ~ space will be devoted to experimental results.

5. COMPUTER SIMULATION OF GROWTH There is no doubt that fractals are formed in nature by random processes, and the modern computer with its high speed and large memory capacity is a powerful tool to simulate randomness. Since randomness can be implemented in literally an infinite number of ways, a large variety of models has been devised. The models are all embedded on a lattice, but when the system is sufficiently large, the structure of the underlying lattice is unimportant. One of the earliest models is the lattice version of the Eden process for the growth of tumor or cell culture.29At the beginning one puts a cell at the origin of a square lattice. The second cell can grow in any one of the four neighboring positions of the first. The two-celled animal now has six nearest-neighbor positions to allow the third cell to grow. This growth algorithm is repeated a large number of times to form a cluster of cells. Such a cluster is found to be compact, namely that their Hausdorff dimension is equal to the space dimension of the lattice. The cluster contains regions of unoccupied sites, but they tend to situate near the periphery.30 A related model has been provided by Sawada et al. to simulate dielectric b r e a k d ~ w n .The ~ ~ model assumes that, when a branch is formed, the probability of growth at the tip is higher than that on the side by a factor R. For R close to 1 the resulting cluster is compact. For R >> 1, the structure of the cluster looks different on different length scales. On a coarse scale the structure is compact, but on a fine scale it has a fractal dimension of around 1.5. The random animal model is a variation of the percolation p r ~ b l e m . ~ ~ . ~ ~ H. E. Stanley, J . Phys. Soc. Japan 52, (Suppl.), 151 (1983). L. M. Sander, in “Kinetics of Aggregation and Gelation” (F. Family and D. P. Landau, eds.), p. 13. North-Holland Publ., Amsterdam, 1984. L. M. Sander, Proc. NATO Adv. Study Inst. Scaling Phenomena Disordered Syst., 1985. 29 M. Eden, Proc. Berkeley Symp. Math. Statist. IV, 233 (1961). ’O P. A. Rikvold, Phys. Rev. A 26,647 (1982). 3 1 Y. Sawada, S. Ohta, M. Yamazaki, and H. Honjo, Phys. Rev. A 26, 3557 (1982); see also L. Pietronero and H. J. Wiesmann, J . Stat. Phys. 36,909 (1984). 32 D. Stauffer, Phys. Rev. Lett. 41, 1333 (1978). 33 H. P. Peters, D. Stauffer, H. P. Holters, and K. Loewenich, Z . Phys. B 34,399 (1979). 26

*’

FRACTALS AND CONDENSED MATTER PHYSICS

229

Consider the site-percolation case. When the sites are being filled in at random, some occupied sites form clusters in the sense that they are connected by nearest-neighbor bonds. These clusters are called lattice animals, and they have the dimension of 1.54 in 2D and 1.89 in 3D.33These numbers are arrived at by a variety of techniques, e.g., Monte Carlo simulation, series expansion, and the renormalization-group method. Gould et al. shed some light on the fractal dimension of the clusters of some growth models by the real-space renormalization-group technique.34 The method is based on the fact that fractals are self-similar. If one assigns a fugacity or weight K to a cluster of a certain size, and finds that a cluster b times larger in linear size has the fugacity K’, where K‘ is related to K by a transformation R ( K ) , then the fractal dimension of the cluster is given by d = In &/ln b, where AK = ( a K ’ / d K )evaluated at the fixed point of iteration K*, which is the solution of K * = R ( K * ) .The method is completely rigorous provided that one starts with a sufficiently large cluster. However, in practice one must start from a single lattice point and derive the transformation R by considering a rather small cell, say b = 2 or 3. This puts a severe limitation on the accuracy of the method. Nevertheless,the authors showed that one can use the method to classify the various growth models into different universality classes, and thus obtain insight into why they result in clusters with different fractal dimensions. The diffusion-limited aggregation (DLA) model was proposed by Witten and Sander to simulate soot particles.35At the start a seed particle is placed at the center of a square grid. A second particle is introduced into the grid and allowed to random walk on the grid. If it visits a site next to the seed particle, it is stuck to the seed to form a two-particle cluster. If, however, it touches the boundary first, it is discarded and another particle is introduced. A third random-walking particle may join the cluster of two by visiting any one of the neighboring sites of the two particles. By repeating the process a large number of times, the resulting cluster has the form shown in Fig. 13. This object has a fractal dimension of 1.66. In 3D the diffusion-limited aggregate has the fractal dimension of 2.5.36The open structure of the cluster comes about as a result of screening, i.e., when a branch is formed, the next particle has a higher probability to be caught near the tip than the side of the branch. Witten and Sander also pointed out that the DLA model is the discrete version of the Langer model for dendritic g r ~ w t h . ~Through ’ the concept of screening, one can perceive a similarity between the DLA cluster and the model of Sawada et al. for dielectric breakdown. H. Gould, F. Family, and H. E. Stanley, Phys. Rev. Lett. 50,686 (1983). T. A. Witten, Jr. and L. M. Sander, Phys. Rev. Lett. 47, 1400 (1981). 36 P. Meakin, Phys. Rev. A 27,604 (1983). 37 J. S. Langer and H. Muller-Krumbhaar, Acta Metall. 26,1681,1689,1697 (1978). 34

35

230

S. H. LIU

FIG. 13. A diffusion-limitedaggregate of 1500 particles as simulated on a 2D squaregrid. The aggregate has d = 1.66 (Witten and Sander35).

The cluster-cluster aggregation model, proposed by Meakin and Kolb et aL3’ independently, is a variation of the DLA model. In this model there are many particles undergoing random walk simultaneously in the lattice, and any two may stick if they find themselves on nearest-neighbor sites. Similarly, a third particle may join a cluster, and clusters may stick and form larger clusters. The result, shown in Fig. 14, is a more open structure with smaller Hausdorff dimension than the DLA cluster, 1.4 for 2D and 1.8 for 3D. The reason for the more open struoture is also due to screening, because a cluster is more effectivelyscreened by the existing structure than a single particle, and is less likely to reach the compact region of the structure. 38

39

P. Meakin, Phys. Rev. Lett. 51, 1119 (1983). M. Kolb, R. Botet, and R. Jullien, Phys. Rev. Lett. 51, 1123 (1983).

FRACTALS AND CONDENSED MATTER PHYSICS

a

< b

23 1

400 LATTICE UNITS

400 LATTICE UNITS

>

FIG.14. Cluster-cluster aggregates as simulated a 2D square grids. The number of particles are (a) N = 10,000;(b) 15,000; (c)20,000;and (d) 25,000. The aggregateshave d = 1.4 (Meakin3').

An analogous model simulates dendritic electrodepo~ition.~~ Consider a square grid. The initial growth sites occupy the horizontal axis, and the random-walking particles are introduced in the upper half-plane at a fixed density. The motion of the particles is restricted so that no two particles can 40

R.F. Voss, J.Stat. Phys. 36,861 (1984);R.F.Voss and M. Tomkiewin, J . Electrochem. SOC.132, 371 (1985).

232

S. H. LIU

< C

< d

400 LATTICE UNITS

400 LATTICE UNITS FIG. 14. (Continued).

>

>

occupy the same site simultaneously. Unlike the cluster-cluster model the moving particles do not stick together. Upon reaching the growth site, a particle may stick with a finite probability. The density of the particles and the sticking probability determine the result of the growth. With high sticking probability the growth resembles a deciduous forest in winter. Reducing the

FRACTALS AND CONDENSED MATTER PHYSICS

233

A

C

FIG. 15. Simulation of dendritic electrodeposition in the low-concentration regime with high (top) to low (bottom)sticking probabilities(Voss4').

sticking probability enables the particles to pack at higher density so that the branches acquire an increasing amount of foliage, as shown in Fig. 15. High concentration and high sticking probability produces a growth that resembles a wire mesh, as in Fig. 16. The effect of reducing the sticking probability is quite similar. The ballistic model of growth is similar to the DLA model except that the particles rain down from above from random position^.^',^^ The cluster is compact on the large length scale, but has the structure of DLA on finer scales. It also has the general appearance of a fan. The angle at the apex of the fan is 41

42

P. Meakin, Phys. Reo. B. 28,5221 (1983). P. Ramanlal and L. M. Sander, Phys. Rev. Lett. 54, 1828 (1985).

234

S. H. LIU

A

B

C

FIG.16. Simulation of dendritic electrodepositionin the high-concentration regime with high (top) to low (bottom) sticking probabilities ( V O S S ~ ~ ) .

approximately 39", and inside the structure there are streaks of unoccupied regions parallel to the edges of the fan. In a recent paper by Ramanlal and Sander42these features of the cluster are explained on the basis of the kinetic theory and stability arguments. This list of models by no means summarizes adequately the present activities in computer simulation. The serious reader is urged to peruse the many articles in the Proceedings of the International Conference on the Kinetics of Aggregation and Gelation. In addition to growth morphology, the growth kinetics has been studied by simulation. However, more insight can be gained from the kinetic-equation approach, which is the subject of the next two subsections.

FRACTALS AND CONDENSED MATTER PHYSICS

235

6. KINETICTHEORY OF AGGREGATION

The kinetic theory of aggregation has a long tradition, owing to its importance in polymer and colloidal chemistry. The mathematical theory was formulated by Smoluchowski in 1916 to study the time evolution when single particles coagulate to form clusters. The process is assumed to be completely irreversible, which means that clusters cannot break up once they are formed. Let Ck(t)be the concentration of clusters of size k, where k = 1,2,3,. . . is the number of particles in the cluster. The time evolution of ck is governed by two processes, the formation of new k clusters by joining together two clusters i and j, i + j = k, and the loss of k clusters when they combine with other clusters. The reaction rate of the channel ij is proportional to the product of the concentrations of the clusters i and j, and the proportionality constant Kij is called the kernel. Thus, the kinetic equation of cluster growth is

where the overdot denotes the time derivative. The equation has been solved for a few simple forms of the kernel, and the solution reflects a variety of interesting physical situations. We will presently discuss the motivations of these simple kernels and their consequences while keeping the mathematics at an elementary level. More details as well as important references can be found in a review by Ziff.43 It is often easier to study the moments of the distribution defined by

The first moment M , represents the total mass of the system. The equation of motion of M , can be found from Eq. (6.1)to be h l ( t ) = - lim

L

m

1 1icicjKij.

L+m i = l j = L

Ordinarily the mass is conserved; i.e., the limit in the above equation vanishes. We will show later that for certain kernels the limit does not vanish, and Eq. (6.3)indicates a loss of mass to infinity. The infinite cluster, which keeps growing, is a large polymer network called the gel. The remaining clusters form the sol. Thus, the nonconservation of mass is the signal for sol-gel transition. In the absence of gel formation, the equations of motion of the higher 43

R. M. Ziff, in “Kinetics and Aggregation and Gelation” (F. Family and D. P. Landau, eds.), p. 191. North-Holland Publ., Amsterdam, 1984.

236

S. H. LIU

moments are

The simplest possible kernel is one with K , equal to a constant independent of i and j . For this case the reaction rate is independent of the sizes of the two clusters. Physically this should be a good approximation for cluster-cluster aggregation in the limit of low concentration and high sticking coefficient so that two clusters join together as soon as they meet. The equations are simplified by a scaling of the time so that K , = 1. One introduces the generating function m

The differential equation for f can be found to be

j

= $f*,

whose solution is f ( x , t ) = f ( X 9 O)/P - (t/2)f(x, O)].

(6.7)

For initial conditions we consider the monodispersive case where only single particles are present; i.e., cl(0) = 1 and all other c’s are zero.Thus we obtain f ( x , 0) = ex - 1. Re-expansion of Eq. (6.7) allows us to find

+ t/2)k+l,

C&(t) = (t/2)&-1/(1

(6.8)

which gives the complete description of the time evolution of the concentration of clusters of various sizes. One finds that for k ) 1 the concentration is initially zero, increases until reaching a maximum at t = k - 1, and decays afterwards. The first few moments can be found by substituting Eq. (6.8) into Eq. (6.2): MOW

=

1/(1

Mi(t) = 1,

M&)

=

1

+ t/a, (6.9)

+ t,

etc. The monotonic decrease of the zeroth moment with time indicates the continued formation of increasingly large clusters. But the fact that the first moment is time independent means that no infinite cluster, the gel, is formed in a finite amount of time. The average number of particles in the clusters is found by ( N ) = M l ( t ) / M o ( t )= 1

+ t/2.

(6.10)

FRACTALS AND CONDENSED MATTER PHYSICS

237

For large t the average cluster grows linearly with time. If the cluster is a fractal of dimension d, we deduce that its linear size r grows according to the power law

r a t'ld.

(6.11)

It will be discussed later how this law has been verified experimentally. A second kernel which has been analyzed extensively is one which is proportional to the sum of the sizes of the two clusters. With a proper scaling of time we can write K , = i + j. For this problem it is much easier to solve the equations for the moments. These are found to be

ho = -MOMl, h 2

= 2M1M2,

(6.12)

etc. For the monodispersive initial condition we have M,(O) = 1 for all n. Under the assumption that no gel-sol transition occurs, we put M,(t) = 1 and obtain the solutions of Eq. (6.12) Mo(t) = e-', M2(t)= e",

(6.13)

etc. The physics is very similar to the constant kernel case, but the linear dependence on time is replaced by an exponential dependence. Physically this kernel is realized approximately when the stickingcoefficientis low so that two clusters do not stick until they have sampled a large number of possible sticking sites. Consider the cluster in Fig. 14. When it encounters another cluster of the same kind, the chances are high that the tip of a branch of one will come close to a site somewhere on a branch of the other. Under low sticking probability, the two clusters are more likely to join if there is a large number of sites on each for the tip of the other to seek a favorable place to form a bond. Consequently, the reaction rate is higher for larger clusters, and the rate of cluster formation is slow at the beginning but accelerates when the average cluster grows sufficiently large. Recent experiments have confirmed ~*~~ this prediction. The concentrations c&) have also been c a l c ~ l a t e d , 4with the results Ck(t) = (1 - ~ ) ( k u ) l - ' e - ~ " / k ! ,

(6.14)

where u = 1 - e-'. The peak of Ck(t)occurs at t = In k/2. The third kernel is interesting because it leads to the formation of a gel. The kernel K , is proportional to the product of i and j, and we set the M. Golovin, Izu. Geofir. Ser. 783 (1963);Bull. Acad USSR Geophys. Ser. (5), 482 (1963). W. T. Scott, J . Atoms. Sci. 25,54 (1968).

44 A. 45

238

S. H. LIU

proportionality constant equal to one. The kernel applies under the condition that the two clusters can form a bond between every possible pair of sites, and this can happen only if the clusters are highly flexible, like polymer molecules. Again, it is easier to analyze the equations of motion for the moments (6.15) etc. In the pre-gel stage we impose the conservation of mass and set M , = 1 for the monodispersive initial condition. This gives Mo(t) = 1 - t/2,

M J t ) = (1 - t)-',

(6.16)

etc. The second moment diverges at t = 1, and at this time a gel is formed. The concentrations in the pre-gel stage can be solved by the generating-function with the result Ck(t) = kk-2tk-'e-k'/k!.

(6.17)

Using this result we can calculate the net flux to the infinite cluster (6.18) where the limit of L + co is taken after t = 1. For t ) 1, the solution for ck is47,48 ck(t) = ck(t = l)/t.

(6.19)

Thus, all finite clusters diminish in number and are in the end devoured by the gel. Other forms of the kernel, such as i" j" and (ij)", have been analyzed in recent years with regard to the formation of gel.49The subject is outside the scope of this review. From the condensed matter point of view, the gel-sol transition is another form of the percolation problem. The critical mechanical rigidity can be compared with the critical conductivity of the percolation network. This topic will be discussed in a later section.

+

J. B. McLeod, Q. J. Math. Oxford 13, 119, 193 (1962). R. M. Ziff, J. Stat. Phys. 23,241 (1980). 48 F. Leyvraz and H. R. Tschudi, J. Phys. A 14,3389 (1981); 15,1951 (1982). 49 F. Leyvraz, in "Kinetics of Aggregation and Gelation" (F. Family and D. P. Landau, eds.), p. 201. North-Holland Publ., Amsterdam, 1984.

46 47

FRACTALS AND CONDENSED MATTER PHYSICS

239

7. KINETICTHEORY OF DENDRITIC GROWTH

A kinetic theory for the growth of dendrites has been developed over the past four decades. The most modern version was given by Langer and MiillerK r ~ m b h a a r . ~In’ this review we will give a simplified but physical account of this highly mathematical theory. The system under consideration consists of a melt in contact with a solid. The solid is supercooled so that it attracts particles from the melt to condense on its surface. The motion of the particles in the melt is controlled by diffusion. The interface between the solid and the melt moves as the solid grows in size. It has been known for a long time that the system has a steady-state solution,sO*slwhich describes a paraboloidal growth front moving at a constant speed in the melt. Later Mullens and Seserkas2 and Kotler and Tillers3 demonstrated that the steady growth front is inherently unstable such that any small fluctuation in the accretion of particles causes the development of higher harmonics in the growth front. The basic physics is best illustrated by a very simple two-dimensional model of Witten and Sander.54The density of particles in the melt is represented by a diffusion field u(r, t ) which satisfies the diffusion equation

au

- = qv2u. at

where q is the diffusion constant. On the boundary between the solid and the melt the boundary condition is u = 0. In addition the growth velocity of the interface is given by the equation

v, = q i i . vul,.

(7.2)

In the two-dimensional model we consider the growth front to be a circle of radius R. For slow diffusion Eq. (6.19) may be approximated by the Laplace equation, whose solution with circular symmetry is u = Alnr

- AlnR,

(7.3)

where the boundary condition u(R) = 0 is satisfied. Since V , = dR/dt, we find

G . P. Ivantzov, Dokl. Akad. Nauk SSSR 58,567 (1947). D. E. Temkin, Dokl. Akad. Nauk SSSR 132, 1307 (1960); Sou. Phys. Dokl. 132,609 (1960). 5 2 W. W. Mullens and R. F. Seserka, J . Appl. Phys. 34,323 (1963). ” G. R. Kotler and W. A. Tiller, J . Cryst. Growth 2,287 (1968). T. A. Witten and L. M. Sander, Phys. Rev. B 27,5686 (1983). 51

’‘

240

S. H. LIU

from Eq. (7.2) the relation

_-

d R - qA dt

R

*

(7.4)

This equation shows that smaller disks grow more rapidly than larger ones. Next, we consider a small distortion of the growth front as a result of some random perturbation r =R

+ 6cos(mO),

(7.5)

with 6 << R. The solution of the Laplace equation containing the distortion is u = A In r

+ B + C cos(mO)/rm.

(74 The boundary condition u = 0 on the surface gives B = -A In R and C = - A 6 R m - ' . The radial gradient of u on the surface is

au

-=

ar

A/r - mC cos(rnO)/rm''

=A/R

+ (m - 1)A6cos(mO)/R2.

(7.7) Putting this into the equation of the growth velocity yields Eq. (7.4) in the zeroth order of 6 and d6 dt

- = (m -

6 1)qAR

in the first order. Together with Eq. (7.4) we find that ( 8 / 6 ) / ( R / R )= m - 1.

(7.9)

This says that any distortion with m > 1 grows faster than the growth circle. Thus, the circle is inherently unstable, and undulations inevitably develop as a result of fluctuation. In three dimensions we expand distortions in spherical harmonics and find that any distortion with 1 > 1 grows faster than the spherical growth front.52 Witten and Sander argue that whenever the lefthand side of Eq. (7.9) is greater than 1, the growth morphology will be a fractal.5 4 Langer and Muller-Krumbhaar solved the more realistic problem of paraboloidal growth front, which moves uniformly in the direction of the symmetry axis.37 They showed that, due to surface tension, the growth velocity is inversed proportional to the square root of the radius of curvature of the growth tip. Similar to the circular disk model, the growth front is unstable and undulations tend to form due to random fluctuations in particle flow. The result is the growth of side branches on the paraboloid. Since the branches have sharper growth tips, they grow faster than the main front and

FRACTALS AND CONDENSED MATTER PHYSICS

24 1

eventuallybecome unstable themselves. In this manner, new branches grow on the branches to give the solid the dendritic structure. 8. EXPERIMENTAL STUDIES

In this section we will capsulize a series of experimental studies of the growth kinetics of gold colloid aggregates carried out by Weitz and his cow o r k e r ~ . ~ ~The - ~ gold ’ colloid has a fascinating history. As told by Weitz and Huang in a recent review,s5it was first made by alchemists in the Middle Ages. Faraday studied the stability of the colloid against aggregation in 1857.58 The unaggregated colloid has the red-wine color, and to account for this phenomenon Mie formulated the first theory of the scattering of electromagnetic radiation from small dielectric particles.” Smoluchowski’s kinetic theory of cluster aggregation, reviewed in Part 111,6,was stimulated by early works on gold colloids. The system was among the first to be studied by transmission electron microscopy,6’ and Turkevitch discovered that the gold particles are highly uniform in size. The colloid is also useful in medicine and biomedical sciences. Recently the aggregates have been employed as substrates for surface enhanced Raman scattering studies.61 The colloid is made by reducing the gold salt NaAuC1, with sodium citrate. When first made the system consists of spherical particles which have an average diameter of 14.5 nm and a root-mean-square deviation in diameter of about 10% of the mean. The particles are covered with negative citrate ions, which give them a mutual Coulomb repulsion and stability against aggregation. Aggregation is initiated by adding a small amount of pyridine to the solution. The neutral pyridine molecules replace the citrate ions as the adsorbed species on the particle surface, thus allowing the particles to come into contact with one another. Initially the particles are probably bound together by Van der Waals interaction, but eventually metallic bonding occurs at the point where the spheres touch. The cluster formation process is entirely irreversible. By varying the amount of pyridine added to the solution, the aggregation rate can be controlled reproducibly over many orders of A. Weitz and J. S. Huang, in “Kinetics of Aggregation and Gelation” (F. Family and D. P. Landau, eds.), p. 19. North-Holland Publ., Amsterdam, 1984. ” D. A. Weitz, J. S. Huang, M. Y. Lin, and J. Sung, Phys. Rev. Lett. 53, 1657 (1984). ” D. A. Weitz, J. S . Huang, M. Y. Lin, and J. Sung, Phys. Rev. Lett. 54,1416 (1985). ” M. Faraday, Philos. Trans. R. SOC.London Ser. A 147, 145 (1857). 59 G. Mie, Ann Phys. 25,377 (1908). 6o J. Turkevitch, P. C. Stevenson, and J. Hillier, Trans. Faraday SOC.Discuss. 11,55 (1951). J. A. Creighton, C. G. Blatchford, and M. G. Albrecht, J . Chem. SOC.Faraday Trans. I1 75,790 ” D.

(1979).

242

S. H. LIU

magnitude. This makes the system ideal for the purpose of investigating the growth kinetics. Previous studies of the structure of the clusters have revealed their selfsimilar property, with a Hausdorff dimension d = 1.75.21To study the growth kinetics Weitz et al. used the dynamic-light-scattering technique to measure the autocorrelation function of the scattered light.62 The principle of the measurement is quite similar to small-angle x-ray scattering except that the scattered photons are counted at two different times and their correlation function is calculated. This way the time variation of the scattering cross section is detected. There are three time scales one must consider. The longest time scale is the growth rate of the aggregate, next is the diffusion rate of the clusters in the liquid, and the shortest is the rotation and vibration of the individual clusters. In the correlation measurement the time interval is adjusted to be sensitive to the diffusion mode. During the measurement the sizes of the clusters change by negligible amounts, and their internal motions are averaged over. The photon autocorrelation function is proportional to the square of the spatial Fourier transform of the two time density-density correlation function p2(r, t ) =

s

p(r,O)p(r

+ r’, t)d3r’.

(8.1)

In case the motion is dominated by diffusion, the Fourier transform of the above expression can be found to be t ) = 4 4 ) exp( - 9 q 2 t ) ,

(8.2) where 4 4 ) is the scattering cross section in Eq. (4.13).Since the system consists of clusters of all sizes, we must add up their individual contributions to obtain

where ci is the concentration of the clusters of size i. The correlation function is no longer a simple exponential function of time, and this is borne out by the experiments.The data are analyzed by taking the initial logarithmic derivative of M, t )

where (9)is the average diffusion constant weighted by the scattering cross section. The latter is a measure of the mass of the cluster. The 42 dependence has been found to hold for a range of momentum transfer, and the slope of K, B. J. Berne and R. Pecora, “Dynamic Light Scattering.” Wiley, New York, 1976.

FRACTALS AND CONDENSED MATTER PHYSICS

243

FIG.17. Log-log plot showing the scaling of the mean radius of gold colloidal clusters with (Weitz et time. The solid line corresponds to R x

versus q 2 gives the average diffusion constant of the clusters.56 Since the diffusion constant is inversely proportional to the radius of gyration R, the analysis gives the average radius of the clusters. The measurement is repeated over a period of time during which the clusters are allowed to grow. The growth rate is slow enough so that the clusters are unchanging in size during each measurement. In this manner, the time development of R is measured and found to be well represented by the power law R a as shown in Fig. 17.56 Recall that for the constant kernel the kinetic theory predicts such a power law with the exponent l/d, we find that the theoretical value for the exponent is 0.57, in agreement with the measured value. The authors also discussed the dependence of the average cluster size on the initial concentration. According to Eq. (6.10) the average mass grows linearly in time, and the proportionality constant is the initial concentration. Therefore, the radius is expected to be proportional to the concentration raised to the power l/d. This relation is well verified. As mentioned earlier, the aggregation rate can be controlled by the amount of pyridine added to the solution. With a minute amount of pyridine the growth rate is found to depend exponentially on time, R cc eC‘, as shown in Fig. 18.” This behavior is predicted for the sum kernel in Part III,6. Thus, the small sticking coefficient is the dominating factor in the growth process. The clusters are more compact, with a fractal dimension around 2.05. For intermediate amounts of pyridine, the growth starts exponentially with time and turns over to a power law at a later time. The fractal dimension lies between the limits 1.75 for diffusion-limited aggregates and 2.05 for reactionlimited aggregates. The results also suggest that the silica aggregates studied by Schaefer et a1.24 and Sinha et aL2’ are grown under reaction-limited conditions.

244

S. H. LIU

TIME

(HRS)

FIG. 18. Increase of characteristic cluster size with time for (a)reaction-rate-limited, (b) crossover, and (c) diffusion-limited kinetics. The insert shows the initial behavior on an expanded scale. The rate-limited data are fit to R a e112.24 (Weitz et ~ 2 1 . ~ ' ) .

The two-dimensional diffusion-limited aggregation model has been realized experimentally by Matsushita et al. by electrolytic deposition of zinc leaves at the interface between two immiscible liquids.62aThe electrolytic cell is a vat containing 2M ZnSO, aqueous solution as the electrolyte. Above the electrolyte is a layer of n-butyl acetate (CH3COO(CH2),CH3).A carbon cathode is inserted so that its tip is just at the interface between the two fluids. Electrodeposition is initiated by a dc voltage between the carbon cathode and a ring-shaped zinc anode. The zinc metal leaf grows two dimensionally along the interface from the tip of the cathode outward. For bias voltages below a critical value, roughly 8 V, the fractal dimension of the metal cluster is a . ~ ~ the constant very close to 1.66,as predicted by computer ~ i m u l a t i o nAbove critical voltage the fractal dimension increases steadily. The electric current is a measure of the growth rate. Below the critical voltage the current increases almost linearly in time after an initial steep rise. This shows that a compact center core is first formed, and subsequently the growth process can be described by the Smoluchowski equation with the constant kernel. The growth process above the critical voltage is not sufficiently understood. IV. Diffusion on Fractal Networks

The problem of diffusion on fractal networks originated from attempts to understand the conductivity threshold of the percolation cluster. Early efforts were based on scaling. Recently exactly soluble models have been devised to check the validity of the scaling arguments. Experimental results are sparse so M. Matsushita, M. Sano, Y. Hayakawa, H. Honjo, and Y. Sawada, Phys. Rev. Lett. 53, 286

(1984).

FRACTALS AND CONDENSED MATTER PHYSICS

245

that comparison between theory and experiment is qualitative at best. The major problem lies in the difficulty in realizing a percolation network in the laboratory. 9. ANOMALOUS DIFFUSION

The question of anomalous diffusion was raised by de G e n n e ~and , ~ ~the scaling theory for the phenomenon was worked out by Alexander and O r b a ~ hGefen , ~ ~ et al.,65and Rammal and Toulouse.66 The problem has been reviewed recently by Aharony6' and Orbach.68 The conductivity of a percolation network is related to the density of charge carriers N and their diffusion coefficient 9 by the relation

(9.1) Since only the electrons on the occupied sites can participate in the conduction process, the quantity N is proportional to P,, the probability that a site belongs to the infinite network. Using the known relations G a (p - pC)' a t-'/' and P, a (p - P , ) ~cc t-@/', we find 9 a <-('-B)Iv = 5-8. (9.2) G a N%

For a part of the network or a small cluster the diffusion coefficient 9(1) can be , f(x) derived from the scaling argument. We assume that 9 ( r ) = r y f ( r / t ) where is constant for small x but has the power-law form for large x. For r >> 5, Eq. (9.2) implies that f ( x ) a x'. Consequently, y = -8 and 9 ( r ) a r-e. The physical meaning of the anomalous diffusion coefficient is as follows. In the normal diffusion problem the mean-square displacement is given by

( 9 )= 9 t ,

(9.3) where the time interval t is proportional to the number of steps taken by the particle. For anomalous diffusion we find that the time is related to the In analogy with the random-walk and the self-avoiding distance by t a random-walk problems, the exponent in the steps-versus-distance relation gives the fractal dimension of the path of the walker. Thus, the path of a random walk on the percolation network has the fractal dimension d, = 8 + 2, which has been estimated to have the value 3.0 for 2D and 3.5 for 3D. P.G. de Gennes, Recherche 7,919 (1976). S. Alexander and R. Orbach, J . Phys. Lett. (Paris) 43, L-625 (1982). " Y. Gefen, A. Aharony, and S. Alexander, Phys. Rev. Lett. 50,77 (1983). 66 R. Rammal and G . Toulouse, J . Phys. Lett. (Paris) 44, L-13 (1983). 'b A. Aharony, Proc. N A T O Adv. Study Inst. Scaling Phenomena Disordered Syst., 1985. R. Orbach, Proc. N A T O Adv. Study Inst. Scaling Phenomena Disordered Syst., 1985. b3

64

246

S. H. LIU

In the time interval t the diffusion distance is given by R a t”(e+2).The volume enclosed within R is V aR d a t d l ( e + 2 ) . (9.4) The density of diffusion modes is proportional to the reciprocal of V, or t - d ’ ( e + 2 ) . In the energy space the mode density is given by the Fourier transform of the time dependence, which is

N ( o ) a ox,

(9.5)

with = 1 - d/(e

+ 2).

(9.6)

The value of x has been estimated to be 0.32 for 2D and 0.30 for 3D. The conductance of the percolation network can be thought of as the input conductance of a resistive network which has a unit resistor connecting every pair of occupied sites. The input and output terminals are separated by a large distance measured by the correlation length 5. As pointed out by Stephen:’ the problem can be generalized by adding a unit capacitance connecting every site to a common ground. The frequency dependence of the conductivity of such a network was studied by Gefen et al? by using the scaling argument. It is assumed that the only time scale in the problem is t a te+*.At low frequencies the ac conductivity reduces to the dc value G ( p ,o)a t-””a ~ ~ ‘ ( o t ) - ~ ’ ,

where x’ = p/v(O + 2). At high frequencies the conductivity has the scaling form G ( p , o ) = o”’g(ot),and the scaling function g is a constant. Thus, the high-frequency behavior is given by the exponent x‘. Using the known values for the critical exponents, one estimates that x’ = 0.32 for 2D and 0.60 for 3D. An experimental attempt to verify this frequency law was not entirely s u c ~ e s s f u lIn . ~this ~ experiment the specimens were thin Au films deposited on insulator substrates. The transmission electron microscope image of a typical film is shown in Fig. 19. The metal part is in the form of interconnected clusters,and by varying the amount of the metal one can make Au films above or below the metal-insulator transition. The ac resistance data of several films are shown in Fig. 20. The five lower curves are for specimens in the metallic regime, and the two upper curves are for insulating ones. At low frequencies the resistance of a metallic film approaches a constant value, but above a characteristic frequency the resistance follows a power law with an exponent 0.95 & 0.05. For an insulating film a power law with the same exponent is seen 69 ’O

M.J. Stephen, Phys. Rev. B 17,4444 (1978). R. B.Laibowitz and Y.Gefen, Phys. Rev. Lett. 53,380 (1984).

FRACTALS AND CONDENSED MATTER PHYSICS

247

FIG.19. Transmissionelectron microscopeimageof a thin gold film deposited on an insulating substrate. This sample is just below the percolation threshold.The typical cluster width is about 10 nm (Laibowitz and Gefen”).

for many decades of frequency down to the lowest frequency of the measurement. There is qualitative agreement between the data and the scalinglaw predictions, but the measured exponent is quite different from the expected value, apparently because the specimens were not accurate realization of the model network. There is also a theoretical problem with the scaling argument. In a diffusion system an ac signal of higher frequency propagates over a smaller volume than a signal of lower frequency. It is not correct to scale the conductivity, which leaves out the volume factor. Instead, one should scale the total admittance in the diffusion volume. In the dc limit the admittance is equal to the conductivity

248

S. H. LIU

lo7,

=

I x

los

h

0 x

-

.-

10’

c

+

A

A

51 F)

+

+

B-la

A

+ + :

0

A

A

A

0

0

0

0

0

0

0

0

v (r

1

0

~

~

0

1

Au Clusters

8-53

+ +

I

0 x

A-I3

-

10‘

I

I

A-43

~

0

e

A A !

R

8-51 10’

-

t 11 10’

o

o

0

0

E-13 I

10’

1 10‘

f(W

0

I 10’

0

0

0

I 1o6

0

0

I 10’

FIG.20. The real part of the ac impedence plotted versus frequency for a number of gold films in the metallic and insulating side. The labels identify the individual samples (Laibowitz and Gefen”).

multiplied by the cross-sectional area tD-’of the diffusion volume and divided by the length 5. This changes the exponent of the frequency dependence to x’ = (p/v - D + 2)/(0 + 2). The correction makes no difference for a 2D system, but changes the frequency exponent for all other values of D. The expression for x’ can be simplified to x’ = 1 - d / ( 0 + 2), which is the same exponent which enters the density of diffusion modes, Eq. (9.6). In the future we will denote the exponent in both cases by x. The fact that the admittance rather than the conductivity obeys the scaling relation will be borne out later after analyzing an exactly soluble model. Unfortunately the disagreement between theory and experiment persists. Consider the electrical network envisioned by Stephen. If the resistors connecting the occupied sites are replaced by inductors, the network will have a spectrum of resonant frequencies. The density of these resonant modes, called fractons by Alexander and O r b a ~ h is, ~related ~ to the density of diffusion modes. One merely needs to replace the frequency by its square in Eq. (9.6) and multiply the result by the frequency. Hence, the density of

FRACTALS AND CONDENSED MATTER PHYSICS

249

+

fractons is Nf(o)a d - ’ , where df = 2d/(B 2) is called the fracton or spectral dimension. The estimated values for df are around 1.4 for D = 2 to 5 and approaching 4 for D = Derrida et al. showed by the effectivemedium theory that the fracton concept only applies when the wavelength of the excitation is less than the correlation length of the cluster.’l The lowerfrequency excitations are phononlike. The concept of fractons was proposed by Alexander and Orbach concerning a series of experiments by Stapleton and co-workers.l0*” In these experiments the spin-lattice-relaxation rates were measured for many protein molecules in the temperature range of 4-20 K. For crystalline solids of dimension D the relaxation rate has been shown to have a power-law temperature dependence with the exponent equal to 3 + 2D.72The dimension of the lattice enters the result through an integration over the square of the frequency spectrum. Stapleton et al. argued that for protein molecules, which are closely simulated by the paths of self-avoiding random walks, the dimension D should be replaced by the fractal dimension d. Indeed this gives values of d in good agreement with direct measurements. Alexander and Orbach pointed out that the spectrum of the elastic vibrations of a fractal system should have the dimension df, which is in general less than d.64*68On the other hand, the observed temperature exponent certainly does not agree with the formula 3 + 2d,. To make the matter more complicated, a number of authors have shown recently that the elastic-vibration problem does not map onto the diffusion problem, and as a result another spectral dimension comes into the picture. We will come back to these questions in the next section. 10. DIFFUSION ON THE SIERPINSKI GASKET

As discussed earlier experimental verification of the multitude of ideas in the last subsection are fraught with problems. This makes exactly soluble models very useful as testing grounds of the predictions. The Sierpinski gasket (SG) and its extension in D-dimensional space have been studied extensively by many investigators. This section is devoted to the diffusion problem on these model systems. We will show first that the diffusion problem is exactly mapped onto a resistor-capacitor (RC) network. The one-dimensional problem suffices to illustrate the basic principle. Consider the RC network in Fig. 21, for which all resistors have the value R and all capacitors the value C. The equations for the ”

’*

B. Derrida, R. Orbach, and K. -W. Yu, Phys. Rev. B 29,6645 (1984). R. Orbach and H. J. Stapleton, in “ElectronParamagnetic Resonance”(S.Gschwind,ed.),Ch. 2. Plenum, New York, 1972.

250

S. H. LIU

FIG.21. The RC equivalent network for a 1D diffusion problem.

voltage V , at the node n and the current I between the nodes n and n

+ 1 are

V,+, - V , = -I,,R, dV, In- In-1= -c-. dt

(10.1)

In the continuum limit these equations reduce to the diffusion equation for both V and I with the diffusion coefficient 9 = l-’RC, 1 being the distance between adjacent nodes. The difference equations are solved by

V , = V, exp( - an - i o t ) and similarly for I,,.The quantity a satisfies the equation cosha = 1 - ioz/2,

(10.2)

where T = RC is the time constant of the network. For oz << 1, one obtains a = ( - i o z ) - ’ ” , The continuum limit is realized when o7 << 1 , which is the limit of low frequency. The characteristic admittance Y(m)is the ratio between Inand V,

Y(o) = R-’(-~ox)’’’

(10.3)

in the low-frequency limit. The same considerations apply to the diffusion problem in an arbitrary dimensional space. On a fractal network the equivalent network is the RC network discussed in the last section, i.e., one with R between nearest-neighbor sites and C connecting every site to the common ground. The SG network is constructed by the procedure shown in Fig. 2. The unit triangle has a resistor R on every side and a capacitor C / 2 between every vertex and the common ground. When three of these units are linked together to form the next stage, the inner vertices have the capacitor C but the three outer ones have C / 2 .For a sufficientlylarge unit all except the three outermost vertices will be connected to ground through a capacitor C .

25 1

FRACTALS AND CONDENSED MATTER PHYSICS

We use the method of decimation, formulated for this problem by Clerc et al., to calculate the ac The voltages at the three vertices of the basic unit are denoted by F,i = 1-3, and the currents flowing into the vertices are li.These variables can be arranged as column vectors, and the voltage and current vectors are related by the admittance matrix

I = YV,

(10.4)

where Y has the simple form

(10.5) with Y, = 2 / R - iwC and Yz = - 1/R. The sinusoidal time dependence e-'"' is assumed. In the next step of the buildup process three units with vertices (1,2,3), (1', 2', 3'), and (l", 2", 3") are linked together by joining the pairs of vertices (2,1'), (3,1"), and (3', 2"). For every pair of linked vertices the voltages are equal and the currents add up to zero. By using these relations, we can eliminate six sets of voltage and current variables from three sets of relations in Eq. (4.10). The result is three linear relations between the voltage and current variables of the three outer vertices (1,2', 3"), and these can be arranged in the form of a matrix equation like Eq. (10.4):

I ' = Y'V',

(10.6)

where the new admittance matrix has the same form as in Eq. (10.5) with Y1

=

YI(Y1

Y, = -Y;(Y,

-

Yd(2Y1 -

+ 3 Y M Y 1 + Yd(2Y1 - Yd,

Y,)/(Y,

+ Y,)(2Y,

-

Y,).

(10.7)

The procedure is repeated for larger networks. In every step the form of the admittance matrix is preserved, and the new components are related to the old ones by the same relations as in Eq. (10.7). Thus, these are the iteration relations for increasingly large SG networks. One finds that for the nth-order The dc problem is trivial to work SG the admittance matrix elements are Yp) = (2/R)(3)"-' and Y$" = -(l/R)(+y-'. We use two of the three terminals as the input and output, then the admittance between the terminals is r, = Yp) - Y$') = (3/R)(3)"-'. If the length of the basic equilaterial triangle is taken as one unit, then for the nthorder SG the length is L, = 2"-'. We therefore can relate Y. cc L;", where a = p / v = (In 5 - In 3)/ln 2. The conductivity G, is found by multiplying the 72a

J. P. Clerc, A. -M. S . Tremblay, G .Albinet, and C. D. Mitescu,J . Phys. (Paris)45, L-913(1984).

73

Y.Gefen, A. Aharony, B. B. Mandelbrot, and S. Kirkpatrick, Phys. Reo. Lett. 47,1771 (1981).

252

S. H. LIU

admittance by the length and dividing the result by the cross-sectional area. For a 2D system the two factors are both L,, so G,, a L;". The carrier density is the number of vertices divided by the area of the nth-order equilateral triangle, with the result N, a L i b , where b = p/v = 2 - ln3/ln2. From these we deduce the exponent for the anomalous-diffusion coefficient 8 = In 5/ln 2 - 2. At low frequencies o << 1/RC, we can carry out the iteration to the first order of o,with the result Y';)(w)= (1/R)(3)"-'(2 - i5"-'w7).

(10.8)

This allows us to relate Yy)(o) = (;)Yy+')(op).

(10.9)

For the infinite network both Y y ) and Y(;+')approach to the same limit. Therefore, Eq. (10.9) determines the frequency exponent of the infinite network Y ( o )a (- iw)x,where x = 1 - In 3/ln 5 = 0.32. This is in complete agreement with the prediction of the anomalous diffusion theory x = 1 - d / ( 8 + 2). Having verified the theory, we can use it to calculate the fracton dimension d , = 2 In 5/ln 3. The nature of the eigenvalue spectrum and eigenstates of the fractons have been discussed in detail by Ramma174,75and Tremblay and Southern.76 For the finite network the ac response has different behavior in three frequency regions. Define 7,, = RC(5)"-' and z = RC. The conductivity is a constant for 07, << 1, follows the power law ( - io)xfor w7, >> 1 << 07, and is - ioC/2 for oz >> 1. The last region is irrelevant because the discrete-network model for diffusion breaks down. The ac conductivity for the finite SG may be compared with that of a percolation network above the threshold, because both have a connected resistive path between the input and output terminals. The conductivity of the infinite SG corresponds to that for a cluster below the threshold. The frequency exponent x = 0.32 is close to the value for the 2D percolation cluster. The same calculation can be carried out for the generalized SG in D dimensions. The exponent relating to L, is a

=

[ln(D

+ 3) - ln(D + l)]/ln2.

(10.10)

The cross-sectional area is L:-', and the length is L , . Thus the exponent for the conductivity is a' = a + D - 2. The exponent b can be found to be D In(D + l)/ln 2. These give the exponent for the anomalous-diffusion coefR. Rammal, Phys. Rev. B 28,4871 (1983). R. Rammal, J . Phys. (Paris)45, 191 (1984). 76 A. -M. S. Tremblay and 9. W. Southern, J . Phys. Lett. (Paris)44, L-843 (1983).

74

75

FRACTALS AND CONDENSED MATTER PHYSICS

ficient 6, = a' - b of Eq. (10.8) is

= ln(D

+ 3)/ln 2 - 2.

253

The D-dimensional generalization [ D - i(D

+ 3)"-'wz].

(10.11)

From this result we verify the scaling relation for the diffusion time 7,, K ( D 3Y-l = L:+2. The exponent for the frequency dependence of the ac conductivity can be calculated directly from the generalization of Eq. (10.9)

+

(10.12)

+

with the result x = 1 - ln(D + l)/ln(D 3) = 1 - d / ( 8 + 2), where d = ln(D l)/ln2 is the fractal dimension of the system. This result agrees with that of the corrected dynamical scaling theory, Eq. (9.6). In fact, Eq. (10.12) shows that it is the admittance, rather than the conductivity, that obeys the scaling relation. This argument was used in the last section to correct the dynamical scaling theory for anomalous diffusion. The SG has been used by Stephen to model the diamagnetism of percolation clusters composed of superconducting grains.77It was predicted on this basis that the diamagnetic susceptibility of the system does not diverge at the percolation threshold. In a similar problem Rammal and Toulouse7' and Alexander 7 9 have investigated the motion of tight-binding electrons on an SG in the presence of a magnetic field. The properties of the eigenvalue spectrum and eigenstates of the system have been discussed.

+

1 1. THEac RESPONSE OF ROUGHINTERFACES The electrical response of the interface between a metal electrode and a solid ionic or liquid electrolyte has been studied for a long time owing to its importance in electrochemistry. The standard theory may be summarized as follows. Any potential difference between the electrode and the bulk of the electrolyte is screened by a space-charge layer in the electrolyte." This layer and the surface of the metal electrode form a parallel-plate capacitor. The ac signal encounters an impedance which consists of the ohmic resistance in the electrolyte and the interfacial capacitance. The resistance in the metal is negligible in comparison. The frequency dependence of the impedance is quite M. J. Stephen, Phys. Lett. W A , 67(1981). R. Rammal and G . Toulouse, Phys. Rev. Lett. 49, 1194 (1982). 79 S. Alexander, Phys. Rev. B 29,5504 (1984). G . H. Wannier, "Statistical Physics," Ch. 17. Wiley, New York, 1966. 77

78

254

S. H. LIU

simple, that the real part is independent of the frequency while the imaginary part is inversely proportional to the frequency. It has been known since 1926 that real systems do not behave in this simple manner.” One invariably finds that both the real and the imaginary parts of the interfacial impedance follow the same power-law frequency dependence, and their magnitudes are so related that the total impedance has the form

Z(o)a (-

io)-x,

(11.1)

where the exponent x lies between 0 and 1. In the language of electrochemistry this impedance of unknown origin is called the constant-phase-angle element or constant-phase element. Within the past decade many investigators have demonstrated that the value of x is intimately related to the roughness of the interface, with x approaching 1 when the surface is made increasingly ~ m o o t h . ~ ~Under - ~ ’ microscopy even well-polished surfaces show long lines of scratches with jagged faces. Efforts have been made to model these grooves by transmission lines, but these models invariably predict x = 4. Lyden et al. observed a similar effect at the interface between a polycrystalline semiconductor and an aqueous electrolyte.86Cohen and Tomkiewicz have suggested that a new percolation process is taking place in the surface layer of the semicond~ctor.~’ Since the effect is also seen for a metal electrode, this shows that diffusion in the electrolyte is an important contributing factor. Le Mehaute and Crepy have suggested a connection between the phenomenon and the fractal nature of the interface, and have argued that x = l/d, where d is the fractal dimension of the projection of the rough interface, 1 < d < 2.88 Their reasoning is rather obscure at places, and it is difficult to assess its validity. A soluble model for the ac response of a rough interface has been proposed by Liu89and refined by Kaplan and Gray.” The model is based on the Cantor bar in Fig. 4. If the sections at the different stages of construction are linked together as shown in Fig. 22, we obtain a fractal model for the grooves on the electrode, or protrusions in the electrolyte. Each groove has the self-similar structure in that it subdivides into two branches, and each branch is similar to I. Wolfe, Phys. Rev. 27,755 (1926). See the review by P. H. Bottelberghs, in “Solid Electrolytes”(P.Hagenmuller and W. Van Gool, eds.),Ch. 10. Academic Press, New York, 1978. 83 R. de Levie, Electrochim Acta 10,113 (1965). 84 P. H. Bottelberghs and G. H. J. Broers, J . Electroanal. Chem. 67, 155 (1976). ” R. D. Armstrong and R. A. Burnham, J . Electroanal. Chem. 72,257 (1976). J. K. Lyden, M. H. Cohen, and M. Tomkiewicz, Phys. Rev. Lett. 47,961 (1981). M. H. Cohen and M. Tomkiewicz, Phys. Rev. B 26,7097 (1982). ” A. Le Mehaute and G . Crepy, Solid State Ionics 9/10,17 (1983). 89 S. H. Liu, Phys. Rev. Lett. 55,529 (1985). 90 T. Kaplan and L. J. Gray, Phys. Rev. B 32,7360 (1985). 82

’’

255

FRACTALS AND CONDENSED MATTER PHYSICS

FIG.22. A Cantor bar model for a rough electrode-electrolyte interface (LiuE9).

the whole groove when magnified by a factor a, a > 2. The dimension of the Cantor bar is d = ln2/lna c 1. The projection of the interface has the dimension d + 1 and the surface itself has the dimension d + 2. The electric circuit analog of the interface is shown in Fig. 23. The resistance R in the electrolyte increases by the ratio a at every stage of branching because of the reduction in cross-sectional area. The capacitance C, which is assumed to be the same at every stage, represents the interfacial capacitance of the two lateral faces of the branch. The interfacial capacitance in the dip between branches becomes decreasingly important at higher branching stages, and can be shown to be an irrelevant parameter; i.e., it does not affect the exponent of the frequency dependence in the low-frequency limit. The common ground is the electrode. The impedance of the network in Fig. 23 between the input terminal and the ground has the form of a continued fraction: Z(O) = R

+ - i o C1 +

2 ---...

1

2 a R - ioC+ a 2 R -

’

(11.2)

where we have used the condensed form of writing the fraction?’ One can show that Z ( w ) satisfies the following frequency scaling relationg0

(11.3) In the low-frequency limit this reduces to (11.4)

It is straightforward to solve for the frequency exponent from this relation, 91

See “Handbook of Mathematical Functions,” p. 19. National Bureau of Standards Applied Mathematics Series 55 (1964).

256

S. H. LIU

R

c+

a3R

*

FIG.23. The RC equivalent network for the rough interface in Fig. 22 (Liu").

with the result x = 1 - In 2/ln a = 1 - d. The model can be generalized to N branches with the ratio a > N. The frequency exponent for this case is x = 1 -In N/ln a = 1 - d. A smooth surface has few grooves (small N) and each branch has a large reduction ratio in area (large a).Consequently d is small and x is close to 1, as seen experimentally. The fractal dimension of the surface is d, = d + 2, so we can relate x to d, by x = 3 - d,. Real interfaces do not have the regularity of the Cantor bar model. Kaplan and Gray have generalized the model such that both the number of branches and the area reduction ratio a are random variables. They have shown that the impedance of the system has a power-law frequency dependence with the

257

FRACTALS AND CONDENSED MATTER PHYSICS

exponent x given by the implicit equation 1 = (a"-')(N).

(11.5)

The fractal dimension of the surface d, is also given by an implicit equation very similar to Eq. (1 1.5). In fact, they showed that the relation x = 3 - d, continues to hold. Since x and d, can be measured independently, this result suggests a possible way to test the theory. One can establish a connection between the result of this model and the anomalous-diffusion theory. The crucial quantity is the parameter 8, which can be directly calculated. Consider a random walker starting from a point in the middle of the network in Fig. 23. If we assume that the walker may take a step to the left or to the right with equal probability, then after t steps there is a net drift to the right by the amount R = t/3 because there are twice as many paths in this direction. If we assume that the probability of each step is proportional to the cross-sectional area of the path, then the walker will drift to the left by the amount R = -t(a - 2)/(a 2). In either case we find 8 = - 1. Putting this result into Eq. (9.6),we obtain x = 1 - d, as obtained earlier. The Cayley tree is a special case of the model network with a = 1. However, the Cayley-tree model gives a fixed x = +.It is necessary for the circuit to have the fractal property in order to obtain values of x other than or 1. We conclude this Part by observing that the theory of anomalous diffusion on fractal networks seems to be on firm ground. The difficulty at the present time is to realize a fractal network in the laboratory, because specimens of percolation systems consist of metal clusters of irregular shapes imbedded in an insulating matrix or deposited on an insulating substrate. It is not clear how the theory needs to be modified for these real systems. We anticipate rapid development in this area in the near future.

+

4

V. Elastic Properties of Fractal Networks This Part deals with two topics, the critical elastic threshold of percolation networks and the spectral dimension of elastic vibrations of fractal systems. The theory is not as well developed as the diffusion theory, and there is much room for careful studies. The problems now under investigation touch upon the general questions of elastic properties of disordered solids. For this reason the field deserves attention by physicists as well as materials scientists. 12. CRITICAL ELASTIC THRESHOLD OF PERCOLATION NETWORKS It was pointed out earlier that the gel-sol transition involves the formation of an infinite cluster, and this has been compared with the percolation problem

258

S. H. LIU

at the threshold. de Gennes argued that the elastic rigidity of the gel is analogous to the critical electrical conductivity of the percolation network.92 The equivalence of the two problems can be established if one uses an elastic model for which the elastic potential of a pair of points on the network is proportional to the square of the relative displacement vector of the two points. Denote the coordinates of the two points by Ri and Rj and their displacements by ri and rj. Then the elastic potential is

Fj = (K/2)(ri - rj)2 ,

(12.1)

where K is the elastic constant. This model for the elastic potential is known as the isotropic force model because Fj is an isotropic function of the relative displacement vector. The model is unrealistic in one respect. Consider a rigid rotation through a small angle 8.The relative displacement vector is given by Ir, - rjl = IRi - Rile. One can see from Eq. (12.1) that the system has a finite resistance against rigid rotation. Feng and Sen showed that the familiar mass and spring model is free from this and the appropriate elastic potential is

yj = (K/2)[(ri - rj) .(Ri - Rj)I2.

(12.2)

This model has been called the central-force model. One can easily see that in this model elastic interactions can not propagate around a right-angle bend, and consequently a percolation network on a square lattice is unstable against shear. To circumvent this problem the authors studied two lattices that are stabilized by central forces alone, the 2D triangular lattice and the 3D facedcentered-cubic(fcc) lattice. They studied numerically the bulk modulus K, and the shear modulus p, as functions of the concentration of bonds, and found that both moduli exhibit the power-law singularity at the threshold K.3 Pe

(P- P C Y ,

(12.3)

where the threshold pc = 0.35 and t = 1.2 for the 2D lattice with isotropic elastic forces and pc = 0.58 and t = 2.4 for the same lattice with central forces. The elastic threshold is much higher for the central-force model, which reflects the fact that the system is less rigid. The larger critical exponent means that the elastic modulus increasesmore slowly, again indicating a softer lattice. Similar results were obtained for the 3D lattice. The authors conclude that the two kinds of forces belong to two distinct universality classes. In other words, the elastic problem is not equivalent to the diffusion problem. Bergman and Kantor reached similar conclusions by studying an elastic Sierpinski gasket, for which every link is a simple spring (see Fig. 24a).94The P. G . de Gennes, J . Phys. Lett. (Paris)37, L-l(l976). S. Feng and P. N. Sen, Phys. Rev. Lett. 52,216 (1984). 94 D. J. Bergman and Y. Kantor, Phys. Rev. Lett. 53,511 (1984).

92

93

FRACTALS AND CONDENSED MATTER PHYSICS

259

basic unit of three equal springs arranged in an equilateral triangle can be stretched in three linearly independent ways, as shown in Fig. 24b. The amplitude of the distortion is infinitesimal. If the same amount of distortion is imposed on the second stage, as shown in Fig. 24c, the three inner springs are not stretched in the first order. As far as the elastic response is concerned, the inner springs may be ignored. At the same time the outer springs are doubled in length; so the spring constant is reduced to one-half of that for the single spring. The argument can be generalized to higher stages of the construction. The effective spring constant is scaled inversely with the length of the unit L. The same consideration applies to the D-dimensional extension of the SG. The elastic modulus is obtained from the spring constant by the ratio of length to cross-sectional area. Thus, the scaling relation for the elastic modulus is K , cc L 1 - D ,

L.4

(12.4)

-

$\ A

Y',

FIG.24. (a) Sierpinski gasket with L = 4 unit cells to a side. Each side is a simple spring with spring constant K. (b) and (c) Three possible deformationsof the L = 1 and 2 gaskets.Note that to the lowest order of deformation only the external bonds change their lengths (Bergman and Kantor94).

260

S. H. LIU

For a percolation cluster the correlation length is given by 5 a ( p - pJ", and the elastic modulus is K, a (p - pJ, we can correlate Eq. (12.4) with K , cc t-''" and write t/v = D - 1.

(12.5)

For the isotropic force model t is the same as p for the conductivity, given under Eq. (10.10)

p/v

=D -2

+ [ln(D + 3) - In(D + 1)]/ln2.

( 12.6)

Comparing these results, one finds

tlv

=

plv

= 0.737,1.585,

1,2,

for D = 2,3, for D = 2,3.

The critical exponents given by the two models are not as different as that found by Feng and Sen, apparently because the Sierpinski gasket, being built from a succession of tripods, is more rigid than a percolation cluster on a triangular lattice. Kantor and Webman realized that the percolation cluster on a square lattice can be stabilized by the addition of bond-bending elastic forces.'' This new elastic force tends to preserve the angle between a pair of bonds emanating from the same vertex. Based on this elastic-force model, they showed that the elastic threshold is identical to the conductivity threshold but the critical exponent for the elastic modulus is much higher than that for the conductivity. Consider a linear chain of N identical springs on the x axis. The nodes on the chain are at the points X I = la, where a is the length of each spring and 1 = O-N. We fixed the end at X , and apply a force F at XN in the x direction. Each spring is distorted by an amount 6 = FIK, where K is the spring constant. The displacement of the endpoint is N6; so the effective spring constant of the chain is KIN. A chain like this has no rigidity against bending unless additional restoring potential is assumed. If we denote the lateral displacement of the nodes by y,, the added potential has the form KIN-1

' = TI 1 (YI-1 = 1

-2y,+y,+A2.

(12.7)

If we hold the first link rigid and apply a shear force F' at X I in the y direction, one can find by elementary means that the distortion of the lth node is y, = (F'/K')(13- 1)/6(N - 1). The distortion at the loaded end is y N = (F'/K')N(N + 1)/6, and the effective spring constant is 6K'/N2 for large N. 95

Y.Kantor and 1. Webman, Phys. Rev. Lett. 52, 1891 (1984).

FRACTALS AND CONDENSED MATTER PHYSICS

26 1

Therefore, a sufficiently long chain is always relatively soft against transverse stress compared with longitudinal stress. The same consideration can be extended to a random chain formed by N bonds {bi} each of length a. The spring part of the elastic energy depends on the change in length {6bi} of the bonds, and the bond-bending part of the energy depends on the change of relative angle {a&}, where &i is the angle between the bonds bi and bi - 1. The total energy has the expression V = 3K’a2C6+Z i

+ i K Ci 6 b f .

(12.8)

When a force F is applied to the end of the chain, the relative changes in the lengths and orientations of the bonds can be found by minimizing W = V F (RL - RN), where the quantity in the parentheses is the displacement of the end of the chain. The component of this displacement in the direction of the force is given by

-

where 4 is a unit vector perpendicular to the plane. The minimization of W leads to

-

6bi = F bi/Ka, where Ri denotes the original position of the endpoint of the vector bi. The minimized elastic energy is given by W = F2NR:J2K’a2

+ F 2 L l/2Ka, l

(12.11)

where R: is the squared radius of gyration of the projection of the sites Ri on the direction of F x 5 (1 2.12) and (12.13) For very long chains the second term in Eq. (12.11) is very small compared with the first one. Thus, the rigidity of the twisted chain is dominated by anglebending elastic forces. The percolation network consists of multiply connected blobs linked together by singly connected chains. The chains are the weakest

262

S. H.

LIU

links in the cluster, and their stability determines the elastic threshold. Since the angle-bendingforces give elastic rigidity to the chains, the elastic threshold must be identical to the conductivity threshold. Keeping only the bond-bending part of the elastic energy, we find that the effective elastic constant, which relates the elastic energy to the displacement of the end of the chain, is given by k, = K ' a 2 / N R : .

(12.14)

Coniglio showed that the number of bonds N in the singly connected part of the percolation backbone diverges near the threshold as (p - pC)-' = 51/v.96 The radius of gyration behaves like the correlation length 5. Therefore, the bulk elastic constant is given by K , a 52-Dk, oc (-(D+l/V). (12.15) From this one finds t = Dv + 1. Putting in the known values of v, one obtains t = 3.6 in 2D and 3.55 in 3D. For D 2 6 the network consists mainly of singly connected chains. For this case v = N cc t 2 ,and R : a t2,and these give t = 4. The conductivity exponent is p = 3. Feng et al. verified some of these theoretical results by computer simulat i ~ n . The ~ ' numerical estimate for a 2D network is t = 3.3 f 0.5 for a range of ratios of K ' / K . When the ratio is large they observed a cross-over to a lower value t v for smaller systems. This result was predicted by Bergman and Kantor for pure central forces.94 The fact that the elastic critical exponent t is much larger than the conductivity exponent p has been confirmed experimentallyby B e n g ~ i g u iA. ~ ~ set of holes is punched into a metal (copper or aluminum) foil at random on the sites of a square lattice. The percolation threshold is reached from above by adding more and more holes. The conductivity and the elastic constants of the perforated sheet are measured. When the data are plotted versus ($c - $), where $ is the fraction of the surface which is removed and $c = 0.6, the critical exponents are found to be t = 3.7 f 0.2 and p = 1.2 f 0.1 for Cu and t = 3.3 f 0.2 and p = 1.1 0.1 for Al. These are in good agreement with the theory. Halperin et al. showed that if, in the above-mentioned experimental setup, the holes are punched such that the centers are not restricted on a square grid, the critical exponents will be further increased.99This is because the threshold property is determined by the weakest links in the system, which are the necks

4,

N

A. Coniglio, Phys. Rev. Lett. 46,250 (1981). S. Feng, P. N. Sen, B. I. Halperin, and C. J. Lobb, Phys. Rev. B 30,5386 (1984). 98 L. Benguigui, Phys. Rev. Lett. 53,2028 (1984). 99 B. I. Halperin, S. Feng, and P. N. Sen, Phys. Rev. Lett. 54,2391 (1985).

96

97

FRACTALS AND CONDENSED MATTER PHYSICS

263

between nearest-neighbor holes. In the lattice-percolation case the necks are of the same width, but in the continuum case there is a distribution of widths, with the narrow ones dominating the elastic critical exponent. The net result is that the elastic critical exponent is increased by 3 in 2D and 3in 3D, while the conductivity exponent remains the same in 2D but increases by 5 in 3D. The authors also pointed out that in the continuum percolation case the exponents are dependent on the model, because different models yield different probability distributions for the narrow necks.

13.

SPECTRAL

DIMENSION OF ELASTIC VIBRATIONS

Since realistic elastic forces give different values of the critical exponent from the isotropic forces, the spectral dimension of the elastic vibrations of a fractal network should also be different from the fracton dimension derived in the last Part. In this Part we will use the term fracton dimension for isotropic force systems, such as RL networks, and the term spectral dimension for elastic systems. We illustrate the problem by solving the normal modes of the linear chain with bond-stretching and bond-bending forces. A mass M is attached at every node. The longitudinal vibrations have the dispersion relation w = (K/M)l/Zak, in the region of small wave vector k. The density of states is a constant and the spectral dimension is 1. The equation of motion of the transverse modes is

(13.1) The dispersion relation is found to be w = (K‘/M)’12(ak)’. The spectral dimension of these modes is i, and the density of states has a singularity (o)-”~ at low frequencies. One can go one step further and solve the vibrational problem of a zigzag chain with right-angle bends. The longitudinal and transverse distortions are coupled, but there remain two branches of normal modes with linear and quadratic dispersion relations, respectively. Liu studied the spectral dimension of a Sierpinski gasket.”’ This system is ideal for studying the compressional waves because it is stabilized by spring forces alone. The basic building block of the system is an equilateral triangle with masses M/2 at every vertex and a spring with spring constant K on every side. When three such units are joined together to form the next stage, the interior vertices will acquire the mass M. For a large system all vertices will loo

S. H. Liu, Phys. Rev. B 30,4045 (1984); in “Kinetics of Aggregation and Gelation” (F. Family and D. P. Landau, eds.), p. 169. North-Holland Publ., Amsterdam, 1984.

264

S. H. LIU

have the same mass except the three outermost vertices, and the mass defects there have little effect on the spectrum. The elastic response of the basic unit is defined by a set of linear equations relating the external forces fi exerted on the vertices and the displacements ri of the vertices, i = 1-3. For a 2D system we arrange fi and ri as six-component column vectors, then the matrix relating them is the 6 x 6 dynamical matrix D'''(o), where D'"(o)

=

r? - ( M u 2 / 2 ) T ,

(13.2)

and r? is the elastic matrix and Fthe unit matrix. Choosing the coordinates of the vertices as (0.3). ( - 4.0). and (4.0). one finds -

1

?

0

i?=K

-a -$I4 -a - $14

-a

0

f

-$I4

-J3/4

-$ $14 3

-5

2

$14 -1

-$I4

:

_-

-3

$14

$14

3

0

0 ( 1 3.3) 0'

0

4

2

-$I4

-

-53

-1

0 0

$14

-$I4

3 4-

In the build-up process three units are linked together. The dynamical matrix of the larger unit is found by the same decimation process for the diffusion problem except that we must cancel vector forces and equate vector displacements at the decimated vertices. In the end the forces at the three new vertices are linked to the displacements of the vertices by a new dynamical , is related to the old dynamical matrix D ' O ) ( o )by a matrix matrix D ( ' ) ( o ) which equation. The mathematics is straightforward but tedious; so we will only discuss the results here. For frequencies very low compared with the characteristic frequency of the basic unit, i.e., ( K / M ) " 2 , the new dynamical matrix has approximately the same form as in Eq. (12.16) but the elastic matrix k? and the mass M are renormalized. The renormalization of k is very simple in that every element is reduced by a factor of 2. The physical reason for this scaling relation has been discussed by Bergman and K a n t ~ rThe . ~ mass ~ scaling is not simple at first. But if we iterate on the relation between D'O) and D"),we can calculate the approximate dynamical matrices for larger units in the lowfrequency limit. It is found that after the first two or three stages the effective mass scales up by a factor of 3 at every iteration. This shows that the lowfrequency excitations are collective vibrations involving all atoms in the system. The frequency spectrum scales according to ( K / M ) ' / ' , so that at every stage the frequency of a normal mode scales down by the factor 6lI2. The scaling factor for the density of states is derived from a consideration of the number of modes in a frequency interval do.At each step all modes in do

FRACTALS AND CONDENSED MATTER PHYSICS

265

are compressed into an interval dw', where w' = w/& O n the other hand, the number of atoms goes up threefold. Thus the conservation of the number of modes gives 3Nn+,(w')do' = Nn(w)dw.

(13.4)

where N,,(o) is the density of states per atom in the nth stage of construction. In the limit of large n the spectral density should be independent of n. As a result, Eq. ( 1 3.4) gives N ( w ) a ode-where d , = 2 In 3/ln 6, which is smaller than the fracton dimension d , = 2 In 3/ln 5. Consequently, the density of low-frequency modes is higher in the real elastic system, and this is consistent with the finding that the elastic network is softer than what the isotropic force model predicts. The same consideration can be applied to Sierpinski gaskets in D-dimensional space. The elastic constant scales with the length of an edge as K a L-', the mass scales as M a Ld, where d is the fractal dimension d = ln(D + l)/ln2. The frequency of a normal mode scales like o a L-(d+l)/z,and finally d , = 2d/(d 1). An analogous discussion for the tortional modes has been done by Webman and Grest'"" and applied to the percolation network. The elastic modulus of the predominantly bending mode is related to the correlation length by K a ( - ( f / v - D + Z ) , while the mass M a Yd. Thus, the frequency of a normal mode behaves like <-A, where A = ( t / v - D + d + 2)/2. We define CT = t / v D + d = ( t - P ) / v . This allows us to write, in close analogy with the diffusion problem, that

+

d , = 2d/(a

+ 2).

(13.5)

Using the estimated values of the exponents, Webman and Grest estimated that d , = 0.8 and 0.9, for percolation networks in 2D and 3D, respectively. The spectral density of these modes is singular at low frequencies. At higher frequencies the modes are predominantly compressional. A large volume of a percolation cluster behaves like a simple spring such that K cc t-'. Following this scaling relation, we find d,

= 2d/(d

+ 1).

(13.6)

For D 2 6 the compressional spectral dimension of the percolation network attains the limit 4, which agrees with the speculation of Alexander and Orba~h.~~ As pointed out in the last Part, the problem of spectral dimension was stimulated by the spin-lattice-relaxation measurements of Stapleton et al. on protein molecules.'0*' l To explain the temperature dependence of the relaxation rate, the authors proposed that the spectral dimension of a fractal looa

I. Webman and G . S. Grest, Phys. Rev. B 31, 1689 (1985).

266

S. H. LIU

system should be equal to its fractal dimension. Replacing the fractal dimension by the fracton dimension spoils the agreement between theory and experiment. Since the spectral dimension is even smaller than the fracton dimension, the formula 3 + 2 4 deviates further from the observed temperature exponent. We are left with the unsatisfactory situation that a large body ' that the of experimental data is not understood. Helman et ~ 1 . ' ~argued fractal dimension was underestimated by Stapleton et al. because they only counted the structure of the backbone. The vibrational spectrum is strongly modified by the crosslinks, which have the effect of increasing the fracton dimension, thus bringing the theory into agreement with experiment. However, the crosslinks are created by hydrogen bonds, which are an order of magnitude weaker than the carbon bonds in the backbone. Thus, it is not proper to consider the two kinds of bonds on the same footing."' Liu pointed out that the weakest spot of the entire theory is in the determination of the coupling strength between the spin and the lattice strain.'02a The standard theory for crystalline solids treats the elastic vibration as propagating waves, and the same result was carried over to fractal systems. This cannot be correct because in fractal systems the translational symmetry is absent and all modes are localized. This introduces a complicated position and energy dependence to the coupling matrix element such that it is not at all obvious whether a scaling relation is possible. Therefore,we now have a better appreciation of the complexity of the problem, but the experimental results remain unexplained. Orbach has applied the idea of fractional spectral dimension to disordered solids such as glasses and epoxy resin.68 These materials are not fractals, because their mass distribution is quite uniform once the length scales exceed the nearest-neighbor distance. Nevertheless, nonintegral spectral dimension may be possible if the modes are localized. Strictly speaking, these modes are elastic vibrations, not the fractons envisioned by Alexander and O r b a ~ h . ~ ~ However, their estimated fracton dimension 1.33 is not too far from the observed value 1.4 for epoxy resin68and the value 1.5 0.2 obtained by He and Thorpe by computer s i m ~ l a t i o n . 'In ~ ~the latter work the authors started with a uniform crystal and then randomly removed bonds. The elastic potential consists of bond-stretching and angle-bending terms of the remaining bonds. Much work needs to be done to elucidate the relation between the elastic vibrations in fractals and in glassy materials.

J. S. Helman, A. Coniglio, and C. Tsallis, Phys. Rev. Lett. 53, 1195 (1984). H. J. Stapleton, Phys. Rev. Lett. 54, 1734 (1985). S. H. Liu, Phys. Rev. B 32,6094 (1985). l o 3 H. He and M. F. Thorpe, Phys. Rev. Lett. 54,2107 (1985).

267

FRACTALS AND CONDENSED MATTER PHYSICS

VI. Magnetic Ordering on Fractal Networks There are many motivations to study the magnetic ordering of spin systems on fractal networks. From the practical point of view this is one way to approach a theory of dilute magnetic systems, i.e., a lattice in which only a fraction of sites are occupied by magnetic atoms. A more esoteric motivation was to implement the nonintegral dimension that enters the modern theory of phase transition.”* At the present time neither program has been carried to fruition, but the accumulation of theoretical insights may form an important cornerstone of the theory of random magnetic systems. 14.

MAGNETIC ORDERING OF ISING SYSTEMS

In a series of publications Gefen et al. have made an exhaustive study of magnetic ordering transitions on a number of fractal l a t t i ~ e s . ’ ~They ~-~~~ demonstrated that, unlike systems with translational symmetry, the critical behavior of Ising systems on fractals depends not only on the dimension of the fractal but also on other geometrical properties as well. We will discuss the new concepts involved with illustrative examples. The basic mathematical tool for this class of Ising problems is real-space renormalization, which is a form of the decimation process used in connection with diffusion and elastic-vibration problems. The principle of the method may be illustrated by considering the well-known 1D Ising problem. Consider two Ising spins coupled ferromagnetically with the coupling strength J. The Hamiltonian of the system is

H

= -J C T ~ C T , .

(14.1)

The system has four configurations, namely (TT), (TJ), (It),and (11).The first and last configurations have the two spins in parallel, and their probability weight is Pl = exp K, where K = J/kT. The two middle configurations have antiparallel spins, and they have the probability weight P2 = exp( - K). The partition function of the system is

Z

= 2(P1

+ Pz) = 4cosh K,

(14.2)

S. K. Ma, “Modern Theory of Critical Phenomena.” Benjamin, Reading, MA, 1976. Y. Gefen, B. B. Mandelbrot, and A. Aharony, Phys. Rev. Lett. 45,855 (1980). lo6 Y. Gefen, A. Aharony, and B. B. Mandelbrot, J. Phys. A 16, 1267 (1983). lo’ Y. Gefen, A. Aharony, Y. Shapir, and B. B. Mandelbrot, J. Phys. A 17,435 (1984). Y. Gefen, A. Aharony, and B. B. Mandelbrot, J. Phys. A 17, 1277 (1984). log Y. Gefen, Y. Meir, B. B. Mandelbrot, and A. Aharony, Phys. Rev. Lett. 50, 145 (1983). lo’

268

S. H. LIU

and the two-spin correlation function is u = (0102) = 2(P1 - P2)/Z = tanhK.

(14.3)

The quantity u is less than or equal to one, the equal sign holds only at zero temperature. We now construct a two-link chain with three spins. The Hamiltonian of the system is

H

=

-J(a102

+

02~3).

( 14.4)

The long-range-order property of the system is determined by the two-spin correlation function of the two end spins 1 and 3, and this quantity is calculated as follows. Consider the configuration for which both spins 1 and 3 are in the up state. This can be realized in two ways, with spin 2 either up or down. Therefore, the probability weight of this arrangement is P i = PI + Pi. Similarly, the configuration with both 1 and 3 in the down state has the same weight. For the two configurations with 1 and 3 in antiparallel alignment, the spin 2 must be parallel to one and antiparallel to the other. Thus, these two configurations have the weights P ; = 2P1P2. The partition function of the three spin systems is Z' = 2(P; + P ; ) = 8 cosh' K, and the two-spin correlation function of the end spins is u' = (0103) = tanh2 K = u2.

(14.5)

The procedure may be repeated to obtain the two-spin correlation function of the end spins of longer chains. In the nth step the correlation function u, is related to that of the (n - 1)th step by 2 u, = u,-1.

(14.6)

The length of the chain after n steps is L, = 2"-'. We can relate u to L, by u, = exp( - L,/<),where is the correlation length defined by

<

(14.7) < = - l/ln tanh K z exp(2K)/2. Long-range order is attained when < + and this occurs at T = 0. A low00,

temperature scaling field can be defined as t = exp( - 2K), then the scaling relation for the correlation length is 5 a t-', with v = 1. The fact that the 1D Ising system does not order at finite temperature can be understood on the basis of a simple argument due to Peierls."' Consider a ferromagnetically ordered chain of length N. The spin order may be interrupted by reversing the moments to the right of any one of the N - 1 points along the chain. The gain in internal energy is J and the gain in entropy is kln(N - 1). At the temperature T the change in free energy is A F = J R. Peierls, Proc. Cambridge Philos. SOC.32,477 (1936).

FRACTALS AND CONDENSED MATTER PHYSICS

269

k T ln(N - l), and it follows that the disordered state is favored as long as T > J / k ln(N - 1). For an infinite chain the ordering temperature approaches zero. Gefen et al. have analyzed a few Ising problems based on the Koch curve with and without branchings.lo6 The simplest case is the Ising problem with nearest-neighbor interactions along the Koch curve. The thermodynamics of the system is identical to the 1D case except that the correlation length is defined as the straight-line distance between the two end spins, rather than the length of the zigzag path along the curve. The critical temperature remains zero, and the critical exponent v is now In 3/ln 4 = l/d. Similar qualitative results are obtained for branching models, and recently Aharony et al. discussed the magnetic correlation on the percolation cluster based on these results." The Ising problem on the Sierpinski gasket is more intere~ting.'~' The basic building block of the system is an equilateral triangle with one Ising spin at every corner. There is ferromagnetic interaction between every pair of spins; so the Hamiltonian of the system is H = -J(cT~o,

+

+

0 2 ~ 3

0301).

(14.8)

There are eight spin configurations, two of which have all spins in the up or down state and the remaining six have one spin in the opposite direction from the other two. The probability weights of the first two configurations are Pl = exp(3K), while the rest have P2 = exp( - K).The partition function is Z

= 2P1

+ 6P2,

(14.9)

and the two-spin correlation function for any pair is u = (01a2) = 2(P1 - P2)/Z = (e3K - e-K)/(e3K

+ 3e -K 1.

(14.10)

Just as in the 1D case, u is less than 1 at any finite temperature. Now we build the next step by joining together three units with spins (010203), (a',a;o;), and (ayo;d!) such that a; = a2, o; = a;, and a3 = a:. We take as new spin variables the three corner spins (o,a;ay), and sum over the directions of the inner spins to obtain the probability weights of the configurations of the larger unit. It is easy to find that P; P;

'I1

+ 3PiPi + 4 P i , = P:P2 + 4PIP$ + 3P3. = P:

(14.11)

A. Aharony, Y. Gefen, and Y. Kantor, Proc. N A T O Ado. Study Inst. Scaling Phenomena Disordered Syst., 1985.

270

S. H. LIU

The two-spin correlation function of the larger unit is U' = ( c ~ c ; = ) 2(P; - P;)/(2P;

= u2/(1

- u + UZ),

+ 6P;) (14.12)

after a little algebra. The procedure can be repeated indefinitely, and we find Eq. (14.12)to be the recurrent relation for the two-spin correlation function of ever larger units. Similar to the 1D case the quantity u, for the nth unit decreases steadily with increasing n at a finite temperature so that long-range order exists only at T = 0. The fact that long-range spin order does not exist at finite temperature can be understood by a simple extension of the Peierls argument. A large portion of the SG can be isolated from the rest by cutting a few as four bonds. A lattice with this property is said to be finitely ramified. By reversing the spins of the isolated portion, the system can gain entropy without a proportionate increase in internal energy. A 3D version of the SG, or the skewed web, is a finitely ramified lattice with the fractal dimension of 2. A similar calculation shows that the critical temperature of the system is zero, even though 2 is the lower critical dimension for the Ising system. In a system without translational symmetry, the ramification of the lattice is a more important factor in determining whether a phase transition may occur at a finite temperature. The percolation cluster is finitely ramified; so some of the findings here may apply to dilute magnetic systems near the percolation limit. The SG possesses the unusual property that it is on the verge of a phase transition. To appreciate this point we rewrite Eq. (14.12)in the form un -1

= 1 - u n- l- 1

We define t = e-4K and write

+ u,!,.

(14.13)

+ 4t + 4t2 + o(t3).

(14.14)

+ 4t + [16(n - 1) + 4]t2 + 0(t3).

(14.15)

u1

=1

Iterating n times, we obtain u, = 1

The small-t expansion breaks down when the t and t Z terms in Eq. (14.15) become comparable; i.e., (n - 1) r 1/4t. The length of the side is L, = 2"-'. Thus, we can define a correlation length by

< = exp[$ In 2 exp(4K)],

<

(14.16)

in the sense that at a finite temperature T, all spins within are ordered. This result is very unusual in that at a temperature T = J/2k or K = 2, the correlation length is Assuming that the nearest-neighbor distance is 3 A, which is the typical interatomic distance, we find that the size of the correlated region is many orders of magnitude larger than the limit of the observable

FRACTALS AND CONDENSED MATTER PHYSICS 1.5

I

I

0

1

I

I

2

3

27 1

1.o

0.5

0

T/ J

FIG.25. The magnetizationcurves of an king model on the Sierpinski gasket. The number by each curve indicates the order of the gasket. The inset illustrates the small tail at the lower end (Liu"*).

universe. At a temperature around T = J/k a system of laboratory size becomes ordered without undergoing a phase transition. The magnetization curve of the system can be calculated by applying an infinitesimal magnetic field and determining the response."' The results for a set of SG of various sizes are plotted in Fig. 25. It can be seen that for systems of laboratory dimensions the magnetization grows quite abruptly at a finite temperature. The only sign that no phase transition occurs is a small tail at the lower end. The specific heat, however, is dominated by short-range effects and exhibits only a Schottky anomaly around T = J/k. The Sierpinski gasket is the simplest of the so-called hierarchical lattices invented by Berker and Ostlund."3 These systems have the common property that real-space renormalization is exact. In a series of publications, Kaufman and Griffiths have studied the topological properties of these lattices in relation to the existence of phase transition, the convexity of free energy, and the possibility of continuously varying critical exponents.' These 'I2

'I5 'I6

S. H. Liu, Phys. Rev. B 32,5804 (1985). A. N. Berker and S. Ostlund, J . Phys. C 12,4961 (1979). M. Kaufman and R. B. Griffiths, Phys. Rev. B 24,496 (1981). R. B. Griffiths and M. Kaufman, Phys. Rev. B 26,522 (1982). M. Kaufman and R. B. Griffiths, Phys. Rev. B 28,3864 (1983).

212

S. H. LIU

problems are too sophisticated for the scope of this review. Equally intriguing is the finding by McKay and Berker that the renormalization-group trajectories for the hierarchical Ising lattice shows chaotic behavior, and the authors suggested this as a model for spin-glasses.' '7 9 1 A more conventional approach to the spin-glass problem was taken by Banavar and Cieplak, who assumed a variety of random interactions between Ising spins on a Sierpinski gasket.' l 9 As pointed out earlier, even without random interactions the Ising system on SG exhibits an abrupt magnetic transition and a broad specific-heat peak, and these properties are shared by spin-glass systems. More work is needed before a definitive conclusion can be reached on whether this line of inquiry will lead to a viable model for dilute magnetic alloys. Infinitely ramified fractal spin systems have been investigated by Gefen et al. on the basis of a lattice called the Sierpinski carpet.'" In the beginning stage the lattice is a square. It is subdivided into b2 squares of which 1' are cut out. The procedure is repeated to the smaller squares until the subdivisions reach microscopic dimensions. The fractal dimension of such an object is d = ln(b2 - lZ)/lnb, which is between 1 and 2. The number of bonds that must be cut to isolate a part of the object increases with the size of the carpet without bound. As a result the spin systems are expected to order with T, # 0, even though the fractal dimension of the system is smaller than the lower critical dimension for the translationally invariant Ising system. The authors have shown that the critical property of the spin system depends on two more geometrical properties, the connectivity and the lacunarity. A part of the carpet can be isolated by making a cut across the carpet. The cut may encounter some holes, and as a result it has the general appearance of a Cantor bar. The smallest fractal dimension of such a cut is defined as the connectivity. In the construction procedure the squares that are cut out may be distributed in different ways, for instance, connected or separated. This results in holes which are either large but few or small but numerous. The quantity lacunarity is defined to describe this distribution of holes, and it measures the failure of the system to have translational invariance. An approximate renormalizationgroup scheme has been constructed and used to find the dependence of the critical temperature on these geometrical factors. In another publication Gefen et al. demonstrated that a hypothetical (D= 1 + +dimensional Ising model may be implemented on Sierpinski carpets of low l a c ~ n a r i t y . To '~~ what extent this program can be carried to higher dimensions is an open question.

'

I

S . R. McKay and A. N. Berker, J . Appl. Phys. 55, 1649 (1984). For a review of the properties of spin glasses see J. A. Mydosh, J . Magn. Magn. Mater. 7,237 (1978).

'I9

J. R. Banavar and M. Cieplak, Phys. Reo. B 28,3813 (1983).

FRACTALS AND CONDENSED MATTER PHYSICS

273

So far all the published works on magnetic problems are concerned with thermodynamical properties. The dynamical response of magnetic moments on fractal lattices should be of great interest in view of the fact that neutron scattering and spin-relaxation data are available for random magnetic systems, and these data are in general poorly understood. At the present level of research activity, one can expect that these challenges will soon be met. ACKNOWLEDGMENTS The author wishes to thank Professor H. Ehrenreich for his encouragement and support during the course of writing this review. He is indebted to participants of his short course at the Oak Ridge National Laboratory for their stimulating questions and for bringing to his attention numerous publications in this area. Informative discussions with R. Orbach, B. B. Mandelbrot, Y. Gefen, S. Alexander, L. Sander, I. Webman, D. Weitz, J. S. Huang, S. K. Sinha, J. B. Bates, S. Spooner, J. S. Lin, M. Rasolt, T. Kaplan, J. C. Wang, and K. R. Subbaswamy are gratefully acknowledged. Work sponsored by the Division of Materials Sciences, U. S. Department of Energy under contract DE-AC05-840R21400 with Martin Marietta Energy Systems, Inc.