\documentstyle[graphicx,aps,epsf,multicol]{revtex}
\begin{document}
\title{Segregation in a one-dimensional model of interacting species}
\author{L.~Frachebourg$^\ast$, P.~L.~Krapivsky$^{\ast}\dag$, and 
E.~Ben-Naim$\ddag$} 
\address{$^\ast$Center for Polymer Studies and Department of Physics,
Boston University, Boston, MA 02215}
\address{$\dag$Courant Institute of Mathematical Sciences, 
New York University, New York, NY 10012-1185}
\address{$\ddag$The James Franck Institute, The University of Chicago, 
Chicago, IL 60637}
\maketitle
\begin{abstract}
We investigate segregation and spatial organization in a one-dimensional
system of $N$ competing species forming a cyclic food chain. For
$N<5$, the system organizes into single-species domains, with an
algebraically growing average size.  For $N=3$ and $N=4$, the domains
are correlated and they organize into ``superdomains'' which are
characterized by an additional length scale. We present scaling
arguments as well as numerical simulations for the leading asymptotic
behavior of the density of interfaces separating neighboring domains.
We also discuss statistical properties of the system such as the
mutation distribution and present an exact solution for the case
$N=3$.

{PACS numbers: 02.50.Ga, 05.70.Ln, 05.40.+j}
\end{abstract}
\begin{multicols}{2}

Coarsening underlies numerous natural processes including phase
separation, grain growth, soap bubbles, and species segregation.  It
is generally believed that coarsening systems exhibit dynamical
scaling\cite{bray}, {\it i.e.}, the typical domain size grows
algebraically with time, $\ell(t)\sim t^\alpha$.  The exponent
$\alpha$ is usually independent of many details of the system such as 
the spatial dimension.  However, much less is known on coarsening in
systems with more than two equilibrium phases.  In this Letter, we
show that species segregation can exhibit two-length scaling rather
than ordinary single-length scaling.

The Lotka-Volterra model of interacting populations ``living'' on a
one-dimensional lattice is a simple system which exhibits species
segregation. The case where $N$ species form a food chain is
especially well suited for studying segregation.  We assume that every
species plays the role of prey and predator simultaneously.  The food
chain is arranged in a cyclic manner. For example, when $N=3$, $A$
eats $B$, $B$ eats $C$, and $C$ eats $A$.  ``Eating'' events involve
nearest neighbors and lead to duplication of the winner and
elimination of the loser, corresponding to the following reaction
scheme
\begin{equation} 
A+B\to 2A, \quad B+C\to 2B, \quad C+A\to 2C.
\label{Eq-1} 
\end{equation}
Here and throughout this study we restrict ourselves to random and
symmetric initial conditions, where the average initial species
densities are all equal $1/N$. Despite the nonconserving nature of the
process, the average densities remain constant in the thermodynamic
limit.

For a large number of species, most pairs of species do not interact and
the system quickly reaches a frozen state.  Previous studies
\cite{bramson,fisch} have mainly concentrated on establishing the upper
bound for $N$ above which the system does not coarsen.  It has been
proved rigorously that the marginal chain length is $N_c=5$
\cite{bramson,fisch}. For $N\geq N_c$ each site quickly reaches a final
frozen state, while for $N<N_c$, the state of each site changes an
infinite number of
times. However, theoretical understanding of the kinetic behavior and
the coarsening properties of the system is still incomplete
\cite{fisch,tainaka}. In this study, we illuminate the rich kinetic
behavior of the system by analyzing the density of interfaces separating
different single-species domains.


For $N=2$, this system is equivalent to the voter model\cite{lig,elp}
which can be solved exactly \cite{glauber}. In terms of interfaces,
the $N=2$ model is equivalent to an ensemble of annihilating random
walks.  The system separates into single species domains.  The average
domain size $\ell$ exhibits a diffusive growth law
$\langle\ell(t)\rangle\sim \sqrt{t}$.

