# A reverse duality for the ASEP with open boundaries

G.M. Schütz

###### Abstract

We prove a duality between the asymmetric simple exclusion process (ASEP) with
non-conservative open boundary conditions and an asymmetric exclusion process
with particle-dependent hopping rates and conservative reflecting boundaries.
This is a reverse duality in the sense that the duality function relates the
measures of the dual processes rather than expectations. Specifically, for a
certain parameter manifold of the boundary parameters of the open ASEP this
duality expresses the time evolution of a family of shock product measures
with $N$ microscopic shocks in terms of the time evolution of $N$ particles in
the dual process. The reverse duality also elucidates some so far poorly
understood properties of the stationary matrix product measures of the open
ASEP given by finite-dimensional matrices.
AMS 2010 subject classifications: 82C22. Secondary: 82C23, 82C26, 60K35

Departamento de Matemática, Instituto Superior Técnico, Universidade de
Lisboa, Av. Rovisco Pais, 1049-001 Lisbon, Portugal

##  1 Introduction

Duality expresses time-dependent expectations for one Markov process in terms
of time-dependent expectations of another Markov process and thus is an
important tool for the study of stochastic interacting particle systems [50,
51, 40], particularly if the dual process has some simple properties that are
“hidden” in the original process. The key object is the duality function $D$:
Processes $\boldsymbol{\eta}(t)$ with state space $\Omega$ and
$\boldsymbol{\xi}(t)$ with state space $\Xi$ respectively are said to be dual
w.r.t. a bounded measurable function $D:\Omega\times\Xi\to{\mathbb{R}}$ if and
only if for all $\boldsymbol{\eta}\in\Omega,\boldsymbol{\xi}\in\Xi$ and $t\geq
0$ one has
${\mathbf{E}}_{\boldsymbol{\eta}}D(\boldsymbol{\eta}(t),\boldsymbol{\xi})={\mathbf{E}}_{\boldsymbol{\xi}}D(\boldsymbol{\eta},\boldsymbol{\xi}(t))$.
The subscript denotes expectation with the initial states
$\boldsymbol{\eta}=\boldsymbol{\eta}(0)$ and
$\boldsymbol{\xi}=\boldsymbol{\xi}(0)$ resp. Considering the duality function
$D(\boldsymbol{\eta},\boldsymbol{\xi})$ to represent a family of functions
$f^{\boldsymbol{\xi}}(\boldsymbol{\eta})=D(\boldsymbol{\eta},\boldsymbol{\xi})$
indexed by $\boldsymbol{\xi}$ then duality expresses for all members of this
family the expectation of this function at time $t$ in terms of the
expectation at time $t$ of the family of functions
$g^{\boldsymbol{\eta}}(\boldsymbol{\xi})=D(\boldsymbol{\eta},\boldsymbol{\xi})$
indexed by $\boldsymbol{\eta}$ of the dual process.

In the last decade a very considerable body of work on applications of duality
to interacting particle systems defined on the one-dimensional integer lattice
has emerged see e.g. [10, 45, 52, 11, 18, 16, 46, 47] and references therein
specifically for very recent work on duality for processes closely related to
the ASEP which is studied in the present paper. To push this field further we
introduce the concept of reverse duality which expresses the time evolution of
probability measures of one process in terms of the dual process rather than
relating time-dependent expectations. Reverse duality can, however, be
rephrased as a conventional duality for the time-reversed processes and thus
extends the framework of conventional duality. We use it to derive a duality
that has so far remained elusive, viz., for the asymmetric simple exclusion
process (ASEP) on a one-dimensional integer lattice of $L$ sites with two non-
conservative open boundaries [49, 44, 25, 63]. On a certain manifold of
boundary rates for injection and absorption the open ASEP is shown to be
reverse-dual to an ASEP with particle-dependent jump rates on the same lattice
but with two conservative reflecting boundaries. Thus we go beyond a recent
proof of self-duality for one open and one reflective boundary condition [46,
47] and an earlier failure to construct a non-trivial duality for two open
boundaries [57].

The reverse duality for the open ASEP constructed here allows for expressing
the time evolution of a certain family of Bernoulli shock measures with $N$
shocks in terms of the evolution of the probability measure of the dual ASEP
with $N$ particles and thus reduces the complexity of the time evolution of
the system with $2^{L}$ states to $L\choose N$ states. Specifically, the time
evolution of a single shock in the many-particle ASEP turns into the simple
problem of a single biased random walker. The evolution of initial measures
with more than one shock can be studied through the complete integrability
[37, 27] which provides an explicit system of Bethe equations for the spectrum
of the generator [55, 56, 22, 23, 70] and for which results on the
eigenfunctions are available [21]. Along the way one also obtains insight into
some of the properties of the stationary matrix product measures of the open
ASEP introduced in [25] which on the parameter manifold studied here is given
by finite-dimensional matrices [54].

The paper is organized as follows. In Sec. 2 we introduce reverse duality in
general terms, and define the ASEP with open boundaries, shock measures for
the ASEP (which give rise to the duality function), and the dual shock ASEP
with particle-dependent hopping rates (which turns out to the reverse dual).
The results, in particular the reverse duality, are stated in precise form as
theorems in Sec. 3, with comments on the significance of these results. The
theorems are formally proved in Sec. 4. In Appendix A we show how the
generator of the open ASEP in the parametrization of the present work is
related to the integrable Heisenberg quantum spin chain with non-hermitian
boundary fields of [55, 56] and in Appendix B we discuss some properties of
the boundary parameters.

##  2 Definitions

###  2.1 General definitions and conventions

We consider the finite integer lattice
$\Lambda_{L}=[L_{-},L_{+}]\cap{\mathbb{Z}}$ with $L=L_{+}-L_{-}+1\geq 2$
sites. The Kronecker-$\delta$ is defined by

