Percolation threshold estimation

August 30, 2013

The first assignment for Algorithms 1 is an estimation for the percolation threshold in an n-by-n array of sites; the system is said to percolate whenever there exists an uninterrupted series of connections between the top “edge” and the bottom one.

We need a fast lookup for the neighborhood of any given site; the neighbors() method checks whether a site is in the middle or on the sides or corners of the grid:

def neighbors(ii, n):
    "returns list of neighbors of ii in n-by-n grid"
    l = []
    if mod(ii, n) > 0:
        l.append(ii-1)
    if mod(ii-n+1, n) > 0 or ii == 0:
        l.append(ii + 1)
    if ii > n-1:
        l.append(ii - n)
    if ii < n*(n-1):
        l.append(ii + n)
    return l

I use two lists: L is a pointer array as in the previous post and O contains the state tags (‘#’ for closed, ‘.’ for open site). There are two extra sites, serving as roots for the top, resp. bottom sites.

def percolation_frac(L, O, sites):
    for ii, isite in enumerate(sites):
        O[isite] = '.'
        if isite < n:
            L[isite] = ridxTop
        elif n*(n-1) < isite <= N:
            L[isite] = ridxBot
        neighb = neighbors(isite, n)
        for inb in neighb:
            if O[inb]=='.':
                L = wqupc2(L, isite, inb)                
        if L[-1]==L[-2]: 
            """ percolation: 
            top and bottom virtual sites 
            belong to same component
            """
            pf = float(ii) / N
            return pf

As there is no analytical means of predicting when the percolation phase transition is going to occur, we just run the experiment in Monte Carlo mode (many times with randomized initialization), and average the per-run percolation threshold estimate.
Theory tells us that in a square lattice, this value lies around 0.5927 .
My implementation bumps into a few false positives, possibly due to the order in which the union operation is performed, thereby skewing the average more toward 0.6something. Need to look into this.

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: