# Particle Systems with Stochastic Passing

## Abstract

We study a system of particles moving on a line in the same direction. Passing is allowed and when a fast particle overtakes a slow particle, it acquires a new velocity drawn from a distribution $`P_0(v)`$, while the slow particle remains unaffected. We show that the system reaches a steady state if $`P_0(v)`$ vanishes at its lower cutoff; otherwise, the system evolves indefinitely.

PACS numbers: 02.50-r, 05.40.+j, 89.40+k, 05.20.Dd

We are interested in behavior of a system of interacting particles moving on the real line in one direction, say to the right. The system is endowed with the following simple dynamics: (i) Particles move freely between “collisions”; (ii) After a collision, or “passing” event, the velocity of the slow particle remains the same, $`v_{\mathrm{slow}}`$=const, while the fast particle instantaneously acquires some new velocity, $`v_{\mathrm{new}}>v_{\mathrm{slow}}`$, drawn from the intrinsic velocity distribution $`P_0(v)`$. We would like to answer the following basic questions about the behavior of the system: Does the velocity distribution $`P(v,t)`$ reach a steady state or the system continues to evolve indefinitely? How does the average velocity depend on time? etc.

Our motivation is primarily conceptual, as we want to understand non-equilibrium infinite particle systems with two-body interactions. Thus we have chosen the simplest dynamics – interactions occur only upon colliding, and only one particle is affected. The appealing simplicity of the model suggests that it might show itself in different natural phenomena, and indeed, we originally arrived to this model in an attempt to mimic traffic on one-lane roads. Somewhat related dynamics were already used in modeling voting systems, force fluctuations in bead packs, asset exchange processes, combinatorial processes, continuous asymmetric exclusion processes, granular gases, and aggregation-fragmentation processes.

Let us first consider discrete velocity distributions. Specifically, we assume that both initial velocities and new velocities are drawn from the same intrinsic distribution $`P_0(v)=p_j\delta (vv_j)`$. For the binary distribution, the system does not evolve at all, so the first non-trivial case is the ternary intrinsic distribution when the system contains slow, moderate, and fast particles. Initially,

$$P_0(v)=p_1\delta (vv_1)+p_2\delta (vv_2)+p_3\delta (vv_3).$$
(1)   
We set

$$p_1+p_2+p_3=1,v_1<v_2<v_3,$$
(2)   
without a loss of generality. When the steady state is reached, the velocity distribution remains ternary,

$$P_{\mathrm{eq}}(v)=p_1\delta (vv_1)+q_2\delta (vv_2)+q_3\delta (vv_3).$$
(3)   
The density of slow particles does not change, while the densities $`q_2`$ and $`q_3`$ of moderate and fast particles differ from the initial values. The final densities are found from a simple probabilistic argument based on the requirement of stationarity. For moderate particles we get

$$p_3(v_2v_1)q_2=p_2(v_3v_1)q_3.$$
(4)   
The left-hand side of Eq. (4) gives the loss in $`q_2`$ which happens when a moderate particle overtakes a slow particle and becomes a fast particle. The right-hand side gives the gain in $`q_2`$ which takes place when a fast particle overtakes a slow particle and converts into a moderate particle. Solving Eq. (4) together with the normalization condition, $`p_1+q_2+q_3=1`$, we find

$$q_2=p_2\frac{p_2+p_3}{p_2+\nu p_3},q_3=\nu p_3\frac{p_2+p_3}{p_2+\nu p_3}.$$
(5)   
where $`\nu =(v_2v_1)/(v_3v_1)`$. Since $`\nu <1`$, we have $`q_2>p_2`$ and $`q_3<p_3`$. Thus the density of moderate particles increases while the density of fast particles decreases. Similarly, one can analyze discrete velocity distributions with more than three particle species. In all cases (i) the system reaches a steady state; (ii) the average velocity decreases and eventually reaches some finite value; (iii) the density of the most slow particle species remains unchanged.

Now we consider a continuous velocity distributions. Let $`[v_{\mathrm{min}},v_{\mathrm{max}}]`$ be a support of the intrinsic velocity distribution $`P_0(v)`$. By Galilean transform, we can set $`v_{\mathrm{min}}=0`$ without loss of generality. We consider unbounded distributions, $`v_{\mathrm{max}}=\mathrm{}`$, although main results equally apply to the cases with finite $`v_{\mathrm{max}}`$.

The passing rule asserts that when a fast particle overtakes a slow particle moving with a velocity $`v_{\mathrm{slow}}`$, the assignment of the new velocity $`v`$ occurs with probability