$\delta_{x,y}:=\left\\{\begin{array}[]{ll}1&\mbox{ if }x=y\\\ 0&\mbox{ else
}\end{array}\right.$ (2.1)

for $x,y$ from any set.

For a finite or countably infinite vector space column vectors with components
$\phi_{i}$ are denoted by the quantum mechanical ket symbol
$|\,{\phi}\,\rangle$ while row vectors denoted by the bra symbol
$\langle\,{\phi}\,|$ have complex conjugate components $\bar{\phi}_{i}$. The
standard scalar product is written in the Dirac convention

$\mbox{$\langle\,{\phi}\,|\,{\psi}\,\rangle$}=\sum_{i}\bar{\phi}_{i}\psi_{i}$
(2.2)

Notice the complex conjugation of the left argument.

######  Definition 2.1 (Kronecker product)

Let $A$ and $B$ be two finite-dimensional matrices with $m_{A}\geq 1$
($m_{B}\geq 1$) rows and $n_{A}\geq 1$ ($n_{B}\geq 1$) columns with matrix
elements $A_{ij}$ and $B_{kl}$ respectively. The Kronecker product $A\otimes
B$ is a $m_{A}m_{B}\times n_{A}n_{B}$-matrix $C$ with matrix elements

$C_{(i-1)m_{B}+k,(j-1)n_{B}+l}=A_{ij}B_{kl}$ (2.3)

with
$i\in\\{1,\dots,m_{A}\\},j\in\\{1,\dots,n_{A}\\},k\in\\{1,\dots,m_{B}\\},l\in\\{1,\dots,n_{B}\\}$.
The $n$-fold Kronecker product of a matrix $A$ with itself is denoted by
$A^{\otimes n}$ and by convention $A^{\otimes 0}=1$ and $A^{\otimes 1}=A$.

######  Definition 2.2 (Local and bond operators)

(i) With the two-dimensional unit matrix denoted by $\mathds{1}$ and an
arbitrary two-dimensional matrix $\mathsf{a}$ the Kronecker product

$\mathsf{a}_{k}:=\mathds{1}^{\otimes(k-L_{-})}\otimes\mathsf{a}\otimes\mathds{1}^{\otimes(L_{+}-k)},\quad
k\in\Lambda_{L}$ (2.4)

is called a local operator acting on site $k\in\Lambda_{L}$.
(ii) For any four-dimensional matrix $\mathsf{b}$, bond operators
$\mathsf{b}_{k,k+1}$ with $L_{-}\leq k\leq L_{+}-1$ are defined by the
Kronecker product

$\mathsf{b}_{k,k+1}=\mathds{1}^{\otimes(k-L_{-})}\otimes\mathsf{b}\otimes\mathds{1}^{\otimes(L_{+}-1-k)}.$
(2.5)

###  2.2 Duality and reverse duality

For many applications it is convenient to express the generator of a Markov
process $\boldsymbol{\eta}(t)$ with countable state space in terms of the
intensity matrix $W$ whose off-diagonal matrix elements
$W_{\boldsymbol{\eta},\boldsymbol{\eta}^{\prime}}$ are the transition rates
from configuration $\boldsymbol{\eta}$ to $\boldsymbol{\eta}^{\prime}$ and
with diagonal elements given by the negative sum
$W_{\boldsymbol{\eta}\boldsymbol{\eta}}=-\sum_{\boldsymbol{\eta}^{\prime}}W_{\boldsymbol{\eta},\boldsymbol{\eta}^{\prime}}$
of outgoing transition rates from state $\boldsymbol{\eta}$. In terms of
intensity matrices $W$ for the process $\boldsymbol{\eta}(t)$ and $Q$ for the
process $\boldsymbol{\xi}(t)$, and with the duality matrix $D$ with the matrix
elements $D_{\boldsymbol{\eta}\xi}=D(\boldsymbol{\eta},\boldsymbol{\xi})$
duality is expressed in matrix form as the relation [72]

$WD=DQ^{T}.$ (2.6)

This notion of duality has the following natural generalization.

###### Definition 2.3

Let $\boldsymbol{\eta}(t)$ be a Markov process with intensity matrix $W$. A
Markov process $\boldsymbol{\xi}(t)$ with intensity matrix $Q$ satisfying

$RW=Q^{T}R$ (2.7)

with a duality matrix $R$ with matrix elements
$R_{\boldsymbol{\xi}\boldsymbol{\eta}}$ is called reverse dual to
$\boldsymbol{\eta}(t)$ w.r.t. the duality function
$R(\boldsymbol{\xi},\boldsymbol{\eta})=R_{\boldsymbol{\xi}\boldsymbol{\eta}}$.

Obviously, reverse duality is a mutual relation between two processes as can
be seen by transposing both sides of the defining equation (2.7) which then
says that $\boldsymbol{\eta}(t)$ is a reverse dual to $\boldsymbol{\xi}(t)$
with duality matrix $R^{T}$. To give some motivation for the notion of revsere
duality assume that the processes $\boldsymbol{\eta}(t)$ and
$\boldsymbol{\xi}(t)$ have invertible invariant measures $\pi_{W}$ and
$\pi_{Q}$ resp., i.e., $\pi_{W}\neq 0$ for all $\boldsymbol{\eta}$ and
$\pi_{Q}\neq 0$ for all $\boldsymbol{\xi}$. This implies that that the
diagonal matrices $\hat{\pi}_{W}^{-1}$ with diagonal entries
$\pi_{W}^{-1}(\boldsymbol{\eta})$ and $\hat{\pi}_{Q}^{-1}$ with diagonal
entries $\pi_{Q}^{-1}(\boldsymbol{\xi})$ exist. By definition, the intensity
matrices for the reversed processes are then given by
$W_{rev}=\hat{\pi}_{W}^{-1}W^{T}\hat{\pi}_{W}$ and
$Q_{rev}=\hat{\pi}_{Q}^{-1}Q^{T}\hat{\pi}_{Q}$, resp., so that
$W^{T}=\hat{\pi}_{W}W_{rev}\hat{\pi}_{W}^{-1}$ and
$Q=\hat{\pi}_{Q}^{-1}Q_{rev}^{T}\hat{\pi}_{Q}$. After transposition and with
the duality matrix $D=\hat{\pi}_{W}^{-1}R^{T}\hat{\pi}_{Q}^{-1}$ the defining
relation (2.7) for reverse duality turns into the conventional duality
relation $W_{rev}D=DQ_{rev}^{T}$ for the reversed processes. In other words,
for reversible processes reverse duality implies a conventional duality for
the reversed processes and vice versa.

While conventional duality relates expectation values of the dual processes,
reverse duality can be used to relate the measures for all times $t\geq 0$. To
see this, we denote by $\mu(t)$ the probability of the configuration
$\boldsymbol{\eta}(t)$ for some unspecified initial measure $\mu(0)$ and for
the process $\boldsymbol{\xi}(t)$ we denote by
$P(\boldsymbol{\xi},t|\boldsymbol{\xi}^{\prime},0)$ the transition probability
at time $t$ from a configuration $\boldsymbol{\xi}^{\prime}$ to a
configuration $\boldsymbol{\xi}$.

###### Theorem 2.4

For a countable set $\Xi$ let $\mathcal{F}_{\Xi}$ be a family of probability
measures for the process $\boldsymbol{\eta}(t)$ with countable state space
$\Omega$ and let the process $\boldsymbol{\xi}(t)$ with state space $\Xi$ be
reverse dual to $\boldsymbol{\eta}(t)$ with duality function
$R(\boldsymbol{\xi},\boldsymbol{\eta})=\mu^{\boldsymbol{\xi}}_{\boldsymbol{\eta}}\in\mathcal{F}_{\Xi}$.
Then at time $t\geq 0$ the measure $\mu^{\boldsymbol{\xi}}(t)$ is given by the
convex combination

$\mu^{\boldsymbol{\xi}}(t)=\sum_{\boldsymbol{\xi}^{\prime}}P(\boldsymbol{\xi},t|\boldsymbol{\xi}^{\prime},0)\mu^{\boldsymbol{\xi}^{\prime}}(0)$
(2.8)

of initial measures $\mu^{\boldsymbol{\xi}^{\prime}}_{0}\in\mathcal{F}_{\Xi}$.

Proof: The time evolution property (2.8) of the measure follows directly from
the semigroup property of the Markov processes $\boldsymbol{\eta}(t)$ and
$\boldsymbol{\xi}(t)$ and Definition 2.3 of reverse duality. $\Box$

Thus if the duality function $R(\boldsymbol{\xi},\boldsymbol{\eta})$
represents a family of probability measures
$\mu^{\boldsymbol{\xi}}_{\boldsymbol{\eta}}$ indexed by $\boldsymbol{\xi}$
then reverse duality expresses the Markovian time evolution of all members of
this family of measures in terms of the probability measure of the reverse
dual process starting from configuation $\boldsymbol{\xi}$.111Duality between
measures is discussed in a different setting also in [19]. In contrast,
considering the duality function of conventional duality to be member of a
measurable family of functions
$f^{\boldsymbol{\xi}}(\boldsymbol{\eta}):=D(\boldsymbol{\eta},\boldsymbol{\xi})$
of the states $\boldsymbol{\eta}$ and denoting the expecation of this function
at time $t\geq 0$ w.r.t. some unspecified initial measure by
${\mathbf{E}}f^{\boldsymbol{\xi}}(t)$, conventional duality reads

${\mathbf{E}}f^{\boldsymbol{\xi}}(t)=\sum_{\boldsymbol{\xi}^{\prime}}P(\boldsymbol{\xi},t|\boldsymbol{\xi}^{\prime},0){\mathbf{E}}f^{\boldsymbol{\xi}^{\prime}}(0).$
(2.9)

This highlights the difference between reverse duality and conventional
duality in terms of the role of the transition probability of the (reverse)
dual process.

For reversible processes $\boldsymbol{\xi}(t)$ with reversible measure $\pi$
Theorem 2.4 can be stated in a slightly modified form.

###### Theorem 2.5

For a countable set $\Xi$ let $\mathcal{F}_{\Xi}$ be a family of probability
measures for the process $\boldsymbol{\eta}(t)$ with countable state space
$\Omega$ and let the process $\boldsymbol{\xi}(t)$ with state space $\Xi$ be
(a) reversible with reversible measure $\pi$ and (b) reverse dual to
$\boldsymbol{\eta}(t)$ with duality function
$R(\boldsymbol{\xi},\boldsymbol{\eta})=\pi(\boldsymbol{\xi})\mu^{\boldsymbol{\xi}}_{\boldsymbol{\eta}}\in\mathcal{F}_{\Xi}$.
Then at time $t\geq 0$ the measure $\mu^{\boldsymbol{\xi}}(t)$ is given by the
convex combination

$\mu^{\boldsymbol{\xi}}(t)=\sum_{\boldsymbol{\xi}^{\prime}}P(\boldsymbol{\xi}^{\prime},t|\boldsymbol{\xi},0)\mu^{\boldsymbol{\xi}^{\prime}}(0)$
(2.10)

of initial measures $\mu^{\boldsymbol{\xi}^{\prime}}(0)\in\mathcal{F}_{\Xi}$.

Proof: With $S:=\hat{\pi}^{-1}R$ reverse duality as defined in (2.7) implies
for $Q=Q^{rev}$ the intertwining relation

$SW=QS.$ (2.11)

The time evolution property (2.10) then follows from the semigroup property of
the Markov processes $\boldsymbol{\eta}(t)$ and $\boldsymbol{\xi}(t)$ and
reversibility of $\boldsymbol{\xi}(t)$. $\Box$

###### Remark 2.6

The difference to Theorem 2.4 appears in the time-reversal of the transition
probability of the reverse dual process and in the relationship between the
duality function $R(\boldsymbol{\xi},\boldsymbol{\eta})$ and the measures
$\mu^{\boldsymbol{\xi}}$.

###### Remark 2.7

For self-dual processes with $W=Q$ the intertwining relation (2.11) expresses
a symmetry of the generator [64, 35, 15] and for invertible $S$ it constitutes
a similarity transformation between processes with intensity matrices $W,Q$
[36].

###  2.3 Intensity matrices and quantum Hamiltonians

We recall that by construction the sum of all matrix elements in each row of
the intensity matrices add up to zero. It is well-known that the transposed
intensity matrix often takes the form of a Hamiltonian operator $H=-W^{T}$ for
quantum mechanical spin chains [1, 53, 68] and concepts from quantum mechanics
can be imported to solve probabilistic problems even though in general $H$ is
not a Hermitian matrix. To illustrate the point in the context of duality we
mention that an invariant measure is a ground state vector of the quantum
system and non-trivial self-dualities arise from non-Abelian symmetries of the
generator [65, 66, 35, 15]. An general relation between duality and spectral
properties of the generator - which play a crucial role in quantum mechanics
and can be computed exactly for integrable systems - was recently worked out
in [58]. Complete integrability can also be used to directly construct novel
dualities [69].

With this in mind we express below duality and reverse duality in terms of the
negative transpose $H$ of the intensity matrix and in slight abuse of language
we also call $H$ the generator of the process. Time-invariance of a measure
$\pi$ means that the vector $|\,{\pi}\,\rangle$ is a right eigenvector of $H$
with eigenvalue $0$. Defining the diagonal matrix $\hat{\pi}$ with the
components $\pi(\cdot)$ on the diagonal, the generator of the time-reversed
process is given by

$H_{rev}:=\hat{\pi}H^{T}\hat{\pi}^{-1}$ (2.12)

and has the same invariant measure $\pi$. With $H=-W^{T}$ and $G=-Q^{T}$
duality as defined in (2.6) then reads

$D^{T}H=G^{T}D^{T}$ (2.13)

with transposed duality matrix $D^{T}$ with matrix elements
$D^{T}_{\boldsymbol{\xi}\boldsymbol{\eta}}=D(\boldsymbol{\eta},\boldsymbol{\xi})$.
Reversed duality as defined in (2.7) becomes

$HR^{T}=R^{T}G^{T}$ (2.14)

and for reversible $G$ this yields the intertwining relation

$HS^{T}=S^{T}G$ (2.15)

with $S^{T}=R^{T}\hat{\pi}^{-1}$.

To rephrase Theorem 2.5 we use quantum mechanical bra-ket notation for vectors
and the scalar product (2.2). Then a probability measure for the process
$\boldsymbol{\eta}(t)$ is a ket vector $|\,{\mu(t)}\,\rangle$ whose components
$\mu_{\boldsymbol{\eta}}(t)$ are the probabilities of configuration
$\boldsymbol{\eta}(t)$. The Markov semigroup property yields
$\mbox{$|\,{\mu(t)}\,\rangle$}=[\mbox{$\langle\,{\mu(0)}\,|$}\exp{(Wt)}]^{T}=\exp{(-Ht)}\mbox{$|\,{\mu(0)}\,\rangle$}$
so that the transition probability
$\mu(\boldsymbol{\eta}^{\prime},t|\boldsymbol{\eta},0)$ from a configuration
$\boldsymbol{\eta}$ to a configuration $\boldsymbol{\eta}^{\prime}$ is given
by

$\mu(\boldsymbol{\eta}^{\prime},t|\boldsymbol{\eta},0)=\mbox{$\langle\,{\boldsymbol{\eta}}\,|$}\mathrm{e}^{Wt}\mbox{$|\,{\boldsymbol{\eta}^{\prime}}\,\rangle$}=\mbox{$\langle\,{\boldsymbol{\eta}^{\prime}}\,|$}\mathrm{e}^{-Ht}\mbox{$|\,{\boldsymbol{\eta}}\,\rangle$}.$
(2.16)

and similarly for the process $\boldsymbol{\xi}(t)$ generated by $G$. The
evolution formula (2.10) of Theorem 2.5 reads in vector form

$\mbox{$|\,{\mu^{\boldsymbol{\xi}}(t)}\,\rangle$}=\sum_{\boldsymbol{\xi}^{\prime}}P(\boldsymbol{\xi}^{\prime},t|\boldsymbol{\xi},0)\mbox{$|\,{\mu^{\boldsymbol{\xi}^{\prime}}(0)}\,\rangle$}.$
(2.17)

with the transition probability
$P(\boldsymbol{\xi}^{\prime},t|\boldsymbol{\xi},0)$ of the dual process.

###### Remark 2.8

Under the conditions of Theorem 2.5 the family $\mathcal{F}_{\Xi}$ of
probability vectors belongs to a right-invariant subspace of the generator $H$
with dimension $d_{\Xi}\leq|\mathcal{F}_{\Xi}|$.

###  2.4 ASEP with open boundaries

####  2.4.1 Generator

In the ASEP each lattice site $k\in\Lambda_{L}$ can be occupied or empty as
expressed by the occupation numbers $\eta_{k}\in\\{0,1\\}$. The particle
configuration on the whole lattice $\lambda_{L}$ is denoted by
$\boldsymbol{\eta}:=\\{\eta_{L_{-}},\eta_{L_{-}+1},\dots,\eta_{L_{+}}\\}$. For
a given configuration $\boldsymbol{\eta}$ we also define the flipped
configuration $\boldsymbol{\eta}^{k}$ with occupation numbers

$\eta^{k}_{l}=\left\\{\begin{array}[]{ll}1-\eta_{k}&\mbox{if }l=k\\\
\eta_{l}&\mbox{else}.\end{array}\right.$ (2.18)

and the swapped configuration $\boldsymbol{\eta}^{kk+1}$ with occupation
numbers

$\eta^{kk+1}_{l}=\left\\{\begin{array}[]{ll}\eta_{k+1}&\mbox{if }l=k\\\
\eta_{k}&\mbox{if }l=k+1\\\ \eta_{l}&\mbox{else}.\end{array}\right.$ (2.19)

Expectations are denoted by $\langle\,{\cdot}\,\rangle$. In particular,
$\langle\,{\eta_{k}}\,\rangle$ is called the local density.

Informally, the Markovian dynamics of the open ASEP can be described as
follows. A particle on a bulk site $k\in[L_{-}+1,L_{+}-1]$ jumps with rate $r$
($\ell$) independently of the other particles to right (left) neighbouring
site $k+1$ ($k-1$), provided the target site is empty. Otherwise the jump
attempt is rejected. On the left boundary site $L_{-}$ a particle jumps with
rate $r$ to the right (if site $L_{-}+1$ is empty) or is annihilated with rate
$\gamma$. If site $L_{-}$ is empty then a particle is created with rate
$\alpha$. On the right boundary site $L_{+}$ a particle jumps with rate $\ell$
to the left (if site $L_{+}-1$ is empty) or is annihilated with rate $\beta$.
If site $L_{+}$ is empty then a particle is created with rate $\delta$.

With the local jump rate

$w_{k}(\boldsymbol{\eta})=r\eta_{k}(1-\eta_{k+1})+\ell(1-\eta_{k})\eta_{k+1}$
(2.20)

and the boundary rates

$b^{-}(\boldsymbol{\eta})=\alpha(1-\eta_{L_{-}})+\gamma\eta_{L_{-}},\quad
b^{+}(\boldsymbol{\eta})=\delta(1-\eta_{L_{+}})+\beta\eta_{L_{+}}$ (2.21)

the generator $\mathcal{L}$ of this process is thus given by

$\mathcal{L}=\mathcal{L}^{bulk}+\mathcal{L}^{-}+\mathcal{L}^{+}$ (2.22)

with the bulk generator

$\mathcal{L}^{bulk}f(\boldsymbol{\eta})=\sum_{k=L_{-}}^{L_{+}-1}w_{k}(\boldsymbol{\eta})[f(\boldsymbol{\eta}^{kk+1})-f(\boldsymbol{\eta})]$
(2.23)

and the boundary generators

$\mathcal{L}^{-}f(\boldsymbol{\eta})=b^{-}(\boldsymbol{\eta})[f(\boldsymbol{\eta}^{L_{-}})-f(\boldsymbol{\eta})],\quad\mathcal{L}^{+}f(\boldsymbol{\eta})=b^{+}(\boldsymbol{\eta})[f(\boldsymbol{\eta}^{L_{+}})-f(\boldsymbol{\eta})].$
(2.24)

Throughout this work it is assumed that all rates
$r,\ell,\alpha,\beta,\gamma,\delta$ are strictly positive and that
$r\neq\ell$.

To write the corresponding intensity matrix we denote the two-dimensional unit
matrix by $\mathds{1}$ and define the matrices

$\displaystyle w^{b}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}0&0&0&0\\\ 0&-\ell&\ell&0\\\
0&r&-r&0\\\ 0&0&0&0\end{array}\right)$ (2.29) $\displaystyle w^{-}$
$\displaystyle=$ $\displaystyle\left(\begin{array}[]{cc}-\alpha&\alpha\\\
\gamma&-\gamma\end{array}\right)$ (2.32) $\displaystyle w^{+}$
$\displaystyle=$ $\displaystyle\left(\begin{array}[]{cc}-\delta&\delta\\\
\beta&-\beta\end{array}\right).$ (2.35)

and the Kronecker products $w^{-}_{L_{-}}:=w^{-}\otimes\mathds{1}^{\otimes
L_{+}-L_{-}}$, $w^{+}_{L_{+}}:=\mathds{1}^{\otimes L_{+}-L_{-}}\otimes w^{+}$,
and $w^{b}_{k,k+1}:=\mathds{1}^{\otimes-L_{-}+k}\otimes
w^{b}\otimes\mathds{1}^{\otimes L_{+}-1-k}$. In terms of these matrices the
intensity matrix for the open ASEP is given by

$W=W_{bulk}+W_{edge}$ (2.36)

with

$W_{bulk}=\sum_{k=L_{-}}^{L_{+}-1}w^{b}_{k,k+1},\quad
W_{edge}=w^{-}_{L_{-}}+w^{+}_{L_{+}}.$ (2.37)

In $H:=-W^{T}$ one recognizes a non-hermitian version of the XXZ quantum spin
chain with non-diagonal boundary fields [55, 56], see the appendix.

To make contact with the well-known properties of the XXZ chain, the
computations below are done using $H$ rather than $W$. Defining the particle
number operator

$\hat{n}:=\left(\begin{array}[]{cc}0&0\\\ 0&1\end{array}\right)$ (2.38)

and the Kronecker product $\hat{n}_{k}:=\mathds{1}^{\otimes-
L_{-}+k}\otimes\hat{n}\otimes\mathds{1}^{\otimes L_{+}-k}$ one gets

$H=H_{bulk}+H_{edge}$ (2.39)

with

$H_{bulk}=\sum_{k=1}^{L-1}\tilde{h}_{k,k+1},\quad
H_{edge}=\tilde{h}^{-}_{L_{-}}+\tilde{h}^{+}_{L_{+}}$ (2.40)

where

$\displaystyle\tilde{h}_{k,k+1}$ $\displaystyle=$
$\displaystyle-(w^{b}_{k,k+1})^{T}+(r-\ell)\left(\hat{n}_{k+1}-\hat{n}_{k}\right)\,=\,\left(\begin{array}[]{cccc}0&0&0&0\\\
0&r&-r&0\\\ 0&-\ell&\ell&0\\\ 0&0&0&0\end{array}\right)_{k,k+1},$ (2.45)
$\displaystyle\tilde{h}^{-}_{L_{-}}$ $\displaystyle=$
$\displaystyle-(w^{-}_{L_{-}})^{T}+(r-\ell)\hat{n}_{L_{-}}\,=\,\left(\begin{array}[]{cc}\alpha&-\gamma\\\
-\alpha&\gamma+r-\ell\end{array}\right)_{L_{-}},$ (2.48)
$\displaystyle\tilde{h}^{+}_{L_{+}}$ $\displaystyle=$
$\displaystyle-(w^{+}_{L_{+}})^{T}+(r-\ell)\hat{n}_{L_{+}}\,=\,\left(\begin{array}[]{cc}\delta&-\beta\\\
-\delta&\beta-(r-\ell)\end{array}\right)_{L_{+}}.$ (2.51)

Notice the addition of the discrete gradient
$(r-\ell)\left(\hat{n}_{k+1}-\hat{n}_{k}\right)$ to the local bulk generators
which cancel in (2.39) with the corresponding contributions to the boundary
generators in (2.48) and (2.51).

####  2.4.2 Particle current

Particle conservation in the bulk is reflected in the discrete continuity
equation

$\mathcal{L}\eta_{k}=j_{k-1}-j_{k}$ (2.52)

where

$j_{k}=r\eta_{k}(1-\eta_{k+1})-\ell(1-\eta_{k})\eta_{k+1},\quad L_{-}\leq
k<L_{+}$ (2.53)

is the instantaneous current across bond $(k,k+1)$. At the boundaries one has

$\mathcal{L}\eta_{L_{-}}=j^{-}-j_{L_{-}},\quad\mathcal{L}\eta_{L_{+}}=j_{L_{+}-1}-j^{+}$
(2.54)

where

$j^{-}=\alpha(1-\eta_{L_{-}})-\gamma\eta_{L_{-}},\quad
j^{+}=\beta\eta_{L_{+}}-\delta(1-\eta_{L_{-}})$ (2.55)

are the boundary currents with a source term a sink term each. Stationarity
implies constant expectation

$j:=\mbox{$\langle\,{j^{-}}\,\rangle$}=\mbox{$\langle\,{j_{k}}\,\rangle$}=\mbox{$\langle\,{j^{+}}\,\rangle$}$
(2.56)

in the invariant measure which was first computed for the TASEP
($\gamma=\delta=\ell=0$) in two different ways, viz., with a matrix product
ansatz in [25] and simultaneously by generating function techniques in [63].
The general case was solved subsequently in [61].

####  2.4.3 Invariant matrix product measures

To discuss some properties of the invariant measure it is convenient to
introduce the jump asymmetry $q$ and the time scale $w$ by

$q:=\sqrt{\frac{r}{\ell}},\quad w:=\sqrt{r\ell}$ (2.57)

and to parametrize the boundary rates as

$\displaystyle\alpha=(r+\omega_{-})\rho_{-},\quad\gamma=(\ell+\omega_{-})(1-\rho_{-})$
(2.58)
$\displaystyle\beta=(r+\omega_{+})(1-\rho_{+}),\quad\delta=(\ell+\omega_{+})\rho_{+}$
(2.59)

where the parameters $\omega_{\pm}$ may be interpreted as boundary jump
barriers and, as seen below, the parameters $\rho_{\pm}$ play the role of
boundary densities. It is convenient to also consider the fugacity

$z\equiv z(\rho)=\frac{\rho}{1-\rho}$ (2.60)

as a function of the density $\rho$. For an indexed density $\rho_{i}$ with
$i$ from any index set we shall use the notation

$z_{i}\equiv z(\rho_{i})$ (2.61)

as shorthand.

The invariant measure of the open ASEP is conveniently characterized in terms
of the functions [61]

$\kappa_{\pm}(x,y):=\frac{1}{2x}(y-x+r-\ell\pm\sqrt{(y-x+r-\ell))^{2}+4xy})$
(2.62)

which are the roots of the quadratic equation

$x\kappa^{2}-(y-x+r-\ell)\kappa-y=0$ (2.63)

and which are related by

$\kappa_{+}(x,y)\kappa_{-}(x,y)=-\frac{y}{x}.$ (2.64)

As shown in [54] on the parameter manifold $\mathcal{B}_{N}$ defined by the
relation

$\kappa_{+}(\alpha,\gamma)\kappa_{+}(\beta,\delta)=q^{2N}$ (2.65)

the invariant measure for the open ASEP can be expressed in terms of a matrix
product measure (MPM) with $(N+1)$-dimensional representation matrices of a
quadratic algebra satisfied by the matrices of the matrix product ansatz
developed in [25]. On the submanifold
$\mathcal{B}^{M}_{N}\subset\mathcal{B}_{N}$ defined by the further constraint

$\kappa_{-}(\alpha,\gamma)\kappa_{-}(\beta,\delta)=q^{-2M},\quad 1\leq M\leq
N$ (2.66)

the MPM has the special property that it exists only for system sizes
$L>N-M+1$ [54] whereas away from this submanifold it exists for all system
sizes [29]. Boundary parameters satisfying (2.65) and (2.66) for $1\leq N\leq
L$ define the setting that the main results proved below apply to.

To get further insight into the significance of the parameter manifold
$\mathcal{B}_{N}$ we note that for $N=0$ the condition (2.65) reduces to
$z_{+}=z_{-}$. As shown in [61], with this constraint on the boundary
parameters the invariant measure is a Bernoulli product measure with density
$\rho_{+}=\rho_{-}$ and the manifolds (B.3) and (B.4) in Appendix B turn into
the stationarity condition (2.56) for the expectation of the current.
Following [49] this observation may be understood as follows. In a semi-
infinite system with $L_{+}=+\infty$ the left boundary rates (2.58) allow for
an invariant Bernoulli product measure with stationary density $\rho_{-}\leq
1/2$ and stationary current given by (B.3) (corresponding to
$j=\mbox{$\langle\,{j^{-}}\,\rangle$}$) while in a semi-infinite system with
$L_{-}=-\infty$ the right boundary rates (2.59) allow for an invariant
Bernoulli product measure with a stationary density $\rho_{+}\geq 1/2$ and
stationary current given by (B.4) (corresponding to
$j=\mbox{$\langle\,{j^{+}}\,\rangle$}$). When $\rho_{+}=\rho_{-}$ or
$\rho_{+}=1-\rho_{-}$ the two currents $\langle\,{j^{-}}\,\rangle$ and
$\langle\,{j^{+}}\,\rangle$ match. Then stationarity of the Bernoulli product
measure holds even in the finite lattice $\Lambda_{L}$ for any $L_{\pm}$ and
any density $\rho:=\rho_{+}=\rho_{-}$.

