# 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/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 '16 at 3:09
• @Wjx How would that solve the problem tho? – Casper Jun 22 '16 at 3:12
• I've already got an elegant solution, I simply need a few coding time~ – Wjx Jun 22 '16 at 3:25
• There's elegant mathematical property in this figure! – Wjx Jun 22 '16 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 – Casper Jun 22 '16 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/2}}, AspectRatio -> Sqrt/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. – Casper Jun 22 '16 at 15:42
• yep~That's even faster than my original solution! – Wjx Jun 22 '16 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 '16 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 '16 at 23:38
• A slight simplification: ListLinePlot[NestList[move, start, steps].{{1, 0}, {1, Sqrt}/2}/2^n, AspectRatio -> Automatic, PlotRange -> {{0, 1}, {0, Sqrt/2}}]`. – J. M.'s technical difficulties Jun 23 '16 at 1:06