Consider now the $N=3$ case. There are two types of interfaces: right
moving ($AB$, $BC$, and $CA$) and left moving ($BA$, $CB$, and $AC$),
denoted by $R$ and $L$, respectively. The interface dynamics and
consequently, the coarsening kinetics are sensitive to the microscopic
realization of the reaction process. For parallel dynamics (bonds
updated simultaneously) opposite moving interfaces annihilate, $R+L\to
\emptyset$, while for sequential dynamics (bonds updated one at a
time) interfaces moving in the same direction react as well, $R+L\to
\emptyset$, $R+R\to L$, and $L+L\to R$.  Hence, for the 3-species
model with parallel dynamics the interface reaction process is
equivalent to the well-known two-velocities ballistic annihilation
process\cite{Elskens} and the interface density, $M(t)$, decays as
$t^{-1/2}$. This behavior can be understood by arguing that in a
linear region of size ${\cal L}$, the imbalance between the number of
left and right moving interfaces is of order $\Delta\sim\sqrt{c_0{\cal
L}}$. After a time $t={\cal L}/v_0$ only this residual fluctuation
remains and as a result the concentration decay $M(t)\sim \Delta/{\cal
L}\sim (c_0/v_0t)^{1/2}$ follows.

The above heuristic picture suggests a special domain pattern.  The
system organizes into ballistically growing superdomains.  Each
superdomain contains interfaces moving in the same direction, while
neighboring superdomains are separated by opposite moving interfaces.
Domains inside each superdomain are arranged cyclically ($ABCABC$ or
$CBACBA$).  In addition to the average size of superdomains, there is an
additional length scale corresponding to the average distance between two
adjacent similar velocity interfaces.  We define these two 
length scales using an illustrative configuration:

\begin{equation} 
B\overbrace{AABBB\underbrace{CCCC}_{\ell}AAABBCCC}^{\cal L}B.
\label{Eq-2}
\end{equation}
The corresponding coarsening exponents, $\alpha$ and $\beta$, are
defined via $\langle\ell(t)\rangle\sim t^{\alpha}$ and $\langle{\cal
L}(t)\rangle\sim t^{\beta}$, respectively.  For $N=3$ with parallel
dynamics we thus have $\alpha=1/2$ and $\beta=1$.

In the complementary sequential dynamics case, interfaces perform a
biased random walk and thus, the ballistic motion is now supplemented
by diffusion.  The system again organizes into domains of right and
left moving interfaces.  Inside a domain, interfaces moving in the
same direction can now annihilate via a diffusive mechanism, unlike
the case of parallel dynamics (more precisely, collision of say two
right moving interfaces gives birth to a left moving interface which
is soon annihilated with the nearest right moving interface).  Similar
single-species annihilation with convective-diffusive transport has
been investigated in Ref.\cite{brk} where the concentration decay
$M(t)\sim t^{-3/4}$ has been established.  This prediction is
consistent with numerical simulations. The simulations also indicate
that the system slowly approaches the asymptotic behavior
$\langle\ell(t)\rangle\sim t^{3/4}$ \cite{lpe}.

The resulting spatial structure is thus similar to the parallel case,
Eq.~(\ref{Eq-2}). However, while the larger length scale remains
unchanged, $\langle{\cal L}(t)\rangle\sim t$, the smaller length scale
is now a geometric average of a diffusive and a ballistic scale. We
conclude that the coarsening patterns are characterized by two length
scales, and the coarsening kinetics are sensitive to the details of
the dynamics.


In the $N=4$ case, there are static interfaces denoted by $S$ ($AC$,
$BD$, $CA$, and $DB$), in addition to the right and left moving
interfaces, ($AB$, $BC$, $CD$, $DA$) and ($BA$, $CB$, $DC$, $AD$),
respectively.  For sequential dynamics, interfaces react upon collision
according to $R+L\to \emptyset$, $R+S\to L$, $R+R\to S$, $L+L\to S$, and
$S+L\to R$.  Under the assumption that neighboring interfaces are
uncorrelated, the interface densities evolve according to the following
rate equations

\begin{eqnarray}    
\label{Eq-3}   
\dot R=&-&2R^2-2RL-RS+SL, \nonumber\\
\dot L=&-&2L^2-2RL-SL+RS, \\
\dot S=&R&^2+L^2-RS-SL.   \nonumber
\end{eqnarray}
Solving these equations subject to the initial conditions 
$R(0)=L(0)=S(0)=1/4$ gives
\begin{equation}    
M(t)={1\over 4+4t}, \quad
S(t)={1\over \sqrt{4+4t}}-{1\over 4+4t}.
\label{Eq-4}
\end{equation}
In the above, $M(t)=R(t)=L(t)$ is the density of moving interfaces.
According to the rate equation theory, the average distance between
two static interfaces, $\langle\ell(t)\rangle \sim t^{1/2}$, grows
slower than the average distance between two moving interfaces,
$\langle{\cal L}(t)\rangle\sim t$.  A nontrivial spatial organization
occurs in which large ``superdomains'' contain many domains of
alternating noninteracting ($AC$ or $BD$) species.  Similar to the
$N=3$ case, there are two relevant growing length scales as in the
following illustration
\begin{equation}
B\overbrace{AACCCAAA\underbrace{CCCC}_{\ell}
AACCAAACCC}^{{\cal L}}D. 
\label{Eq-5}
\end{equation} 
Numerical simulations agree qualitatively with this picture.
However, the quantitative predictions for the coarsening exponents
fail.  

