I am having some difficulties implementing a Monte Carlo simulation in order to compare the results to analytical calculations. The following code simulates a case with one random walker and a stationary trap. The results from the code below match the analytical calculations. The average walk-length for a 5x5 lattice turns out to be about 31.66, which is correct.

stepTypes = {{0, 1}, {1, 0}, {0, -1}, {-1, 0}};
pos1 = RandomInteger[{0, 4}, 2];
While[pos1 != {2, 2}, (*{2,2} is the stationary trap*)
 pos1 = Mod[pos1 + RandomChoice[stepTypes], 5];
Print["Walk length is:", index];
index = 0;

But when analyzing the same case, only replacing the trap with another random walker. The results dont match at all, they tend to be higher than the analytical and the difference increases hyperbolically as the lattice size grows.

stepTypes = {{0, 1}, {1, 0}, {0, -1}, {-1, 0}};
pos1 = RandomInteger[{0, 4}, 2];
pos2 = RandomInteger[{0, 4}, 2];
While[pos1 != pos2,
  pos1 = Mod[pos1 + RandomChoice[stepTypes], 5];
  pos2 = Mod[pos2 + RandomChoice[stepTypes], 5];

Print["Walk length is:", index];
index = 0;

By simulating the first case I was able to tell that my boundary conditions and the random number generation work correct. Not sure if the problem is in the way I synchronize the two walkers' movement? While loop tests the condition after both walkers have moved and increments the walk-length counter once both have moved as well. Any ideas of why the second case does not work? Perhaps it has to do with a subsequent random number generation?

Note: The original code runs this simulation N number of times and computes the average.


N     Walker/Trap     Two Walkers
3     8.99            8.00
4     18.33           22.23
5     31.66           26.10
6     49.24           52.67
7     71.82           56.61
  • 2
    $\begingroup$ Just skimming, but: if #1 at +5 and #2 at +4, can you end up having #1 go left and #2 go right, so that #1 at +4 and #2 at +5, but #1 and #2 not equal, so your program doesn't realize the two have collided? $\endgroup$ – user1722 Apr 18 '16 at 1:39
  • $\begingroup$ Very Good point. Not sure if that case is supposed to count as a collision, but I will check as soon as I get home. Any simple way to check for it ? $\endgroup$ – Casper Apr 18 '16 at 2:00
  • $\begingroup$ You can do two checks for equality, one after each step left/right, or you can keep the directions as variables and make sure that the pre and post position of the two walkers is the same (ie, #1 left of #2 in both cases or #1 right of #2 in both cases). If you don't count that as a collision, #1 could pass #2 and continue going to infinity, no? $\endgroup$ – user1722 Apr 18 '16 at 2:02
  • 1
    $\begingroup$ The problem is more subtle than I first thought. There only seems to be trouble with square lattices with sides having an even number of points. Limiting your code to odd lattices, I get {{5, 31.09}, {7, 78.28}, {9, 128.01}, {11, 227.43}, {13, 254.77}, {15, 288.32}} with the means being taken over 100 iterations. How does that compare to answers you expect? There was no exponential growth of run time as the lattice size increased. $\endgroup$ – m_goldberg Apr 18 '16 at 2:29
  • $\begingroup$ Yes, i knew about his issue. And the results you got are correct (with a little more iterations they would match the analytic exactly). When computing even sized lattices I ran into infinite loops due to some polarity problems, so the two walkers never meet - not sure how to fix this. Perhaps implementing a check for the case @barrycarter mentioned would fix it. $\endgroup$ – Casper Apr 18 '16 at 2:50

This is not answer, but an extremely long comment.

I find this problem very interesting, but haven't been able to solve it. In my attempts, I developed a tool to visualize the the two-walker random walk. I am posting this tool because I think it might be useful to the OP or anyone else looking this problem for exploring what's going on.

steps = {{0, 1}, {1, 0}, {0, -1}, {-1, 0}};

simulation[n_, max_] :=
  Module[{index, w1, w2, pos1, pos2},
    index = 0;
    pos1 = RandomInteger[{0, n - 1}, 2];
    pos2 = RandomInteger[{0, n - 1}, 2];
    w1 = {pos1}; w2 = {pos2};
    While[pos1 != pos2 && index <= max,
      AppendTo[w1, pos1 = Mod[pos1 + RandomChoice[steps], n]];
      AppendTo[w2, pos2 = Mod[pos2 + RandomChoice[steps], n]];
    {w1, w2}]

path[pts : {{_, _} ...}, color_?ColorQ, r_Real] := 
   {Opacity[.5], Line[Subsequences[pts, {2}]]},
   {Thick, Circle[pts[[1]], r]}, (* start marker *)
   Circle[pts[[-1]], 1.6 r]}     (* end marker *)

visualize[n_, max_: 500] :=
  Module[{w1, w2, r},
    r = .016 (n - 1);
    {w1, w2} = simulation[n, max];
         Point @ Flatten[CoordinateBoundsArray[{{0, n - 1}, {0, n - 1}}], 
       1]},                  (* lattice points *)
      {path[w1, Red, 1.4 r], (* red walker path *)
       path[w2, Blue, r]}},  (* blue walker path *)
      PlotRangePadding -> Scaled[.05]]]