$$P_0\left(v|v_{\mathrm{slow}}\right)=P_0(v)\frac{\theta (vv_{\mathrm{slow}})}{_{v_{\mathrm{slow}}}^{\mathrm{}}𝑑v^{}P_0(v^{})}.$$
(6)   
Eq. (6) guarantees that $`v>v_{\mathrm{slow}}`$ and that the normalization requirement, $`𝑑vP_0\left(v|v_{\mathrm{slow}}\right)=1`$, is obeyed.

Now we can write a Boltzmann equation for the velocity distribution $`P(v,t)`$:

$`{\displaystyle \frac{P(v,t)}{t}}=P(v,t){\displaystyle _0^v}𝑑v^{}(vv^{})P(v^{},t)`$ (7)   
$`+{\displaystyle _0^v}𝑑v_2P_0(v|v_2){\displaystyle _{v_2}^{\mathrm{}}}𝑑v_1(v_1v_2)P(v_1,t)P(v_2,t).`$ (8)   
The first term on the right-hand side of Eq. (7) describes loss in $`P(v,t)`$ due to collisions with slow particles: Collisions occur with rate proportional to velocity difference, and the integration limits ensure that only collisions with slower particles are taken into account. The second, gain, term accounts for the increase of $`P(v,t)`$ due to a random assignment of velocity $`v`$ after collision.

We could not solve Eq. (7) in the general case of an arbitrary intrinsic velocity distribution $`P_0(v)`$. Attempts to find a solution even for some particularly simple $`P_0(v)`$, e.g., linear, exponential, or uniform, turned out to be fruitless as well. Thus we proceed by employing asymptotic, approximate, and numerical techniques.

We start by looking at the asymptotic behavior of $`P(v)`$ in the small velocity limit. Let $`vu(t)`$, where $`u(t)`$ is the average velocity,

$`u(t)v={\displaystyle _0^{\mathrm{}}}𝑑vvP(v,t).`$ (9)   
Then Eq. (7) simplifies to

$`{\displaystyle \frac{P(v,t)}{t}}`$ $`=`$ $`P_0(v)u(t){\displaystyle _0^v}𝑑v_2P(v_2,t)`$ (10)   
$``$ $`P(v,t){\displaystyle _0^v}𝑑v^{}(vv^{})P(v^{},t).`$ (11)   
To probe the small $`v`$ behavior, we need to know $`P_0(v)`$ at $`v0`$. Let us consider a family of intrinsic velocity distributions that behave algebraically:

$$P_0(v)Av^\mu \mathrm{when}v0.$$
(12)   
Now assume that the system reaches the steady state: $`P(v,t)P_{\mathrm{eq}}(v)`$ and $`u(t)u_{\mathrm{eq}}`$. Plugging these and Eq. (12) into Eq. (10), we find that the velocity distribution also behaves algebraically in the small velocity limit,

$$P_{\mathrm{eq}}(v)(\mu +1)Au_{\mathrm{eq}}v^{\mu 1}\mathrm{when}v0.$$
(13)   
In other words, the steady state velocity distribution scales as $`v^1P_0(v)`$. Recalling the normalization requirement, $`_0^{\mathrm{}}𝑑vP(v)=1`$, we see that it is possible only when $`\mu >0`$. Thus our assumption that the system reaches a steady state is certainly wrong when $`\mu 0`$. In this region, we anticipate that the system will evolve indefinitely. Note that both the exponential and uniform intrinsic distributions belong to the borderline case of $`\mu =0`$ that separates stationary and evolutionary regimes; for them an anomalously slow kinetics is anticipated.

To probe the behavior of evolving systems, we assume that in the long-time limit there is a very small fraction of “active” particles that move with velocities $`v1`$ and the vast majority of “creeping” particles that hardly move at all. We ignore collisions between active particles since their density is very low. We also ignore collisions between creeping particles since their relative velocity is very small. This picture suggests that only collisions between active and creeping particles matter. Thence, the velocity distribution of active particles obeys

$$\frac{P(v,t)}{t}=P_0(v)u(t)vP(v,t).$$
(14)   
Eq. (14) may at best describe the evolution process in the long-time limit. However, for the sake of tractability, we apply it to the whole time range and use the natural initial condition $`P(v,0)=P_0(v)`$. Eq. (14) is an inhomogeneous linear differential equation which is easily solved to give

$$P(v,t)=P_0(v)e^{vt}\left[1+_0^t𝑑t^{}u(t^{})e^{vt^{}}\right].$$
(15)   
This solution implies

$$P(v,t)u(t)v^1P_0(v)\mathrm{for}v\frac{1}{t},$$
(16)   
which resembles Eq. (13).