\begin{figure}
\vspace{-.1in}
\narrowtext
%\epsfxsize=\hsize
%\epsfbox{fig1.eps}
\centerline{\includegraphics[width=8cm]{fig1}}
\vspace{-.1in}
\caption{The concentrations of stationary (filled symbols) and moving
(open symbols) interfaces as a function of time for the 4-species
model with sequential (diamonds) or parallel (circles) dynamics in a
log-log plot.  Lines of slope $-1/3$ (solid) and $-2/3$ (dotted) are
shown as a reference.  An average over 100 systems of size $10^6$ was
taken.
\label{fig1}}
\end{figure}
\vspace{-.2in}


In the following, we use heuristic arguments to obtain the exponent
values.  Numerical simulations indicate that parallel and
sequential dynamics are asymptotically equivalent and thus, we restrict 
ourselves to the former simpler case.  The annihilation reaction
$R+L\to \emptyset$ is supplemented by the exchange reaction $R+S\to L$
and $L+S\to R$.  According to the rate theory as well as the
simulations $M(t)\ll S(t)$, and thus, we assume an alternating spatial
structure of ``empty'' regions (with no more than one moving
interface) and ``stationary'' regions (with many stationary interfaces
inside any such region).  If the interface densities obey scaling,
then the size of the empty and the stationary regions should be
comparable.  The average size of an empty or a stationary region is
therefore of the order of $M^{-1}$.  The number of stationary
interfaces inside a stationary region is of the order of $S/M$.  The
evolution proceeds as follows: a moving interface hits the least
stationary particle and bounces back.  Then this interface hits the
least stationary particle of the neighboring stationary region, and
bounces back again.  This ``zig-zag'' process continues until one of
these stationary regions ``melts'', thereby giving birth to a larger
empty region.  If there is a moving particle inside merging empty
region, the two moving particles quickly annihilate. Otherwise,  
the moving particle continues to eliminate stationary
interfaces.  The typical time $\tau$ for a stationary region to melt
is $\tau=M^{-1}\times S/M=S/M^2$.  This melting time $\tau$ is also
the typical time for annihilation of a moving interface and thus,
\begin{equation}    
\dot M\sim -{M\over \tau}\sim -{M^3\over S}.
\label{Eq-6}
\end{equation} 
Using $M(t)\sim \langle{\cal L}(t)\rangle^{-1}\sim t^{\alpha}$ and
$S(t)\sim \langle{\ell(t)}\rangle^{-1}\sim t^{\beta}$, the exponent
relation $2\beta-\alpha=1$ emerges. A second independent exponent
relation, $\alpha+\beta=1$, will be presented in the discussion of the
mutation distribution below.  Combining these two relations we find
that $\alpha=1/3$ and $\beta=2/3$. These values are in good agreement
with parallel as well as sequential simulations. It is seen that for
$N=4$, the coarsening kinetics are independent of the details of the
dynamics, in contrast with the $N=3$ behavior.

For the 5-species cyclic Lotka-Volterra model, it is known that
the system approaches a frozen state\cite{bramson,fisch}.
Nevertheless, it is useful to consider the interface dynamics for the
$N=5$ case, where there are two types of stationary interfaces, $S_R$
($AC,BD,CE,DA,EB$) and $S_L$ ($AD,BE,CA,DB,EC$), in addition to the
right and left moving interfaces, $R$ ($AB,BC,CD,DE,EA$) and $L$
($BA,CB,DC,AD,AE$).  The reaction process is symbolized by $R+L\to
\emptyset$, $R+S_L\to L$, $R+S_R\to S_L$, $S_R+L\to R$, $S_L+L\to
S_R$, $R+R\to S_R$, and $L+L\to S_L$. It is straightforward to
generalize the rate equations (\ref{Eq-3}) to this case as well, and
we merely quote the results.  According to these equations, the static
interfaces approach a final nonzero value $S(t)\to S_{\infty}$, and
the mobile interfaces decay exponentially, $M(t)\sim
\exp(-S_{\infty}t)$, as $t\to\infty$.  Interestingly, the rate equations
correctly predict the marginal number of species $N_c=5$. 

