
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

Given an 



M
×
N


{\displaystyle M\times N}

 rectangular array of cells numbered 




c
e
l
l

[
0..
M
−
1
,
0..
N
−
1
]


{\displaystyle \mathrm {cell} [0..M-1,0..N-1]}

assume 



M


{\displaystyle M}

 is horizontal and 



N


{\displaystyle N}

 is downwards.

Assume that the probability of any cell being filled is a constant 



p


{\displaystyle p}

 where

 



0.0
≤
p
≤
1.0


{\displaystyle 0.0\leq p\leq 1.0}


The task

Simulate creating the array of cells with probability 



p


{\displaystyle p}

 and then
testing if there is a route through adjacent filled cells from any on row 



0


{\displaystyle 0}

 to any on row 



N


{\displaystyle N}

, i.e. testing for site percolation.

Given 



p


{\displaystyle p}

 repeat the percolation 



t


{\displaystyle t}

 times to estimate the proportion of times that the fluid can percolate to the bottom for any given 



p


{\displaystyle p}

.

Show how the probability of percolating through the random grid changes with 



p


{\displaystyle p}

 going from 



0.0


{\displaystyle 0.0}

 to 



1.0


{\displaystyle 1.0}

 in 



0.1


{\displaystyle 0.1}

 increments and with the number of repetitions to estimate the fraction at any given 



p


{\displaystyle p}

 as 



t
>=
100


{\displaystyle t>=100}

.

Use an 



M
=
15
,
N
=
15


{\displaystyle M=15,N=15}

 grid of cells for all cases.

Optionally depict a percolation through a cell grid graphically.

Show all output on this page.

 
! loosely translated from python.
! compilation: gfortran -Wall -std=f2008 thisfile.f08
 
!$ a=site && gfortran -o $a -g -O0 -Wall -std=f2008 $a.f08 && $a
!100 trials per
!Fill Fraction goal(%)    simulated through paths(%)
!           0                          0
!          10                          0
!          20                          0
!          30                          0
!          40                          0
!          50                          6
!
!
!    b b b   b   h   j     m m m
!  b b   b b b   h h   m m m m m
!    b b b       h h h     m    
!    b     h h h   h h h h      
!  b b   h h   h h h h   h h h  
!  b b b   h h h   h h h h h h h
!  b b   @   h   h   h h h h h  
!      @ @       h h h h h h h h
!    @ @ @ @       h h   h   h  
!      @ @ @ @       h h h h h h
!      @ @ @   h h h h     h h h
!  @ @ @     h h   h   h     h h
!    @       h         h h h h h
!  @     h h   h     h h h     h
!  @ @   h h h h h h h   h h   h
!          60                         59
!          70                         97
!          80                        100
!          90                        100
!         100                        100
 
program percolation_site
  implicit none
  integer, parameter :: m=15,n=15,t=100
  !integer, parameter :: m=2,n=2,t=8
  integer(kind=1), dimension(m, n) :: grid
  real :: p
  integer :: i, ip, trial, successes
  logical :: success, unseen, q
  data unseen/.true./
  write(6,'(i3,a11)') t,' trials per'
  write(6,'(a21,a30)') 'Fill Fraction goal(%)','simulated through paths(%)'
  do ip=0, 10
     p = ip/10.0
     successes = 0
     do trial = 1, t
        call newgrid(grid, p)
        success = .false.
        do i=1, m
           q = walk(grid, i)    ! deliberately compute all paths
           success = success .or. q
        end do
        if ((ip == 6) .and. unseen) then
           call display(grid)
           unseen = .false.
        end if
        successes = successes + merge(1, 0, success)
     end do
     write(6,'(9x,i3,24x,i3)')ip*10,nint(100*real(successes)/real(t))
  end do
 
contains
 
  logical function walk(grid, start)
    integer(kind=1), dimension(m,n), intent(inout) :: grid
    integer, intent(in) :: start
    walk = rwalk(grid, 1, start, int(start+1,1))
  end function walk
 
  recursive function rwalk(grid, i, j, k) result(through)
    logical :: through
    integer(kind=1), dimension(m,n), intent(inout) :: grid
    integer, intent(in) :: i, j
    integer(kind=1), intent(in) :: k
    logical, dimension(4) :: q
    !out of bounds
    through = .false.
    if (i < 1) return
    if (m < i) return
    if (j < 1) return
    if (n < j) return
    !visited or non-pore
    if (1_1 /= grid(i, j)) return
    !update grid and recurse with neighbors.  deny 'shortcircuit' evaluation
    grid(i, j) = k
    q(1) = rwalk(grid,i+0,j+1,k)
    q(2) = rwalk(grid,i+0,j-1,k)
    q(3) = rwalk(grid,i+1,j+0,k)
    q(4) = rwalk(grid,i-1,j+0,k)
    !newly discovered outlet
    through = (i == m) .or. any(q)
  end function rwalk
 
  subroutine newgrid(grid, probability)
    implicit none
    real :: probability
    integer(kind=1), dimension(m,n), intent(out) :: grid
    real, dimension(m,n) :: harvest
    call random_number(harvest)
    grid = merge(1_1, 0_1, harvest < probability)
  end subroutine newgrid
 
  subroutine display(grid)
    integer(kind=1), dimension(m,n), intent(in) :: grid
    integer :: i, j, k, L
    character(len=n*2) :: lineout
    write(6,'(/)')
    lineout = ' '
    do i=1,m
       do j=1,n
          k = j+j
          L = grid(i,j)+1
          lineout(k:k) = ' @abcdefghijklmnopqrstuvwxyz'(L:L)
       end do
       write(6,*) lineout
    end do
  end subroutine display
 
end program percolation_site
 