Irreversible Markov Dynamics and Hydrodynamics for KPZ States in the Stochastic Six Vertex Model

We introduce a family of Markov growth processes on discrete height functions defined on the 2-dimensional square lattice. Each height function corresponds to a configuration of the six vertex model on the infinite square lattice. We focus on the stochastic six vertex model corresponding to a particular two-parameter family of weights within the ferroelectric regime. It is believed (and partially proven, see Aggarwal, arXiv:2004.13272) that the stochastic six vertex model displays nontrivial pure (i.e., translation invariant and ergodic) Gibbs states of two types, KPZ and liquid. These phases have very different long-range correlation structure. The Markov processes we construct preserve the KPZ pure states in the full plane. We also show that the same processes put on the torus preserve arbitrary Gibbs measures for generic six vertex weights (not necessarily in the ferroelectric regime). Our dynamics arise naturally from the Yang-Baxter equation for the six vertex model via its bijectivisation, a technique first used in Bufetov-Petrov (arXiv:1712.04584). The dynamics we construct are irreversible; in particular, the height function has a nonzero average drift. In each KPZ pure state, we explicitly compute the average drift (also known as the current) as a function of the slope. We use this to analyze the hydrodynamics of a non-stationary version of our process acting on quarter plane stochastic six vertex configurations. The fixed-time limit shapes in the quarter plane model were obtained in Borodin-Corwin-Gorin (arXiv:1407.6729).