As in the $N=4$ case the qualitative predications of the rate
equations are correct, but the quantitative predictions fail.  Since
the density of mobile interfaces rapidly decreases while the density
of stationary interfaces remains finite we can ignore collisions
between mobile interfaces.  We should estimate the survival
probability of a mobile interface in a sea of stationary ones.  There
are two reactions in which moving interfaces survive although they
change their type, $R+S_L\to L$ and $L+S_R\to R$.  Thus, a 
moving interface is long lived in the following environment: $\cdots
S_RS_RS_RS_RMS_LS_LS_LS_L\cdots$.  Clearly, in such configurations the
zig-zag reaction process takes place.  The moving interface travels to
the right during a time $t_0=(c_0v_0)^{-1}$, eliminates a stationary
interface and travels to the left a time of order $2t_0$, eliminates
an interface and travels back to the right, {\it etc}. Thus, to
eliminate $N_s$ interfaces, the moving interface should spend a time
of order $t\simeq t_0\sum_{i=1}^{N_s}i=t_0N_s(N_s+1)/2$.  Therefore,
the number of stationary interfaces $N_s(t)$ eliminated by a moving
interface scales with time as $N_s(t)\sim \sqrt{c_0v_0t}$.  Special
configuration of length $N_s$ are encountered with probability
$\propto \exp(-N_s)$, and thus, the density of moving
interfaces exhibits a stretched exponential decay, $M(t)\propto
\exp(-\sqrt{c_0v_0t})$.  Hence, the approach towards
the frozen state is slowed down due to spatial fluctuations.


It was pointed out recently that nontrivial behavior underlies 
low-activity or persistent sites in coarsening systems
\cite{derrida}.  It is useful to consider the mutation distribution,
$P_n(t)$, defined as the fraction of sites that have mutated (changed
their state) exactly $n$ times during the time interval $[0:t]$.  Mutation
kinetics and coarsening kinetics are intimately related.  Let $\langle
n(t)\rangle=\sum_n nP_n(t)\sim t^{\nu}$ be the average number of
mutations.  Since every motion of an interface contributes to an
increase in the number of mutations in one site, the mutation rate
equals the density of moving interfaces, $d\langle
n(t)\rangle/dt=M(t)$. Using $M(t)\sim t^{-\mu}$ one has $\nu=1-\mu$.
In the closely related voter model \cite {lig,elp,glauber}
(corresponding to $N=2$), it was found that the mutation distribution
$P_n(t)$ obeys scaling \cite{elp}
\begin{equation}
P_n(t)={1\over \langle n(t)\rangle}\,
\Phi\left({n\over \langle n(t)\rangle}\right). 
\label{Eq-8}
\end{equation} 
The scaling function has the following limiting behaviors
\begin{equation}
\Phi(z)\sim\cases{
z^{\gamma}&$z\ll 1$;\cr
\exp(-z^{\delta})&$z\gg 1$}
\label{Eq-9}
\end{equation}
The behavior in the small argument limit reflects the decay of
persistent sites. For the voter model $P_0(t)\sim t^{-\theta}$ with
$\theta=3/8$ \cite{derrida}.  If this power-law decay holds generally,
then the exponent relation $\theta=\nu(\gamma+1)$ should be satisfied.
The large $z$ limit describes ultra-active sites. A convenient way to
estimate the fraction of such sites is to consider sites which make of
the order of one mutations per unit time.  At time $t$, the fraction
of these rapidly mutating sites is exponentially suppressed,
$P_t(t)\propto \exp(-t)$.  It is therefore natural to assume the
exponential form $\Phi(z)\sim \exp(-z^{\delta})$ for the tail of the
scaling distribution, thereby implying an additional exponent relation
$\delta=1/(1-\nu)$.

