Lattice gas with nearest-neighbor exclusion in a shear-like field

Abstract

We present Monte Carlo simulations of the lattice gas with nearest-neighbor exclusion and Kawasaki (hopping) dynamics (hard square lattice gas), under the influence of a nonuniform drive, on the square lattice. The drive, which favors motion along the +x axis and inhibits motion in the opposite direction, varies linearly in the y direction. Our lattice has rigid walls at the end points in the y direction and periodic boundaries along the drive. We find that this model has transition to a sublattice-ordered phase at a density of about 0.298, lower than in equilibrium (rhoc <FONT FACE=Symbol>@</FONT> 0.37), but somewhat higher than in the uniformly driven case at maximal bias (rhoc <FONT FACE=Symbol>@</FONT> 0.272). For smaller global densities (r < 0.33), the ordering occurs with particle accumulation in the low-drive region. Above this density we observe a surprising reversal in the density profile, with particles migrating to the high-drive region.

Nearest-neighbor exclusion; Monte Carlo; Non-equilibrium lattice gas; Shear drive

Lattice gas with nearest-neighbor exclusion in a shear-like field

Fabricio Q. Potiguar; Ronald Dickman

Departamento de Física, ICEx, Universidade Federal de Minas Gerais, 30123-970, Belo Horizonte, MG, Brazil

ABSTRACT

We present Monte Carlo simulations of the lattice gas with nearest-neighbor exclusion and Kawasaki (hopping) dynamics (hard square lattice gas), under the influence of a nonuniform drive, on the square lattice. The drive, which favors motion along the +x axis and inhibits motion in the opposite direction, varies linearly in the y direction. Our lattice has rigid walls at the end points in the y direction and periodic boundaries along the drive. We find that this model has transition to a sublattice-ordered phase at a density of about 0.298, lower than in equilibrium (r_c@ 0.37), but somewhat higher than in the uniformly driven case at maximal bias (r_c@ 0.272). For smaller global densities (r < 0.33), the ordering occurs with particle accumulation in the low-drive region. Above this density we observe a surprising reversal in the density profile, with particles migrating to the high-drive region.

Keywords: Nearest-neighbor exclusion; Monte Carlo; Non-equilibrium lattice gas; Shear drive

I. INTRODUCTION

The lattice gas with nearest-neighbor exclusion (NNE) is the infinite repulsion limit of the ordinary repulsive lattice gas. Here, particles are forbidden to occupy the same or neighboring sites. Its equilibrium version, in the square lattice, has a continuous order-disorder phase transition at the critical density r_c» 0.37 [1-5] which belongs to the Ising universality class. This model is also known as the hard square lattice gas. The square lattice, like other bipartite lattices, may be divided into two sublattices, A and B, such that the neighbors of all sites in A belong to B and vice-versa. For densities above r_c the two sublattices have unequal densities.

More recently, nonequilibrium versions of the NNE model were studied [6,7]. It was found [6] that the critical density varies with drive intensity: the higher the drive, the lower the critical density. The transition is continuous for low bias but becomes first order if the bias strength is > 0.75. Above the transition density, the system separates into regions of low and high local density, with the high-density region essentially frozen. Szolnoki and Szabó [7] extended the dynamics to include next-nearest-neighbor (diagonal) hops, and observed a similar variation of the critical density with drive strength, but with a homogeneous stationary state. Continuous phase transitions in this version of the model fall in the Ising class, as the equilibrium case. These driven lattice gases are known as driven diffusive models (DDS) [8,9].

Here we consider a hard-core DDS with nearest-neighbor hopping dynamics, in which the drive is nonuniform. Specifically, the probability for a particle at y to attempt a jump to the right (x ® x+1) is given by:

for y = 1, 2,...,L, on a square lattice of L² sites (a jump to the left occurs with probability P_l(y) = 1/2-P_r(y)). Jumps in the ±y direction occur with probability 1/4. We present numerical results from Monte Carlo (MC) simulations of the model.

Our main objective is to obtain the critical properties of the model. To determine the critical density, we study the behavior of the order parameter and the stationary current as functions of global density. Of particular interest are the current and the density profiles as functions of y which show how the system organizes under the influence of the drive. In the following section (II) we detail the model and simulation procedure. In Section III, we present numerical results and discussions. Final considerations are reserved for section IV.

