Toom's NEC rule

In the midst of a casual conversation sometime in 1984, Charlie Bennett mentioned to me that Andrei Toom [1] had constructed a probabilistic cellular automaton (CA) [2] with a rather amazing steady-state phase diagram. Toom had proved rigorously that the phase diagram exhibited “generic bistability,” that is, two genuinely stable phases coexisting over a nonzero fraction of the parameter space. Now, almost twenty years later, this innocuous-sounding result remains one of the more remarkable surprises I have experienced in thirty-odd years of work in statistical mechanics.

Before exploring the mechanism and implications of this phenomenon, let us summarize more precisely what Toom's so-called “NEC model” or “NEC rule” is and does. The CA is defined on a two-dimensional square lattice, each site i of which is occupied by a binary variable, S_i, which assumes the values ±1. We often refer to these variables as “spins,” which point either “up” (+1), or “down” (-1), employing the terminology of the Ising model [3], whose analogy to the NEC cellular automaton is most revealing, as we shall see. The dynamical rule under which the {S_i } evolve in time can be described as follows:

1. Any given spin S_i, denoted C (for center), evolves in discrete time according to a triangular neighborhood of spins consisting of C itself, its neighbor to the north (denoted N), and its neighbor to the east (denoted E); hence the name “NEC.” All spins are updated simultaneously.
2. The rule consists of two parts, the first (a) deterministic, and the second (b) probabilistic: (a) The first part is a simple majority rule: If two or more of the spins N, E, and C point up (down) at time t, then the value, denoted F (for future) of spin C at time t + 1 is up (down). (b) If the result of applying (a) is that F points up, then flip F down with probability q; if the result of (a) is that F points down, then flip F up with probability p.

