Table Of Content

Particle current fluctuations in a variant of asymmetric Glauber model S. R. Masharian,1,∗ P. Torkaman,2,† and F. H. Jafarpour2,‡ 1Islamic Azad University, Hamedan Branch, Hamedan, Iran 2Physics Department, Bu-Ali Sina University, 65174-4161 Hamedan, Iran We study the total particle current fluctuations in a one-dimensional stochastic system of clas- sical particles consisting of branching and death processes which is a variant of asymmetric zero- temperature Glauber dynamics. The full spectrum of a modified Hamiltonian, whose minimum eigenvalue generates thelarge deviation function for the total particle current fluctuations through 4 a Legendre-Fenchel transformation, is obtained analytically. Three examples are presented and 1 0 numerically exact results are compared to ouranalytical calculations. 2 PACSnumbers: 05.40.-a,05.70.Ln,05.20.-y n Keywords: driven-diffusivesystemsinone-dimension,particlecurrentfluctuations,Glauber model a J 7 I. INTRODUCTION tracted from the first and the last sites of the lattice respectively. In the bulk of the lattice, on the other hand, ] theparticlesaresubjectedtototallyasymmetricbranch- h Mostofone-dimensionalnon-equilibriumsystemswith c stochastic dynamics show unique collective behaviors in ing and death processes. The steady-state of this sys- e temhasalreadybeencalculatedusingthematrixproduct their steady-states which usually can not be found in m method [14]. It is known that this steady-state can be their equilibrium counterparts. Non-equilibrium phase - written as a linear combination of product shock mea- transitionandshockformationaretwoexamplesofthese t a remarkable behaviors. Non-equilibrium systems are also sures with one shock front and that these shock fronts t have simple random walk dynamics [15]. s interesting to study from a mathematical point of view. . Wearespeciallyinterestedinthetotalparticlecurrent t These systems have opened up new horizon of research a fluctuations in this system. As we mentioned, the large in the field of exactly solvable systems [1, 2]. m deviation function for the particle current fluctuations Duringthelastdecadeseveralnon-equilibriumexactly - can be obtained from the Legendre-Fenchel transforma- d solvable systems have been introduced and studied in tion of the minimum eigenvalue of a modified Hamil- n related literature. On the other hand, different math- tonian. This modified Hamiltonian can be constructed o ematical techniques have been developed to study their c steady-state properties. A matrix method, knownas the from the stochastic time evolution operator of the sys- [ tem (sometimes called the Hamiltonian). This can be matrixproductmethod, isintroducedandusedtocalcu- done by multiplying the non-diagonal elements of the 1 late the steady-state and the average value of the phys- Hamiltonianofthesystembyanexponentialfactorwhich v ical quantities in the steady-state of these systems [3]. 5 The excitations, which give the relaxation times, can counts the particle jumps contributing to the total par- 7 also be obtained using the Bethe Ansatz [4]. Recently, ticle current, see for example [8–10]. It turns out that 2 the modified Hamiltonian associated with the total par- therehasbeenattemptstoestablishconnectionsbetween 1 ticle current fluctuations can be fully diagonalized. The Bethe Ansatz and matrix product method [5]. Other in- . 1 teresting quantities include the large deviation function key point is to change the basis ofthe vector space inan 0 appropriatewaybyintroducingaproductshockmeasure for the probability distribution of fluctuating quantities, 4 withmultipleshockfronts. Inthisnewbasisthemodified such as the particle current, in the steady-state of these 1 Hamiltonian becomes an upper block-bidiagonal matrix systems [6, 7]. The large deviation function can be ob- : v tained,throughaLegendre-Fencheltransformation,from whichis mucheasierto workwith, because we onlyneed i to diagonalize the diagonal blocks. X the minimum eigenvalue of a modified Hamiltonian. We Our analytical investigations reveal that for the small are then basically left with finding the minimum eigen- r a value of a matrix. The number of systems for which this particle current fluctuations (smaller than the average particle current in the steady-state to be more precise) quantity can be calculated exactly is very limited. the eigenvector associated with the minimum eigenvalue In this paper we consider a stochastic system of clas- of the modified Hamiltonian should be written as a lin- sical particles in which the particles interact with each earcombinationofproductshockmeasureswitha single other according to a variant of the zero-temperature shock front. In contrast, for the large particle current Glauberdynamicsonalatticewithopenboundaries[11– fluctuations (larger than the average particle current in 13]. More precisely the particles are injected and ex- the steady-state to be more precise) it should be written as a linear combination of product shock measures with more than one shock front. The validity of our analytical calculations is checked by comparing the analytical ∗Electronicaddress: [email protected] †Electronicaddress: [email protected] results with those obtained from numerical diagonaliza- ‡Electronicaddress: [email protected] tion of the modified Hamiltonian. 2 This paper is organized as follows. In the second sec- of ,A are given by {∅ } tionwewillreviewtheknownresultsonthe steady-state 0 0 ω 0 properties of the system. The total particle current is − 1 0 0 0 0 also introduced and its averagevalue in the steady-state h=  , 0 0 ω +ω 0 is calculated. In the third section we will briefly review 1 2 0 0 ω 0 the basics of the particle current fluctuations. The forth  − 2  section is devoted to the diagonalization of the modified α 0 Hamiltonian. The minimum eigenvalue of the modified h(1) = , (cid:18) α 0(cid:19) Hamiltonian will be discussed in the fifth section. We − will compare the analytical and numerical results in the 0 β sixth section. The concluding remarks are also given in h(L) = − . the last section. (cid:18)0 β (cid:19) TherighteigenvectorwithvanishingeigenvalueofHamil- tonian H gives the steady-state probability distribution II. STEADY-STATE vectorof the system. It is knownthat this vectorcan be written as a linear combination of product shock mea- sureswitha singleshock front[15]. It turns outthatthe Let’s consider a lattice of length L. We assume that dynamics of the position of a shock front is similar to each lattice site can be occupied by at most one particle thatofabiasedrandomwalkermovingonafinite lattice oravacancy. Thereactionrulesbetweentwoconsecutive withreflectingboundaries. This steady-stateprobability sites k and k+1 on the lattice are as follow distributionvectorhasalsobeenobtainedusingamatrix product method in [14]. By associating the operators E A with the rate ω ∅ −→ ∅∅ 1 (1) and D to the presence of a vacancy and a particle in a A AA with the rate ω . ∅ −→ 2 givenlattice site, the steady-state weight of any configuration τ , ,τ is proportionalto 1 L inwhichaparticle(vacancy)islabeledwithA( ). Apar- { ··· } ∅ ticle can enter the system from the left boundary of the L lattice with the rate α. A particle canalsoleavethe sys- W (τkD+(1 τk)E)V (4) hh | − | ii temfromtherightboundarywiththerateβ. Thismodel kY=1 is an asymmetric variant of zero-temperature Glauber in which τ = 0 if the lattice-site k is empty and τ = 1 k k dynamics [11–13]. The time evolution of the probabil- if it is occupied by a particle. In (4), V and W ity distribution vector P(t) is given by a master equa- are two auxiliary vectors. It has been sh|owinithat thhhese| | i tion [2] twooperatorsandvectorshaveatwo-dimensionalmatrix representation given by [14] d P(t) = H P(t) (2) 0 0 1 0 dt| i − | i D =(cid:18)d ω2 (cid:19), E =(cid:18) d 0(cid:19), ω1 − where the Hamiltonian H is an stochastic square matrix H = L−1 ⊗(k−1) h ⊗(L−k−1) |Vii=(cid:18) (ω2−−ω1β1+ωβ2)dω1 (cid:19), (5) k=1 I ⊗ ⊗I P (cid:0) (cid:1) + h(1) ⊗(L−1) (3) W = (ω1−ω2+α)d 1 ⊗I hh | (cid:16) α (cid:17) + ⊗(L−1) h(L) inwhichdisafreeparameter. Using(4)and(5)onecan I ⊗ easily calculate the weight of any configuration in the in which is a 2 2 identity matrix and that we have steady-state and also the average value of the physical I × defined quantities, such as the particle current, in the long-time limit. V⊗N V V V . Letuscall ρk (t)theaveragelocaldensityofparticles ≡ ⊗ ···⊗ at the latticehsitei k at time t. Using (1) and considering N times theinjectionandextractionofparticlesattheboundaries | {z } the time evolution of this quantity is given by Introducing the basis kets d ρ (t) = α 1 ρ ω ρ (1 ρ ) , 1 0 dth 1i h − 1i− 1h 1 − 2 i = , A = |∅i (cid:18)0(cid:19) | i (cid:18)1(cid:19) d ρ (t) = ω ρ (1 ρ ) ω ρ (1 ρ ) ,(6) k 2 k−1 k 1 k k+1 dth i h − i− h − i the matrix representation of h in the basis of d , A,A ,AA and that of h(1) and h(L) in the basis ρL (t) = ω2 ρL−1(1 ρL) β ρL . {∅∅ ∅ ∅ } dth i h − i− h i 3 in which k = 2, ,L 1. The average local density of and also the source terms which are defined as follows ··· − particlesisrelatedtotheaverageparticlecurrentthrough a continuity equation S = α 1 ρ +ω ρ (1 ρ ) , 1 1 2 1 2 h − i h − i d S = ω ρ (1 ρ ) ω ρ (1 ρ ) , (10) ρk (t)= Jk−1 (t) Jk (t)+Sk(t) (7) k 2h k − k+1 i− 1h k−1 − k i dth i h i −h i S = β ρ ω ρ (1 ρ ) . L L 1 L−1 L − h i− h − i for k = 1, ,L. We define J (t) as the average local k ··· h i particle current from the lattice site k to k+1 at time in which k = 2, ,L 1. The average particle current t. Sk(t) is also called a source term. In (7) we have also in(9)canbeund·e·r·stoo−dbyinvestigating(1). Thesecond assumed that J0 (t) = JL (t) = 0. In the steady-state dynamicalrulein(1)clearlyincreasetheparticlecurrent. h i h i the time dependency of the quantities will be dropped The first dynamical rule in (1) can be considered as a and we find backward movement of a vacancy. This is equivalent to a forwardmovementofa particle which, again,increases S = J J for k=1, ,L. (8) k h ki−h k−1i ··· theparticlecurrent. Similarexamplesofparticlecurrent in the presence of source terms can be found in [16]. Comparing (6) and (7) one finds the following relation for the average particle current in the steady-state The average local particle current in the steady-state canbecalculatedusingthematrixproductmethod. The J =(ω +ω ) ρ (1 ρ ) for k =1, ,L 1 (9) result is k 1 2 k k+1 h i h − i ··· − J = (ω +ω ) ρ (1 ρ ) k 1 2 k k+1 h i h − i W k−1(τ D+(1 τ )E)DE L (τ D+(1 τ )E)V = (ω +ω ) {τ}hh | i=1 i − i i=k+2 i − i | ii 1 2 P Q L Q W (τ D+(1 τ )E)V {τ}hh | i=1 i − i | ii W Ck−1DECPL−k−1 V Q = (ω +ω )hh | | ii 1 2 W CL V hh | | ii (ω2 ω2)αβ(ω2)k = 2 − 1 ω1 αω (β ω +ω )(ω2)L βω (α+ω ω ) 1 − 1 2 ω1 − 2 1− 2 in which we have defined C =D+E. Defining the total current J in the system is given by particle current as P(J) e−tI(J) (13) L−1 ≈ J = J (11) k h i h i Xk=1 which is valid for t . The large deviation function → ∞ I(J) measures the rate at which the total particle cur- it is easy to see that in the limit of L we have →∞ rent deviates from its average value. It is known that β(ω1+ω2) for ω <ω , the large deviation function I(J) is, according to the β−ω1+ω2 1 2 Gärtner-Ellis theorem, the Legendre-Fenchel transform  hJi= 2ω1 for ω1 =ω2, (12) oH˜f,tdheenomteindimbyumΛ∗e(iλg)e,n[v6a]lue of a modified Hamiltonian  αα(+ωω11+−ωω22) for ω1 >ω2. I(J)=mλax(Λ∗(λ)−Jλ). (14) This indicates that the systemundergoesa phasetransi- tion at ω = ω . The phase ω > ω (ω < ω ) is called 1 2 1 2 1 2 The modified Hamiltonian is defined as follows the low-density (high-density) phase. H˜ = L−1 ⊗(k−1) ˜h ⊗(L−k−1) k=1 I ⊗ ⊗I III. PARTICLE CURRENT FLUCTUATIONS P (cid:0) (cid:1) + h˜(1) ⊗(L−1) (15) ⊗I Assuming that the large deviation principle holds, the probability distribution for observing the total particle + ⊗(L−1) h˜(L) I ⊗ 4 ρ = 1 ρ = 0 01 k1 n1 k2 n2 k3 kf+1 L FIG. 1: A simple sketch of a base vector definedin (17). in which respect to λ at λ = 0 gives the total current defined in (11) 0 0 ω e−λ 0 1 − h˜ =0 0 0 0  , dΛ∗ 0 0 ω +ω 0 J = . (16) 0 0 1ω2e−λ2 0  h i dλ (cid:12)(cid:12)λ=0  −  (cid:12) In the following section we will show that H˜, defined h˜(1) = α 0 , in (15), can be diagonalized exactly. (cid:18) α 0(cid:19) − 0 β h˜(L) = − . IV. DIAGONALIZATION OF H˜ (cid:18)0 β (cid:19) NotethatsincethemodifiedHamiltoniandefinedin(15) In order to diagonalize the modified Hamiltonian H˜ becomes equal to the stochastic time evolution operator defined in (15) we start with redefining the basis of the given in (3) at λ = 0 then we have Λ∗(λ = 0) = 0. vector space by introducing the following product shock Finally, the first derivative of the minimum eigenvalue measure with N =2f +1 shock fronts k1,n1,k2,n2, ,nf,kf+1 = A ⊗k1 ⊗(n1−k1) A ⊗(k2−n1) A ⊗(kf+1−nf) ⊗(L−kf+1) (17) | ··· i | i ⊗|∅i ⊗| i ⊗···⊗| i ⊗|∅i in which N = 1,3,5, ,L+1 for an even L and N = the modified Hamiltonian (15) has the following upper ··· 1,3,5, ,L for an odd L and that 0 k < n < block-bidiagonalmatrix representation 1 1 ··· ≤ k < < n < k L. A simple sketch of such a2prod·u·c·t shofck mefa+su1re≤with multiple shock fronts is A1 B3 0 0 ··· 0 0 0 A B 0 0 0 given in Fig. 1.  3 5 ···  0 0 A B 0 0 5 7 For a given N the number of these vectors is sim- H˜ = ... ... ·.·..· ...  (20) ply given by a Binomial coefficient C (L +   L+1,N ≡  0 0 0 0 B 0  1)!/(N!(L+1 N)!). Nowthedimensionalityofthevec-  L−1  −  0 0 0 0 ··· A B  torspaceconstructedwiththese vectorscanbe obtained  ··· L−1 L+1   0 0 0 0 0 A  as  ··· L+1  inwhichA forN =1,3, ,L+1isaC C N L+1,N L+1,N ··· × matrix and B for N =3,5, ,L+1 is a C CL+1,N =2L (18) C matriNx. The matrix·e·l·ements of thesLe+m1,aNt−ri2ce×s XN arLe+g1i,vNen explicitly in the Appendix A. Although the modified Hamiltonian (20) is not a regardless of whether L is even or odd. The vec- stochasticmatrix;however,its matrixstructure suggests tors (17) make a complete orthonormalbasis for our 2L- the following picture which will be discussedin more de- dimensional vector space. Assuming L is an even num- tailintheforthcomingsections. Acting(20)on k1 with | i ber [20], in the basis k1 =0,1, ,L,whichisaproductshockmeasurewitha ··· single shock front, gives a linear combination of product shock measures with a single shock front. These series k , k ,n ,k , , k ,n , ,n ,k (19) of evolution equations are quite similar to the evolution 1 1 1 2 1 1 f f+1 {| i | i ··· | ··· i} 5 equations for a particle at the lattice site k1 performing A. Diagonalization of A1 a biased random walk on a one-dimensional lattice with dreeflfiencetianngibnovuanridaanrtiesse.ctIonr,owthheirchwowridlls,btehcealvleecdtor1s,{in|kt1hi}e inAth1eisbaasCisL(+117,)1,×thCeLf+ol1l,o1wtirnidgiamgaotnraixlmreaptrreisxenwthaitcihonhas, sensethatactingH˜ onanymemberofthissectSorgivesa linearcombinationofthevectorsinthesamesector. The α ω e−λ 0 0 0 0 1 − ··· matrix elements of A1 in (20) determine the coefficients  α ω1+ω2 ω1e−λ 0 0 0  − − ··· of these linear expansions. 0 ω e−λ ω +ω 0 0 0 2 1 2  − ···  On the other hand, acting (20) on k1,n1,k2 with A1 = ... ... ... ... ... ... ...  | i   m0 e≤asukr1e w<ithn1th<reek2sh≤ockLf,rownhtsic,hgiivsesaaplriondeuacrtcoshmobcik-  00 00 00 ··· ω1ω+e−ωλ2 −ωω+1e−ωλ 0β  nation of the product shock measures with a single or  0 0 0 ··· − 20 1ω e−λ2 −β  2 three shock fronts. These series of evolution equations  ··· −  are quite similar to those of two random walkers at the The structure of this matrix reminds us of the evolu- lattice sites k and k besides an obstacle at the lattice tion operator for a biased random walk moving on a 1 2 site n which does not have any dynamics. The reason one-dimensional lattice of length L + 1 with reflecting 1 that the shock front at the lattice site n (or the obsta- boundaries –although the reader should note that the 1 cle)doesnothaveanydynamicscanbeeasilyunderstood evolution operator is not a stochastic matrix. In an ap- by looking at (1). In fact, the position of a shock front propriatebasis 0 , 1 , , L theevolutionequations {| i | i ··· | i} of type 0 01 1 is not affected by these dynamical for the position of the random walker can be formally ··· ··· rules. As long as the random walkers are more than a written as follows single lattice site away from the obstacle, they perform A 0 = α1 +α0 , biasedrandomwalksonthelattice. Therandomwalkers 1| i − | i | i A k = ω e−λ k 1 ω e−λ k +1 +(ω +ω )k , also reflect from the boundaries of the lattice whenever 1| 1i − 1 | 1− i− 2 | 1 i 1 2 | 1i A L = β L 1 +β L they reach to the boundaries. When one of the random 1| i − | − i | i walkers arrives at an obstacle, the random walker and with k = 1, ,L 1. Now the eigenvectors and also 1 the obstacle both disappear; however, the other random ··· − the eigenvalues of A can be obtained using the same 1 walker continues to perform a biased random walk on approachemployed in [17] by writing the lattice. No new random walker or obstacle will be created once they disappear. The matrices A and B A Λ =Λ (λ)Λ (21) 3 3 1 1 1 1 | i | i are responsible for the dynamics of these random walkers. The above argumentsuggests that the vectors in and considering 1 S besides the vectors k ,n ,k define an invariant sec- {| 1 1 2i} L toofrt,hwishsiecchtowrilglivbeescaallilnedeaSr3c.omAbcitninagtioHñoofnthaenyvemcteomrsbienr |Λ1i= Ck1|k1i. (22) kX1=0 the same sector. Substituting (22) in (21) and using the evolution equa- The next invariant sector, which will be called , is 5 tions for the shock front one can calculate C ’s by ap- defined by the vectors in and k ,n ,k ,n ,kS in k1 S3 {| 1 1 2 2 3i} plying a plane wave Ansatz. Defining which 0 k < n <k <n < k L. In this case we 1 1 2 2 3 ≤ ≤ havethreerandomwalkerswhichareseparatedfromeach ω α β η 2, ζ 1 eλ, ξ 1 eλ other by two obstacles. Once a random walker meets an ≡rω ≡ − ω ≡ − ω 1 2 1 obstacle at two consecutive lattice sites, they disappear. The matrices A and B generates the dynamics of the and 5 5 randomwalkersandtheirinteractionswiththeobstacles. F(x,y,z) e−λ y+x+(x−1 y)z (y+y−1) ≡ − − This procedure can be continued to see that there are (cid:0) (cid:1) L/2+1 invariant sectors. we find tonIniaonr(d2e0r)twoeficnadntfhoelleoiwgetnwvoaleuqeusiovafltehnetmscoedniafireiods.HFarmoiml- Ck1 =ηk1(a1(z1)ζz)1kδ1k1+,0(a1(z1−1ξ))zδk1−1k,L1 (23) − − one hand, the eigenvectors of this matrix can be written in which as a linear combination of the vectors in each invariant sector. This helps us find all of its eigenvalues. On the a(z ) F(z ,η,ζ) F(z−1,η−1,ξ) other hand, since the eigenvalues of the modified Hamil- 1 = 1 = z−2L 1 . tonianareequaltothoseofA ’sforN =1,3, ,L+1, a(z1−1) −F(z1−1,η,ζ) − 1 F(z1,η−1,ξ) N ··· one can diagonalize each AN separately to calculate the It also turns out that the eigenvalues of A are given by 1 eigenvaluesof (20). We will employ the secondapproach which is the subject of the forthcoming sections. Λ (λ)=(ω +ω ) e−λ√ω ω (z +z−1). (24) 1 1 2 − 1 2 1 1 6 The equation governing z is also given by hand, for a given n , A(n1) is a n (L n )-dimensional 1 1 3 1 − 1 block tridiagonal matrix with the following structure F(z−1,η,ζ)F(z−1,η−1,ξ) z2L = 1 1 . (25) 1 F(z1,η,ζ)F(z1,η−1,ξ) A˜ B˜ 0 0 0 0 1 1 ··· Itcanbeseenthattheequation(25)has2L+4solutions.  B˜0 A˜2 B˜1 ··· 0 0 0  0 B˜ A˜ 0 0 0 Two of these solutions i.e. z = 1 have to be excluded 2 2 1  ···  sinceforthesevaluesofz1 theco±rrespondingeigenvector A3(n1) = ... ... ... ... ... ... ...  vanishes. Ontheotherhand,ifz1issolutionfor(25)then  0 0 0 A˜ B˜ 0  z−1 is also a solution. This means that the remaining  ··· 2 1  1  0 0 0 B˜ A˜ B˜  2L+2 solutions result in L+1 eigenvalues. For λ = 0  ··· 2 2 1   0 0 0 0 B˜ A˜  the pair z1 = η±1 corresponds to the eigenvalue Λ1 = inwhichA˜ andA˜ are(L n )··(·L n )2squa2rematrices 0. Finally, it can be shown that the solutions of the 1 2 − 1 × − 1 whose matrix representations are given in Appendix B. equation (25) are either phases i.e. z = 1 or they are real numbers. For the phase solution|s1|z1 = eiθ and the Onω teh−eλothanerdhBãnd=, weωhea−vλe dewfihneerdeB˜0is=a−(LαI, nB˜1)= smallest eigenvalue in this case is given by − 1 I 2 − 2 I I − 1 × (L n ) identity matrix. Noting that 1 − Λphase(λ)=(ω +ω ) 2e−λ√ω ω . (26) 1 1 2 − 1 2 L−1 The real solutions of the equation (25) are much easier n1(L n1)=CL+1,3 − to be found in the thermodynamic limit L . Let us nX1=1 →∞ restrict the real solutions to z > 1. It turns out that 1 | | the equation (25) has two real solutions in the thermo- it is easy to check the dimensionality of A . 3 dynamic limit Investigatingthe structureofA(n1) fora givenn sug- 3 1 gests that it can be regarded as a non-stochastic evolu- z(1) =G(η−1,ζ), z(2) =G(η,ξ) (27) 1 1 tion operator for two biased random walkers, moving on a one-dimensional lattice of length L+1 with reflecting where G(x,y) is defined as follows boundaries, which are separated by an obstacle. Let us denote the position of the first and the second random G(x,y) = 1 eλ(1+x2)+(y 1) 2x(cid:16) − walker on the lattice by k1 and k2 respectively. The obstacle is at the lattice site n . For a fixed n the random 1 1 + (eλ(1+x2)+(y 1))2 4x2y . walkers can only hop into the lattice sites which satisfy p − − (cid:17) the condition 0 k1 < n1 < k2 L. In terms of the matrix elements≤of A(n1), both ra≤ndom walkers reflect Substituting (27) in (24) gives the corresponding eigen- 3 valueswhichwillbedenotedbyΛ(1) andΛ(2). Whenever from the boundaries and also the obstacle. One should 1 1 note that the matrix A , in contrast to H˜, does not al- these eigenvalues exist, they will be definitely smaller 3 lowtherandomwalkersto mergewiththe obstacles. We than Λphase(λ). The conditions under which Λ(1) and 1 1 remind the reader that B’s in (20) were responsible for Λ(2) exist, will be discussed later. disappearance of the random walkers and the obstacles. 1 For a given n let us introduce an appropriate n (L 1 1 − n )-dimensional basis k ,k with 0 k n 1 1 1 2 1 1 B. Diagonalization of A3 and n + 1 k {|L. Wi}e arrange≤thes≤e vect−ors 1 2 ≤ ≤ as 0,n + 1 , 0,n + 2 , 0,L , 1,n + 1 , 1,n + 1 1 1 1 bloIcnkthdeiagboansiasl(m17a)trtihxewmitahtrLix A13 bisloackCs.L+T1h,3e×seCbLlo+c1k,3s 12,i,L·{·|·.|1I,nLtih,·is·i·b,|a|snis1−th1e,nev1io+·lu·1·tii|,o|nn1ei−qu|1a,tnio1n+s2fio,ri··t|·h,e|nra1n−- willbe calledA(n1) withn =1−, ,L 1. Onthe other domi}walkers can be written as follows 3 1 ··· − 7 A(n1) k ,k = ω e−λ k 1,k ω e−λ(1 δ )k +1,k 3 | 1 2i − 1 | 1− 2i− 2 − k1,n1−1 | 1 2i ω e−λ(1 δ )k ,k 1 ω e−λ k ,k +1 − 1 − n1+1,k2 | 1 2− i− 2 | 1 2 i + 2(ω +ω )k ,k for 1 k n 1,n +1 k L 1, 1 2 1 2 1 1 1 2 | i ≤ ≤ − ≤ ≤ − A(n1) 0,k = α(1 δ )1,k ω e−λ(1 δ )0,k 1 3 | 2i − − n1,1 | 2i− 1 − n1+1,k2 | 2− i ω e−λ 0,k +1 +(α+ω +ω )0,k for n +1 k L 1 , 2 2 1 2 2 1 2 − | i | i ≤ ≤ − A(n1) k ,L = β(1 δ )k ,L 1 ω e−λ(1 δ )k +1,L 3 | 1 i − − n1,L−1 | 1 − i− 2 − k1,n1−1 | 1 i ω e−λ k 1,L +(β+ω +ω )k ,L for 1 k n 1, (28) 1 1 1 2 1 1 1 − | − i | i ≤ ≤ − A(n1) 0,L = α(1 δ )1,L β(1 δ )0,L 1 +(α+β)0,L . 3 | i − − n1,1 | i− − n1,L−1 | − i | i These equations can be used to find the eigenvalues and (z ,z ) is a solution then the pairs (z−1,z−1), (z ,z−1) 1 2 1 2 1 2 eigenvectors of A3(n1). The reader can easily convince and (z1−1,z2) are also the solution, one finds n1(L−n1) himself that in the above mentioned basis, the matri- solutions (or eigenvalues) by mixing the solutions of the ces A˜ ’s are responsible for moving the position of the equations (32). 1,2 second random walker while B˜ ’s are responsible for In summary, for each n1 = 1, ,L 1 one solves 0,1,2 ··· − the equations (32) to find z and z . Substituting these moving the position of the first random walker. For a 1 2 into (31) gives the corresponding eigenvalues. The total given n , the eigenvalue equation 1 number of eigenvalues of A obtainedin this way will be 3 A(n1) Λ =Λ (λ)Λ (29) L(L2 1)/6. 3 | 3i 3 | 3i − can now be solved by using (28) and introducing C. Diagonalization of AN n1−1 L Λ = C k ,k . (30) | 3i k1,k2| 1 2i Inthe basis (17)the matrix AN is a CL+1,N CL+1,N kX1=0k2=Xn1+1 block diagonal matrix. Our procedure in the ×preceding sectionscanbecontinuedtoseethateachblockofA for By considering a plane wave ansatz and after some N N = 1,3, ,L+1 can be regarded as a non-stochastic straightforward calculations one finds that, for a given ··· evolution operator for f + 1 = (N + 1)/2 biased ran- n , the coefficients C are given by 1 k1,k2 dom walkers at the positions k = k ,k , ,k 1 2 f+1 { } { ··· } which are separated by f obstacles at the positions 2i=1ηki ai(zi)ziki +ai(zi−1)zi−ki n = n ,n , ,n given that 0 k < n < Ck1,k2 = Q (1(cid:16) ζ)δk1,0(1 ξ)δk2,L (cid:17) k{2}< n2{<1 2<··n· f <f}kf+1 L. In o≤ther1words,1for − − ··· ≤ i=1,2, ,f +1 the positions of the obstacles and the for 0 k1 n1 1 and n1+1 k2 L in which random·w·a·lkers should satisfy the following constraints ≤ ≤ − ≤ ≤ aa11((zz1−11)) =−FF((zz1−11,η,η,ζ,ζ)) =−z1−2n1 , ni−1+2≤ni ≤L−(N −2i), (33) a2(z2) = z−2LF(z2−1,η−1,ξ) = z−2n1 . ni−1+1≤ki ≤ni−1 a2(z2−1) − 2 F(z2,η−1,ξ) − 2 with n 1 and n L+1. For a given N each block 0 f ≡− ≡ The eigenvalues of A3(n1) are also given by ofAN,whichwillbecalledA(N{n}),isaD{n}-dimensional square matrix where D is given by {n} Λ (λ)=2(ω +ω ) e−λ√ω ω (z +z−1+z +z−1) (31) 3 1 2 − 1 2 1 1 2 2 D =n (n n 1) (n n 1)(L n ) {n} 1 2 1 f f−1 f − − ··· − − − in which the equations governing z and z are 1 2 with the following property z2n1 = F(z1−1,η,ζ) , 1 F(z1,η,ζ) D{n} =CL+1,N . (32) X{n} z2(L−n1) = F(z2−1,η−1,ξ) . 2 F(z2,η−1,ξ) InordertodiagonalizeA({n}) weconsideranappropriate N The firstequationin(32)has2n +2solutionswhile the D -dimensional basis k ,k , ,k and write 1 {n} 1 2 f+1 {| ··· i} secondequationhas2(L n )+2solutions. Excludingthe 1 solutionsz = 1andz−= 1andnotingthatifthepair A({n}) Λ =Λ (λ)Λ 1 ± 2 ± N | Ni N | Ni 8 in which the eigenvectors of AN are written as follows A. The case λ≥0 Λ = C k ,k , ,k . (34) The formula (35) suggests that as λ + the mini- | Ni X{k} k1,k2,···,kf+1| 1 2 ··· f+1i mum eigenvalue of H˜ should come from→the∞eigenvalues of A with the least number of random walkers, which N The coefficients can be calculated using a plane wave is in this case f = 0 i.e. the eigenvalues of A . Let us 1 ansatz and one finds work in the thermodynamic limit L . In this limit, →∞ A has two discrete eigenvalues Λ(1) and Λ(2) which can 1 1 1 fi=+11ηki ai(zi)ziki +ai(zi−1)zi−ki be calculated by substituting (27) in (24). These two Ck1,k2,···,kf+1 = Q (1 (cid:16)ζ)δk1,0(1 ξ)δkf+1,L (cid:17) λeigen0vathlueesmginoimtouzmereoigaesnλva→lue0.ofWtheehmavoedfiofiuenddHtahmatiltfoor- − − ≥ nian is either Λ(1) or Λ(2) of A depending on the values in which 1 1 1 of the microscopic reaction rates ω and ω and also on 1 2 λ. Defining a1(z1) = F(z1,η,ζ) = z−2n1 , a1(z1−1) −F(z1−1,η,ζ) − 1 (αω βω )2+(α β)2ω ω 1 2 1 2 λ ln[ − − ] aaii(z(zi−i1)) =−zi−2ni−1 =−zi−2ni fori=2,··· ,f , c ≡ (α−β)[(αω1−βω2)(ω1+ω2)+αβ(ω2−ω1)] we bring a summery of the results in the following: af+1(zf+1) = z−2LF(zf−+11,η−1,ξ) = z−2nf . af+1(zf−+11) − f+1F(zf+1,η−1,ξ) − f+1 • ωIn1t>hiωs2phase, for 0 λ λ , we have Λ∗ = Λ(1). ≤ ≤ c 1 The eigenvalues are also given by For λ λ the minimum eigenvalue can be deter- c ≥ minedusingtheleftpanelofFig.2. Ascanbeseen f+1 in the region III (the shaded area) the minimum ΛN(λ)=(f+1)(ω1+ω2)−e−λ√ω1ω2 (zi+zi−1) (35) eigenvalue is givenby Λ(12) while it is givenby Λ(11) Xi=1 intheregionsIandII.Thedifferencebetweenthese regionsisrelatedtothe asymptoticbehaviorofthe in which zi’s satisfy the following equations minimum eigenvalue Λ∗ when λ + which is → ∞ given by z2n1 = F(z1−1,η,ζ) , 1 F(z1,η,ζ) α in region I, z22(n2−n1) =z32(n3−n2) =···=zf2(nf−nf−1) =1, (36) λ→lim+∞Λ∗(λ)= βω1+ω2 iinn rreeggiioonn IIIII,. (37) zf2(+L1−nf) = FF((zzff−++111,,ηη−−11,,ξξ)) . • ωIn2t>hiωs1phase, for 0 λ λ , we have Λ∗ = Λ(2). ≤ ≤ c 1 For λ λ the minimum eigenvalue can be deter- In summary, the eigenvalues of AN can be calculated ≥ c mined using the right panel of Fig. 2. As can be as follows: we first fix the position of the obstacles n { } seenintheregionI(theshadedarea)theminimum which should satisfy the first relation in (33). We will then solve the equations (36) and substitute their solu- eigenvalue is givenby Λ(11) while it is givenby Λ(12) tions in (35) which gives the corresponding eigenvalues. inthe regionsIIandIII. Theasymptoticbehaviors For each set of n one find D eigenvalues. of the minimum eigenvalue Λ∗ in different regions {n} { } is given by (37). Now that allof the eigenvalues of the modified Hamil- tonian(15)areknown,wediscussaboutthesmallestone Usingtheabovedescriptionoftheminimumeigenvalue intheforthcomingsectionfromwhichthelargedeviation Λ∗ and (16) one can easily reproduce the results in (12). function for the total particle current can be calculated. It is worth mentioning here that the first derivative of the minimum eigenvalue Λ∗ is not continuous at λ . c This means that the large deviation function for the to- V. MINIMUM EIGENVALUE OF H˜ tal particle current fluctuations, which can be obtained using (14), is a linear function of J for J J J [6] a b ≤ ≤ Itisclearthattheminimumeigenvalueofthemodified I(J)=Λ(1)(λ ) λ J =Λ(2)(λ ) λ J (38) HamiltonianH˜ dependsonboththemicroscopicreaction 1 c − c 1 c − c rates and also λ. For each value of λ there are 2L eigen- where for ω >ω 1 2 values. Thissectionisdividedintotwoparts. Inthefirst part we will consider the case λ 0. The second part is dΛ(2) dΛ(1) ≥ J = 1 and J = 1 , devoted in the case λ≤0. a dλ (cid:12)λ=λc b dλ (cid:12)λ=λc (cid:12) (cid:12) (cid:12) (cid:12) 9 β β I II I II ∗ (1) ∗ (2) Λ =Λ1 Λ =Λ1 ) ) ω2 (1Λ1 ω2 (1Λ1 + + = = 1 1 ω ∗ ω ∗ Λ Λ ∗ (2) ∗ (2) Λ =Λ1 Λ =Λ1 III III ω1+ω2 α ω1+ω2 α FIG.2: TheminimumeigenvalueofthemodifiedHamiltonianH˜ forλ>λc isgivenbydifferentexpressionsindifferentregions (see text). and that for ω >ω is the number of random walkers which can be obtained 2 1 from (40). We are actually looking for the value of X dΛ(1) dΛ(2) for which (35) is minimum. One can easily see that as J = 1 and J = 1 . a dλ (cid:12)λ=λc b dλ (cid:12)λ=λc λ→−∞ the equations (36) become (cid:12) (cid:12) The minimum eig(cid:12)envalue of the modified(cid:12) Hamiltonian z12n1 ≃1, for λ 0 generates, using (14), the large deviation z2(n2−n1) = =z2(nX−1−nX−2) =1, functio≥n for the total particle current fluctuations for 2 ··· X−1 0 J J . z2(L−nX−1) 1 ≤ ≤h i X ≃ in which n ’s should satisfy (33). It turns out that these i equations generate the minimum eigenvalue of (35) pro- B. The case λ≤0 videdthatthedistributionoftheobstaclesonthe lattice is uniform i.e. For λ 0 the situation is quite different. Let us first ≤ L define Min(ΛN) as the smallest eigenvalue of AN. We n1 =n2 n1 = =L nX−1 have found that for a given finite L the minimum eigen- − ··· − ≃ X value of the modified Hamiltonian Λ∗ is given by Λ∗ = which means z = z = = z = eiθ with θ Xπ. Min(Λ2i−1) for λci ≤ λ ≤ λci−1 in which i = 1,2,··· ,X It is now clear1that2for v·e·r·y smaXll and negative≃valuLes bydefining λc0 =0andλcX =−∞. Here X is adiscrete of λ, (35) takes its minimum value provided that the parameter which maximizes expression Xπ X X Xπ Xcos( ). (39) (z +z−1)=2 cosθ =2Xcos( ) L i i L Xi=1 Xi=1 and that it can take one of the following values becomes maximum. 1,2,3, ,L+2 for an even L, We havenot been able to find exactanalyticalexpres- ··· 2 sionsforλ ’sfori=1,2, ,X 1;however,itispossi- X = (40) ci ··· −  1,2,3, ,L+1 for an odd L. ble tofindλci numericallybysolvingthe followingequa- ··· 2 tion  This means that in order to calculate the large devia- Min(Λ ) =Min(Λ ) . (41) 2i+1 2i−1 toifonlenfugtnhctLio,nwfeoronthlyenteoetdaltpoakrtniocwle tchuerrmenitniimnuamsyeisgteenm- (cid:12)(cid:12)λ=λci (cid:12)(cid:12)λ=λci Ourexactnumerical(cid:12)calculationsshowth(cid:12)atalthoughthe values of A ’s for i=1, ,X. This has been shown 2i−1 ··· minimum eigenvalue of the modified Hamiltonian Λ∗ is schematically in Fig. 3. The inset of this figure shows continuous at λ ’s, its first derivative is not continuous X as a function of L. For instance, for a system of ci at these points which results in length 21 L 23 we only need to know Min(Λ ) up 1 ≤ ≤ to Min(Λ ). This can be understood as follows. Let us 11 I(J) = Min(Λ ) λ J assume that as λ→−∞ the minimum eigenvalue of the 2i−1 (cid:12)(cid:12)λ=λci − ci modifiedHamiltonianisgivenbyMin(Λ2X−1)whichcan = Min(Λ )(cid:12) λ J be obtained from (35) by substituting N = 2X −1. X 2i+1 (cid:12)(cid:12)λ=λci − ci (cid:12) 10 -2 -1 L* 1 2 Λ Λ Λ Λ cX-1 c3 c2 c1 -2 Λ -4 L-1 L L L L 8 L2X LH5 LH3 LH1 LH1 76 HMin Min Min Min-6 Min X 45 3 2 1 1 6 1013 17 2124 28 -8 L FIG. 3: Schematic of the structure of the minimum eigenvalue of the modified Hamiltonian Λ∗ as a function of λ. The inset showsthemaximumnumberoftherandomwalkers X contributingintheminimumeigenvalueΛ∗ for agiven LuptoL=30. For more information see the text. for J J J where where z∗ is the real solution of (25) for λ = λ∗ which is 2i−1 ≤ ≤ 2i 1 givenbyG(η−1,ζ)λ=λ∗ andG(η,ξ)λ=λ∗ forω1 >ω2and d ω >ω respective|ly. Ontheotherh|and,asλapproaches J = Min(Λ ) , 2 1 2i−1 dλ 2i−1 (cid:12)(cid:12)λ=λci to λ∗ from below we have d (cid:12) J2i = dλMin(Λ2i+1)(cid:12)(cid:12)λ=λci . Min(Λ2X−1)(cid:12)(cid:12)λ=λ∗ ≃ (ω1+ω2)(cid:16)X (44) We havealso found that, for a give(cid:12)nfinite L, J2i J2i−1 (cid:12) 1 X (z∗+z∗−1) decreases as i increases. This will be discussed−in the −2Xi=1 i i (cid:17) next section in terms of three examples. Inwhatfollowswewillconsiderthelarge-Llimitwhich in which z∗’s for i = 1, ,X are the solutions of the seemstobe mucheasiertomanage. Itturnsoutthatfor equations (i36) at λ = λ∗··f·or f = X 1. We have nu- − L >> 1, λ λ for i = 2, ,X 1 drops to zero merically checked that for ω > ω the solutions of (36) | ci − ci−1| ··· − 1 2 as iL−2. On the other hand, it can be shown that in the satisfy thermodynamic limit L we have →∞ X λ∗ λ =λ = =λ ln2√ω1ω2 . z1∗ ≃G(η−1,ζ)|λ=λ∗ , (zi∗+zi∗−1)≃2(X−1) ≡ c1 c2 ··· cX−1 ≃ ω +ω Xi=2 1 2 while for ω >ω This means that in the large-L limit, we only need to 2 1 work with the minimum eigenvalues of A and A 1 2X−1 X−1 i.e. zX∗ ≃G(η,ξ)|λ=λ∗ , (zi∗+zi∗−1)≃2(X−1). Min(Λ ) for λ∗ λ + , Xi=1 1 Λ∗ = ≤ ≤ ∞ Replacingtheseinto(44)andcomparingitwith(43)con-  Min(Λ ) for λ λ∗ . 2X−1 −∞≤ ≤ firms (42).  For λ << λ∗ and in the large-L limit, Min(Λ ) is At λ=λ∗ we have found that 2X−1 approximately given by Min(Λ1) = Min(Λ2X−1) , Xπ d (cid:12)(cid:12)(cid:12)λ=λ∗ d (cid:12)(cid:12)(cid:12)λ=λ∗ Min(Λ2X−1)≃X(cid:16)(ω1+ω2)−2√ω1ω2cos( L )e−λ(cid:17). Min(Λ ) = Min(Λ ) (42) 1 2X−1 dλ (cid:12)λ=λ∗ dλ (cid:12)λ=λ∗ The minimum eigenvalue of the modified Hamiltonian (cid:12) (cid:12) (cid:12) (cid:12) for λ 0 generates, using (14), the large deviation func- which can be explained as follows. For L >> 1, as λ ≤ tionforthetotalparticlecurrentfluctuationsforJ J . approaches to λ∗ from above, we find using (24) that ≥h i In the next section we will check the validity of the 1 above mentioned results by studying three different ex- Min(Λ1)(cid:12)λ=λ∗ ≃(ω1+ω2)(1− 2(z1∗+z1∗−1)) (43) amples. (cid:12) (cid:12)