Spontaneous symmetry breaking in a bridge model fed by junctions

Vladislav Popkov^1,2, Martin R Evans³ and David Mukamel⁴

¹ Dipartimento di Fisica E R Caianiello', and Consorzio Nazionale Interuniversitario per le Scienze Fisiche della Materia (CNISM), Università di Salerno, Baronissi, Italy
² Interdisziplinäres Zentrum fur Komplexe Systeme, Römerstrasse 164, D-53117 Bonn, Germany
³ SUPA, School of Physics, University of Edinburgh, Mayfield Road, Edinburgh EH9 3JZ, UK
⁴ Department of Physics of Complex Systems, The Weizmann Institute of Science, Rehovot 76100, Israel

E-mail: popkov@sa.infn.it, m.evans@ed.ac.uk and david.mukamel@weizmann.ac.il

Received 29 July 2008, in final form 11 September 2008
Published 1 October 2008

1. Introduction and model definition

Nonequilibrium stationary states (NESS), in which probability currents are supported and invariant measures are not generally of Gibbs–Boltzmann form, offer many surprising phenomena not seen in equilibrium systems. Examples include boundary-induced phase transitions and spontaneous symmetry breaking (SSB) in one-dimensional systems. The minimal models for NESS are driven particle models such as the totally asymmetric exclusion process (TASEP) and its relatives. Some exact solutions for these models have been found, for example the TASEP, leading to an understanding of boundary-induced phase transitions [1]. However, important questions remain as to the nature of SSB in driven one-dimensional systems.

The SSB phenomenon was first observed in what is often referred to as the bridge model [2], and subsequently in some other models [3, 4]. The bridge model comprises a one-dimensional lattice on which positive particles (pluses) move to the right and negative particles (minuses) move to the left. At the left boundary of the lattice, pluses may enter the lattice and minuses may leave; at the right boundary of the lattice minuses may enter the lattice and pluses may leave. When the exit rate at which particles may leave the lattice is lowered, the stationary state changes from a symmetric one, in which the currents and bulk densities of pluses and minuses are equal, to a symmetry-broken state where there is a majority species of particle with larger current and higher bulk density than the minority species. However, the symmetry-broken state has not been solved exactly except in the limit in which the exit rate tends to zero where it has been rigorously proven that SSB occurs [7]. Also a model with deterministic bulk dynamics exhibiting a symmetry-broken phase has been solved exactly [5, 6].

The transition to this symmetry-broken state remains a subject of debate. A mean-field theory originally predicted that the transition should occur via an intermediate, weakly symmetry-broken phase [8]. Monte Carlo simulations have shown that the mean-field theory does not correctly predict the position of the transition [9]; however, at least on finite systems, an intermediate phase is seen and the transition from symmetric to the strongly symmetry-broken phase occurs through a sequence of transitions [10]. This sequence has also been observed in a related model [4]. However, it has been suggested that the region of parameter space occupied by the intermediate phase and over which the sequence of transitions occurs disappears in the infinite system size limit [11].

The failure of the mean-field approximation to exactly predict the phase diagram can be traced back to the boundary conditions of the bridge model which do not correspond to particle reservoirs at fixed density. Instead there are effective impurities at the boundaries [12]. Also the hydrodynamic limit of the bridge model is not well defined over the whole phase space. In particular, it breaks down at the SSB transition in the limit of small input and exit rates [12].

In this paper we introduce a new class of bridge models demonstrating the SSB phenomenon. In these models the input and output of the pluses and minuses are governed by TASEPs. Thus one can think of the ends of the bridge as junctions where TASEPs for the pluses and minuses merge. The input and output rates at the bridge are not external parameters as in TASEP; rather, they are determined self-consistently by the dynamics of the bridge and its feeding segments. The models exhibit complex phase diagrams including a phase co-existence region showing hysteresis and quite special triple phase co-existence point, which appear to be correctly predicted by a mean-field theory.

Our model is defined as follows. We consider a chain of the length L occupied by two species of particles, `plus' particles moving to the right and `minus' particles moving to the left, with hardcore exclusion and random sequential update. The chain is region II of figure 1 and we shall refer to it as the bridge. At each end of the bridge, there are junctions where the chain splits into two parallel segments, one containing only plus particles and holes and the other only minus particles and holes (sections I and III in figure 1). At the external boundaries of the parallel segments, usual TASEP rules are applied. Namely, plus particles are injected into their segment of section I at the left with rate α, if the first site is empty, and are removed with rate β from the right end (upper segment of section III in figure 1), if the last site is occupied. Likewise minus particles are injected into their segment of section III with rate α and removed with rate β from the left end of their segment of section I. In sections I and III pluses hop to the right inside their segments with rate 1. Plus particles enter the bridge at the left junction if the first site of section II is empty and leave it at the other junction if the first site of the plus segment of section III is empty, both with rate K. Inside the bridge, particles exchange with empty sites and with the minuses with the same rate 1. Similarly, inside the lattice minus particles hop to the left with rate 1 in sections I and III, enter and leave the bridge if the entrance or exit site is empty with rate K, and inside the bridge exchange with both pluses and empty sites with rate 1. The model is symmetric with respect to simultaneous charge inversion and left–right reflection.