To close the solution of Eq. (15), we must determine $`u(t)`$. It is possible to plug (15) into the definition of the average velocity, Eq. (9), and get an integral equation for $`u(t)`$. In the following we use another approach, which is technically simpler. Note that the density of active particles, $`𝑑vP(v,t)`$, is manifestly conserved by Eq. (14). After integration over velocity, Eq. (15) becomes

$$1=\widehat{P}_0(t)+_0^t𝑑t^{}u(t^{})\widehat{P}_0(tt^{}),$$
(17)   
where $`\widehat{P}_0`$ is the Laplace transform of the intrinsic velocity distribution,

$$\widehat{P}_0(t)=_0^{\mathrm{}}𝑑vP_0(v)e^{vt}.$$
(18)   
One can guess the long time behavior of $`P(v)`$ without actually solving Eq. (17). Let us assume that the average velocity varies slowly with $`t`$. Then the integral on the right-hand side of Eq. (17) can be estimated as $`u(t)_0^t𝑑t^{}\widehat{P}_0(t^{})`$, which implies

$$u(t)\left[_0^t𝑑t^{}\widehat{P}_0(t^{})\right]^1.$$
(19)   
For an intrinsic velocity distribution with an algebraic behavior (12) in the small-$`v`$ limit, we have $`\widehat{P}_0(t)t^{1\mu }`$ for large $`t`$. Hence $`_0^t𝑑t^{}\widehat{P}_0(t^{})t^\mu `$ for $`\mu <0`$, and it follows from Eq. (19) that $`u(t)t^\mu `$. The above derivation is quite careless, though the final result is correct. Now we derive this result in a more rigorous way.

The convolution form of the integral in Eq. (17) suggests to apply the Laplace transform once more. It yields

$$\frac{1}{s}=_0^{\mathrm{}}𝑑v\frac{P_0(v)}{s+v}+\widehat{u}(s)_0^{\mathrm{}}𝑑v\frac{P_0(v)}{s+v}.$$
(20)   
Here

$$\widehat{u}(s)=_0^{\mathrm{}}𝑑tu(t)e^{st}$$
(21)   
is the Laplace transform of the average velocity. The double Laplace transform of the intrinsic velocity distribution has been simplified:

$$_0^{\mathrm{}}𝑑te^{st}_0^{\mathrm{}}𝑑vP_0(v)e^{vt}=_0^{\mathrm{}}𝑑v\frac{P_0(v)}{s+v}.$$
(22)   
Thus the Laplace transform of the average velocity is

$$\widehat{u}(s)=1+\left[s_0^{\mathrm{}}𝑑v\frac{P_0(v)}{s+v}\right]^1.$$
(23)   
Generally, one cannot obtain more explicit results. Given that the above approach describes only the long-time asymptotics, let us focus on this regime. To probe the long-time behavior, one should determine the small $`s`$ asymptotics of $`\widehat{u}(s)`$. For algebraic intrinsic velocity distributions (12), the asymptotics of (22) reads

$`{\displaystyle _0^{\mathrm{}}}𝑑v{\displaystyle \frac{P_0(v)}{s+v}}As^\mu {\displaystyle _0^{\mathrm{}}}𝑑w{\displaystyle \frac{w^\mu }{w+1}}={\displaystyle \frac{A\pi }{\mathrm{sin}(\pi \mu )}}s^\mu .`$ (24)   
This applies for $`1<\mu <0`$ (the lower bound comes from the normalization requirement, $`𝑑vP_0(v)=1`$). Plugging (24) into (23) yields

$$\widehat{u}(s)\frac{\mathrm{sin}(\pi \mu )}{A\pi }s^{1\mu }\mathrm{for}s0,$$
(25)   
and by making the inverse Laplace transform, we finally arrive at

$$u(t)\frac{\mathrm{sin}(\pi \mu )}{A\pi \mathrm{\Gamma }(1+\mu )}t^\mu \mathrm{for}t\mathrm{}.$$
(26)   
This result agrees with the asymptotics we naively derived earlier.

A special consideration is required for the borderline case of $`\mu =0`$. For concreteness, consider the exponential intrinsic distribution, $`P_0(v)=\mathrm{exp}(v)`$. Its double Laplace transform reads

$`{\displaystyle _0^{\mathrm{}}}𝑑v{\displaystyle \frac{e^v}{s+v}}=e^sE_1(s),`$   
where

$`E_1(s)={\displaystyle _1^{\mathrm{}}}𝑑x{\displaystyle \frac{e^{xs}}{x}}`$   
is the exponential integral. As a result, Eq. (23) becomes