Overview
The study of models for surface growth is pervasive in many areas of science and engineering. Some physical examples include front propagation in combustion and crystal growth [DRS16].
A large collection of surface growth models falls under the umbrella of the KPZ universality class, named after the Kardar-Parisi-Zhang stochastic partial differential equation [KPZ86]. Models in the KPZ class are formulated as randomly growing height functions h t (x), where x ∈ R d is the space variable, and t ∈ R ≥0 is time. The models in the KPZ class have the following three characteristic properties [Cor12], [HT15] [QS15]: • Smoothing. The dynamics tends to force large fluctuations in the height function back towards the mean.
• Slope dependent growth speed. The average velocity of growth of the height function at a point only depends on its average slope around that point.
• Space-time uncorrelated noise. The randomness in the model comes from a random environment which is space-time uncorrelated, such as the white noise in the KPZ equation.
For models in the KPZ class, and for growth models in general, it is often intractable to obtain exact expressions for various observables. However, there exist integrable models in which algebraic structure allows for a precise analysis of observables. In this work, we construct a new integrable random growth model in two space dimensions, based on the six vertex model. The latter is well-studied in statistical mechanics by techniques of quantum integrability (in particular, Bethe Ansatz). We refer to the book [Bax07] for an introduction, and also to [Res10] for a more recent survey of the six vertex model. There are two basic questions in the study of a growth model: the identification of the translation invariant stationary measures (that is, measures which are invariant both under space translation and under the stochastic dynamics), and the computation of the current J(∇h). The current is the velocity of the height function growth as a function of a particular interface slope ∇h t (x), where ∇ is the gradient in space.
In the case of higher dimensions, the mathematical study of KPZ growth models is much more limited (for example, see [HT15] for a discussion of applicable numerical and renormalization group methods). In the (2+1)-dimensional situation (the case d = 2), the sign of the Hessian of the current distinguishes between two different subclasses, isotropic (positive Hessian) and anisotropic (negative Hessian), of the KPZ class, with very different large scale behavior. While in the isotropic case there exist only numerical estimates of the (nonzero) scaling exponents, for anisotropic KPZ models (AKPZ, for short) it is expected that α = β = 0. In particular, the average fluctuation of the height function grows at a rate slower than any polynomial in t. See also the recent results [CET20], [CET21] on scaling exponents and renormalization in the anisotropic KPZ partial differential equation.
The first (2+1)-dimensional anisotropic KPZ model which was studied rigorously is related to the dimer model on the hexagonal lattice (we refer to [Ken09] for details on this dimer model). This Markov dynamics was introduced in [BF14] in the non-stationary regime (that is, for a particular densely packed initial condition). After that, Toninelli [Ton17] showed that the Markov dynamics on Gibbs measures on the full plane is well-defined (that is, almost surely it does not make infinitely many jumps in finite time in a finite region) and preserves pure Gibbs states of an arbitrary slope (s, t). A pure state is, by definition, a translation invariant and ergodic Gibbs measure, and s and t may be interpreted as average densities of dimers of any two orientations (out of three orientations on the hexagonal lattice). The existence of the dynamics and the preservation of pure Gibbs states follows from a coupling of the dynamics on a large torus T L with that on Z 2 , and requires nontrivial probabilistic arguments.
An explicit formula for the current J(s, t) for this dynamics on Z 2 was conjectured in [BF14], [Ton17] based on particular cases, and then proven in [CF17]. The negative Hessian of J(s, t) puts the model into the anisotropic KPZ class. Under this dynamics, the height function fluctuations should grow at most as √ log t as t → +∞. This was shown in [BF14] for a particular initial condition, and in [Ton17] under some conditions on the slopes. In a subsequent work on the AKPZ class [CFT19], an analogous dynamics on dimer coverings of the square lattice is studied. The current for this process is computed explicitly, and it has negative Hessian. Moreover, [CFT19] presented a new simplified argument for the √ log t fluctuations for both square and hexagonal cases, which in particular removes the technical assumption of [Ton17].
It is known that the pure Gibbs states (with (s, t) away from a finite number of points in the polygon of allowed slopes) for the hexagonal dimer model are all liquid [KOS06], that is, the variance of the height function difference h t (x) − h t (y) grows logarithmically with distance |x − y| (with time t fixed and x, y ∈ R 2 ). Moreover, the limiting fluctuation field is the conformally invariant (massless) Gaussian Free Field [Ken08], [Pet15]. In other words, the Markov dynamics on the hexagonal dimer model indeed displays scaling exponents α = β = 0.
We remark that the celebrated domino shuffling of [EKLP92], [Pro03] also fits into the framework of [BF14] and thus belongs to the anisotropic KPZ class. We refer to [BF15], [CT19], and [BT18] for details.
The stochastic six vertex model introduced in [GS92] may be viewed as a certain model of interacting dimers. The presence of interaction introduces a new type of pure Gibbs states, the KPZ pure states π(s), arising for a special one-parameter family of slopes (s, t) with t = ϕ(s), where ϕ(s) = ϕ(s | u) := s s + u − su . (1.1) Here u ∈ (0, 1) is a parameter of the six vertex weights (i.e., it parametrizes the Gibbs property which gives rise to various pure states) which we will usually omit in the notation. The KPZ pure states exhibit scaling exponent α = 1 3 along a single direction in the plane, and 1 2 along all other directions [Agg18]. In particular, the limiting fluctuations cannot be described by a conformally invariant field. The name "KPZ state" comes from the fact that this stochastic six vertex model configuration in the plane may be viewed as a trajectory of a stationary Markov chain on particle configurations on Z (coming from the vertex model's transfer matrix which is a stochastic matrix) belonging to the (1+1)-dimensional KPZ universality class [BCG16], [CGST20].
It is believed (but not yet proven) that the stochastic six vertex model also displays liquid pure Gibbs states arising for generic slopes (s, t) outside a certain region H in the two-dimensional space of slopes (see Figure 6 for an illustration). The slopes on the boundary of H correspond to KPZ pure states, and there do not exist pure states with slopes inside H. This was conjectured in [SB94], [BS95], and proven recently in [Agg20b].
The main results of our paper are: • We construct an irreversible continuous time Markov dynamics C(t), t ∈ R ≥0 , on six vertex configurations in the full plane Z 2 . The dynamics preserves the KPZ pure Gibbs state π(s).
More precisely, we show that the full plane Markov dynamics is well-defined (i.e., almost surely does not make infinitely many jumps in finite time in a finite region) when started from π(s) (Theorem 5.7), and that π(s) is its stationary distribution (Theorem 5.8).
• We show that the current (that is, the average drift of the height function) under the dynamics C started from π(s) is given by an explicit formula (Theorem 7.1), (1.2) • In Section 8 we present an analogue of the continuous time dynamics C(t) which lives on six vertex configurations on a torus. This construction is independent of the full plane one. Moreover, we are able to extend our Markov processes to the case of the most general six vertex weights (traditionally denoted by a 1 , a 2 , b 1 , b 2 , c 1 , c 2 ), and for pure states on the torus with arbitrary average slopes (s, t) (in particular, away from the one-parameter family of KPZ pure states). Our torus dynamics is structurally similar to the one of [BB17], but they are different.
Let us make some general remarks on these results, and then in Section 1.2 below we formulate them in detail.
Algebraically, both the construction of the full plane dynamics and the dynamics on the torus come from the bijectivisation of the Yang-Baxter equation [BP19], which turns the equation into a Markov map. Its action maps the random six vertex configuration into a random six vertex configuration with swapped spectral parameters. Composing these Markov maps and passing to a Poisson type continuous time limit following the approach in [PS21] produces the jump rates in the full plane dynamics. In the torus case, we deform the six vertex graph on the torus by a certain twist to define the Markov maps, and then observe that in the continuous time limit the twist trivializes with probability going to one. We also present another proof of the fact that the dynamics on the torus preserves the Gibbs distributions using an adaptation of [BB17], which allows to generalize our construction to arbitrary six vertex weights and arbitrary average slopes.
For the full plane model, as in the dimer case considered in [Ton17], we employ a nontrivial probabilistic argument to prove that the Markov dynamics C(t) is well-defined. There are two major obstacles which we overcome compared to the dimer case. First, instead of coupling the full plane dynamics to the Hammersley process [Ham72], [AD95] to upper bound the probability of long jumps, we need to consider a new particle system in which particles can jump and annihilate one another. For this Annihilation-Jump process we are able to produce the required estimates. Furthermore, unlike in the dimer case, the stationary measure π(s) does not possess a determinantal structure, so we have to restrict our analysis to KPZ pure states admitting an explicit description [Agg18] as trajectories of a (1+1)-dimensional system, see Section 2.3 below.
We compute the current (1.2) for all KPZ pure states of the stochastic six vertex model. One in principle could employ our dynamics on the torus to get the current for general six vertex weights and arbitrary slopes, but this computation seems out of reach with existing techniques.
We believe that our Markov processes on the torus (extended to the general six vertex weights) should belong to the anisotropic KPZ class, at least for the weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 which are "sufficiently close" to the dimer situation (cf. fluctuation universality for small perturbations away from the dimer case in [GMT20]). The six vertex model corresponds to a dimer model for several choices of parameters. Let us mention two cases: • Setting a 1 = a 2 = b 1 = b 2 = 1 and c 1 = c 2 = √ 2 maps six vertex configurations into domino tilings [EKLP92], [ZJ00], [FS06a].
• Setting a 2 = 0 and b 1 b 2 = c 1 c 2 forbids the vertex (1, 1; 1, 1) and maps six vertex configurations into Gibbs ensembles of nonintersecting lattice paths. Here "Gibbs" means that the measure is invariant under uniform resampling. Such nonintersecting paths are readily identified with lozenge tilings, for example, see [BS10,Figures 2,3]. In this case we relate our processes on the torus to the ones studied in [BF14], [Ton17], [CFT19], see Section 8.3.
Both these degenerations to dimer models, however, cannot be obtained by directly specializing parameters in the stochastic six vertex model. For the latter, we know the current J(s, ϕ(s)) only along a one-dimensional curve of slopes, and so we cannot access the full Hessian of J(s, t) to determine its sign.
Above we have outlined and briefly discussed our main results. In Section 1.2 below we define the Markov dynamics C(t) on six vertex configurations and formulate the main theorems.
We define a Markov process whose state space is formed by collections of up-right lattice paths in Z 2 which can meet at a vertex but cannot cross or share an edge. Allowed path configurations are in bijection with height functions defined (up to a constant) on the faces of the lattice such that around each vertex the differences of the height functions are of one of six types displayed in Figure 1, top. An example of a state is given in Figure 1, bottom.
The Markov process C(t), t ∈ R ≥0 , is driven by a collection of independent Poisson clocks. Each vertical edge of Z 2 has a Poisson clock, and the rate of the clock at an edge (x, y)−(x, y +1), x, y ∈ Z, is a function of the local path configuration around that edge, as shown in the table in  1   1 δ  When the local configuration at some (x, y) ∈ Z 2 corresponds to an up jump (Figure 2, left) and the clock rings, then the path passing through the edge (x, y) − (x + 1, y) jumps up by 1. Instantaneously, each path through (x + i, y) − (x + i + 1, y) jumps up as well, for i = 1, 2, . . . , c, where c is the largest integer such that for each i = 1, 2, . . . , c, (x+i, y)−(x+i+1, y) is occupied and the horizontal edge immediately above is unoccupied. Then (x, y) − (x, y + 1) becomes occupied, and we also remove the vertical path from the vertical edge (x + c + 1, y) − (x + c + 1, y + 1) which was previously occupied. Similarly, when the local configuration at (x, y) corresponds to a down jump, then the path occupying the horizontal edge (x, y + 1) − (x + 1, y + 1) jumps down, and again a maximal sequence of adjacent occupied horizontal edges also jumps down. In words, both for up and down jumps we restore the path configuration to the allowed state using the most possible horizontal edges to the right of (x, y). See also Figure 11 in the text for an illustration of a jump.
The KPZ pure Gibbs state π(s) is a translation invariant ergodic Gibbs measure on full plane path configurations under which the average density of occupied vertical edges is s ∈ (0, 1), and the average density of occupied horizontal edges is ϕ(s).
Theorem (Theorems 5.7 and 5.8 in the text). Let 0 < s < 1. For a set Ω of initial states (Definition 5.2 below) which is probability 1 under the Gibbs measure π(s), there exists a Markov process C(t) whose generator acts according to the jump rates defined above. Furthermore, π(s) Up jump is stationary under C(t).
We can define the height change at a face under the Markov chain C(t) as follows. If a path jumps up past a face, then the height at this face height decreases by 1, and if it jumps down past a face, then the height increases by 1. Denote the resulting randomly evolving height function by h t (x, y). Define the current in the KPZ pure state π(s) by where the dynamics starts in stationarity from π(s). Due to stationarity, the ratio in the righthand side of (1.3) is independent of t, so J(s, ϕ(s)) is well-defined.
Theorem (Theorem 7.1 in the text). The current is given by formula (1.2).
Formula (1.2) for the current allows to predict a hydrodynamic limit equation (in two space dimensions) for the height function in a non-stationary version of our dynamics acting in the quadrant Z ≥0 × Z ≥1 . (In fact, we defined the full plane dynamics after looking at the bulk behavior of the quadrant dynamics.) Under this non-stationary dynamics which we denote by Q(τ ), the edge Poisson clocks have rates depending on the lattice coordinate y and time τ , but otherwise the dynamics is the same as C(t), see Section 4.3 below. The action of Q(τ ) changes (in distribution) the Gibbs property of the stochastic six vertex model by continuously increasing the spectral parameter from u to some terminal value u + η ∈ (0, 1), see Theorem 4.9 in the text. We match the current (1.2) to a (heuristic) hydrodynamic limit of Q(τ ). The latter should continuously transform the explicit limit shapes (obtained in [BCG16], see also [RS18], [Agg20a]) of the stochastic six vertex model in the quadrant with empty bottom and fully packed left boundary conditions. See Figure 14, left, for a simulation of the stochastic six vertex model with such boundary conditions.

Outline
In Section 2 we define the stochastic six vertex model, describe Gibbs pure states, and review the Yang-Baxter equation for the six stochastic six vertex model. In Section 3 we employ bijectivisation to turn the Yang-Baxter equation into a Markov map, and introduce the resulting discrete time Markov dynamics on six vertex configurations in the quadrant. In Section 4 we look at its continuous time Poisson type limit leading to a Markov dynamics Q(τ ) in the quadrant, and state the measure mapping property. In Section 5 we define the full plane dynamics C(t) and show that it exists and preserves the KPZ pure state π(s). The proof of the main estimate is postponed to Section 6, where we describe the coupling with the Annihilation-Jump particle system. In Section 7 we calculate the current and discuss the hydrodynamic limit of the dynamics Q(τ ) in the quadrant. In Section 8 we present the construction of the dynamics on the torus, first using bijectivisation, and then using an adaptation of the proof in [BB17] to generalize our Markov chain to general six vertex weights and arbitrary average slopes.
For an event or condition A, 1 A stands for the indicator of A.

Acknowledgments
We

Stochastic six vertex model
Here we define our main "static" objects -certain families of probability measures on configurations of up-right lattice paths on edges of Z 2 given by the stochastic six vertex model.

Vertex weights
Throughout the paper we fix the "quantization" parameter q ∈ [0, 1). The six vertex model is a probability measure on configurations of up-right paths on the twodimensional discrete lattice Z 2 , such that there is at most one path allowed per edge. The paths are allowed to meet at a vertex. See Figure 4 below for an example of a path configuration.
For a single vertex, let i 1 , j 1 , i 2 , j 2 ∈ {0, 1} be the number of paths entering and exiting the vertex as shown in Figure 3, left. Due to the path preservation, it must be i 1 + j 1 = i 2 + j 2 , which leaves six possibilities for a vertex. We assign the following weights to these six possible vertices: w u (0, 0; 0, 0) = w u (1, 1; 1, 1) = 1, (2.1) Here u ∈ (0, 1) is the spectral parameter (also called the rapidity) which may change from one vertex to another.
We will use the parametrization of the vertex weights by (q, u) as in (2.1) interchangeably with another one. This other parametrization involves two parameters (δ 1 , δ 2 ) with 0 ≤ δ 1 < δ 2 ≤ 1, see Figure 3 for an illustration. The two parametrizations are related by In particular, the ratio δ 1 /δ 2 is fixed throughout. The weights (2.1) satisfy the stochasticity property at each vertex, namely, Thus, our model is a particular case of the general asymmetric six vertex model with weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 specialized as a 1 = a 2 = 1, b 1 = 1 − c 1 = δ 1 , and b 2 = 1 − c 2 = δ 2 (cf. Figure 3). The stochastic case of the six vertex model was introduced by Gwa and Spohn [GS92] who studied its hydrodynamics in 1+1 dimension. Asymptotic Tracy-Widom GUE fluctuations in this model were obtained in Borodin-Corwin-Gorin [BCG16].

Gibbs property
Let us fix a finite rectangle Λ = {x, x + 1, . . . , x + R} × {y, y + 1, . . . , y + R } ⊂ Z 2 , and specify boundary conditions, that is, the positions of entering paths at the bottom and the left boundaries, as well as the positions of exiting paths at the upper and right boundaries (see Figure 4 for an example). Observe that the set of all up-right path configurations in Λ with specified incoming and outgoing boundary conditions is finite. We consider a probability measure on these path configurations depending on a sequence of spectral parameters u = (u y , u y+1 , . . . , u y+R ) (with u l ∈ (0, 1) for all l) given by . (2.4) The product in (2.4) is taken over all vertices (k, l) in the lattice, and for each vertex we pick one of the weights from (2.1) (with the corresponding spectral parameter u = u l ) depending on the occupation numbers i 1 , j 1 , i 2 , j 2 of the edges adjacent to that vertex. Here Z is the partition function, that is, the probability normalizing constant which ensures that the right-hand side of (2.4) sums to 1 over all possible path configurations with the given boundary conditions. One can also consider a probability measure P free u on path configurations in the rectangle Λ with free outgoing boundary conditions, for which the locations and numbers of outgoing paths on the left and top boundaries of the rectangle are not specified. Due to the stochasticity condition (2.3), the partition function of P free u is simply equal to 1, and the measure P free u can be sampled by running a row-to-row Markov chain based on the vertex weights w u , see Figure 5 for an illustration. The conditional distributions of P free u with fixed locations of all outgoing paths are precisely the measures P rect u . We use the measures P rect u (2.4) as a building block for random ensembles of up-right paths in infinite regions of Z 2 .
Definition 2.1 (Stochastic six vertex u-Gibbs measures). Let u = {u l ∈ (0, 1) : l ∈ Z} be a sequence of spectral parameters, and Λ ⊆ Z 2 a finite or infinite rectangular region. A probability measure P on configurations of up-right paths in Λ is called u-Gibbs if for any finite rectangle Λ = {x, . . . , x + R} × {y, . . . , y + R } ⊂ Λ, the conditional distribution of the up-right paths in Λ with arbitrary fixed incoming and outgoing boundary conditions on all four sides of Λ is P rect u| Λ . Here u| Λ := (u y , u y+1 , . . . , u y+R ) is the corresponding restriction of the spectral parameter sequence.
w uy w uy+1 w uy+2 Figure 5: Sampling of P free u by the row-to-row Markov chain. At each step, we perform the sequential update (from left to right) in a single row, using the vertex weights w u y+i , i = 0, 1, 2, and the incoming paths from the left and from the row below. For example, after one step the probability of getting the displayed configuration is equal to w uy (1, 0; 1, 0)w uy (1, 0; 0, 1)w uy (0, 1; 1, 0).
Remark 2.2. In Definition 2.1 it suffices to consider only the choices of boundary conditions for Λ having nonzero P-probability.
Let Λ be finite. In this case the u-Gibbs measures of Definition 2.1 on up-right path ensembles in Λ are mixtures of the measures P rect u corresponding to taking random boundary conditions on all four sides of the finite rectangle Λ. Since the set of all possible boundary conditions is finite, the mixture is also finite. The measure P free u is a particular example of a u-Gibbs measure with random boundary conditions (on the right and top boundaries only).
Definition 2.3 (Homogeneous stochastic six vertex Gibbs measures). The above Definition 2.1 deals with inhomogeneous (in the vertical direction) Gibbs measures for the stochastic six vertex model. Setting u l ≡ u ∈ (0, 1) for all l ∈ Z leads to the important subclass of homogeneous stochastic six vertex Gibbs measures. The homogeneous Gibbs property is indexed by two parameters (q, u) (equivalently, (δ 1 , δ 2 ), see (2.2)).
Remark 2.4. The space of homogeneous Gibbs measures coincides with the space of the Gibbs measures for the symmetric ferroelectric six vertex model. The symmetry means that a 1 = a 2 = a, b 1 = b 2 = b, and c 1 = c 2 = c in the notation as in Figure 3, and the ferroelectric condition reads ∆ = (a 2 + b 2 − c 2 )/(2ab) > 1. We refer to [Agg18, Appendix A.1] for details.
Let us define a family of u-Gibbs measures on up-right paths in the quadrant Z ≥0 ×Z ≥1 which are analogues of P free u and can also be sampled by running a row-to-row Markov chain: Definition 2.5 (Stochastic six vertex model in the quadrant). Fix a sequence of spectral parameters u l ∈ (0, 1), l = 1, 2, . . ., and take Λ to be the quadrant Z ≥0 × Z ≥1 . The stochastic six vertex in the quadrant with step (also called half domain wall ) boundary conditions [GS92], [BCG16] has empty bottom boundary and an incoming path at each horizontal edge (−1, l) − (0, l), l ≥ 1, along the left boundary. It is the unique u-Gibbs measure for which for any h ≥ 0, l ≥ 1 the marginal distribution of the configuration in a finite rectangle of the form {0, 1, . . . , h}×{1, . . . , l} is P free u (with empty and packed incoming conditions at the bottom, respectively, the left boundary). The uniqueness of thus described measure follows in a standard way from the Kolmogorov extension theorem, since the marginals P free u in finite rectangles are compatible. More generally, if the bottom boundary has paths incoming at locations encoded by a (finite or infinite) subset λ ⊆ Z ≥0 , and the left boundary is packed, then we call these the step-λ boundary conditions. If the left boundary is empty and the bottom boundary is encoded by λ, we refer to this as the empty-λ boundary conditions. We denote by P s6v u the distribution of the stochastic six vertex model in the quadrant (with step-λ or empty-λ boundary conditions, which is specified in each case separately).

Pure states
Here we discuss homogeneous Gibbs measures on path ensembles in the full plane satisfying certain natural assumptions.
Definition 2.6. Consider the homogeneous stochastic six vertex Gibbs property (Definition 2.3) with some parameters (q, u). A translation invariant, ergodic Gibbs measure (also called a pure state, for short) is a Gibbs probability measure on configurations of up-right paths in the whole plane Z 2 which satisfies two additional properties: • The distribution of the path ensemble does not change under shifts of the underlying lattice Z 2 by arbitrary elements of Z 2 ; • The measure is ergodic, which, by definition, means that the probability of any translation invariant event (from the σ-algebra associated with path configurations on Z 2 ) is either 0 or 1.
In particular, each pure state admits a slope (s, t) ∈ In other words, s and t are the average densities of the vertical and horizontal edges under the pure state. If a pure state has slope (s, t), we denote it by π s,t . Let us recall known results and predictions about pure states π s,t . We refer to [Agg20b] for details, and follow the notation of that paper. Let We usually omit the spectral parameter u in the notation. Define the following subsets of [0, 1] 2 : Let H be the open subset of [0, 1] 2 between h 1 and h 2 , so that ∂H = h 1 ∪ h 2 . See Figure 6 for an illustration. By [Agg20b, Theorem 1.2], for (s, t) ∈ H, there are no pure states with this slope, and for each (s, t) ∈ ∂H, there is exactly one pure state with this slope. Since the phase diagram in Figure 6 is symmetric in (s, t), it suffices to consider only π(s) := π s,ϕ(s) . For path configurations, the symmetry (s, t) ↔ (t, s) is realized by the reflection across the line y = x.
The measure π(s) is precisely the trajectory of the stationary stochastic six vertex model from [Agg18]. Namely, under π(s) the joint distribution of the locations of the occupied entering vertical and horizontal edges along the boundary of any quadrant {x, x + 1, . . . } × {y, y + 1, . . . } ⊂ Z 2 is given by the Bernoulli product measure with density s for the vertical and ϕ(s) for the horizontal edges, respectively (that is, each edge is occupied independently with probability s for vertical and ϕ(s) for horizontal entering edges). Then given a boundary condition, we use the stochastic weights depending on (q, u) to sample the configuration in the quadrant according to the measure P s6v u . In other words, the row-to-row transfer matrix in the stochastic six vertex model defines an interacting particle system in (1+1) dimensions (that is, a discrete time Markov process on configurations in Z). The path configuration in Z 2 under π(s) is the trajectory of the evolution of the Bernoulli measure of density s (on a horizontal slice) under this interacting particle system.
We refer to the part ∂H of the phase diagram as the KPZ phase of the stochastic six vertex model, since in this phase the fluctuations of the height function live on scale N 1/3 (in domains of size N ) along the characteristic direction, and are described by the Baik-Rains distribution [Agg18] (the distribution was introduced in [BR00]). This asymptotic fluctuation behavior is typical for stationary models in the KPZ universality class [Cor12].
We discussed pure states corresponding to slopes (s, t) from H := H ∪ ∂H. For (s, t) not on the boundary of [0, 1] 2 but outside of H, the pure states conjecturally exist and exhibit liquid phase behavior in the sense of [KOS06]. That is, the height fluctuations in a domain of size N grow logarithmically with N . The liquid phase behavior for the stochastic six vertex model is a major open problem. In the present paper we mostly deal with pure states in the KPZ phase, and it would be very interesting to extend at least some of our results (besides the existence of the irreversible dynamics on the torus) to the liquid phase.
(2.8) See Figure 7 for an illustration of the sums involved in both sides of (2.8).
Proof. There are finitely many choices of the boundary conditions i 1 , i 2 , i 3 , j 1 , j 2 , j 3 , and for each such choice the Yang-Baxter equation (2.8) (an identity between rational functions in u, v, and q) is verified in a straightforward way. Figure 7: An illustration of the Yang-Baxter equation. The straight vertex in the lighter shade has weight w u , and in the darker shade has weight w v . The cross vertex has weight X u,v .
Remark 2.8. One can readily check that for any choice of i 1 , i 2 , i 3 , j 1 , j 2 , j 3 the number of summands in each part of the Yang-Baxter equation (2.8) is at most two.

Row operators
Let V = C 2 with the canonical orthonormal basis e 0 , e 1 . We think that the space V is associated with the vertical direction at a vertex, and use the vertex weights w u (2.1) to define four operators A u , B u , C u , and D u in V . The operators A u and D u act diagonally: A u e 0 = w u (0, 0; 0, 0)e 0 , A u e 1 = w u (1, 0; 1, 0)e 1 D u e 0 = w u (0, 1; 0, 1)e 0 , D u e 1 = w u (1, 1; 1, 1)e 1 , and the other two operators act as B u e 0 = w u (0, 1; 1, 0)e 1 , B u e 1 = 0, C u e 1 = w u (1, 0; 0, 1)e 0 , C u e 0 = 0.
Stacking the vertices horizontally, one can define the action of the operators A u , B u , C u , D u in finite tensor powers V ⊗k . This can be done inductively: Here Multiplying two operators like A u A v (acting in V ) corresponds to stacking two vertices vertically, with the spectral parameters u and v in the lower and the upper vertex, respectively. Thanks to the Yang-Baxter equation (Proposition 2.7), the operators satisfy a number of quadratic relations. We do not need all of the relations, so let us list some of them which are used below in the paper: (2.13) By summing (2.10) and (2.11), we see that Note relations (2.9)-(2.13) hold not only when acting in V (in a single-vertex situation, which is Proposition 2.7), but also in each finite tensor power V ⊗k . This is due to the fact that we can iterate the Yang-Baxter equation and move the cross vertex horizontally through a row of adjacent vertices.

Row operators in the half-infinite tensor product
We need the half-infinite tensor product V [0,∞) of the spaces V with the fixed vector e 0 . By definition, the space V [0,∞) has the orthonormal basis indexed by finite subsets λ ⊂ Z ≥0 (recall notation from Section 1.4), where We see that all but finitely many of the m i 's are equal to 0. We do not need the Hilbert space completion of V [0,∞) since all our expressions involving V [0,∞) below are in terms of matrix elements. Let us define how the operators A u and B u act in Indeed, each matrix element A u e λ , e µ with (λ) = (µ) is the product of the vertex weights w u (2.1) over all vertices indexed by Z ≥0 , and in this product all but finitely many of the vertices have weight w u (0, 0; 0, 0) = 1. This implies that the action of A u (in terms of matrix elements) is well-defined, and similarly for B u e λ , e ν , where (ν) = (λ) + 1.
Note that the other two operators, C u and D u , do not act in the infinite tensor product V [0,∞) due to the presence of infinitely many vertex weights w u (0, 1; 0, 1) = 1 in the corresponding products for their matrix elements.
Remark 2.9. The operators in V [0,∞) allow to express probabilities in the stochastic six vertex model in the quadrant (Definition 2.5) as matrix elements. For example, for the step boundary conditions, the probability to observe a configuration λ of occupied vertical edges at height l ≥ 1 from the bottom is given by This probability is nonzero if and only if (λ) = l.

Bijectivisation of summation identities
We recall the notion of bijectivisation from Bufetov-Petrov [BP19, Section 2]. Suppose we have two disjoint finite sets A, B, so that each element x ∈ A ∪ B is equipped with a positive weight (3.1) Identity (3.1) defines probability distributions on A and B with probability weights proportional to {w(a)} a∈A and {w(b)} b∈B , respectively. A bijectivisation is a coupling between these two probability distributions, expressed via conditional probabilities: Definition 3.1 (Bijectivisation). A bijectivisation is a family of forward and backward transition probabilities • Sum to one property: (3.2) • Reversibility condition: For general A and B, a bijectivisation is not unique.
Remark 3.2. In the special case when |A| or |B| is equal to 1, there is at most one bijectivisation. More precisely, in this case the solution {p fwd (a → b), p bwd (b → a)} a∈A, b∈B to the family of linear equations (3.2)-(3.3) is unique. When, moreover, this solution is nonnegative, then it determines a bona fide bijectivisation (i.e., a stochastic Markov map).
The row label is the cross vertex state on the left hand side, and the column label is the cross vertex state on the right hand side. The cross vertex states determine i 1 , i 2 , k 1 , k 2 and j 1 , j 2 , k 1 , k 2 (recall the labels in Figure 7). Throughout most of the table, the remaining states of the edges are determined uniquely (otherwise I, J are incompatible). However, in four cells additional information is required, and it is provided in terms of the indicators 1 k 3 =0 and 1 k 3 =1 .

Bijectivisation of the Yang-Baxter equation
Here we apply the general definition of the bijectivisation from the previous Section 3.1 to the Yang-Baxter equation of Proposition 2.7.
Denote I = {i 1 , i 2 , i 3 } and J = {j 1 , j 2 , j 3 }. The combination of I and J encodes a choice of the boundary conditions in (2.8). For each such choice, let A I,J and B I,J index the nonzero summands in the left, respectively, the right-hand side of (2.8). One can check that for any I, J, we have one of the following three possibilities: • Either the sums in both sides of (2.8) are empty For example, this happens when i 1 +i 2 +i 3 = j 1 + j 2 + j 3 . Call this the incompatible case.
• Or one or both of |A I,J | and |B I,J | is equal to 1. Call this the one-to-two case.
• We have |A I,J | = |B I,J | = 2, but the Yang-Baxter equation (2.8) looks as w(a 1 ) + w(a 2 ) = w(b 1 ) + w(b 2 ), where w(a 1 ) = w(b 1 ) and w(a 2 ) = w(b 2 ). In other words, we can always match each of the two terms in the left-hand side to the corresponding term in the right-hand side which has the same weight. Call this the two-to-two case.
For the one-to-two case, there is at most one nonnegative bijectivisation of the Yang-Baxter equation (by Remark 3.2). In the two-to-two case, we make the natural choice and assign the bijectivisation to be deterministic. That is, in the notation of the two-to-two case, set p fwd (a 1 → b 1 ) = p fwd (a 2 → b 2 ) = 1, and all the other forward probabilities to zero (and similarly to the backward probabilities). In this way, for any choice of the boundary conditions I, J, we have outlined a single solution to the linear equations (3.2)-(3.3). Denote this solution by (when I and J are incompatible, by agreement, we set all these probabilities equal to 0). Here u, v is the spectral parameter of the lower (resp. the upper) vertex in the left-hand side of the Yang-Baxter equation (2.8). Denote Finally, conditions 0 < u < v < 1 together with 0 < q < 1 ensure that α, β, γ belong to (0, 1), which guarantees the nonnegativity of the forward probabilities in Figure 8. The fact that the backward probabilities are then also nonnegative follows from the reversibility.
Depending on the left boundary, the cross is fully empty or fully occupied. In both cases its weight is equal to 1, see (2.7) and (2.1).
• For each i = 0, 1, 2, . . . , drag the cross past the i th vertical edge in the lattice. This is a random operation involving the transition probabilities p fwd u,v from Figure 8. Since κ, µ, λ are finite subsets, this random procedure eventually becomes deterministic far to the right because all edges are eventually empty. This corresponds to the forward transition probabilities becoming trivial: The resulting state of the cross vertex is empty, and it has weight 1.
Denote the law of the resulting random ν by U u,v (µ → ν | κ, λ). We call U u,v the two-row Markov transition operator which takes µ to a random ν given the boundary conditions encoded by κ, λ.

Action on two-row Gibbs measures
Within the setting of the previous Section 3.3, let us show how the operator U u,v (Definition 3.4) swaps the spectral parameters u ↔ v in Gibbs measures on path configurations in the two-row lattice Z ≥0 × {1, 2}.
Assume that the old triple (κ, µ, λ) with (λ) = (µ) + 1 = (κ) + 2 has the (u, v)-Gibbs distribution. This means that the conditional distribution of µ is where Z = Z(u, v) is the normalizing constant (the other case (λ) = (µ) = (κ) is similar, but the operators B should be replaced with A). The Yang-Baxter equation (Proposition 2.7) implies that Z(u, v) is symmetric in u, v. Indeed, after summing over µ to get Z(u, v), we see that the symmetry follows from the commutation of the operators B u and B v , see (2.12)-(2.13). The Markov operator U u,v extends this symmetry to the level of probability distributions: The same identity also holds with the operators B replaced everywhere by A.
Proof. We prove the statement with B, the one with A is analogous. The action of U u,v starts by adding the full cross (having weight 1) on the left. Then it proceeds by dragging the cross through the lattice, and sampling µ given κ,μ, λ by a sequence of p fwd u,v steps. Thus, we can write the left-hand side of (3.6) as the product of local vertex weights and the local transition probabilities p fwd u,v . The dragging of the cross and the summation overμ in the left-hand side of (3.6) allows to repeatedly apply the local relation which follows from the properties (3.2)-(3.3) of the bijectivisation (recall that ) in the right-hand side of the local relation signifies the appearance of the local (v, u)-Gibbs distribution. The additional cross vertex weight X u,v (k 2 , k 1 ; j 2 , j 1 ) participates in the next local relation in the process of dragging the cross. Ultimately, one arrives at the configuration (κ, µ, λ) whose weight has swapped spectral parameters. The eventual state of the cross vertex is empty because it is to the right of the rightmost path. Since the weight of an empty cross is 1, we can remove the cross vertex weight, and we are left with the right hand side of the desired identity (3.6).

Action on infinite configurations
In the previous Sections 3.3 and 3.4 we defined the two-row Markov operator U u,v (μ → µ | κ, λ) acting on finite configurations of vertical arrows. Here we extend this operator to possibly infinite subsets of Z ≥0 . Recall that for µ ⊆ Z ≥0 and h ≥ 1, the truncation is Pick κ, λ ⊆ Z ≥0 , and assume that µ ⊆ Z ≥0 has a (u, v)-Gibbs distribution on the two-row lattice Z ≥0 × {1, 2}. That is, for any h ≥ 1 and any choice of the occupations of the horizontal edges (h − 1, 1) − (h, 1) and (h − 1, 2) − (h, 2) on the right boundary, the conditional distribution of µ [<h] is given by P rect (u,v) with the corresponding boundary conditions. For h ≥ 1, denote by U ) the probability that by randomly dragging the cross to the right past the horizontal coordinate h − 1 (using the probabilities p fwd u,v from Figure 8), the state of the internal edges at 0, 1, . . . , h − 1 is given by ν [<h] . By the same computation as in the proof of Proposition 3.5, the distribution of ν [<h] is (v, u)-Gibbs. Moreover, as h → +∞, the distributions of ν [<h] are compatible. Thus, by the Kolmogorov extension, we arrive at a random subset ν ⊆ Z ≥0 which has a (v, u)-Gibbs distribution.
Let us now take the stochastic six vertex model in the quadrant as in Definition 2.5, with the spectral parameters u = (u 1 , u 2 , . . .) and step-λ or empty-λ boundary conditions. The full configuration of the up-right paths in the quadrant can be encoded by a sequence λ (i) , i ≥ 0, of subsets of Z ≥0 , where λ (0) = λ (the bottom boundary condition), and each λ (i) represents which vertical edges among Definition 3.6. Fix k ≥ 1 and assume that u k < u k+1 . Let L k,u be the Markov transition operator acting on the sequence (λ (i) ) i≥0 by randomly changing λ (k) to ν (k) according to the two-row transition probability Proposition 3.7. We have where P s6v u L k,u is the action of a Markov operator on a probability measure, s k is the k-th elementary permutation, and the boundary conditions of both stochastic six vertex models in the quadrant in (3.7) are the same.
, and ν = ν (k) . The subsets κ, µ, ρ encode the configuration of the stochastic six vertex model before the application of L k,u , and ν is the result of this application.
Under the step-λ boundary conditions, for any truncation h ≥ 1 the conditional distribution of µ [<h] given For the empty-λ boundary conditions, both operators B + D should be replaced by A + C with the same spectral parameters, and the argument is analogous. Applying the operator U [<h] u,v which maps µ [<h] to a random subset ν [<h] , and arguing as in the proof of Proposition 3.5, we see that ). Indeed, this follows from the commutation of B u + D u with B v + D v , which in turn is a consequence of the Yang-Baxter equation, see (2.12)-(2.13). This completes the proof.