II. SIMULATIONS

We consider a square lattice of length L (L² sites) with N particles (N < L²/2). (Most of our studies use L = 100.) The initial configuration is prepared via random sequential adsorption (RSA) [10,11] of particles, always respecting the excluded-volume condition. The dynamics consists of random selection of particles and assignment of a new (trial) position at one of the nearest neighbor sites with the probabilities introduced above. If the trial position does not violate the exclusion constraint, the move is accepted. Each MC time unit corresponds to N attempted moves. We perform averages over 10⁶ MC time units in the stationary state. Twenty values of the density r = N/L², ranging from r = 0.01 to r = 0.37 are studied.

III. RESULTS

The order parameter is defined as the difference in sublattice occupancies per particle,

where N_A(B) is the number of particles in sublattice A(B) and á ñ indicates an average over realizations. The stationary current is defined as the difference between the number of jumps along the drive less the number contrary to it, per site and unit time.

In Fig. 1 we show the average value of the order parameter and, in the inset, its corresponding fluctuation. The rapid rise in the order parameter and the associated maximum in its variance signal a continuous transition to an ordered phase. The data suggest a critical density of r_c = 0.298, although a precise determination will require studies of larger systems. In any case, it is clear that the critical density is lower than that for the equilibrium transition [2]. The apparently continuous transition in the presence of a nonuniform drive is likely to be due the fact that the system does not order all at once.

In Fig. 2, we show the stationary current for several system sizes. This quantity displays the same behavior as in the uniformly driven case: it increases at small densities (reflecting the increasing number of carriers) and decreases for larger densities (due to the reduction in available space for movement). The maximum value of áj(t)ñ falls at roughly in the same density as in the uniform drive case [6]. Interestingly, the phase transition near r = 0.298 is associated with a plateau in the current (see inset of fig. 2).

The order parameter and the stationary current present strong fluctuations for densities above r = 0.32. The evolution of these quantities typically displays sudden jumps between the ordered and the disordered state, a fact already observed in the uniformly driven case. In these cases, the relaxation time to the stationary state is on the order of 10⁷ MC steps. When the global density approaches the maximal value we studied here, this slow relaxation is suppressed due to the quick formation of a jammed region in the high-field portion of the lattice (see below the results for the profiles and the configuration snapshots). The drive provokes formation of organized structures while their thermal motion provides a mechanism for breaking such clusters. We will see that such behavior may find a parallel in actual physical systems.

Of interest is how the system organizes as the global density increases. To understand this, we determine the stationary density and current profiles, r(y) and j(y), respectively. These quantities are shown in Figs. 3 and 4.

Particles concentrate first, as the global density increases, in the low field region (this concentration is observed at densities as low as 0.20). This is surprising given the finding that a strong drive favors order [6]. On this basis one might expect particles to concentrate first in the high-drive region. In fact just the opposite occurs: for a global densities r < 0.32, the density profile r(y) (Fig. 3) is highly skewed to the region around y = 1. Note that the local density is > 0.37 in this region (for r between 0.30 and 0.32), that is, greater than or equal to r_c in equilibrium. The density profile decays monotonically with increasing y (except for small density oscillations induced by the wall at y = 1). For r = 0.33, the local density instead increases with y, reaching a peak near y = 47, after which it decays in an approximately linear fashion until y = L.

The density profiles show that particles tend to accumulate in the layers in contact with the rigid walls (y = 1 and y = L). This accumulation is due to excluded volume effects [12] that entropically favor enhanced densities at rigid surfaces. The effect of this accumulation shows up in the current profile (Fig. 4), as a sharp peak at y = L.

The observations made in the density profiles have reflections at the current profiles (Fig. 4). The current is much smaller in the small-y region (low-bias) and increases monotonically with y. It clearly shows that ordering occurs in the low-drive region for densities between 0.30 and 0.32, as described above. Note also that the current in the high-drive region is close to zero for r > 0.35, signaling a greater accumulation of particles in this part of the lattice.