$$\widehat{u}(s)=1+\frac{1}{se^sE_1(s)}.$$
(27)   
Using the well-known asymptotics of the the exponential integral, $`E_1(s)=\mathrm{ln}s\gamma +𝒪(s)`$, (where $`\gamma 0.5772`$ is Euler’s constant), we transform Eq. (27) into

$$\widehat{u}(s)=\frac{1}{s(\mathrm{ln}s+\gamma )}+𝒪\left(\frac{1}{\mathrm{ln}s}\right).$$
(28)   
Performing the inverse Laplace transform gives

$$u(t)\frac{1}{\mathrm{ln}t}\mathrm{for}t\mathrm{}.$$
(29)   
To summarize, for the family of intrinsic velocity distribution with algebraic behavior in the small $`v`$ limit (12), our predictions for the long-time asymptotics of the average velocity $`u(t)`$ are:

$$u(t)\{\begin{array}{cc}\mathrm{const},\hfill & \text{for }\mu >0\text{;}\hfill \\ (\mathrm{ln}t)^1,\hfill & \text{for }\mu =0\text{;}\hfill \\ t^\mu ,\hfill & \text{for }1<\mu <0\text{.}\hfill \end{array}$$
(30)   
To check the validity of asymptotic predictions and, more generally, to see if the mean-field theory is applicable at all, we perform molecular dynamics simulations and solve the Boltzmann equation (7) numerically. To sample distinct regimes predicted in (30) we consider the intrinsic velocity distribution

$`P_0(v)={\displaystyle \frac{v^\mu e^v}{\mathrm{\Gamma }(\mu +1)}}`$ (31)   
with $`\mu =1,0,1/2`$.

In molecular dynamics simulations, we place $`N`$ particles onto the ring of length $`L=N`$ so that the average density is equal to one. Most of our simulations are performed for $`N=510^4`$ particles, but we also simulated twice larger system and found no appreciable difference. Initially, particle velocities are randomly drawn from the distribution $`P_0(v)`$. The model is updated according to the collision-time-list algorithm suggested in Ref..

To solve the Eq. (7) numerically, we use Euler’s time update with both uniform and non-uniform grid; in the latter case, we take $`v_N=(N/N_{max})^4v_{max}`$ velocity grid with $`v_{max}=15`$ and $`N_{max}=300500`$. Integrals on the right-hand side of Eq. (7) are calculated using the trapezoid rule; time increment $`\delta t=0.1`$ was found to be suitable for all three $`P_0(v)`$.

The main conclusion is that results of molecular dynamics simulations and numerical solutions of the mean-field equation are virtually identical (see, e.g., Fig. 1). It confirms our assumption that the system remains well-stirred and no appreciable spatial correlations develop. Additionally, Fig. 1 shows that for $`P_0(v)=ve^v`$, the approach of $`P(v,t)`$ to the steady state is non-uniform in velocity. This is caused by the obvious fact that for any finite time, the velocity distribution $`P(v,t)`$ must still vanish at the lower cutoff as $`P_0(v)`$ does. In other words, the steady state (13) is reached outside the “boundary layer”, $`vv_{}(t)`$, while within the boundary layer, $`vv_{}(t)`$, the velocity distribution continues to evolve. The threshold velocity $`v_{}(t)`$ is estimated by evaluating the first, leading term in the right-hand side of Eq. (10): $`t^1Pv_{}^{\mu +1}P`$, which implies $`v_{}t^{1/(\mu +1)}`$. The width of the boundary layer shrinks with time but the boundary layer still exists ad infinitum.

Figs. 2–4 plot the average velocity vs. time for the intrinsic velocity distributions (31) with $`\mu =1,0,1/2`$, respectively. We find good agreement with the theoretical prediction of Eq. (30) when $`\mu 0`$. On Fig. 5 the plot of the local exponent $`\alpha (t)d\mathrm{ln}[u(t)]/d\mathrm{ln}[t]`$ vs. $`t^{1/2}`$ is shown for $`\mu =1/2`$. The results of extrapolation of $`\alpha (t)`$ to the $`t\mathrm{}`$ limit are not contradicting the theoretical prediction, $`\alpha =1/2`$.

In summary, we have shown that the fate of the system of passing particles is determined by the behavior of the intrinsic velocity distribution near its lower cutoff: If $`P_0(v)`$ vanishes in this limit, the system reaches a steady state; otherwise, the evolution continues forever. Comparison between solutions of the mean-field Boltzmann equation and results of molecular dynamics simulations suggests that the mean-field theory description is exact. It will be interesting to confirm this result rigorously.

We gratefully acknowledge partial support from NSF and ARO.