Two-row Markov operator via coin flips
Here we present another description of the Markov operator U u,v which is adapted to taking the continuous time limit in Section 4 below. Fix spectral parameters 0 < u < v < 1, and recall that we are dragging a cross through the two-row lattice Z ≥0 × {1, 2}, with spectral parameters u and v on rows 1 and 2, respectively.
Let us argue in the setting of finite subsets. Assume that κ, µ, λ with either (λ) = (µ) + 1 = (κ) + 2 or (λ) = (µ) = (κ) encode the configurations of occupied vertical edges as in Figure 9. We aim to present an algorithm sampling ν from the distribution U u,v (µ → ν | κ, λ). This algorithm involves independent random coin flips to initiate jumps at each horizontal coordinate i. If a jump is initiated, then it propagates to the right according to certain rules. A jump can be either up or down, let us describe them.
Denote the horizontal edges (i, 1) − (i + 1, 1) and (i, 2) − (i + 1, 2) by e i,1 and e i,2 , respectively. The occupation of these edges changes after initiating a jump (and the vertical edge's occupation changes accordingly). Note that the occupation of the previous horizontal edges e i−1,1 , e i−1,2 does not change after a jump at i is initiated.
Initiating up jump. An up jump can be initiated at i if the local path configuration is one of , , . (3.8) After initiating an up jump, the path occupying e i,1 jumps up from e i,1 to e i,2 , and the vertical edge (i, 1) − (i, 2) becomes occupied. For example, in the first case in (3.8) we have → .
Propagation of up jump. The up jump propagates as follows. Find the largest c = c(i) = c(i, κ, λ, µ) ≥ 0 such that e i+j,1 is occupied and e i+j,2 is unoccupied for all j = 0, 1, . . . , c, and for all these edges we swap the occupation status of e i+j,1 and e i+j,2 . We also change the occupation status of (i + c + 1, 1) − (i + c + 1, 2) from occupied to unoccupied. This results in a well defined configuration. See Figure 11 for an illustration.
Initiating down jump. A down jump can be initiated at i if the local configuration is one of , , .
Then the path occupying e i,2 jumps down from e i,2 to e i,1 , and the vertical edge (i, 1) − (i, 2) becomes unoccupied. For example, → .
Propagation of down jump. The down jump propagates as follows. Find the largest c = c(i) = c(i, κ, λ, µ) ≥ 0 such that e i+j,1 is unoccupied and e i+j,2 is occupied for all j = 0, 1, . . . , c. For all these edges, swap the occupation status of e i+j,1 and e i+j,2 . Change the state of the last vertical edge (i + c + 1, 1) − (i + c + 1, 2) from unoccupied to occupied.
Sampling procedure. Having defined initiation and propagation of the up and down jumps, we are now in a position to describe how the random finite subset ν with distribution U u,v is sampled. Start by setting i = 0. While i ≤ λ 1 (where λ 1 is the position of the rightmost occupied vertical edge), sequentially repeat the following steps: 1. If we can initiate an up jump at i, do so according to an independent coin flip with probability 1 − γ, 1 − α, or 1 − β (see (3.5)) depending on the local configuration as in the table in Figure 10, left.
2. If we can initiate a down jump at i, do so according to an independent coin flip with probability 1 − γ, 1 − α, or 1 − β (see (3.5)) depending on the local configuration as in the table in Figure 10, right.
3. If we initiated an up or down jump, let c(i) be as described above, and propagate the jump to the configuration at positions i + 1, . . . , i + c(i). Then set i = i + c(i) + 1.
4. If we cannot initiate any jump at i, set i = i + 1.
Proof. From the table of the forward transition probabilities p fwd u,v in Figure 8 we see that the probability to change the occupation state of a vertical edge by dragging the cross is equal to 1 − α, 1 − β, or 1 − γ depending on the state of the horizontal and cross edges around. Indeed, the state of the cross before and after the dragging is determined uniquely by the initiated jump. The jump is not initiated with the complementary probability α, β, or γ, respectively, which corresponds to another state of the cross vertex after the dragging. Next, the rules for propagation of up or down jumps come from the parts of the table in Figure 8 where p fwd u,v = 1. This completes the proof. * Figure 11: A generic up jump with propagation. In this case an up jump was initiated at the vertical edge with an asterisk. Note that the propagation can continue past vertical paths which go straight up.
Let us extend Proposition 3.7 to infinite subsets κ, µ, λ (the setup of Section 3.5). For each h ≥ 1, we may truncate the sampling algorithm by allowing the initiation of jumps only at i = 0, 1, . . . , h−1. Reading the resulting vertical path configuration after all the jumps, we clearly ). For h → +∞, these truncations of the sampling algorithm are consistent. Therefore, Proposition 3.7 holds when the Markov operator L k,u is defined as above in this subsection using jumps initiation and propagation.
Remark 3.9. When u = v, from (3.5) we see that α = β = γ = 1. This means that no up or down jumps can be initiated, so the Markov map U u,u is simply the identity operator. In this section we construct a continuous time Markov chain on the homogeneous stochastic six vertex model configurations in the quadrant Z ≥0 × Z ≥1 . The continuous time Markov chain is a Poisson type limit of the discrete time Markov chain obtain by a repeated application of the operators L k,u (Definition 3.6), as the inhomogeneous parameters u i become equal. The Taylor expansion of the transition probabilities around the identity (cf. Remark 3.9) leads to the desired continuous time dynamics. This Poisson type limit is analogous to the one in [PS21]. In particular, we also get space-inhomogeneous jump rates linearly depending on the y coordinate.
Definition 4.1. Denote by L u the one-step Markov operator which is the result of the application of the infinite sequence of Markov operators L 1,u , L 2,s 1 u , L 3,s 2 s 1 u , . . . (in this order). Here each s k is the elementary permutation (k, k + 1), the first operator L 1,u involves the spectral parameters u 1 , u 2 , the next operator L 2,s 1 u involves u 1 , u 3 , and so on.
Let (λ (y) ) y≥0 be the random sequence encoding the result of the application of L u to the sequence (λ (y) ) y≥0 we started with. The new random sequence is well-defined since for any finite y 0 , the first several layers (λ (y) ) 0≤y≤y 0 are obtained from (λ (y) ) y≥0 by finitely many Markov operators. Proposition 4.2. The operator L u acts on the measure P u s6v as P s6v u L u = P s6v Su .
See Figure 12 for an illustration.
Proof of Proposition 4.2. Immediately follows by iterating Proposition 3.7.

