q-Zero Range has Random Walking Shocks

...but no other surprises in the zero range world. We check all nearest neighbour 1-dimensional asymmetric zero range processes for random walking product shock measures as demonstrated already for a few cases in the literature. We find the totally asymmetric version of the celebrated q-zero range process as the only new example besides an already known model of doubly infinite occupation numbers and exponentially increasing jump rates. We also examine the interaction of shocks, which appears somewhat more involved for q-zero range than in the already known cases.


Introduction
Since the early investigations of their hydrodynamic scaling limits (e.g., Rezakhanlou [27]), it has been clear that several 1-dimensional asymmetric interacting particle systems develop shocks; moving discontinuities of the rescaled density profile. Such shocks are very well understood on the hydrodynamic PDE level, but less so in the stochastic model itself. This microscopic phenomenon has been investigated by several authors, we give a brief summary of some of the important steps in the area.
Second class particles are coupling objects on the level of particles that trace the difference between two coupled instances of a particle model with slightly differing initial conditions. We describe the precise coupled dynamics for our models in Sect. 2. Briefly, second class particles use up the difference of jump rates between the coupled pair of models. Their connections to shocks have a long history of exploration, we only mention a few instances of this here. Derrida et al. [15] investigated stationary distributions in the asymmetric simple exclusion process (ASEP) as seen by a second class particle, and realised that their description greatly simplifies under certain condition on the jump of the densities in the shock. Namely, the stationary distribution becomes a product of Bernoulli distributions with different constant densities on the two sides of the second class particle. Balázs [3] reproduced this result in the exponential bricklayers process (EBLP). Next, Belitsky and Schütz [10] were able to transfer the result in a fixed frame of reference by observing that under the same condition on shock densities, the ASEP evolves a shock product distribution into translations of this same distribution, moreover the coefficients are the transition probabilities of a drifted simple random walk. We therefore call this phenomenon random walking shocks, and will give its precise formulation in Sect. 3.
Belitsky and Schütz [10] could also handle multiple shocks, exploring the interaction between them. Random walking shocks of ASEP interact like another ASEP, except the jump rates of these shock-walkers change with the rank as we consider them from left to right. This results in an attraction between the shock-random walkers, making them appear as one large shock after rescaling. Balázs [4] again reproduced this result in the EBLP and noted that, rather than excluding each other, in this model the random walkers almost act independently except for their rank-based jump rates. This again makes the walkers attract each other and a group of microscopic shocks appear as one large shock in the hydrodynamic scale. A precise statement for multiple random walking shocks follows in Sect. 3.
The random walking shock results had not used the second class particle until Balázs et al. [7] reproduced both random walking shock results (for ASEP and EBLP) with a second class particle added to the position where the densities jump. This on one hand simplified the arguments, on the other hand explained both the stationarity results from the viewpoint of the second class particle, and the random walking shock measures that first came without the second class particle. They also added a generalised exponential zero range process; a doubly infinite extension of zero range processes with an exponential form of the jump rates. This process came naturally as it has half the dynamics of EBLP. Certain other processes (e.g., [7,24]) are known to possess the random walking shock property, to the best knowledge of the authors none of which belong to the zero range family.
Stochastic duality, if discovered for a model, is a very powerful tool for the analysis of Markov processes. Given a process X (t) with state space Ω X , to establish stochastic duality one needs to find a dual process Y (t) on a state space Ω Y and a duality function Here E x refers to expectation w.r.t. the dynamics of process X (·) started from state X (0) = x, and the meaning of E y is analogous. Often the model Y (·) turns out to be simpler than X (·), in which case various properties of this latter can be deduced from such a duality. A general account on duality for interacting particle systems can be found in Liggett's book [22]. The function H and the process Y (·) are in general not easy to find. One of the early works discussing these in an asymmetric scenario is due to Schütz [28]. Such duality properties made an important role in the proofs of [10]. Recently, systematic work on duality for particle systems has been performed, we cite Groenevelt [17], Redig and Sau [25,26], and references therein. Multispecies models (i.e., with second class particles) also possess duality properties as demonstrated by Belitsky and Schütz [11] and later generalized by Kuan [20]. The importance of these developments for our context is highlighted in Belitsky and Schütz [12], where it becomes apparent that random walking shocks in ASEP, including the second class particle, are directly related to stochastic duality functions. The key observation is that certain duality functions can act as Radon-Nikodym derivatives between well-known stationary distributions and random walking shock measures. The dual process is (the timereversal of) the random walkers themselves.
In this work we perform a systematic search within 1-dimensional nearest neighbour attractive zero range processes and find the totally asymmetric q-zero range (q-TAZRP) as the only new model that has the random walking shocks property. This model is remarkable for its algebraic (and now, probabilistic) properties. We do not use dualities, however in light of recent developments our results give strong indications regarding the nature of dualities behind q-TAZRP. Indeed, some details already emerged in private communications with Frank Redig-further investigations will be part of future work.
We believe that our result and purely probabilistic approach will also nicely complement the more analytic research that has been done on q-TAZRP, some of which we briefly mention next.
Some of the early mentions of the model is by Povolotsky [23], where it is found as one suitable for Bethe ansatz. Remarkably, q-TAZRP is one of the examples that satisfy the microscopic concavity property, hence Kardar-Parisi-Zhang (KPZ)-scaling of current fluctuations could be proved using purely probabilistic methods in [8]. It also arises as the inter-particle distance process for q-TASEP, distortion of totally asymmetric simple exclusion and a very celebrated particle system in the integrable probability literature. Duality and determinantal formulas allow exact calculations in Borodin et al. [14] that show connections of q-TASEP to the celebrated KPZ equation. It has also been shown that q-TASEP arises as certain projection of McDonald processes, see Borodin and Corwin [13]. Bethe ansatz and integrable probability techniques allowed Korhonen and Lee [19] to derive exact formulas for certain transition probabilities in q-TAZRP. Ferrari and Vető [16] and Barraquand [9] proved Tracy-Widom limits for current fluctuations in q-TASEP, while recently Imamura and Sasamoto [18] proved KPZ-related scaling limits of the tagged particle distribution. Finally, we mention Barraquand [9] and Lee, Wang [21] as investigations for the non-homogeneous version of the dynamics.
In this article we stick to simple probabilistic methods, to be more precise we apply the arguments of [7] to carry out our search. For this reason at several points we simply refer to [7] rather than repeat all arguments from therein. Our notation fully conforms [7]. We start by recalling parts of the zero range framework in Sect. 2, then state and prove our results in Sects. 3 and 4.