On the other hand, when $\rho_{+}=1-\rho_{-}$ stationary can only be attained
with a non-homogeneous density profile which has been argued to be a shock
profile [63, 67] analogous to a domain wall in driven non-equilibrium systems
[30, 48, 20, 41, 17] but with diffusive fluctuations of the domain wall
position [28, 62, 23] proved rigorously for the ASEP on the infinite integer
lattice in [32]. Any mismatch of the boundary currents in the range
$\rho_{-}<1/2$ and $\rho_{+}>1/2$ is then expected to be still realized by a
shock profile, but with a shock localized close to the left boundary for
$\mbox{$\langle\,{j^{+}}\,\rangle$}<\mbox{$\langle\,{j^{-}}\,\rangle$}$ or
close to the right boundary for
$\mbox{$\langle\,{j^{+}}\,\rangle$}>\mbox{$\langle\,{j^{-}}\,\rangle$}$, thus
explaining the origin of the phase diagram of the open ASEP in terms of domain
wall dynamics [42]. From a mathematical perspective, this reasoning is a
conjecture based on microscopic shock stability proved only in the
hydrodynamic limit [2]. As will be seen below, the reverse duality constructed
below provides a rigorous proof of this conjecture on microscopic scale for
the parameter manifold $\mathcal{B}^{1}_{N}$ defined by (2.65) and (2.66).

Finally we note that for system size $L=N$ the invariant measure on the
parameter manifold $\mathcal{B}^{1}_{N}$ is a zero-current blocking measure
restricted to the finite lattice [12]. This blocking measure, originally
defined for the infinite system [50], is a product measure with strictly
increasing marginal fugacities $z_{k}\propto q^{2(k-L_{-})}$,
$k\in\Lambda_{L}$, is also an invariant measure of the ASEP with reflecting
boundaries [60] where $\alpha=\beta=\gamma=\delta=0$.

###  2.5 Shocks in the ASEP

####  2.5.1 Bernoulli shock measures

In the spirit of [6, 10] we extend the notion of Bernoulli shock measures to
the finite lattice. To this end and to allow for more compact notation
auxiliary boundary sites are introduced.

###### Definition 2.9

(Bernoulli shock measures) For auxiliary boundary sites $x_{0}:=L--1$ and
$x_{N+1}:=L_{+}+1$ and a nonempty ordered set
$\mathbf{x}:=\\{x_{1},\dots,x_{N})\\}$ of $N$ lattice sites the product
measure

$\mu^{\mathbf{x}}_{\boldsymbol{\eta}}=\prod_{k=L_{-}}^{L_{+}}p^{\mathbf{x}}_{\eta_{k}}$
(2.67)

with marginals

$p^{\mathbf{x}}_{\eta_{k}}=\left\\{\begin{array}[]{ll}(1-\rho^{\star}_{i})(1-\eta_{k})+\rho^{\star}_{i}\eta_{k}&k=x_{i},\quad
1\leq i\leq N\\\
(1-\rho_{i})(1-\eta_{k})+\rho_{i}\eta_{k}&x_{i}<k<x_{i+1},\quad 0\leq i\leq
N\\\ \end{array}\right.$ (2.68)

is called a Bernoulli shock measure with $N$ microscopic shocks at positions
$x_{i}\in\Lambda_{L}$ and left boundary density $\rho_{0}$, right boundary
density $\rho_{N}$, bulk densities $\rho_{i}$ for $1\leq i\leq N-1$, and shock
densities $\rho^{\star}_{i}$ for $1\leq i\leq N$.

This definition is motivated by the fact that seen from a shock position the
local particle configuration of an infinite system becomes uncorrelated at
large distances [31, 24]. The Bernoulli shock measures extend this property to
all lattice sites which was found in [26] for specific shock densities whose
fugacities satisfy the relation

$\frac{z_{i}}{z_{i-1}}=q^{2}$ (2.69)

and second-class particles at the shock positions.

The vector representation of a shock measure with $N$ shocks is given by

$\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}:=\sum_{\boldsymbol{\eta}}\mu^{\mathbf{x}}_{\boldsymbol{\eta}}\mbox{$|\,{\boldsymbol{\eta}}\,\rangle$}$
(2.70)

and with the column vector
$\mbox{$|\,{\rho_{\cdot}}\,\rangle$}:=(1-\rho_{\cdot},\rho_{\cdot})^{T}$ one
gets the Kronecker product

$\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}=\mbox{$|\,{\rho_{0}}\,\rangle$}^{\otimes(x_{1}-L_{-})}\otimes\mbox{$|\,{\rho^{\star}_{1}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}^{\otimes(x_{2}-x_{1}-1)}\otimes\dots\otimes\mbox{$|\,{\rho^{\star}_{N}}\,\rangle$}\otimes\mbox{$|\,{\rho_{N}}\,\rangle$}^{\otimes(L_{+}-x_{N})}$
(2.71)

which reflects the factorized form of the shock measure.

####  2.5.2 Microscopic and macoscopic shocks

For a Bernoulli shock measure the expectation of the particle current (2.53)
inside a domain with density $\rho_{i}$ is given by

$j_{i}\equiv j(\rho_{i})=(r-\ell)\rho_{i}(1-\rho_{i}).$ (2.72)

As shown in [6] for shock densities $\rho^{\star}_{i}=1$ the quantities

$\displaystyle d^{\ell}_{i}$ $\displaystyle=$
$\displaystyle\frac{j_{i-1}}{\rho_{i}-\rho_{i-1}}$ (2.73) $\displaystyle
d^{r}_{i}$ $\displaystyle=$ $\displaystyle\frac{j_{i}}{\rho_{i}-\rho_{i-1}}$
(2.74)

are the shock hopping rates of a microscopically stable shock with bulk
densities $\rho_{i-1},\rho_{i}$ satisfying (2.69) which we therefore call the
microscopic shock stability condition. In [43, 4, 9, 10] the same microscopic
shock jump rates were found with a second-class particle [31] at the shock
position which is equivalent to $\rho^{\star}_{i}=0$. The average shock
velocity of shock $i$ between bulk densities $\rho_{i}$ and $\rho_{i+1}$ is
the difference

$v_{i}=d^{r}_{i}-d^{\ell}_{i}$ (2.75)

of shock jump rates and the variance of the shock position is described by the
diffusion coefficient

$D_{i}=\frac{1}{2}(d^{r}_{i}+d^{\ell}_{i})$ (2.76)

which appear also on macroscopic scale in the random motion of a
macroscopically stable shock in the ASEP [31, 32]. Notice, however, that
because of the convexity of the current-density relation
$j\propto\rho(1-\rho)$ only the milder condition Rankine Hugoniot condition
[59] is required for the macroscopic stability of a shock which for the ASEP
only implies $\rho_{i+1}>\rho_{i}$ rather than the strict condition (2.69).
For consecutive shocks this convexity also yields $v_{i}>v_{i+1}$ for all $i$
so that after finite microscopic time consecutive shocks become very close and
form a bound state for all subsequent times [6]. On macroscopic scale this
phenomenon corresponds to a coalescence of shocks [33]. Also other particle-
conserving models with random-walking shocks are known [43, 3, 39, 7, 5]

####  2.5.3 Shock random walks and MPMs for the open ASEP

A shock measure with $N$ stable shocks which all satisfy the microscopic
stability criterion (2.69) then has the property

$\frac{z_{N}}{z_{0}}=q^{2N}.$ (2.77)

If the boundary parameters of the open ASEP satisfy (2.65) with $z_{-}=z_{0}$
then $z_{+}=z_{N}$ then also (2.77) is satisfied and the $(N+1)$-dimensional
representations of [54] of the stationary matrix product algebra of [25] are
convex combinations of homogenous Bernoulli product measures with fugacities
$z_{0}$ and $z_{N}$ and Bernoulli shock measures with $N$ shocks satisfying
(2.69) and (2.65) [38]. Therefore we refer to (2.65) as boundary shock
stability condition.

In fact, it was suggested even earlier in [43] that the shocks whose invariant
distribution in the finite system is described by the MPM remain stable during
the stochastic time evolution and perform a random walk dynamics. Here we
focus on boundary rates satisfying also the condition (2.77) and use reverse
duality to prove this conjecture and make it precise. In particular, we
identify the shock densities $\rho^{\star}_{i}$ which allows us to elucidate
special properties of MPMs satisfying (2.66).

###  2.6 Shock ASEP with particle-dependent hopping rates

Using duality it is proved below that on microscopic scale shocks in the open
ASEP on $\Lambda$, or more precisely the microscopic shock positions perform a
random motion that we introduced in more general form in [10] under the name
shock exclusion process in the context of self-duality for a multi-species
ASEP. Here we go beyond [10] by considering non-conservative open boundaries
instead of the more straightforward case of conservative dynamics on
${\mathbb{Z}}$. However, we focus on only one species of particles.

For a single species of exclusion particles the shock exclusion process, has
particles, labelled by $i$ with $1\leq i\leq N$ that are located on lattice
sites $x_{i}$ with the single-file condition
$x_{0}<x_{1}<x_{2}\dots<x_{N}<x_{N+1}$ with the fixed auxiliary boundary
coordinates $x_{0}:=L_{-}-1$ and $x_{N+1}=L_{+}+1$. Thus the configuration of
all particles is represented by the coordinate vector
$\mathbf{x}:=(x_{1},x_{2},\dots,x_{N})$ and each particle has its individual
jump rate

$w_{i}(\mathbf{x})=w^{\ell}_{i}(\mathbf{x})+w^{r}_{i}(\mathbf{x})$ (2.78)

with left- and right jump rates

$\displaystyle w^{\ell}_{i}(\mathbf{x})$ $\displaystyle=$ $\displaystyle
d^{\ell}_{i}(1-\delta_{x_{i},x_{i-1}+1})(1-\delta_{x_{i},x_{i+1}})$ (2.79)
$\displaystyle w^{r}_{i}(\mathbf{x})$ $\displaystyle=$ $\displaystyle
d^{r}_{i}(1-\delta_{x_{i},x_{i+1}-1})(1-\delta_{x_{i},x_{i-1}})$ (2.80)

given by the shock jump rates (2.73) and (2.74). The terms with the
Kroneckersymbols entail the single-file exclusion principle and the boundary
terms $x_{0}$ and $x_{N+1}$ appearing in the rates $w_{1}(\mathbf{x})$ and
$w_{N}(\mathbf{x})$ express the particle number conservation at the reflecting
boundaries.

The generator of the shock ASEP can be defined more formally by introducing
for a given coordinate vector $\mathbf{x}=(x_{1},\dots,x_{N})$ the locally
shifted coordinate vectors $\mathbf{x}^{i\pm}$ with the coordinates

$x^{i\pm}_{j}=x_{j}\pm\delta_{i,j}.$ (2.81)

###### Definition 2.10

(Shock exclusion process) The $N$-particle shock exclusion process of
particles $i$, $1\leq i\leq N$ that are located on lattice sites
$x_{i}\in\Lambda_{L}$ with the single-file condition
$x_{0}<x_{1}<x_{2}\dots<x_{N}<x_{N+1}$ and with jump rates (2.79), (2.80) is
defined by the generator

$\mathcal{M}f(\mathbf{x})=\sum_{i=1}^{N}\mathcal{M}_{i}f(\mathbf{x})$ (2.82)

with the single-particle hopping generators

$\mathcal{M}_{i}f(\mathbf{x})=w_{i}^{\ell}(\mathbf{x})(f(\mathbf{x}^{i-})-f(\mathbf{x}))+w_{i}^{r}(\mathbf{x})(f(\mathbf{x}^{i+})-f(\mathbf{x}))$
(2.83)

for fixed parameters $\rho_{i}\in[0,1]$, $0\leq i\leq N$.

###### Remark 2.11

For $N=1$ the shock exclusion process reduces to a biased random walk on
$\Lambda_{L}$ with generator

$\mathcal{M}f(x)=w^{\ell}(x)(f(x-1)-f(x))+w^{r}(x)(f(x+1)-f(x))$ (2.84)

and with left- and right jump rates

$\displaystyle w^{\ell}(x)$ $\displaystyle=$ $\displaystyle
d^{\ell}_{1}(1-\delta_{x,L_{-}})(1-\delta_{x,L_{+}+1})$ (2.85) $\displaystyle
w^{r}(x)$ $\displaystyle=$ $\displaystyle
d^{r}_{1}(1-\delta_{x,L_{+}})(1-\delta_{x,L_{-}-1})$ (2.86)

for any fixed pair of densities $\rho_{0},\rho_{1}\in[0,1]$. The intensity
matrix $Q$ is the tridiagonal Toeplitz matrix with matrix elements

$\displaystyle Q_{xy}$ $\displaystyle=$ $\displaystyle
d^{\ell}_{1}(1-\delta_{x,L_{-}})(\delta_{y,x-1}-\delta_{y,x})+d^{r}_{1}(1-\delta_{x,L_{+}})(\delta_{y,x+1}-\delta_{y,x})$
(2.87)

for $x,y\in\Lambda_{L}$.

##  3 Reverse duality for the open ASEP

The main result is the reverse duality of the non-conservative open ASEP and
the conservative shock ASEP with reflecting boundaries. The duality function
is given by the Bernoulli shock measures of Definition 2.9 with $N$ shocks and
a specific choice of shock densities determined by condition 2.66. Since the
main ideas in the proof for general $N$ appear already in the simpler case of
a single shock, we separate the single-shock case (which also allows for more
explicit information on the time evolution of the shock measure) from the
general case. Before stating these results we assert reversibility of the
shock exclusion process.

###### Proposition 3.1

The $N$-particle shock exclusion process 2.10 on $\Lambda_{L}$ with reflecting
boundaries is reversible w.r.t. the unnormalized measure

$\pi(\mathbf{x})=\prod_{i=1}^{N}d_{i}^{2x_{i}}$ (3.1)

where

$d_{i}:=\sqrt{\frac{d^{r}_{i}}{d^{\ell}_{i}}}$ (3.2)

is the hopping asymmetry of particle $i$.

Now we are in a position to present the main results.

###### Theorem 3.2

Let $W$ be the intensity matrix (2.36) of the open ASEP with boundary rates
(2.58) and (2.59) and for parameters $\rho_{0}$ and $\rho_{1}$ let $Q$ be the
intensity matrix of a simple biased random walk with generator (2.84).
Further, let $\mu^{x}$ be the shock measure (2.67) with left bulk density
$\rho_{0}=\rho_{-}$ and shock density

$\rho_{1}^{\star}=\frac{\alpha}{\alpha+\gamma}.$ (3.3)

The intensity matrices $W$ and $Q$ satisfy the reverse-duality relation

$RW=Q^{T}R$ (3.4)

w.r.t. the duality matrix $R$ with matrix elements
$R_{x,\boldsymbol{\eta}}=d_{1}^{2x}\mu^{x}(\boldsymbol{\eta})$ if and only if
the following three conditions are satisfied:
(i) The bulk shock stability condition (2.69) is satisfied for $i=N=1$,
(ii) The boundary rates are on the manifold $\mathcal{B}^{1}_{1}$ .

###### Theorem 3.3

Denote by $\mu^{x}_{t}$ the distribution at time $t$ of the open ASEP,
starting at $t=0$ from a shock measure $\mu^{x}_{0}=\mu^{x}$ and let
Conditions (i) - (iii) of Theorem 3.2 be satisfied. Then, for any
$x\in\Lambda_{L}$

$\mu^{x}(t)=\sum_{y=L_{-}}^{L_{+}}P(y,t|x,0)\,\mu^{y}(0)$ (3.5)

where

$P(y,t|x,0)=\frac{d_{1}^{2}-1}{d_{1}^{2L}-1}d_{1}^{2(y-L_{-})}+\frac{2}{L}\sum_{p=1}^{L-1}d_{1}^{y-x}\psi_{p}(x)\psi_{p}(y)\frac{w}{\epsilon_{p}}\mathrm{e}^{-\epsilon_{p}t}$
(3.6)

with

$\displaystyle\psi_{p}(y)$ $\displaystyle:=$ $\displaystyle
d_{1}\sin{\left(\frac{\pi p}{L}(y+1-L_{-})\right)}-\sin{\left(\frac{\pi
p}{L}(y-L_{-})\right)}$ (3.7) $\displaystyle\epsilon_{p}$ $\displaystyle=$
$\displaystyle w\left[d_{1}+d^{-1}_{1}-2\cos{\left(\frac{\pi
p}{L}\right)}\right]$ (3.8)

is the transition probability of the biased random walk generated by the
intensity matrix $Q$ and starting at time $t=0$ from $x$. The limit
$\mu^{\ast}_{1}:=\lim_{t\to\infty}\mu^{x}(t)$ is the unique invariant measure
of the open ASEP with boundary parameters as specified by Conditions (ii) and
(iii) and is given by the convex combination

$\mu^{\ast}_{1}=\frac{d_{1}^{2}-1}{d_{1}^{2L}-1}\sum_{y=L_{-}}^{L_{+}}d_{1}^{2(y-L_{-})}\mu^{y}.$
(3.9)

of shock measures $\mu^{y}$.

###### Remark 3.4

For asymmetric exclusion processes with open boundaries similar random walk
dynamics were found earlier for a single anti-shock with bulk densities
related by the inverse stability condition $z_{1}=q^{-2}z_{0}$ under
conditioning the process on a special atypical value of the time-integrated
particle current [8]. Boundary parameters in that work corresponded to
$\omega_{+}=\omega_{-}=0$. In [39] random walk dynamics of a single shock were
proved for a discrete-time totally asymmetric simple exclusion process with
deterministic sublattice update and stochastic open boundaries. In this model
there is no constraint on the boundary parameters.

###### Remark 3.5

According to condition (iii) in Theorem 3.2 the invariant measure (3.9) can be
expressed by the two-dimensional representation of the stationary matrix
product algebra characterized in the Appendices of [54] which exists for any
$L\geq 2$ and which is analyzed further in [12]. The homogeneous Bernoulli
measure with densities $\rho_{0}$ and $\rho_{1}$ generically contribute to the
two-dimensional MPM [38] but are absent in the convex combination (3.9). It
would be interesting to explore whether two-dimensional time-dependent
matrices can capture the dynamics of some subspace as has been proved for the
SSEP with open boundaries [68] and also for a non-conservative annihilating
random walks with pair deposition [65, 69].

###### Remark 3.6

The spectrum of the generator given by the eigenvalues (3.8) yields a subset
of eigenvalues of the generator of the open ASEP and is in agreement with the
picture of spectral properties arising from a shock random walk explored in
[28, 62, 23].

As pointed out above, Theorem 3.2 is a special case of the following result
for duality functions given by Bernoulli shock measures with $N>1$ shocks.

###### Theorem 3.7

Let $W$ be the intensity matrix (2.36) of the open ASEP with boundary rates
(2.58) and (2.59) and for parameters $\rho_{0},\dots,\rho_{N}$ let $Q$ be the
intensity matrix of the $N$-particle shock exclusion process of Definition
2.10 and reversible measure (3.1). Further, let $\mu^{\mathbf{x}}$ be the
shock measure (2.67) with left boundary density $\rho_{0}=\rho_{-}$ and shock
fugacities

$z^{\star}_{i}=\frac{\alpha}{\gamma}q^{2(i-1)}$ (3.10)

for $1\leq i\leq N$. The intensity matrices $W$ and $Q$ satisfy the reverse-
duality relation

$RW=Q^{T}R$ (3.11)

w.r.t. the duality matrix $R$ with matrix elements
$R_{x,\boldsymbol{\eta}}=\pi(\mathbf{x})\mu^{\mathbf{x}}_{\boldsymbol{\eta}}$
if and only if the following three conditions are satisfied:
(i) The bulk shock stability condition (2.69) is satisfied for all
$i\in\\{1,\dots,N\\}$,
(ii) The boundary rates are on the manifold $\mathcal{B}^{1}_{N}$ for $1\leq
N\leq L$.

###### Remark 3.8

The conservative reflective boundaries of the reverse dual are in contrast to
the conventional duality for the open symmetric simple exclusion process
(SSEP) which is dual to the SSEP with nonconservative absorbing boundaries
noticed in [71] and later fully developed in [13, 34], including the
generalization to symmetric partial exclusion [64]. An open question is
whether a similar conservative reverse duality can be constructed for the
conservative asymmetric partial exclusion process studied in [14] when non-
conservative open boundaries are added.

Since the shock exclusion process with reflecting boundaries is reversible for
any $N$ also Theorem 3.3 has a natural generalization for $N>1$.

###### Theorem 3.9

Let $\mu^{\mathbf{x}}(t)$ denote the distribution at time $t$ of the ASEP of
Theorem 3.7, starting from an $N$-shock measure
$\mu^{\mathbf{x}}=\mu^{\mathbf{x}}(0)$. Then

$\mu^{\mathbf{x}}(t)=\sum_{\mathbf{y}}P(\mathbf{y},t|\mathbf{x},0)\,\mu^{\mathbf{y}}$
(3.12)

where $P(\mathbf{y},t|\mathbf{x},0)$ is the transition probability of the
shock exclusion process 2.10. The limit
$\mu^{\ast}_{N}:=\lim_{t\to\infty}\mu^{\mathbf{x}}(t)$ is the unique invariant
measure of the open ASEP with boundary parameters as specified by Conditions
(ii) and (iii) and is given by the unnormalized convex combination

$\tilde{\mu}^{\ast}_{N}=\frac{j_{N}^{2y_{N}}}{j_{1}^{2y_{1}}}\sum_{\mathbf{y}}\prod_{i=1}^{N-1}j_{i}^{2(y_{i}-y_{i+1})}\mu^{\mathbf{y}}$
(3.13)

of shock measures with $N$ shocks.

###### Remark 3.10

The spectrum of the generator of the $N$-particle shock exclusion process can
be calculated from Bethe ansatz [55, 56, 22, 23, 70]. However, there is no
simple expression in closed form. For results on the eigenfunctions see [21].
Further progress along these lines has been made recently for the SSEP with
open boundaries by using duality [34]. It would be most interesting to
generalize this approach to the open ASEP to cover the full manifold of
boundary parameters.

##  4 Proofs

We prove reverse duality and its consequences for the time evolution of shock
measures by studying the action of the quantum Hamiltonian on the vector
representation of the shock measures which are shown to form an invariant
subspace. The coefficients appearing in these computations yield the reverse
duality with the shock exclusion process. Sufficient conditions for a similar
property were conjectured earlier in [43] for second-class particles at the
shock positions and in [38] for a single shock. Here we obtain rigorously the
necessary and sufficient conditions on the parameters of the ASEP for the
reverse duality for a specific choice of shock densities for $N\geq 1$
shocks..

###  4.1 Preliminaries

Everywhere below a vector denoted by $|\,{abc\dots}\,\rangle$ is defined to be
the Kronecker product
$\mbox{$|\,{a}\,\rangle$}\otimes\mbox{$|\,{b}\,\rangle$}\otimes\mbox{$|\,{c}\,\rangle$}\otimes\dots$
of two-dimensional vectors
$\mbox{$|\,{a}\,\rangle$},\mbox{$|\,{b}\,\rangle$},\mbox{$|\,{c}\,\rangle$}\dots\in{\mathbb{R}}^{2}$
with components $(a_{0},a_{1}),(b_{0},b_{1}),(c_{0},c_{1}),\dots$. For vectors
with components of the form $a_{0}=1-\rho$, $a_{1}=\rho$ and $0\leq\rho\leq 1$
we use greek letters and write $|\,{\rho}\,\rangle$ and call $\rho$ a density.
In slight abuse of wording we call for any vector $|\,{a}\,\rangle$ with
$a_{0}\neq 0$ the ratio

$z(a):=\frac{a_{1}}{a_{0}}$ (4.1)

the fugacity. The functions

$\epsilon_{-}(\cdot):=\alpha-\gamma
z(\cdot),\quad\epsilon_{+}(\cdot):=\delta-\beta z(\cdot)$ (4.2)

play a role as eigenvalues of the boundary operators $h^{\pm}$.

For two vectors $|\,{a}\,\rangle$, $|\,{b}\,\rangle$ we introduce the
determinant

$\Delta(a,b):=a_{0}b_{1}-a_{1}b_{0}$ (4.3)

and recall that linear dependence of $|\,{a}\,\rangle$ and $|\,{b}\,\rangle$
is equivalent to $\Delta(a,b)=0$. For $z(b)=q^{2}z(a)$ one has

$\Delta(a,b)=(r-\ell)\frac{a_{0}b_{1}}{r}=(r-\ell)\frac{a_{1}b_{0}}{\ell}.$
(4.4)

For linearly independent vectors we also define the functions

$\displaystyle
d(a,b):=(r-\ell)\frac{a_{1}a_{0}}{\Delta(a,b)},\quad\tilde{d}(a,b):=(r-\ell)\frac{a_{1}b_{0}}{\Delta(a,b)}.$
(4.5) $\displaystyle
d_{\pm}(a,b):=\frac{\epsilon_{\pm}(a)a_{0}(a_{0}+a_{1})}{\Delta(a,b)},\quad\tilde{d}_{\pm}(a,b):=\frac{\epsilon_{\pm}(a)a_{0}(b_{0}+b_{1})}{\Delta(a,b)}.$
(4.6)

which play a role for the transition rates of the dual process. We recall that
in all assertions and in all the proofs below it is tacitly assumed that all
boundary rates and bulk jump rates of the open ASEP are strictly positive and
that $q\neq 1$.

######  Lemma 4.1 (Eigenvectors)

(i) For any two vectors $|\,{a}\,\rangle$, $|\,{b}\,\rangle$ the bulk
eigenvalue property

$h\mbox{$|\,{ab}\,\rangle$}=0$ (4.7)

is satisfied if and only if $\Delta(a,b)=0$.
(ii) For vectors $|\,{a^{\pm}}\,\rangle$ with components such that
$a^{\pm}_{0}a^{\pm}_{1}\neq 0$ and $a^{\pm}_{0}+a^{\pm}_{1}\neq 0$ the two
boundary eigenvalue equations

$h^{\pm}\mbox{$|\,{a^{\pm}}\,\rangle$}=\epsilon_{\pm}(a^{\pm})\mbox{$|\,{a^{\pm}}\,\rangle$}$
(4.8)

have a unique solution with strictly positive fugacities given by

$z(a^{-})=\kappa^{-1}(\alpha,\gamma),\quad z(a^{+})=\kappa_{+}(\beta,\delta)$
(4.9)

with $\kappa_{+}(\cdot,\cdot)$ defined in (2.62). The eigenvalues are

$\displaystyle\epsilon_{-}(\cdot)$ $\displaystyle=$
$\displaystyle\alpha-\gamma z(\cdot)\,=\,(r-\ell)\frac{z(\cdot)}{1+z(\cdot)}$
(4.10) $\displaystyle\epsilon_{+}(\cdot)$ $\displaystyle=$
$\displaystyle\delta-\beta z(\cdot)\,=\,-(r-\ell)\frac{z(\cdot)}{1+z(\cdot)}.$
(4.11)

For $a^{\pm}_{0}a^{\pm}_{1}=0$ or $a^{\pm}_{0}+a^{\pm}_{1}=0$ the eigenvalue
equation (4.8) has no solution for strictly positive boundary rates and
$r\neq\ell$.

Proof: Consider two vectors $|\,{a}\,\rangle$ with components $a_{0},a_{1}$
and $|\,{b}\,\rangle$ with components $b_{0},b_{1}$ and determinant
$\Delta_{ab}$ defined in (4.3).
(i) For the product vector one gets from the definition (2.45)

$\displaystyle h\mbox{$|\,{ab}\,\rangle$}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}0&0&0&0\\\ 0&r&-r&0\\\
0&-\ell&\ell&0\\\
0&0&0&0\end{array}\right)\left(\begin{array}[]{c}a_{0}b_{0}\\\ a_{0}b_{1}\\\
a_{1}b_{0}\\\
a_{1}b_{1}\end{array}\right)=\Delta_{ab}\left(\begin{array}[]{c}0\\\ r\\\
-\ell\\\ 0\end{array}\right)$ (4.24)

