This notebook contains additional material to the Mathematica GuideBooks  (Springer-Verlag, New York, 2004) http://www.MathematicaGuideBooks.org . It uses the stylesheet of the Mathematica GuideBooks, which is embedded.

Graphics Volume, Chapter 2

Consider 1000 random walkers on a 3D cubic lattice, each making 1000 steps. Carry out these walks and visualize the lattice points that have been visited and the ones that have not been visited.

sitesVisitedByNRandomWalkers[n_, steps_]  :=
Module[{pos, visitedQ,
        dirs = {{-1, 0, 0}, {1, 0, 0}, {0, -1, 0}, {0, 1, 0}, 
                {0, 0, -1}, {0, 0, 1}}},
       Do[pos = {0, 0, 0};
          Do[pos = pos + dirs[[Random[Integer, {1, 6}]]];
             visitedQ[pos] = True, {steps}], {n}]; 
          #[[1, 1, 1]]& /@ DownValues[visitedQ]]

We start with 100 walkers, each making 100 steps.

SeedRandom[1];
Length[visitedSites = sitesVisitedByNRandomWalkers[100, 100]] // Timing

The function makeCube returns the faces of the cube sitting at a lattice point {x,y,z}. (To avoid rational numbers, we multiply all coordinates by 2.)

makeCube[{x_, y_, z_}] :=  Map[(2{x, y, z} + #)&, 
  {{#8, #4, #2, #6}, {#8, #6, #5, #7}, {#8, #7, #3, #4}, 
   {#4, #3, #1, #2}, {#1, #3, #7, #5}, {#2, #1, #5, #6}}&[
  {-1, -1, -1}, {-1, -1, 1}, {-1, 1, -1}, {-1, 1, 1}, 
  { 1, -1, -1}, { 1, -1, 1}, { 1, 1, -1}, { 1, 1, 1}], {2}]

The outermost and the inner boundaries of the lattice points visited by the walkers are the cube faces occurring exactly once in the list of all cube faces.

makeBoundaryCubeFaces[visitedSites_] := #[[1, 2]]& /@ 
 Select[Split[Sort[{Plus @@ #/4, #}& /@ 
              Flatten[makeCube /@ visitedSites, 1]],
           #1[[1]] === #2[[1]]&], Length[#] === 1&];

boundaryCubeFaces = makeBoundaryCubeFaces[visitedSites];

To find the connected components of the surfaces, we associate the cube faces boundaryCubeFaces to all of its edges.

edges[face_] := Sort /@ Partition[Append[face, face[[1]]], 2, 1]

makeEdgesAndFaces[boundaryCubeFaces_] := 
({#[[1, 1]], Last /@ #}& /@ Split[
  Sort[Flatten[Function[f, {#, f}& /@ edges[f]] /@ 
     Take[boundaryCubeFaces, All], 1]], #1[[1]] === #2[[1]]&]);

edgesAndFaces = makeEdgesAndFaces[boundaryCubeFaces];

The function neighborFaces[edge] returns a list of all faces sharing the edge edge.

Function[{edge, faceList}, neighborFaces[edge] = faceList] @@@ edgesAndFaces;

Next, we separate the surface into its connected surface parts. The function makeConnectedSurfaces does this by selecting recursively all faces that are connected through a common edge.

makeConnectedSurfaces[] := 
Module[{faceBags = {}, counter = 0, allEdges, xMax, startEdge,
        currentFaces, newFaceBag, currentEdges, newFaces},
Off[Unset::norep];
(* as long as there are unused edges *)
While[DownValues[neighborFaces] =!= {}, 
(* all remaining edges *)
allEdges = #[[1, 1, 1]]& /@ DownValues[neighborFaces];
(* outermost edge *)
xMax = Max[Transpose[Flatten[allEdges, 1]][[1]]];
startEdge = Select[allEdges, (#[[1, 1]] === xMax || 
                              #[[2, 1]] === xMax)&, 1][[1]];
(* neighbor faces of startEdge *)
currentFaces = neighborFaces[startEdge];
neighborFaces[startEdge] =.;
newFaceBag = {};
(* extend current surface part as long as possible *)
While[
 (* the edges of the current faces *)
 currentEdges = Flatten[edges /@ currentFaces, 1] // Union;
 (* new faces on the outside *)
 newFaces = DeleteCases[Flatten[neighborFaces /@ currentEdges, 1],
                        _neighborFaces] // Union;
 newFaces =!= {},
 (* keep outer layer of faces *)
 newFaceBag = {newFaceBag, newFaces};
 (neighborFaces[#]=.)& /@ currentEdges;
 currentFaces = newFaces];
 (* keep new component *)
 faceBags = {faceBags, Union[Level[newFaceBag, {-3}]]}];
 (* return connected components *) Level[faceBags, {-4}]]

For the example under consideration, we have four components. The outermost one is clearly dominating in size.

walkSurfaces = makeConnectedSurfaces[];
{First[#], Length[#]}& /@ Split[Length /@ walkSurfaces]

Here are two graphics. The left graphic shows the outer boundary of all sites visited and the right graphic shows the not visited sites that are inside the left graphic.

(* color polygon of one surface component 
   according to distance to the origin *)
inOutColoredPolys[polys_] := 
Module[{polyDistances = {Sqrt[#.#]&[Plus @@ #/4], #}& /@ N[polys], 
        min, max},
       {min, max} = {Min[#], Max[#]}&[First /@ polyDistances];
       {SurfaceColor[#, #, 2.2]&[Hue[0.78 - 0.78 (#[[1]] - min)/(max - min)]], 
                                     Polygon[#[[2]]]}& /@ polyDistances];
                                              
(* randomly colored polygons *)
randomlyColorPolys[polys_] := 
{SurfaceColor[Hue[Random[]], Hue[Random[]], 3 Random[]], Polygon /@ polys}

showWalkSurfaces[walkSurfaces_] := 
Show[GraphicsArray[
{(* color large component according to distance to the origin *)
 Graphics3D[{EdgeForm[], inOutColoredPolys[First[walkSurfaces]]}, 
           PlotRange -> All, Boxed -> False],
 (* color inner components randomly *)
 Graphics3D[{EdgeForm[], randomlyColorPolys /@ Rest[walkSurfaces]}, 
            PlotRange -> All, Boxed -> False]}]];

showWalkSurfaces[walkSurfaces];

Now, we repeat the above calculations with 1000 walkers, each taking 1000 steps. The walks can be carried out in less than half a minute on a 2 GHz computer.

SeedRandom[1];
Length[visitedSites = sitesVisitedByNRandomWalkers[1000, 1000]] // Timing

Repeating the steps from above, we obtain now 806 connected surface components.

(* all cube faces *)
boundaryCubeFaces = makeBoundaryCubeFaces[visitedSites];
(* edges and faces *)
edgesAndFaces = makeEdgesAndFaces[boundaryCubeFaces];
(* make definitions for neighborFaces *)
Function[{edge, faceList}, neighborFaces[edge] = faceList] @@@ edgesAndFaces;
(* no longer needed *) Clear[edgesAndFaces];
walkSurfaces = makeConnectedSurfaces[];
Length[walkSurfaces]

Again, we have one large component representing the outermost visited sites and many small not visited sites in the inside. Here is a plot of the logarithm of the distribution.

ListPlot[{Length[#], Log @ First[#]}& /@ Split[Sort[Length /@ walkSurfaces]], 
         PlotRange -> All, Frame -> True, Axes -> False]

Here are the connected components.

showWalkSurfaces[walkSurfaces]

By calculating the volume and the surface area of the clusters for various number of walkers (=steps), we see that the surface roughness is approximately constant.

volumeSurfaceData = 
Table[{o, Length[#], Length[makeBoundaryCubeFaces[#]]}&[
                            sitesVisitedByNRandomWalkers[o, o]], 
      {o, 100, 1000, 100}]
      
N[Log[#3, #2]& @@@ volumeSurfaceData]

For many results on many random walkers, see [✶].