Here's some work being done by Ric Sherlock and Devon McCormick to come up with a good Rosetta Code entry for the "Brownian Tree" problem - another name for Diffusion-Limited Aggregation.

Here's an initial attempt with a few helper functions:

```NB.* brownianTree3.ijs: make diffusion-limited aggregation.

NB.* brownianTree3: aggregate points in cluster by Brownian motion.
brownianTree3=: 3 : 0
'world npts'=. y           NB. Boolean mat world, # pts to insert.
if. 0=#world do. sz=. 100  NB. Default to 100x100 world
peri=. 2#3-~sz         NB. Track periphery of cluster
world=. 1 (<<.-:peri)} (2#sz)\$0 end.  NB. Seed center
peri=. 2#3-~sz=. #world
for_pt. i. npts do.        NB. Random pt on periphery,
xy=. >:?peri           NB.  not on populated cell or edge.
while. -.+./world{~<"1 nn=. getNbrs xy do. NB. No populated neighbor?
xy=. nn {~ ?#nn                        NB. Move to random neighbor
if. +./ xy (1&>@[ , >) peri do.        NB.  near periphery.
xy=. >:?peri
end.
end.
world=. (1) (<xy)}world
end.
NB.EG brownianTree3 '';100
)

NB.* getNbrs: get neighbors of a cell.
getNbrs=: (0 0 -.~ <: 3 3#: i.9) +"1 ]
NB.* pad0s: pad rows and cols of numeric mat with x zeroes in each direction.
NB.* minmaxrow: indexes of min and max occupied (=1) cell of boolean mat.
minmaxrow=: [:(([:<./0{"1]) , [:>./1{"1]) [:>(<_ __)-.~ [:(<./,>./) &.>[:I.&.><"1
NB.* minmaxrc: min and max occupied (=1) row and column of boolean mat.
minmaxrc=: minmaxrow,minmaxrow@:|:
NB.* trim0s: trim all-zero rows and columns from boolean mat.
trim0s=: [: (] #"1~ 0 +./ .~:])] #~ 0 +./ .~:"1 ]
NB.* density: portion of occupied cells in boolean mat.
density=: ([: +/,) % [:#,```

One problem is that this is very inefficient code for a large, empty matrix. That's why the arguments are an existing matrix and a number of points: this allows us to start with a small world so that short random walks will encounter the cluster, then increase the size of the matrix by padding the borders as the cluster grows. This means a typical use of this code will be something like the following.

```   sz=. +/,w3 [ tm=. 6!:2 'w3=. brownianTree3 '''';100'  NB. Start tracking number of cells and run-times.
sz,:tm
99
2.6285197
viewmat w3```

This gives something like this:

Next, aggregate more cells to the existing cluster:

```   sz=. sz,+/,w3 [ tm=. tm,6!:2 'w3=. brownianTree3 w3;500'
sz,:tm
99       568
2.4815702 3.0286292```

Notice how much more efficient is the second round since we have a larger target, hence shorter random walks. This gives us something like this:

Are we close to the borders of our world yet?

```   minmaxrc w3
21 82
19 81```

No, this looks like a sufficient margin to continue adding more points. Let's incorporate this check as part of our routine:

```   minmaxrc w3 [ sz=. sz,+/,w3 [ tm=. tm,6!:2 'w3=. brownianTree3 w3;1000'
8 88
10 90
sz,:tm
99       568      1392
2.4815702 3.0286292 1.0529608```

Now it looks like this:

Since we're approaching the borders of our world, let's make it larger the next time we add more points:

```   minmaxrc w3 [ sz=. sz,+/,w3 [ tm=. tm,6!:2 'w3=. brownianTree3 (50 pad0s w3);2000'
sz,:tm
33 165
29 168
sz,:tm
99       568      1392      3187
2.4815702 3.0286292 1.0529608 21.477968```

We notice a drop in efficiency as we increased the empty area by padding the borders to give us a 200x200 world.

We have enough free border space to add another bunch of points without further padding for now:

```   minmaxrc w3 [ sz=. sz,+/,w3 [ tm=. tm,6!:2 'w3=. brownianTree3 w3;2000'
21 178
20 177
sz,:tm
99       568      1392      3187      4882
2.4815702 3.0286292 1.0529608 21.477968 6.5040638```

The aggregate now looks like this:

If we continue growing this cluster in this fashion, we'll see some aesthetic flaws in our eventual result. Here's a sample of a 500x500 world:

This looks too "filled-in" in the center. Also, it's assuming a suspiciously square outline - possibly because random, independent (X,Y) selection fails to pick points uniformly on a circle? Or maybe it has do with some other aspect of how we choose the random "release points" for the new particles.

DevonMcCormick/DLA03 (last edited 2010-05-05 17:02:04 by DevonMcCormick)