The Gray-Scott Model: A Turing Pattern Cellular Automaton
Adding reactions to our diffusion automaton
Now that we have established a cellular automaton for coarse-grained particle diffusion, we will add to it the three reactions that we introduced in the previous lesson, which are reproduced below.
- A “feed” reaction in which new A particles are fed into the system at a constant rate.
- A “death” reaction in which B particles are removed from the system at a rate proportional to their current concentration.
- A “reproduction” reaction A + 2B → 3B.
STOP: How might we incorporate these reactions into our automaton?
First, we have the feed reaction, which takes place at a rate f. It is tempting to simply add some constant value f to the concentration of A in each cell in each time step. However, if [A] were close to 1, then adding f to it could cause [A] to exceed 1, which we may wish to avoid.
Instead, if a cell has current concentration [A], then we will add f(1-[A]) to this cell’s concentration of A particles. For example, if [A] is equal to 0.01, then we will add 0.99f to the cell because the current concentration is low. If [A] is equal to 0.8, then we will only add 0.2f to the current concentration of A particles.
Second, we consider the death reaction of B particles, which takes place at rate k. Recall that k is proportional to the current concentration of B particles. As a result, we will subtract k · [B] from the current concentration of B particles.
Third, we have the reproduction reaction A + 2B → 3B, which takes place at a rate r. The higher the concentration of A and B, the more this reaction will take place. Furthermore, because we need two B particles in order for the collision to occur, the reaction should be more rare if we have a low concentration of B than if we have a low concentration of A. To model this situation, if a given cell is represented by the concentrations ([A], [B]), then we will subtract r · [A] · [B]2 from the concentration of A and add r · [A] · [B]2 to the concentration of B in the next time step.
We now just need to combine these reactions with diffusion. Say that as the result of diffusion, the change in its concentrations are ΔA and ΔB, where a negative number represents particles leaving the cell, and a positive number represents particles entering the cell. Then in the next time step, the particle concentrations [A]new and [B]new are given by the following equations:
[A]new = [A] + ΔA + f(1-[A]) - r · [A] · [B]2;
[B]new = [B] + ΔB - k · [B] + r · [A] · [B]2.
Applying these reaction-diffusion computations over all cells in parallel and over many generations constitutes a cellular automaton called the Gray-Scott model.1
Before continuing, let us consider an example of how a single cell might update its concentration of both particle types as a result of reaction and diffusion. Say that we have the following hypothetical parameter values:
- dA = 0.2;
- dB = 0.1;
- f = 0.3;
- k = 0.4;
- r = 1 (the value typically always used in the Gray-Scott model).
Furthermore, say that our cell has the current concentrations ([A], [B]) = (0.7, 0.5). Then as a result of diffusion, the cell’s concentration of A will decrease by 0.7 · dA = 0.14, and its concentration of B will decrease by 0.5 · dB = 0.05. It will also receive particles from neighboring cells; for example, say that it receives an increase to its concentration of A by 0.08 and an increase to its concentration of B by 0.06 as the result of diffusion from neighbors. Therefore, the net concentration changes due to diffusion are ΔA = 0.08 - 0.14 = -0.06, and ΔB = 0.06 - 0.05 = 0.01.
Now we will consider the three reactions. The feed reaction will cause the cell’s concentration of A to increase by (1 - [A]) · f = 0.09. The kill reaction will cause its concentration of B to decrease by k · [B] = 0.2. And the reproduction reaction will mean that the concentration of A decreases by [A] · [B]2 = 0.175, with the concentration of B increasing by the same amount.
As the result of all these processes, we update the concentrations of A and B to the following values ([A]new, [B]new) in the next time step according to our equations above.
[A]new = 0.7 - 0.06 + 0.09 - 0.175 = 0.555
[B]new = 0.5 + 0.01 - 0.2 + 0.175 = 0.485
We are now ready to implement the Gray-Scott model in the following tutorial. The question is: even though we have built a coarser-grained simulation than the previous lesson, will we still see Turing patterns?
Reflection on the Gray-Scott model
In contrast to the particle-based simulator introduced earlier, the Gray-Scott model produced an animation in under a minute on a laptop. We show the results of this model in the videos that follow. Throughout these animations, we use the parameters dA = 1.0, dB = 0.5, and r = 1, and we color each cell according to its value of [B]/([A]+[B]) using the “Spectral” color map.
Our first video shows an animation of the Gray-Scott model using the parameters f = 0.034 and k = 0.095. We use a comparable initial configuration of the automaton as in the diffusion example, in which a cluster of B particles are found in a board full of A particles.
If we expand the size of the simulation and add multiple clusters of B particles to the automaton, then the patterns become more complex as waves of B particles collide.
If we keep the feed rate constant and increase the kill rate slightly to k = 0.097, then the patterns change significantly into spots.
If we make the A particles a little happier as well, increasing f to 0.038 and k to 0.099, then we have a different striped pattern.
And if we increase f to 0.042 and k to 0.101, then again we see spots.
The point is that very slight changes in our model’s parameters can produce drastically different results in terms of the patterns that we witness. In this prologue’s conclusion, we will connect this observation back to our original motivation of identifying the cause for animal skin patterns.
-
P. Gray and S.K. Scott, Autocatalytic reactions in the isothermal, continuous stirred tank reactor: isolas and other forms of multistability, Chemical Engineering Science 38 (1983) 29-43. ↩