Here are two interesting paths (from a code developer's point-of-view).



The above walk is interesting because the walk ended at the starting point of the blue walker and demonstrates why the marker circles are scaled the way the are.



This one is interesting because it is one of the relatively rare even lattice simulations to come to a successful conclusion.

I really hope someone will come up with an answer to this question soon because I spending too much time on and not getting much out it.

  • $\begingroup$ Thank you for this comment. I have found a solution already tho, take a look at the comments above. $\endgroup$ – Casper Apr 19 '16 at 4:42
  • $\begingroup$ I find this problem interesting, is it possible to post details of your solution @FEARxxx ? $\endgroup$ – thils Apr 20 '16 at 20:45
  • $\begingroup$ @thils the discrepancy had to do with the way the walkers were subsequently moved. You have to either move synchronously (on separately threads perhaps) or you have to make it check for a specific case where the walkers are excatly one step away and their next step is towards each other. In other words, if they swap positions then there was a collision. $\endgroup$ – Casper Apr 21 '16 at 3:00

The discrepancy had to do with the way the walkers were subsequently moved. You have to either move synchronously (on separate threads perhaps) or you have to make it check for a specific case where the walkers are exactly one step away and their next step is towards each other. In other words, if they swap positions then there was a collision.

Thanks everyone.


Here is my code for a single walker:

singleWalker[] := Module[{
    stepTypes = {{0, 1}, {1, 0}, {0, -1}, {-1, 0}},
    pos1 = RandomInteger[{0, 4}, 2]
    Mod[# + RandomChoice[stepTypes], 5] &,
    # != {2, 2} &

Doing walks = (singleWalker[] & /@ Range[1000]); and measuring N@Mean[Length /@ walks] gives (in my case) 31.013. I'm not confident I'm measuring exactly the right length (I took Rest to count the number of steps, rather than the number of places the walker stood, which is different by 1).

You can also look at Histogram[Length /@ walks].

Here's my (edited to handle the switcheroo/synchronization case as discussed in the comments) code for two walkers:

twoWalkers[] := 
  Module[{stepTypes = {{0, 1}, {1, 0}, {0, -1}, {-1, 0}}, 
          pos1 = RandomInteger[{0, 4}, 2], pos2 = RandomInteger[{0, 4}, 2]}, 
        Mod[# + {RandomChoice[stepTypes], RandomChoice[stepTypes]}, 5] &,
        {pos1, pos2}, 
        Not@Or[#2[[1]] == #2[[2]], And[#1 == Reverse[#2]]] &, 2];

Then, doing walks2 = (twoWalkers[] & /@ Range[1000]); and N@Mean[Length /@ walks2] gives about 35.47 25.8.

I think the issue is that you are measuring far too few N. Monte Carlo often requires hundreds or thousands of measurements to get a good answer, and you need to convince yourself that you're in the "infinite statistics" limit. The dramatic hopping around that you get for small N suggests that you're not taking a reliable average (ie. that your ensemble is too small).

This is related to the fact that the Histograms are very long-tailed.

  • $\begingroup$ You should have gotten 26.10 and i would say 5000 iterations would be enough. Monte Carlo sims do require a huge amount of statistical data (for it to be accurate) but as that number approaches infinity it gets closer to a specific value, so if after say 5000 iterations you are getting the same mean value, then there is no need to do more. btw as you can read in the comments above, the problem was related to the fact that the walkers didnt move synchronously. Thank you for your answer. $\endgroup$ – Casper Apr 21 '16 at 3:06
  • $\begingroup$ I just ran 100000 trajectories and found that two walkers on a 5x5 grid again to be 35.4851. $\endgroup$ – evanb Apr 21 '16 at 4:22
  • $\begingroup$ The same synchronization problem i had at the very beginning. Look at my answer. Your program does not count for some collisions. $\endgroup$ – Casper Apr 21 '16 at 12:24
  • $\begingroup$ Ah, I see. A switch-places step is a collision. $\endgroup$ – evanb Apr 21 '16 at 17:27
  • 1
    $\begingroup$ I implemented the switching case, and got just about 26. $\endgroup$ – evanb Apr 21 '16 at 17:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.