May 4, 2013

Random Walk

I use Mathematica to play with 2D Brownian Motion a little bit. According to Smoluchowski and Einstein's theory, Brownian Motion is described by random walk (RW). There are several methods to generate random walk in Mathematica, such as the use of Accumulate + RandomChoice, which is an example in the documentation:
Graphics[
  Line[
    Accumulate[
      RandomChoice[{{-1, 0}, {1, 0}, {0, 1}, {0, -1}}, 1000]
    ]
  ]
]
Fig. 1,  Random walk on a 2D square lattice
Fig. 2, a random walk with 10^4 steps on the square lattice

Fig. 3, calculation of the fractal dimension using box counting method. $a$ is the size of the boxes, the vertical values counts the number of boxes that cover the random walk points. It gives a fractal dimension of ~1.602.
The mean distance after $N$ steps is proportional to $\sqrt{N}$. This is known as the Einstein's law. Fig. 4 shows the test of this law.

Fig. 4, the mean distance. The nonlinear fit shows that the mean distance agree with the Einstein's law.


Similarly, we can work out a random walk on a triangular lattice.
Graphics[
  Line[Accumulate[
    RandomChoice[{{-1, 0}, {1, 0}, {1/2, Sqrt[3]/2}, {-1/2, 
       Sqrt[3]/2}, {1/2, -Sqrt[3]/2}, {-1/2, -Sqrt[3]/2}}, 1000]]]]

Fig. 5, Random walk on a 2D triangular lattice


In order to incorporate more sophisticated random walks, I use Nest and NestList.

 
 v = 1;
 rw[r_] := r + v * Through[{Cos,Sin}[RandomReal[{0,2*Pi}]]];
 urw[r_] := r + RandomReal[{0,2*v}] * Through[{Cos,Sin}[RandomReal[{0,2*Pi}]]];
 mbrw[r_] := r + RandomVariate[MaxwellDistribution[v]]* Through[{Cos,Sin}[RandomReal[{0,2*Pi}]]];

 n = 10^5;
 Show[
  Graphics[{Thickness[Large],Circle[{0, 0}, Sqrt[n]], 
            Arrow[{{0, 0}, Sqrt[n] {Cos[Pi/4], Sin[Pi/4]}}], 
            Text[Style[Sqrt[N], 20], 0.5 Sqrt[n] {Cos[Pi/4], Sin[Pi/4]}],
            PointSize[0.01], Point[{0, 0}],
            Red, Opacity[0.5], Thin, Line[NestList[rw,{0,0},n]],
            Green, Opacity[0.5], Thin, Line[NestList[urw,{0,0},n]],
            Yellow, Opacity[0.5],Thin, Line[NestList[mbrw,{0,0},n]]}
  ],
 ImageSize -> 400]
 

Fig. 6, Colorlines: Red, random walk; Green, uniform distribution random walk; Yellow, Maxwell-Boltzmann random walk.

Fig. 7, RW in a circle.

Fig. 8, RW in a circle, accelerated in radial direction.
Now one can study the quantitative behavior of Brownian Motion, following the historical path of Langevin, Smoluchowski and Einstein. I'll put that part aside.

The next question is what will happen if self-interactions are added. One of these self-interactions is self-avoid (SA). In SARW, the new walk try to avoid its history. Depending on how to interpret the "history", one may avoid either the "path" or the "points". I am more interesting in the points here, because points may represent the "system configurations" (system states/phase space points/field configuration/functionals), whereas the path represents another "field" such as the gauge field/gauge links etc. Avoiding points is the simplest self-avoid self-interaction.