A possible explanation for the surprising reversal of the density and current profiles with increasing global density is related to the formation of a jammed region, as observed in the uniformly driven system [6]. We term a region "jammed" when all particles in this region have their movement severely reduced due to aggregation induced by the field, leading to a vanishing diffusion coefficient and mobility. When the global density is too low for such a region to form, particles tend to collect in the low-drive region because a strong drive tends to destroy the local correlations needed for particles to pack to high density, even if such packing does not result in long-range order. The depletion of the high-drive region appears to be the reason for the plateau in the current observed around r = 0.30 in Fig. 2. Fig. 3 shows that as the global density increases, the local density in the high-drive region remains nearly constant, so the global current hardly varies, since its main contribution comes from particles in the high-field region.

When, on the other hand, the global density is sufficiently high, an irreversible accumulation of particles happens in the high-field region, so that the low-drive region has fewer particles than at lower global densities, for which there is no jammed region.

We may now identify two factors leading to the continuous variation of the order parameter with density shown in Fig. 1. One reason is that ordering begins in the low-bias region. Studies of the uniformly driven system show that the transition is continuous under a weak bias [6]. The second reason is that the width of the ordered region grows continuously with increasing density.

To illustrate these ideas, we show in Fig. 5 a configuration for r = 0.31 and L = 100. As expected, the low-drive region is very dense and contains few mobile particles. In the uniformly driven system (at maximum bias) one observes, at this density, formation of "herringbone" pattern of diagonal stripes, pointing along the drive, with particles in this structure essentially frozen. In the present case the low-drive region is highly ordered, with almost all particles in the same sublattice, but there is no sign of the herringbone pattern. The high-drive region is disordered, permitting the high currents and lower densities reported above. Several clusters of particles exist in the high-drive region, but they are not large enough to cause jamming.

We show in fig. 6 a typical configuration at global density r = 0.35. Evidently the long diagonal line of particles at the upper right is associated with jamming in the high-drive region. The empty triangular region implies a decrease in the local density with increasing y. Particles are not free to enter this region since all particles along the diagonal edge are immobile. These observations are supported by the the density and current profiles (Figs. 3 and 4). The density is roughly constant in the middle portion of the lattice and begins to decrease near y = 78, where the empty triangular region begins. The current is only appreciably different from zero in the lower portion of the lattice, as signalled in Fig. 6 by the presence of particles both sublattices. The diagonal edges observed in configurations at this density (always in the high-drive region), are extremely long-lived structures, since only the particle at the tip of the line can move without violating the exclusion constraint. In this configuration, such particles are blocked by others, providing a virtually infinite lifetime to this structure. The large number of voids in the high-bias region of fig. 6 results from the formation of the jammed structure during the early stages of the dynamics (confirmed by the histogram of displacements in y during thermalization, not shown here, [13]).

To study correlations between the particles, we determine the radial distribution function, g(r) in the high- and low-drive regions (Fig. 7). This function is proportional to the probability of finding a pair of particles separated by a distance r, and is normalized so that g ® 1 as r ® ¥. For purposes of determining g(r), the low-drive region is taken as the strip 6 < y < 14, while the high-drive region comprises 86 < y < 94.

The g(r) curves for global density r = 0.31 show that the high- and low-drive structures are markedly different. In the low-bias region the peaks are much larger due to the sublattice ordering associated with packing of particles (as evidenced by the configuration in Fig. 5), and is compatible with the existence of long-range order. The high-drive region shows little structure; the oscillations in g(r) decay rapidly with distance. The picture for r = 0.35 is quite different. The sharpness of the peaks in the curve for r = 0.35 in the high-bias region reflects the very different sublattice densities, as does the fact that g @ 0 for r = 3, and . The radial distribution function in the low-drive region, for this density, exhibits less structure, indicating the more equal sublattice occupancies.

A final question is whether all the observed ordering is somehow induced by the rigid walls. We studied a system with periodic boundaries in y direction and a symmetric drive profile and found essentially the same results for the current and density profiles [13].

IV. CONCLUSIONS

We studied a lattice gas with nearest-neighbor exclusion driven by a nonuniform, shear-like drive, on the square lattice, under nearest-neighbor hopping dynamics. The problem is of interest as an example of the surprising behavior to be found in a simple nonequilibrium system. We find that the model undergoes a continuous order-disorder transition at a critical density of about r_c = 0.298. This is unlike the uniformly driven model, in which the transition is discontinuous for a bias > 0.75. The stationary current follows roughly the same trends as in the uniformly driven case, but exhibits a plateau in the neighborhood of the phase transition.

