Lecture 26: Introduction to spatial models through lattices (4 days)

(previous, next)



Percolation

Thought experiment: sand in a funnel, with water on top. (dribble castles at the beach)

Flow-rate depends on how much free space there is left...

Site Percolation on a random square lattice

Simplest version of site percolation: Suppose we have square tiles, randomly black and white.

In canopy...

p = 0.2
n = 20
A = floor(rand(n,n) + p)
spy(A)

We can use this kind of code to create a menagery of random media by tuning the occupency probability \(p\) and the size \(n\).

[Show code]
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
#!/usr/bin/env python2

from scipy import *
from matplotlib.pyplot import *

for n in [16,32]:
    figure(1)
    clf()
    i = 1
    for p in [0.1, 0.3, 0.4, 0.5, 0.6, 0.9]:
        A = floor(rand(n,n)+p)
        subplot(2,3, i)
        i += 1
        spy(A)
        xlabel('p=%.2f'%p, fontsize=18)
        gca().get_xaxis().set_ticks([])
        gca().get_yaxis().set_ticks([])
    savefig('siteperc_n=%02d.jpg'%n, bbox_inches = 'tight')

n=64
i=1
for p in [0.3, 0.4, 0.5]:
    A = floor(rand(n,n)+p)
    subplot(1,3, i)
    i += 1
    spy(A)
    xlabel('p=%.2f'%p, fontsize=18)
    gca().get_xaxis().set_ticks([])
    gca().get_yaxis().set_ticks([])

savefig('siteperc_n=%02d.jpg'%n, bbox_inches = 'tight')

n = 16 cells

n = 32 cells

n = 64 cells

How far can you walk on white tiles without jumping? We can do this by hand, or with a recursive algorithm.

[Show code]
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
#!/usr/bin/env python

from scipy import *

def printbinarymatrix(A):
    m,n = A.shape
    s = '\n'.join([ \
            ''.join([( ' ' if i < 0.5 else '*' ) for i in row]) for row in A])
    s = '\n'.join(['-'*(n+2),s,'-'*(n+2)])
    print s
    return

def calc_percolation_depth(A, i0=None):
    """
    Given a matrix A of {0,1}'s, find the farthest across the
    matrix we can walk on 0-entries, when
    using a von Neumann neighborhood, if we start
    at the left side in row i0.

    If i0 == None, return the maximum over all initial rows.
    """

    upper_bound = A.shape[1]
    if i0 == None:
        return max([calc_percolation_depth(A,i) for i in xrange(upper_bound)])

    neighbors = [(1,0), (-1,0), (0,1), (0,-1)]

    # to track which sites we have visited ...
    visited = 0*A
    def next(i,j):
        if i < 0 or i >= A.shape[0]:
            return 0
        if j >= upper_bound:
            return upper_bound
        if j < 0 or j > A.shape[1]:
            return 0
        if visited[i,j] == 1:
            return 0
        visited[i,j] = 1
        if A[i,j] > 0.5:
            return j
        return max(j,max([next(i+dx,j+dy) for (dx,dy) in neighbors]))
    return next(i0,0)

def test():
    p = 0.45
    N = 12
    A = floor(rand(N,N) + p)
    assert set([0,1]) == set(map(int,A.flatten()))
    print "Here we have a %dx%d matrix, with fill fraction %f"%(N,N,p)
    printbinarymatrix(A)
    print "If we start on the left side, how far toward the right can we move on blanks?"
    for i in range(N):
        d = calc_percolation_depth(A,i)
        print "%2d steps on row %d"%(d,i)
    print "Max depth: ", calc_percolation_depth(A)
    print "Now, double-check with your own eyes.  Are all of these correct?"