Poisson type limit and jump rates
Here we employ the description of the Markov operator L k,u in terms of independent coin flips (Proposition 3.8) to obtain a Poisson type limit of the transition probabilities.
We utilize the expansion from Lemma 4.3 together with the iterated moving of the bottom spectral parameter up to infinity, as defined in the previous Section 4.1. In this subsection we outline the main expansions, and in the next Section 4.2 we define the generator of the continuous time dynamics, and show the existence of the dynamics. Define (The order of the composition of Markov operators means that L u is applied first.) By Proposition 4.2, we have P s6v u Q m = P s6v S m u . Let us take spectral parameters close to each other. Fix real parameters u and η > 0 such that u, u + η ∈ (0, 1). Let ε > 0 be sufficiently small, and define Denote by τ ∈ R ≥0 the rescaled time. Then we have As ε → 0, the parameters u i (4.2) become all equal to u, and u[τ ] i become all equal to u + (1 − e −τ )η. The difference of the spectral parameters u[τ ] 1 and u[τ ] k+1 (which are exchanged at horizontal layer k at the step L S τ /ε −1 u in the chain Q τ /ε ) has the form Therefore, the ε → 0 limit of the Markov transition operators Q τ /ε should lead to a continuous time Markov chain with the transition semigroup (Q(τ )) τ ∈R ≥0 , which acts on the homogeneous stochastic six vertex model in the quadrant as P s6v, hom u Q(τ ) = P s6v, hom u+(1−e −τ )η . We see that from τ = 0 to τ = +∞, the chain Q(τ ) continuously increases the spectral parameter u to u + η. The definition of Q(τ ) employs the probabilities a(u), b(u), c(u) (4.1), and is given in the next Section 4.3.

Continuous time chain in the quadrant
Let us define the continuous time Markov semigroup Q(τ ) in terms of its generator G quad u,η,τ . The generator depends on u, η, and also on the time variable τ . The latter means that the continuous time Markov chain is time-inhomogeneous. First, recall a basic definition: For short, we say that the arrivals τ i in the Poisson process occur according to an exponential clock with time-dependent rate r(τ ).
Definition 4.5. Let v 1 , v 2 be a pair of vertices at vertically adjacent positions (j, k), (j, k + 1) in the quadrant. If the local configuration of the paths around v 1 , v 2 is one of the six configurations in Figure 10, we call the pair (v 1 , v 2 ) a seed pair.
We attach to each seed pair an independent exponential clock with the time-dependent rate where k is the y-coordinate of v 1 , and R v 1 ,v 2 (u) is given in Figure 13. The rate (4.4) is the coefficient by ε in the expansion of 1 − α, 1 − β, or 1 − γ as in Lemma 4.3, where we took into account the inhomogeneity coming from exchanging the spectral parameters, see (4.3).
When the clock at a seed pair (v 1 , v 2 ) rings, this generates an up or down jump of the horizontal path, as illustrated in Figure 10. This jump then instantaneously (at the same time moment, without any waiting) propagates to the right according to the rules given in Section 3.6. Thus defined jumps lead to the following infinitesimal generator of a time-inhomogeneous continuous time Markov chain: Definition 4.6 (Generator). Let σ denote a configuration of up-right paths in the quadrant. If (v 1 , v 2 ) is a seed pair for σ, denote by σ (v 1 ,v 2 ) the result of initiating an (up or down) jump of the horizontal path at (v 1 , v 2 ), and the propagation of this jump according to the rules in Section 3.6. Let f (σ) be a cylindric function. That is, f depends on σ only through the restriction of σ to a finite window inside the quadrant (this window depends on f ). The generator of the dynamics on the stochastic six vertex model in the quadrant is, by definition, the operator acting as Here the sum is over all seed pairs (v 1 , v 2 ) of σ. While the number of seed pairs may be infinite, the action of G quad u,η,τ (4.5) is well-defined on cylindric functions.
We aim to define a Markov semigroup Q(τ ) with generator (4.5) which can start from configurations belonging to a certain space of regular initial configurations Ω quad . This makes sure that Q(τ ) does not make infinitely many jumps through a finite space in finite time. Let us define the space of initial configurations, and then prove that the semigroup Q(τ ) and the corresponding Markov process σ(τ ) exists.  1; 1, 1)).
Proof. In the step-λ case, the configuration in the region [0, R ]×[R, ∞) is stochastically monotone in λ. That is, when one adds an extra occupied initial edge to λ, the probability that the configuration in [0, R ] × [R, ∞) is fully packed does not decrease. This follows by considering the stochastic six vertex weights (Figure 3), and observing that for fixed j 1 , the probability that i 2 = 1 increases when i 1 = 1. When λ = ∅, the probability of Ω quad is equal to 1 thanks to the Law of Large Numbers established in [BCG16]. Indeed, the latter states that the bottom boundary of the fully occupied region in the quadrant is linear, see Figure 14, left, for a simulation.
Similarly, for empty-λ boundary conditions, the probability that [0, R ] × [R, ∞) is empty of paths does not increase when adding extra occupied edges to λ. When λ = Z ≥0 , the probability of Ω quad is also equal to 1. Indeed, since δ 1 (the probability of going up) is smaller than δ 2 , the whole region which is slightly above the diagonal of the quadrant is empty, see Figure 14, right, for a simulation.
Theorem 4.9. Assume that σ(0) ∈ Ω quad . There exists a continuous time Markov chain σ(τ ), τ ∈ R ≥0 , on configurations of up-right paths whose generator at time τ is G quad u,η,τ given by (4.5) (here we use the convention around time-dependent rates, cf. Definition 4.4). Moreover, the transition operator Q(τ ) of this Markov chain acts on the homogeneous six vertex model (with step-λ or empty-λ boundary conditions for arbitrary fixed λ ⊂ Z ≥0 ) as follows: Proof. Recall the truncation of the bijectivisation Markov operator U to U [<h] , see Section 3.5.
Since the jump rates in our generator G quad u,η,τ are Poisson limits of the ones in the operators U , the generator is compatible with these truncations, too. That is, G quad u,η,τ preserves the space of cylindric functions f (σ) which depend on σ only through the restriction of σ to a finite vertical strip {0, 1, . . . , R } × Z ≥1 (where R > 0 is fixed). In other words, due to the very construction of the jump rates in G quad u,η,τ , the dependence does not propagate from {R + 1, R + 2, . . .} × Z ≥1 into {0, 1, . . . , R } × Z ≥1 . Therefore, it suffices to show existence of the process σ(τ ) restricted to an arbitrary vertical strip {0, 1, . . . , R } × Z ≥1 , and the full quadrant process then arises by Kolmogorov extension.
Consider the following collection of independent time-inhomogeneous Poisson processes: Here (j, k) are the lattice coordinates, and the Poisson process W l j,k = W l j,k (τ ) has time-dependent rate k ηe −τ l (u + (1 − e −τ )η), where l is one of the letters a, b, or c, see (4.1). We think of the Poisson process W l j,k as attached to the vertical edge (k, j) − (k, j + 1). Fix an initial condition σ(0) ∈ Ω quad , and R , We can define the evolution σ(τ ) restricted to {0, 1, . . . , R } × Z ≥1 as a function of all the Poisson processes W l j,k , in the spirit of the Harris graphical construction [Har78]. That is, for each seed pair (v 1 , v 2 ), if there is an arrival in one of the Poisson processes attached to the edge v 1 − v 2 , then it triggers the corresponding jump of the horizontal path up or down, which propagates to the right. We call such an arrival in a Poisson process the seed arrival.
To complete the construction of the process in {0, 1, . . . , R } × Z ≥1 , we need to show that in any finite time interval, there are almost surely finitely many seed arrivals in the Poisson processes W l j,k , where 0 ≤ j ≤ R , k ≥ 1. Note that this is not obvious since the jump rates depend linearly on the vertical coordinate k, and hence are unbounded.
Since the configuration in [0, R ]×[R, ∞) is initially empty, infinitely many seed arrivals might arise only when up-right paths perform infinitely many up jumps. More precisely, the top of the paths crossing the vertical line with x-coordinate R must jump up infinitely many times. See Figure 15 for an illustration.
The process of up jumps of the top path is bounded from above by the pure birth process with rate(k → k + 1) = Ck, for some constant C > 0. (This process is also called the Yule process.) It is well-known (e.g., [KM57]) that this pure birth process does not run off to infinity in finite time because the sum of its inverse rates diverges: ∞ k=1 (Ck) −1 = +∞. This implies that the desired process σ(τ ) restricted to {0, 1, . . . , R } × Z ≥1 does not perform infinitely many jumps in finite time, and hence completes the construction of the dynamics σ(τ ) (with the Markov generator Q(τ )) in the quadrant.
Formula (4.6) for the action of Q(τ ) on the homogeneous stochastic six vertex model follows as a Poisson limit of Proposition 4.2, as explained in Section 4.2.
5 Markov process preserving full plane Gibbs measures

Bulk limit of the quadrant dynamics. Heuristics
In this section we discuss the full plane continuous time dynamics arising in the bulk of the process Q(τ ) constructed in Section 4 above. Consider running Q(τ ) with an initial configuration sampled from the stochastic six vertex model P s6v, hom u with, say, step-∅ boundary conditions. Let ε > 0 be small, and consider a rectangular part of the lattice around a point with coordinates ( x/ε , y/ε ). Assume that the limit ε → 0 preserves the lattice scale, that is, the rectangular part of the lattice turns into the full plane Z 2 as ε → 0.
The local statistics of the path configuration around ( x/ε , y/ε ) are described [Agg20a] by the pure state of a slope (s, t) belonging to either the KPZ phase or the frozen phase, in the terminology from Section 2.3. See also Figure 14, left, for a simulation. In the rest of the discussion we ignore the frozen part and focus on the KPZ phase.
Slowing down Q(τ ) so that it runs at speed ε, the transition rates (for initiating jumps) around ( x/ε , y/ε ) in a finite time interval [0, τ 0 ] have the form Here O(ε) may depend on τ 0 , but τ 0 is fixed. Since y is also fixed, we see that as ε → 0, the jump rates converge to finite values proportional to R v 1 ,v 2 (u). In other words, we arrive at a (so far, hypothetical) full plane continuous time Markov chain with homogeneous rates const · R v 1 ,v 2 (u). Throughout this section we denote the full plane chain by C(t), t ∈ R ≥0 . Therefore, if there is reasonable locality in the original process Q(τ ) in the quadrant (more precisely, if we can turn off jumps outside of a large enough box around ( x/ε , y/ε ) without affecting the process in a smaller box), then in the bulk limit regime the mapping of the measures (4.6) turns into the statement that C(t) should preserve the local distribution, the pure Gibbs state with slope (s, t) in the KPZ phase.
In Sections 5 and 6 we prove the existence of the process C(t) together with the preservation of the KPZ pure states. We obtain C(t) directly from the jump rates, and not as a bulk limit of the dynamics in the quadrant. However, we employ the dynamics in the quadrant to show that the full plane process preserves the KPZ pure Gibbs state.