Our results show that this transition is due to the concentration of particles at the low-bias region, for global densities between r = 0.30 and 0.32. Remarkably, the nonuniform drive induces a highly nonuniform density profile, expelling particles from the high-bias region. The effect is sufficiently strong to induce sublattice ordering in the low-bias region. Thus the drive favors a class of configurations that, on the basis of entropy maximization, are extremely unlikely. Note that at these densities there is no jamming, i.e., the system is ergodic. Migration of particles to the low-bias region appears to derive from the destruction of short-range correlations (required for efficient packing), by the drive.

For higher densities, we observe a completely inverted picture, with formation of jammed structures in the high-drive region, while particles outside this region remain mobile. The jammed region is characterized by a dense (r > 0.37) strip of particles; at higher global densities this region displays long diagonal chains of particles associated with voids. Other variants of the model, with different drive geometries [14] show similar effects.

The parallel between our lattice gas and real systems is built by the observation of the slow relaxation at high densities. Our particles have (infinite) repulsive interactions and thermal motion (provided by the Monte Carlo algorithm). A thermal system composed of particles that have highly repulsive interactions, albeit not infinite ones, is a colloidal suspension. In fact, the experiments of Bertrand et al. [15] show a similar effect to the slow relaxation observed here. They studied suspensions at several densities under the influence of shear. The suspension presents shear-thickening at intermediate densities, where a small vibration drives the system back to a fluid state. At higher densities, after shear is applied, the suspension forms a paste, becoming trapped in this jammed state. In our case, thermal motion, or vibration, is always present, as well as shear, so we do not observe the breakdown of the jammed state at intermediary densities, but a slow relaxation towards a denser state. In a granular system, which is, by definition, athermal, effects analogous to those of thermal agitation can be produced by shaking. This raises the possibility that the behavior identified in the sheared lattice gas might also be observed in a sheared packing if, besides the shear drive, continuous shaking were applied to the grains. This suggests that the model studied here can be extended to study the dynamics of certain complex fluids, a subject we intend to explore in future work.

Acknowledgments

We thank the Brazilian agencies CAPES, CNPq and Fapemig for support.

[14] Studies of a system with half zero and half maximum drive, and a driven system with variable density are in progress.

Received on 7 September, 2005

• [1] D. S. Gaunt and M. E. Fisher, J. Chem. Phys. 43, 2840 (1965).
• [2] F. H. Ree and D. A. Chesnut, J. Chem. Phys. 45, 3983 (1966).
• [3] L. K. Runnels and L. L. Combs, J. Chem. Phys. 45, 2482 (1966).
• [4] D. S. Gaunt and M. E. Fisher, J. Chem. Phys. 46, 3237 (1967).
• [5] R. J. Baxter, I. G. Enting, and S. K. Tsang, J. Stat. Phys. 22, 465 (1980).
• [6] R. Dickman, Phyr. Rev. E 64, 16124 (2001).
• [7] A. Szolnoky and G. Szabó, Phys. Rev. E 65, 47101 (2002).
• [8] B. Schmittmann and R. K. P. Zia, Statistical Mechanics of Driven Diffusive Systems, vol. 17, Phase Transitions and Critical Phenomena (Academic Press, London, 1995).
• [9] J. Marro and R. Dickman, Nonequilibrium phase transitions in lattice models (Cambridge University Press, Cambridge, 1999).
• [10] P. Meakin, J. L. Cardy, E. Loh, and D. J. Scalapino, J. Chem. Phys. 86, 2380 (1987).
• [11] R. Dickman, J.-S. Wang, and I. Jensen, J. Chem. Phys. 94, 8252 (1991).
• [12] R. Dickman, P. Attard, and V. Simonian, J. Chem. Phys. 107 , 205 (1997).
• [13] F. Q. Potiguar, R. Dickman, submitted to Phys. Rev. E (cond-mat/0507114).
• [15] E. Bertrand, J. Bibette, V. Schmitt, Phys. Rev. E 66, 60401 (2002).

Publication Dates

• Publication in this collection
  16 Oct 2006
• Date of issue
  Sept 2006

History

• Received
  07 Sept 2005