Zero Range Process
We consider 1-dimensional asymmetric zero range processes with site occupation numbers ω i in the range either I = Z or I = Z + = {0, 1, 2, . . .}. Due to the zero range-type dynamics, these numbers cannot be bounded from above, while any finite lower bound can be transformed into 0 by a shift in the interpretation and the rate function of the jumps. Particles jump from site i with rate g(ω i ) choosing the right or left neighbouring site according to probabilities 0 < p ≤ 1 and 1 − p. The rate function g is non-negative and we assume it to be non-decreasing and not the constant function throughout Z. To preserve the state space, we also assume g(0) = 0 if I = Z + . The state space of the process is Ω = I Z , and the formal infinitesimal generator is Under appropriate initial data, the dynamics of this family has been constructed for g's of bounded increments by Andjel [1], and of the very similar bricklayers processes in the totally asymmetric case (i.e., p = 1) up to exponentially growing g's in Balázs et al. [6]. Non-decreasing g's make the model attractive, and it becomes a member of the family of misanthrope processes as described in [7] with f = g from therein. The product measure μ θ : = i∈Z μ θ with marginal on site i ∈ Z is stationary for θ < θ <θ with θ : = lim z→∞ ln(g(z)), and This is the range that makes the sum of μ θ (z) convergent. Here we used the usual definition Expectation w.r.t. this marginal will be denoted by E θ . We briefly comment on hydrodynamics of the model. Define the density of particles under the product measure μ θ as ρ(θ) : = E θ ω i . It is readily checked from the definition that this is a strictly increasing function; its inverse is denoted by θ(ρ). With the hydrodynamics flux, which is the expected jump rate in stationarity at density ρ: the density profile is expected to satisfy with the rescaled time and space variables T = εt and X = εi as ε → 0. This has in fact been proved for several cases of the rates g, see e.g. Rezakhanlou [27] or Bahadoran, Guiol, Ravishankar and Saada [2] for details. It has also been shown in [5] that H is convex, resp. concave, if g is so. In both cases, similarly to the inviscid Burgers equation, (3) develops discontinuous shock solutions. These are increasing jumps for concave H and decreasing jumps for convex H . Shocks of the PDE (3) have a well-defined velocity, described by the Rankine-Hugoniot formula for the case of densities λ and ρ on the two sides of the shock. The subject of this paper is the microscopic structure of shocks, and the way we handle this is via inserting second class particles as markers for the shock positions. These are objects arising from basic coupling of the models. We briefly summarize this classical object for zero range processes. We consider two realisations ω(t) and ζ (t) of the same zero range dynamics with initial states that satisfy ω i (0) ≤ ζ i (0) for each i ∈ Z. The formal infinitesimal generator for the joint evolution of this pair makes sure that -the order ω i (t) ≤ ζ i (t) is a.s. kept for any time t ≥ 0, -marginally ω(t) follows the zero range dynamics, -marginally ζ (t) also follows the zero range dynamics.
This generator on suitable functions ϕ of the pair (ω, ζ) is given by The rates in this generator are non-negative due to g being non-decreasing. We say that there are ζ i (t) − ω i (t) many second class particles at site i at time t. It is easy to see that the coupled dynamics conserves the number of second class particles.