Admissible configurations of up-right paths
We aim to construct a continuous time Markov process C(t) on up-right path configurations in the full plane Z 2 with the following dynamics. Recall that to each path configuration we associate its seed pairs (v 1 , v 2 ) of vertically adjacent vertices, see Definition 4.5.
Definition 5.1 (Jumps in the process C(t)). Each seed pair (v 1 , v 2 ) has an exponential clock with rate R v 1 ,v 2 (u), see Figure 13. When the clock at (v 1 , v 2 ) rings, the horizontal path passing through this seed pair jumps up or down. This jump then instantaneously propagates to the right until it finds a vertical edge where it can stop (see Figure 11 for an illustration, and Section 3.6 for the definition of jump propagation). Note that in contrast with the dynamics in the quadrant, in the full plane case all rates are time-independent and are homogeneous throughout the plane.
Since the jumps may potentially propagate very far to the right, it is not immediately clear how to define the process C(t) with the jump rates from Definition 5.1 even locally. Indeed, a new jump can be initiated anywhere along a very long horizontal path, and close to its right end there could be infinitely many propagated jumps in finite time. In other words, defining the generator of this process formally as similarly to (4.5), where the sum is over all seed pairs (v 1 , v 2 ) of σ, may lead to divergence for some configurations σ. Therefore, let us first define the space of admissible configurations of up-right paths: Definition 5.2. Fix real R, A > 0. Let Λ R := [−R, R] 2 ∩ Z 2 be the finite square with side 2R around the origin. Let E R,A be the event (i.e., a subspace of up-right path configurations) such that there exists a horizontal sequence of vertices (x 0 , y), (x 0 + 1, y), . . . , (x 0 + n, y) =: (x 1 , y) where the path configuration around (x 0 + j, y) is not equal to for any 0 ≤ j ≤ n, where n > A and (x 1 , y) ∈ Λ R (here is a shorthand for the vertex (0, 1; 1, 0)). In words, configurations in E R,A have horizontal strings of occupied edges of length > A in an R-neighborhood of the origin. Clearly, for B > A we have E R,B ⊂ E R,A .
Define the set of admissible configurations to be In words, configurations in Ω are such that for some p ≥ 1 and all large enough R, every horizontal string of adjacent vertices of length larger than (log R) p ending in Λ R contains .
Next, let us show that under pure states in the KPZ phase, almost surely the configuration of up-right paths in admissible. Fix the parameters of the model q, u ∈ (0, 1), and the density of vertical occupied edges s = ρ ∈ (0, 1). Let π(ρ) = π ρ,ϕ(ρ) be the corresponding KPZ pure state (Section 2.3). Denote this measure by π, for short.
Lemma 5.3. Let ζ > 0. There exists a constant C ζ > 0 such that for all R > 0, Proof. It suffices to show that for some constant 0 < θ < 1, the probability under π of a particular sequence of vertices (x, y), (x+1, y), . . . , (x+n, y) not having a vertex equal to is upper bounded by θ n . Indeed, then we have by taking union bound: π(E R,C ζ log R ) ≤ const · R 2 θ C ζ log R = const · R 2−C ζ log(1/θ) , which may be made less than C ζ R −ζ by a choice of C ζ . Now, conditioning on the configuration up to the vertex (x+i−1, y), we have two possibilities (recall the vertex weights in Figure 3): • The horizontal edge exiting from (x + i − 1, y) is occupied. Conditioning on this, the probability of at (x + i, y) is at least (1 − ρ)(1 − δ 2 ).
), leads to the upper bound of the probability that no vertex in the sequence is by θ n/2 . This completes the proof.
Proof. Follows from Lemma 5.3 by a Borel-Cantelli type argument. Indeed, taking ζ > 1 in Lemma 5.3, we have Since E R,C ζ log R ⊃ E R,(log R) 2 for all R large enough, we see that π(∪ R 0 ∩ R≥R 0 E R,(log R) 2 ) = 0. This implies that π(Ω c ) = 0 as Ω c is the intersection for all p ≥ 1, and the part with p = 2 already has zero probability.

Formulations
We construct the full plane dynamics C(t) as a limit of the truncated processes, which exists for admissible configurations.
Definition 5.5 (Truncated process). Let R > 0, and define C R (t), t ∈ R ≥0 , to be a continuous time Markov process on configurations of up-right paths with the jumps described in Definition 5.1, with the modification that new jumps can be initiated only by seed pairs (v 1 , v 2 ) inside the finite square Λ R .
The truncated processes C R are clearly well-defined.
Remark 5.6. For a fixed initial configuration C R (0) = C 0 the distributions of {C R } t∈R ≥0 can be coupled for various values of R. Indeed, a natural coupling of C R and C R , R > R, corresponds to using the same Poisson clocks for both C R and C R inside Λ R , and turning on additional Poisson clocks for all seed pairs in Λ R \ Λ R for C R .
We aim to show that when the initial configuration C 0 belongs to Ω (5.2), then the restrictions of the processes {C R (t)} 0≤t≤t 0 to Λ K for fixed K, t 0 stabilize as R → +∞. This would lead to the two main results of this section which we now formulate.
Theorem 5.7 (Existence). For all C 0 ∈ Ω, there exists a Markov process {C(t)} t∈R ≥0 started from C(0) = C 0 , which evolves according to the jumps described in Definition 5.1. That is, the action of the generator (G u f )(σ) (5.1) is well-defined for σ ∈ Ω, and corresponds to a continuous time Markov process. Moreover, for all t we have C(t) ∈ Ω.
Recall that π is the KPZ pure Gibbs state with density ρ of the occupied vertical edges, where 0 < ρ < 1.
Theorem 5.8 (Preservation of KPZ pure states). The Markov process C(t) from Theorem 5.7 preserves the measure π. More precisely, for any cylindric function f we have where E C 0 is the expectation with respect to the Markov chain C(t) started from C 0 , and π[f ] is the integral of f over the measure π.
Note that in (5.3) we could integrate over the set of all up-right path configurations instead of admissible configurations Ω, but this is the same since π(Ω) = 1 by Proposition 5.4.
The proofs of Theorems 5.7 and 5.8 occupy the rest of Section 5 and also Section 6.

Proof of Theorem 5.7
The following is the first and main lemma in our argument. Recall the sets E R,A from Definition 5.2 in which there are long horizontal paths leading to far propagation of jumps. Throughout the rest of the section we assume the following lemma.
Lemma 5.9 (Main estimate). Fix real numbers ζ > 0, p 0 > 0, C > 0, T > 0, and R > 0. There is a constantC =C ζ such that the following holds for allR ≥ R > 0. If we start the truncated dynamics CR(t) from a configuration C 0 / ∈ E R,C(log R) p 0 , then for some p ≥ p 0 , In words, if there is initially a bound on the length of jump propagation, then at each finite time a slightly worse bound holds with high probability. The proof of Lemma 5.9 utilizes a nontrivial coupling and is postponed till the next Section 6.
Given a path configuration C, for e an edge in Z 2 , denote by δ e = δ e (C) the path occupation indicator variable for that edge. For a trajectory of the truncated process C R (t), denote the corresponding edge indicator at time t by δ R e (t). The random variables δ R e (t) are naturally coupled for various values of R, see Remark 5.6.
Lemma 5.10. Let ζ, T > 0, and e be a fixed edge of Z 2 . Let the initial condition C 0 belong to Ω. Then there exists a constant C = C ζ > 0 such that for any sufficiently large R , R with R > R we have In words, the probability that two truncations diverge on a fixed edge is small.
Proof of Lemma 5.10. We closely follow the proof of Proposition 7.6 of [Ton17], adapting it to our setting. Suppose that e is the vertical edge ((0, 0), (0, 1)). This does not restrict the generality since a change in a horizontal edge is accompanied by a change of a vertical edge which is sufficiently close (since C 0 ∈ Ω).
Plug in the constant p 0 from the fact that C 0 ∈ Ω and C = 1 into Lemma 5.9. Let p,C be the constants from Lemma 5.9. Let E 1 denote the event that either C R (t ) ∈ E R,C ζ (log R) p or C R (t ) ∈ E R,C ζ (log R) p for some time t ∈ [0, T ]. In words, under E 1 there exists a time t when in C R (t ) or C R (t ) we see a horizontal string of vertices of lengthC ζ (log R) p , none of which are . By Lemma 5.9, P(E 1 ) ≤CR −ζ , which means that in the rest of the proof we may assume that E 1 did not happen.
Suppose δ R e (t) = δ R e (t) for some t ∈ [0, T ] and E 1 did not happen. Denote this event by E 2 . Let t 1 ≤ T be the first time at which δ R e (t + 1 ) = δ R e (t + 1 ). (Throughout the proof, t ± mean one-sided limits.) Then, there must have been a clock that rang at time t 1 which caused e to change in one but not the other of the configurations C R and C R . Note that in both processes C R and C R the jumps cannot propagate by more thanC(log R) p . Thus, there must be an edge e 1 touching the rectangle [−C(log R) p , 0] × [0, 1] such that δ R e 1 (t − 1 ) = δ R e 1 (t − 1 ). Let t 2 < t 1 be the first time at which δ R e 1 (t + 2 ) = δ R e 1 (t + 2 ), and there must be another edge e 2 to the left of e 1 such that δ R e 2 (t − 2 ) = δ R e 2 (t − 2 ), and the first time t 3 when the occupations of e 2 diverged in two processes. We continue this argument, and obtain a sequence of edges and times, which terminates with an edge e n outside of Λ R . We see that C R could not ever change the state of e n , and it might have changed under C R .
To summarize, in E 2 there exists a sequence of clocks that ring in Λ R , at times 0 < t n < t n−1 < · · · < t 2 < t 1 and positions (x 1 , y 1 ), (x 2 , y 2 ), . . . , (x n , y n ), such that |y i − y i−1 | ≤ 1 and x i − x i−1 ≤C(log R) p . This implies that n ≥ R/(C(log R) p ). Let us bound the probability of E 2 from above. We use two more observations: • We choose n locations in [0, R] where the clocks must ring such that the distance between two consecutive locations is ≤C(log R) p . Therefore, each next location chooses from at most 3C(log R) p available locations, and thus the total number of ways to choose the locations is upper bounded by (3C(log R) p ) n .
• There are at least n clock rings in [0, T ], and we know that the rate of each clock ringing is bounded by a constant. Therefore, the total number of clock rings during [0, T ] is stochastically dominated by a Poisson random variable with mean θT for some 0 < θ < +∞. Therefore, P(Poisson(θT ) ≥ n) ≤ const · (eθT ) n /(n!).
Therefore, we have For fixed R, the sum over n is, up to a constant, bounded from above by its first term (C ) n 0 (log R) pn 0 /(n 0 !), where n 0 = R/(C(log R) p ). As R → +∞, one readily checks that this first term goes to zero faster than any power of R, and so is of order o(R −ζ ). Combining the bounds on the probabilities of E 1 and E 2 yields the result.
Lemma 5.11. Fix an admissible initial configuration C 0 ∈ Ω, and an edge e. Recall that δ R e (t) is the occupation of e under the R-truncated process. With probability 1, the limit δ e (t) := lim Since there are countably many edges, the almost sure limit in (5.4) also exists simultaneously for all edges e.
Proof of Lemma 5.11. Let ∆t > 0 be fixed. If lim R→∞ (δ R e (t), t ∈ [0, ∆t]) does not exist in the sense of uniform convergence, then for infinitely many positive integers R, there exists a time t ∈ [0, ∆t] for which δ R e (t) = δ R+1 e (t). However, by Lemma 5.10, for arbitrary ζ > 1 and a constantC > 0 depending on ζ but not on R, we have By the Borel-Cantelli lemma, for any ∆t, we get, on a probability 1 event, the uniform convergence of the paths {(δ R e (t), t ∈ [0, ∆t])} as R → ∞. Taking the intersection of these events as ∆t ranges over positive integers yields the result. Lemma 5.9 and the proof technique of Lemma 5.10 imply the following bound on the propagation speed: Lemma 5.12. Let R 1 and R 2 > 0 be truncation radii, let R > 0 be an integer with R < min(R 1 , R 2 ), and let k < R be finite (i.e., not thought of as large). Fix time ∆t > 0.
Let C 1 (0) = C 1 , C 2 (0) = C 2 ∈ Ω be configurations that agree inside Λ R . Then there exists a coupling of the trajectories of C R 1 1 (t) and C R 2 2 (t), under which they agree inside Λ k for all t ∈ [0, ∆t] with probability at least 1 − CR −ζ for C independent of R, R 1 , R 2 .
Let us now show that the edge random variables δ e (t) indeed lead to a Markov process, which, moreover, stays admissible throughout the whole time.
Lemma 5.13. For any admissible C(0) = C 0 ∈ Ω, the joint distribution of the edge occupation trajectories {δ e (t)} t≥0 over all edges e defines a Markov process C(t) with values in Ω started from Proof. The fact that the stochastic process C(t) (coming from the limit in Lemma 5.11) stays in the space Ω of admissible configurations once started in Ω follows asR → ∞ limit of the estimate in Lemma 5.9.
Let us show the Markov property of the stochastic process C(t). Let F {C(t)} t∈[0,s] = F {δ e 1 (t), . . . , δ e k (t)} t∈[0,s] be a bounded continuous functional on the space of right continuous paths in Ω, which depends only on the evolution of finitely many edges e 1 , . . . , e k , and only on their values up to some time s. Let F t 0 be the σ-algebra generated by the random variables δ e (t) with t ≤ t 0 . Then we have (5.5) In the final expression above, C R (t 0 ) is a random configuration obtained from running the truncated dynamics started from C 0 for time t 0 . The first equality is by the bounded convergence theorem, and the second is by the Markov property for the truncated dynamics. We claim that the final expression in (5.5) is equal to First, with high probability we have C R (t 0 ) = C(t 0 ) on a large square Λ R 0 for large but fixed R 0 < R, as R → ∞. Indeed, by Lemma 5.10, for a constant C R 0 ,ζ independent of R we have As a result, almost surely lim Now let us apply Lemma 5.12, or, more precisely, its limiting version as R 1 → +∞. Conditioned on the fact that the configurations C(t 0 ) and C R (t 0 ) agree inside Λ R 0 , we can upper bound the probability that the edge indicators for e 1 , . . . , e k are different in C(t) (which is the limit of C R 1 (t) as R 1 → +∞) and C R (t) for some time t ∈ [0, s]. Namely, for some C 0 independent of R 0 . Thus, with probability 1 we have Taking R 0 → ∞ through the positive integers gives the Markov property.
Proof of Theorem 5.7. To finalize the proof of Theorem 5.7, it remains to verify the correct jump rates in the generator G u (5.1), as the existence of the process in the space Ω of admissible configurations follows from Lemma 5.13. Take a finite box Λ k and a large integer R > k. We approximate the dynamics of {δ e (t)} e∈Λ k by the truncated dynamics {δ R e (t)} e∈Λ k . Keeping track of the error terms in Lemmas 5.9 and 5.10, we see that the probability that {δ e (t)} e∈Λ k = {δ R e (t)} e∈Λ k for some t ∈ [0, ∆t] is upper bounded by C∆tR −ζ for any ζ > 0 and for some C depending on ζ.
Let f be a cylindric function depending only on the occupation of the edges in Λ k . Denote by G R u the generator as in (5.1), but with with Poisson clocks outside of Λ R turned off. Let f 0 := f (C 0 ) be the value of f at the initial configuration C 0 ∈ Ω. We have Then, taking ∆t → 0 we have lim sup for some constant C. Sending R → ∞ gives the desired generator.

