# Percolation on a 2d grid graph

using EasyABM

## Step 1: Create Model

In this model we will work solely with the graph and won't need agents. We initially create a grid graph, and then create our model as follows.

graph = square_grid_graph(20,20);
model = create_graph_model(graph, open_nodes=Int[], percolated=false, prob=0.0)

## Step 2: Initialise the model

In the initialiser! we make the node at the upper right corner of the grid graph to be open and also give it color cl"blue" to indicate presence of water at this node. All other nodes are kept closed and their color is set to cl"black".

function initialiser!(model)
model.properties.num_nodes=length(vertices(model.graph))
for node in vertices(model.graph)
if node==model.properties.num_nodes #last node is open & has water
model.graph.nodesprops[node].open=true
model.graph.nodesprops[node].color=cl"blue"
else #all other nodes are closed
model.graph.nodesprops[node].open=false
model.graph.nodesprops[node].color=cl"black"
end
end
push!(model.properties.open_nodes, model.properties.num_nodes)
model.properties.percolated=false
end

init_model!(model, initialiser= initialiser!, props_to_record = Dict("nodes"=>Set([:color])))

## Step 3: Defining the step_rule! and running the model

In this step we implement the step logic of the percolation model in the step_rule! function and run the model for 50 steps. At each step of the simulation we make 10 randomly chosen nodes open and let water flow to them (i.e. assign them blue color) if connected to any source node (i.e. blue colored node).

function let_water_out!(node, model)
water_source_nodes=[node]
loop_condition = true
n = model.properties.num_nodes
while length(water_source_nodes)>0 && loop_condition
nd = pop!(water_source_nodes)
nbrs = neighbor_nodes(nd, model)
for nbr in nbrs
if (model.graph.nodesprops[nbr].open)&&(model.graph.nodesprops[nbr].color==cl"white")
model.graph.nodesprops[nbr].color = cl"blue"
if (nbr<=sqrt(n))
model.properties.percolated=true
model.properties.prob = length(model.properties.open_nodes)/n
loop_condition=false
break
end
push!(water_source_nodes, nbr)
end
end
end
end

function let_water_in!(node, model)
nbrs = neighbor_nodes(node, model)
n=model.properties.num_nodes
for nbr in nbrs
if model.graph.nodesprops[nbr].open && (model.graph.nodesprops[nbr].color==cl"blue")
model.graph.nodesprops[node].color = cl"blue"
if (node<=sqrt(n))
model.properties.percolated=true
model.properties.prob = length(model.properties.open_nodes)/n
return
end
break
end
end
end

function step_rule!(model)
num=0
n=model.properties.num_nodes
percolated = model.properties.percolated
#open 10 nodes per step
while (num<10) && (length(model.properties.open_nodes)<n) && (!percolated)
node = rand(1:n)
if !(model.graph.nodesprops[node].open)
num+=1
model.graph.nodesprops[node].open = true
push!(model.properties.open_nodes, node)
if (node<=n)&&(node>=(n-sqrt(n)+1))#upper most row
model.graph.nodesprops[node].color=cl"blue"
else
model.graph.nodesprops[node].color=cl"white"
let_water_in!(node, model)
end
if (model.graph.nodesprops[node].color == cl"blue") && !(model.properties.percolated)
let_water_out!(node, model)
end
percolated = model.properties.percolated
end
end
end
run_model!(model, steps = 50, step_rule = step_rule!)

## Step 4: Visualisation

In order to draw the model at a specific frame, say 4th, one can use draw_frame(model, frame = 4). If one wants to see the animation of the model run, it can be done as

animate_sim(model)

The code below calculates the average probability at which percolation occurs for an ensemble of model runs.

function calculate_percolation_probability(;grid_graph_size=20, attempts=10)
n=grid_graph_size
frames=Int(ceil(n^2/10)) # every step 10 random nodes are opened, so n^2/10 steps must be enough
percolation_probs = Float64[]
for _ in 1:attempts
graph = square_grid_graph(n,n)
model = create_graph_model(graph, graphics=false, open_nodes=Int[], percolated=false, prob=0.0) #don't need visualisation so set graphics=false
init_model!(model, initialiser= initialiser!)
run_model!(model, steps=frames, step_rule = step_rule! )
if model.properties.percolated
push!(percolation_probs, model.properties.prob)
end
end
percolation_probs
end

probs=calculate_percolation_probability(grid_graph_size=20, attempts=80);
sum(probs)/length(probs)
0.5919687499999999

## References

1. https://en.wikipedia.org/wiki/Percolation_theory