Random Walking Shocks
We follow the ideas and methods in [7]; for completeness, we repeat the definitions in our context. For comparison, the rates from [7] now read p(y, z) = pg(y) and q(y, z) = (1 − p)g(z).
We start with a single shock, where we only have one second class particle in the system. The shock measure we consider in this case is a product (for sites) probability distribution on the pair (ω, ζ): where the indices represent sites of Z, ν θ is the stationary measure (2) concentrated on the diagonal of the single-site coupling state space I × I : ν is an arbitrary probability measure on I × I with exactly one second class particle: and θ , σ are the parameters of the shock measure on the left-and right-hand sides, respectively. This distribution describes a two-sided stationary distribution for sites i = j with constant densities and no second class particles, and one second class particle in a yet unknown distributionμ for ω j = ζ j − 1 at site j. We say that this distribution satisfies the random walking shock property, if its evolution according to the coupled dynamics (5) obeys with non-negative numbers P and Q that can depend on parameters θ and σ . This says that ν j evolves into a linear combination of translations of this same measure, and the coefficients in this linear combination agree with the transition probabilities of a simple random walk of right jump rate P and left jump rate Q.
Here is our first result of a systematic search within the class of zero range processes above for random walking shocks.

Theorem 1
The random walking shock property holds for zero range processes as described above in exactly the following cases.
We now comment on these cases.
1. This is the case of independent walkers. A second class particle in this model is yet another independent particle in the system with no interactions whatsoever, and it indeed sees a flat μ θ =μ = μ σ scenario stationary. No shocks occur in this description, P and Q are the constant rates of a single walker. 2. This is q-TAZRP, its rate can be rewritten as (a constant multiple of) g(y) = 1 − q y , and the shock condition is θ − σ = ln q < 0. The jump rate is concave, hence so is the hydrodynamic flux [5]. A shock therefore jumps upwards in density as we cross it from left to right. 3. This is the convex counterpart of q-TAZRP, only differs from the next case in having a lower bound on particle occupations. Shocks are negative jumps, accordingly. 4. This doubly infinite generalised exponential zero range has been known to have random walking shocks since [7], in fact its bricklayers counterpart, without the second class particles, was already investigated in [4]. A slight novelty here is the constant b which has not been noted before in this context. Indeed, without the constant b the marginals have the remarkable property that μ θ is a shift by one of μ σ . This property is lost with the introduction of b = 0 (it is obviously false in the above cases when the marginals are both concentrated on I = Z + ). Nevertheless, the random walking shocks structure survives this addition of b in the rates.
Next we discuss the case of multiple shocks, each with a second class particle. The shock measure now takes the form ν σ ,m : = i∈Z where the vector m is of non-negative integer components m i , and is a probability distribution on the one-site coupling space I × I . It forces exactly m i second class particles at site i, and features a yet unknown distributionμ σ i ,m i for ω i = ζ i − m i . We require that μ σ,0 = μ σ ; the distribution on a site without second class particles agrees with (2) of some parameter θ < σ i <θ . We say that ν σ ,m satisfies the random walking shocks property if, according to (5), its time evolution satisfies d dt with positive parameters P and Q that now can depend on the number m i and m i+1 of second class particles on the sites affected by the change, as well as the parameters σ i and σ i+1 on these sites. In m i,i+1 , the number of second class particles is changed analogously to (1), while σ i+1,i is left to be defined in the theorem below. We examine this phenomenon for cases 2, 3 and 4 from above.

