Bidirectional Non-Markovian Exclusion Processes
aa r X i v : . [ n li n . C G ] J a n
20 January 2020
Bidirectional Non-Markovian Exclusion Processes
Robin Jose, Chikashi Arita and Ludger Santen
Department of Theoretical Physics & Center for Biophysics, Saarland University,66123 Saarbr¨ucken, GermanyE-mail: [email protected]
Abstract.
Bidirectional transport in (quasi) one-dimensional systems genericallyleads to cluster-formation and small particle currents. This kind of transport canbe described by the asymmetric simple exclusion process (ASEP) with two species ofparticles. In this work, we consider the effect of non-Markovian site exchange timesbetween particles. Different realizations of the exchange process can be considered:The exchange times can be assigned to the lattice bonds or each particle. In the lattercase we specify additionally which of the two exchange times is executed, the earlierone (minimum rule) or the later one (maximum rule). In a combined numerical andanalytical approach we find evidence that we recover the same asymptotic behavior asfor unidirectional transport for most realizations of the exchange process. Differencesin the asymptotic behavior of the system have been found for the minimum rule whichis more efficient for fast decaying exchange time distributions.
1. Introduction
One of the broadly investigated fields in non-equilibrium physics is actively driventransport. These processes can be found in different topics such as pedestrian dynamics[1, 2, 3, 4], vehicle traffic [5, 6, 7, 3], and intracellular transport of molecular motorsalong cellular filaments [8, 9, 10, 11, 12, 13, 14].A common tool to model active transport is a lattice gas [15, 5]. These stochasticprocesses are defined in a very simple way, but lead to many interesting phenomena [16].Particles hop stochastically to their nearest neighbor sites on the lattice, but the hoppingrates are spatially biased, and this asymmetry causes a non-vanishing flow of particles ina specific direction. The particular case where particles are allowed to unidirectionallyhop in one dimension is called the totally asymmetric exclusion process (TASEP) [17].One of the most basic assumptions in the TASEP is that each site of the lattice is eitheroccupied by a particle or empty. Due to this exclusion principle, hopping is prohibitedif the target site is already occupied by another particle. Therefore particles behave asan obstacle for each other, in other words, particles themselves serve as an environmentand influence motility.In bidirectional transport, however, particles have to be transported in oppositedirections. Adaptation of the TASEP can be done by distinguishing two different speciesof particles with opposing directions on the same lattice. By introducing an positionexchange rate, several intriguing scenarios have been reported, such as spontaneoussymmetry breaking [18, 19]. In case of slow position exchange interactions, the particleflux is determined by the exchange times of particles from different species similar toa defect in the unidirectional TASEP [20]. Bidirectional TASEP models have beenmodified by introducing a second dimension to describe pedestrian dynamics [21, 22]or intracellular traffic on polar filaments by adding additional lanes [23, 24, 9]. Despitesituations where symmetry is broken and the system organizes into lanes [22], particleinteractions often lead to cluster formation and is therefore a limiting factor for transport[25]. In this paper, the focus is on active particles in confinement. Here, not only theaspect of non-equilibrium drive but also crowded and confined environments is expectedto heavily influence transport processes. In the field of glass theory for example, apopular approach is to describe a particle trapped inside a cage , denoting the potentialcreated by its neighbor particles [26]. Such interactions can affect the waiting timedistribution of particle movements as it can be seen in a trap model by Bouchaud etal. [27]. Here, a particle falls into a trap of potential depth E which is exponentiallydistributed and escape from it following a Poisson process with a rate depending on theenergy E . This combination leads to algebraically (power law) distributed waitingtimes for particles to escape from traps. It is therefore not guaranteed that in acomplex environment, properties from an exponential distribution, resulting in constantrates which are typically used in TASEP models. However, the choice of waiting timedistribution can be crucial for the systems phenomena because heavy tails induce ahigher statistical weight for extreme values such as for the scale free family of Paretodistributions [28].Recent studies by Concannon et al. [29] and Khoromskaia et al. [30] took a stepforward to investigate transport behavior in the framework of unidirectional exclusionprocesses for non-Markovian waiting time distributions. It was found that, beside a fluidphase, the particles form condensates which are complete in space and time and hence aflux depending on the system size. This phenomenon differs from typical condensationsappearing beside a stable current flowing out of clusters which is seen in models relatedto Markovian processes [25, 31].The influence of crowded environments reflected in non-Poissonian waiting timeson bidirectional transport is however not fully understood. Combining the aspects ofbidirectional lattice gases and non-Poissonian waiting time distributions for exchangeprocesses, here we will investigate two-species non-Markovian TASEPs.This work is organized as follows. In section 2, we develop the model including threesub-versions for realizing exchange process between particles. We then show analyticalestimations and simulation results for single and many particle dynamics in section 3and compare them no simulation results. Finally we discuss our results in section 4.
2. The model
In order to mimic bidirectional transport on a track we first introduce plus- and minus-directed particles on a one-dimensional lattice of L sites, each of which is either occupiedby a plus (“+”) particle, occupied by a minus (“ − ”) particle, or empty. Empty sides aredenoted as holes if the particle density ρ <
1. Particles are identical up to their direction.We have three types of stochastic, microscopic events between two neighboring sites,i.e. +0 ⇒ , (1)0 − ⇒ − , (2)+ − ⇒ − + . (3)This two-species TASEP has two conserved quantities, i.e. the numbers N ± of plus andminus particles, under periodic boundary conditions.In a crowded environment, the exchange process (3) can be very different fromfree hopping events and turn out to be the major criterium for estimating the particleflux similar to bottlenecks in unidirectional exclusion processes [20]. In a first approachwe therefore impose Poissonian stochasticity on the normal jumps (1) and (2), but aheavy tailed non-Markovian waiting time distribution on exchange processes (3) in thissection. It has been observed that the particle flux is heavily influenced by particlecondensates, both in a two-species Markovian TASEP [18] and in a one-species non-Markovian TASEP [29]. We see a similar phenomenon in first results of our two-speciesnon-Markovian model, (figure 1 (a)). Macroscopic clusters block particle flow for a longtime interval interrupted by short boosts of particles hopping in a free space outsidethe clusters. This blockage is the major inhibitor of particle flux so that the waitingtime distribution for exchange processes mainly controls the transport property of oursystem. We therefore focus on the cluster region indicated by the green ellipse in figure1 (a) leaving only the two states “+” and “ − ” for a site in the lattice in the following. By concentrating on clusters and therefore neglecting holes, we consider a system fullyoccupied by plus and minus particles without holes N + + N − = L , for which an exemplaryconfiguration is given in figure 1 (b).We remark that usually the two-species TASEP with N + + N − = L is consideredas standard TASEP below. We will explain the difference to the standard one-speciesTASEP in detail.We first consider the standard one-species TASEP, where each site i takes states η i = 1 (occupied by a particle) or η i = 0 (empty). In most of the cases of the TASEP, theexponential distribution p ( t ) = α − e − αt ( t >
0) is used to generate waiting times betweentwo consecutive stochastic events, hence the TASEP is usually a Markov process. The space timetime space a) b)c) pb sb t i t i t i+1 t i+1 t i+1 t i t i+2 t i+2 t i+2 t i+3 t i+3 t i t i t i t i+1 t i+1 t i+1 t i+3 Figure 1: (a)
Kymograph section of a simulation of our two-species TASEP with holes,color-coded as plus particles in blue, minus particles in red and holes in white. Hoppingtimes (1) and (2) are distributed exponentially, exchange times (3) by a power law. Thesystem size is L = 1000 with N + = N − = 200, the shown time interval is approximately15000 time units. In this paper, we focus on the cluster, indicated by a green ellipsein the plot, leading to the model without holes defined in chapter 2.1. (b) Schemeof the two-species model without holes. Triangles mimic particles, the tip indicatesthe direction. Exchange processes are only allowed for the configuration (+ − ) (3). (c) Comparison of the two update schemes particle-based (left) and site-based (right)regarding particle trajectories (red lines). Timelines (dashed lines) for hopping events(blue bars) are carried by particles, hence identical to the red trajectory line in the leftscheme or are fixed to the sites (right). Hopping events can be rejected due to exclusion(red crosses).parameter α stands for the hopping rate and is independent of time or the currentsystem state. A common way to assign waiting times to the Markovian TASEP is touse Gillespies direct method or first reaction method [32]. In the latter one, times aregenerated from the exponential distribution for every possible hopping transition butonly the smallest one determines the process which is executed and the other timesare not used in the update mechanism anymore. In order to increase computationalefficiency, a modified waiting time algorithm called next reaction method [33] is storingthe assigned times for every process and realizing them successively if the transition isallowed.Using the next reaction method, we distinguish two ways to assign waiting times tothe Markovian TASEP, illustrated in figure 1 (c). The first one is a site-based update.We generate and store a series of times { t i , t i , . . . } for each bond between sites i and i + 1. The difference t k +1 i − t ki with t i = 0 must obey the exponential distribution. Attime t = t ki , one checks whether the local configuration is appropriate for hopping, i.e. η i = 1, η i +1 = 0. If this is the case, we move the particle from i to i + 1. If not, themove is rejected. The second one is a particle-based update. We give a series of times { t i , t i , . . . } to each particle labeled by i . Again the difference t k +1 i − t ki should followan exponential distribution. The particle i attempts to hop to its right neighbor site attime t = t ki , but the jump is again allowed only when the target site is empty. Thesetwo update schemes eventually yield equivalent dynamics to the particles, as long as weuse an exponential distribution.However, for a probability density function (PDF) p ( t ) which belongs to a powerlaw p ( t ) = ( < t < γ − t − γ t >
1, (4)with an exponent γ > { t i , t i , . . . } for each bond, such that now thedifference t k +1 i − t ki obeys the algebraic distribution (4). If we find η i = + and η i +1 = − at time t = t ki , we simply exchange the positions, i.e. we get a new configuration with η i = − and η i +1 = +.On the other hand, the particle-based update becomes a little complicated, and weneed to further divide it into sub-schemes. The simplest one is considered as follows. Weonly assign the time series to the + particles. At time t = t ki , the + particle labeled by i wants to hop rightward. This is allowed only if there is a − particle on the target site. Inother words, − particles are regarded as holes, and the systems is completely identicalto the one-species non-Markovian TASEP that was introduced by Concannon et al. Thus, we name this rule particle-based- asymmetrical update, because this case does nothold plus-minus (or particle-hole after the identification) symmetry for Non-Markovianprocesses.Now we wish to look for rules which do not break the plus-minus symmetry todefine two-species bidirectional transport using the particle-based approach. We assigna time-series to each plus and minus particle { t j , t j , . . . } . At time t = t ij the particle isactivated if a “+ − ” configuration is given. We now consider two different local updateschemes for a plus particle i with waiting time t ki and a minus particle with waiting time t mj which are located at neighboring sites l (plus particle) and l + 1 (minus particle).1) The minimum rule : An exchange between a pair of “+ − ” particles is executedif one of the two particles is active. This means that the minimum of the two waitingtimes determines the particle exchange.2) The maximum rule : An exchange of the particles is executed if both particles areactive, i.e. a particle exchange happens after max { t ki , t mj } . In both cases, the particlesbecome passive after exchanging positions.Let us summarize the four types of the update rules in the following: site-basedparticle-based asymmetricsymmetric ( maximumminimum (5)
3. Results
In this section, we show simulation results of our non-Markovian two-species TASEP,with a completely filled system, i.e. N + + N − = L . We first discuss the site-basedupdate from [30] as an algorithm for bidirectional transport in section 3.1 and then thepreviously described three types of particle-based rules in section 3.2. For each of themwe evaluate dynamics by measuring the PDFs for the effective waiting time, i.e. theduration between two realized exchange processes of a particle. We discuss the effectivewaiting time distribution, in the following two situations; the system with only oneplus particle and L − single particle dynamics (spd), and theequally-filled case where N + = N − called many particle dynamics (mpd). We discussthe spd as a reference in order to highlight the collective effects modifying the originalPDF p ( t ) to the PDF for effective waiting times of exchange processes. With manyparticle dynamics, we investigate the influence of exclusion on the effective PDF. Incase of mpd, we also investigate transport efficiency by measuring particle flows. Wecompare the flows for the three symmetrical model rules in section 3.3. The first bidirectional update we deal with is the site-based model inspired by [30]. Westart with spd, i.e. the situation where there is only one plus particle, i.e. N + = 1and N − = L −
1, and we probe its motion in the environment of minus particles.We give analytical estimates for and measure by simulations the effective waiting timedistribution for site-exchange events, and compare the tail exponents of this quantity.For the site-based update scheme, analytical estimates can be deduced from renewaltheory as in [30]. Let us assume that the plus particle arrives at a site i , where an internalclock is already running since the last event on that site at time t ki , giving the clock an a) sb sb (cid:1) =2.5 (cid:0) =2.5 b) spd mpd tt P ( t ) -8 -6 -4 -2 -1 eff. time p(t)r(t) 10 -8 -6 -4 -2 -1 eff. timep(t)r(t) Figure 2:
Site-based rule:
PDFs of the effective exchange time for a system of L = 100for the site-based rule. (a) Single-particle dynamics, i.e. N + = 1 and N − = L −
1. (b)Many-particle dynamics, i.e. N + = N − = L/ t − t ki . The particle has to wait for the remaining time until t k +1 i to execute thefollowing step. This remaining time is called the residual waiting time distributed bythe PDF r ( t ) which in general is different from the original PDF p ( t ). In the case where p ( t ) is given by equation (4), one finds the residual waiting time PDF in the form r ( t ) = ( γ − γ − < t < γ − γ − t − γ t ≥ r ( t ) > < t <
1, while theoriginal p ( t ) is zero. More remarkably, the residual waiting time has a tail exponentshifted by 1, meaning that very large values have a higher statistical weight than theoriginal PDF p ( t ). In figure 2 (a), our simulation results for spd agree to the residualwaiting time r ( t ).As a next step, we introduce N + = L/ r ( t ) calculatedfor a single particle. This similar result is expected since the waiting times are renewedin both cases for every site in the lattice, no matter if a plus particle or a minus particleis occupying the site.After discussing the influence of the used update method and the particle interactionon the PDF of exchange times, we will continue by examining transport properties in abidirectional, two-species system. As discussed of effective exchange time distributions,we will continue by examining how these effects are reflected in the transport propertiesin the bidirectional two-species TASEP with N + = N − = L/ J . In [29], it was shown that for the unidirectional non-Markovian TASEP, systemsize has an impact on J in the clustering phase. Since we find similar shifts of the PDFfor the site-based update, we also expect a transition from a size-dependent to a size- sb -4 -3 -2 -1 J L γ = γ (cid:2) γ (cid:3) γ (cid:4) γ (cid:5) (cid:6) .0 λ (cid:7) (cid:8) Figure 3:
Site-based rule:
Length dependency of the particle flux J for differentvalues of γ in the PDF (4) and mpd. The black line serves as a comparison to the fluxgenerated by an exponential distributed waiting times with exponent λ .independent regime with growing γ . We test this transition by plotting the particle fluxversus the length of the system in figure 3 for different exponents γ . Up to a finite-sizeeffect, no significant dependency of J is observed when using a γ > L , so that no significant size-dependency is observed for L > J is decreasing with the system size for γ <
3. At the same time, theerror bars are larger and it is difficult to judge the limit for infinite system size fromnumerical results. It is, therefore, necessary to argue with additional information abouta transition from length dependency to constant fluxes. From the residual waiting time,we know that the average effective waiting time diverges at γ = 3. So we expect the fluxto vanish below this critical value γ c = 3 for infinite systems similar as pointed out in[30]. Below γ c , the transport efficiency is determined by the asymptotics of the effectivewaiting time distribution. This effect induces a strong size dependence of the flow and γ =2.5 γ =2.7 γ =3.0 γ =4.0 γ =5.0a) b) -8 -6 -4 -2 -1 P ( t ) t L= 5L= 10L= 40L= 100p(t)r(t)
Figure 4: (a) Kymographes for the site-based model with runtime: 8 · system size L = 1000 and equal particle numbers N + = N − . For γ < γ = 2 . L = 1000 which cover a time period of 8 · time units for differentexponents in figure 4 (a). Homogeneous spatiotemporal structures are obtained for γ ≥
4, while stable clusters emerge for γ ≤
3. In this regime, the lifetime of the clustersare comparable to the simulation time, which make it difficult to reach the stationarystate of the system numerically. The slow relaxation of the system is caused by rareextreme values for the waiting times which develop in an aging phase inside a clusterand block other exchange events such as in [29]. We estimate the uncertainty of J byusing the partial time averages J n = (cid:16)P t n +∆ t = t n J t (cid:17) / (cid:16)P t n +∆ t = t n t (cid:17) where ∆ is the completetime of measurement divided by the number of partial blocks and t n = n ∆.For the site-based model, we expect that the residual waiting time PDF (6)describes exactly the effective exchange time distribution. The results of our simulations,however, slightly differ from this prediction. We expect that the small deviation can beattributed to a finite-size effect, since local waiting times may exceed the typical timea particle needs for a complete tour in a finite periodic system. This finite-size effect isobservable in figure 4 (b) where tails differ from r ( t ) for small system sizes. As a reference system for particle-based updates, we startwith the asymmetrical rule, which is identical to the model of Concannon et al. [29]. Inthe paper, it was pointed out that interaction via exclusion leads to a shift of one forexponents in the hopping time PDF for a unidirectional many particle system. Again,we will show PDFs for single-particle and many-particle dynamics.The spd under the particle-based, asymmetric update rule is much easier than for a(cid:9) pb asym pb asym (cid:10) =2.5 (cid:11) =2.5 b) spd mpd tt P ( t ) -8 -6 -4 -2 -1 eff. timep(t)r(t) 10 -8 -6 -4 -2 -1 eff. timep(t)r(t)1 Figure 5:
Particle-based asymmetrical rule:
PDFs of the effective exchangetime for a system of L = 100 for the asymmetrical particle-based rule. (a) Single-particle dynamics, i.e. N + = 1 and N − = L −
1. (b) Many-particle dynamics, i.e. N + = N − = L/ { t , t , . . . } contains the times for which the plusparticle moves, since it is never blocked by other plus particles. In other words, theeffective waiting time distribution is p ( t ) itself, and we can regard the plus particlesimply as a non-Markovian random walker while minus particles serve as passive holes.The motion of the particle is completely determined by { t , t , . . . } . Our simulationresults for this scenario are shown in figure 5.Now turning to the mpd realization of the asymmetric update, we get back thescenario discussed in [29], since minus particles are passive and correspond to holes ofthe one-species TASEP. For completeness, the results are shown in figure 5 (b). Asexpected, we observe that for small times near t = 1, the simulation data points followthe original density p ( t ) but then the exponent of PDF changes to γ − γ ≤ γ c has already been discussed in [29], hence we continue with the symmetric rulesfor particle-based updates. Let us turn to the first particle-based-symmetrical update, themaximum rule, which was introduced in section 2. This rule does not break thesymmetry between plus and minus particles because the exchange process is triggeredby both particles (with index i and j ) that have to be activated for an exchange processfirst. Hence, one has to choose the maximum from the two next event times in thetime-series t ki , t ℓj , in order to determine when the plus particle actually exchanges withits neighbor. We estimate PDFs for the symmetrical maximum update by calculatingthe density of the maximum of two random variables X and Y with density p ( t ) and p ( t ).Because all particles follow to their own time series { t i , t i , . . . } and { t j , t j , . . . } ,different situations appear for neighboring pairs of plus and minus particles. There isalways one particle inducing the exchange process, e.g. the one with the later time forthe maximum update rule. However, the process can be induced by a plus particle or aminus particle.Let us assume that a plus particle induced an exchange k and then becominginvolved in a new exchange process k + 1 with a new partner. Here, the plus particle’sclock has no age immediately after the last executed exchange process. In contrast, thenew neighbor minus particle already is located on its position for some time meaning aclock with a non zero age. Hence, the plus particle follows the density p plus ( t ) = p ( t )but we assume that the minus particle’s residual waiting time is rather described by p minus ( t ) = r ( t ) because it is standing in the queue of other minus particles. In thesingle particle case, the maximum of two random times is taken from a time X thatfollows p ( t ) and a time Y that follows r ( t ) as max ( X p ( t ) , Y r ( t ) ) which we call a mixedscenario in the following.However, for symmetric rules, minus particles also can introduce an exchangeprocess. Let us follow a plus particle again, but this time the exchange k was inducedby its former neighboring particle. The plus particle has a new exchange partner for the1process k + 1 again, but this time also a non zero age, just as its new neighbor minusparticle which was not involved in the former process k . In this case, we assume bothresidual waiting times are distributed by p plus ( t ) = p minus ( t ) = r ( t ). The maximum oftwo random times X and Y is now chosen as max ( X r ( t ) , Y r ( t ) ), called the pure scenario .First, we give the cumulative distribution function (CFD) for both p ( t ) and r ( t ): P ( t ) = ( < t < − t − γ t ≥
1, (7) R ( t ) = ( γ − γ − t < t < − γ − t − γ t ≥
1. (8)With these, we can now calculate the density of the maximum X and Y , i.e. f mixmax ( t ) = ddt [ P ( t ) R ( t )] (9)= ( < t < γ − t − γ + γ − γ − t − γ + − γγ − t − γ ) t ≥
1, (10) f puremax ( t ) = ddt [ R ( t ) R ( t )] (11)= (cid:16) γ − γ − (cid:17) t < t < γ − γ − t − γ − γ − γ − t − γ t ≥
1. (12)Both, the pure and mixed scenarios have the leading exponent 1 − γ which is equal tothe exponent of r ( t ). In the pure scenario, we get an estimate for effective waiting timessmaller than 1. This expression is expected to overestimate the weight of the waitingtimes because also the mixed scenario is contributing to the exchange processes, alwayswith times larger than one. However, we will use the pure scenario as an estimate for a) pb max (cid:12) =2.5 b) spd -8 -6 -4 -2 -1 P ( t ) t eff. timer(t) fm(cid:13)x (t) pb max (cid:14) =2.5mpd - (cid:15) (cid:16) (cid:17)2 (cid:18) (cid:19) t e ff . (cid:20)(cid:21)(cid:22)(cid:23) r(t) (cid:24)(cid:25)(cid:26)(cid:27) (t) Figure 6:
Particle-based maximum rule:
PDFs of the effective exchange time for asystem of L = 100 for the particle-based maximum rule. (a) Single-particle dynamics,i.e. N + = 1 and N − = L −
1. (b) Many-particle dynamics, i.e. N + = N − = L/ f max ( t ) = f puremax ( t ) = (cid:16) γ − γ − (cid:17) t < t < γ − γ − t − γ − γ − γ − t − γ t ≥
1. (13)We compare these estimates with simulation results in figure 6 (a) for spd. Here,the tail behavior of our simulation results are in good agreement with f max ( t ). Also,short time behavior is well approximated when relating the sharp increase at t ≈ t < t > γ = 2 . f max ( t ) to simulation results of the mpd in figure 6 (b). Again wesee that the tail behavior is well described by f max ( t ) and r ( t ). The short time behavioris still close to the estimate but the tail has a higher statistical weight comparing tospd. The shift in the exponent that is seen when comparing the original waiting timePDF and the effective exchange time PDF is also consistent with the results for the fluxin the maximum rule which is shown in figure 8 (b). We observe no significant changesin the flux for system sizes larger than L = 100 if γ > γ c but a flux that vanishes with L for γ < γ c , similar to the results found for the site-based model. The second particle-based-symmetrical update is the minimumrule, also introduced in section 2. Also, this rule does not break the symmetry betweenplus and minus particles. The exchange process is triggered by the first particle that isactivated for an exchange process at the minimum time of the two next event times ineach particle’s time-series t ki , t ℓj . We calculate an estimates for PDFs from the minimumof two random variables X and Y with density p ( t ) and p ( t ).As in the maximum model, we have plus induced and minus induced exchanges.In our analytical estimate we use again the assumption that the particle with an agedresidual waiting time is distributed by r ( t ) so that we have to calculate min ( X p ( t ) , Y r ( t ) )for the mixed scenario and min ( X r ( t ) , Y r ( t ) ) for the pure scenario.We again use the CFDs of equations (7) and (8) to calculate the minimum densityin both scenarios: f mixmin ( t ) = p ( t ) [1 − R ( t )] + r ( t ) [1 − P ( t )] (14)= ( γ − γ − < t < γ − γ − t − γ t ≥
1, (15) f puremin ( t ) = 2 r ( t ) [1 − R ( t )] (16)= γ − γ − (cid:16) − γ − γ − t (cid:17) < t < γ − γ − t − γ t ≥
1. (17)We realize that the largest tail is 3 − γ which we get for the pure scenario. Sinceboth scenarios can be observed in the exclusion process, each of them contributes to the3effective exchange time but the tail behavior is determined by the slower process thatcan block transport completely. We therefore set f min ( t ) = f puremin ( t ) = γ − γ − (cid:16) − γ − γ − t (cid:17) < t < γ − γ − t − γ t ≥
1. (18)Note that f min ( t ) leads to a increase in transport efficiency for γ > f min ( t ) becomes larger than − γ if γ >
3. This counter-intuitiveresult follows from the assumption that both waiting times of interfacing particlesare distributed by r ( t ) instead of p ( t ). However, the prediction would mean aslower exchange than in the asymmetric particle-based rule where minus particles arecompletely passive. We will see that our estimates actually describes the simulationresults for spd only for γ > γ = 2 .
5. Here,the tail behavior is well represented by p ( t ).In order to understand the origin of this difference in tail exponents, we measure a) pb min (cid:28) =2.5 b) mpdpb min pb min (cid:29) =4.0 (cid:30) =4.0spd mpd c) d) pb min (cid:31) =2.5spd ! " % & ’ ( t ) ( ff ) *+,/ p(t)r(t) (t)10 ; < > ( t ) t ? ff @ ABCD p(t)r(t) EFGH (t) 10 I K MN O PQ ff R STUV r(t)
WXYZ (t) [ ff \ ]^_‘ r(t)c b cdgh (t)10 j k lo p q t Figure 7:
Particle-based minimum rule:
PDFs of the effective exchange time for asystem of L = 100 for the particle-based minimum rule. (a) Single-particle dynamics,i.e. N + = 1 and N − = L − γ = 2 .
5. (b) Many-particle dynamics,i.e. N + = N − = L/ γ = 2 .
5. (c) Single-particle dynamics for exponent γ = 4 .
0. (d) Many-particle dynamics for exponent γ = 4 .
0. The dashed line shows thetail behavior of f min ( t ) close to the simulation data.4 -4 -3 -2 -1 L γ = rs γ = tu γ = 4 v -4 - w - y -1 z { L γ = |} γ = ~(cid:127) γ = 4 (cid:128) (cid:129)(cid:130) (cid:131)(cid:132) min max (cid:133) Figure 8: (a) Minimum and (b) maximum rule:
Length dependency of the particleflux J for different values of the exponent γ in the PDF (4) and mpd.the residual waiting time carried by the plus particle, which we call ˜ p ( t ). We find thatin the tail ˜ p ( t ) ≈ r ( t ) for γ > γ < γ = −
2. Calculating the minimum with such an exponent from ˜ p would leadto an f min ( t ) with exponents ≤ − γ , i.e. p ( t ) serves as a upper limit (see Appendix Bfor details).We can understand this deviation by realizing that the plus induced dynamics isgetting more important for if γ is below the critical value γ c . The influence of the tailin the residual waiting time is important for inducing events by minus particles that didstand in the queue for a long time. In contrast, the plus particle is more often responsiblefor inducing the events and consequently determines the effective exchange time. Weshow that the ratio of plus induced events is growing in this regime in Appendix B. Thetime average in the calculation of the residual waiting time in Appendix A is not validdue to temporal correlations in this scenario.For mpd, similar behavior is observed. The estimate f min is well suited to thesimulation result if γ >
3, which is shown in figure 7 (d). In results for small exponents γ < r ( t ).The relevance of the residual waiting time is caused by the dominance of long exchangetimes in long queues. Furthermore, passively exchanged particles keep their event timeafter an exchange process which leads to long range correlations of particles exchangetimes.Similar to the results of the site-based model and the maximum model, the fluxin the minimum rule does not show significant changes with L in the fast decayingregime above γ c , which is shown in figure 8 (a). For γ < γ c , again a dependency on thesystem length in the data supports the qualitative difference between the regimes foundin results from effective exchange times above. However, the flux clearly is higher inthe minimum rule than in the maximum rule (panel b)). We want to further study thedifference of the applied model rules on the particle flux in the next section.5 a) b) c) γ maxminsb γ maxminsb γ (cid:134) (cid:135) (cid:136) J min L= 100min L= 1000 Figure 9: Particle flux depending on the exponent γ for the three different update rules(site-based green, minimum orange, maximum blue). Error bars are drawn from semvalues of 10 realizations. (a) The system size is L = 100, (b) L = 1000. (c) Section for γ ≤ J for L = 100 (blue) and L = 1000 (red). The exponents for effective exchange times found for the different model rules aresummarized in table 1. We now compare the flux which is generated by the threesymmetrical updates in figure 9. Even though the asymptotic behavior of the threerules are similar (except for γ > γ c in the minimum model), the short exchange timesinfluence the value of the flux. This leads to significant differences between the site-basedmodel, the particle-based minimum model and the particle-based maximum model for γ > γ c . As expected, the maximum rule is really slower than the site-based model andthe minimum rule can enhance the transport. For γ < γ c however, the flux generatedin the maximum rule is close to the flux by the site-based rule. In the minimum rule,we still measure higher fluxes, both for the system size of L = 100 (a) and L = 1000(b). Note that we observe finite-size effects in these results which is shown in panelc) where the data points deviate for the different system sizes. By the analysis of tailexponents, we expect J = 0 for the infinite system in the stationary state such as in theother models.Table 1: The resulting exponent seen in the effective waiting time PDFs in the differentupdate rules for single-particle dynamics and many-particles dynamics.single-particle dynamics many-particle dynamicssite-based 1 − γ − γ particle-based asymmetrical − γ − γ particle-based maximum 1 − γ − γ particle-based minimum − γ − γ for γ < γ > − γ − γ for γ < γ >
4. Conclusion
In our contribution, we analyzed different bidirectional variants of the TASEP withnon-Markovian exchange dynamics. These models are relevant for one-dimensionaltransport problems in crowded environments, where the high density of particle clustersleads to small effective exchange rates of particle positions. A possible realization ofthe bidirectional transport model would include two oppositely moving particle speciesand holes. A model of this kind would combine a Markovian particle-dynamics, whichwould be applied when the particles move toward an empty site and a non-Markovianparticle-dynamics, which governs the particle-exchange. Simulation results show strongcondensation of the particles, which implies that the bidirectional transport capacity isdetermined by the efficiency of the exchange processes rather than by the time spent inthe low density area. Therefore, we restricted our analysis to symmetric and fully-filledsystems. This choice reduces considerably the corrections to scaling for small systemsizes.Modeling bidirectional transport of active particles with lattice gases allowsassigning the exchange times to the particles as well as to the lattice. In the latter case,we can map the problem to the uni-directional process, since pairs of oppositely movingparticles behave as particles and holes in the uni-directional case. This is even true forspd which correspond to a uni-directional system with a single hole where the dynamicsof the particle is governed by the residual waiting time. Significant differences to theuni-directional case exist if the exchange times are assigned to the particles. This canbe realized in a symmetric or in an asymmetric way, wherein the latter case the reactiontimes are assigned to only one-particle species. The asymmetric case implies that oneparticle species can be assumed to be passive which is why we find the asymptoticsof the uni-directional non-Markovian TASEP for the single-particle as well as for themany-particle dynamics.Assigning exchange times symmetrically to both particle species implies that twoexchange times are given for a pair of oppositely moving particles. Therefore, onehas to define an additional selection rule. In this work we have chosen two extremecases that preserve the symmetry between the two types of particles, i.e. either theminimum or the maximum of the two waiting-times will be selected. In case of themaximum rule, residual waiting time and f max ( t ) show the same asymptotics. Thereforethe maximum rule modifies indeed the effective exchange time distribution but to theexponent compared to the residual waiting time of the site-based model. Significantdifferences exist only in the γ > γ c regime for the minimum rule. In this case, we findthat the asymptotics is governed by the distribution of independent residual times wherethe asymptotics of effective exchange time distribution is given by 3 − γ . For γ < γ c this is not the leading contribution. Here, the asymptotic behavior is in accordancewith the asymmetric particle based model which is given by the residual waiting timefor mpd. This effect is the result of passively moving particles which keep the assignedexchange time. For small values of γ the dynamic is, as for the other cases, governed7by pairs of particles with long residual waiting times.Our results underline the universality of the findings which have been discussedfor the uni-directional non-Markovian TASEP [29, 30]. In this class of models, manyparticle effects generically lead to a dominant contribution of the residual waiting timefor γ < γ c . Here, the configurations are characterized by large particle clusters and a sizedependent flow of particles. The transport capacity of large systems in this parameterregime is extremely low compared to their Markovian counterparts. For γ > γ c however,we observe homogeneous particle configurations and size-independent values of the flowwhich differ for the different implementations of the dynamics.Our findings can be relevant for bidirectional flows under strong confinement asfor example in narrow escape problems in pedestrian dynamics [5] or intracellulartransport in axons and dendrites [35] where the effective exchange dynamics can benon-Markovian. Acknowledgments
This work was funded by the Deutsche Forschungsgemeinschaft (DFG) through Collab-orative Research Center SFB 1027 (Project A8).8
Appendix A. Calculation of residual waiting times
The waiting time PDF for a single particle site exchange in the site-based model iscalculated by using renewal theory following [36]. In particular, we will determine theresidual waiting time until the next exchange event occurs if two particles are in thelocal (+ − ) configuration. We start with a renewal process for renewal waiting times X n distributed by Eq. (4), which are given to a site in the lattice. The N -th renewal of thewaiting time on this site occurs at time S N = N X n =1 X n , (A.1)i.e. we can count the number of passed renewal events N ( t ) at each time t .For a time t > S N ( t ) , the waiting times of this site are called the duration ofrenewal time intervals ˜ X ( t ) = X N ( t )+1 = S N ( t )+1 − S N ( t ) (see figure A1). The renewalprocess also has an age Z ( t ) = t − S N ( t ) as well as a residual life (residual waiting time) Y ( t ) = S N ( t )+1 − t until the next renewal event takes place at time S N ( t )+1 . The residualwaiting time is therefore also written as Y ( t ) = ˜ X ( t ) − Z ( t ).We now calculate the time averaged CDF of the residual waiting time Y ( t ), i.e. F Y ( y ) = Pr { Y ( t ) ≤ y } that gives the fraction of time that the residual waiting time issmaller than a given y . We can invent an indicator reward function R ( t ) to determine S S S S N-1 S N S N+1 t0 Y(t)Z(t)X X X X N X N+1 timetimeS S S S N-1 S N S N+1 X X X X N X N+1 Y ( t ) Figure A1:
Top:
The renewal process is determined by the time series X n , n ∈ N ,build from the algebraic waiting time PDF in Eq. (4). Summing up these waiting times S , S , ... S N gives the time for the next event at S N +1 . For a time S n ≤ t ≤ S N +1 , thecurrent process has the age Z and the residual life Y . Bottom:
The residual life Y ( t )is a step wise function of time, decaying from X n to 0 during the time in the intervalbetween S n − and S n .9if the residual waiting time is actually smaller or not, i.e. R ( t ) = R ( Z ( t ) , ˜ X ( t )) = ( X ( t ) − Z ( t ) ≤ y Y ( t ) can then be found by using the time averageover the indicator function F Y ( y ) = lim t →∞ t Z t R ( τ )d τ = 1 X Z x = yx =0 Pr { X > x } d x, (A.3)where X = γ − γ − denotes the mean value of the renewal event duration for the PDF inEq. (4).We now use this framework to determine the residual waiting time for the renewalprocess with algebraic waiting time PDF Eq. (4), i.e. f X ( x ) = ( < x < γ − x − γ x > F X ( x ) = ( x < − x − γ x ≥
1. (A.4)We use equation A.3 to determine the CDF F Y ( y ) = 1 X Z y (cid:0) − (1 − x − γ )Θ( x − (cid:1) d x (A.5)= γ − γ − y y < γ − γ − (cid:18) γ − γ − − γ y − γ (cid:19) y ≥
1, (A.6)and finally arrive at the result Y ( y ) = γ − γ − y < γ − γ − y − γ y ≥ Y ( y ).0 Appendix B. Additional measurements for the minimum rule insingle-particle dynamics
In this appendix, we further examine results in the single-particle dynamics, minimumupdate rule. As we have seen in figure 7, the effective waiting time distribution followsthe estimation f puremin only for exponents γ > f mixmin , which is expected since the plus particle always has a non-aged waiting time. This result is in contrast to the result of the bulk text, which is alsoshown in figure B1 (b) for comparability. For the minimum rule, effective waiting timesdo not follow f mixmin , hence the age of the plus particle plays a role. However, neither dothey follow f puremin which is expected for the minimum of to random variables distributedby the residual time r ( t ).In a second step, we show simulation results for residual waiting times of the plusparticle in figure B2, which we call ˜ p ( t ) in the following. In the fast regime of γ >
3, themeasurements follow the theoretical estimate r ( t ). However, this changes for exponents γ <
3. The simulation result do not follow r ( t ) anymore but rather stay close tothe asymptotic of γ = 2. This slope for ˜ p ( t ) is consistent with the effective exchangetime observed in figure 7 (a), where simulation results follow p ( t ) ≈ t − γ in the tail,when considering the minimum out of a random variable distributed by r ( t ) for minusparticles and ˜ p ( t ) for plus.We see that the time-averaged estimate r ( t ) is not valid anymore for the residualwaiting time of a single plus particle in the minimum model for γ <
3. We further givean argument for the break down of validity by showing that the fraction of exchange a) b) delete plus times pb min γ =2.4 γ =2.4spd spd -8 -6 -4 -2 -1 P ( t ) t eff. timep(t) (cid:137)(cid:138)(cid:139)(cid:140) (t) (cid:141)(cid:142)(cid:143)(cid:144)(cid:145)(cid:146)(cid:147)(cid:148)(cid:149)(cid:150)(cid:151) (t) 10 -8 -6 -4 - (cid:152) -1 (cid:153) (cid:154) t p(t) (cid:155)(cid:156)(cid:157)(cid:158) (t) (cid:159)(cid:160)¡¢£⁄¥ƒ§¤' (t) “ ff « ‹›fifl Figure B1: a) Single-particle dynamics exclusion process where residual waiting timesof the plus particle is deleted after each exchange process, independent whether the plusparticle was active or passive in the exchange. b) The normal spd minimum rule fromthe bulk text for comparison.1 a) pb min γ =3.5 b) spd pb min pb min γ =2.7 γ (cid:176)–†‡ spd spd c) d) pb min γ =4.2spd pb min γ =2.2spd ·(cid:181) ¶• -8 -6 -4 - ‚ -1 „ ” ( t ) p »…‰ (cid:190) t ¿ p(t)r(t)t - (cid:192) -8 -6 -4 - ` ´ˆ˜ ¯ t ˘ p(t)r(t)t - ˙ -8 -6 -4 - ¨ -1 (cid:201) ˚ ( t ) p ¸(cid:204)˝ ˛ t ˇ p(t)r(t)t - — -8 -6 -4 - (cid:209) -1 (cid:210) (cid:211) p (cid:212)(cid:213)(cid:214) (cid:215) t (cid:216) p(t)r(t)t - (cid:217) -8 -6 -4 - (cid:218) -1 (cid:219) (cid:220) ( t ) t p (cid:221)(cid:222)(cid:223) (cid:224) t Æ p(t)r(t)t - (cid:226) -1 ª (cid:228) γ r a t i o p l u s i ndu c ed t > 0t > 1 t > 30t > 100 Figure B2:
Particle-based minimum rule
PDFs of the residual waiting time of aplus particle ˜ p for a system of L = 100 filled by N + = 1 and N − = L − γ = 4 .
2, (b) γ = 3 .
5, (c) γ = 3 .
1, (d) γ = 2 .
7, (e) γ = 2 .
2. Blue data points show theresidual waiting times of the plus particle, red p ( t ) original waiting time distribution,green r ( t ) residual waiting time from renewal theory, yellow a constant function t − as areference line. f) The measured ratio of exchange processes which have been induced byan active plus particle in the spd minimum case. Statistics over all exchange events arecolored blue, exchange processes with an effective waiting time of at least 1 are orange,30 green and at least 100 red. Data points are missing if no such high waiting timeshave been observed in the simulation, errors bars show the sem.2events induced by the plus particle increases for γ < Appendix C. Update scheme
For the Markovian TASEP: We also remark that, thanks to the memoryless property ofthe Markovian TASEP, one can practically generate the next time t k +1 i at every t = t ki .The Markov property does not hold for the algebraic distribution p ( t ) in equation (4).To evolve the system in time, we use a modified waiting time algorithm ( next reactionmethod [33]) similar to [29, 34] .Times for all events (particle-based or site-based) are initialized at the beginningof the simulation ( τ = 0). The shortest time t α is then chosen from a list of all waitingtimes to be the absolute time for the next event α . In the realization of the event, thesystem time is increased up to this point in time τ i +1 = τ i + t α . Anyhow, the processis only executed if the local particle configuration is appropriate. After the realization,the waiting time of the event is renewed by taking a new time t new from the distribution p ( t ) added to the current system time t α = τ + t new . This time is then placed into thelist for the event α and the procedure is repeated.3 [1] Mohcine Chraibi, Antoine Tordeux, Andreas Schadschneider, and Armin Seyfried. Modellingof pedestrian and evacuation dynamics. In Complex Dynamics of Traffic Management , pages649–669. Springer US, New York, NY, 2019.[2] Livio Gibelli and Nicola Bellomo.
Crowd Dynamics, Volume 1: Theory, Models, and SafetyProblems . Springer, 2019.[3] Dirk Helbing. Collective phenomena and states in traffic and self-driven many-particle systems.
Selected papers of the Twelfth International Workshop on Computational Materials Science(CMS2002) , 30(1):180–187, May 2004.[4] Gunnar G. Lvs. Modeling and simulation of pedestrian traffic flow.
Transportation Research PartB: Methodological , 28(6):429–443, 1994.[5] Andreas Schadschneider, Debashish Chowdhury, and Katsuhiro Nishinari.
Stochastic transport incomplex systems: from molecules to vehicles . Elsevier, 2010.[6] Debashish Chowdhury, Ludger Santen, and Andreas Schadschneider. Statistical physics ofvehicular traffic and some related systems.
Physics Reports , 329(4):199–329, 2000.[7] Boris S. Kerner. Complex dynamics of traffic management, introduction to. In
Encyclopedia ofComplexity and Systems Science , pages 1177–1180. Springer New York, New York, NY, 2009.[8] Mithila Burute and Lukas C. Kapitein. Cellular logistics: Unraveling the interplay betweenmicrotubule organization and intracellular transport.
Annu. Rev. Cell Dev. Biol. , 35(1):29–54, October 2019.[9] Cecile Appert-Rolland, Maximilian Ebbinghaus, and Ludger Santen. Intracellular transport drivenby cytoskeletal motors: General mechanisms and defects.
Physics Reports , 593:1–59, 2015.[10] Debashish Chowdhury. Stochastic mechano-chemical kinetics of molecular motors: Amultidisciplinary enterprise from a physicists perspective.
Stochastic mechano-chemical kineticsof molecular motors: A multidisciplinary enterprise from a physicists perspective , 529(1):1–197,August 2013.[11] William O. Hancock. Bidirectional cargo transport: moving beyond tug of war.
Nature ReviewsMolecular Cell Biology , 15(9):615–628, September 2014.[12] Nobutaka Hirokawa, Shinsuke Niwa, and Yosuke Tanaka. Molecular motors in neurons: Transportmechanisms and roles in brain function, development, and disease.
Neuron , 68(4):610–638,November 2010.[13] Nobutaka Hirokawa and Reiko Takemura. Molecular motors and mechanisms of directionaltransport in neurons.
Nature Reviews Neuroscience , 6:201, February 2005.[14] Anthony Brown. Axonal transport of membranous and nonmembranous cargoes.
J Cell Biol ,160(6):817, March 2003.[15] T. Chou, K. Mallick, and R. K. P. Zia. Non-equilibrium statistical mechanics: from a paradigmaticmodel to biological transport.
Reports on progress in physics , 74(11):116601, 2011.[16] David Mukamel. Phase transitions in nonequilibrium systems.
Soft and fragile matter , page 237,2000.[17] Bernard Derrida. An exactly soluble non-equilibrium system: the asymmetric simple exclusionprocess.
Physics Reports , 301(1-3):65–83, 1998.[18] Martin R. Evans, Damien P. Foster, Claude Godrche, and David Mukamel. Spontaneous symmetrybreaking in a one dimensional driven diffusive system.
Physical review letters , 74(2):208, 1995.[19] Peter F. Arndt, Thomas Heinzel, and Vladimir Rittenberg. Spontaneous breaking of translationalinvariance in one-dimensional stationary states on a ring.
Journal of Physics A: Mathematicaland General , 31(2):L45, 1998.[20] Philip Greulich and Andreas Schadschneider. Single-bottleneck approximation for driven latticegases with disorder and open boundary conditions.
Journal of Statistical Mechanics: Theoryand Experiment , 2008(04):P04009, 2008.[21] Carsten Burstedde, Kai Klauck, Andreas Schadschneider, and Johannes Zittartz. Simulationof pedestrian dynamics using a two-dimensional cellular automaton.
Physica A: StatisticalMechanics and its Applications , 295(3-4):507–525, 2001. [22] Stefan Nowak and Andreas Schadschneider. Quantitative analysis of pedestrian counterflow in acellular automaton model. Physical Review E , 85(6):066128, 2012.[23] Stefan Klumpp and Reinhard Lipowsky. Phase transitions in systems with two species of molecularmotors.
EPL (Europhysics Letters) , 66(1):90, 2004.[24] Maximilian Ebbinghaus and Ludger Santen. A model for bidirectional traffic of cytoskeletalmotors.
Journal of Statistical Mechanics: Theory and Experiment , 2009(03):P03030, 2009.[25] M. Ebbinghaus, C. Appert-Rolland, and L. Santen. Particle interactions and lattice dynamics:scenarios for efficient bidirectional stochastic transport?
Journal of Statistical Mechanics:Theory and Experiment , 2011(07):P07004, 2011.[26] Jean-Philippe Bouchaud. Weak ergodicity breaking and aging in disordered systems.
Journal dePhysique I , 2(9):1705–1713, 1992.[27] Ccile Monthus and Jean-Philippe Bouchaud. Models of traps and glass phenomenology.
Journalof Physics A: Mathematical and General , 29(14):3847, 1996.[28] Mark E. J. Newman. Power laws, pareto distributions and zipf’s law.
Contemporary physics ,46(5):323–351, 2005.[29] Robert J. Concannon and Richard A. Blythe. Spatiotemporally complete condensation in a non-poissonian exclusion process.
Physical review letters , 112(5):050603, 2014.[30] Diana Khoromskaia, Rosemary J. Harris, and Stefan Grosskinsky. Dynamics of non-markovianexclusion processes.
Journal of Statistical Mechanics: Theory and Experiment , 2014(12):P12013,2014.[31] Rui Jiang, Katsuhiro Nishinari, Mao-Bin Hu, Yong-Hong Wu, and Qing-Song Wu. Phaseseparation in a bidirectional two-lane asymmetric exclusion process.
Journal of StatisticalPhysics , 136(1):73–88, July 2009.[32] Daniel T. Gillespie. A general method for numerically simulating the stochastic time evolution ofcoupled chemical reactions.
Journal of computational physics , 22(4):403–434, 1976.[33] Michael A. Gibson and Jehoshua Bruck. Efficient exact stochastic simulation of chemical systemswith many species and many channels.
The journal of physical chemistry A , 104(9):1876–1889,2000.[34] Mieke Gorissen and Carlo Vanderzande. Ribosome dwell times and the protein copy numberdistribution.
Journal of Statistical Physics , 148(4):628–636, 2012.[35] O. A. Shemesh and M. E. Spira. Paclitaxel induces axonal microtubules polar reconfigurationand impaired organelle transport: implications for the pathogenesis of paclitaxel-inducedpolyneuropathy.
Acta Neuropathologica , 119(2):235–248, February 2010.[36] Robert G. Gallager.