
Percolation Simulation

This is a simulation of aspects of mathematical percolation theory.

For other percolation simulations, see Category:Percolation Simulations, or:
1D finite grid simulation
Mean run density 
2D finite grid simulations
Site percolation | Bond percolation | Mean cluster density

Let 



c


{\displaystyle c}

 be a 2D boolean square matrix of 



n
×
n


{\displaystyle n\times n}

 values of either 1 or 0 where the
probability of any value being 1 is 



p


{\displaystyle p}

, (and of 0 is therefore 



1
−
p


{\displaystyle 1-p}

).
We define a cluster of 1's as being a group of 1's connected vertically or
horizontally (i.e., using the Von Neumann neighborhood rule) and bounded by either 



0


{\displaystyle 0}

 or by the limits of the matrix.
Let the number of such clusters in such a randomly constructed matrix be 




C

n




{\displaystyle C_{n}}

.

Percolation theory states that 



K
(
p
)


{\displaystyle K(p)}

 (the mean cluster density) will satisfy 



K
(
p
)
=

C

n



/


n

2




{\displaystyle K(p)=C_{n}/n^{2}}

 as 



n


{\displaystyle n}

 tends to infinity. For 



p
=
0.5


{\displaystyle p=0.5}

, 



K
(
p
)


{\displaystyle K(p)}

 is found numerically to approximate 



0.065770


{\displaystyle 0.065770}

...

Task
Show the effect of varying 



n


{\displaystyle n}

 on the accuracy of simulated 



K
(
p
)


{\displaystyle K(p)}

 for 



p
=
0.5


{\displaystyle p=0.5}

 and
for values of 



n


{\displaystyle n}

 up to at least 



1000


{\displaystyle 1000}

.
Any calculation of 




C

n




{\displaystyle C_{n}}

 for finite 



n


{\displaystyle n}

 is subject to randomness, so an approximation should be
computed as the average of 



t


{\displaystyle t}

 runs, where 



t


{\displaystyle t}

 ≥ 



5


{\displaystyle 5}

.

For extra credit, graphically show clusters in a 



15
×
15


{\displaystyle 15\times 15}

, 



p
=
0.5


{\displaystyle p=0.5}

 grid.

Show your output here.

See also
 s-Cluster on Wolfram mathworld.
// version 1.2.10
 
import java.util.Random
 
val rand = Random()
const val RAND_MAX = 32767
 
lateinit var map: IntArray
var w = 0
var ww = 0
 
const val ALPHA = "+.ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz"
const val ALEN = ALPHA.length - 3
 
fun makeMap(p: Double) {
    val thresh = (p * RAND_MAX).toInt()
    ww = w * w
    var i = ww
    map = IntArray(i)
    while (i-- != 0) {
        val r = rand.nextInt(RAND_MAX + 1)
        if (r < thresh) map[i] = -1
    }
}
 
fun showCluster() {
    var k = 0
    for (i in 0 until w) {
        for (j in 0 until w) {
            val s = map[k++]
            val c = if (s < ALEN) ALPHA[1 + s] else '?'
            print(" $c")
        }
        println()
    }
}
 
fun recur(x: Int, v: Int) {
    if ((x in 0 until ww) && map[x] == -1) {
        map[x] = v
        recur(x - w, v)
        recur(x - 1, v)
        recur(x + 1, v)
        recur(x + w, v)
    }
}
 
fun countClusters(): Int {
    var cls = 0
    for (i in 0 until ww) {
        if (map[i] != -1) continue
        recur(i, ++cls)
    }
    return cls
}
 
fun tests(n: Int, p: Double): Double {
    var k = 0.0
    for (i in 0 until n) {
        makeMap(p)
        k += countClusters().toDouble() / ww
    }
    return k / n
}
 
fun main(args: Array<String>) {
    w = 15
    makeMap(0.5)
    val cls = countClusters()
    println("width = 15, p = 0.5, $cls clusters:")
    showCluster()
 
    println("\np = 0.5, iter = 5:")
    w = 1 shl 2
    while (w <= 1 shl 13) {
        val t = tests(5, 0.5)
        println("%5d %9.6f".format(w, t))
        w = w shl 1
    }
}