Theorem 2 The random walking shocks property holds in each of the cases 2, 3 and 4 of Theorem 1 if
Then, Notice that all σ i 's must be bounded from above by ln a in case 2, and from below by ln b in case 4 for the distributions to make sense.
The interpretation is again that the product distribution (8) evolves into a linear combination of similar product measures with changed σ and m parameters, and the coefficients can be given an interacting random walks interpretation. To see this, notice that the steps m i,i+1 and σ i,i+1 keep (9) for all times. The walkers are the second class particles themselves, jumping from i to i + 1 with rate P(m i , m i+1 , σ i , σ i+1 ). This is clearly a generalization of the random walking shock property (7) for a single shock. The rates P give insights on the interaction between shocks. Notice first that the form we obtained for P is rather similar to the rates g of the model itself. This is something one has seen in [7] to some extent in both ASEP and EBLP, and it repeats here again. Multiple walkers are allowed to occupy the same site. In this case, we see stronger interactions than with EBLP before: except for b = 0 in case 4, the rate P of jump from a site i will not coincide with the sum of the rates of individual shock-walkers. This is due to the convexity of the exponential function.
To see that shock-walkers stay in tight distance to each other, we can repeat the argument for Lemma 4.2 of [4]. Namely, first label the walkers in increasing order from left to right (and run the label dynamics to keep this order). Then notice that two neighbouring walkers, if both are alone on their sites i < j, have higher rate to decrease their distance then to increase it. Adding higher-labeled walkers on site j blocks the front walker due to the labeling dynamics, while adding lower-labeled ones on site i makes m i , hence the rate of decreasing the distance, even higher. Thus, each inter-walker distance is stochastically bounded by independent Geometric distributions.
Thinking of hydrodynamics, it is natural to expect that our shock-walkers have the Rankine-Hugoniot velocity (4) as their mean velocity. We verify this directly in certain cases:

Proposition 3
The Rankine-Hugoniot formula (4) holds for the velocity in single shock cases of Theorem 1, and the two-shock cases of Theorem 2.
Due to the more complicated interaction of shocks, we do not see an easy way to verify this for more than two shocks, as has been done for the examples of [7].
There is a well established connection between exclusion and zero range processes that works via mapping inter-particle distances of exclusion into zero range models. The classical ASEP maps to constant rate asymmetric zero range, and q-TASEP maps to q-TAZRP this way. Notice that random walking shocks in the form discussed here are not known to appear in constant rate zero range, nor in q-TASEP. The mapping of course gives something, but this seems rather unnatural as the second class particles of the respective models do not get maped into each other either. The random walking shocks of ASEP, with [7] or without [10] the second class particle, give a strange statement for constant rate zero range that the authors can best imagine by simply undoing the mapping back to ASEP, which does not do much good to the state of the art. The random walking shocks discovered here for q-TAZRP also mean something in q-TASEP, but this also seems to be best imagined by undoing the mapping back to q-TAZRP. In fact, it is via the mapping that already flat, constant density stationary distributions of q-TASEP are best understood.

Proofs
We build on [7], and cite various steps and formulas from there across our proofs.