To implement SA, the whole list has to be stored.

 greedymbsarw[rn_] := Module[{new , eps = 1, v = 1},

  (* WriteString["stdout", Length[rn], " "];*)
  new = rn[[-1]] + RandomVariate[MaxwellDistribution[v]]* Through[{Cos, Sin}[RandomReal[{0, 2 Pi}]]];
  While[
   Or @@ Table[Norm[new - pt] < eps, {pt, rn}],
   new = rn[[-1]] + 
      RandomVariate[MaxwellDistribution[v]]* Through[{Cos, Sin}[RandomReal[{0, 2 Pi}]]];
   ];
  Append[rn, new]
  ];

  lazymbsarw[rn_] := 
    Module[{new , eps = 1, v = 1},
      new = rn[[-1]] + RandomVariate[MaxwellDistribution[v]]* Through[{Cos, Sin}[RandomReal[{0, 2 Pi}]]];
      If[ And @@ Table[Norm[new - pt] > eps, {pt, rn}], Append[rn, new], rn]  
    ];

 n = 10^5;
 Show[
  Graphics[{Thickness[Large],Circle[{0, 0}, Sqrt[n]], 
            Arrow[{{0, 0}, Sqrt[n] {Cos[Pi/4], Sin[Pi/4]}}], 
            Text[Style[Sqrt[N], 20], 0.5 Sqrt[n] {Cos[Pi/4], Sin[Pi/4]}],
            PointSize[0.01], Point[{0, 0}],
            Red, Opacity[0.5], Thin, Line[NestList[mbrw,{0,0},n]],
            Green, Opacity[0.5], Thin, Line[Nest[lazymbsarw,{{0,0}},n]],
            Blue, Opacity[0.5], Thin, Line[Nest[greedymbsarw, {{0, 0}}, n]]}
  ],
 ImageSize -> 400]

In these two implementations, greedy version will loop till a new legal walk is found; whereas the lazy version stations if new legal walk if found at the first proposal. It is important to use Maxwell-Boltzmann (or other global ) RW to avoid dead-lock. Note that for the greedy version may take a very long time. One thus may also add a upper bound for trials.
Colorlines: Red, Maxwell-Boltzmann RW; Green, lazy SARW; Blue, greedy SARW.
Again, I'll leave the discussion of quantitative feature later.

We can play with more interactions.





update:
OK, I spent sometime to workout the random self-avoid walk, with respect to avoiding the path instead of the points. Essentially, one only need to compute if two segments cross. The function to do that is,

IntersectQ[x1_, y1_, x2_, y2_, x1p_, y1p_, x2p_, y2p_] := 
 Module[{j1, j2, det}, 
  j1 = (x1 + x2 - x1p - x2p)*(y1p - y2p) - (y1 + y2 - y1p - 
       y2p)*(x1p - x2p);
  j2 = (x1 + x2 - x1p - x2p)*(y1 - y2) - (y1 + y2 - y1p - y2p)*(x1 - 
       x2);
  det = (x1 - x2)*(y1p - y2p) - (y1 - y2)*(x1p - x2p);
  Abs[j1] < Abs[det] && Abs[j2] < Abs[det]]

As one may check, $t=J_1/\det, s=J_2/\det$ is the solution (parameters) of the intersection of two lines. The explicit coordinate is $(x, y) = (\frac{x_1+x_2}{2}+t\frac{x_2-x_1}{2}, \frac{y_1+y_2}{2}+t\frac{y_2-y_1}{2})$. If $-1 < t < 1, -1 < s < 1$ then the segments intersect.

lazyrsaw[rn_] := Module[{new, eps = 1, l}, 
   l = Length[rn];
   new = rn[[-1]] + 
     RandomVariate[MaxwellDistribution[1]]*
      Through[{Cos, Sin}[RandomReal[{0, 2 Pi}]]];
   If[And @@ 
     Table[! (IntersectQ[rn[[-1, 1]], rn[[-1, 2]], new[[1]], new[[2]], 
         rn[[i, 1]], rn[[i, 2]], rn[[i + 1, 1]], rn[[i + 1, 2]]]), {i,
        1, l - 1}], Append[rn, new], rn
    ]
 ];

random SAW with lazy algorithm
It's clear that the RW easily gets stuck in the corner, though possibility exist for it to jump out of it. It's better to take a look at the typical behavior of random SAW:

$\bar{N}$ is the average steps of the RSAWs.




.

No comments:

Post a Comment