Untangling Complex Systems, page 59
tions are always present and they may be emphasized in far-from-equilibrium conditions, each RD
process creates a unique pattern. In this respect, the RD phenomena are a form of artistic expression
that has two authors: the person who selects the reagents, prepare the hydrogel stamp and the dry
gel, and nature that lays down the colors and creates the ultimate structures.
9.10.2 reacTion-diffusion Processes in Technology
The Reaction-Diffusion processes are also useful for the technological development in two fields,
at least. First, in the field of material science. In fact, RD processes originate macroscopic patterns
that have ordered features spanning many spatial scales up to the micro- or even the nano-level
(Grzybowski 2009). What is remarkable is the possibility of generating structures having dimensions
significantly smaller than those characterizing the initial distribution of reagents. Moreover, the tiny
features of the structures may be finely shaped by changing macroscopic initial conditions, such as
the concentrations of the reagents, their gradients, and the identity of the hosting medium.
Another appealing technological application of the RD phenomena is in the development of sys-
tems of artificial intelligence (Adamatzky and De Lacy Costello 2012; Gentili 2013). The artificial
RD systems, described in this chapter, are smart because they imitate the computational strategies
of their biological counterparts tightly. The next generation of computing machines might be just
Reaction-Diffusion Computers. An RD computer would be a spatially extended chemical system,
which would process information using interacting growing patterns, excitable, sub-excitable and
diffusive waves. In RD processors, the elementary computing elements are micro-volumes, working
in parallel and communicating with the closest neighbors through diffusion. Information is encoded
as the concentration of chemicals, and the computation is performed via the interaction of wave
fronts. The RD processors could be implemented in either geometrically constrained architectures
or free space or encapsulated in elastic membrane or hybrid designs. On the other hand, stationary
patterns could work as memory elements of the futuristic Chemical Computers wherein processors
and memory will not be physically separated like in our current electronic computers that are based
on the Von Neumann architecture.
Presumably, in the next future, we will witness a great contribution of the RD processes to the
development of materials science and the unconventional information technology.
9.11 KEY QUESTIONS
• What happens if an autocatalytic reaction is carried out in an unstirred container?
• Describe the two-variables model for the formation of Turing patterns.
• Which are the conditions required to observe Turing patterns?
• What was the fortuitous event that allowed discovering Turing patterns in a chemical labo-
ratory? Which was the trick to achieve three-dimensional Turing patterns?
• Repeat the procedure for nondimensionalizing a system of differential equations.
• Which are the principal processes responsible for the embryological development?
The Emergence of Order in Space
291
• Make examples of the application of the Turing’s Reaction-Diffusion model.
• Which are the roles of cytoskeleton within a cell?
• Which are the transport processes that can compete with diffusion in pattern formation
within a cell?
• Describe how a molecular motor works.
• What distinguishes a chemical wave from a physical wave?
• Present the Fisher-Kolmogorov equation.
• Present some examples of chemical waves in biology.
• What is a Liesegang pattern and how does it form?
• Which are possible applications of Reaction-Diffusion patterns?
9.12 KEY WORDS
Balance equation; Turing’s Reaction-Diffusion model; Positional Information; Molecular motor;
Advection; Cell polarity; Trigger waves; Phase waves; Periodic precipitations.
9.13 HINTS FOR FURTHER READING
• For deepening the subject of self-organization and evolution in biological systems, I sug-
gest reading the books The Origins of Order by Kauffmann (1993), On Growth and Form
by Thompson (1917), which allows contemplating the forms of life, Life’s Ratchet. How
Molecular Machines Extract Order from Chaos by Hoffmann (2012), and the review by
Lander (2011).
• More information about chemical waves can be found in the book Chemical Waves and
Pattern edited by Kapral and Showalter (1995), and in the book Chemical Oscillations,
Waves and Turbulence by Kuramoto (2003). In biology, there are many examples of chem-
ical waves as shown in the review by Deneke and Di Talia (2018), and in that by Volpert
and Petrovskii (2009). The chemical processes responsible for intracellular calcium waves
are described in the review by Berridge et al. (2000). Strogatz (2004) dedicates two chap-
ters of his intriguing book, titled Sync, to the theme of brain waves.
9.14 EXERCISES
9.1. If a reaction with an autocatalytic step is left unstirred, exciting phenomena can be
observed. So, go to your wet laboratory, wear a white coat, gloves, and safety glasses.
Prepare the following solutions using deionized water as the solvent.
• Solution A: H SO (sulfuric acid) 0.76 M.
2
4
• Solution B: KBrO (potassium bromate) 0.35 M.
3
• Solution C: CH (COOH) (malonic acid) 1.2 M.
2
2
• Solution D: Ce(SO ) (cerium(IV) sulphate) 0.005 M in an aqueous solution containing
4 2
H SO 0.76 M.
2
4
• Solution E: Ferroin (tris(1,10-phenanthroline)iron(III) sulphate) 0.025 M.
In a test tube having a diameter of about 1.6 cm, introduce the following amounts of the stock
solutions: 2 mL of A, 2 mL of B, 2 mL of C, and 2 mL of D. Stir the final solution with a glass
rod, efficiently and quickly. Then, drop 200 μL of E, vigorously. Do not shake the solution,
but fix the test tube to a laboratory stand through a clamp. Then, wait and stare at the solution.
a. What do you see?
b. Try to determine the wavelength of the periodic spatial structure you observe.
c. Apply equation [9.10] to determine the rate constant of the autocatalytic
step responsible for the pattern you see. Remember that, according to the
292
Untangling Complex Systems
FKN model, the autocatalytic step in the mechanism of the BZ reaction is
HBrO
−
+
3
+
4
+
2 + BrO3 + H
3
+ Ce
2
→ Ce
2
+ 2HBrO2 + H2O. In equation [9.10], b repre-
sents [BrO−3], D is the diffusion coefficient of HBrO
M
2, which is equal to 2 × 10−5 cm2s−1.
d. Repeat the experiment by changing only the amount of KBrO . Instead of 2 mL, add
3
200 μL of solution B plus 1.8 mL of water in the test tube. What happens to the wave-
length of the periodic spatial structure?
9.2. In Figure 9.29, there are three spatial distributions of the concentration u in a mono-dimensional “box” of dimension 1.
The first case (trace labeled as 1) refers to a situation of uniform distribution of u along
x: in fact, u = 1 everywhere along x in the interval [0, 1]. Such spatial distribution of u comes from a perfect mixing. On the other hand, trace 2 refers to the case of a linear
increase of the concentration of u from 0 at x = 0, up to u = 2 at x = 1. Finally, trace 3 refers to a step function: u = 0 in the first half of the x box, whereas it becomes equal to 2 in the
second half. In any case, the average value of u within the range [0, 1], is 1.
Calculate the average rates of hypothetical reactions that are governed by the following
kinetic laws: (a) rate ∝ u; (b) rate ∝ u 2; (c) rate ∝ u 3; (d) rate e( u
∝
− )
1 ; (e) rate
( u
∝
− )
10 1 .
Makes final considerations about the importance of mixing, especially when the kinetic
law is non-linear.
9.3. Find the minimum of det( J )′, defined in equation [9.20], with respect to K.
9.4. Try to find an analogy of the “local self-activation and lateral inhibition” that gives rise to
Turing structures in nature and/or in our daily life.
9.5. The conditions required to have diffusion-driven instability are:
I.
tr ( J ) = R
( )
x X
R Y
+ ( )
s
y
< 0
s
II.
det ( J ) = ( R
( )
)
x X R
Y
R Y R X
( )
> 0
s
y
− ( )
s
x
( )
s
y
s
III.
( D ( )
)
X
Ry Y
D R X
+
( ) > 0
s
Y
x
s
/
1 2
IV.
( D
) 2
([
]
)
X R
y Y
( )
D R ( X )
> D D
R ( X ) R Y
( )
R Y
( )
+
[
]
R
y( X )
s
Y
x
s
X
Y
x
s y
−[
]
s
x
s
s
Apply these conditions to the following set of Reaction-Diffusion equations for the station-
ary state ( x
)
s , ys = ( ,
0 0):
∂ x =α x(1− r 2
2
1 y ) + y (1− r 2 x ) + ( Dδ )∇ x
∂ t
∂ y
r 1
= β
α
y 1+
xy + x(γ + r
2
)
2 y + δ ∇ y
∂ t
β
2.0
3
1.5
2
u 1.0
1
0.5
0.0
0.0
0.2
0.4
0.6
0.8
1.0
x
FIGURE 9.29 Three trends of the concentration u versus the x-coordinate.
The Emergence of Order in Space
293
Such model has been obtained by expanding the functions R( X) and R( Y) in a Taylor series about the homogeneous stationary state ( x
)
s , ys , neglecting terms of order higher than
cubic (Barrio et al. 1999).
Moreover, determine the critical wavenumber of the Turing structure as a function of
the diffusion-reaction parameters, using equation [9.24]. If the length of the domain is L,
how many waves appear in the pattern?
9.6. Consider the reaction-diffusion equations of the exercise 9.5. The parameters have the fol-
lowing values: D = 0.516; δ = 0.0021; α = 0.899; β = −0.91; r
1 = 3.5; r 2 = 0; γ = −α = −0.899.
What do you expect if the system is in a mono-dimensional space?
Let us assume that the spatial domain is L = [−1, 1] and that we have zero-flux bound-
ary conditions. Note that the Laplacian operator for diffusion must be calculated by a
finite difference method that approximates the derivatives in each direction. Consult any
textbook regarding the numerical solutions of partial differential equations. For instance,
the textbook by Morton and Mayers (2005).
9.7. Consider the following reaction-diffusion equations:
∂ u
2
2
2
u
u
= α u(1− rv 1 ) + v(1− r 2 u) + ( Dδ ) ∂
+ ∂
∂ t
∂ x 2 ∂ 2
y
∂ v
2
2
= β
α r 1
∂ v
v
v 1+
uv + u(γ + r )
2 v + δ
+ ∂
2
2
∂ t
β
x
∂
y
∂
These equations are formally identical to those presented in exercises 9.5 and 9.6, but now
we use different symbols for the variables for a matter of clarity.
Solve the system of partial differential equations on a square domain with a side
included between [−1, +1], over time and in three different situations that differ only in the
values of the r and parameters. Use the Euler’s method (see Appendix A).
1
r 2
The
parameters
r and affect the interaction of the activator and inhibitor. First case:
1
r 2
r
1 = 3.5 and r 2 = 0.
Second
case:
r
1 = 0.02 and r 2 = 0.2.
Third
case:
r
1 = 3.5 and r 2 = 0.2.
The values of the other parameters are always the same in the three cases. They are:
D = 0.516, δ = 0.0021, α = 0.899, β = −0.91, γ = − α.
Note
that
D = ( D
)
U DV is the ratio between the diffusion coefficients of the activa-
tor and inhibitor, respectively. D must be less than 1 since the diffusion of the inhibitor
must be greater than that of the activator. The value of δ acts to scale the diffusion
compared to the chemical reactions. In every case, assume to have zero-flux boundary
conditions.
9.8. On the webpage of Prof. Shigeru Kondo’s research group (http://www.fbs.osaka-u.ac.jp/
labs/skondo/) it is possible to download a Reaction-Diffusion simulator. An easy guide to the simulator can be found on the Supporting Material of the paper published by Kondo
and Miura in 2010. The simulator solves the following system of partial linear differential
equations of the Turing model:
∂ u = a
2
uu − buv + cu − duu + Du∇ u
∂ t
∂ v = a
2
vu + bvv + cv − dvv + Dv∇ v
∂ t
It is possible to play with the values of the parameters and see how they affect the final pat-
tern. Try the combinations of the parameters values proposed by the authors of the simulator.
294
Untangling Complex Systems
In the end, try your combination of the parameters values, after checking that the four
conditions required to generate a Turing structure hold.
9.9. Determine the steady-state solution for the CIMA reaction in the homogeneous case
(equation [9.35]). Define the values of the a and b parameters that guarantee a stable
steady state.
9.10. In the paper by Lengyel and Epstein (1991) there are the values of the kinetic con-
stants for the simplified kinetic model of the CDIMA reaction [9.31 and 9.32]. They are:
k
3
−
1
−
5
−
3
1
−
1
−
3
−
1
−
a
1 = 7 5
. ×10 s , k b 1 = 5×10 M, k 2 = 6×10 M s , k 3 b = 2 6
. 5×10 s , α = × −
1 10 14 2
M .
If the initial concentrations of the reactants are [CH
3
]
−
3
[ ]
−
2 (COOH)2
= 10 M, I2 =10 M,
0
0
[ClO ]
4
2
= 6× −
10 M, and they are maintained constant, determine how large should be the
0
ratio d = ( D
)
Y DX to observe the formation of Turing patterns.
9.11. In the Schnackenberg model, the reaction steps are:
A k 1
→
X
X
k 2
→
P
2 X + Y
k 3
→
3 X
k
B
4
