02+Lecture+Balance+Cond

=Global balance and incompressible flow= An alternative formulation of the global balance condition derived last week, now expressed in terms of "flow"

math P(a \to a) + P (a \to b) + P(a \to c) = \pi(a) \quad \text{"Everything has to enter"} math

math P(a \to a) + P(b\to a) + P(c\to a) = \pi(a)\quad\text{"Everything has to exit"} math

with the flow defined as

math P(a \to b) = \pi(a) p(a \to b) math    //π// ( 1 ) = 1  /   2      //π// ( 2 ) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1  <span class="mo" style="font-family: MathJax_Main;">/   <span class="mn" style="font-family: MathJax_Main;">2    <span style="clip: rect(1.712em,1000em,3.081em,-0.535em); left: 0em; position: absolute; top: -4.668em;">  <span class="mi" style="font-family: MathJax_Math;">//π// <span class="mo" style="font-family: MathJax_Main;">( <span class="mn" style="font-family: MathJax_Main;">4 <span class="mo" style="font-family: MathJax_Main;">) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1  <span class="mo" style="font-family: MathJax_Main;">/   <span class="mn" style="font-family: MathJax_Main;">2 <span class="mo" style="font-family: MathJax_Main;">)    <span style="clip: rect(1.712em,1000em,3.081em,-0.535em); left: 0em; position: absolute; top: -3.268em;">  <span class="mi" style="font-family: MathJax_Math;">//π// <span class="mo" style="font-family: MathJax_Main;">( <span class="mn" style="font-family: MathJax_Main;">5 <span class="mo" style="font-family: MathJax_Main;">) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1    <span style="clip: rect(1.712em,1000em,3.081em,-0.535em); left: 0em; position: absolute; top: -1.868em;">  <span class="mi" style="font-family: MathJax_Math;">//π// <span class="mo" style="font-family: MathJax_Main;">( <span class="mn" style="font-family: MathJax_Main;">6 <span class="mo" style="font-family: MathJax_Main;">) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1  <span class="mo" style="font-family: MathJax_Main;">/   <span class="mn" style="font-family: MathJax_Main;">2    <span style="clip: rect(1.712em,1000em,3.081em,-0.535em); left: 0em; position: absolute; top: -0.468em;">  <span class="mi" style="font-family: MathJax_Math;">//π// <span class="mo" style="font-family: MathJax_Main;">( <span class="mn" style="font-family: MathJax_Main;">8 <span class="mo" style="font-family: MathJax_Main;">) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1  <span class="mo" style="font-family: MathJax_Main;">/   <span class="mn" style="font-family: MathJax_Main;">2    <span style="clip: rect(1.712em,1000em,3.081em,-0.535em); left: 0em; position: absolute; top: 0.932em;">  <span class="mi" style="font-family: MathJax_Math;">//π// <span class="mo" style="font-family: MathJax_Main;">( <span class="mn" style="font-family: MathJax_Main;">9 <span class="mo" style="font-family: MathJax_Main;">) <span class="mo" style="font-family: MathJax_Main; padding-left: 0.278em;">= <span class="mn" style="font-family: MathJax_Main; padding-left: 0.278em;">1  <span class="mo" style="font-family: MathJax_Main;">/   <span class="mn" style="font-family: MathJax_Main;">2    <span style="clip: rect(1.768em,1000em,2.853em,-0.554em); left: 0em; position: absolute; top: 2.276em;">  <span class="mtext" style="font-family: MathJax_Main;">all others

Application of the balance condition (schema)
Here is the construction of a valid Monte Carlo algorithm for the homogeneous 3x3 game (each arrow weighs 1/2).

We will show later how to make this ad-hoc method into a systematic method, for the inhomogeneous 3x3 pebble game.

Application of the balance condition (Python code and graphics output)
code format="python" import random, pylab

moves = {1: [(2, 0.5)], 2: [(3, 1.0)], 3: [(6, 1.0)], 4: [(1, 0.5), (5, 1.0)], 5: [(2, 0.5), (8, 1.0)], 6: [(5, 0.5), (9, 1.0)], 7: [(4, 1.0)], 8: [(7, 1.0)], 9: [(8, 0.5)]} all_pos = [] pos = 9 for iter in range(5000000): Upsilon = random.random for k in moves[pos]: if Upsilon < k[1]: pos = k[0] break all_pos.append(pos) pylab.figure(1) pylab.hist(all_pos, bins=9, range=(0.5,9.5), normed=True) x_vec = range(1, 10) y_vec = [1.0/9.0 for k in range(9)] pylab.plot(x_vec,y_vec,'ro') pylab.title('incompressible 3x3_homogeneous') pylab.savefig('incompressible_3x3.png') pylab.show code