The mutation distribution can be exactly calculated for the $N=3$ case
with parallel dynamics. Since the initial $\pm v_0$ interface
velocities are uncorrelated, the interface ballistic annihilation
process can be mapped into a random walk problem.  As a result, the
interface density is found from first passage properties of a random
walk.  Similarly, the mutation distribution can be shown to be
equivalent to the probability that the minimum of a $t$-step random
walk is exactly $n$ (for details, see\cite{lpe}).  Using the
definition of Eq.~(\ref{Eq-8}), and the average number of mutations
$\langle n(t)\rangle\sim\sqrt{4t/3}$, the exact scaling distribution
can be written
\begin{equation}
\Phi(z)={4\over \sqrt{\pi}}\,e^{-z^2}{\rm Erf}(z),
\label{Eq-10}
\end{equation}
with $z=n/\langle n(t)\rangle$. The limiting behaviors of this scaling
function agree with the predictions of Eq.~(\ref{Eq-9}), and the
scaling exponents $\theta=1$, $\nu=1/2$, $\delta=2$, and $\gamma=1$,
satisfy the predicted scaling relations. 

We turn now to the $N=4$ case where according to the above discussion
zig-zag reactions $R+S\to L$ and $L+S\to R$ dominate over the
annihilation reaction $R+L\to \emptyset$ in the long-time limit.  We
therefore consider a simpler solvable case where a single mobile
interface is placed in a regular sea of static interfaces to evaluate
the scaling function $\Phi(z)$. This
interface moves one site to the right, two to the left, three to the
right {\it etc.}  Similar to the above discussion on the survival
probability in the $N=5$ case, at time $t$ this interface has
eliminated $N_s\sim (t/t_0)^{1/2}$ static interfaces, with
$t_0=(c_0v_0)^{-1}$. The origin is visited $N_s t_0$ times, site $1$
is visited $(N_s-1)t_0$, site $-1$ is visited $(N_s-2)t_0$, {\it etc.}
This implies that the mutation distribution is $P_n(t)=\langle
n\rangle^{-1}\Phi(n/\langle n\rangle)$, with $\langle n\rangle\sim N_s
t_0$ and $\Phi(z)=1$ for $z<1$ and $\Phi(z)=0$ for $z>1$. Therefore,
$\gamma=0$.  This approximation is inappropriate for predicting the
tail of $\Phi(z)$ which is sensitive to annihilation of the moving
interfaces.  However, in the small $z$ limit the annihilation process
should be negligible, and thus $\gamma=0$.  The fraction of unvisited
sites is equivalent asymptotically to the survival probability of a
stationary interface and thus $\theta=\alpha$.  Using the previously
established relations $\nu=1-\mu$, $\theta=\nu(\gamma+1)$, and
$\mu=\beta$, we obtain the second independent exponent relation
$\alpha+\beta=1$ which was used to obtain the asymptotic behavior of
the mobile and the static interfaces.

So far, we investigated the cyclic Lotka-Volterra model with {\it
asymmetric} interactions.  To learn how the coarsening kinetics depend
on the interaction rules it is also useful to consider a {\it
symmetric} interaction rule where both of the reaction channels
$A+B\to 2A$ and $A+B\to 2B$ are allowed.  In this case, asymptotically
exact results can be obtained.  We discuss only the $N=4$ case since
the $N=2$ and $N=3$ cases reduce to the voter model (see\cite{elp}).
Denote the static interfaces ($AC$, $CA$, $BD$, $DB$) by $S$ and
moving interfaces by $M$.  The symmetric ``eating'' rule implies that
moving interfaces perform simple random walks.  Interfaces react
according to $M+M\to S$, and $M+S\to M$ or $M+S\to S$ depending on the
local environment.  Moving interfaces undergo diffusive annihilation
and thus, $M(t)\sim t^{-1/2}$.  The fraction of surviving stationary
interfaces is proportional asymptotically to the fraction of sites
which have not been visited by mobile interfaces up to time $t$,
$S(t)\sim P_0(t)\sim t^{-3/8}$.  We should also take into account
creation of stationary interfaces by annihilation of moving
interfaces. This process produces new stationary interfaces with rate
of the order $-dM/dt$. Thus, the stationary interface density
satisfies the rate equation $\dot S=\dot P_0-\dot M$.  Combining this
equation with $P_0(t)\sim t^{-3/8}$\cite{derrida} and $M(t)\sim
t^{-1/2}$, we find that the surviving interfaces provide the dominant
contribution while those created in the process $M+M\rightarrow S$
contribute only to a correction of the order $t^{-1/8}$.  Thus a
two-scale structure similar to Eq.~(\ref{Eq-5}), with
$\langle\ell(t)\rangle\sim t^{3/8}$ and $\langle{\cal L}(t)\rangle\sim
t^{1/2}$, emerges. Hence, it is seen that the coarsening exponents can
be sensitive to the details of the interaction rules.


