
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 



v


{\displaystyle v}

 be a vector of 



n


{\displaystyle n}

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



p


{\displaystyle p}

; the probability of a value being 0 is therefore 



1
−
p


{\displaystyle 1-p}

.
Define a run of 1s as being a group of consecutive 1s in the vector bounded
either by the limits of the vector or by a 0. Let the number of such runs in a given
vector of length 



n


{\displaystyle n}

 be 




R

n




{\displaystyle R_{n}}

.

For example, the following vector has 




R

10


=
3


{\displaystyle R_{10}=3}



[1 1 0 0 0 1 0 1 1 1]
 ^^^       ^   ^^^^^

Percolation theory states that





K
(
p
)
=

lim

n
→
∞



R

n



/

n
=
p
(
1
−
p
)


{\displaystyle K(p)=\lim _{n\to \infty }R_{n}/n=p(1-p)}


Task
Any calculation of 




R

n



/

n


{\displaystyle R_{n}/n}

 for finite 



n


{\displaystyle n}

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



t


{\displaystyle t}

 runs, where 



t
≥
100


{\displaystyle t\geq 100}

.

For values of 



p


{\displaystyle p}

 of 0.1, 0.3, 0.5, 0.7, and 0.9, show the effect of varying 



n


{\displaystyle n}


on the accuracy of simulated 



K
(
p
)


{\displaystyle K(p)}

.

Show your output here.

See also
 s-Run on Wolfram mathworld.
#lang racket
(require racket/fixnum)
(define t (make-parameter 100))
 
(define (Rn v)
  (define (inner-Rn rv idx b-1)
    (define b (fxvector-ref v idx))
    (define rv+ (if (and (= b 1) (= b-1 0)) (add1 rv) rv))
    (if (zero? idx) rv+ (inner-Rn rv+ (sub1 idx) b)))
  (inner-Rn 0 (sub1 (fxvector-length v)) 0))
 
(define ((make-random-bit-vector p) n)
  (for/fxvector
   #:length n ((i n))
   (if (<= (random) p) 1 0)))
 
(define (Rn/n l->p n) (/ (Rn (l->p n)) n))
 
(for ((p (in-list '(1/10 3/10 1/2 7/10 9/10))))
  (define l->p (make-random-bit-vector p))
  (define Kp (* p (- 1 p)))
  (printf "p = ~a\tK(p) =\t~a\t~a~%" p Kp (real->decimal-string Kp 4))
  (for ((n (in-list '(10 100 1000 10000))))
    (define sum-Rn/n (for/sum ((i (in-range (t)))) (Rn/n l->p n)))
    (define sum-Rn/n/t (/ sum-Rn/n (t)))
    (printf "mean(R_~a/~a) =\t~a\t~a~%"
            n n sum-Rn/n/t (real->decimal-string sum-Rn/n/t 4)))
  (newline))
 
(module+ test
  (require rackunit)
  (check-eq? (Rn (fxvector 1 1 0 0 0 1 0 1 1 1)) 3))