Balance condition (shaving algorithm, schema)
To apply the incompressible-flow approach to a more complicated system, as shown below, we may incrementally take out (shave off) loops one after the other. For example, we may take out a first loop along the sites 1-2-3-6-9-8-7-1, which brings us to a square with decremented external weights (as shown in the lecture).

Balance condition (shaving algorithm, Python code and graphics output)
code format="python" import random, pylab

moves = {1: [(2, 1.0/4.0)], 2: [(3, 1.0)], 3: [(6, 7.0/11.0)], 4: [(1, 1.0)], 5: [(2, 1.0/3.0), (6, 1.0)], 6:[(9, 4.0/9.0), (3, 1.0)],      7: [(4, 1.0/17.0), (8, 4.0/17.0)], 8: [(7, 4.0/7.0), (5, 1.0)],       9: [(8, 4.0/6.0)]} all_pos = [] pos = 9 for iter in range(500000): Upsilon = random.random for k in moves[pos]: if Upsilon < k[1]: pos = k[0] break all_pos.append(pos) pylab.figure(1) pylab.hist(all_pos,bins=9,range=(0.5,9.5),normed=True) x_vec=range(1, 10) y_vec=[4.0, 2.0, 11.0, 1.0, 3.0, 9.0, 17.0, 7.0, 6.0] y_vec=[a / 60. for a in y_vec] pylab.plot(x_vec, y_vec, 'ro') pylab.title('Histogram inhomogeneous global 3x3) pylab.savefig('inhomogeneous_global_3x3.png') pylab.show code



=Lifting=

"Lifting", an extremely attractive idea that takes the concept of loops and incompressible flows one step farther: One draws the loops, but tries to remain on them as much as allowed (by the global balance condition). In order to lift the loops one from the other, one introduces an additional variable in the definition of a Markov chain. This is trivial to implement on a ring of //L// sites with uniform probabilities. The gain is spectacular, as one can now visit all sites in a time linear in //L//.

Lifting in one dimension (schema)
To implement the lifting algorithm on an //L//-site ring (with general probabilities), we need to respect the incompressibility condition at each site. This is illustrated with the below schema

Lifting in one dimension (Python code and graphics output)
Here, an implemented Python program with the lifting idea. The probability distribution is taken as a random variable. Try to comment out the seed function (the line with random.seed('smac')), in order to understand its purpose. code format="python" import random, pylab, math

random.seed('smac') Ntime = 100000000 L = 8 Pi = [random.random for i in range(L)] all_pos = [] k = 0 dir = 1 for itime in range(10000000): k_p = (k + 1) % L   p_p = min(1, Pi[k_p] / Pi[k]) k_m = (k - 1) % L   p_m = min(1, Pi[k_m] / Pi[k]) Upsilon = random.random if dir == 1 : if Upsilon < p_p: k = k_p elif Upsilon < p_m: dir = -1 elif dir == -1 : if Upsilon < p_m: k = k_m elif Upsilon < p_p: dir = 1 if itime % 100 == 0 : all_pos.append(k) pylab.hist(all_pos, bins=L,range=(-0.5,L-0.5),normed=True) pylab.title('Lifting algorithm in one d ' ) pylab.xlabel('Site') pylab.ylabel('Occupation probability') Z = sum(Pi[k] for k in range(L)) grid = [k for k in range(L)] Boltzmann = [Pi[k] / Z for k in range(L)] pylab.plot(grid, Boltzmann, 'ro') pylab.savefig('Lifting_1d.png') pylab.show code

Metropolis-Hastings algorithm
The Metropolis-Hastings algorithm is again concerned with the detailed balance condition, but explicitly splitting the probability to move from configuration i to configuration j into an a priori probability to propose the move and a probability to accept them. We could have used this concept already in the discussion of the 3x3 pebble game, where is equal to for neighboring sites and zero otherwise. Here, we discuss it on the heliport.

See SMAC sect. 1.1.6. Understanding of a priori probabilities will be essential in later lectures. =References=