which proves part (i) of the Lemma.
(ii) To make notation less heavy we drop the superscripts $\pm$ in the vectors
$|\,{a^{\pm}}\,\rangle$ in the following computations. From the definition
(2.48) of the boundary matrix one obtains

$h^{-}\mbox{$|\,{a}\,\rangle$}=\left(\begin{array}[]{c}\alpha a_{0}-\gamma
a_{1}\\\ -\alpha a_{0}+\gamma a_{1}+(r-\ell)a_{1}\end{array}\right)$ (4.25)

from which one reads off immediately that for strictly positive rates
$\alpha,\gamma$ there exists no eigenvector with $a_{0}=0$ or $a_{1}=0$, i.e.,
for $a_{0}a_{1}=0$. To proceed with the proof of the eigenvalue properties for
a vector $|\,{a}\,\rangle$ such that $a_{0}a_{1}\neq 0$ we define

$g_{\pm}(a):=\mp\frac{r-\ell}{\epsilon_{-}(a)}-\frac{1}{z(a)}$ (4.26)

and note that for any linearly independent vector $|\,{b}\,\rangle$ the
equality

$\displaystyle h^{-}\mbox{$|\,{a}\,\rangle$}$ $\displaystyle=$
$\displaystyle\epsilon_{-}(a)\left(\begin{array}[]{c}a_{0}\\\
g_{-}(a)a_{1}\end{array}\right)$ (4.29) $\displaystyle=$
$\displaystyle\left(\tilde{d}_{-}(a,b)-\tilde{d}(a,b)\right)\mbox{$|\,{a}\,\rangle$}-\left(d_{-}(a,b)-d(a,b)\right)\mbox{$|\,{b}\,\rangle$}$
(4.30)

holds without further constraints on the vectors $|\,{a}\,\rangle$ and
$|\,{b}\,\rangle$.

If an eigenvector exists (which requires $d_{-}(a,b)=d(a,b)$ or equivalently
$g_{-}(a)=1$) then the eigenvalue is of the form $\epsilon_{-}(a)$ as defined
in (4.2) and one deduces that there exists no eigenvector with $a_{1}=-a_{0}$
since then $z=-1$ and therefore $g(a)\neq 1$ for $r\neq\ell$. From
$g_{-}(a)=1$ one also obtains the second equality in (4.10) and together with
the first equality one arrives at the quadratic equation

$z^{2}+\frac{\gamma-\alpha+r-\ell}{\gamma}z-\frac{\alpha}{\gamma}=0$ (4.31)

for the fugacity $z\equiv z(a)$. The positive root is given by
$z=\kappa^{-1}(\alpha,\gamma)$ which yields (up to an arbitrary multiplicative
factor) the components of the eigenvector $|\,{a}\,\rangle$ in terms of the
boundary rates as in (4.9).

For the right boundary matrix one gets

$\displaystyle h^{+}\mbox{$|\,{a}\,\rangle$}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}\delta a_{0}-\beta a_{1}\\\ -\delta
a_{0}+\beta a_{1}-(r-\ell)a_{1}\end{array}\right)$ (4.34) $\displaystyle=$
$\displaystyle\epsilon_{+}(a)\left(\begin{array}[]{c}a_{0}\\\
g_{+}(a)a_{1}\end{array}\right)$ (4.37) $\displaystyle=$
$\displaystyle\left(\tilde{d}_{+}(a,b)+\tilde{d}(a,b)\right)\mbox{$|\,{a}\,\rangle$}-\left(d_{+}(a,b)+d(a,b)\right)\mbox{$|\,{b}\,\rangle$}.$
(4.38)

with a linearly independent vector $|\,{b}\,\rangle$. Eigenvectors with
$a_{0}a_{1}=0$ or $a_{0}+a_{1}=0$ cannot exist as can be seen by arguments
analogous to the case treated above. For an eigenvector one immediately reads
off the first equality in (4.11). Moreover, an eigenvector satisfies
$d_{+}(a,b)=-d(a,b)$ which gives the second equality in (4.11) and together
with the first equality in (4.11) the quadratic equation

$z^{2}-\frac{\delta-\beta+r-\ell}{\beta}z-\frac{\delta}{\beta}=0$ (4.39)

for the fugacity. The positive root is given by $z=\kappa_{+}(\beta,\delta)$
which yields (4.9) for the right boundary. $\Box$

###### Corollary 4.2

A probability vector $|\,{\pi}\,\rangle$ with factorized boundary marginals
such that
$\mbox{$|\,{\pi}\,\rangle$}=\mbox{$|\,{\rho_{-}}\,\rangle$}\otimes\mbox{$|\,{\tilde{\pi}}\,\rangle$}\otimes\mbox{$|\,{\rho_{+}}\,\rangle$}$
and with $z_{+}=q^{2N}z_{-}$ is an eigenvector of the boundary operator
$h^{-}_{L_{-}}+h^{+}_{L_{+}}$ if and only if (2.65) holds.

######  Lemma 4.3 (Projection)

(i) For vectors $|\,{a^{(i)}}\,\rangle$, $|\,{\tilde{a}^{(i)}}\,\rangle$,
$i\in\\{1,2\\}$, such that $a^{(i)}_{0}a^{(i)}_{1}\neq 0$ and
$\tilde{a}^{(i)}_{0}\tilde{a}^{(i)}_{1}\neq 0$ and a pair of linearly
independent vectors ($|\,{b^{(i)}}\,\rangle$, $|\,{\tilde{b}^{(i)}}\,\rangle$)
such that $|\,{b^{(1)}}\,\rangle$ is linearly independent of
$|\,{\tilde{a}^{(1)}}\,\rangle$ and $|\,{b^{(2)}}\,\rangle$ is linearly
independent of $|\,{\tilde{a}^{(2)}}\,\rangle$ the two-site jump generator $h$
has the three-dimensional projection properties

$\displaystyle(A)\qquad h\mbox{$|\,{a^{(1)}b^{(1)}}\,\rangle$}$
$\displaystyle=$ $\displaystyle
u^{(1)}\mbox{$|\,{a^{(1)}b^{(1)}}\,\rangle$}-v^{(1)}\mbox{$|\,{b^{(1)}\tilde{a}^{(1)}}\,\rangle$}-w^{(1)}\mbox{$|\,{a^{(1)}\tilde{b}^{(1)}}\,\rangle$}$
(4.40)

with coefficients

$u^{(1)}=\ell-\tilde{d}(b^{(1)},\tilde{b}^{(1)}),\quad
v^{(1)}=d(a^{(1)},\tilde{a}^{(1)}),\quad w^{(1)}=-d(b^{(1)},\tilde{b}^{(1)})$
(4.41)

and

$\displaystyle(B)\qquad h\mbox{$|\,{b^{(2)}\tilde{a}^{(2)}}\,\rangle$}$
$\displaystyle=$ $\displaystyle
u^{(2)}\mbox{$|\,{b^{(2)}\tilde{a}^{(2)}}\,\rangle$}-v^{(2)}\mbox{$|\,{a^{(2)}b^{(2)}}\,\rangle$}-w^{(2)}\mbox{$|\,{\tilde{b}^{(2)}\tilde{a}^{(2)}}\,\rangle$}$
(4.42)

with coefficients

$u^{(2)}=r+\tilde{d}(b^{(2)},\tilde{b}^{(2)}),\quad
v^{(2)}=-d(\tilde{a}^{(2)},a^{(2)}),\quad w^{(2)}=d(b^{(2)},\tilde{b}^{(2)})$
(4.43)

if and only if the condition

$\frac{z(\tilde{a}^{(i)})}{z(a^{(i)})}=q^{2},\quad$ (4.44)

on the vectors $|\,{a^{(i)}}\,\rangle$, $|\,{\tilde{a}^{(i)}}\,\rangle$ are
satisfied.
(ii) For a pair of linearly independent vectors ($|\,{a}\,\rangle$,
$|\,{\tilde{a}}\,\rangle$) such that $a_{0}a_{1}\neq 0$ and
$\tilde{a}_{0}\tilde{a}_{1}\neq 0$ and vectors $|\,{\tilde{c}}\,\rangle$
linearly independent of $|\,{a}\,\rangle$, $|\,{c}\,\rangle$ linearly
independent of $|\,{\tilde{a}}\,\rangle$ the two-site jump generator $h$ has
the three-dimensional projection property

