Phase Coexistence in Driven One Dimensional Transport

4 Pages·0.24 MB·English

by  A. Parmeggiani

Phase Coexistence in Driven One Dimensional Transport
A. Parmeggiani (1), T. Franosch (1), and E. Frey (1,2)
(1) Hahn-Meitner Institut, Abteilung Theorie, Glienicker Str. 100, D-14109 Berlin, Germany
(2) Fachbereich Physik, Freie Universit¨at Berlin, Arnimallee 14, D-14195 Berlin, Germany
(Dated: February 2, 2008)
3
0 We study a one-dimensional totally asymmetric exclusion process with random particle attach-
0
mentsanddetachmentsinthebulk. Theresultingdynamicsleadstounexpectedstationaryregimes
2
for large butfinite systems. Such regimes are characterized by aphase coexistence of low and high

n densityregionsseparatedbydomainwalls. Weuseamean-fieldapproachtointerpretthenumerical
a resultsobtainedbyMonte-Carlosimulationsandwepredictthephasediagramofthisnon-conserved
J dynamicsin thethermodynamic limit.

4
2 PACSnumbers: 02.50.Ey,05.40.-a,64.60.-i,72.70.+m

]
h Even some of the simplest driven diffusive systems in allowing particle attachment and detachment. We are
c
one dimension show surprisingly rich and complex be- interested in the limit where the kinetic rates are such
e
m havior which is rather unexpected when looked at with that the incoming and outgoing fluxes of particles at the
experience gained from equilibrium phenomena [1]. A boundaries and in the bulk are comparable. This im-
-
t particularly illuminating example are boundary-induced pliesthataparticle,injectedattheboundaryorattached
a
t phase transitions in driven one-dimensional (1D) trans- somewhere in the bulk, remains long enough on the lat-
s
port processes, such as the Totally Asymmetric Simple ticetomoveafinitefractionofthetotalsystemsize. New
.
t
a Exclusion Process (TASEP). The model, originally pro- phenomenaareexpectedintheregimeofcompetitionbe-
m posed in [2], consists of particles hopping unidirection- tween TASEP and LK for a large but finite system. Of
- ally with hard-core exclusion along a 1D lattice. Due course, the dynamics in an infinitely large system would
d to conservation of the particle current in the bulk, the becompletelydominatedbythebulkadsorptionanddes-
n
ratesofincomingoroutgoingparticlesattheboundaries orptionrates. Itturnsoutthatthepresenceofthekinetic
o
drive the systemto non-trivialstationarystates [3]. The ratessignificantly change the picture of TASEP,produc-
c
[ resulting phase diagram shows continuous and discon- ing a completely reorganized phase diagram. We shall

tinuous transitions of the average density of particles in showbycomputersimulationsandmean-fieldarguments
1
the limit of large system sizes. These results were ob- that,inthisnon-conserveddynamics,onecanhavephase
v
5 tainedfirstinmean-fieldtheoryandthenextendedwhen coexistence where low and high density phases are sepa-
7 a complete analytical solution was presented solving ex- rated by stable discontinuities in the density profile.
4 plicitely the recursion relations of the model or using a The model we discuss here is directly inspired by the
1
matrix product ansatz technique [4]. unidirectional motion of many motor proteins along cy-
0
toskeletal filaments [6]. Motors advance along the fil-
3 The TASEP is one out of many examples for driven
0 systems with stationary non-equilibrium states, which ament while attachment and detachment of motors be-
/ tweenthecytoplasmandthefilamentoccur[7]. Recently,
t cannotbedescribedintermsofBoltzmannweights. This
a ithasbeenshownthatsuchdynamicscanberelevantfor
has tobe contrastedwithprocesseslikethe bulk adsorp-
m
modelingthefilopodgrowthineukaryoticcellsproduced
tion/desorption kinetics of particles on a lattice coupled
- bymotorproteinsinteractingwithinactinfilaments[8].
d toareservoir(“LangmuirKinetics”,LK),whosestation-
n ary state is well described within standard concepts of
o equilibrium statistical mechanics. Here, particles adsorb particle vacancy
c
at an empty site or desorb from an occupied one with
: wD wA
v fixed respective kinetic rates obeying detailed balance. //
Xi The bulk density profile at equilibrium is described by a b
a Langmuir isotherm, determined solely by the ratio of FIG. 1: TASEPscheme with bulk attachment/detachment.
r
a thetwokinetic rates[5],asgivenbytheGibbs ensemble.
We consider a 1D lattice composed of sites i=1,...,N
Due to the presence of the particle reservoir there is no
(Fig.1). Theconfigurationsaredescribedintermsofoc-
conservation of particles and no net particle current in
cupation numbers n =1 for a site occupied by a particle
i
the bulk. It is interesting to ask what can be expected
and n =0 for an empty site (vacancy). The dynamics
i
in coupling two processes which have genuinely different
is determined by a master equation for the probabilities
dynamics and stationary states, like TASEP with open
to find a particular configuration {n }. We apply the
i
boundaries and LK.
following dynamical rules. For each time step, a site i
In this letter, we relax the constraint that the con- is chosen at random. A particle at site i can jump to
served dynamics in the bulk imposes to the TASEP by site i+1 if unoccupied (we fix units of time by putting
2
this rate equal to unity). In the bulk i=2,..,N −1, a
particle can also leave the lattice with site-independent 0 0.2 0.4 0.6 0.8 1
detachmentrateω orfillthe site (if empty)with arate
D (a)
ω by attachment. At the boundaries, a particle can fill
A
0.8 W D =10
a vacancywith arate α atsite i=1,ora vacancycanbe
formedbyremovingaparticlefromthelatticewitharate
β atsitei=N. Werefrainfromgivingexplicitlythemas- 0.6 W D =0.1
ter equation for the probabilities. Correlations induced ni r(x)
intothemany-particleproblemcanbeconvenientlystud- 0.4
ied within an operator representation in Fock space [9].
W D =0.001
Then the equations of the bulk dynamics read: 0.2
dn
i
=ni−1(1−ni)−ni(1−ni+1)+ωA(1−ni)−ωDni, 0
dt (b)
(1a)
while at the boundaries one obtains: 0.8
5
dn1/dt=α(1−n1)−n1(1−n2), (1b) 0.6 43
dnN/dt=nN−1(1−nN)−βnN. ni k=2
0.4
By taking averages [10] one observes that in order to
compute the time evolution of hn (t)i one needs the cor-
i
0.2
responding averages of higher order correlations. In or-
der to obtainan exact solution, elaborate techniques are
necessary. WerestrictthediscussiontoMonte-Carlosim- 0
(c)
ulations (MCS) and a mean-field approximation (MFA)
which we shall apply below. 00..88 W D =10
0.75
The system exhibits a particle-hole symmetry in the
following sense. A jump of a particle to the right corre- 00..66
sponds to a vacancy move by one step to the left. Sim- r(x) 1.0 0.2 0.1 0.08
0.06
ilarly, a particle entering the system at the left bound- 00..44 0.4
ary can be interpreted as a vacancy leaving the lattice, 0.04
and vice versa for the right boundary. Attachment and
00..22
detachmentofparticlesinthe bulk ismappedto detach- W D =0.001
ment and attachment of vacancies, respectively.
00
We are interested in large system sizes (N ≫ 1) and, 00 00..22 00..44 00..66 00..88 11
x
eventually, in the “thermodynamic limit” N → ∞. In
FIG. 2: (a) Average density profile hnii computed by MCS
this case,the studyofthe competitionbetweenbulk and
(continuous line) and average density profile ρ(x) computed
boundarydynamicsneedsthatthe kineticratesdecrease
by numerical integration of MFA stationary state equations
simultaneously with the system size. More precisely, we (2) (dashed line) in therescaled variable x=i/N for N=103
define the“reduced”ratesΩA andΩD asΩA=ωAN and with α=0.2, β=0.6, K=3 and different kinetics rates ΩD
Ω =ω N, keeping Ω ,Ω ,α,β fixed as N→∞. Note indicatedinthegraph.(b)MCSaveragedensityprofilefordif-
D D A D
thatthebindingconstantK=ω /ω remainsunchanged ferentsystemsizes,sameα,β,K asbefore,andΩD=0.1. The
A D
width of thesteep rise decreases with increasing system sizes
whenpassingtothethermodynamiclimit. Moreover,for
N=10kwithk=2,3,4,5indicatedinthegraph. (c)MFAav-
ωA=ωD=0, one arrives back at the TASEP respecting eragedensityprofilefor ε=10−3,sameα,β,K asbefore,and
the same particle-hole symmetry described above.
different kinetic rates ΩD indicated in the graph. The hor-
We have performed extensive computer simula- izontal dashed line for ρ(x)=0.75 represents the Langmuir
tions [17] to obtain the average density profile in the isotherm for K=3.
stationarystate. Weillustratetypicalphenomenabyfol-
lowingapathinparameterspacealongcurveswithfixed
by the well-knownratio K/(1+K) of Langmuir equilib-
α,β, and K while increasing Ω = Ω /K. Fig. 2(a)
D A
riumdensity[18]. Anewfeatureappearsforintermediate
shows the density profile for three different values of the
rates, Ω =0.1, precisely when bulk and boundary dy-
kinetic rates. At small kinetic rates, Ω ,Ω ≪ α,β, D
A D
namics compete. The density exhibits a non-monotonic
the average density hn i in the bulk is practically con-
i
structure in bulk, characterized by a region of low and
stantand closeto the low-densityvalue predicted by the
high density connected by a steep rise.
TASEP, hn i=α. Conversely, at high kinetic rates the
i
bulk profile is structureless and essentially determined Fig. 2(b) shows the density profile for different sys-
3
tem sizes. One observes a decrease of the width of the findgoodagreementofMFAcomparedwithMCSforthe
transitionregionasthenumberofsitesisincreased. The full range of kinetic rates in the limit of large N.
simulationsuggestsadiscontinuityoftheprofileinterms In analogy with fluid dynamics, to describe these re-
of the rescaled variable x=i/N upon approaching the sults one considers an effective current density which for
infinite system limit. A preliminary finite size scaling ourproblemreadsj=−(ε/2)∂ ρ+ρ(1−ρ). Abbreviating
x
analysis is compatible with a rescaled transition width thefluxesfromandtothereservoirbyF =Ω (1−ρ)and
A A
which scales as N−ν with ν ≃ 0.5. This is very differ- F =Ω ρ, Eq. (2) can be read as a balance equation:
D D
ent from the mean-field result, νMFA=1 [19]. Thus we ∂xj=FA−FD. Since there aretwoboundaryconditions
have identified an intermediate parameter range where one has to be carefulwhen discardingthe secondderiva-
lowandhighdensityphasescoexistseparatedbyasharp tive in (2) for a small prefactor ε. The average profile is
domain wall (DW). We also find that the discontinuity then governed by similar physics as the Burgers’ equa-
inthedensityseemstobestableoratleastlocalizedina tion in the inviscid limit [11]. Generically one expects
small region compared to the system size [20]. This has shocks (here, “domainwalls”,DW) in the bulk and den-
to be contrasted with the domain wall (“shock”) found sitylayersattheboundaries(“boundarylayers”). Cross-
in the TASEP right at the phase boundary between the ing a DW, the current j remains continuous in the limit
high and low density phase (α=β <1/2) which is de- ε→0, while boundary layers form whenever the density
localized and moves as a random walker once it is far associate to the bulk current does not fit the boundary
from the system boundaries [9]. condition. To better understand these features, we have
Phenomena like phase separation/coexistence have exploredthedependenceofthedensityprofileρ(x)onΩ
D
previously been observed in non-homogeneous systems for fixed α,β and K (see Fig. 2(c)). For small kinetic
with open boundaries like TASEP with isolated local- rates, Ω =10−3, the profile is close to the one expected
D
ized defects [11, 12]. The location of the domain walls from TASEP, with a boundary layer bridging the bulk
areexpectedandfoundtobeidenticaltothedefectposi- density up to ρ=1−β (not resolved in Fig. 2(c)). In-
tions. Incontrast,thelocationoftheDWinourhomoge- creasing Ω the slope of the bulk density increases. For
D
nousmodelisself-tunedanddeterminedbythe valuesof Ω >0.05, MFA connects a region of low density (LD),
D
the kinetic rates (see below). In systems with periodic i.e. ρ(x)<1/2 to a high density region(HD), ρ(x)>1/2,
boundary conditions (which are not the subject of this by a DW. Whereas the solution close to the left bound-
letter)phaseseparationhasbeenfoundinTASEPwitha ary is smooth, one finds a boundary layer at the right
blockage[13], quenched disorder[14], orin homogeneous end bridging densities ρ=1/2 down to ρ=1−β. For
systemswithmulti-speciesparticledynamics[15](seefor larger Ω the DW moves to the left, while the slope of
D
a general criterion [16]). the LD region increases and the HD profile flattens ap-
To rationalize all these findings we have developed a proaching the Langmuir density value K/(K +1). For
mean-field theory. Defining ρ = hn i, the MFA con- Ω >1 the DW remains practically localized at the left
i i D
sistsoftakingtheaverageofEqs.(1a,1b)andfactorizing boundary. Note that the DW location strongly depends
the two-site correlations, hnini+1i=ρiρi+1. Then Eqs. on typical values of the bulk kinetic rates when they are
(1a,1b) display the same form provided that the binary comparable with the boundary rates α and β.
occupationnumberni isreplacedbythecontinuousvari- In the inviscidlimit ε→0 the complete phase diagram
ableρiwith0≤ρi≤1. Theequationsarenowinterpreted can be obtained analytically within MFA, up to some
as ordinary differential equations. treatment of the density discontinuities. Interestingly,
Toobtainananalyticallytractablesystemofequations the solution found is never given by either a constant
we have coarse-grained the discrete lattice with lattice low/high density profiles as in TASEP or the Langmuir
constant ε=L/N to a continuum. For fixed total length isotherm, but by a completely new set of solutions [22].
L=1andN→∞,ε→0onegetsthenonlineardifferential The mean-field analytical solution allows to draw the
equation for the average profile in the stationary state, phase diagram and compare it to TASEP. Fig. 3 repre-
sents a cut through the phase diagram for Ω =0.1 and
ε D
2∂x2ρ+(2ρ−1)∂xρ+ΩA(1−ρ)−ΩDρ=0, (2) K =3 with α and β used as control parameters. One
finds an extended LD-HD coexistence region separating
where positions are measured by the rescaled variable a LD and a HD phase. At the boundaries of the coexis-
x = i/N, 0 ≤ x ≤ 1. Equations (1b) translate now to tence region, the DW between the low and high density
boundary conditions for the density field, ρ(0)=α and phases are located in the proximity of the open ends of
ρ(1)=1−β.OneobservesthatMFArespectstheparticle- the1Dlattice. Forsmallαthe DWdevelopsatthe right
hole symmetry mentioned above, provided that when end, x=1 and moves to the left as α increases. At the
ρ(x) 7→ 1−ρ(1−x) one interchanges α↔β, Ω ↔Ω . phase boundary between the coexistence region and the
A D
Duetothispropertywecanrestrictthediscussiontothe HDphase,theDWislocatedattheleftendofthelattice,
case Ω >Ω [21]. The numerical mean-field solutions x=0. In both cases, when the DW enters and leaves the
A D
are included in Fig. 2(a) for different values of Ω . We lattice, it matches with a boundary layer connecting the
D
4
ht
1 g
b ei b=0.2
h
0.8 0all .4 b=0.4
w
LD 0.6 D 0main .2 b->0.5 [1] Ban.dSCchrmitiictatmlPahnennoamndenRa,.KV.oPl..1Z7i,a,Edins.PCh.aDseomTbraannsditiJo.nLs.
H o Lebowitz (Academic Press, London, 1995).
=1 D- d 00 0.1 0.2 0.3 [2] C.T. MacDonald, J.H. Gibbs and A.C. Pipkin, Biopoly-
0.4 x L mers 6, 1 (1968).
a [3] J. Krug, Phys.Rev.Lett. 67, 1882 (1991).
0
= [4] B. Derrida, E. Domany and D. Mukamel, J. Stat. Phys.
0.2 x HD
69, 667 (1992); G.M. Schu¨tz and E. Domany, J. Stat.
Phys. 72, 277 (1993); B. Derrida et al., J. Phys. A 26,
0 1493 (1993).
0 0.2 0.4 0.6 0.8 1
[5] R.H. Fowler, Statistical Mechanics, (Cambridge, 1936).
a
FIG. 3: Phase diagrams obtained by the exact solution of [6] A.Albertsetal.,TheMolecularBiologyoftheCell,(Gar-
the stationary mean-field equation (2) in the inviscid limit land, New York,1994).
for K=3 and ΩD=0.1. The inset shows the dependence of [7] Our model could also be relevant for studies of surface-
theDW amplitude on α for different values of β. adsorption and growth in presence of biased diffusion or
totrafficmodelswithbulkon-offramps.Seee.g.V.Priv-
man (Ed.) Nonequilibrium Statistical Mechanics in One
bulk density with the density given by the correspond-
Dimension (Cambridge Univ. Press, Cambridge, 1997)
ing boundary condition. (In this case, boundary layers
and D.Chowdury et al.,Phys. Rep.329, 199 (2000).
appear only for α>1/2 and/or β >1/2). MFA shows [8] K. Kruse and K. Sekimoto, Phys. Rev. E 66, 031904
that the bulk density profile becomes independent of β (2002).
for β > 1/2. Hence, in this regime for a given α only [9] Forareviewseee.g.: G.M. Schu¨tz,inPhase Transitions
the magnitude of the boundary layers changes, but not and Critical Phenomena,vol.19,Eds.C.DombandJ.L.
the profile of the bulk density. This explains the vertical Lebowitz (Academic Press, San Diego, 2001).
[10] The notion of averaging is non-trivalsince the system is
phase boundaries of the coexistence region for β≥1/2.
far from equilibrium. Nevertheless, one expects the ini-
In addition to its location the DW is also character-
tial state to relax to a stationary ensemble over which
ized by its height (see Fig. 3). For β<1/2 we find that
averages should betaken.
the height discontinuously jumps to a finite value upon [11] S.A.Janowsky andJ.L. Lebowitz, Phys.Rev.A 45, 618
entering the coexistence regionfrom the LD phase. This (1992).
has to be contrasted with the case β ≥ 1/2 where the [12] A.B. Kolomeisky, J. Phys. A 31, 1153 (1998).
[13] G.M. Schu¨tz, J. Stat.Phys. 71, 471 (1993).
phase transition fromthe LD to the coexistence phase is
[14] G. Tripathy and M. Barma, Phys. Rev. E 58, 1911
characterized by a continuous increase in the height of
(1998).
the DW (comparethe dashedline inFig.3). InMFA we
[15] M.R. Evanset al., Phys.Rev.E 58, 2764 (1998).
findthatbothDWamplitudeandpositionexhibitpower
[16] Y. Kafri et al.,Phys. Rev.Lett.89, 035702-1 (2002).
law behavior with (α−αc)1/3 and (α−αc)2/3, respec- [17] MCS were performed with random sequential updating
tively. At the phase boundary to the HD phase the DW in discrete time: at each time step, a site is chosen at
height always jumps to zero discontinuously. random and the dynamical rules are applied. Time and
WorkingoutthecompleterichscenariofordifferentΩ sample averages havebeenevaluated.Theresultingpro-
A
files coincide in both averaging procedures for given pa-
and Ω in the limit N ≫1 needs a detailed analysis of
D rameters and different system sizes. In simulations for
thephasediagram[22]. Wejustmentionthattheoriginal
Figs.2(a)and2(b)stationaryprofileshavebeenobtained
maximal current phase of the TASEP appears for ΩA= over 105 time averages (with a typical time interval of
ΩD only, where the Langmuir density is valued to 1/2. ∼10 N between each average step).
This can be proved numerically as well as analytically. [18] Inthiscase,duringthedetachmentandattachmentpro-
Conversely,assoonasΩ 6=Ω , MFApredictsthatsuch cesses, we consider a constant particle concentration in
A D
a phase continuously disappears in favor of a HD phase the reservoir cr (i.e. cr=1).
[19] Upon introducing x′=x/ε, a scaling analysis of Eq. (2)
if Ω >Ω (respectively LD phase if Ω <Ω ).
A D A D showsthatthedomainwallwidthscalesasνMFA=1[22].
The authors thank P. Benetatos and J.E. Santos for
[20] As a simple test, we artificially shifted the front of the
useful discussions and comments. This work was par-
discontinuityand let thesystem evolvein time.Thedis-
tiallyfundedbytheDFGundercontractnos. FR850/4-1 continuity always moved back to its initial position.
and SFB 413. A.P. is supported by a Marie-Curie Fel- [21] The case ΩA=ΩD introduces new symmetries [22].
lowship no. HPMF-CT-2002-01529. [22] A. Parmeggiani, T. Franosch and E. Frey,unpublished.