Proof of Theorem 5.8
We now show that the Markov process C(t) preserves the KPZ pure Gibbs state π = π(ρ) with density ρ of the occupied vertical edges, where 0 < ρ < 1. For this, we utilize the dynamics in the quadrant constructed in Section 4.
Definition 5.14. Fix a large positive integer N . Let ν N,u be the stochastic six vertex measure with parameters u, q (cf. Section 2.1; and we suppress the dependence on q) in the quadrant with the empty-λ boundary conditions (Definition 2.5), where λ ⊂ Z ≥0 is a random subset defined as follows. For each N 2 ≤ i ≤ 3N 2 , independently toss a coin with probability of success ρ. When the coin is a success, include i into λ. There are no incoming vertical paths at the bottom boundary outside the interval N 2 , 3N 2 .
Let M = M (N ) be such that the ratio M/N converges to a positive constant ≤ 1−δ 2 20 , where δ 2 is in a small neighborhood of δ 2 (u) (2.2). For any k < M , define the k × k box Λ For a cylindric function f on full plane configurations such that f ∈ σ(Λ k ) (i.e., depending only on the path configuration in Λ k = [−k, k] 2 ), we denote by f N its pullback under the shift by (−N, −M ). That is, f N is a function depending only on the configuration in the shifted k × k box Λ N k .
Denote by C N (t) the process in the quadrant (Definition 4.6 and Theorem 4.9) started from ν N,u and run at the (slowed down) speed 1/M (N ). Take R < M , and let C R N (t) be the truncation of the process C N (t) corresponding to turning off the Poisson clocks outside of the shifted lattice square [N − R, N + R] × [M − R, M + R]. We assume that C R N (t) also starts from ν N,u . Recall that we denote by C(t) and C R (t) the usual and the truncated full plane processes, where for C R (t) we turn off the Poisson clocks outside the lattice square Λ R = [−R, R] 2 . Let C(t) and C R (t) start from the pure state π.
To establish Theorem 5.8, it suffices to show that E π [f (C(t))] = π[f ] for any bounded cylindric function f ∈ σ(Λ k ). Recall that E π is the expectation with respect to the process started from the pure state π, and π[f ] is the integral of the function f against π.
We will approximate f (C(t)) by f N (C N (t)) employing a sequence of couplings coming from the bounds on information propagation (established in Section 5.4), and also from a coupling lemma and the local statistics theorem of [Agg20a] applied to ν N,u . The result will then follow from the mapping of the stochastic six vertex measures under the dynamics on the quadrant (Theorem 4.9), and continuity of ν N,u in u.
We proceed by establishing several lemmas.
Proof. For large enough C we know that π(C(0) ∈ E R,C log R ) ≤ CR −ζ , see Lemma 5.3. Furthermore, given an initial configuration C(0) such that C(0) / ∈ E R,C log R , under the standard coupling of the full plane process C(t) with its truncated counterpart C R (t), the configurations C(t) and C R (t) agree on Λ k (the subset determining the values of f ) except on an event with probability at most O(R −ζ ).
Proof. This is a statement about finite state space Markov chains. It suffices to see that the jump rates of C R N converge to those of C R , and that ν N,u | Λ N R converges to π| Λ R (if we identify the two boxes). The former claim is established by the computation in Section 5.1. The latter claim follows, for example, from an application of Proposition 2.17 of [Agg20a].
Note that the events E R,A (Definition 5.2) make sense for quarter plane configurations, so long as we replace the box Λ R centered at (0, 0) with the shifted box Λ N R , where R is fixed, and N is sufficiently large so that R < min(N, M (N )). Denote the corresponding shifted event centered at (N, M ) by E N R,A .
Lemma 5.17. There exists C such that for all positive integers R < min(N, M (N )) we have where C N (t) is the quarter plane process started from ν N,u .
Proof. This is established similarly to Lemma 5.9, see also Lemma 5.12.
Lemma 5.18. We have Proof. This follows similarly to Lemma 5.10. Indeed, one can bound the probability that, at any time in some compact time interval, either C R N (t) or C N (t) develop long sequences of consecutive horizontally adjacent vertices which are not . These bounds are provided by Lemma 5.17 and an analogue of Lemma 5.9 (proved in the same manner by coupling, see Section 6 below).
Proof. Here we use the fact that the map ϕ(ρ) (2.6) is continuous in u. Namely, for each N > 0 let M = M (N ) = 1−δ 2 −ε 20 N for some ε > 0 small enough. Note that we can couple the boundary conditions of a six vertex configuration C N sampled from ν N,u with those of C sampled from π such they agree with probability 1 on [N/2, 3N/2] × {0}. Therefore by Proposition 2.17 of [Agg20a], for some constant c independent of N , there is a coupling of the configurations C N and C such that they agree on the set of edges in [N − M, N + M ] × [0, 2M ] with probability at least 1 − c −1 e −cM . We apply this with u replaced by u + t/M and π replaced by the Gibbs measurẽ π M with x slope ρ and spectral parameter u + t/M . Then we note thatπ M | Λ N k can in turn be coupled with π| Λ k (the marginal of the Gibbs measure with spectral parameter u on Λ k , where the edges in Λ k are identified with those of Λ N k ), such that they agree with probability 1 − ε M , with ε M → 0 as M → ∞. Therefore, we get a coupling of ν N,u+t/M | Λ N k with π| Λ k such that the configurations agree with probability lower bounded by 1 − c −1 e −cM − ε M . Since M → ∞ as N → ∞, this implies the claim. Now we can finish the proof of Theorem 5.8. Using Lemmas 5.15, 5.16 and 5.18, the measure mapping property of the quarter plane dynamics (Theorem 4.9), and then Lemma 5.19, we get where ε N → 0 as N → ∞ with R fixed, and ε R → 0 as R → ∞. So taking N → ∞ first, then R → ∞, yields the result.

Proof of Lemma 5.9 via coupling
Recall that E R,A is the set of up-right path configurations which have horizontal strings of occupied edges (i.e., no vertices of occupation type ) of length > A, in an R-neighborhood of the origin. Under the truncated full plane dynamics C R (t), configurations in E R,A might lead to jump propagation of length > A. In this section we prove Lemma 5.9 which states that if the initial configuration was not in E R,C(log R) p 0 , then with high probability it will never be in E R,C(log R) p with some p ≥ p 0 under the dynamics CR(t), up to time t ≤ T . HereR ≥ R is an arbitrary integer (and we keep the notationR throughout the section for consistency with Lemma 5.9). We achieve this bound via a monotone coupling of the dynamics CR(t) on a given horizontal slice with a particle system which is easier to analyze.

Annihilation-Jump particle system
We start by defining the jump rates of the Annihilation-Jump particle system (AJ). Its state space consists of a sequence of ⊕ particles a i ∈ Z ≤0 , and a sequence of particles b i ∈ Z ≤0 , which satisfy either · · · < b 2 < a 2 < b 1 < a 1 ≤ 0, or · · · < a 2 < b 2 < a 1 < b 1 ≤ 0. Pick m > 0, this parameter is the overall time scaling factor in the AJ dynamics. The possible jumps and their rates are as follows: • Any particle at position z, with the closest particle to its left at position y < z, jumps to position x at rate m for any y < x < z. This rule does not distinguish ⊕ and particles, and is like in the Hammersley process [Ham72], [AD95].
• Take any particle b i and let the closest ⊕ particle to its right be a j (so j = i or i −1). The pair (b i , a j ) disappears at rate 2m, and the sequences of remaining particles {a k } k =j , {b k } k =i are relabeled. If i = 1 and b 1 is the rightmost particle, then b 1 simply disappears.
• Take any ⊕ particle a i and let the closest particle to its right be b j (so j = i or i − 1). The pair (a i , b j ) disappears at rate 2m, and the remaining particles are relabeled. If i = 1 and a 1 is the rightmost particle, then a 1 simply disappears.
Similarly to the full plane dynamics C(t), a priori it is not clear that the AJ particle system is well-defined. Indeed, there may be initial configurations leading to infinitely many jumps through a finite space in finite time. However, we only need to analyze a truncated version of AJ.
Definition 6.1 (Truncated AJ particle system). LetR, N 0 > 0 be positive integers which determine the truncation. Given an initial configuration, define a i (t) = aR i (t), b i (t) = bR i (t) to evolve according the following rules. Put an independent rate m Poisson clock P z (t) at each lattice site z ∈ Z ≤0 . If P z , |z| ≤R, rings at time t, then we take the leftmost particle to the right of z (regardless of whether it is ⊕ or ) and place it into z. If the clock P z rings at a particle z = a i or z = b j , ignore this ring.
Also take a collection of rate 2m Poisson clocks P a i (t) and P b i (t) for each integer i ≥ 1. For these clocks, i rings at time t, and i ≤ N 0 , then b i (t − ) and the ⊕ particle to its right, a j (t − ), annihilate each other and disappear, and we relabel the particles. If i = 1 and b 1 is the rightmost particle, then b 1 simply disappears.
• If P a i rings at time t, and i ≤ N 0 , then a i (t − ) and the particle to its right, b j (t − ), annihilate each other and disappear, and we relabel the particles. If i = 1 and a 1 is the rightmost particle, then a 1 simply disappears.
In the truncated dynamics we ignore all other clock rings.

Stochastic domination
Here we describe in which sense the AJ system dominates the full plane dynamics.
Let λ (0) , λ (1) be two infinite interlacing signatures, that is, λ i+1 for all i ∈ Z. This pair represents a six vertex path configuration on a one-row lattice Z × {1}, by encoding the positions of occupied vertical edges entering (λ (0) ) and leaving (λ (1) ) this row. By analogy with the AJ system, below we also refer to occupied vertical edges as particles.
as follows: In other words, the sequence A = (A 1 > A 2 > · · · ) indexes the positions of the particles in λ (1) with no particle at the same position in λ (0) . Similarly, In other words, the sequence B = (B 1 > B 2 > · · · ) indexes the positions of the particles in λ (0) with no particle at the same position in λ (1) .
For two consecutive horizontal rows (at the 0-th and the first slice) of a full plane six vertex model configuration evolving under some Markov dynamics, denote by A(t), B(t) the timedependent random variables corresponding to this dynamics.
If s = (s 1 > s 2 > s 3 > · · · ), r = (r 1 > r 2 > r 3 > · · · ) are two decreasing sequences of numbers, then we say s ≤ r if s i ≤ r i for all i = 1, 2, 3, . . . . We are now in a position to formulate the lemma on stochastic domination. Lemma 6.3. Let CR(t) denote the trajectory of the full plane six vertex model under the truncated dynamics, which gives rise to the quantities A(t), B(t). We can couple CR(t) with a truncated AJ particle system (a(t), b(t)) with initial configuration a(0) = A(0), b(0) = B(0) and suitableR, N 0 , and m, in such a way that with probability 1, for all t ≥ 0 we have In words, the AJ particles are always further to the left compared to their six vertex model counterparts. We say that (6.1) means that the AJ system stochastically dominates this slice of the full plane dynamics. See Figure 16 for an illustration.
Proof. We couple the Poisson clocks in the truncated full plane six vertex model dynamics CR(t) with the clocks in a suitably truncated AJ particle system. For the truncation in the AJ system we take the sameR, and let let N 0 be the maximum integer such that A N 0 (0), B N 0 (0) ≥ −R.  Figure 16: Putting the six vertex and the AJ configurations together. In the figure, the configurations are ordered as in (6.1).
in the full plane dynamics defined in (4.1). When a clock at an edge rings (almost surely, there is at most one ring in finite time in the truncated process), depending on the current configuration and possibly on the outcome of an independent coin flip, this ring initiates a jump in CR, or we ignore it. Here coin flips are needed to model smaller jump rates. For example, if m = a > b, then we model rings at rate b from the Poisson clock of rate m and the coin with probability of success b/m. At time t − , denote the six vertex configuration on rows 0 and 1 by λ (0) = λ (0) (t − ), λ (1) = λ (1) (t − ), and the AJ system configuration by a = a(t − ), b = b(t − ). We first define an auxiliary AJ like particle system which is a function of the Poisson clock rings in CR. It then will be evident that this auxiliary dynamics is dominated by the AJ system in the same sense as in (6.1), which will lead to the desired statement. We refer to this auxiliary particle system as the "AJ system" to simplify the wording.
Pick z ∈ Z with |z| ≤R, and consider the configuration of occupied edges around z. There are six possible configurations corresponding to the allowed configurations under the six vertex model (cf. Figure 3). We denote vertical edge locations corresponding to the configurations at horizontal levels 0 and 1 by (z, 0) and (z, 1), respectively. In the six cases below we assume that a clock rings at one of these edges, and define how the particle system changes as a result of this ring.
1. Let the path configuration around z be (0, 0; 0, 0), which means that λ (1) Figure 13 we see that only the clock at the bottom vertical edge (z, 0) may ring. In this case, we let the leftmost AJ particle to the right of z (if it exists) jump into z, provided that z is not occupied by another AJ particle. If no such AJ particle exists, ignore this ring. From now on, we will simply say "the particle to the right of z attempts to jump into z".
On the six vertex side, after this clock ring the edge (z, 0) might stay empty if λ (0) , λ (−1) locally do not look like the configurations in Figure 13, left, or if the coin flips do not produce an actual jump initiation in CR. In the remaining situation when the edge (z, 0) becomes occupied, one of the vertical occupied edges at (z , 0), z > z, must instantaneously become empty due to the jump propagation. If z is not one of the B j 's, then this corresponds to a right jump in the six vertex model (namely, a "creation" of a new pair A k , B k to the right of z, and relabeling of A, B). Alternatively, z could be equal to the leftmost of the B j 's which are greater than z, and this is the furthest left jump that may occur under CR. Clearly, in all these cases the domination (6.1) is preserved by the jumps in the two systems.
Here is an illustration of the moves in the six vertex model and the AJ system, with the furthest possible left jump of (A, B) under CR: The remaining five cases are considered similarly, and we discuss them in less detail. One readily sees that in each of the remaining five cases, the domination (6.1) is preserved.
2. Let the path configuration around z be (1, 1; 1, 1), which means that z = λ i−1 . In this case, only the clock at the bottom edge (z, 0) may ring, and in the AJ system we define that the particle to the right of z attempts to jump into z. In the illustration below we display the furthest possible left jump of (A, B) under CR, which arises when the edge (z, 0) becomes empty, and the furthest possible edge to the right becomes occupied: 3. Let the path configuration around z be (0, 1; 0, 1), which means that λ i . In this case, only the clock at the top edge (z, 1) may ring, and in the AJ system we define that the particle to the right of z attempts to jump into z. In the illustration below we display the furthest possible left jump of (A, B) under CR: 4. Let the path configuration around z be (1, 0; 1, 0), which means that z = λ (1) i . In this case, only the clock at the top edge (z, 1) may ring, and in the AJ system we define that the particle to the right of z attempts to jump into z. In the illustration below we display the furthest possible left jump of (A, B) under CR, which arises when the edge (z, 1) becomes empty, and the furthest possible edge to the right becomes occupied: 5. Let the path configuration around z be (0, 1; 1, 0), which means that λ i−1 . Then, by definition, z = A k for some k ≤ i. Here the clock may ring at either (z, 0) or (z, 1) (which corresponds to the double rate 2m of annihilation in the AJ system). In both situations, we define that a k and the particle to its right annihilate (or a 1 disappears if it was the rightmost of the AJ particles). Here are the illustrations of both cases: In CR, both clock rings might lead to no change, or to a right jump of A k , or to an annihilation of A k and B k−1 (and relabeling of the remaining particles in A, B). In all these cases, the domination (6.1) is preserved.
6. Let the path configuration around z be (1, 0; 0, 1), which means that λ (1) i . Then, by definition, z = B k for some k ≤ i. Here the clock may ring at either (z, 0) or (z, 1) (again, this corresponds to the fact that the annihilation rate is 2m). In both situations, we define that b k and the particle to its right annihilate (or b 1 disappears if it was the rightmost of the AJ particles). Here are the illustrations of both cases, and we similarly see that the domination (6.1) is preserved: We have thus constructed the auxiliary AJ like particle system on Z which is a function of the rate m Poisson clock rings in CR, and which dominates the dynamics of (A, B). In the auxiliary system, particles indeed annihilate at rate 2m per pair (a, b) or (b, a), but the auxiliary system lacks some of the AJ system's left jumps. Namely, if a clock from CR at z / ∈ {a i } i≥1 ∪ {b i } i≥1 contributes to an annihilation event, then it did not produce a jumping event in the auxiliary system. Therefore, adding extra independent jump events at rate m produces the full AJ system, which clearly jumps to the left more often than the auxiliary system. We see that the coupling we constructed indeed satisfies the domination (6.1), which completes the proof of Lemma 6.3.