In conclusion, we investigated coarsening in a one-dimensional model
of competing species.  Interestingly, the coarsening properties of
this model may depend on the details of the dynamics.  The spatial
patterns are characterized by the existence of two characteristic
length scales, the average length of the single-species domains,
$\langle\ell(t)\rangle\sim t^{\alpha}$, and the average length of
superdomains, $\langle{\cal L}(t)\rangle\sim t^{\beta}$ (the
corresponding coarsening and mutation exponents are summarized in
Table 1).  This unusual behavior is an example of scaling violation in
a system with a scalar order parameter (similar behavior has been
observed in a few one-dimensional systems with vector order parameter
\cite{bray,xy}).  The above results contrasts the mean-field infinite
dimension limit, indicating that fluctuations play an important role
in sufficiently low spatial dimensions. It will be interesting to
establish whether similar behavior underlies other competing
population systems as well.

\vspace{.2in}
\centerline{
\begin{tabular}{|l|c|c|c|c|}
\hline
N&$\alpha$&$\beta$&$\nu$&$\theta$\\ 
\hline
3 (parallel)  &\,1/2\,&\,1\,&\,1/2\,&\,1\,\\
3 (sequential)&\,3/4\,&\,1\,&\,1/4\,&\,1\,\\
4             &\,1/3\,&\,2/3\,&\,1/3\,&\,1/3\,\\
4 (symmetric)  &\,3/8\,&\,1/2\,&\,1/2\,&\,3/8\,\\
\hline
\end{tabular}}
\vspace{.1in} {\small {\bf Table 1}: Coarsening
($\langle\ell(t)\rangle\sim t^{\alpha}$ and $\langle{\cal L}(t)\rangle\sim
t^{\beta}$) and persistence ($\langle n(t)\rangle\sim t^{\nu}$ and 
$P_0(t)\sim t^{\theta}$) exponents in 1D.}
\vspace{.2in}
  
We thank S.~Ispolatov, G.~Mazenko, J.~Percus, and
S.~Redner for discussions.  L.~F. is supported by the Swiss NSF,
P.~L.~K. is supported in part by a grant from NSF.  
E.~B.~ is supported in part by NSF Award Number 92-08527, and by
the MRSEC Program of the NSF under Award Number DMR-9400379.

\begin{thebibliography}{99}

\bibitem{bray}
      A.~J.~Bray, {\sl Adv. Phys.} {\bf 43}, 357 (1994).

\bibitem{bramson} 
 M.~Bramson and D.~Griffeath, {\sl Ann. Prob.} {\bf 17}, 26
      (1989).

\bibitem{fisch}
      R.~Fisch, {\sl Physica D}{\bf 45}, 19 (1990); 
      {\sl Ann. Prob.} {\bf 20}, 1528 (1992).

\bibitem{tainaka}
      K.~Tainaka, {\sl J. Phys. Soc. Jpn.} {\bf 57}, 2588 (1988);
      {\sl Phys. Rev. Lett.} {\bf 63}, 2688 (1989); {\sl Phys. Rev. E} 
      {\bf 50}, 3401 (1994).

\bibitem{lig}    
      T.~M.~Liggett, {\it Interacting Particle Systems}
      (Springer, New York, 1985).

\bibitem{elp}    
      E.~Ben-Naim, L.~Frachebourg, and P.~L.~Krapivsky, {\sl Phys. Rev. E} 
      {\bf 53}, 3078 (1996). 

\bibitem{glauber} 
      R.~J.~Glauber, {\sl J. Math Phys.} {\bf 4}, 294 (1963). 


\bibitem{Elskens}   
      Y.~Elskens and H.~L.~Frisch, {\sl Phys. Rev. A} {\bf 31}, 3812 (1985).

\bibitem{brk}
      E.~Ben-Naim, S.~Redner, and P.~L.~Krapivsky, {\sl J.~Phys.~A}, in press.

\bibitem{lpe}    
      L.~Frachebourg, P.~L.~Krapivsky, and E.~Ben-Naim, unpublished.

\bibitem{derrida}
      B.~Derrida, {\sl J.~Phys.~A} {\bf 28}, 1481 (1995);
      B.~Derrida, V.~Hakim, and V.~Pasquier, {\sl Phys. Rev. Lett.} 
      {\bf 75}, 751 (1995).

\bibitem{xy}
      A.~D.~Rutenberg and A.~J.~Bray,  {\sl Phys. Rev. Lett.} 
      {\bf 74}, 3836 (1995).


\end{thebibliography}

\end{multicols} 
\end{document}
