Evacuation simulation using Cellular Automata (part 1)

In this image, we have a neural network with one layer and 3 neurons. The neural network gets the data from the input layer and transforms the data according to the structure defined.

Each of these arrows and neurons represents an operation in the network. Usually, these operations are multiplication and sum between the data and weights.

These weights are calculated using an optimization algorithm called backpropagation. Backpropagation is a method that exploits the chain rule of calculus to obtain the best neural network weights.

In this article, we will optimize our neural network without backpropagation. Instead, we will apply a bio-inspired algorithm Particle Swarm Optimization.

Cellular Automata — brief overview

A cellular automaton consists of a regular grid of cells, such that each cell can be in one of a finite number of states, and the grid itself can be in any finite number of dimensions.

For each cell, a set of cells called its neighborhood is defined relative to the specified cell. An initial state (time t = 0) is selected by assigning a state for each cell.

A new generation is created (advancing t by 1), according to some update rule that determines the new state of each cell in terms of the current state of the cell and the states of the cells in its neighborhood.

The rule for updating the state of cells can be either deterministic or based on some probability distribution, and in most cases, it stays the same for the entire simulation.

If the update rule is not deterministic, we have a stochastic cellular automaton. The order in which the update rule is applied is extremely important, and is most cases it is applied to the whole grid simultaneously, which is known as a parallel or synchronous update. 

In this case, the state of every cell in the model is updated together, before any of the new states influence other cells. Every cell at step t+1 sees the same “frozen picture” of the grid state at the previous step t.

In contrast, an asynchronous update is able to update individual cells independently, in such a way that the new state of a cell can immediately affect the calculation of states in its neighbouring cells. It is important not to confuse update rule with update order. 

The update rule is the function that takes as input the current cell state and the current cell neighborhood, and outputs a new state to update the current cell. The update order is just the cell order in which this function is applied (synchronously or using some asynchronous strategy).

Given two identical CAs (having the same dimension, update rule and initial state), the choice of a synchronous or asynchronous update, in most cases, will produce a completely different CA behavior.

Which one to choose will depend on the simulation goals and assumptions. The synchronous approach assumes the presence of a global clock to ensure all cells are updated together. While convenient for some systems, this might be an unrealistic assumption if the model is intended to represent, for example, a living system where there is no evidence of the presence of such a device.

Another point of attention when implementing a CA is what happens at the boundaries of the grid. If, for example, we are using a Moore neighborhood, what is the neighborhood for cells at the edges?

One way to solve this is to let these cells have a constant state. Another method is to define neighborhoods differently for these cells. As in the case of the update order, the method chosen will affect the CA behavior.

Cellular automata is a huge topic, with plenty of books and papers on the subject. If you are new to CAs, a good starting point is the Wikipedia article and this Stanford page on the subject.

CA evacuation model

The obvious choice here is a 2D CA, which represents a room where each cell can be either empty, occupied by a person, occupied by an obstacle, or be an exit cell. He have just defined the set of possible states for a cell and the dimension of our CA. Now comes the hard part: how to define the update rule and update order?

Since we are dealing with people moving inside a room, by “update rule” we really mean, how people move in respect to each other and possible obstacles? Let’s first start with the easier rules, which are very suitable for our simulation purposes:

• Obstacle cells never change their state and cannot be occupied.
• Exit cells never change their state.
• Each empty or exit cell can be occupied by at most one individual.

The rules for the actual motion of people are a bit more involved and depend on the update order chosen. For our simulation we are going to use a parallel (aka synchronous) update rule.

This means that at very iteration, people may try to move to the same location, which goes against rule 3 above. To solve this, we will resolve these conflicts in a clever way.

Now, how do people choose their move? In an evacuation situation, the goal is always the closest exit. So our move algorithm for a person occupying cell C is as follows:

In other words, we look around the current cell and choose only the cells that we can move to (empty ones or some exit) that are closer or at a equal distance to some exit than the cell we currently are.

Then we just choose the move that gets us closest to some exit. What this means is that people will try to move to a better spot if they can, but will also try to move even when their distance to some exit stays the same with the movement. Movements that get a person farther away from some exit are never considered.

To calculate the distance cost from cell A to the nearest exit E we should take obstacles into account. The Dijkstra algorithm will handle that for us. Each cell will have an associated cost value indicating how farther away it is from the nearest exit. This map is sometimes called a floor field.

Two things are important to note: first, diagonal movements have a cost of 2, horizontal or vertical movements cost 1. Second, we are not considering people as obstacles. 