Proof of Theorem 1
We plug (6) in (49) and (50) of [7], with marginals (2) but no assumption onμ. Also, we rename x : = ω −1 , y : = ω 0 , z : = ω 1 from there. Then checking the existence of random walking shocks amounts to ensuring that are each constant as functions of y and x, resp. z. As p > 0 in all cases, the second line tells us thatμ = μ σ .
If p < 1 then this implies in the first line that g is an affine function, and since it cannot go negative, the range of occupations cannot be doubly infinite, that is I = Z + . Furthermore, this line also tells us that μ θ =μ, which in turn is equal to μ σ , and this gives case 1 of the theorem.
If p = 1 then the first line of the display is empty, and in the second linê The general solution is where the only restriction so far is that g(0) = 0 if I = Z + and g(y) ≥ 0 for all y ∈ I must hold. For I = Z this meansâ > 0 and b ≥ 0, and attractivity requires θ > σ as well. Moreover, θ > σ > θ = ln b. This is case 4. When I = Z + ,â > 0 requires again θ > σ, and we need b = −â for g(0) = 0. This is case 3. However, we can also haveâ < 0 with θ < σ <θ = ln b and b = −â. This is case 2 after a sign change in the constant, a = −â.
We next check that A of (52) in [7] is also independent of x, y and z. We only indicate how it is done in the linear case 1, and concentrate on the other cases. For these latter, + g(y + 1) · μ σ (y + 1) · e −σ · g(z) μ σ (y) − g(y + 1) + g(y) − g(z) = g(ω b+1 ) − g(ω a−1 ) + g(y) · e θ −σ − g(y) + g(z) − g(y + 1) + g(y) − g(z) = g(ω b+1 ) − g(ω a−1 ) + g(y) · e θ −σ − g(y + 1), For the linear case 1 the above covers the right jumps part of the dynamics. The exponential vanishes and the increment is constant. The left jump parts work similarly. For all other cases we now use (10): This being a constant, (10) with p = 1 does give models with random walking shocks. The left and right jump rates of the shock are, respectively (cf. [7]) This latter is further simplified by where we last used a general property of expectations within the whole family of misanthrope processes. Thus, which gives the rates as seen in the respective cases of the theorem. Notice that this also coincides with the negative of the expectation of (11) w.r.t. ω a−1 and ω b+1 .
Proof of Proposition 3 (single shock case) We start with the analogue of (12) for expectations. Notice that the function ρ(θ) can be given in a somewhat explicit expression using q-series or just summing up the below, but the next calculation is sufficient for our purposes. All shock cases of Theorem 1 are captured by (10), which we therefore use in the third step.
The change y = z + 1 of summation variable was possible due to the first term being zero if I = Z + , and summability when I = Z. With this result, the Rankine-Hugoniot velocity reads (13), as required.

Proof of Theorem 2
We need to check [7, (54), (55), (56)] with the assumptions we have. Equation (55) is trivially zero. Checking (56), we'll use ∓ as − for case 2 and + for cases 3 and 4 . With this notation we rewrite g(y) = ∓a · e ∓β y ∓ b which now conforms the statement of Theorem 1, with b = −a for cases 2 and 3 . We have, from [7, (56)], Repeating (12) in this scenario, This gives the P rates in the statement of the theorem, as soon as we can see that A i of [7, (54)] does not depend on ω a , . . . , ω b and its mean w.r.t. ω a−1 and ω b+1 exactly cancels the previous display under the summation of [7, (53)]. To this order, we proceed with calculating [7, (54)] for our setting.
Performing the sum in [7, (53)], the first two terms telescope: At this point we see that indeed the sum does not depend on ω a , . . . , ω b , as required. By assumption, see [7], all second class particles are in the volume a . . . b of dependence for the test function ϕ. Hence m a−1 = m b+1 = 0. By E θ g(ω) = e θ , we can calculate the mean All that is left is to check that b i=a−1 (C i +E i ) cancels with this. Notice that σ i ∓βm i = σ i−1 . Thus parts of C i + E i rewrite into a telescopic expression: Adding the two previous displays results in But this is zero due to m a−1 = m b+1 = 0 and (9).

Proof of Proposition 3 (two shocks)
The distance between the two walkers is a birth and death chain with rates, to be specified later, When the chain is at 0, the two walkers coincide, and we have a single rate with m = 2 for one of them to move. This coincides with r ↑ 0 . Otherwise both will have rates with m = 1 to move; the one behind has rate r ↓ while the one ahead has rate r ↑ + . The center of mass of the two walkers moves 1 2 units with these rates, hence has velocity The Rankine-Hugoniot velocity reads V = H ρ(θ + 2β) − H ρ(θ) ρ(θ + 2β) − ρ(θ) = e θ +2β − e θ ρ(θ + 2β) − ρ(θ + β) + ρ(θ + β) − ρ(θ) .
We now return to (14): e θ +2β b − e θ +β + e θ +β b − e θ +2β = C as required once b = a is substituted as was done for case 2. For case 3, change the sign of a, while for case 4, turn a into b. In both cases change the sign of β as well. In the last step for case 3, substitute b = −a as was done in that case. Otherwise the proof is word for word identical.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.