$\displaystyle(C)\qquad h\mbox{$|\,{a\tilde{a}}\,\rangle$}$ $\displaystyle=$
$\displaystyle
u\mbox{$|\,{a\tilde{a}}\,\rangle$}-w\mbox{$|\,{\tilde{c}\tilde{a}}\,\rangle$}-v\mbox{$|\,{ac}\,\rangle$}$
(4.45)

with coefficients

$u=\tilde{d}(a,\tilde{c})-\tilde{d}(\tilde{a},c),\quad v=-d(\tilde{a},c),\quad
w=d(a,\tilde{c})$ (4.46)

if and only if the condition (4.44) is satisfied.
(iii) If $a_{0}a_{1}=0$ or $\tilde{a}_{0}\tilde{a}_{1}=0$ the projection
equations (4.40), (4.42), and (4.45) have no solution for strictly positive
boundary rates and $r\neq\ell$.

Proof: We focus on the projection property (A) and drop the superscript $(1)$
to make notation less heavy. The proof of the projection properties (B) and
(C) is analogous. In simplied notation (4.40) reads
$h\mbox{$|\,{ab}\,\rangle$}=u\mbox{$|\,{ab}\,\rangle$}-v\mbox{$|\,{b\tilde{a}}\,\rangle$}-w\mbox{$|\,{a\tilde{b}}\,\rangle$}$
from which using (4.24) one obtains are the four linear equations

$\displaystyle 0$ $\displaystyle=$ $\displaystyle
ua_{0}b_{0}-vb_{0}\tilde{a}_{0}-wa_{0}\tilde{b}_{0}$ (4.47) $\displaystyle
r(a_{0}b_{1}-a_{1}b_{0})$ $\displaystyle=$ $\displaystyle
ua_{0}b_{1}-vb_{0}\tilde{a}_{1}-wa_{0}\tilde{b}_{1}$ (4.48)
$\displaystyle-\ell(a_{0}b_{1}-a_{1}b_{0})$ $\displaystyle=$ $\displaystyle
ua_{1}b_{0}-vb_{1}\tilde{a}_{0}-wa_{1}\tilde{b}_{0}$ (4.49) $\displaystyle 0$
$\displaystyle=$ $\displaystyle
ua_{1}b_{1}-vb_{1}\tilde{a}_{1}-wa_{1}\tilde{b}_{1}$ (4.50)

for the three coefficients $u,v,w$. (a) Consider $a_{0}=0$ which implies by
linear independence $b_{0}\neq 0$. Then the first two equations yield $v=r=0$
in contradiction to $r>0$. Hence there exists no solution for $a_{0}=0$.
(b) Consider $a_{0}\neq 0$ and $\tilde{a}_{0}=0$. The first and third equation
reduce to

$\displaystyle 0$ $\displaystyle=$ $\displaystyle ub_{0}-w\tilde{b}_{0}$
(4.51) $\displaystyle-\ell(a_{0}b_{1}-a_{1}b_{0})$ $\displaystyle=$
$\displaystyle a_{1}\left(ub_{0}-w\tilde{b}_{0}\right)$ (4.52)

which due to linear independence yields $\ell=0$, in contradiction to the
assumption $\ell>0$. Hence there exists no solution for $\tilde{a}_{0}=0$. The
proof that there exists no solution for $a_{1}=0$ or $\tilde{a}_{1}=0$ is
analogous.
(c) To proceed with the generic case we introduce the fugacities

$\tilde{z}:=\frac{\tilde{a}_{1}}{\tilde{a}_{0}},\quad z:=\frac{a_{1}}{a_{0}}$
(4.53)

and consider next $b_{0}=0$. Using linear independence the four equations
yield immediately the unique solution

$\displaystyle u=r,\quad v=\ell\frac{a_{0}}{\tilde{a}_{0}},\quad w=0$ (4.54)

if and only if $\tilde{z}=q^{2}z$ in agreement with the assertion of the
theorem. For $b_{0}\neq 0$ it is useful to define

$\tilde{v}:=v\frac{\tilde{a}_{0}}{a_{0}},\quad\tilde{w}:=w\frac{\tilde{b}_{0}}{b_{0}}.$
(4.55)

The four equations turn into

$\displaystyle u$ $\displaystyle=$ $\displaystyle\ell+\tilde{w}$ (4.56)
$\displaystyle\tilde{v}$ $\displaystyle=$ $\displaystyle\ell$ (4.57)

and

$\displaystyle r(a_{0}b_{1}-a_{1}b_{0})$ $\displaystyle=$ $\displaystyle
a_{0}\left[(\ell+\tilde{w})b_{1}-\tilde{w}\frac{b_{0}}{\tilde{b}_{0}}\tilde{b}_{1}\right]-\ell\frac{a_{0}}{\tilde{a}_{0}}b_{0}\tilde{a}_{1}$
(4.58) $\displaystyle 0$ $\displaystyle=$ $\displaystyle
a_{1}\left[(\ell+\tilde{w})b_{1}-\tilde{w}\frac{b_{0}}{\tilde{b}_{0}}\tilde{b}_{1}\right]-\ell\frac{a_{0}}{\tilde{a}_{0}}b_{1}\tilde{a}_{1}.$
(4.59)

Inserting the second into the first equation then yields

$\displaystyle\tilde{w}$ $\displaystyle=$
$\displaystyle\frac{(r-\ell)b_{1}\tilde{b}_{0}}{b_{1}\tilde{b}_{0}-b_{0}\tilde{b}_{1}}\,=\,-\frac{(r-\ell)b_{1}\tilde{b}_{0}}{\Delta_{b\tilde{b}}}$
(4.60)

and this solution exists if and only if $\tilde{z}=q^{2}z$ as the second
equation shows. Thus with (4.3) one arrives at

$\displaystyle u$ $\displaystyle=$
$\displaystyle\ell-(r-\ell)\frac{b_{1}\tilde{b}_{0}}{\Delta_{b\tilde{b}}}$
(4.61) $\displaystyle v$ $\displaystyle=$
$\displaystyle\ell\frac{a_{0}}{\tilde{a}_{0}}\,=\,(r-\ell)\frac{a_{0}a_{1}}{\Delta_{a\tilde{a}}}$
(4.62) $\displaystyle w$ $\displaystyle=$
$\displaystyle-(r-\ell)\frac{b_{0}b_{1}}{\Delta_{b\tilde{b}}}$ (4.63)

and with the functions $d(\cdot,\cdot)$ and $\tilde{d}(\cdot,\cdot)$ defined
in (4.5) one obtains the coefficients (4.41). $\Box$

###### Corollary 4.4

For shocks satisfying the shock stability criterion (2.69) the shock jump
rates $d_{ir,\ell}$ satisfy the identities

$\displaystyle d(\rho_{i-1},\rho_{i})$ $\displaystyle=$
$\displaystyle\ell\frac{1-\rho_{i-1}}{1-\rho_{i}}\,=\,r\frac{\rho_{i-1}}{\rho_{i}}\,=\,\ell(1-\rho_{i-1})+r\rho_{i-1}\,=\,d^{\ell}_{i}$
(4.64) $\displaystyle-d(\rho_{i},\rho_{i-1})$ $\displaystyle=$
$\displaystyle\ell\frac{\rho_{i}}{\rho_{i-1}}\,=\,r\frac{1-\rho_{i}}{1-\rho_{i-1}}\,=\,r(1-\rho_{i})+\ell\rho_{i}\,=\,d^{r}_{i}$
(4.67) $\displaystyle\epsilon_{+}(\rho_{i})\,=\,d^{r}_{i}-r$
$\displaystyle\epsilon_{-}(\rho_{i-1})\,=\,d^{\ell}_{i}-\ell.$

which also imply $d^{r}_{i}d^{\ell}_{i}=r\ell$ independently of $i$.

###  4.2 Proof of Proposition 3.1

The proof of reversibility of the shock exclusion process) is essentially
trivial. According to Definition 2.10 the intensity matrix of the $N$-particle
shock exclusion process is given by the matrix elements

$Q_{\mathbf{x}\mathbf{x}^{\prime}}=\sum_{i=1}^{N}\left[w^{\ell}_{i}(\mathbf{x})(\delta_{\mathbf{x}^{\prime},\mathbf{x}^{i-}}-\delta_{\mathbf{x}^{\prime},\mathbf{x}})+w^{r}_{i}(\mathbf{x})(\delta_{\mathbf{x}^{\prime},\mathbf{x}^{i+}}-\delta_{\mathbf{x}^{\prime},\mathbf{x}})\right].$
(4.68)

Using
$\delta_{\mathbf{y}^{i\pm},\mathbf{x}}=\delta_{\mathbf{y},\mathbf{x}^{i\mp}}$
one obtains

$\displaystyle
d_{i}^{2}w^{\ell}_{i}(\mathbf{x}^{\prime})\delta_{\mathbf{x},(\mathbf{x}^{\prime})^{i-}}$
$\displaystyle=$ $\displaystyle
w^{r}_{i}(\mathbf{x})\delta_{\mathbf{x}^{\prime},\mathbf{x}^{i+}}$ (4.69)
$\displaystyle
d_{i}^{-2}w^{r}_{i}(\mathbf{x}^{\prime})\delta_{\mathbf{x},(\mathbf{x}^{\prime})^{i+}}$
$\displaystyle=$ $\displaystyle
w^{\ell}_{i}(\mathbf{x})\delta_{\mathbf{x}^{\prime},\mathbf{x}^{i-}}$ (4.70)

and $\pi(\mathbf{x}^{i\pm})/\pi(\mathbf{x})=d_{i}^{\pm 2}$ yields
$Q_{\mathbf{y}\mathbf{x}}=Q_{\mathbf{x}\mathbf{y}}\pi(\mathbf{x})/\pi(\mathbf{y})$
which is the definition of reversibility, viz.,
$Q^{T}=\hat{\pi}Q\hat{\pi}^{-1}$, in terms of the matrix elements. $\Box$

###  4.3 Reverse duality involving one shock

The modified duality matrix $S$ with elements
$S_{x,\boldsymbol{\eta}}=\mbox{$\langle\,{\boldsymbol{\eta}}\,|\,{\mu^{x}}\,\rangle$}=\mu^{x}(\boldsymbol{\eta})$
is the matrix with the column vectors $|\,{\mu^{x}}\,\rangle$ as rows. For one
shock the vector representation (2.71) of the shock measure of Theorem 2.4 is
given by

$\mbox{$|\,{\mu^{x}}\,\rangle$}=\mbox{$|\,{\rho_{0}}\,\rangle$}^{\otimes(x-L_{-})}\otimes\mbox{$|\,{\rho^{\star}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}^{\otimes(L_{+}-x)},\quad
L_{-}\leq x\leq L_{+}$ (4.71)

with $\rho^{\star}$ given by (3.3).

####  4.3.1 Proof of Theorem 3.2

One has to find the necessary and sufficient conditions under which the shock
measures $|\,{\mu^{x}}\,\rangle$ form an invariant subspace under the action
of $H$, i.e.,

$H\mbox{$|\,{\mu^{x}}\,\rangle$}=\sum_{y=L_{-}}^{L_{+}}G_{yx}\mbox{$|\,{\mu^{y}}\,\rangle$}$
(4.72)

with a matrix $G$ and then show that $G=-Q^{T}$ for the intensity matrix $Q$
of the biased random walk defined in Remark 2.11. To do so, we first employ
the projection lemma 4.3 together with the eigenvector properties detailed in
Lemma 4.1 and its Corollary 4.2 to prove (4.72) if and only if (2.65), (2.66)
and (2.77) hold with $N=m+1=1$.

##### Shock in the range $L_{-}<x<L_{+}$:

From the product form of the shock measure, the eigenvalue properties (4.7)
and (4.8) in Lemma 4.1 and Corollary 4.2 one concludes

$\displaystyle H\mbox{$|\,{\mu^{x}}\,\rangle$}$ $\displaystyle=$
$\displaystyle(\tilde{h}^{-}_{L_{-}}+\tilde{h}_{x-1,x}+\tilde{h}_{x,x+1}+\tilde{h}^{+}_{L_{+}})\mbox{$|\,{\mu^{x}}\,\rangle$}$
(4.73) $\displaystyle=$
$\displaystyle\mbox{$|\,{\rho_{0}}\,\rangle$}^{\otimes(x-L_{-}-1)}\otimes[(\tilde{h}_{12}+\tilde{h}_{23})\mbox{$|\,{\rho_{0}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}^{\star}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}]\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}^{\otimes(L_{+}-x-1)}$
$\displaystyle+(\epsilon_{-}(\rho_{0})+\epsilon_{+}(\rho_{1}))\mbox{$|\,{\mu^{x}}\,\rangle$}.$

Taking $a^{(1)}=a^{(2)}=\rho_{0}$, $\tilde{a}^{(1)}=\tilde{a}^{(2)}=\rho_{1}$,
$b^{(1)}=b^{(2)}=\rho_{1}^{\star}$ in (4.40) and (4.42) of Lemma 4.3 yields
after taking Kronecker products on the right in (4.40) and on the left in
(4.42)

$\displaystyle\tilde{h}_{12}\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(\ell-\tilde{d}(\rho_{1}^{\star},\tilde{b}^{(1)})\right)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle-d(\rho_{0},\rho_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}\rho_{1}}\,\rangle$}+d(\rho_{1}^{\star},\tilde{b}^{(1)})\mbox{$|\,{\rho_{0}\tilde{b}^{(1)}\rho_{1}}\,\rangle$}$
$\displaystyle\tilde{h}_{23}\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(r+\tilde{d}(\rho_{1}^{\star},\tilde{b}^{(2)})\right)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
(4.74)
$\displaystyle+d(\rho_{1},\rho_{0})\mbox{$|\,{\rho_{0}\rho_{0}\rho_{1}^{\star}}\,\rangle$}-d(\rho_{1}^{\star},\tilde{b}^{(2)})\mbox{$|\,{\rho_{0}\tilde{b}^{(2)}\rho_{1}}\,\rangle$}.$

Since $\tilde{b}^{(1)}$ and $\tilde{b}^{(2)}$ are arbitrary (except for the
immaterial inequality with $\rho_{1}^{\star}$) we can without loss of
generality take $\tilde{b}^{(1)}=\tilde{b}^{(2)}\neq\rho^{\star}$. Thus the
projection lemma asserts that for $q>1$, $\rho_{1}^{\star}\in[0,1]$ and any
pair of densities $\rho_{0},\rho_{1}\in(0,1)$ such that
$0<\rho_{0}<\rho_{1}<1$ and $z_{1}=q^{2}z_{0}$ the three-site jump generator
$\tilde{h}_{12}+\tilde{h}_{23}$ projects the vector
$\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}:=\mbox{$|\,{\rho_{0}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}^{\star}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}$
on the three-dimensional subspace of spanned by
$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$,
$\mbox{$|\,{\rho_{1}^{\star}\rho_{1}\rho_{1}}\,\rangle$}:=\mbox{$|\,{\rho_{1}^{\star}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}^{\otimes
2}$, and
$\mbox{$|\,{\rho_{0}\rho_{0}\rho_{1}^{\star}}\,\rangle$}:=\mbox{$|\,{\rho_{0}}\,\rangle$}^{\otimes
2}\otimes\mbox{$|\,{\rho_{1}^{\star}}\,\rangle$}$ as

$(\tilde{h}_{12}+\tilde{h}_{23})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}=(r+\ell)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}\rho_{1}}\,\rangle$}-v_{01}\mbox{$|\,{\rho_{1}^{\star}\rho_{1}\rho_{1}}\,\rangle$}-w_{01}\mbox{$|\,{\rho_{0}\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
(4.75)

with constants

$u_{01}=r+\ell,\quad w_{01}=-d(\rho_{1},\rho_{0}),\quad
v_{01}=d(\rho_{0},\rho_{1})$ (4.76)

that do not depend on $\rho_{1}^{\star}$. The identities (4.64) - (4.67) then
yield for arbitrary shock density $\rho_{1}^{\star}$ the bulk closure relation

$\displaystyle H\mbox{$|\,{\mu^{x}}\,\rangle$}$ $\displaystyle=$
$\displaystyle[r+\ell+\epsilon_{-}(\rho_{0})+\epsilon_{+}(\rho_{1})]\mbox{$|\,{\mu^{x}}\,\rangle$}-d^{\ell}_{1}\mbox{$|\,{\mu^{x-1}}\,\rangle$}-d^{r}_{1}\mbox{$|\,{\mu^{x+1}}\,\rangle$}$
(4.77) $\displaystyle=$
$\displaystyle(d^{\ell}_{1}+d^{r}_{1})\mbox{$|\,{\mu^{x}}\,\rangle$}-d^{\ell}_{1}\mbox{$|\,{\mu^{x-1}}\,\rangle$}-d^{r}_{1}\mbox{$|\,{\mu^{x+1}}\,\rangle$}$

if and only if the conditions (2.65) and (2.77) are satisfied with $N=1$.

##### Shock at the left boundary:

From the eigenvalue properties in Lemma 4.1 one gets

$\displaystyle H\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle(h^{-}_{L_{-}}+\tilde{h}_{L_{-},L_{-}+1}+h^{+}_{L_{+}})\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}$
(4.78) $\displaystyle=$
$\displaystyle[(h^{-}_{1}+h)\mbox{$|\,{\rho_{1}^{\star}}\,\rangle$}\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}]\otimes\mbox{$|\,{\rho_{1}}\,\rangle$}^{\otimes{L-2}}+\epsilon_{+}(\rho_{1})\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}.$

To handle the left boundary term in (4.78) use case B in the projection lemma
4.3 by choosing in (4.42) $a^{(2)}=\rho_{0}$, $b^{(2)}=\rho_{1}^{\star}$, and
$\tilde{a}^{(2)}=\rho_{1}$ and by choosing in (4.30) $a=\rho_{1}^{\star}$
taking the Kronecker product with $|\,{\rho_{1}}\,\rangle$ on the right. This
yields