In this work we consider only K  =  1, leaving α and β as the model parameters. We calculate the phase diagram of the model using a mean-field approximation and direct simulation of the dynamics.

2. Phase diagram and phase transition lines

Before discussing the phase diagram of our model, it is useful to recall the phase diagram of the TASEP (for plus particles moving to the right). The phases are distinguished by the large system size limits of the expressions for ρ the bulk density of particles far from the boundaries and j  =  ρ(1 – ρ) the current of particles. When α < 1/2 and β > α the low-density (LD) phase occurs where the density in the bulk is equal to α and is determined by the left boundary; when β < 1/2 and α > β the high-density (HD) phase occurs where the density in the bulk is equal to 1 – β and is determined by the right boundary; when α > 1/2 and β > 1/2 the maximal current (MC) phase occurs where the bulk density is 1/2.

We now present the phase diagram of the bridge model fed by junctions obtained from a mean-field analysis detailed in the following section. The resulting α–β phase diagram is given in figure 2; it contains three different phases:

The phase diagram (figure 2) is reminiscent of the TASEP phase diagram discussed above, but with several important differences. First the LDS2 phase replaces the maximal current phase and the transition lines to the LDS2 phase are at α  =  1/3, β  =  1/3 rather than α  =  1/2, β  =  1/2. Second, the SSB phase replaces the high-density phase of the TASEP. Clearly a high-density symmetric phase is not possible in the present model since in that phase the densities of both species on the bridge section have to be greater than 1/2. Third, there is a co-existence region (1/3 < β < 2/5, α > β) where the system may be in either of the LDS2 and SSB phases; this produces an interesting kink in the phase boundary of the SSB phase.

We now discuss the phase transitions and transition lines separating the phases described above. For the TASEP, as well as other nonequilibrium steady states with a non-zero conserved current, the `order' of a phase transition is determined by the non-analyticity of the current at the transition, i.e. the order of the derivative of the current at which a singularity appears [1]. In the present model we find some novel features for the nonequilibrium phase transitions. Notably, across all the transition lines some of the bulk densities change discontinuously, as can be seen by comparing stationary densities in different phases. Indeed for the LDS1/SSB transition there is even a discontinuity in the current of the minority species. We now describe the behaviors at the transitions.

2.1.  LDS1/LDS2 transition line: α  =  1/3, β > α

The bulk densities in the bridge are continuous across the transition and the currents are continuous. In section I the plus density jumps discontinuously across the transition from 1/3 to 2/3, whereas the minus density is continuous. Similarly in section III, the minus density jumps discontinuously from 1/3 to 2/3. On the transition line one finds a shock in section I separating regions of plus densities 1/3 and 2/3 and a shock in section III separating regions of minus densities 2/3 and 1/3.

2.2.  LDS1/SSB transition line: α  =  β < 1/3

The bulk densities in all sections are discontinuous across the transition: the bulk density of the majority species jumps from α to 1 – α in all sections; the bulk density of the minority species (taken as minus) jumps from α to α/2 in sections I and II and from α to 1 – α/2 in section III. The majority current is continuous whereas the minority current jumps from α(1 – α) to α/2(1 – α/2). On the transition line two types of shock configurations are in fact observed. First there is symmetric configuration where the plus density is α in sections I and II, and in section III there is a shock separating regions of plus densities α and 1 – α. Similarly the minus density is α in sections II and III, and in section I there is a shock separating regions of minus densities 1 – α and α. Second, there is an asymmetric shock configuration where, taking the majority species to be plus, the bulk plus density is 1 – α in sections II and III and in section I there is a shock between regions of plus density α and 1 – α. The minority species density, taken to be minus, is α/2 in sections I and II and 1 – α/2 in section III. Here the currents of plus and minus are unequal.

2.3.  LDS2/SSB co-existence region: α > β, 1/3 ≤ β ≤ 2/5

In this region, denoted `SSB + LDS2' in the phase diagram (figure 2), the system can be in either of the two symmetry-related SSB states or the LDS2 phase. In the stochastic model of infinite size, the three phases are stable. In finite systems flips between the phases take place, with typical time between the flips growing exponentially with the system size L (this issue will be addressed in section 4). Note that the spatial co-existence of the phases (as in [3], for example) is not possible since they carry different currents.