Completing the proof of Lemma 5.9 via an estimate in the AJ system
To finalize the proof of Lemma 5.9 we need to show that, given that the initial configuration was not in E R,C(log R) p 0 for some R, C (that is, it did not have long strings without vertices of type ), then with high probability the configuration will not be in E R,C(log R) p up to time t ≤ T , where p ≥ p 0 andR ≥ R. Due to the Z 2 translation invariance of the dynamics, we may consider the event that a long string of vertices without starts at (0, 1), and then take a union bound over all vertices in Λ R (there are const · R 2 of them).
By the stochastic domination shown in Lemma 6.3, it suffices to upper bound the following probability under the AJ system: (6.2) Here the AJ particle system is truncated atR, N 0 , and starts from the configuration a(0) = A(0), b(0) = B(0), as described before Lemma 6.3. Note that since the AJ system jumps only to the left, (6.2) also is an upper bound for P AJ (∃ t ∈ [0, T ] : a 1 (t) ≤ −C(log R) p ). The latter quantity is an upper bound for the corresponding quantity in the dynamics CR. The constants C > 0 and p ≥ p 0 > 0 in (6.2) will be determined later to make the probability (6.2) sufficiently small, namely, of order R −ζ for any ζ > 0 (as R grows).
Let n be such that a n (0) ≤ −C(log R) p . Observe that n is large for large R because the initial configuration is not in E R,C(log R) p 0 , namely, Moreover, a n (0) cannot be too large in the absolute value because the initial configuration is not in E R,C(log R) p 0 , namely, |a n (0)| ≤ Cn(log R) p 0 . (6.4) To determine if the event in (6.2) occurred, we may only look at the behavior of the particles started in (a n (0), 0] up to time T . Moreover, by taking a n (0) sufficiently small, we may assume that the event in (6.2) is due entirely to a combination of annihilations, and jumps caused by Poisson clocks in the interval (a n (0), 0] up to time T . Moreover, a 1 (T ) ≤ −C(log R) p implies that there are no particles left at time T between −C(log R) p and 0. We treat separately the cases when this absence of particles is mainly due to annihilations or mainly due to particles jumping out.
Denote by K(n) the (random) number of particles out of the ones started at time t = 0 in (a n (0), 0], and which did not get annihilated up to time t = T . Fix r ∈ 1 2 , 1 . If a 1 (T ) ≤ −C(log R) p and K(n) < n r , then there were many annihilations, while if K(n) ≥ n r , then there should have been many jumps. First, we estimate the probability of many annihilations: Lemma 6.4. For some c > 0 we have P AJ (K(n) < n r ) ≤ exp{−cn r } for all large enough n.
By (6.3), exp {−cn r } decays to zero faster than any power of R as R → +∞ provided that (p − p 0 )r > 1 (which is ensured by choosing p large enough).
Proof of Lemma 6.4. We need to bound the probability that at least n− n r /2 particles disappear during time T . This is upper bounded by P n i= n r /2 where {E i } are independent exponential random variables, where the rate of E i is equal to 4i (the rate is 4 since particles can annihilate with their left or right neighbors). This expression is bounded using [Jan18, Theorem 5.1.(iii)] with a * = 4 n r /2 (minimum of the rates of the E i 's), µ = 1−r 4 log n + O(1) (mean of the sum), λ = T /µ ≤ 1, and the bound is of the form e −a * µ(λ−1−log λ) . The dominating term in the exponent (going fastest to −∞ as n → +∞) is −a * µ log(1/λ) ∼ −const · n r (log n)(log log n), which leads to an estimate ≤ e −cn r , as desired.
Lemma 6.4 bounds the number of annihilations. Let us now consider the case K(n) ≥ n r . Then the main contribution to the event a 1 (T ) ≤ −C(log R) p comes from having many jumps. Without annihilations, the AJ system is simply the Hammersley process on Z, a natural analogue of the Hammersley process on R introduced in [Ham72], see also [AD95], [Sep01]. The Hammersley process on Z (which is essentially the same as the PushTASEP via the particle-hole involution) is coupled to the semi-discrete Poisson last passage percolation [Ton17, Section 5]. We thus observe that when K(n) ≥ n r , the event a 1 (T ) ≤ −C(log R) p may occur only if there exists an up-right path of percolation length at least n r in the (rate 1) semi-discrete Poisson last passage percolation in the space-time region (a n (0), 0] × [0, T ] ⊂ Z × R (see Figure 17 for an illustration). Denote by E n the event that such a path exists. We upper bound its probability as follows: T 0 a n (0) The last passage percolation length is the maximal number of Poisson points collected by an up-right path from (a n (0) + 1, 0) to (0, T ), where we count at most one point per lattice site. In terms of jumps in the AJ system, the Poisson points collected by the maximal up-right path are the times when particles jump into the corresponding locations. We count at most ony point per lattice site since several clock rings at one site might not lead to actual jumps in the AJ system. Lemma 6.5. Let n be such that a n (0) ≤ −C(log R) p . Then we have for somec > 0: Proof. We first estimate the probability that the maximal (in the sense of last passage percolation) up-right path has length M , where M ≥ n r . The existence of such a path means that there are M Poisson points picked up by the path (see Figure 17 for an illustration). The number of configurations of such points is equal to the number of sequences x 1 < x 2 < . . . < x M with a n (0) < x 1 and x M ≤ 0, which is upper bounded by |an (0) We see that these quantities decay in M faster than the terms of a geometric series, so their sum over M ≥ n r (coming from the union bound) is bounded by a constant times the first term. This produces the desired estimate.
To finalize the proof of Lemma 5.9, we pick the constantsC and p to bound the desired probability P AJ (a 1 (T ) ≤ −C(log R) p ) from above. We use Lemmas 6.4 and 6.5 and a union bound over n ≥ n 1 := C C (log R) p−p 0 , so that a n (0) ≤ −C(log R) p : (6.5) The first series in (6.5) is bounded from above by (power of n 1 ) · exp (−const · n r 1 ) = (power of log R) · exp −const · (log R) r(p−p 0 ) . (6.6) Indeed, one can bound the tail n≥n 1 e −cn r by e −cn r 1 + ∞ n 1 e −cx r dx. The integral is equal to a constant times an incomplete Gamma integral of the form ∞ z t a−1 e −t dt. For fixed a and large z (which is our case), the behavior is of the form e −z times a power of z (e.g., see [DLMF, 8.11(i)]). We thus see that by (6.6), the first series in (6.5) decays faster than any power of R as R → +∞ as long as (p − p 0 )r > 1.
For the second summand in (6.5) we have, using (6.3) and (6.4): log c|a n (0)eT | n r ( n r !) 2 ≤ const · n r log |a n (0)| + n r log(eT ) − 2rn r log n + 2n r + O(log n) ≤ const · n r (1 − 2r)(p − p 0 ) log(log R) + p 0 log(log R) + lower order terms Since r > 1 2 , by taking p large enough we may make the first term dominate as R → +∞. This leads to an overall decay of the second sum in (6.5) as R → +∞, faster than any power of R.
To conclude, the desired probability (6.5) is O(R −ζ ) for any ζ > 0. Taking the final union bound over all vertices in Λ R , as discussed in the beginning of this Section 6.3, multiplies our estimate by const · R 2 . With this factor the probability still decays faster than any power of R, and so we arrive at the estimate in Lemma 5.9. 7 Current and hydrodynamics 7.1 Computing the current Recall that with each path configuration of the six vertex model in Z 2 we associate the height function h(x, y) (defined up to a constant), see Section 1.2 in the Introduction. The full plane dynamics C(t) gives rise to the time-dependent random function h t (x, y). For the KPZ pure state π(s) (we recall its definition in Section 2.3), we define the corresponding current (average change of the height function) by J(s, ϕ(s)) := 1 t E π(s) (h t (0, 0) − h 0 (0, 0)) , (7.1) where the initial height function h 0 corresponds to the configuration distributed as π(s). Instead of (0, 0), we could take an arbitrary face of the lattice (by translation invariance of the measure and the dynamics). The right-hand side of (7.1) is independent of t, so we may send t → 0, and write J(s, ϕ(s)) = ∂ ∂t E π(s) h t (0, 0) t=0 .
Therefore, we may compute the current by looking at the Markov generator G u of C(t) given by (5.1). Namely, we have J(s, ϕ(s)) = where ∆ v 1 ,v 2 h 0 (0, 0) is the (signed) change of the height function at (0, 0) triggered by initiating the jump at the vertical edge v 1 − v 2 . In taking the expectation, we assume that h 0 is distributed according to π(s). Recall R v 1 ,v 2 (u) is equal to either of the quantities a(u), b(u), or c(u) (4.1) depending on the surrounding paths, see Figure 13.
Proof of Theorem 7.1. We use (7.2) and the description of π(s) as a trajectory of the stationary stochastic six vertex model, as discussed in Section 2.3. Throughout this proof, we denote t := ϕ(s), for short, and use the notation δ 1 , δ 2 for vertex weights, see (2.2). We compute the average change of height at the face (0, 0) by assuming that this change comes from a jump of a horizontal path of a specified structure to the left of (0, 0). Namely, to the left of (0, 0) we find the rightmost pair of vertices of one of the following six kinds (here we use the traditional names for the six vertices, see Figure 3): (c 1 , a 1 ), (b 2 , c 2 ), (c 1 , c 2 ) (for up jumps, which contribute −1 height change); (b 1 , c 1 ), (c 2 , a 2 ), (b 1 , a 2 ) (for down jumps, which contribute +1 height change), which is at position ((x, 0), (x, 1)) for some x ≤ 0. (Note that (c 1 , c 2 ) is not a seed pair, but the others are). Our horizontal path starts from this initial pair of vertices, crosses n ≥ 0 vertical paths (i.e., formed by three occupied vertical edges in a single column), and contains n + 1 uninterrupted horizontal strings of b 2 vertices of lengths k 0 , k 1 , . . . , k n ≥ 0 in between the vertical paths. The contribution from vertices in the region {x + 1, x + 2, . . . , 0} × {0, 1} to the probability of such a two-layer path configuration in Z × {0, 1} is (sδ 1 ) n ((1 − s)δ 2 ) |k| , where we denote |k| = k 0 + k 1 + . . . + k n . Note that this contribution is the same in the two cases when the horizontal path goes through the bottom or the top layer. The rate of the height change contains a term involving the vertex weights of the pair at (x, 0), (x, 1) times the rate of initiating jump at this pair (this term is equal to zero for (c 1 , c 2 ) because it is not a seed pair). Moreover, for up jumps, we also need to add a term accounting for |k| extra seed pairs along the horizontal path where a jump may also be initiated.
Overall, we obtain the following expression for the current: First one can compute the sum over k, using k=(k 0 ,k 1 ,...,kn)∈Z n+1 Employing these identities leaves only the summation over n, which is readily computed. After necessary simplifications, we arrive at the desired formula (7.3).

Heuristic hydrodynamics in the quadrant
This subsection presents a heuristic discussion of some hydrodynamic Burgers type equations in one and two space dimensions related to the stochastic six vertex model in the quadrant.
Recall the Markov dynamics Q(τ ) with the infinitesimal generator G quad u,η,τ (4.5). The dynamics Q(τ ) acts on path configurations in the quadrant Z ≥0 × Z ≥1 , with step-λ or empty-λ boundary conditions. Recall that the subset λ ⊂ Z ≥0 encodes the locations of incoming vertical arrows along the bottom boundary, and this subset stays fixed throughout the dynamics. By Theorem 4.9, Q(τ ) changes (in distribution) the Gibbs property of the stochastic six vertex model by continuously increasing the spectral parameter from u to some terminal value u + η ∈ (0, 1). Namely, we have P s6v, hom u Q(τ ) = P s6v, hom u+(1−e −τ ) η , where P s6v, hom u denotes the homogeneous stochastic six vertex model with spectral parameter u and our fixed step-λ or empty-λ boundary conditions. Denote the evolving spectral parameter by u(τ ) := u + (1 − e −τ )η. (7.4) Now consider the limit when the lattice coordinates are ( x/ε , y/ε ) for some (x, y) ∈ R 2 ≥0 , and ε 0. Assume that the subset λ = λ (ε) depends on ε and behaves regularly in the sense that its height function is for all ε and x ∈ R ≥0 , whereλ is a fixed nondecreasing function with slope ≤ 1. For each τ , consider the (random) height function h τ (k, l), where k ≥ 0, l ≥ 1, of the stochastic six vertex model P s6v, hom , which is defined as the (signed) number of up-right paths crossed between (k, l) and (0, 0). The height function increases by crossing a path going north or west, and decreases otherwise. See Figure 1 from the Introduction for an illustration.
From [Agg20a, Theorem 1.1] we know that the random height function h τ (k, l) admits a limit shape lim with convergence in probability. Here H(τ, x, y) is a nonrandom function with boundary conditions for all τ : H(τ, x, 0) =λ(x), H(τ, 0, y) = 0, empty-λ boundary conditions; y, step-λ boundary conditions. Also denote ρ(τ, x, y) := − ∂ ∂x H(τ, x, y), (7.5) this is the density of the occupied vertical edges near the global location (x, y). As our discussion in the current Section 7.2 is heuristic, we assume that the derivative (7.5) exists in a suitable sense (and similarly for all other derivatives below). There are two types of differential equations the function H(τ, x, y) should satisfy: • For each fixed τ , the density ρ(τ, x, y) should satisfy a version of the Burgers equation in (1+1) dimensions [Agg20a, Theorem 1.1]: This equation corresponds to the slice by slice evolution under the transfer matrix of the stochastic six vertex model. Here y plays the role of time, and ϕ(ρ | u(τ )) is the particle current in stationarity on Z at density ρ.
• The height function should satisfy: This equation corresponds to the fact that the H's are the limit shapes of the random height functions h τ (k, l). Indeed, the latter are obtained from h 0 (k, l) (by means of increasing τ ) using the Markov dynamics Q(τ ). Finally, the average velocity of the height function in the bulk around (x, y) under Q(τ ) is e −τ η y J(ρ, ϕ(ρ | u(τ ))), where ρ = ρ(τ, x, y) is the density of the occupied vertical edges, and the factor e −τ η y is due to the inhomogeneity of the edge Poisson clocks, cf. (4.5).
(7.13) Substituting the expression for ϕ (2.6) and using the Burgers equation (7.6) at time τ to express ρ y through ρ x , we see that (7.13) reduces tõ as desired.

Dynamics on the torus
In this section we define an analogue of the full plane dynamics C(t) constructed in Section 5 which acts on up-right path configurations on the torus. Our dynamics preserve the six vertex model Gibbs measures with arbitrary slopes (s, t). We present two proofs. The first immediately follows from the Yang-Baxter equation and its bijectivisation, and involves a certain discrete twist of the six vertex graph on the torus. For simplicity, in the first proof we only consider the stochastic six vertex weights for which we already have explicit jump rates. For the second proof, we mimic the argument of [BB17] involving symmetry of jump rates, and this allows to generalize our torus dynamics to arbitrary six vertex weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 .

Bijectivisation on the torus
Suppose we have an M × N torus Z/M Z × Z/N Z denoted by T = T M,N . Consider the set S k 1 ,k 2 of configurations of up-right paths on the torus with fixed overall height change k 1 and k 2 in the x and y directions, respectively. Let µ k 1 ,k 2 denote the Gibbs measure on S k 1 ,k 2 given by a choice of six vertex Boltzmann weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 (see Figure 3 for an illustration of the weights). We index the vertices by (x, y) where x ∈ {0, . . . , M − 1}, y ∈ {0, . . . , N − 1}. If x or y is outside of this range, we reduce these coordinates modulo M or N , respectively. First, we consider the special stochastic case a 1 = a 2 = 1, b 1 = 1 − c 1 = δ 1 , b 2 = 1 − c 2 = δ 2 , where δ 1 , δ 2 depend on q and u, see (2.2). Let us define a continuous time Markov chain on up-right path configurations on the torus. We employ the notion of seed pairs (Definition 4.5) and the jump rates R v 1 ,v 2 (u) from Figure 13. Definition 8.1 (Continuous time Markov dynamics on the torus). Each vertical edge (v 1 , v 2 ) which is a seed pair has an exponential clock with rate R v 1 ,v 2 (u). When the clock at (v 1 , v 2 ) rings, the horizontal path passing through this seed pair jumps up or down depending on the local path configuration around the edge v 1 − v 2 . This jump then instantaneously propagates to the right. The jump propagation may stop in two ways: • Either there exists a configuration to the right where the jump may stop in the same way as in the full plane dynamics, according to the rules described in Section 3.6. See Figure 11 for an illustration.
• Or there is no stopping configuration, and the jump propagation has to make a loop around the torus, leading to a jump of a full straight horizontal path. See Figure 18 for an illustration.
We denote by L(τ ) the continuous time Markov semigroup of thus defined process.
* Figure 18: An example of a jump where the propagation goes all the way around the torus, and the full straight horizontal path jumps up. The identification of the horizontal edges is also shown.
We will show that in the case of the stochastic six vertex weights, the Markov chain L(τ ) preserves the Gibbs measure µ k 1 ,k 2 for any k 1 , k 2 . We achieve this by constructing L(τ ) as a Poisson type continuous time limit of a discrete time Markov chain coming from the bijectivisation of the Yang-Baxter equation for the stochastic six vertex model defined in Section 3. In the torus case it turn out to be very convenient to define the discrete Markov chain not on the straight torus T, but on its suitably twisted version.
Definition 8.2. Let the (discretely) twisted torusT =T M,N be the graph displayed in Figure 19. We associate the spectral parameters u, u + ε, . . . , u + ε with the horizontal strands as shown in this figure, where 0 < u < 1 and 0 < ε < 1 − u.
Since the graphT is embedded into the usual torus, onT the height function is still welldefined (up to a constant) on the faces of the fundamental domain. The horizontal and the vertical height change along a generator of each homology class are also well defined. See Figure 20 for an example of a path configuration and the height function. Denote by S twist k 1 ,k 2 the six vertex configurations onT where the horizontal height change is k 1 , and the vertical height change is k 2 . Denote by µ ε,k 1 ,k 2 the Gibbs measure restricted to six vertex configurations in S twist k 1 ,k 2 , where we use stochastic vertex weights parameterized by the spectral parameters u on the bottom row and u+ε on every other row. For cross vertices, we take their weights equal to X u,u+ε = w u/(u+ε) (2.7), which are the weights entering the Yang-Baxter equation. Note that all these vertex weights oñ T are nonnegative. u u + ε u + ε u + ε u + ε u + ε u u + ε u + ε u + ε u + ε u + ε Figure 19: The twisted six vertex graph on the torus with M = 4 and N = 6. The identified edges are indicated in bold, and note that on the very left we also identify the pieces of a diagonal cross. The spectral parameters associated to the horizontal strands are u, u + ε, u + ε, . . . , u + ε. Definition 8.3. For s, s ∈ S twist k 1 ,k 2 we define the Markov transition probability L ε (s → s ) by starting from s and performing the following sequence of N random updates: • First, drag the cross vertex between rows 0 and 1 through the lattice until it is to the right of x = M − 1. Each step of dragging the cross is a random update coming from the bijectivisation of the Yang-Baxter equation, see Figure 8. After this, the spectral parameters on rows 0, . . . , N − 1 are u + ε, u, u + ε, . . . , u + ε.
• Then drag the cross between rows 1 and 2 through the lattice using bijectivisation until it is to the right of x = M −1. After this, the spectral parameters are u+ε, u+ε, u, u+ε, . . . , u+ε. . . .
• At the last step, drag the cross between rows N − 1 and 0 through the lattice. After this step, the lattice returns to the original state, and the spectral parameters are back to u, u + ε, u + ε, . . . , u + ε.
The following statement ensures that L ε is well-defined: Lemma 8.4. The random updates described in Definition 8.3 preserve S twist k 1 ,k 2 . Proof. At each step of dragging the cross the path configuration and the associated height function only change locally, so the height change over the whole torus is preserved.
Proof. This follows from the fact that each step of dragging the cross through the whole torus maps the current Gibbs measure into a Gibbs measure with swapped spectral parameters (see Propositions 3.5 and 3.7). After all N steps, the spectral parameters are back to the original sequence u, u + ε, u + ε, . . . , u + ε, and hence the Gibbs measure is preserved.
Let us now discuss the limit as ε → 0.
Proposition 8.6. In the limit as ε → 0, the Gibbs measure µ ε,k 1 ,k 2 on the twisted graph can be identified with the measure µ k 1 ,k 2 on the straight torus T.
Proof. As ε → 0, the cross vertex weights become X u,u = w 1 , and place zero weight onto the vertices (1, 0; 1, 0) and (0, 1; 0, 1), see (2.1). With this restriction on the cross vertex types, we may identify path configurations onT with those on T, see Figure 21 for an illustration. Therefore, the measure µ ε,k 1 ,k 2 for ε = 0 is determined only by the usual vertices and not the cross vertices, and thus coincides with µ k 1 ,k 2 . Figure 21: A path configuration under µ ε,k 1 ,k 2 for ε = 0. In this case, the state of a cross vertex is completely determined by the paths to the right of it. Cutting the cross vertices out and identifying the strands as shown by the dashed lines leads to the straight torus graph T.
Consider the Poisson type continuous time limit as ε → 0 of the iterated Markov transition operators L τ /ε ε , where τ ∈ R ≥0 . We know from Section 4 that muliple dragging of the crosses should be replaced by jumps of the horizontal paths initiated by Poisson clocks placed onto vertical edges. One readily sees that in this limit we have where L(τ ) is the Markov transition operator (over time τ ) from Definition 8.1. This convergence together with Propositions 8.5 and 8.6 implies the following result: Theorem 8.7. The continuous time Markov process L(τ ) on stochastic six vertex configurations on the torus T preserves the measure µ k 1 ,k 2 for arbitrary k 1 , k 2 .
Remark 8.8 (Comparison with a similar process from [BB17]). In [BB17], the authors introduce a continuous time Markov process (denote it byL) preserving the measure µ k 1 ,k 2 on the torus. RotatingL by π/2 counterclockwise and reflecting along the x direction makes the jumps inL move the horizontal paths by jumps triggered by vertical edges. That is, we may describe both L andL in similar terms. Moreover, after a scalar time renormalization, some of the jumps and their jump rates exactly coincide in L andL. However, this transformation does not make all the rates in both processes equal, which shows that the two processes are not the same.
The processes L andL share another common feature, namely, that one can prove the preservation of the measure µ k 1 ,k 2 under both using a certain symmetry of the jump rates. We present such an argument for our processes L in Section 8.2 below.

Dynamics on the torus for general six vertex model
Up to an overall constant η > 0, the jump rates a(u), b(u), and c(u) (4.1) in the Markov chain L(τ ) from Definition 8.1 can be written in terms of the general six vertex weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 (recall Figure 3) as follows: Extending the definition, let L(τ ) denote the Markov process on configurations on the torus with the jump rates depending on the generic six vertex weights a 1 , a 2 , b 1 , b 2 , c 1 , c 2 as in (8.1). Similarly, let µ k 1 ,k 2 be the six vertex model on the torus T determined by these generic vertex weights (and horizontal and vertical height changes k 1 , k 2 ). We can extend Theorem 8.7 to the general weights.
Proof. For brevity, we will not reproduce here the details from [BB17], and simply follow the argument and notation of that paper. The main ingredient is to check that [BB17, Lemma 2] (a symmetry of the jump rates under the flip transformation) applies to our Markov process L(τ ). This is indeed the case, as can be seen by a direct inspection of all the jump rates. Then the preservation of µ k 1 ,k 2 under L(τ ) follows similarly to [BB17,Theorem 5]. We need to show s µ k 1 ,k 2 (s)Rate(s → s 0 ) − µ k 1 ,k 2 (s 0 ) s 2 Rate(s 0 → s 2 ) = µ k 1 ,k 2 (s 0 ) s Rate(s 0 → s) − Rate(s 0 → s) .

Degeneration to five vertex model and lozenge tilings
Let us set the weight a 2 of the vertex (1, 1; 1, 1) to zero. This turns the six vertex model into the five vertex model, which may be viewed as a certain model of nonintersecting (but interacting) paths, or, equivalently, lozenge tilings on the triangular lattice (with interacting lozenges). The five vertex model admits a more detailed asymptotic analysis than the general six vertex one by means of the Bethe Ansatz, for example, see the recent work [dGKW21].
Further letting b 1 b 2 = c 1 c 2 makes the five vertex model free fermion by removing the interaction. The free fermion five vertex model is equivalent to a dimer model, and may be analyzed asymptotically through determinantal point processes, see, for example, [ABPW21].
Let us consider the degeneration of our Markov dynamics on the torus at a 2 = 0. Under a suitable renormalization, the rates (8.1) for a 2 = 0 reduce to Because a 2 = 0 and b = 0, out of six possible seed pairs in Figure 13 only two may lead to a jump. Both these seed pairs correspond to up jumps, so our Markov process becomes totally asymmetric. By Theorem 8.9, the dynamics on the torus with rates (8.2) preserves the five vertex model (in particular, the one considered in [dGKW21]). Determining the particle current of the dynamics preserving the five vertex model could be simpler than in the general six vertex case, but this is outside the scope of the present work.
In the free fermion case b 1 b 2 = c 1 c 2 we see that c = a. After mapping nonintersecting paths of the free fermion five vertex model to lozenge tilings, one readily sees that our dynamics reduces to the totally asymmetric case of the interacting Hammersley processes studied in [Ton17], [CF17] [CFT19].