diffusion limited aggregation simulation (a re-posting)
- To: mathgroup at yoda.physics.unc.edu
- Subject: diffusion limited aggregation simulation (a re-posting)
- From: gaylord at ux1.cso.uiuc.edu
- Date: Sat, 20 Jun 1992 09:12:39 -0500
apologies to everyone. my repost of my post didn't have a subject so this
is a repost of the repost.
==================
i've heard of some interest in monte carlo simulations in physics and i
have a bunch which i am finally beginning to publish in various places.
below is the code for (a) diffusion-limited aggregation and (b) determining
its fractal dimensions.
DLA refers to the the growth of clusters by the sticking-together of
colliding particles.
DLA's are a model for crystal growth and snow flake formation.
the full paper on DLA with graphics output and complete (extensive) text
is in the "Mathematica in Education" newsletter, volume 1, number 3 Spring
1992
[contact wellin at sonoma.edu for information on this publication which costs
$15 for 4 issues (good value)].
i can post this or even the entire paper to Math Source if someone will
tell me how to do it without two tons of paperwork. (btw -has MathSource
set up its ftp site yet?).
other models I will be putting out on the net over the course of the
summer:
spreading phenomena [epidemics, forest fires, percolation, fluid flow
through porous media, gelation].
sandpile dynamics {really neat-looking}.
etc., etc.
richard
btw -
some of this code is a little old and so while i have tried to make it
'functional' it still retains procedural vestiges (eg., Do) and outdated M
stuff (Block should be Module).
comments and suggestions are welcomed.
oh, one other thing: my friend Bill Tyndall wrote some of the programs with
me (he's co-author of the full article). if there's any part of the code
you don't like, Bill probably wrote it.
============================
=======================
Diffusion-Limited Aggregation (* see Scientific American, 256(1), 94-100
(1987) article by Sanders *)
aggregate growth starts from a single occupied 'seed' site.
a particle begins a rw from a randomly chosen location along the perimeter
of a circle of radius, r, centered on the seed.
the rw continues until it is too far outside of the circle (in which case
the walk is terminated) or it lands on a site next to the aggregate, in
which case the particle joins the structure.
after a given particle completes its trajectory, another particle is
emitted from the circle and moves by the same rules.
the process continues until a large dla containing n sites is created.
the program:
notes:
a two-dimensional square lattice is used.
the seed is located at {0, 0}.
each rw'er starts from the lattice location closest to the randomly
selected site on the circumference of a circle of radiur, r
(notes: r can be made to vary with the size of the dla and the rw's step
length can be proportional to its distance from the perimeter).
the outer 'terminating' boundary is a square of length (r + m).
the graphics output uses Rectangle, colored so that a site added to the
aggregate at a later time are lighter than a site added earlier (note: Text
can also be used to number each Rectangle in the order in which it is
added)
Clear[dLA]
dLA[r_, n_Integer, m_Integer] :=
Block[{angle, loc, particleCount = 0,
stepChoices = {{1, 0}, {0, 1}, {-1, 0}, {0, -1}}, structure = {}},
occSites = {{0,0}};
While[Length[occSites] < n,
angle = Random[Real, {0, N[2 Pi]}];
loc = Floor[{r Cos[angle], r Sin[angle]}];
++particleCount;
While[Abs[loc[[1]]] < r + m && Abs[loc[[2]]] < r + m &&
!MemberQ[occSites, loc] &&
Length[Intersection[occSites,
Map[(# + loc)&, stepChoices]]] == 0,
loc += stepChoices[[Random[Integer, {1, 4}]]] ];
If[Length[Intersection[occSites,
Map[# + loc &, stepChoices]]] =!= 0,
AppendTo[occSites, loc] ] ];
Print["The number of particles released was ", particleCount];
Print["The size of the dla is ", Length[occSites]];
Do[AppendTo[structure,
{RGBColor[0.1, 0.9 i/Length[occSites], 0.52],
Rectangle[occSites[[i]] - {0.5, 0.5},
occSites[[i]] + {0.5, 0.5}],
RGBColor[1, 1, 1],
Text[i, occSites[[i]]]}],
{i, 1, Length[occSites]}];
Show[Graphics[structure], Axes->None, AspectRatio ->1,
PlotRange->{{-(r +1), (r +1)}, {-(r + 1), (r + 1)}}]
]
----------------------------------
SeedRandom[0]
Timing[dLA[5, 35, 5]]
----------------------------
the fractal dimension of an object is a measure of the "compactness" of its
structure (ie., how much space is filled by the object).
as growth of a DLA proceeds and the object becomes larger, an initially
smooth aggregate or cluster becomes more irregularly shaped and hence its
fractal dimension decreases (eg., in a dla, there is a screening effect
which increases the likelihood that a rw'er contacts an exposed exterior
portion of the dla before it penetrates into the more shielded interior
portions).
the program:
a list of the occupied sites comprising the object is used.
a list "fracList" is constructed consisting of ordered pairs, each having
a second element which is the length of a square, 2L and a first element
which is the number density of occupied sites within that square (ie., the
number of occupied sites within a square running between -L and L in each
direction, divided by the total number of sites in the square).
the fractal dimension, df, of the object is determined using <<Graphics.m
and Fit[N[Map[Log, fractalDataList]], {1, x}, x]
(note: df varies with the size of the object so that df should be
calculated by fitting various parts of fracList (using Drop[fracList,...])
which are chosen by a visual inspection of LogLogListPlot[fracList]).
<<Graphics.m
Clear[fractalDim]
fractalDim[occSites_List] :=
Block[{latticeDensity},
latDensity[t_Integer] :=
N[Count[occSites, {x_?(Abs[#] <= t &), y_?(Abs[#] <= t &)}]]/(2t +1)^2;
fracList = Table[{2s, latDensity[s]}, {s, Max[Abs[Flatten[occSites]]]];
fracdim = Fit[ N[Map[Log, fractalDataList]], {1, x}, x];
Print["The fractal dimension of is ", Coefficient[fracdim, x]];
LogLogListPlot[fractalDataList, PlotStyle->{PointSize[.045],
RGBColor[0.763, 0.244, 0.251]}] ]
fractalDim[occSites]
"if you're not programming functionally, then you must be programming
dysfunctionally"
"is there a pattern to this madness ?"
richard j. gaylord, university of illinois, gaylord at ux1.cso.uiuc.edu
"if you're not programming functionally, then you must be programming
dysfunctionally"
"is there a pattern to this madness ?"
richard j. gaylord, university of illinois, gaylord at ux1.cso.uiuc.edu