# Random walk on a Sierpinski gasket

I am trying to simulate a random walk on a Sierpinski gasket. The best strategy i could come up with is to use Nearest point function to determine the next possible step of my walker. But this creates an issue at n = 2 generation where the are two unwanted steps the walker shouldn't be able to go to. In other words, it shouldn't be able to cross the white triangle in the middle. There should always be 4 possible steps (except for the corners where there is only 2). Any ideas on how this case can be handled?

sierpinski[{a_, b_, c_}] := With[{ab = (a + b)/2, bc = (b + c)/2,
ca = (a + c)/2}, {{a, ab, ca}, {ab, b, bc}, {ca, bc, c}}];
pts = {{0, 0}, {1, 0}, {1/2, Sqrt[3]/2}};
n = 2;
d = Nest[Join @@ sierpinski /@ # &, {pts}, n];
list = Flatten[d, 1];
list = DeleteDuplicates[list];
point = RandomChoice[list];
next = Nearest[DeleteCases[list, point], point];
Graphics[{EdgeForm@Black, Polygon@d, Red, PointSize[0.025], Point[point], Blue, Point[next]}]

• Why not try making it a graph?
– Wjx
Jun 22, 2016 at 3:09
• @Wjx How would that solve the problem tho? Jun 22, 2016 at 3:12
• I've already got an elegant solution, I simply need a few coding time~
– Wjx
Jun 22, 2016 at 3:25
• There's elegant mathematical property in this figure!
– Wjx
Jun 22, 2016 at 3:55
• I think i got it. G = GraphData[{"Sierpinski", 3}]; and then walkers position is simply n = VertexList[NeighborhoodGraph[G, pos]]; so position = RandomChoice from that list Jun 22, 2016 at 4:27

Well, facinating question! Code first:

Clear["*"];

n = 5;
steps = 1000;
dir = {{0, 0}, {1, 0}, {0, 1}};

start = FromDigits[#, 2] & /@ Transpose@Append[RandomChoice[dir, n - 1], RandomInteger[1, 2]];

move[pt_] :=With[{c = RandomChoice@Pick[-dir, BitAnd @@ (# + pt) & /@ -dir, 0]}, pt + c +RandomChoice@DeleteCases[dir, -c]];
move[{2^n, 0}] := {2^n - 1, RandomInteger@1};
move[{0, 2^n}] := {RandomInteger@1, 2^n - 1};
move[{0, 0}] := RandomChoice@Rest@dir;

ListLinePlot[{#1 + #2/2, #2*Sqrt@3/2} & @@@NestList[move, start, steps],
PlotRange -> {{0, 2^n}, {0, 2^n*Sqrt[3]/2}}, AspectRatio -> Sqrt[3]/2]


This code will generate the following figure:

But how can this code be so elegant?(Well, is it?------Let me assume it is~) The key is how this pattern can be generated.

First is about points! Points, if exist on this figure, its base 2 expression of x and y coordinate after throw away the last digit must satisfy the following:

For one corresponding digit, there shouldn't be two 1s.

So the starting point could be determined.

Similarly, you can find the pattern of moving action------Simply calculate the mid point of the move and the rule will be the same!

So, the basic idea will be using BitAnd and a few simple list manipulation will do the rest!

Hope this can help~

n = 5;
steps = 10000;

G = GraphData[{"Sierpinski", n}];
start = RandomChoice@VertexList@G;
proc[pos_] :=
RandomChoice@DeleteCases[VertexList@NeighborhoodGraph[G, pos], pos];
Nest[proc, start, steps]; // AbsoluteTiming


This method can generate a result for n<6 because other wise the GraphData won't tell you the graph~

# Edit of Method 2

Try change this part:

DeleteCases[VertexList@NeighborhoodGraph[G, pos], pos]


to

VertexInComponents[G,pos,1]


Then it will be even faster than my first solution if n is small~

BTW: I checked about the efficiency of the two methods. The latter one takes much longer time(10*or even greater), I think that's because NeighbourhoodGraph is time consuming and Graph approch to this problem, afterall, is not that direct.

• yeah, Union[VertexInComponent[G, pos, 1], VertexOutComponent[G, pos, 1]]; returns identical results but much faster. Jun 22, 2016 at 15:42
• yep~That's even faster than my original solution!
– Wjx
Jun 22, 2016 at 15:47
• I think just VertexInComponent[G, pos, 1] is enough but Union[VertexInComponent[G, pos, 1], VertexOutComponent[G, pos, 1]],or maybe I miss something.
– yode
Jun 22, 2016 at 17:49
• @yode Well, I didn't check this graph before I made the change. Apologize~ If this graph is without direction, yes, this will be enough.
– Wjx
Jun 22, 2016 at 23:38
• A slight simplification: ListLinePlot[NestList[move, start, steps].{{1, 0}, {1, Sqrt[3]}/2}/2^n, AspectRatio -> Automatic, PlotRange -> {{0, 1}, {0, Sqrt[3]/2}}]`. Jun 23, 2016 at 1:06