$\displaystyle h\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(r+\tilde{d}(\rho_{1}^{\star},\tilde{b}^{(2)})\right)\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
(4.79)
$\displaystyle+d(\rho_{1},\rho_{0})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}-d(\rho_{1}^{\star},\tilde{b}^{(2)})\mbox{$|\,{\tilde{b}^{(2)}\rho_{1}}\,\rangle$}$

and

$\displaystyle\tilde{h}_{1}\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(\tilde{d}_{-}(\rho_{1}^{\star},b)-\tilde{d}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}-\left(d_{-}(\rho_{1}^{\star},b)-d(\rho_{1}^{\star},b)\right)\mbox{$|\,{b\rho_{1}}\,\rangle$}$
(4.80)

Taking without loss of generality $b=\tilde{b}^{(2)}\neq\rho_{1}^{\star}$ one
finds

$\displaystyle(h+\tilde{h}_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(r+\tilde{d}_{-}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
(4.81)
$\displaystyle+d(\rho_{1},\rho_{0})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}-d_{-}(\rho_{1}^{\star},b)\mbox{$|\,{b\rho_{1}}\,\rangle$}.$

Observe next that for $\rho_{1}^{\star}$ as specified in the Theorem one has
$\epsilon_{-}(\rho_{1}^{\star})=0$ so that
$\tilde{d}_{-}(\rho_{1}^{\star},b)=d_{-}(\rho_{1}^{\star},b)=0$ for any vector
$b\neq\rho_{1}^{\star}$. Therefore

$\displaystyle(h+\tilde{h}_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}$
$\displaystyle=$ $\displaystyle
r\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}+d(\rho_{1},\rho_{0})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}.$
(4.82)

Using (4.67) and (4.67) the conclusion from these computations is the left
closure relation

$H\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}=[r+\epsilon_{+}(\rho_{1})]\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}+d(\rho_{1},\rho_{0})\mbox{$|\,{\mu^{L_{-}+1}}\,\rangle$}=d^{r}_{1}\left(\mbox{$|\,{\mu^{L_{-}}}\,\rangle$}-\mbox{$|\,{\mu^{L_{-}+1}}\,\rangle$}\right)$
(4.83)

if and only if $\rho_{1}^{\star}=\alpha/(\alpha+\gamma)$ as specified in the
Theorem.

##### Shock at the right boundary:

From Lemma 4.1 one gets in a similar fashion

$\displaystyle H\mbox{$|\,{\mu^{L_{+}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle(h^{-}_{L_{-}}+\tilde{h}_{L_{+}-1,L_{+}}+h^{+}_{L_{+}})\mbox{$|\,{\mu^{L_{+}}}\,\rangle$}$
(4.84) $\displaystyle=$
$\displaystyle\epsilon_{-}(\rho_{0})\mbox{$|\,{\mu^{L_{+}}}\,\rangle$}+\mbox{$|\,{\rho_{0}}\,\rangle$}^{\otimes(L-2)}\otimes[(h+h^{+}_{2})\mbox{$|\,{\rho_{0}}\,\rangle$}\otimes\mbox{$|\,{\rho^{\star}}\,\rangle$}]$

To treat (4.84) we use the projection formula (4.40) in Lemma 4.3 with the
replacements $a^{(1)}\to\rho_{0}$, $b^{(1)}\to\rho_{1}^{\star}$,
$\tilde{a}^{(1)}\to\rho_{1}$ and take $a=\rho_{1}^{\star}$ in (4.38). After
left Kronecker multiplication with $|\,{\rho_{0}}\,\rangle$ this yields

$\displaystyle h\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(\ell-\tilde{d}(\rho_{1}^{\star},\tilde{b}^{(1)})\right)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
(4.85)
$\displaystyle-d(\rho_{0},\rho_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}+d(\rho_{1}^{\star},\tilde{b}^{(1)})\mbox{$|\,{\rho_{0}\tilde{b}^{(1)}}\,\rangle$}$

and

$h^{+}_{2}\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}=\left(\tilde{d}(\rho_{1}^{\star},b)+d_{+}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{0}\rho^{\star}}\,\rangle$}-\left(d(\rho_{1}^{\star},b)+d_{+}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{0}b}\,\rangle$}$
(4.86)

With $b=\tilde{b}^{(1)}\neq\rho_{0}$ one thus finds the projection property

$\displaystyle(h+h^{+}_{2})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(\ell-\tilde{d}(\rho_{1}^{\star},\tilde{b}^{(1)})+\tilde{d}(\rho_{1}^{\star},b)+d_{+}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
(4.87)
$\displaystyle-d(\rho_{0},\rho_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}+d(\rho_{1}^{\star},\tilde{b}^{(1)})\mbox{$|\,{\rho_{0}\tilde{b}^{(1)}}\,\rangle$}$
$\displaystyle-\left(d(\rho_{1}^{\star},b)+d_{+}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{0}b}\,\rangle$}$

which with (4.67) and (4.67) and by choosing without loss of generality
$b=\tilde{b}^{(1)}\neq\rho_{1}^{\star}$ leads to

$\displaystyle(h+h^{+}_{2})\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(\ell+d_{+}(\rho_{1}^{\star},b)\right)\mbox{$|\,{\rho_{0}\rho_{1}^{\star}}\,\rangle$}$
(4.88)
$\displaystyle-d(\rho_{0},\rho_{1})\mbox{$|\,{\rho_{1}^{\star}\rho_{1}}\,\rangle$}-d_{+}(\rho_{1}^{\star},b)\mbox{$|\,{\rho_{0}b}\,\rangle$}$

Since by assumption 2.66 holds one has $\epsilon_{+}(\rho_{1})=0$ so that
$d_{+}(\rho_{1}^{\star},b)=0$ for all $b\neq\rho_{1}^{\star}$. With (4.64) we
arrive at the right closure relation