The NEC rule on a finite lattice is fully probabilistic in that there is a nonzero probability of any state, {S_i }, of the system evolving into any other state, {S_i^' }, in a single time step. The rule obviously depends on the two parameters, p and q, whose rough effect is clear: When p > q, up spins are favored over down, since the rule flips more down spins to up than the reverse; when p < q, down spins are favored over up; and when p = q up and down spins are treated equivalently, implying a global up–down symmetry. It follows that the “bias,” i.e., the quantity h ≡ (p - q)/(p + q), is a rough analogue of the magnetic field in the Ising model. To complete the Ising analogy, one further defines the “noise” in the system, T ≡ p + q, and the expectation value M = <S_i> of any spin, which are respective analogues of the temperature and average magnetization in the Ising model.

Toom proved that in the long-time limit, the NEC model achieves a steady state characterized by the schematic phase diagram of Figure 1, shown in the T–h plane: For sufficiently large values of h and/or T, the model has a unique stable phase with a steady-state magnetization M with the same sign as the bias h. This is exactly what one would expect. Unanticipated, however, is the (roughly triangular) shaded region of the phase diagram near T = h = 0, where the model exhibits two coexisting stable phases, one with M > 0 and the other with M < 0. Along the T axis, the two stable phases have equal and opposite magnetization values, ±M, with M decreasing continuously to zero at the second-order critical point located at h = 0 and T = T_c, the “critical temperature,” or “critical noise” value. At nonzero h, two curved first-order phase boundaries, mirror-symmetric about the h = 0 axis, separate the region of coexisting stable phases from the single-phase region.

Figure 1

Comparison of the NEC and Ising models

To understand how extraordinary this deceptively simple phase diagram really is, let us compare it with that of its much more familiar Ising-model analogue. First recall that the simplest near-neighbor, ferromagnetic equilibrium Ising model on a square lattice in two dimensions is defined by the classical Hamiltonian

Here ∑_<ij> denotes the sum over all nearest-neighbor pairs of sites i and j;¹ the exchange constant J is positive, so that the system lowers its energy by aligning near-neighbor spins; and h, the magnetic field, tends to make all spins point up if h > 0 and point down if h < 0. The stationary equilibrium phase of the Ising model at temperature T is of course described by the Boltzmann distribution, the probability of finding the system in state {S_i } being proportional to
e^{-H({S_i })/k_BT}, where k_B is Boltzmann's constant. The famous phase diagram summarizing the stationary properties of the Ising model is shown in Figure 2, again in the T–h plane. (We use the symbols h and T to represent both the field and temperature in the Ising model and the analogous quantities in the NEC rule.) This phase diagram is strikingly different from that of the NEC rule. In particular, the two-phase coexistence region of Figure 2 occupies only a set of measure zero of the parameter space, being confined to the portion of the T-axis below the critical temperature, T_c, of the Ising model. On this special line, the Ising system indeed exhibits two stable coexisting ferromagnetic phases with equal and opposite magnetizations, consistent with both the up—down symmetry of H({S_i }) at h = 0 and the NEC phase diagram in Figure 1. For any nonzero h, however, the Ising system has a unique stable (paramagnetic) phase, its magnetization M having the same sign as h. For h, = 0 and T > T_c, the Ising model has, like its NEC counterpart, a unique stable paramagnetic phase with M = 0.

Figure 2

The following simple argument demonstrates that the Ising phase diagam is prototypical of equilibrium systems, where two-phase coexistence occurs only along a surface of dimension one less than that of the parameter space: Any equilibrium system that depends on N parameters, ₁, ₂, … , _N, has a free energy F = F(₁, ₂, … , _N). At any point of two-phase coexistence, the free energies, F₁ and F₂, of the two phases must be equal; the resulting equation, F₁(₁, ₂, … , _N) = F₂(₁, ₂, … , _N) can be solved for one of the _i in terms of the other N – 1, producing an (N – 1)-dimensional coexistence surface in the N-dimensional parameter space. This is why coexistence surfaces are called “surfaces” or “curves,” rather than “volumes” or “areas.”

We therefore conclude that the NEC model must be a nonequilibrium or irreversible system. Such systems are usually thought of as open systems, driven externally in some way that prevents them from reaching equilibrium. However, many nonequilibrium toy models such as Toom's NEC rule are defined on closed systems, the effects of the openness and external driving being incorporated in the dynamical rule. To give a fanciful example dear to Charlie's heart, one could imagine one of Maxwell's demons living on a two-dimensional lattice of spins. Armed with a tiny magnetometer, the demon moves around the lattice, flipping spins in such a way as to implement Toom's rule. In describing the resulting system, however, we can safely ignore the demon, i.e., the coupling of the system to the external world, and simply imagine that the system follows the NEC dynamics, without inquiring about the origin of the rules.

In any event, it is easy to infer directly that the NEC model is not an equilibrium system, since it fails to satisfy detailed balance for any energy function. Recall that the dynamics of a noisy, discrete-time system is said to obey detailed balance for the energy function E({S_i }) if, given that the system is in the state {S_i } at time t, the probability, Q({S_i^' },|{S_i }), of its undergoing a transition to the state {S_i^' }, at time t + 1 satisfies the equation

(2)

It is readily apparent that systems obeying this condition have, in the long-time limit t → ∞, a stationary state described by the Boltzmann distribution,² wherein the probability of occurrence of any state {S_i } is proportional to e^{-E({S_i })/k_BT}. To show that the NEC rule cannot obey detailed balance for any conceivable energy function E({S_i }), assume the contrary, i.e., that there exists an E({S_i }) for which the NEC transition probability Q satisfies Equation (2). The fact that the evolution of any site C depends on its neighbor to the north, site N, implies that E({S_i }) contains the term S_CS_N. (To see this, note that an arbitrary nonpathological function of the binary variables S_C and S_N can be written as a linear combination of the four terms: 1, S_N, S_C, and S_CS_N. Of these four, only the S_CS_N term makes the time evolution of site C depend on the value of site N.) But the presence of the term S_CS_N in E({S_i }) implies in turn that the evolution of site N depends on C, its neighbor to the south, which is false, whereupon the reductio ad absurdum argument is complete.

Let us now try to understand on an intuitive level how the nonequilibrium NEC rule manages to achieve generic multistability. A relatively simple explanation can be formulated in terms of droplet dynamics [4]. It is helpful to review how such arguments rationalize the phase diagram of the Ising model before tackling the NEC problem.

Rationalizing generic multistability: Droplet arguments

Let us begin at h = 0 and T < T_c in the Ising model, i.e., on the line of two-phase coexistence. To understand how both of the oppositely magnetized phases can be stable, imagine preparing the system in an initial state in which a straight domain wall separates a domain of up spins from one of down spins. Now at h = 0 this domain wall cannot translate with nonzero velocity, since that would violate the up—down symmetry of Hamiltonian (1), i.e., the equivalence of the up and down spins. Thus, in the thermodynamic limit, the wall can develop local wiggles, but its average position must remain fixed. If this position is denoted by R, the situation is described by the trivial equation R/t = 0.

Next imagine starting the system off in an initial state consisting of a single circular droplet of radius R of up (down) spins immersed in a sea of down (up) spins. Such a droplet can lower its surface energy by shrinking in size, its dynamics being governed by the phenomenological equation [5]

R/t –/R ,

(3)

where is the surface tension of the the droplet. In the limit R → ∞, this equation reduces to R/t = 0, as it must, since the droplet becomes the flat interface considered just above. We conclude that droplets of minority spins tend to shrink and disappear when h = 0.

Now imagine starting from an initial condition wherein all of the spins in the system are down (up). Because of thermal fluctuations, minority up (down) droplets will form spontaneously throughout the system. Since, as we have seen, however, such droplets tend to shrink rather than grow, their total weight will remain limited, provided that the temperature of the system is sufficiently low. Thus, for low enough T—meaning, in practice, T < T_c— a state with an initial positive (negative) magnetization will, in the thermodynamic limit, have M > 0 (M < 0) in perpetuity. This explains heuristically how the two oppositely magnetized phases of the Ising model with h = 0 and T < T_c can be genuinely stable.

Finally, let us next move off the coexistence curve of the Ising model by considering h > 0. Here up spins are favored energetically over down, so a flat interface separating domains of up and down spins will translate with a uniform velocity proportional to h itself, for small h. This situation is represented by the phenomenological equation R/t ~ h. A droplet with radius R of up spins immersed in a sea of down spins is then described by the equation

R/t –/R + h ,

(4)

which reduces to the required R/t = h as R → ∞, and also incorporates the effect of surface tension. This equation provides a heuristic description of the important phenomenon of metastability: Droplets with radii R less than the “critical droplet size” R_c ≡ /h will shrink with time; however, any droplet with a radius R > R_c will grow forever according to Equation (4). Even if one prepares the system in an initial state in which all spins point down, and even if h and T are very small, thermal fluctuations will eventually nucleate a density of up-spin droplets whose radii exceed R_c. These droplets will grow and merge, ultimately producing a stationary state with positive magnetization. This nucleation and growth process explains why for h ≠ 0 the Ising model has a unique stable phase with nonzero magnetization. Note, however, that if /h, and hence R_c, is sufficiently large, it can take considerable time for the nucleation of critical-size droplets. (The time required for the nucleation of such droplets grows exponentially with R_c [6].) Thus, a metastable phase (e.g., a phase with negative magnetization for h > 0), can persist for extremely long times before being destabilized by droplet growth.

To apply this experience with droplet dynamics to Toom's CA, let us return to the NEC rule in the deterministic limit p = q = 0, analogous to h = T = 0 in the Ising model, and imagine preparing the system in an initial condition with a straight interface separating domains of up and down spins (Figure 3). The NEC majority updating rule immediately implies that such interfaces are totally immobile, provided they are oriented with the lattice axes, i.e., in the vertical or horizontal directions [Figures 3(a) and 3(b)]. This mimics the situation in the Ising model. However, the interface oriented at 45° along the NW—SE axis, as in Figure 3(c), translates downward, maintaining its orientation, by one lattice spacing with every time step, since any site C just below the interface has both its N and E neighbors on the opposite side of the boundary. This phenomenon, which holds regardless of whether the up or down domain lies below the domain wall, is a consequence of the peculiar spatial asymmetry of the Toom rule.

Figure 3

The striking ramifications of this behavior become obvious when one considers the evolution of a droplet of either sign, immersed in a sea of the opposite sign. First consider an isosceles right triangular droplet with one vertical side of length R, one horizontal side with length R, and one side oriented from NW to SE (Figure 4): The horizontal and vertical sides remain fixed in position, while the diagonal one translates downward, shrinking the droplet to zero in a time proportional to R.

Figure 4

Since a circular droplet of radius R cannot shrink more slowly than the isosceles right triangular droplet circumscribed around it (Figure 4), it follows that such a circular droplet (of either sign) must also shrink to zero in a time proportional to R. Thus, the dynamics of such a droplet is described by the phenomenological equation R/t ~ –1. Now, the effect of allowing the parameters p and q to be nonzero can only be to add terms of O(p, q) to the right side of this equation, provided p and q are small compared to unity. In this case, the droplet equation therefore takes the form

R/t –1 + O(p, q) ,

(5)

which implies that even for nonzero p and q, droplets of either sign will shrink in a time proportional to their radius, so long as p and q are sufficiently small. In this regime, therefore (which of course corresponds precisely to the two-phase coexistence region of the NEC model), the sign of the initial magnetization of the system will be preserved in the long-time limit, since droplets of either sign are efficiently suppressed. This accounts in a simple intuitive way for the remarkable generic multistability of Toom's rule: Even for p > q, say (meaning up droplets being favored over down), up droplets that nucleate spontaneously in an initial sea of down spins cannot grow to destabilize the system, provided p is sufficiently small. Hence, there is a genuinely stable state with M < 0 in this regime, in addition to the expected one with M > 0. Thus we understand how the irreversibility embodied in the spatial asymmetry of the Toom rule eliminates errors of both signs efficiently enough to generate generic two-phase coexistence.^3,4

It is important to try to identify the minimum conditions required to produce such generic multistability. It has been argued [7] that neither the simultaneous updating nor the discreteness of time, space, and variables inherent in CA is necessary: The only indispensable ingredient seems to be irreversibility in the specific form of spatially asymmetric rules or partial differential equations.

The generic stabilization of complex structures

Having understood heuristically the unconventional phase diagram (Figure 1) of the Toom rule, let us briefly explore its implications for the stabilization of complex structures, or “error correction.” The basic question is inspired by the evolution of living creatures, specifically the fact that the complexity of such creatures seems to be increasing with time. There is little doubt that modern man, for example, is (despite the odd counterexample one is occasionally distressed to encounter), a more complex creature than anything that existed on this planet ten million years ago. This increase in complexity obviously requires nonequilibrium processes working against entropy, and raises an important question of principle: Are complex structures simply metastable, or are there local nonequilibrium rules, operating at nonzero temperature and therefore noisy, that are capable of producing genuinely stable structures of arbitrary complexity under generic conditions, i.e., with no tuning of parameters? At least for strings of binary variables, this last question was answered compellingly in the affirmative through a specific one-dimensional model due to Gacs [8] and through a much simpler three-dimensional model due to Gacs and Reif [9], which is based on the Toom rule. The rest of this paper is devoted to a discussion of this latter construction and its significance.⁵

Imagine for specificity some complex pattern, or at least the information required to construct it, encoded in a long binary sequence. (A binary encoding, in 1's and –1's, of all the information contained in the human genome is a particularly stirring example.) Suppose that the first few bits in that pattern of interest, which we denote S*, or S_i* for i = 1, 2, 3, … , are

Assume, moreover, that the desired sequence is a fixed point of a one-dimensional deterministic, local, dynamical rule—a one-dimensional CA, say—that has been used to create the sequence from some appropriate initial condition, through updates of the system. Since one-dimensional deterministic CA rules capable of universal computation (i.e., capable of carrying out an arbitrary digital computation, given the right initial condition), have been proven [10] to exist, this assumption is not unreasonable.⁶ (In the example above, in which the sequence is supposed to represent a blueprint for the human genome, we can loosely think of the automaton that produces it as mimicking, in a grossly oversimplified way, some of the effects of evolution.) For compactness of terminology, let us refer to the one-dimensional automaton used to produce the sequence of interest as the “creation rule.”

As long as this automaton is deterministic, then once it reaches, after iterations, the desired fixed point above, nothing changes further, and the sequence of interest persists forever. To address the question of the stability of complex patterns under noisy conditions, however, we must give the creation rule some probabilistic character. This is easily accomplished by taking the result of each step of the deterministic creation rule, and then flipping all down spins to up and up spins to down with probabilities r₊ and r_-, respectively, just as in the Toom rule. Of course, this will soon destroy the desired pattern (6) above, since, even for very small r₊ and r_-, nothing prevents the errors from proliferating with time, particularly because a faulty result produced by the noise at time t is likely to cause even the deterministic piece of the rule to generate further faulty results at later times. The probabilistic creation rule is therefore a very poor model for generating stable structures.

Obviously, what is required to rescue the model is a mechanism for employing redundancy to correct errors. Following Gacs and Reif [9], therefore, let us imagine a system consisting not of a single one-dimensional CA but a large number of identical one-dimensional automata, all running the creation rule with the same properly chosen initial condition. Imagine further that all of these one-dimensional CA chains are oriented in the same (say z) direction, and that they are positioned on the sites of a square lattice in the x and y directions, forming a three-dimensional cubic lattice overall (Figure 5). Furthermore, align the one-dimensional automata in the z direction so that the corresponding sites of each one lie in the same XY plane, making the values of all the sites in any given XY plane identical initially. If we now run the deterministic creation rule on each of the chains, the values of all of the sites in any XY plane will remain identical, and each automaton will reach the fixed point S* after time steps.

Figure 5

Now imagine adding to the creation rule the probabilistic component discussed above, with small error rates r_±. Errors will then proliferate randomly along each of the chains, making the values of the sites in any given XY plane nonuniform and rendering each of the automata incapable of achieving the desired fixed point. To suppress this proliferation, Gacs and Reif proposed coupling the equivalent sites of the one-dimensional automata in each XY plane with the Toom rule, producing a fully three-dimensional automaton in which each probabilistic creation-rule update of the one-dimensional chains is followed by s time steps of the two-dimensional Toom rule in every XY plane of the system. The values of p and q for the Toom rule are chosen small enough so that the rule is in its region of two-phase coexistence. The value of s is chosen large enough to bring the average magnetization of each plane to within some desired fraction of the steady-state value of the NEC rule (see footnote 3) with the sign dictated by the result of the previous update of the creation rule.

This algorithm guarantees that after any number, t, of updates of the creation rule [i.e., after t(1 + s) total updates of the system], the average magnetization, M_i(t), of the ith plane has the correct sign, viz., sgn[M_i(t)] = S_i^D(t), where S_i^D(t) is the value of the ith site of the deterministic creation rule after t updates. Thus, for all t ≥ and for every site i, sgn[M_i(t)] = S_i* the desired fixed point value of the deterministic creation rule. In this average sense, therefore—the only sense in which a noisy system such as our three-dimensional CA can hope to reproduce the result of a deterministic calculation—the Gacs—Reif construction succeeds in generating precisely the outcome of the deterministic creation rule.

Gacs and Reif have proven that their three-dimensional CA accomplishes this task for the creation rule (or any other deterministic one-dimensional CA). The physical mechanism underlying the proof is quite easy to understand, however: While each update under the probabilistic creation rule introduces a small fraction of O(r_±) of new sites with the wrong sign (errors) within any plane i, running the NEC rule for sufficient time s eliminates enough of these errors to bring the average magnetization M_i of that plane arbitrarily close to the steady-state value of the NEC rule with the correct sign. Thus errors cannot proliferate in time. After each step of the creation rule and the associated s steps of the NEC rule, the fraction of errors in each plane is ~O(p, q). Hence, for small p and q, the sign of M_i(t) for any i will be identical to the desired result, S_i^D(t), of the deterministic creation rule.

Since this error-correction mechanism works for any sufficiently small values of r_±, p, and q, it is truly generic, requiring no tuning of parameters. Thus, it establishes that arbitrary one-dimensional binary sequences can indeed be stable in nonequilibrium systems under generic conditions and in the presence of noise. Note that the equilibrium Ising model along its coexistence curve, h = 0, also shrinks minority islands of either sign, and so can be used in place of the NEC rule to perform in-plane error correction. However, this requires the tuning of the magnetic field h to zero, and so cannot be done generically. It is a distinct advantage not to need to invoke any mysterious tuning of parameters to stabilize complex structures.

One could of course carry out generic error correction with the Ising model for T < T_c and h very small: h  J say. Since this relies on the metastable phase (the phase with negative magnetization for h > 0, say) being very long-lived in the small-h limit, it will fail when the unique stable phase of the Ising model in nonzero field ultimately manifests itself. However, because the time required for the nucleation and growth of this phase is exponential in J/h, one can readily make it exceed the age of the universe by simply choosing a small h. For all practical purposes, therefore, this metastable generic error correction will work as well⁷ as the real thing. It is nonetheless comforting to know that simple local nonequilibrium rules that cleverly exploit spatial asymmetry can truly stabilize complex structures generically.

Acknowledgments

I thank Charlie Bennett for teaching me many wonderful things about science, for giving and taking hundreds of jokes, and for generous and good-humored collegiality over the past 25 years.

References and notes

Footnotes

¹One obvious difference between the two-dimensional Ising and NEC models is that in the former, the evolution in time of a particular spin depends on all four neighboring spins, rather than just two; the spatial asymmetry in the NEC rule is central to its behavior.
²The existence of a stationary Boltzmann distribution generally implies that the associated stationary phase is the unique stable phase of the system, but there are exceptions which we ignore here.
³Starting from an arbitrary initial condition, the magnetization of the NEC model in its two-phase region relaxes exponentially, with a characteristic relaxation time, to the steady-state value for the appropriate stable phase of the model. The same is true of the Ising model; see [6].
⁴As in the Ising model, two-phase coexistence of the NEC model occurs only in the thermodynamic limit: For any finite system, droplets encompassing the entire system form with nonzero probability. Since the time required to nucleate such droplets increases exponentially with sample size, however, reasonably large systems will exhibit two-phase coexistence over immeasurably long time scales.
⁵We use the term “complex” in a loose, intuitive way here, begging the difficult and contentious issue of providing a careful definition. In this sense, our attitude toward complexity is much like the one expressed toward pornography by former U.S. Supreme Court Justice Potter Stewart, who famously opined: “I know it when I see it.”
⁶It is important to keep in mind, however, that there is, in general, no practical procedure for determining the right initial condition.
⁷Recall from Equations (3) and (5), however, that in the NEC model, droplets of radius R disappear in a time of O(R), more efficiently than in the Ising model, where the time required is of O(R²).