test()
+------------+
|*  * * ** * |
|***   **   *|
|    *  *  * |
|*** *   *** |
|  **    *** |
|*    *   ***|
| *  * * *  *|
|   *** *   *|
|** ******   |
|**   ***   *|
| *** *  *** |
|* ** ***  * |
+------------+
If we start on the left side, how far toward the right can we move on 0's?
 0 steps on row 0
 0 steps on row 1
 6 steps on row 2
 0 steps on row 3
 9 steps on row 4
 0 steps on row 5
 9 steps on row 6
 9 steps on row 7
 0 steps on row 8
 0 steps on row 9
 1 steps on row 10
 0 steps on row 11
Max depth:  9
Now, double-check with your own eyes.  Are all of these correct?

For what p would the water stop flowing?

Percolation threshold for site-percolation on a von-Neumann neighborhood... p_c = 0.5927... (no exact formula is known)

Bond percolation on a square lattice

Schelling's community segregation model

Conceptual model

Visual implementation at Case2014

How do we define a "neighborhood" in two or three dimensions?

Conway's Game of Life

Spatial rewrite-rules

Cellular automata created by Stanislaw Ulam and John von Neumann in the 1940's at LANL.

Wolfram's 1-d automata

Rule enumeration method

compact state space -> must eventually repeat. (what will the period be?)

Example: ??

Classifications - period (all will have finite period) - steady-states * attractors independent of initial condition? * attractors fix parts of initial condition (error-correcting) - amphichiral oscillators (period 2, 4, ...) - shifting patterns (period corresponds to domain harmonics) - long-period acyclics (automata chaos)

Conway's game of life - rewritting in 2 dimensions

square lattice, each cite is on/off, and has 8 neighbors.

Totalistic -- depends only on the number of neighbors, not the configuration

The rules of the game are as follows ...

Discrete-state, discrete-time, like a Markov chain, but deterministic (no randomness!).

Solutions ..

Patterns are arbitrarily complex, in the sense that they are Turing complete?

Code inspection of implementation in python of Conway's "game of life"

[Show code]
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
"""
A simple and fast implementation of conway's
game-of-life using convolutions

Bug: Currently does not play well with Canopy.
    Have not figured out why.

Timothy Reluga
2016.03.25
"""

import numpy as np
import time
from scipy.signal import convolve, fftconvolve
import matplotlib.pyplot as plt

#print plt.get_backend()

def conway(state):
    """
    Conway's game of life state transition
    """
    kernel = np.array([[1,1,1],[1,0,1],[1,1,1]])

    # computes sums around each pixel
    #sums = convolve(state, kernel, mode='same')  # slow normal convolutions
    sums = fftconvolve(state, kernel, mode='same').round() # fast convolutions

    newstate = np.zeros(sums.shape)
    newstate[np.where((sums == 2) & (state == 1))] = 1
    newstate[np.where((sums == 3) & (state == 1))] = 1
    newstate[np.where((sums == 3) & (state == 0))] = 1

    # return new state
    return newstate

def main():
    p = 0.1
    sleeptime = 0.05

    # set up board
    m,n = 150,200
    A = np.floor(p+np.random.random(m*n)).reshape((m, n))

    # plot each frame
    plt.ion()
    plt.figure(1, figsize=(9,9))
    img_plot = plt.imshow(A, interpolation="nearest", cmap = plt.cm.gray)
    plt.show(block=False)
    while True:
        A = conway(A)
        img_plot.set_data(A)
        plt.draw()
        time.sleep(sleeptime)

if __name__ == "__main__":
    main()

SmoothLife video from youtube

Open problem: Chiral rules?

Must be translationally symmetric, but not chirally symmetric (i.e. the direction of rotation matters!)

Lattice--Boltzmann models of electrical signal transmission in the heart

Your browser doesn't support HTML5 Canvas.

Understanding Navier--Stokes in terms of cellular automata

The Navier-Stokes equations are of central importance in applied-mathematics -- such importance that a prize of 1 million dollars has been offered to anybody who can show that they fully describe the motion of a fluid, and never break-down in finite time.

HHP and FHP models of fluids and Navier stokes as lattice gasses ( see paper notes)