$\displaystyle H\mbox{$|\,{\mu^{L_{+}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle
d^{\ell}_{1}\left(\mbox{$|\,{\mu^{L_{+}}}\,\rangle$}-\mbox{$|\,{\mu^{L_{+}-1}}\,\rangle$}\right).$
(4.89)

Conversely, if 2.66 does not hold then closure is violated.

It is thus proved that the set of vectors $|\,{\mu^{x}}\,\rangle$,
$x\in\Lambda_{L}$ forms an invariant subspace under the action of $H$ if and
only of the conditions (i) - (iii) hold. This allows for writing (4.72) with
coefficients

$G_{yx}=\left\\{\begin{array}[]{ll}d^{r}_{1}\left(\delta_{y,L_{-}}-\delta_{y,L_{-}+1}\right)&x=L_{-}\\\
(d^{\ell}_{1}+d^{r}_{1})\delta_{x,y}-d^{\ell}_{1}\delta_{x,y+1}-d^{r}_{1}\delta_{x,y-1}&L_{-}<x<L_{+}\\\
d^{\ell}_{1}\left(\delta_{y,L_{+}}-\delta_{y,L_{+}-1}\right)&x=L_{+}\end{array}\right.$
(4.90)

Comparing with the intensity matrix $Q$ of the biased random walk of 2.11 one
finds $G_{xy}=Q_{xy}$. Thus with $H=-W^{T}$ and
$S_{x\boldsymbol{\eta}}:=\pi^{-1}(x)R(x,\boldsymbol{\eta})=\mu^{x}(\boldsymbol{\eta})$
one gets $SW=QS$ which according to (2.11) in the proof of Theorem 2.5 is
equivalent to (3.4) due to reversibility of the shock exclusion process w.r.t.
the measure $\pi(x)=d_{1}^{x}$ of the biased random walk. $\Box$

####  4.3.2 Proof of Theorem 3.3

The assertion (3.5) as such, i.e., without specifying the transition
probability, is an immediate corollary of the time-reversed duality
established in Theorem (3.2). The specific formula (3.6) follows from the fact
that the intensity matrix for a biased simple random walk with constant rates
is a tridiagonal Toeplitz matrix. For readers not familiar with such matrices
we outline how to arrive at (3.6): (a) The transition probability is by
definition the matrix element $P(y,t|x,0)=(\exp{Gt})_{x,y}$ of the exponential
of the intensity matrix. (b) The non-symmetric intensity matrix can be
symmetrized by the ground state transformation
$\tilde{G}=\hat{\pi}G\hat{\pi}^{-1}$ where $\hat{\pi}$ is the diagonal matrix
with elements $d_{1}^{x-L_{-}}$ on the diagonal. (c) Diagonalizing the
tridiagonal Toeplitz matrix $G$ (which is textbook) then yields the
eigenvalues $\epsilon_{0}=0$ and $\epsilon_{p}$ (3.8) for $1\leq p\leq L-1$
and the associated eigenvectors of $\tilde{G}$ with components that can be
expressed as

$\displaystyle\Psi_{0}(y)$ $\displaystyle=$
$\displaystyle\sqrt{\frac{d_{1}^{2}-1}{d_{1}^{2L}-1}}d_{1}^{y-L_{-}}$ (4.91)
$\displaystyle\Psi_{p}(y)$ $\displaystyle=$
$\displaystyle\sqrt{\frac{2}{L}}\frac{\psi_{p}(y)}{d_{1}-\mathrm{e}^{-\frac{i\pi
p}{L}}}$ (4.92)

and satisfy the orthogonality relations

$\sum_{y=-L_{-}}^{L_{+}}\Psi_{p}(y)\Psi^{\ast}_{q}(y)=\delta_{p,q}\quad\forall
p,q\in\\{0,1,\dots,L-1\\}.$ (4.93)

Orthogonality can be proved by expanding for $d_{1}\neq 1$ the denominator in
the power series

$\frac{1}{d_{1}-\mathrm{e}^{-\frac{i\pi
p}{L}}}=\left\\{\begin{array}[]{ll}\displaystyle\frac{1}{1-d_{1}^{-2L}}\sum_{n=0}^{2L-1}d_{1}^{-n-1}\mathrm{e}^{-\frac{in\pi
p}{L}}&d_{1}>1\\\
\displaystyle-\frac{1}{1-d_{1}^{2L}}\sum_{n=0}^{2L-1}d_{1}^{n}\mathrm{e}^{\frac{i(n+1)\pi
p}{L}}&d_{1}<1\end{array}\right.$ (4.94)

and using the orthogonality relations

$\displaystyle\frac{2}{L}\sum_{y=L_{-}}^{L_{+}}\sin{\frac{\pi
p(y+1-L_{-})}{L}}\sin{\frac{\pi q(y+1-L_{-})}{L}}$ $\displaystyle=$
$\displaystyle\sum_{k\in{\mathbb{Z}}}\left(\delta_{p-q,2kL}-\delta_{p+q,2kL}\right)$
(4.95) $\displaystyle\frac{2}{L}\sum_{p=1}^{L-1}\sin{\frac{\pi
p(x+1-L_{-})}{L}}\sin{\frac{\pi p(y+1-L_{-})}{L}}$ $\displaystyle=$
$\displaystyle\sum_{k\in{\mathbb{Z}}}\left(\delta_{x-y,2kL}-\delta_{x+y,2kL}\right)$
(4.96)

and similar addition formulae involving three trigonometric functions. (d) It
then follows that $P(y,0|x,0)=\delta_{y,x}$ and eigenvector decomposition of
the matrix element $(\exp{Gt})_{x,y}$ leads to (3.6). Finally, we note that
existence of the limit $t\to\infty$ can be read off the fact that the
eigenvalues (3.8) of the intensity matrix are strictly positive for all
$p\in\\{1,\dots,L-1\\}$. Uniqueness of the invariant measure follows from
ergodicity of the process. $\Box$

###  4.4 Reverse duality involving $N$ shocks

The strategy of the proof is analogous to the case of one shock which can be
used as it is for shock measures that do not have neighbouring shock
positions.

####  4.4.1 Proof of Theorem 3.7

##### Isolated shocks

If all shocks are isolated, i.e., there is at least one site with density
$\rho_{i}$ between consecutive shocks with densities $\rho^{\star}_{i}$
$\rho^{\star}_{i+1}$, then according to (A) and (B) of the projection lemma
4.3 the generator $H$ acts on the shock measure by acting on all shocks
separately as in the case of a single shock. In other words, for $N$ shocks at
positions $x_{i}$ such that $x_{i+1}-x_{i}>1$ for all $i$ one gets

$\displaystyle H\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle\left(h^{-}_{L_{-}}+\sum_{i=1}^{N}\left(\tilde{h}_{x_{i-1}}+\tilde{h}_{x_{i}}\right)+h^{+}_{L_{+}}\right)\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$
(4.98) $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\left[d^{\ell}_{i}\left(1-\delta_{x_{1},L_{-}}\right)\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)\right.$
$\displaystyle\left.+d^{r}_{i}\left(1-\delta_{x_{N},L_{+}}\right)\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)\right]$

which implies invariance of the subspace with $N$ shocks with bulk densities
satisfying the microscopic stability condition (2.69) and arbitrary shock
densities provided that no shocks at nearest neighbour sites appear in the sum
over the shock positions.

The new situation for which closure needs to be proved are bunches of
neighbouring shocks. We call $n$ shocks on consecutive sites
$x_{k},x_{k+1},\dots,x_{k+n-1}$ with shock densities
$\rho^{\star}_{i+1},\rho^{\star}_{i+2},\dots,\rho^{\star}_{i+n}$ and bulk
densities $\rho_{i}$ on site $x_{k-1}$ and $\rho_{i+n}$ on site $x_{k+n}$
satisfying (2.77) with $N=n$ a shock bunch of size $n$ at position $x_{k}$. An
isolated shock is called a bunch of size 1. Evidently, for a total of $N$
shocks, the bunch size is in the range $1\leq n\leq N$, the shock index is in
range $0\leq i\leq N-n$, and the bunch position is in the range $L_{-}\leq
k\leq L_{+}-N+1$. The new technical ingredient to prove closure for bunches is
assertion (C) in Lemma 4.3.

##### Bulk bunches of size $n$:

For fixed $n$ we define the product vectors

$|\,{\varXi_{i}}\,\rangle$ $\displaystyle:=$
$|\,{\rho_{i}\rho^{\star}_{i+1}\rho^{\star}_{i+2}\dots\rho^{\star}_{i-1+n}\rho^{\star}_{i+n}\rho_{i+n}}\,\rangle$
(4.99) $|\,{\varOmega_{i}}\,\rangle$ $\displaystyle:=$
$\displaystyle(\tilde{h}_{12}+\sum_{k=2}^{n}\tilde{h}_{kk+1}+\tilde{h}_{n+1,n+2})\mbox{$|\,{\varXi_{i}}\,\rangle$}$
(4.100) $|\,{\varXi^{\ell}_{i+1}}\,\rangle$ $\displaystyle:=$
$|\,{\rho^{\star}_{i+1}\rho_{i+1}\rho^{\star}_{i+2}\dots\rho^{\star}_{i-1+n}\rho^{\star}_{i+n}\rho_{i+n}}\,\rangle$
(4.101) $|\,{\varXi^{r}_{i+n}}\,\rangle$ $\displaystyle:=$
$|\,{\rho_{i}\rho^{\star}_{i+1}\rho^{\star}_{i+2}\dots\rho^{\star}_{i-1+n}\rho_{i-1+n}\rho^{\star}_{i+n}}\,\rangle$
(4.102) $|\,{\varUpsilon_{ij}}\,\rangle$ $\displaystyle:=$
$\displaystyle\mbox{$|\,{\rho_{i}\rho^{\star}_{i+1}\dots\rho^{\star}_{j-1}c_{j}\rho^{\star}_{j+1}\dots\rho^{\star}_{i-1+n}\rho_{i+n}}\,\rangle$},\quad
i<j\leq i+n$ (4.103) $|\,{\tilde{\varUpsilon}_{ij}}\,\rangle$
$\displaystyle:=$
$\displaystyle\mbox{$|\,{\rho_{i}\rho^{\star}_{i+1}\dots\rho^{\star}_{j-1}\tilde{c}_{j}\rho^{\star}_{j+1}\dots\rho^{\star}_{i-1+n}\rho_{i+n}}\,\rangle$},\quad
i<j\leq i-1+n.$ (4.104)

To compute $|\,{\varOmega_{i}}\,\rangle$ we use Lemma 4.3, specifically
assertion (A) for the action of $\tilde{h}_{12}$, assertion (B) for the action
of $\tilde{h}_{n+1,n+2}$, and assertion (A) for the action of the remaining
local bulk jump operators $\tilde{h}_{kk+1}$ in the summation symbol. Setting
$a^{(1)}=\rho_{i}$, $\tilde{a}^{(1)}=\rho_{i+1}$,
$b^{(1)}=\rho^{\star}_{i+1}$, $\tilde{b}^{(1)}=c_{i+1}$ in (A),
$a=\rho^{\star}_{k}$, $\tilde{a}=\rho^{\star}_{k+1}$,
$\tilde{c}=\tilde{c}_{k}$, $c=c_{k+1}$ for $i+1\leq k\leq i+n-1$ in (C), and
$a^{(2)}=\rho_{i+n-1}$, $\tilde{a}^{(2)}=\rho_{i+n}$,
$b^{(2)}=\rho^{\star}_{i+n}$, $\tilde{b}^{(2)}=\tilde{c}_{i+n}$ this yields

$|\,{\varOmega_{i}}\,\rangle$ $\displaystyle=$
$\displaystyle\left(r+\ell+\sum_{k=i+1}^{i+n}\left(\tilde{d}(\rho^{\star}_{k},\tilde{c}_{k})-\tilde{d}(\rho^{\star}_{k},c_{k})\right)\right)\mbox{$|\,{\varXi_{i}}\,\rangle$}$
(4.105)
$\displaystyle-d(\rho_{i},\rho_{i+1})\mbox{$|\,{\varXi^{\ell}_{i+1}}\,\rangle$}+d(\rho_{i+n},\rho_{i+n-1})\mbox{$|\,{\varXi^{r}_{i+n}}\,\rangle$}$
$\displaystyle+\sum_{k=i+1}^{i+n}\left(d(\rho^{\star}_{k},c_{k})\mbox{$|\,{\varUpsilon_{ik}}\,\rangle$}-d(\rho^{\star}_{k},\tilde{c}_{k})\mbox{$|\,{\tilde{\varUpsilon}_{ik}}\,\rangle$}\right)$

Without loss of generality we choose $c_{k}=\tilde{c}_{k}\neq\rho^{\star}_{k}$
so that the sums vanish and using (4.64) and (4.67) one finds

$|\,{\varOmega_{i}}\,\rangle$ $\displaystyle=$
$\displaystyle\left(r+\ell\right)\mbox{$|\,{\varXi_{i}}\,\rangle$}$ (4.106)
$\displaystyle-d^{\ell}_{i+1}\mbox{$|\,{\varXi^{\ell}_{i+1}}\,\rangle$}-d^{r}_{i+n}\mbox{$|\,{\varXi^{r}_{i+n}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left[d^{\ell}_{i+1}+d^{r}_{i+1}+(r-\ell)(\rho_{i}-\rho_{i-1})\right]\mbox{$|\,{\varXi_{i}}\,\rangle$}$
$\displaystyle-d(\rho_{i},\rho_{i+1})\mbox{$|\,{\varXi^{\ell}_{i+1}}\,\rangle$}+d(\rho_{i+n},\rho_{i+n-1})\mbox{$|\,{\varXi^{r}_{i+n}}\,\rangle$}$

if and only if (2.69) which is equivalent to (4.44) holds for all consecutive
shock densities and all consecutive bulk densities.

Applying this relation to a measure $|\,{\mu^{\mathbf{x}}}\,\rangle$ with $N$
shocks at bulk positions $L_{-}<x_{1}<\dots<x_{N}<L_{+}$ arranged in $k$
bunches one thus finds an invariant subspace under the action of the generator
by noting

$\displaystyle H\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle[\epsilon_{-}(\rho_{0})+\epsilon_{+}(\rho_{N})]\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$
(4.107)
$\displaystyle+\sum_{i=1}^{N}(1-\delta_{x_{i-1},x_{i}-1})\left([d^{\ell}_{i}-(r-\ell)\rho_{i-1}]\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-d^{\ell}_{i}\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)$
$\displaystyle+\sum_{i=1}^{N}(1-\delta_{x_{i+1},x_{i}+1})\left([d^{r}_{i}+(r-\ell)\rho_{i}]\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-d^{r}_{i}\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}d^{\ell}_{i}(1-\delta_{x_{i-1},x_{i}-1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)$
$\displaystyle+\sum_{i=1}^{N}d^{r}_{i}(1-\delta_{x_{i+1},x_{i}+1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)$
$\displaystyle+(r-\ell)\left(\delta_{x_{1},L_{-}}\rho_{0}-\delta_{L_{+},x_{N}}\rho_{N}\right)\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$

where in the first equality the bunch property and in the second equality we
have used cancellation of the terms in the telescopic sum over the diagonal
density terms with the boundary eigenvalues. Observing finally that
$\delta_{x_{1},L_{-}}=\delta_{L_{+},x_{N}}=0$ for
$L_{-}<x_{1}<\dots<x_{N}<L_{+}$ we arrive at

$\displaystyle H\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}d^{\ell}_{i}(1-\delta_{x_{i-1},x_{i}-1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)$
(4.108)
$\displaystyle+\sum_{i=1}^{N}d^{r}_{i}(1-\delta_{x_{i+1},x_{i}+1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)$
$\displaystyle=$
$\displaystyle\sum_{i=1}^{N}\left[w^{\ell}_{i}\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)+w^{r}_{i}\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)\right]$
(4.109)

for shocks not at the boundary sites. It remains to prove closure when one or
two shocks are at the boundary sites $L_{\pm}$.

##### Bunch of shocks at the boundaries:

We first consider a bunch of $n$ shocks at the left boundary but no shock at
the right boundary, i.e., $x_{1}=L_{-}$ and $x_{N}<L_{+}$. For fixed $n$ we
define the product vectors

$|\,{\varXi^{-}}\,\rangle$ $\displaystyle:=$
$|\,{\rho^{\star}_{1}\rho^{\star}_{2}\dots\rho^{\star}_{n-1}\rho^{\star}_{n}\rho_{n}}\,\rangle$
(4.110) $|\,{\varOmega^{-}}\,\rangle$ $\displaystyle:=$
$\displaystyle(h^{-}_{1}+\sum_{k=1}^{n-1}\tilde{h}_{kk+1}+\tilde{h}_{n,n+1})\mbox{$|\,{\varXi_{-}}\,\rangle$}$
(4.111) $|\,{\varXi^{-}_{n}}\,\rangle$ $\displaystyle:=$
$|\,{\rho^{\star}_{1}\rho^{\star}_{2}\dots\rho^{\star}_{n-1}\rho_{n-1}\rho^{\star}_{n}}\,\rangle$
(4.112) $|\,{\varUpsilon^{-}_{j}}\,\rangle$ $\displaystyle:=$
$\displaystyle\mbox{$|\,{\rho^{\star}_{1}\dots\rho^{\star}_{j-1}c_{j}\rho^{\star}_{j+1}\dots\rho^{\star}_{n-1}\rho^{\star}_{n}\rho_{n}}\,\rangle$},\quad
1\leq j\leq n$ (4.113) $|\,{\tilde{\varUpsilon}^{-}_{j}}\,\rangle$
$\displaystyle:=$
$\displaystyle\mbox{$|\,{\rho^{\star}_{1}\dots\rho^{\star}_{j-1}\tilde{c}_{j}\rho^{\star}_{j+1}\dots\rho^{\star}_{n-1}\rho^{\star}_{n}\rho_{n}}\,\rangle$},\quad
1\leq j\leq n$ (4.114)

To compute $|\,{\varOmega^{-}}\,\rangle$ we use assertion (B) in Lemma 4.3 for
the action of $\tilde{h}_{n,n+1}$ and assertion (C) for the action of the
remaining local bulk jump operators $\tilde{h}_{kk+1}$ in the summation
symbol. With $a=\rho^{\star}_{1}$ and $b=c_{1}$ for the boundary terms,
$a=\rho^{\star}_{k}$, $\tilde{a}=\rho^{\star}_{k+1}$, $c=c_{k+1}$, and
$\tilde{c}=\tilde{c}_{k}$ for assertion (C), and $b^{(2)}=\rho^{\star}_{n}$,
$\tilde{a}^{(2)}=\rho_{n}$, $a^{(2)}=\rho_{n-1}$,
$\tilde{b}^{(2)}=\tilde{c}_{n}$ for assertion (B) this yields

$|\,{\varOmega^{-}}\,\rangle$ $\displaystyle=$
$\displaystyle\left(\tilde{d}_{-}(\rho^{\star}_{1},c_{1})-\tilde{d}(\rho^{\star}_{1},c_{1})\right)\mbox{$|\,{\varOmega^{-}}\,\rangle$}$
(4.116)
$\displaystyle-\left(d_{-}(\rho^{\star}_{1},c_{1})-d(\rho^{\star}_{1},c_{1})\right)\mbox{$|\,{\varUpsilon^{-}_{1}}\,\rangle$}$
$\displaystyle+\sum_{k=1}^{n-1}\left(\tilde{d}(\rho^{\star}_{k},\tilde{c}_{k})-\tilde{d}(\rho^{\star}_{k+1},c_{k+1})\right)\mbox{$|\,{\varOmega^{-}}\,\rangle$}$
$\displaystyle-\sum_{k=1}^{n-1}d(\rho^{\star}_{k},\tilde{c}_{k})\mbox{$|\,{\tilde{\varUpsilon}^{-}_{k}}\,\rangle$}+\sum_{k=1}^{n-1}d(\rho^{\star}_{k+1},c_{k+1})\mbox{$|\,{\varUpsilon^{-}_{k+1}}\,\rangle$}$
$\displaystyle+\left(r+\tilde{d}(\rho^{\star}_{n},\tilde{c}_{n})\right)\mbox{$|\,{\varOmega^{-}}\,\rangle$}$
$\displaystyle+d(\rho_{n},\rho_{n-1})\mbox{$|\,{\varXi^{-}_{n}}\,\rangle$}-d(\rho^{\star}_{n},\tilde{c}_{n})\mbox{$|\,{\tilde{\varUpsilon}^{-}_{n}}\,\rangle$}$
$\displaystyle=$
$\displaystyle\left(r+\tilde{d}_{-}(\rho^{\star}_{1},c_{1})\right)\mbox{$|\,{\varOmega^{-}}\,\rangle$}-d_{-}(\rho^{\star}_{1},c_{1})\mbox{$|\,{\varUpsilon^{-}_{1}}\,\rangle$}+d(\rho_{n},\rho_{n-1})\mbox{$|\,{\varXi^{-}_{n}}\,\rangle$}$
$\displaystyle+\sum_{k=1}^{n}\left(d(\rho^{\star}_{k},c_{k})\mbox{$|\,{\varUpsilon^{-}_{k}}\,\rangle$}-d(\rho^{\star}_{k},\tilde{c}_{k})\mbox{$|\,{\tilde{\varUpsilon}^{-}_{k}}\,\rangle$}\right).$

Now we choose use without loss of generality $\tilde{c}_{k}=c_{k}$ and observe
that by assumption $\rho^{\star}_{1}=\alpha/(\alpha+\gamma)$ so that
$d_{-}(\rho^{\star}_{1},c_{1})=\tilde{d}_{-}(\rho^{\star}_{1},c_{1})=0$ for
any choice of $c_{1}$. Therefore

$|\,{\varOmega^{-}}\,\rangle$ $\displaystyle=$ $\displaystyle
r\mbox{$|\,{\varOmega^{-}}\,\rangle$}+d(\rho_{n},\rho_{n-1})\mbox{$|\,{\varXi^{-}_{n}}\,\rangle$}$
(4.117) $\displaystyle=$
$\displaystyle[d^{r}_{n}+(r-\ell)\rho_{n}]\mbox{$|\,{\varOmega^{-}}\,\rangle$}-d^{r}_{n}\mbox{$|\,{\varXi^{-}_{n}}\,\rangle$}$

where in the second equality (4.64) was used and one gets for consecutive
bunches with a total of $N$ shocks

$\displaystyle H\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{N}d^{\ell}_{i}(1-\delta_{x_{i-1},x_{i}-1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i-}}}\,\rangle$}\right)$
(4.118)
$\displaystyle+\sum_{i=1}^{N}d^{r}_{i}(1-\delta_{x_{i+1},x_{i}+1})\left(\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}-\mbox{$|\,{\mu^{\mathbf{x}^{i+}}}\,\rangle$}\right)$
$\displaystyle-(r-\ell)\rho_{0}\left(1-\delta_{x_{1},L_{-}}\right)\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}.$

Since for a shock at the left boundary site one has $1-\delta_{x_{1},L_{-}}=0$
one recovers (4.109) also for $x_{1}=L_{-}$. For a bunch of shocks only at the
right boundary or two bunches on both boundaries similar computations using
(2.66) with $N=m+1$ also lead to (4.109) which is therefore valid for any set
of sites $\mathbf{x}$.

Hence invariance of the subspace spanned by the vectors
$|\,{\mu^{\mathbf{x}}}\,\rangle$ is proved. The closure relation (4.109) can
be recast as

$H\mbox{$|\,{\mu^{\mathbf{x}}}\,\rangle$}=-\sum_{\mathbf{y}}Q_{\mathbf{x}\mathbf{y}}\mbox{$|\,{\mu^{\mathbf{y}}}\,\rangle$}$
(4.119)

with the matrix elements $Q_{\mathbf{x}\mathbf{y}}$ (4.68) of the intensity
matrix for $N$-particle shock exclusion process. With the same arguments as at
the end of the proof of Theorem 3.4 one arrives at (3.11). $\Box$

####  4.4.2 Proof of Theorem 3.9

The time evolution formula (3.12) follows from Theorem 2.5. Existence of the
limit (3.13) follows trivially from the finite state space of the ASEP defined
on $\Lambda_{L}$ (which guarantees existence of an invariant measure) and
uniqueness follows from ergodicity which is guaranteed by the strict
positivity of the boundary rates. $\Box$

## Acknowledgements

GMS thanks the Galileo-Galilei Institute in Florence for kind hospitality and
support during the scientific program on “Randomness, Integrability, and
Universality” and particularly C. Giardinà and J. de Gier for inspriring
discussions and pointing out useful references. This work is financially
supported by FCT/Portugal through CAMGSD, IST-ID, Projects UIDB/04459/2020 and
UIDP/04459/2020.

##  Appendix A The XXZ quantum spin chain with non-diagonal boundary fields

We denote by $\mathds{1}$ the two-dimensional unit matrix and by $\mathbf{1}$
the unit matrix of dimension $2^{L}$. We also recall the standard definitions
the Pauli matrices and some of their linear combinations.

###### Definition A.1

The Pauli matrices are the $2\times 2$ matrices

$\sigma^{x}=\left(\begin{array}[]{cc}0&1\\\
1&0\end{array}\right),\quad\sigma^{y}=\left(\begin{array}[]{cc}0&-i\\\
i&0\end{array}\right),\quad\sigma^{z}=\left(\begin{array}[]{cc}1&0\\\
0&-1\end{array}\right)$ (A.1)

and

$\sigma^{\pm}=\frac{1}{2}\left(\sigma^{x}\pm i\sigma^{y}\right)$ (A.2)

are the particle annihilation and creation operators. The diagonal matrices

$\hat{n}=\frac{1}{2}\left(\mathds{1}-\sigma^{z}\right),\quad\hat{v}=\frac{1}{2}\left(\mathds{1}+\sigma^{z}\right)$
(A.3)

are the particle and vacancy projectors respectively.

With these notations the negative transpose $H$ (2.39) of the intensity matrix
$M$ (2.36) is given by

$\displaystyle H$ $\displaystyle=$
$\displaystyle-\sum_{k=L_{-}}^{L_{+}-1}[r(\sigma^{+}_{k}\sigma^{-}_{k+1}-\hat{n}_{k}\hat{v}_{k+1})+\ell(\sigma^{-}_{k}\sigma^{+}_{k+1}-\hat{v}_{k}\hat{n}_{k+1})]$
(A.4)
$\displaystyle-[\alpha(\sigma^{-}_{L_{-}}-\hat{v}_{L_{-}})+\gamma(\sigma^{+}_{L_{-}}-\hat{n}_{L_{-}})+\delta(\sigma^{-}_{L_{+}}-\hat{v}_{L_{+}})+\beta(\sigma^{+}_{L_{+}}-\hat{n}_{L_{+}})]$

A transformation with the diagonal matrix

$\hat{Q}:=\prod_{k=L_{-}}^{L_{+}-1}q^{k\hat{n}_{k}}$ (A.5)

symmetrizes the bulk jump matrices. With the parameters

$\displaystyle\theta:=\ln{q},\quad w:=\sqrt{r\ell}$ (A.6) $\displaystyle
E_{0}:=(L-1)\cosh{(\theta)}+\frac{\alpha+\beta+\gamma+\delta}{w}$ (A.7)

and with the parametrization

$\displaystyle\alpha$ $\displaystyle=$
$\displaystyle\frac{w}{2}\sinh{\theta}\frac{\mathrm{e}^{\phi_{-}-\psi_{-}}}{\sinh{\phi_{-}}\cosh{\psi_{-}}}$
(A.8) $\displaystyle\gamma$ $\displaystyle=$
$\displaystyle\frac{w}{2}\sinh{\theta}\frac{\mathrm{e}^{\psi_{-}-\phi_{-}}}{\sinh{\phi_{-}}\cosh{\psi_{-}}}$
(A.9) $\displaystyle\beta$ $\displaystyle=$
$\displaystyle\frac{w}{2}\sinh{\theta}\frac{\mathrm{e}^{\phi_{+}-\psi_{+}}}{\sinh{\phi_{+}}\cosh{\psi_{+}}}$
(A.10) $\displaystyle\delta$ $\displaystyle=$
$\displaystyle\frac{w}{2}\sinh{\theta}\frac{\mathrm{e}^{-\phi_{+}+\psi_{+}}}{\sinh{\phi_{+}}\cosh{\psi_{+}}}$
(A.11) $\displaystyle\theta_{-}$ $\displaystyle=$
$\displaystyle\psi_{-}-\phi_{-}+\theta L_{-}$ (A.12) $\displaystyle\theta_{+}$
$\displaystyle=$ $\displaystyle\phi_{+}-\psi_{+}+\theta L_{+}$ (A.13)

one gets the quantum Hamiltonian

$\displaystyle H^{XXZ}$ $\displaystyle=$ $\displaystyle\hat{Q}^{-1}H\hat{Q}$
(A.15) $\displaystyle=$
$\displaystyle-\frac{w}{2}\left\\{\sum_{k=L_{-}}^{L_{+}-1}\left[\sigma^{x}_{k}\sigma^{x}_{k+1}+\sigma^{y}_{k}\sigma^{y}_{k+1}+\cosh{(\theta)}\sigma^{z}_{k}\sigma^{z}_{k+1}\right]-E_{0}\mathbf{1}\right.$
$\displaystyle+\frac{\sinh{\theta}}{\sinh{\phi_{-}}\cosh{\psi_{-}}}\left(\mathrm{e}^{-\theta_{-}}\sigma^{-}_{L_{-}}+\mathrm{e}^{\theta_{-}}\sigma^{+}_{L_{-}}+\sinh{\psi_{-}}\cosh{\phi_{-}}\sigma^{z}_{L_{-}}\right)$
$\displaystyle\left.+\frac{\sinh{\theta}}{\sinh{\phi_{+}}\cosh{\psi_{+}}}\left(\mathrm{e}^{-\theta_{+}}\sigma^{-}_{L_{+}}+\mathrm{e}^{\theta_{+}}\sigma^{+}_{L_{+}}-\sinh{\psi_{+}}\cosh{\phi_{+}}\sigma^{z}_{L_{+}}\right)\right\\}.$

of the XXZ spin chain with non-diagonal and non-hermitian boundary fields in
the form presented in [55] for $w=1$ and $E_{0}=0$.

In terms of the parameters (A.6) - (A.13) the boundary functions (2.62) are
given by

$\displaystyle\kappa_{+}(\alpha,\gamma)$ $\displaystyle=$
$\displaystyle\mathrm{e}^{\phi_{-}+\psi_{-}-\theta L_{-}+\theta_{-}}$ (A.16)
$\displaystyle\kappa_{+}(\beta,\delta)$ $\displaystyle=$
$\displaystyle\mathrm{e}^{\phi_{+}+\psi_{+}-\theta_{+}+\theta L_{+}}$ (A.17)

so that the condition (2.65) reads

$\phi_{-}+\psi_{-}+\phi_{+}+\psi_{+}=\theta_{+}-\theta_{-}+(2N-L+1)\theta.$
(A.18)

This is the integrability condition found in [55]. $\Box$

The special condition (2.66) reads

$\displaystyle\frac{\alpha\beta}{\gamma\delta}=\mathrm{e}^{2(\phi_{-}+\phi_{+}-\psi_{-}-\psi_{-})}=\mathrm{e}^{-2\theta(N-M)}$
(A.19)

which can be recast as

$\displaystyle\theta_{+}-\theta_{-}+\theta(N-M-L+1)=0.$ (A.20)

##  Appendix B On the parametrization of the open ASEP

With the parametrization (2.58) and (2.59) one gets

$\displaystyle\kappa_{+}(\alpha,\gamma)=z_{-}^{-1},\quad\kappa_{-}(\alpha,\gamma)=-\frac{\ell+\omega_{-}}{r+\omega_{-}}$
(B.1)
$\displaystyle\kappa_{+}(\beta,\delta)=z_{+},\quad\kappa_{-}(\beta,\delta)=-\frac{\ell+\omega_{+}}{r+\omega_{+}}.$
(B.2)

The independence of the function $\kappa_{+}(\alpha,\gamma)$ of $\omega_{-}$
can be seen by noting that $x/\rho-y/(1-\rho)=r-\ell$ and therefore
$y-x+r-\ell=xz_{-}^{-1}-yz_{-}$ so that the argument of the square root
becomes $(xz_{-}^{-1}+yz_{-})^{2}$. With the role of $\rho$ and $1-\rho$
interchanged one then also realizes that $\kappa_{+}(\beta,\delta)$ is
independent of $\omega_{+}$. Thus on the parameter manifold given by

$(r-\ell)\rho_{-}(1-\rho_{-})=\alpha(1-\rho_{-})-\gamma\rho_{-}$ (B.3)

one has $\kappa_{+}(\alpha,\gamma)=z_{-}^{-1}$ and on the manifold

$(r-\ell)\rho_{+}(1-\rho_{+})=\beta\rho_{+}-\delta(1-\rho_{+})$ (B.4)

one has $\kappa_{+}(\beta,\delta)=z_{+}$.

For $q>1$, boundary densities $\rho_{\pm}\in(0,1)$, and boundary parameters in
${\mathbb{R}}^{+}$ one has $\omega_{\pm}>-\ell$ and correspondingly
$-1<\kappa_{-}(x,y)<0$. Moreover, $\kappa_{+}(x,y)\in{\mathbb{R}}^{+}$. Then
(2.64) implies $\rho_{0}<\sigma^{\star}_{1}<\rho_{1}$. On the other hand, for
$q<1$ where a downward shock with $\rho_{0}>\rho_{1}$ is macroscopically
stable one gets $-\infty<\kappa_{-}(x,y)<-1$ and therefore
$\rho_{1}<\sigma^{\star}_{1}<\rho_{0}$. Hence either way the shock density is
between the boundary densities of the shock.

Alternatively, (2.66) can be written

$\frac{\alpha\beta}{\gamma\delta}=q^{2(N-M)}$ (B.5)

In terms of the parameters $\rho_{\pm}$ the manifold $\mathcal{B}_{N}$ defined
by (2.65) is given by

$\frac{z_{+}}{z_{-}}=q^{2N},\quad N\in{\mathbb{N}}^{+}$ (B.6)

and (2.66) defining the submanifold $\mathcal{B}^{M}_{N}$ can be written as
the constraint

$\frac{(r+\omega_{-})(r+\omega_{+})}{(\ell+\omega_{-})(\ell+\omega_{+})}=q^{2M},\quad
1\leq M\leq N$ (B.7)

on the boundary parameters $\omega_{\pm}$.

## References

  * [1] F.C. Alcaraz, M. Droz, M. Henkel, and V. Rittenberg, Reaction-Diffusion Processes, Critical Dynamics and Quantum ak16ains, Ann. Phys. 230, 250–302 (1994).
  * [2] C. Bahadoran, Hydrodynamics and Hydrostatics for a Class of Asymmetric Particle Systems with Open Boundaries, Commun. Math. Phys. 310(1), 1–24 (2012).
  * [3] Balázs, M.: Multiple shocks in bricklayers’ model. J. Stat. Phys. 117(1–2), 77–98 (2004).
  * [4] Balázs, M., Farkas, G., Kovács, P., and Rákos, A.: Random walk of second class particles in product shock measures. J. Stat. Phys. 139(2), 252–279 (2010)
  * [5] Balázs, M., Duffy, L. and Pantelli, D. q-Zero Range has Random Walking Shocks. J Stat Phys 174, 958–971 (2019).
  * [6] Belitsky, V., Schütz, G.M.: Diffusion and scattering of shocks in the partially asymmetric simple exclusion process. Electron. J. Probab. 7, paper 11, 1-21 (2002)
  * [7] V.Belitsky and G.M. Schütz, Microscopic structure of shocks and antishocks in the ASEP conditioned on low current, J. Stat. Phys. 152, 93–111 (2013).
  * [8] V Belitsky and G M Schütz, Antishocks in the ASEP with open boundaries conditioned on low current J. Phys. A: Math. Theor. 46, 295004 (2013)
  * [9] V. Belitsky and G.M. Schütz, Self-Duality for the Two-Component Asymmetric Simple Exclusion Process. J. Math. Phys. 56, 083302 (2015).
  * [10] V. Belitsky and G. M. Schütz, Self-duality and shock dynamics in the $n$-species priority ASEP, Stoch. Proc. Appl. 128, 1165–1207 (2018).
  * [11] A. Borodin, I. Corwin, Dynamic ASEP, Duality, and Continuous $q^{-1}$-Hermite Polynomials, Int. Math. Res. Notices, Vol. 2020, No. 3, 641–668 (2020).
  * [12] W. Bryc and M. Swieca, On Matrix Product Ansatz for Asymmetric Simple Exclusion Process with Open Boundary in the Singular Case, J. Stat. Phys. 177, 252–284 (2019)
  * [13] G. Carinci, C. Giardina, C. Giberti, and F. Redig, Duality for Stochastic Models of Transport, J. Stat. Phys. 152, 657–697 (2013).
  * [14] G. Carinci, C. Giardinà, F. Redig, and T. Sasamoto, A generalized Asymmetric Exclusion Process with $U_{q}(\mathfrak{sl}_{2})$ stochastic duality, Probab. Theory Relat. Fields 166, 887–933 (2016).
  * [15] G. Carinci, C. Franceschini, C. Giardinà, W. Groenevelt, and F. Redig, Orthogonal Dualities of Markov Processes and Unitary Symmetries, SIGMA 15, 053 (2019)
  * [16] G. Carinci, C. Franceschini, and W. Groenevelt, $q$-Orthogonal dualities for asymmetric particle systems Electron. J. Probab. 26: 1–38 (2021).
  * [17] S. Chakraborty, S. Pal, S. Chatterjee and M. Barma, Large compact clusters and fast dynamics in coupled nonequilibrium systems, Phys. Rev. E 93, 050102(R) (2016).
  * [18] Z. Chen, J. de Gier, and M. Wheeler, Integrable Stochastic Dualities and the Deformed Knizhnik–Zamolodchikov Equation, Int. Math. Res. Notices, Vol. 2020, No. 19, 5872–5925 (2020).
  * [19] Chung, K. L., Walsh, J. B., Markov processes, Brownian motion, and time symmetry. 2nd edition, (Springer, New York, 2005)
  * [20] M. Clincy and M.R. Evans, Phase transition in the ABC model, Phys. Rev. E 67, 066115 (2003).
  * [21] N. Crampé, E. Ragoucy, and D. Simon, Construction of a coordinate Bethe ansatz for the asymmetric simple exclusion process with open boundaries, J. Stat. Mech. 2010, P11038 (2010)
  * [22] J. de Gier and F.H.L Essler, Bethe ansatz solution of the asymmetric exclusion process with open boundaries, Phys. Rev. Lett. 95, 240601 (2005).
  * [23] J. de Gier and F. H. L. Essler, Exact spectral gaps of the asymmetric exclusion process with open boundaries, J. Stat. Mech. 2006, P12011 (2006).
  * [24] B. Derrida, S.A. Janowsky, J.L. Lebowitz, and E.R. Speer, Exact solution of the totally asymmetric simple exclusion process: Shock profiles, J. Stat. Phys. 73, 813–842 (1993).
  * [25] B. Derrida, M.R. Evans, V. Hakim, and V. Pasquier, Exact solution of a 1D asymmetric exclusion model using a matrix formulation, J. Phys. A: Math. Gen. 26, 1493–1517 (1993).
  * [26] B. Derrida, J. L. Lebowitz, E. R. Speer, Shock profiles in the asymmetric simple exclusion process in one dimension, J. Stat. Phys. 89, 135–167 (1997).
  * [27] H.-J. de Vega and A. Gonzalez-Ruiz, Boundary K-matrices for the XYZ, XXZ and XXX spin chains, J. Phys. A: Math. Gen. 27, 6129–6137 (1994).
  * [28] M. Dudziński and G.M. Schütz, Relaxation spectrum of the asymmetric exclusion process with open boundaries. J. Phys. A: Math. Gen. 33, 8351–8364 (2000).
  * [29] F.H.L. Essler and V. Rittenberg, Representations of the quadratic algebra and partially asymmetric diffusion with open boundaries J. Phys. A: Math. Gen. 29, 3375–3407 (1996).
  * [30] M.R. Evans, Y. Kafri, H.M. Koduvely, and D. Mukamel, Phase separation and coarsening in one-dimensional driven diffusive systems: Local dynamics leading to long-range Hamiltonians, Phys. Rev. E 58 2764–2778 (1998).
  * [31] Ferrari, P.A., Kipnis, C., Saada, E.: Microscopic Structure of Travelling Waves in the Asymmetric Simple Exclusion Process. Ann. Probab., 19(1) 226–244 (1991)
  * [32] P.A. Ferrari and L.R.G. Fontes, Shock fluctuations in the asymmetric simple exclusion process, Probab. Theory Relat. Fields 99, 305–319 (1994)
  * [33] P.A. Ferrari, L.R.G. Fontes, and M.E. Vares, The asymmetric simple exclusion model with multiple shocks, Ann. Inst. H. Poincaré Probab. Stat. 36(2), 109–126 (2000)
  * [34] R. Frassek, C. Giardinà, and J. Kurchan, Duality and hidden equilibrium in transport models, SciPost Phys. 9, 054 (2020)
  * [35] C. Giardinà, J. Kurchan, F. Redig, and K. Vafayi, Duality and Hidden Symmetries in Interacting Particle Systems, J. Stat. Phys. 135, 25–55 (2009).
  * [36] M. Henkel, E. Orlandini, and G.M. Schütz, Equivalences between stochastic processes, J. Phys. A: Math. Gen. 28, 6335–6344 (1995).
  * [37] T. Inami and H. Konno, Integrable XYZ Spin Chain with Boundaries, J. Phys. A: Math Gen. 27, L913–L918 (1994).
  * [38] F. H. Jafarpour and S. R. Masharian, Matrix Product Steady States as Superposition of Product Shock Measures in 1D Driven Systems J. Stat. Mech. P10013 (2007)
  * [39] F. H. Jafarpour and S. R. Masharian, Temporal evolution of product shock measures in TASEP with sublattice-parallel update, Phys. Rev. E 79, 051124 (2009)
  * [40] S. Jansen and N. Kurt, On the notion(s) of duality for Markov processes, Prob. Surveys 11, 59–120 (2014).
  * [41] Y. Kafri, E. Levine, D. Mukamel, G.M. Schütz, and R.D. Willmann Phase-separation transition in one-dimensional driven models, Phys. Rev. E 68, 035101(R) (2003)
  * [42] A.B. Kolomeisky, G.M. Schütz, E.B. Kolomeisky and J.P. Straley, Phase diagram of one-dimensional driven lattice gases with open boundaries, J. Phys. A: Math. Gen. 31, 6911–6919 (1998).
  * [43] K. Krebs, F.H. Jafarpour, and G.M. Schütz, Microscopic structure of travelling wave solutions in a class of stochastic interacting particle systems, N. J. Phys. 5, 145.1–145.14 (2003).
  * [44] J Krug, Boundary-induced phase transitions in driven diffusive systems. Phys. Rev. Lett. 67, 1882–1885 (1991).
  * [45] J Kuan, A Multi-species ASEP(q,j) and q-TAZRP with Stochastic Duality, Int. Math. Res. Notices, Vol. 2018, No. 17, 5378–5416.
  * [46] J Kuan, Algebraic Symmetry and Self–Duality of an Open ASEP, Math. Phys. Anal. Geom. 24 12 (2021)
  * [47] J Kuan, Two Dualities: Markov and Schur–Weyl, Int. Math. Res. Notices, Vol. 2022, No. 13, 9633–9662 (2022).
  * [48] R. Lahiri, M. Barma, and S. Ramaswamy, Strong phase separation in a model of sedimenting lattices. Phys. Rev. E 61, 1648–1658 (2000).
  * [49] T.M. Liggett, Coupling the simple exclusion process. Ann. Probab. 4, 339–356 (1976).
  * [50] T.M. Liggett, Interacting particle systems Springer, Berlin, (1985).
  * [51] T. M. Liggett, Stochastic Interacting Systems: Contact, Voter and Exclusion Processes Springer, Berlin (1999).
  * [52] Y. Lin, Markov duality for stochastic six vertex model Electron. Commun. Probab. 24, 1-17 (2019).
  * [53] P. Lloyd, A. Sudbury, and P. Donnelly, Quantum operators in classical probability theory: I. “Quantum spin” techniques and the exclusion model of diffusion, Stoch. Proc. Appl. 61, 205–221 (1996).
  * [54] K. Mallick and S. Sandow, Finite-dimensional representations of the quadratic algebra: Applications to the exclusion process J. Phys. A: Math. Gen. 30, 4513–4526 (1997).
  * [55] R. I. Nepomechie and F. Ravanini, Completeness of the Bethe Ansatz solution of the open XXZ chain with nondiagonal boundary terms, J. Phys. A: Math. Gen. 37, 433–440 (2004).
  * [56] Rafael I. Nepomechie, Bethe Ansatz solution of the open XXZ chain with nondiagonal boundary terms, J. Phys. A: Math. Gen. 37, 433–440 (2004).
  * [57] J Ohkubo, On dualities for SSEP and ASEP with open boundary conditions, J. Phys. A: Math. Theor. 50 095004 (2017).
  * [58] F. Redig, F. Sau, Stochastic duality and eigenfunctions, https://arxiv.org/abs/1805.01318
  * [59] F. Rezakhanlou, Microscopic structure of shocks in one conservation laws Annales de l’I. H. P., section C, tome 12, no 2 (1995), p. 119–153
  * [60] S. Sandow and G.M. Schütz, On $U_{q}[SU(2)]$-Symmetric Driven Diffusion, Europhys. Lett. 27, 7–12 (1994).
  * [61] S. Sandow, Partially asymmetric exclusion process with open boundaries, Phys. Rev. E 50, 2660–2667 (1994).
  * [62] Santen, L., Appert, C. The Asymmetric Exclusion Process Revisited: Fluctuations and Dynamics in the Domain Wall Picture. J. Stat. Phys. 106, 187–199 (2002).
  * [63] G. Schütz and E. Domany, Phase transitions in an exactly soluble one- dimensional asymmetric exclusion model. J. Stat. Phys. 72, 277–296 (1993).
  * [64] G. Schütz, and S. Sandow, Non-abelian symmetries of stochastic processes: derivation of correlation functions for random vertex models and disordered interacting many-particle systems, Phys. Rev. E 49, 2726–2744 (1994).
  * [65] G.M. Schütz, Diffusion-annihilation in the presence of a driving field. J. Phys. A: Math. Gen. 28, 3405–3415 (1995).
  * [66] G.M. Schütz, Duality relations for asymmetric exclusion processes, J. Stat. Phys. 86, 1265–1287 (1997).
  * [67] G.M. Schütz, The Heisenberg chain as a dynamical model for protein synthesis - Some theoretical and experimental results, Int. J. Mod. Phys. B 11, 197–202 (1997).
  * [68] G.M. Schütz, Exactly solvable models for many-body systems far from equilibrium, in: Phase Transitions and Critical Phenomena. Vol. 19, C. Domb and J. Lebowitz (eds.), Academic Press, London (2001).
  * [69] G.M. Schütz, Duality from integrability: annihilating random walks with pair deposition, J. Phys. A: Math. Theor. 53 355003 (2020)
  * [70] D. Simon, Construction of a coordinate Bethe ansatz for the asymmetric simple exclusion process with open boundaries, J. Stat. Mech. 2009, P07017 (2009).
  * [71] H. Spohn, Long-range correlations for stochastic lattice gases in a non-equilibrium steady state. J. Phys. A: Math. Gen. 16, 4275–4291 (1983).
  * [72] A. Sudbury and P. Lloyd, Quantum operators in classical probability theory. II: The concept of duality in interacting particle systems, Ann. Probab. 23(4), 1816–1830 (1995).