That would be like solving an ever-changing maze, which makes no sense. More importantly, this allows us to calculate a Dijkstra map once and reuse it for the entire simulation!

One could argue that when choosing an exit, the congestion around the exits is an important factor. This is certainly true, it’s a kind of a dilemma: do I go for the nearest exit, which is crowded, or do I go towards some farther away exit, which is less crowded?

There are plenty of papers that implement this behavior, but in our model we choose not to to keep it simple. For example, in this paper they use parameter α (in the range [0, 1]) which reflects the degree of impatience of people during the evacuation.

When α = 0 only the distance factor is considered, when α = 1 only congestion is considered, and for 0 < α < 1, these factors are weighted in order to choose a target exit. Let’s recap our rules of movement:

• Obstacle cells never change their state and cannot be occupied.
• Exit cells never change their state.
• Each empty or exit cell can be occupied by at most one individual.
• People will try to move to the nearest empty cell that is closer or at a equal distance from some exit than the cell they currently are.
• If there is no available cell that fits the rules above, no movement will take place.

Let’s analyze two cases below. In both figures we see a yellow cell and its moore neighborhood. Grey cells are obstacles and E marks the exit. Cells with number are empty and free to move.

Suppose the yellow cell is the current cell and there is an individual there trying to move. The numbers are the distance to the exit, and we see that only the values marked in green are less than or equal to 5. In this case, the cell marked with 4 is the best move.

In Figure 2 we see another situation. The red cross is occupied by another person, we cannot move there. In this case the green cell marked with 5, despite not improving our distance towards the exit, will be chosen.

One thing is missing in our implementation: how do we handle conflicts?Since we are doing a parallel update, more than one individual might try to move to the same empty cell.

To resolve this, we simply don’t perform the move right away, instead we keep the “candidate moves” for all cells in a dictionary.

After all candidate moves have been parsed, we iterate this dictionary, and if there is more than one candidate for a cell, we simply randomly choose one that will be eligible to move, and update the grid accordingly.

Right now, our only source of randomness is the way in which we solve move conflicts. To make things a little bit more realistic let’s add a variable called panic_prob, which represents the probability that a person will not move, even if they can.

We could model this variable by drawing from a Beta distribution, so that each person would have its own panic probability, but to simplify things, we will just set it to a constant value. We have all the information we need to implement this, so let’s code!

Click: CA evacuation model code

First, we define a CACell class to represent our cell state. Notice that besides the state, we also have a variable called properties, which we’ll use to add some properties to the cell (like panic probability, for example).

Then we define the CARoom class, which contains all our simulation logic. We define the room dimensions and set each cell to be a square of size 0.4m.

As each cell can be occupied by a at most one person, this is the value most used in the literature, giving a density of approximately 6 persons per m2. We also allocate space for our cells, create a reference to the exit cells and initialize the dijkstra map to None.

Besides the usual get and set stuff, we have some interesting methods. The get_graph method loops through the cells and connects the empty ones, setting their squared distance as the connection weight (the real distance does not matter).

The get_dijkistra_map method uses the graph generated in get_graph and Dijkstra’s algorithm to assign a cost to every empty cell, representing the distance towards the nearest exit.

The step_parallel method is where things happen. It represents a single step of the simulation. We basically loop through each cell and when we find a cell with the state CACell.PERSON_STATE, we try to move to the nearest empty cell that will get us closer to the nearest exit, as explained in the previous section.

Notice that before moving, we flip a coin to check for a panic situation, in which case there is no movement. Also, notice that all moves are stored in a variable called new_state, indexed by the cell.

We later iterate new_state and resolve every cell that has two or more candidate moves, picking one at random. Finally, while resolving these conflicts, we also remove from the simulation every person that has reached an exit cell and return the total number of evacuees in this step.

We end part 1 of this series testing what we have done so far: we create a simple 10m x 10m room with external walls, a centered square-shaped obstacle and a single exit. Then for each empty cell, we add a person with a probability given by full_factor. Then we call step_parallel until everyone is evacuated.

Notice that since each cell is a 0.4m x 0.4m square and we have external walls and a 4 x 4 cells obstacle, there is room for approximately 23 * 23 -16 = 513 people. Given a full_factor of 0.3, we have an expected value of ~154 evacuees. We also implemented a simple __str__ method, so that you can see what is happening in the console.

Our Dijkstra implementation is very short:

To run everything, just add CA.py and dijkstra.py to the same folder and type in:

python3 CA.py

In part 2 we are going to setup some experiments and run them using the multiprocessing module, to speed things up. Also, it is important to be able to visualize the simulations, so pygame will help us with that. See you later!