In the infinite system, the co-existence region SSB + LDS2 entails hysteresis, which is observed when keeping the injection rate α > 2/5 constant and increasing/decreasing adiabatically the extraction rate β from 0 (SSB phase) to 1 (LDS2) and back. Indeed, along the path 0 → 1 the SSB/LDS2 transition happens at β  =  2/5, while on the backward path it will happen at β  =  1/3. The hysteresis is illustrated in figure 4.

A kink point with coordinates α  =  β  =  2/5 on the phase diagram is an unusual one. In phase diagrams of systems at equilibrium, such a kink would imply that there must be another first-order transition line emerging from that point, making it a triple point. This is a result of the Clausius–Clapeyron relation. Our model is out of equilibrium, and the existence of a kink, although unusual, does not violate any rule.

3. Mean-field solution

We now derive the phase diagram through a mean-field approximation commonly used for one-dimensional stochastic systems, wherein two-point correlation functions are replaced by products of one-point correlation functions. For the TASEP this mean-field approximation is known to predict the correct phase diagram and stationary bulk densities [13].

Let us denote the stationary density of plus particles at site k as p_k. Then, the mean-field approximation gives the stationary flux of pluses as

where S  =  I, II, III and p_S is the limiting bulk density far away from the boundaries in segment S. Similarly, the density of minus particles at site k is m_k and we write

In the stationary state the bulk densities p_S, m_S can be different in different segments S but the currents j ₊ , j_– have the same values everywhere. At the boundaries the currents read

and at the junctions,

Following [8, 14], we define effective entrance rates α ⁺ _S, α^–_S and exit rates β ⁺ _S, β^–_S for each segment. For plus particles

The rates α^–_S, β^–_S for minus particles are obtained by the substitutions I ? III and p_k → m_3L–k + 1. For each species, each segment S can be viewed as a TASEP model with, e.g. for the pluses, the effective injection/extraction rates α ⁺ _S, β ⁺ _S. In the large segment length limit, the bulk density of pluses p_S in segment S can be read from the phase diagram of a TASEP [15, 16], namely

In the HD phase the density profile is flat at the right boundary and in the LD phase the density profile is flat at the left boundary. For minuses, the corresponding bulk densities are given by (9)–(11) with substitutions  + →– and p → m. Thus, in the large segment length limit, each segment of the bridge model must exhibit one of the solutions (9)–(11). The possible solutions of the mean-field equations (1)–(6) (MFE) are listed below.

The question of the accuracy of the MF solution remains. Discrepancies between the predicted stationary bulk densities and Monte Carlo simulations lie within numerical error bars, except for the LDS2 phase. In this phase the bulk density from Monte Carlo simulations p_II ≈ 0.339 appears to differ slightly from that predicted by mean field p_II  =  1/3. This discrepancy derives from the fact that while the LDS1 and SSB phases are controlled by uncorrelated inputs from the outer boundaries, the LDS2 phase is controlled by the injection rate from the junctions which could imply correlations in the input. It appears that such correlations result in only a small shift of the LDS2 LDS1 and LDS2 → SSB phase boundaries to α ≊ 0.339 < β and β ≊ 0.339 < α, while other phase boundaries remain as in figure 2. In view of this, we believe that the mean-field predictions are qualitatively correct and are quantitatively very accurate.

4. Mechanisms of flips between symmetry-broken states in the SSB phase

In a finite system, flips between the two symmetry-broken states (majority and minority species interchange) are observed (figure 5(a)). To demonstrate spontaneous symmetry breaking in the SSB phase and the SSB + LDS2 region (figure 2), we need to show that the flipping times between the stable states in finite systems diverge exponentially with the system size L. For the SSB phase the mechanism is best illustrated in the SSB region with the additional condition α < 1/2. In this region the plus profile in section I of figure 3(b) is produced by a domain wall (a shock) between a region adjacent to the boundary with density α and a region with density 1 – β. The domain wall is biased to the left and therefore sticks to the left boundary. (Note that the motion of the plus domain wall does not affect the state of minus particles.) In the SSB configuration the entrance of minuses to the bridge is partially blocked by the high density 1 – β of pluses. To flip the SSB configuration (figure 3(b)), the domain wall of pluses p_I  =  1 – β sticking to the left boundary must retreat (against its bias) back to segment III in order to give an opportunity for the minuses to take over the key middle segment, block the entrance of pluses on it and make the flip. The time to wait for such an improbable event grows exponentially with the system size L. Roughly, one can estimate the flipping time as follows. The shock of pluses (α, 1 – β) is driven to the right with the rate r  =  β(1 – β)/(1 – β – α) and to the left with the rate l  =  α(1 – α)/(1 – β – α). In the SSB region (13) l > r. Let us represent the shock position with a phantom particle which hops to the nearest left/right site with the rates l and r. Introducing the stationary probability η_k for the phantom particle to be at site k, where k  =  0, 1,  ..., 2L we find η_k l  =  η_k–1 r so that η_k  =  η₀(r/l)^k. We can estimate the flipping time by 1/η_2L,^Note1 giving

In the rest of the SSB region (for α > 1/2) the mechanism is similar, but the density associated with the left-hand boundary is 1/2 rather than α. Consequently, in the flipping time estimate α(1 – α) has to be substituted with 1/4 for this case.

For the region SSB + LDS2 with three stable phases, flips between the steady states are also caused by the shock motion, but their mechanism is different. Firstly, there is no direct transition between the two SSB states, only through the intermediate phase SSB1 ⟺ LDS2 ⟺ SSB2, see figure 5(b). The existence of the intermediate LDS2 state makes the transition much easier: indeed, to flip from the LDS2 to SSB state, a shock need only cross one segment as opposed to two segments in direct transitions between the two symmetry-broken states. This explains the big difference in system sizes between figures 5(a) and (b).

5. Conclusions

In this work we have introduced a bridge model fed by junctions where the input and output streams of the bridge are themselves TASEPs. This has allowed the phase diagram predicted by the mean-field approximation to be confirmed by numerical simulations as being qualitatively correct and rather accurate quantitatively. The phase diagram exhibits a number of novel and interesting features that we discussed in section 2. We note here some further points that will be pursued in future work.

A particularly interesting feature of the phase diagram is that we have two types of first-order transitions: one with a co-existence line (α  =  β < 1/3) and the other with a co-existence region (marked LDS2 + SSB in figure 2). These transitions replace the transition to the symmetry-broken phase via an intermediate phase exhibited in the bridge model mean-field phase diagram. In both cases the co-existence is between different possible steady states for the systems. The dynamics of how the system flips between these steady states is rather intricate but can be understood in terms of domain wall dynamics (work in progress).

In the phase diagram (figure 2) a phase resembling the maximal current phase of the TASEP is absent. We also note that none of the solutions of the MFE supports maximal current j ₊   =  j_–  =  1/4 in the system. The reason is that for the case K  =  1 considered here, the junctions act as effective bottlenecks. For larger K (more precisely, for K > 2), a new phase reminiscent of the MC current phase in the TASEP does appear.

Acknowledgments

This work has been supported by DFG within project KR 1123/1-2, the Albert Einstein Minerva Center for Theoretical Physics, and the Israel Science Foundation (ISF). The authors thank The Isaac Newton Institute, Cambridge, where the major part of this work was done, for hospitality during the Principles of the Dynamics of Nonequilibrium Systems programme. VP thanks G Schütz for discussions.

[1] 

Blythe R A and Evans M R 2007 J. Phys. A: Math. Theor. 40 R333–441 
IOPscience

[2] 

Evans M R, Foster D P, Godrèche C and Mukamel D 1995 Phys. Rev. Lett. 74 208 
CrossRefPubMed

[3] 

Popkov V and Peschel I 2001 Phys. Rev. E 64 026126 
CrossRef

[4] 

Pronina E and Kolomeisky A B 2007 J. Phys. A: Math. Theor. 40 2275–86 
IOPscience

[5] 

Willmann R D, Schutz G M and Grosskinsky S 2005 Europhys. Lett. 71 542–8 
IOPscience

[6] 

Grosskinsky S, Schutz G M and Willmann R D 2007 J. Stat. Phys. 128 587–606 
CrossRef

[7] 

Godreche C, Luck J-M, Evans M R, Mukamel D, Sandow S and Speer E R 1995 J. Phys. A: Math. Gen. 28 6039 
IOPscience

[8] 

Evans M R, Foster D P, Godrèche C and Mukamel D 1995 J. Stat. Phys. 80 69 
CrossRef

[9] 

Arndt P F, Heinzel T and Rittenberg V 1998 J. Stat. Phys. 90 783 
CrossRef

[10] 

Clincy M, Evans M R and Mukamel D 2001 J. Phys. A: Math. Gen. 34 9923–37 
IOPscience

[11] 

Erickson D W, Pruessner G, Schmittmann B and Zia R K P 2005 J. Phys. A: Math. Gen. 38 L659–65 
IOPscience

[12] 

Popkov V and Schütz G M 2004 J. Stat. Mech. P12004 
IOPscience

[13] 

Derrida B, Domany E and Mukamel D 1992 J. Stat. Phys 69 667 
CrossRef

[14] 

Pronina E and Kolomeisky A B 2005 J. Stat. Mech. P07010 
IOPscience

[15] 

Derrida B, Evans M R, Hakim V and Pasquier V 1993 J. Phys. A: Math. Gen. 26 1493 
IOPscience

[16] 

Schütz G and Domany E 1993 J. Stat. Phys. 72 277–96 
CrossRef

Notes

Note1

 More accurately, the flipping time is related to the mean first passage time of the phantom particle from site k  =  0 to the site k  =  2L. This will not change the results qualitatively.