Density and current profiles in $U_q(A^{(1)}_2)$ zero range process

TL;DR

This paper analyzes how fixed first-class particles affect the distribution and flow of second-class particles in an integrable zero range process derived from the stochastic $R$ matrix of $U_q(A^{(1)}_2)$, providing exact formulas for density and current.

Contribution

It introduces exact formulas for density and current profiles of second-class particles influenced by fixed first-class particles in an integrable zero range process.

Findings

Exact formulas for local density of second-class particles

Exact formulas for current of second-class particles

Influence of fixed first-class particles on second-class particle distribution

Abstract

The stochastic $R$ matrix for $U_q(A^{(1)}_n)$ introduced recently gives rise to an integrable zero range process of $n$ classes of particles in one dimension. For $n=2$ we investigate how finitely many first class particles fixed as defects influence the grand canonical ensemble of the second class particles. By using the matrix product stationary probabilities involving infinite products of $q$-bosons, exact formulas are derived for the local density and current of the second class particles in the large volume limit.

Equations

$$\mathbb{P}({\sigma }_1,\dots ,{\sigma }_L)=\mathrm{Tr}(X_{{\sigma }_1}\cdots X_{{\sigma }_L}),$$

$$X_{{\alpha }_1,{\alpha }_2}=\frac{(\mu ;q{)}_{{\alpha }_1+{\alpha }_2}}{(q;q{)}_{{\alpha }_1}(q;q{)}_{{\alpha }_2}}\frac{(\mu \mathbf{b};q{)}_{\infty }}{(\mathbf{b};q{)}_{\infty }}{\mathbf{k}}^{{\alpha }_2}{\mathbf{c}}^{{\alpha }_1},\frac{(\mu \mathbf{b};q{)}_{\infty }}{(\mathbf{b};q{)}_{\infty }}=\sum _{j\ge 0}\frac{(\mu ;q{)}_j}{(q;q{)}_j}{\mathbf{b}}^j,$$

$$\mathbf{b}|m⟩$$

$$A_{{\alpha }_1}=\sum _{{\alpha }_2\ge 0}{y}^{{\alpha }_2}{X}_{{\alpha }_1,{\alpha }_2}=\frac{(\mu ;q{)}_{{\alpha }_1}}{(q;q{)}_{{\alpha }_1}}\frac{(\mu \mathbf{b};q{)}_{\infty }}{(\mathbf{b};q{)}_{\infty }}(y\mathbf{k};q{)}_{\infty }^{-1}{\mathbf{c}}^{{\alpha }_1}(\mu y\mathbf{k};q{)}_{\infty }.$$

$$\displaystyle\mathrm{Tr}\bigl{(}A_{d_{1}}\cdots A_{d_{s}}A_{0}^{L-s}\bigr{)}$$

$$\mathcal{S}(\lambda ,\mu )(|\alpha ⟩\otimes |\beta ⟩)=\sum _{\gamma ,\delta \in {\mathbb{Z}}_{\ge 0}^n}\mathcal{S}(\lambda ,\mu {)}_{\alpha ,\beta }^{\gamma ,\delta }|\gamma ⟩\otimes |\delta ⟩,$$

$$\mathcal{S}(\lambda ,\mu {)}_{\alpha ,\beta }^{\gamma ,\delta }=\theta (\gamma +\delta =\alpha +\beta ){\Phi }_q(\gamma |\beta ;\lambda ,\mu ),$$

$${\Phi }_q(\gamma |\beta ;\lambda ,\mu )=q^{\phi (\beta -\gamma ,\gamma )}{\left(\frac{\mu }{\lambda }\right)}^{|\gamma |}\frac{(\lambda {)}_{|\gamma |}(\frac{\mu }{\lambda }{)}_{|\beta |-|\gamma |}}{(\mu {)}_{|\beta |}}\prod _{i=1}^n{\left(\genfrac{}{}{0ex}{}{{\beta }_i}{{\gamma }_i}\right)}_q,\phi (\alpha ,\beta )$$

$${\mathcal{S}}_{1,2}({\nu }_1,{\nu }_2){\mathcal{S}}_{1,3}({\nu }_1,{\nu }_3){\mathcal{S}}_{2,3}({\nu }_2,{\nu }_3)={\mathcal{S}}_{2,3}({\nu }_2,{\nu }_3){\mathcal{S}}_{1,3}({\nu }_1,{\nu }_3){\mathcal{S}}_{1,2}({\nu }_1,{\nu }_2),$$

$$\sum _{\gamma ,\delta \in {\mathbb{Z}}_{\ge 0}^n}\mathcal{S}(\lambda ,\mu {)}_{\alpha ,\beta }^{\gamma ,\delta }=1(\forall \alpha ,\beta \in {\mathbb{Z}}_{\ge 0}^n).$$

$$\sum _{\gamma \in {\mathbb{Z}}_{\ge 0}^n}{\Phi }_q(\gamma |\beta ;\lambda ,\mu )=1(\forall \beta \in {\mathbb{Z}}_{\ge 0}^n),$$

$$T(\lambda |{\mu }_1,\dots ,{\mu }_L)={\mathrm{Tr}}_W\left({\mathcal{S}}_{0,L}(\lambda ,{\mu }_L)\cdots {\mathcal{S}}_{0,1}(\lambda ,{\mu }_1)\right)\in \mathrm{End}(W^{\otimes L})$$

$$T(|{\beta }_1⟩\otimes \cdots \otimes |{\beta }_L⟩)=\sum _{{\alpha }_1,\dots ,{\alpha }_L\in {\mathbb{Z}}_{\ge 0}^n}{T}_{{\beta }_1,\dots ,{\beta }_L}^{{\alpha }_1,\dots ,{\alpha }_L}|{\alpha }_1⟩\otimes \cdots \otimes |{\alpha }_L⟩\in W^{\otimes L}.$$

$$T_{{\beta }_1,\dots ,{\beta }_L}^{{\alpha }_1,\dots ,{\alpha }_L}=0\text{unless}{\alpha }_1+\cdots +{\alpha }_L={\beta }_1+\cdots +{\beta }_L\in {\mathbb{Z}}_{\ge 0}^n.$$

$$[T(\lambda |{\mu }_1,\dots ,{\mu }_L),T({\lambda }^{\prime }|{\mu }_1,\dots ,{\mu }_L)]=0.$$

$$|P(t+1)⟩=T(\lambda |{\mu }_1,\dots ,{\mu }_L)|P(t)⟩\in W^{\otimes L}.$$

$$S(m)=\{({\sigma }_1,\dots ,{\sigma }_L)\in ({\mathbb{Z}}_{\ge 0}^n{)}^L\mid {\sigma }_1+\cdots +{\sigma }_L=m\},$$

$$\frac{d}{dt}|P(t)⟩$$

$$H_{+}$$

$$h_{+}(|\alpha ⟩\otimes |\beta ⟩)$$

$$-ϵ\sum _{i=0}^{|\alpha |-1}\frac{q^i}{1-\mu q^i}|\alpha ⟩\otimes |\beta ⟩,$$

$$h_{-}(|\alpha ⟩\otimes |\beta ⟩)$$

$$-ϵ\sum _{i=0}^{|\beta |-1}\frac{1}{1-\mu q^i}|\alpha ⟩\otimes |\beta ⟩,$$

$$w_{+}(\gamma |\alpha )$$

$$w_{-}(\gamma |\beta )$$

$$|\overline{P}⟩=T(\lambda |{\mu }_1,\dots ,{\mu }_L)|\overline{P}⟩.$$

$$|\overline{P}⟩=\sum _{({\sigma }_1,\dots ,{\sigma }_L)\in S(m)}\mathbb{P}({\sigma }_1,\dots ,{\sigma }_L)|{\sigma }_1,\dots ,{\sigma }_L⟩$$

$$|\overline{P}⟩=$$

$$+(2+q^2)(1+\kappa )(1+q\kappa )|\varnothing ,22,1⟩+(1+2q^2)(1+\kappa )(1+q\kappa )|\varnothing ,1,22⟩$$

$$+(1+q)(1+2q)(1+\kappa )(1+q\kappa )|\varnothing ,12,2⟩+(1+q)(1+q+q^2)(1+\kappa {)}^2|2,2,1⟩+\text{cyclic},$$

$$\begin{array}{rl}|\overline{P}⟩=& 3(1+q\kappa )(2+q+\kappa +2q\kappa )(1+q^2\kappa )|\varnothing ,\varnothing ,112⟩\\ & +(1+q)(1+\kappa )(1+q\kappa )(3+3q+3q^2+2\kappa +2q\kappa +5q^2\kappa )|\varnothing ,1,12⟩\\ & +(1+\kappa )(1+q\kappa )(3+3q+3q^2+\kappa +5q\kappa +2q^2\kappa +q^3\kappa )|\varnothing ,2,11⟩\\ & +(1+q)(1+\kappa )(1+q\kappa )(5+2q+2q^2+3\kappa +3q\kappa +3q^2\kappa )|\varnothing ,12,1⟩\\ & +(1+\kappa )(1+q\kappa )(1+2q+5q^2+q^3+3q\kappa +3q^2\kappa +3q^3\kappa )|\varnothing ,11,2⟩\\ & +(1+q)(1+q+q^2)(1+\kappa {)}^2(2+q+\kappa +2q\kappa )|2,1,1⟩+\text{cyclic}.\end{array}$$

$$\mathbf{k}\mathbf{b}=q\mathbf{b}\mathbf{k},\mathbf{k}\mathbf{c}=q^{-1}\mathbf{c}\mathbf{k},\mathbf{b}\mathbf{c}=1-\mathbf{k},\mathbf{c}\mathbf{b}=1-q\mathbf{k}.$$

$$\begin{array}{rl}\mathbf{b}|m⟩& =|m+1⟩,\mathbf{c}|m⟩=(1-q^m)|m-1⟩,\mathbf{k}|m⟩=q^m|m⟩,\\ ⟨m|\mathbf{c}& =⟨m+1|,⟨m|\mathbf{b}=⟨m-1|(1-q^m),⟨m|\mathbf{k}=⟨m|q^m,\end{array}$$

$$\mathcal{B}$$

$$X_{\alpha }=\frac{(\mu {)}_{{\alpha }_1+{\alpha }_2}}{(q{)}_{{\alpha }_1}(q{)}_{{\alpha }_2}}\frac{(\mu \mathbf{b}{)}_{\infty }}{(\mathbf{b}{)}_{\infty }}{\mathbf{k}}^{{\alpha }_2}{\mathbf{c}}^{{\alpha }_1}\text{for }\alpha =({\alpha }_1,{\alpha }_2)\in {\mathbb{Z}}_{\ge 0}^2.$$

Full text

Density and current profiles

in $\boldsymbol{U_{q}(A^{(1)}_{2})}$ zero range process

A. Kuniba

atsuo.s.kuniba@gmail.com

Institute of Physics, University of Tokyo, Komaba, Tokyo 153-8902, Japan

 and 

V. V. Mangazeev

vladimir.mangazeev@anu.edu.au

Department of Theoretical Physics, Research School of Physics and Engineering, Australian National University, Canberra, ACT 0200, Australia.

Abstract

The stochastic $R$ matrix for $U_{q}(A^{(1)}_{n})$ introduced recently gives rise to an integrable zero range process of $n$ classes of particles in one dimension. For $n=2$ we investigate how finitely many first class particles fixed as defects influence the grand canonical ensemble of the second class particles. By using the matrix product stationary probabilities involving infinite products of $q$-bosons, exact formulas are derived for the local density and current of the second class particles in the large volume limit.

1. Introduction

Zero range processes (ZRPs) [19] are stochastic particle systems on lattice modeling various flows in granules, queuing networks, traffic and so forth. Their characteristic feature is that particles are allowed to share a site and hop over the lattice with the rate that only depends on the occupancy and the list of leaving particles at the departure site111 This slightly generalizes the original terminology in that an arbitrary number of particles are allowed to jump out simultaneously. Such multiple jumps can be suppressed by setting $\mu\rightarrow 0$ in our model. See the explanation after (2.17).. To describe their hydrodynamic limit and the rich behavior like condensation has been an important issue in non-equilibrium statistical mechanics. See for example [5, 8, 9] and references therein.

In the recent work [11], a stochastic $R$ matrix for the quantum affine algebra $U_{q}(A^{(1)}_{n})$ was constructed. It gives rise to discrete and continuous time Markov processes associated with a commuting family of Markov transfer matrices. They are formulated as stochastic dynamics of $n$ classes of particles in one dimension with zero range type interaction. Many integrable Markov processes studied earlier, e.g. [2, 3, 12, 16, 18, 20, 21] can be identified with their special cases as summarized in [10, Fig.1,2]. In this paper we will be concerned with the version of the model introduced in [11, Sec.3.3, 3.4], which will be called the $U_{q}(A^{(1)}_{n})$ ZRP. When $n=1$ it reduces to the zero range chipping model introduced in [16].

Stationary states of the $U_{q}(A^{(1)}_{n})$ ZRP were obtained in [13, 14]. Let $\sigma_{i}=(\sigma_{i,1},\ldots,\sigma_{i,n})\in{\mathbb{Z}}_{\geq 0}^{n}$ be a local state, which means that there are $\sigma_{i,k}$ particles of class $k$ at the lattice site $i$. For the length $L$ periodic chain, the probability of finding the system in a given configuration $(\sigma_{1},\ldots,\sigma_{L})\in({\mathbb{Z}}_{\geq 0}^{n})^{L}$ is expressed, up to normalization, by the matrix product formula

[TABLE]

where $X_{\sigma_{i}}$ is an operator acting on the tensor product $F^{\otimes\frac{n}{2}(n-1)}$ of the $q$-boson Fock space $F=\oplus_{m\geq 0}{\mathbb{C}}(q)|m\rangle$. For $n=1$, the $X_{\sigma_{i}}$ is just a scalar meaning that the stationary measure is factorized. However it is an exceptional feature limited to $n=1$. For instance when $n=2$, the operator $X_{\alpha_{1},\alpha_{2}}$ for the local state $(\alpha_{1},\alpha_{2})\in{\mathbb{Z}}_{\geq 0}^{2}$ reads

[TABLE]

where $\mu$ is another model parameter and the symbol $(z;q)_{m}$ is the $q$-shifted factorial defined in the end of this section. The operators ${\bf b},{\bf c}$ and ${\bf k}$ are the $q$-boson creation, annihilation and the number operators acting on $F$ as

[TABLE]

For $n$ general the operator $X_{\alpha_{1},\ldots,\alpha_{n}}$ for the local state $(\alpha_{1},\ldots,\alpha_{n})\in{\mathbb{Z}}_{\geq 0}^{n}$ possesses a nested structure with respect to the rank $n$ [14]. Thus the $U_{q}(A^{(1)}_{n})$ ZRPs form the first systematic examples of multispecies (or multi-class) ZRPs whose stationary measure on the ring is not factorized. The relevant matrix product operators are also quite distinct from those in the exclusion type processes (cf. [4, 6, 15]) in that they involve quantum dilogarithm type infinite products of $q$-bosons, offering a challenge to extract physics of the model.

With this background in mind we present in this paper a modest analysis of stationary properties of the $U_{q}(A^{(1)}_{2})$ ZRP based on the matrix product formula (1.1)–(1.2). We introduce finitely many first class particles as defects and investigate their influence on the second class particles whose density is kept finite in the infinite volume limit. To motivate this setting, although somewhat technically, note that one must pick an equal number of ${\bf b}$’s and ${\bf c}$’s to get a nonzero contribution to the trace (1.1). Therefore the sum of the expansion index $j$ in (1.2) coming from $X_{\sigma_{1}},\ldots,X_{\sigma_{L}}$ in (1.1) should coincide with the total number of the first class particles. Gathering all such contributions in the limit $L\rightarrow\infty$ is a feasible task at least if the first class particles are kept finite.

The second class particles will be treated in the grand canonical ensemble, which means that $X_{\alpha_{1},\alpha_{2}}$ is effectively replaced by the generating series with respect to $\alpha_{2}$ in the fugacity $y$:

[TABLE]

See (4.1) and (3.3). We stay in the regime $0<q,\mu<1$ where there is no symptom of condensation. See the remarks around (5.9) concerning this point. The equivalence with the canonical ensemble treatment will be argued in Section 6.3. The basic quantity is the probability $P(r,m)$ that exactly $m$ second class particles are found at the site $r$ under the condition that the sites $1,2,\ldots,s$ of the periodic lattice ${\mathbb{Z}}_{L}$ of size $L$ contain $d_{1},\ldots,d_{s}$ first class particles. To avoid the ambiguity we assume $d_{1},d_{s}\geq 1$ but the choice $d_{i}=0$ is still allowed for $0<i<s$. So they form a cluster of defects of size $s$ in general. Up to normalization the conditional probability $P(r,m)$ is given by replacing the $r$ th operator from the left in

[TABLE]

by $y^{m}X_{d_{r},m}$ if $1\leq r\leq s$ and by $y^{m}X_{0,m}$ elsewhere222 Precise treatment involves a regularization as mentioned after (4.6).. Our task is to evaluate it in the infinite volume limit $L\rightarrow\infty$ with $s$ and $d_{1},\ldots,d_{s}$ kept fixed. The limit separates the periodic lattice ${\mathbb{Z}}_{L}$ into the three distinct regions I, II and III, which are the inside ($1\leq r\leq s$), the right ($r>s$), and the left ($r\leq 0$) of the defect cluster, respectively. The ${\mathbb{Z}}_{L}$ periodicity of the lattice implies that the probability $P(r,m)$ for the right region $r>s$ and the left region $r\leq 0$ should match when $r\to\infty$. This is indeed the case as seen from (6.15) and (6.19).

Once the conditional probability is determined, the local density and current of the second class particles are derived at any site $r$. They are physical quantities seen from the defects. The necessary calculations are elementary. The final results are summarized in Theorem 4, 7 and 9 for the regions I, II and III, respectively. They are expressed in terms of the $q$-digamma function and its derivative together with the functions $G_{m,l}(d_{1},\ldots,d_{s})$ (7.10) which incorporate the effect of defects. Curiously the latters are related to the monodromy matrices of the $U_{q}(A^{(1)}_{1})$ ZRP containing the fugacity $y$ as a spectral parameter. In Proposition 10 we will also show that the defects decrease the second class particles in the entire system exactly by their number $d_{1}+\cdots+d_{s}$ compared from the defect-free situation.

We present the profiles of the local density and currents in a number of figures for various values of $q,\mu,\rho$ and the defect pattern $(d_{1},\ldots,d_{s})$, where $\rho$ denotes the average density of the second class particles. The density profiles exhibit a peak and a valley at the left and the right boundaries of the defect cluster, respectively. In the current profiles the peak is not observed but other behavior is more or less similar to the density. The detail is dependent on the pattern $(d_{1},\ldots,d_{s})$ and it is not easy to provide an intuitive explanation in general. However our result captures the mode $\eta_{j}$ (5.3) controlling the correlation length, which shows that the influence of the defects reaches longer distance when the average density $\rho$ is lower. Moreover we provide especially simple profiles in the limits $\rho\rightarrow 0$ and $\rho\rightarrow\infty$ in the case of homogeneous defects in Figure 4, which makes the general case easy to infer.

We remark that analyses of density and current based on the grand canonical ensemble similar to the present paper have been done for a class of two species exclusion type processes. See for example [4, 17] and references therein.

The layout of the paper is as follows. We recall the $U_{q}(A^{(1)}_{n})$ ZRP [11] in Section 2 and the matrix product formula for the stationary probabilities for $n=2$ [13] in Section 3. The local density and current of the second class particles in the grand canonical ensemble are formulated in Section 4. We first deal with the defect-free single species case in Section 5 as a preparatory warm-up partly reproducing known facts in earlier works, e.g. [2, 5, 16]. Our main results in the presence of defects are given in Section 6 with a number of figures showing the density and current profiles. Their derivation are detailed in Section 7. Section 8 contains a brief summary and discussion. Technical lemmas are collected in Appendix A.

Throughout the paper we use the notation $\theta(\mathrm{true})=1,\theta(\mathrm{false})=0$, $(z)_{m}=(z;q)_{m}=\prod_{j=0}^{m-1}(1-zq^{j})$ and the $q$-binomial $\binom{m}{k}_{\!q}=\theta(k\in[0,m])\frac{(q)_{m}}{(q)_{k}(q)_{m-k}}$. The symbols $(z)_{m}$ appearing in this paper always mean $(z;q)_{m}$. For integer arrays $\alpha=(\alpha_{1},\ldots,\alpha_{m}),\beta=(\beta_{1},\ldots,\beta_{m})$ of any length $m$, we write $|\alpha|=\alpha_{1}+\cdots+\alpha_{m}$ and the Kronecker delta $\delta_{\alpha,\beta}=\delta^{\alpha}_{\beta}=\prod_{i=1}^{m}\theta(\alpha_{i}=\beta_{i})$. The relation $\alpha\leq\beta$ or equivalently $\beta\geq\alpha$ is defined by $\beta-\alpha\in{\mathbb{Z}}^{m}_{\geq 0}$.

2. $U_{q}(A^{(1)}_{n})$ zero range process

2.1. Stochastic $R$ matrix

Set $W=\bigoplus_{\alpha=(\alpha_{1},\ldots,\alpha_{n})\in{\mathbb{Z}}_{\geq 0}^{n}}{\mathbb{C}}|\alpha\rangle$. Define the operator $\mathscr{S}(\lambda,\mu)\in\mathrm{End}(W\otimes W)$ depending on the parameters $\lambda$ and $\mu$ by

[TABLE]

where $\Phi_{q}(\gamma|\beta;\lambda,\mu)$ is given by

[TABLE]

The sum (2.1) is finite due to the $\theta$ factor in (2.2). According to the direct sum decomposition $W\otimes W=\bigoplus_{\gamma\in{\mathbb{Z}}_{\geq 0}^{n}}\left(\bigoplus_{\alpha+\beta=\gamma}{\mathbb{C}}|\alpha\rangle\otimes|\beta\rangle\right)$, $\mathscr{S}(\lambda,\mu)$ splits into the corresponding submatrices. Note that $\mathscr{S}(\lambda,\mu)^{\gamma,\delta}_{\alpha,\beta}=0$ unless $\gamma\leq\beta$ hence $\alpha\leq\delta$ as well. The difference property $\mathscr{S}(\lambda,\mu)=\mathscr{S}(c\lambda,c\mu)$ is absent. We call $\mathscr{S}(\lambda,\mu)$ the stochastic $R$ matrix and depict its elements as

[TABLE]

The stochastic $R$ matrix was constructed [11] based on the quantum $R$ matrix of the symmetric tensor representation of the quantum affine algebra $U_{q}(A^{(1)}_{n})$. It satisfies the Yang-Baxter equation and the sum-to-unity condition [11]:

[TABLE]

The latter is a consequence of the sum rule [11]:

[TABLE]

where the summand is zero unless $\gamma\leq\beta$.

2.2. Discrete time

$U_{q}(A^{(1)}_{n})$ ZRP

For a positive integer $L$ we introduce the operator

[TABLE]

depending on the parameters $\mu_{1},\ldots,\mu_{L}$ and $\lambda$. To explain the notation, consider the space $W\otimes W^{\otimes L}$ labeled as $W_{0}\otimes W_{1}\otimes\cdots\otimes W_{L}$ for distinction ($W_{i}=W$). Then $\mathscr{S}_{0,i}(\lambda,\mu_{i})$ acts as the stochastic $R$ matrix $\mathscr{S}(\lambda,\mu_{i})$ (2.1) on $W_{0}\otimes W_{i}$ and as the identity elsewhere. The trace in (2.7) is taken over $W_{0}$ leaving an operator acting on $W_{1}\otimes\cdots\otimes W_{L}=W^{\otimes L}$.

Write the action of $T=T(\lambda|\mu_{1},\ldots,\mu_{L})$ as

[TABLE]

Then the matrix element is depicted by the concatenation of (2.4) as

[TABLE]

which is a customary diagram for row transfer matrices of vertex models on the length $L$ periodic lattice [1]. By the construction it satisfies the weight conservation:

[TABLE]

Thanks to the Yang-Baxter equation (2.5), the matrix (2.7) forms a commuting family (cf. [1]):

[TABLE]

Let $t$ be a time variable and consider the evolution equation

[TABLE]

Due to the weight conservation (2.9) it splits into finite-dimensional subspaces which we call sectors. In terms of the array $m=(m_{1},\ldots,m_{n})\in{\mathbb{Z}}^{n}_{\geq 0}$ and the set

[TABLE]

the corresponding sector, which will also be referred to as $m$, is given by $\oplus_{(\sigma_{1},\ldots,\sigma_{L})\in S(m)}{\mathbb{C}}|\sigma_{1}\rangle\otimes\cdots\otimes|\sigma_{L}\rangle$. The vector $|\sigma_{1}\rangle\otimes\cdots\otimes|\sigma_{L}\rangle\in W^{\otimes L}$ with $\sigma_{i}=(\sigma_{i,1},\ldots,\sigma_{i,n})\in{\mathbb{Z}}_{\geq 0}^{n}$ is interpreted as a state of the system where the $i$ th site from the left is populated with $\sigma_{i,a}$ particles of the $a$ th class or species. The array $m=(m_{1},\ldots,m_{n})$ indicates that there are $m_{a}$ particles of class $a$ in total in the corresponding sector.

In order to interpret (2.11) as the master equation of a discrete time Markov process, the matrix $T=T(\lambda|\mu_{1},\ldots,\mu_{L})$ should fulfill the following conditions:

1. (i)

Non-negativity; all the elements (2.8) belong to ${\mathbb{R}}_{\geq 0}$, 2. (ii)

Sum-to-unity; $\sum_{\alpha_{1},\ldots,\alpha_{L}\in{\mathbb{Z}}_{\geq 0}^{n}}T_{\beta_{1},\ldots,\beta_{L}}^{\alpha_{1},\ldots,\alpha_{L}}=1$ for any $(\beta_{1},\ldots,\beta_{L})\in({\mathbb{Z}}_{\geq 0}^{n})^{L}$.

The property (i) holds if $\Phi_{q}(\gamma|\beta;\lambda,\mu_{i})\geq 0$ for all $i\in{\mathbb{Z}}_{L}$. This is achieved by taking $0<\mu^{\epsilon}_{i}<\lambda^{\epsilon}<1,q^{\epsilon}<1$ for $\epsilon=\pm 1$. The property (ii) assures the total probability conservation and can be shown by using (2.6) in a similar manner to [11, Sec.3.2].

Henceforth we call the $T(\lambda|\mu_{1},\ldots,\mu_{L})$ Markov transfer matrix assuming $0<\mu^{\epsilon}_{i}<\lambda^{\epsilon}<1,q^{\epsilon}<1$ always. The choice of $\epsilon=\pm 1$ specifies one of the two physical regimes of the system. The evolution equation (2.11) describes a stochastic dynamics of $n$ classes of particles hopping to the right periodically via an extra lane (horizontal arrows in (2.8)) which particles get on or get off when they leave or arrive at a site. The rate of such local processes is specified by (2.2), (2.3) and (2.4). For $n=1$ and the homogeneous choice $\mu_{1}=\cdots=\mu_{L}$, it reduces to the model introduced in [16].

2.3. Continuous time

$U_{q}(A^{(1)}_{n})$ ZRP

One can derive continuous time versions of (2.11) from the homogeneous case $\mu_{1}=\cdots=\mu_{L}=\mu$ by taking the logarithmic derivative either at $\lambda=1$ or $\lambda=\mu$ [11, Sec.3.4]. The result is given by

[TABLE]

with $H_{\pm}=\sum_{i\in{\mathbb{Z}}_{L}}h_{\pm,i,i+1}$. The summands $h_{\pm,i,i+1}$ act on the adjacent $(i,i\!+\!1)$ th components of $W^{\otimes L}$ as $h_{\pm}$ and as the identity elsewhere. The pairwise interactions $h_{\pm}$ are specified by

[TABLE]

where $\varphi(\alpha,\beta)$ is defined in (2.3) and $0^{n}$ stands for $(0,\ldots,0)\in{\mathbb{Z}}_{\geq 0}^{n}$. We have included the sign factor $\epsilon=\pm 1$ to cover the two regimes. Denote the action of the matrices ${\mathcal{H}}=H_{+},H_{-}$ on the base vectors as ${\mathcal{H}}(|\beta_{1}\rangle\otimes\cdots\otimes|\beta_{L}\rangle)=\sum_{\alpha_{1},\ldots,\alpha_{L}}{\mathcal{H}}^{\alpha_{1},\ldots,\alpha_{L}}_{\beta_{1},\ldots,\beta_{L}}|\alpha_{1}\rangle\otimes\cdots\otimes|\alpha_{L}\rangle$. The equation (2.12) can be viewed as the master equation of a continuous time Markov process if the following conditions are satisfied:

(i)’ Non-negativity; ${\mathcal{H}}^{\alpha_{1},\ldots,\alpha_{L}}_{\beta_{1},\ldots,\beta_{L}}\in{\mathbb{R}}_{\geq 0}$ for any pair such that $(\alpha_{1},\ldots,\alpha_{L})\neq(\beta_{1},\ldots,\beta_{L})$,

(ii)’ Sum-to-zero; $\sum_{\alpha_{1},\ldots,\alpha_{L}}{\mathcal{H}}^{\alpha_{1},\ldots,\alpha_{L}}_{\beta_{1},\ldots,\beta_{L}}=0$ for any $(\beta_{1},\ldots,\beta_{L})$.

The latter represents the total probability conservation. It was shown in [11] that (i)’ and (ii)’ are satisfied if $0<q^{\epsilon},\mu^{\epsilon}<1$ for $\epsilon=\pm 1$.

The commutativity (2.10) leads to $[H_{+},H_{-}]=0$. Therefore they share the same eigenvectors with the superposition $H=H(a,b,\epsilon,q,\mu)=aH_{+}(\epsilon,q,\mu)+bH_{-}(\epsilon,q,\mu)$ in (2.12). A curious symmetry $H(a,b,-\epsilon,q^{-1},\mu^{-1})=\mathscr{P}H(\mu b,\mu a,\epsilon,q,\mu)\mathscr{P}^{-1}$ is known to hold [11, Rem.9], where $\mathscr{P}=\mathscr{P}^{-1}\in\mathrm{End}(W^{\otimes L})$ is the ‘parity’ operator reversing the sites as $\mathscr{P}(|\alpha_{1}\rangle\otimes\cdots\otimes|\alpha_{L}\rangle)=|\alpha_{L}\rangle\otimes\cdots\otimes|\alpha_{1}\rangle$.

The Markov processes (2.12) is naturally interpreted as the stochastic dynamics of $n$ classes of particles on the ring of length $L$. The base vector $|\alpha_{1}\rangle\otimes\cdots\otimes|\alpha_{L}\rangle$ with $\alpha_{i}=(\alpha_{i,1}\ldots,\alpha_{i,n})\in{\mathbb{Z}}^{n}_{\geq 0}$ represents a state in which there are $\alpha_{i,a}$ class $a$ particles at the $i$ th site. There is no constraint on the number of particles that occupy a site. The matrices $H_{+}$ and $H_{-}$ describe their stochastic hopping to the right and the left nearest neighbor sites, respectively. The transition rate can be read off the first terms on the RHS of (2.14) and (2.15), where the array $\gamma=(\gamma_{1},\ldots,\gamma_{n})$ specifies the numbers of particles that are jumping out. The superposition $H(a,b,\epsilon,q,\mu)$ corresponds to a mixture of such right and left moving dynamics. The rate is determined from the original occupancy ($\alpha$ for $h_{+}$ and $\beta$ for $h_{-}$) and the list of leaving particles ($\gamma$ for $h_{\pm}$) at the departure site and it is independent of the status of the destination site. Thus it defines a ZRP of $n$ classes of particles in a slightly generalized sense in that the rate is allowed to depend on $\gamma$. Here is a snapshot of the system for the $n=2$ case.

\ellipse8020\bullet$$\bullet$$\circ$$\circ\,\circ$$\bullet$$\circ$$\circ$$\circ$$\circ$$\bullet$$\bullet$$\circ$$\bullet$$\bullet$$\circ$$\cdot$$\cdot$$\cdot$$\cdot$$\cdot$$\cdot$$\circ\circ$$\circ\circ$$\circ$$\bullet\bullet$$\bullet$$\bullet$$\circ$$\circ$$\circ$$\bullet$$\bullet$$\circ$$\bullet$$\bullet$$\sigma_{L}$$\sigma_{1}$$\sigma_{i-1}$$\sigma_{i}$$\sigma_{i+1}$$\cdot$$\cdot$$\cdot$$\cdot$$\cdot$$h_{-,i-1,i}$$h_{+,i,i+1}$$\circ\,first class particle$\bullet\,$second class particle$\sigma_{i-1}=(0,1)$,$\sigma_{i}=(3,2)$,$\sigma_{i+1}=(1,2)$.

The local hopping in (2.14) and (2.15) with $n=2$ is depicted as

rate$\displaystyle{aw_{+}(\gamma|\alpha)}:$\drawline(0,5)(0,0) \drawline(0,0)(50,0)\drawline(50,0)(50,5) \drawline(25,0)(25,5) \drawline(12,12)(12,17)\drawline(12,17)(37,17) \overbrace{\circ...\circ}^{\alpha_{1}}\overbrace{\bullet...\bullet}^{\alpha_{2}}$$\circ...\circ\bullet...\bullet$$\overbrace{\circ...\circ}^{\gamma_{1}}\overbrace{\bullet...\bullet}^{\gamma_{2}}rate$\displaystyle{bw_{-}(\gamma|\beta)}:$\drawline(0,5)(0,0) \drawline(0,0)(50,0)\drawline(50,0)(50,5) \drawline(25,0)(25,5) \drawline(37,12)(37,17)\drawline(37,17)(12,17) \overbrace{\circ...\circ}^{\beta_{1}}\overbrace{\bullet...\bullet}^{\beta_{2}}$$\circ...\circ\bullet...\bullet$$\overbrace{\circ...\circ}^{\gamma_{1}}\overbrace{\bullet...\bullet}^{\gamma_{2}}

where the rate $w_{+}(\gamma|\alpha)$ and $w_{-}(\gamma|\beta)$ with $\epsilon=+1$ are given for $|\gamma|>0$ by

[TABLE]

The integrable Markov processes recalled here cover several models studied earlier. When $\epsilon=1,\mu\rightarrow 0$ in $H_{+}$, the multiple jumps $|\gamma|>1$ are suppressed in (2.14). So if $\gamma_{a}=1$ and the other components of $\gamma$ are 0, the rate reduces to $q^{\alpha_{1}+\cdots+\alpha_{a-1}}\frac{1-q^{\alpha_{a}}}{1-q}$. This reproduces the $n$-species $q$-boson process in [21] whose $n=1$ case further goes back to [18]. For $n=1$, there are numerous works including [3, 2, 16, 20] for example. One can overview their interrelation in [10, Fig.1,2]. When $\epsilon=1,(\mu,q)\rightarrow(0,0)$ in $H_{-}$, a kinematic constraint $\varphi(\gamma,\beta-\gamma)=\sum_{1\leq i<j\leq n}\gamma_{i}(\beta_{j}-\gamma_{j})=0$ occurs in (2.15). In order that $\gamma_{a}>0$ happens, the equalities $\gamma_{a+1}=\beta_{a+1},\gamma_{a+2}=\beta_{a+2},\ldots,\gamma_{n}=\beta_{n}$ must hold. It means that larger class particles have the priority to jump out, which precisely reproduces the $n$ class totally asymmetric zero range process explored in [12] after reversing the labeling of the classes $1,2,\ldots,n$ of the particles.

3. Stationary states

3.1. Definition and example

By definition a stationary state of the discrete time $U_{q}(A^{(1)}_{n})$ ZRP (2.11) is a vector $|\overline{P}\rangle\in W^{\otimes L}$ such that

[TABLE]

The stationary state is unique within each sector $m$ up to normalization. Apart from $m$, it depends on $q$ and the inhomogeneity parameters $\mu_{1},\ldots,\mu_{L}$ but not on $\lambda$ thanks to the commutativity (2.10). Sectors $m=(m_{1},\ldots,m_{n})$ such that $\forall m_{a}\geq 1$ are called basic. Non-basic sectors are equivalent to a basic sector of some $n^{\prime}<n$ models with a suitable relabeling of the classes. Henceforth we concentrate on the basic sectors. The coefficient appearing in the expansion

[TABLE]

is the stationary probability if it is normalized as $\sum_{(\sigma_{1},\ldots,\sigma_{L})\in S(m)}{\mathbb{P}}(\sigma_{1},\ldots,\sigma_{L})=1$. In this paper we will abuse the terminology also for the unnormalized states and probabilities.

The stationary states for the continuous Markov process (2.12) are those $|\overline{P}\rangle$ that satisfy $H|\overline{P}\rangle=0$. They are obtained from the discrete time ones just by the specialization to the homogeneous case $\mu_{1}=\cdots=\mu_{L}=\mu$. This is because (2.12) is an infinitesimal version of the commuting time evolutions (2.11) by the construction. In particular the stationary states are independent of $a$ and $b$ in (2.12).

Example 1**.**

Set $(n,L)=(2,3)$ and consider the homogeneous case $\mu_{1}=\mu_{2}=\mu_{3}=-\kappa$. The stationary state in the sector $m=(1,2)$ is given by

[TABLE]

where $|\emptyset,2,12\rangle$ for example is the multiset representation of the base vector $|(0,0)\rangle\otimes|(0,1)\rangle\otimes|(1,1)\rangle$ in the multiplicity representation. The terms “cyclic” are those obtained by the shift $|\sigma_{1},\sigma_{2},\sigma_{3}\rangle\rightarrow|\sigma_{i+1},\sigma_{i+2},\sigma_{i+3}\rangle$ for $i\in{\mathbb{Z}}_{3}\setminus\{0\}$. Similarly the stationary state in the sector $m=(2,1)$ is given by

[TABLE]

It has been conjectured [14, Ex.4] that for any sector $m$ there is a normalization such that $\mathbb{P}(\sigma_{1},\ldots,\sigma_{L})\in{\mathbb{Z}}_{\geq 0}[q,-\mu_{1},\ldots,-\mu_{L}]$ for all $(\sigma_{1},\ldots,\sigma_{L})\in S(m)$.

3.2. Matrix product formula

In [14] a matrix product formula for the stationary probability was obtained for general $n$ and inhomogeneity $\mu_{1},\ldots,\mu_{L}$. In the rest of the paper we shall exclusively deal with the $n=2$ case [13] with the homogeneous choice $\mu_{1}=\cdots=\mu_{L}=\mu$ in the regime $0<q,\mu<1$. In the continuous time setting, it corresponds to $\epsilon=+1$ in (2.13)–(2.15).

To describe the result we need the $q$-boson algebra $\mathcal{B}$ generated by ${\bf b},{\bf c},{\bf k}$ obeying the relations

[TABLE]

Let $F=\bigoplus_{m\geq 0}{\mathbb{C}}(q)|m\rangle$ be the Fock space and $F^{\ast}=\bigoplus_{m\geq 0}{\mathbb{C}}(q)\langle m|$ be its dual on which the $q$-boson operators ${\bf b},{\bf c},{\bf k}$ act as

[TABLE]

where $|\!-\!1\rangle=\langle-1|=0$. They satisfy the defining relations (3.3). The bilinear pairing of $F^{\ast}$ and $F$ is specified as $\langle m|m^{\prime}\rangle=\delta_{m,m^{\prime}}(q)_{m}$. Then $\langle m|(Q|m^{\prime}\rangle)=(\langle m|Q)|m^{\prime}\rangle$ is valid and the trace is given by $\mathrm{Tr}\,Q=\sum_{m\geq 0}\frac{\langle m|Q|m\rangle}{(q)_{m}}$. As a vector space, the $q$-boson algebra $\mathcal{B}$ has the direct sum decomposition

[TABLE]

The trace $\mathrm{Tr}\,Q$ is finite and nonzero only if $Q\in\bigoplus_{r\geq 1}\mathcal{B}^{r}_{0}$ when it is evaluated by $\mathrm{Tr}({\bf k}^{r})=(1-q^{r})^{-1}$. For $Q\in\mathcal{B}$, we let $Q^{\prime}$ denote the projection of $Q$ onto $\bigoplus_{r\geq 1}\mathcal{B}^{r}_{0}$.

We introduce the operator depending on $\mu,q$ which acts on (a completion of) $F^{\ast}$ and $F$:

[TABLE]

We have written $X_{(\alpha_{1},\alpha_{2})}$ as $X_{\alpha_{1},\alpha_{2}}$ for simplicity. The ratio of the infinite product of operators here is to be understood via the series expansion

[TABLE]

Theorem 2**.**

For any sector $m=(m_{1},m_{2})$ with $m_{2}\geq 1$, the stationary probability in (3.1) of the continuous time $U_{q}(A^{(1)}_{2})$ ZRP (2.12) or the discrete time $U_{q}(A^{(1)}_{2})$ ZRP (2.11) with homogeneous parameters $\mu_{1}=\cdots=\mu_{L}=\mu$ is expressed in the matrix product form

[TABLE]

This result was obtained in [13, eq.(42)]. The formula (3.6) depends on $\alpha_{1}$ and $\alpha_{2}$ quite differently, indicating distinct features between the two classes of particles. Interestingly this is a reflection of a tiny asymmetry of the hopping rate (2.16) and (2.17) under the interchange $(\alpha_{1},\beta_{1},\gamma_{1})\leftrightarrow(\alpha_{2},\beta_{2},\gamma_{2})$.

Example 3**.**

Consider the first two terms in (3.2).

[TABLE]

Setting $\mu=-\kappa$ and substituting $\mathrm{Tr}({\bf b}^{2}{\bf k}{\bf c}^{2})=q^{-1}\mathrm{Tr}({\bf b}^{2}{\bf c}{\bf k}{\bf c})=q^{-2}\mathrm{Tr}({\bf k}{\bf b}^{2}{\bf c}^{2})=(1-q^{3})^{-1}$ and $\mathrm{Tr}({\bf b}{\bf c}{\bf b}{\bf k}{\bf c})=(1+q^{2})^{2}\frac{(q)_{1}(q)_{2}}{(q)_{4}}$, we reproduce the first two terms in (3.2) after removing the common factor $(q)_{1}(q)_{2}(q)_{3}/(1+\kappa)^{2}$.

4. Density and currents of second class particles:

Formulation of the problem

Now we come to the main theme of the paper, the stationary quantities in a grand canonical ensemble with respect to the second class particles in the infinite volume limit. The first class particles are kept finite and regarded as defects. In this section we give a general formulation of the problem only deferring the results and their derivation to subsequent sections. We fix the parameters $0<q,\mu<1$ and consider the $\epsilon=+1$ case of (2.13)–(2.15) for the continuous time ZRP (2.12).

Let us introduce the generating series of $X_{m,n}$ (3.6) with respect to the second class particles with fugacity $y$:

[TABLE]

The quantity $g_{m}$ introduced here will be used very frequently in the sequel.

We consider the conditional probability in the stationary states supposing that there are $d_{i}$ first class particles at site $i$ for $i=1,\ldots,s$, and no first class particle is present elsewhere. So they form a fixed cluster of defects whose spatial extension is $s$ and the total number of defect particles is $d_{1}+\cdots+d_{s}$. In terms of the site variable $\sigma_{i}=(\sigma_{i,1},\sigma_{i,2})$, the condition is expressed as $\sigma_{i,1}=\theta(1\leq i\leq s)d_{i}$, allowing $\sigma_{i,2}$ still to fluctuate everywhere. We assume that $d_{1}\geq 1$ and $d_{s}\geq 1$ to fix the location of the defect cluster but allow the choice $d_{i}=0$ for $1<i<s$. The basic quantity is the conditional probability that a given site $r$ contains exactly $n$ second class particles. In the customary notation for the probability $P({\rm A}|{\rm B})$ of the event A under the condition B, we set

[TABLE]

There are three distinguish regions in the infinite volume limit $L\rightarrow\infty$ as

[TABLE]

where $r\leq 0$ should be understood as the site $r+L\in{\mathbb{Z}}_{L}$ in the limit $L\rightarrow\infty$. They correspond to the inside, the right and the left side of the defect cluster, respectively. The quantity (4.2) will be denoted by $P_{\rm I}(r,n),P_{\rm II}(r,n)$ and $P_{\rm III}(r,n)$ accordingly.

From the matrix product formula (3.8) and the definition (4.1), we have

[TABLE]

where $(\cdots)^{\prime}$ is defined after (3.5). In view of the matrix product operators (3.6) and (4.1), the prime restricts the ensemble to those sectors containing at least one second class particle for which the traces become finite. The formulas (4.5)–(4.6) fix the relative weight of such sectors, thereby specify what is meant by the “grand canonical ensemble” with respect to the second class particles with fugacity $y$. Obviously for any $r$ the normalization $\sum_{n\geq 0}P(r,n)=1$ should be fulfilled in each region.

Once the probability $P(r,n)$ is obtained, one can evaluate various physical quantities. In this paper we investigate the expectation number and currents of the second class particles at site $r$. The former is defined by

[TABLE]

According to the regions, it will be denoted by $\rho_{\rm I}(r),\rho_{\rm II}(r)$ and $\rho_{\rm III}(r)$. We call them density for simplicity. As for the current, there are two components $J(r)_{+}$ and $J(r)_{-}$ associated with the local Markov matrices $h_{+}$ (2.14) and $h_{-}$ (2.15) with $\epsilon=+1$ respectively:

[TABLE]

They sum up the contributions from the $l$ hopping second class particles out of total $n+\theta(1\leq r\leq s)d_{r}$ occupants weighted by the probability $P(r,n)$ and the hopping rate $w_{\pm}$ in (2.16) and (2.17). From (2.12) and the picture before (2.16), the total current from the site $r$ to $r+1$ is obtained by the superposition

[TABLE]

Thus it suffices to investigate $J(r)_{+}$ and $J(r)_{-}$ separately. We will also write them as $J_{\rm I}(r)_{\pm},J_{\rm II}(r)_{\pm}$ and $J_{\rm III}(r)_{\pm}$ according to the regions.

5. Defect-free case

First we illustrate the analysis on the defect-free case $s=0$ as a warm-up. It corresponds to the single species model, and some of the contents are well known by earlier works, e.g. [2, 5, 16].

In the absence of defects, the system acquires the ${\mathbb{Z}}_{L}$ translational symmetry. As the result, the probabilities (4.5) and (4.6) are equal and independent of $r$. So we simply denote it by $P(n)$. It is calculated as

[TABLE]

where the prime is defined under (3.5) and we have set

[TABLE]

These quantities will be utilized frequently in the subsequent calculations. We will abbreviate $\eta_{m}(q^{k}y)$ to $\eta_{m}$ if and only if $k=0$.

We suppose $0<y<1$ whose consistency will be confirmed shortly. Then

[TABLE]

holds due to $0<q,\mu<1$. In Appendix A we show that the limit $\lim_{L\rightarrow\infty}$ and the infinite sum $\sum_{m\geq 0}$ in (5.1) may be interchanged, i.e.,

[TABLE]

Since this is just $\sum_{m\geq 0}q^{mn}\delta_{m,0}=1$, we obtain the probability

[TABLE]

The correct normalization $\sum_{n\geq 0}P(n)=1$ has been achieved by virtue of (3.7). This formally agrees with the probability [2, eq.(53)] upon identification of the parameter $\alpha$ there with the fugacity $y$ here.

Next we relate the density $\rho$ of the second class particles to the fugacity $y$. In the grand canonical ensemble under consideration, it is evaluated as

[TABLE]

As mentioned before (5.5), the “partition function” $\mathrm{Tr}(A_{0}^{L})^{\prime}$ here may be replaced with $\Lambda(y)^{L}$ as $L$ goes to infinity, leading to

[TABLE]

The $f(\zeta)$ is a version of the $q$-digamma function. It monotonously grows from [math] to $\infty$ as $\zeta$ changes from [math] to $1$ behaving as $f(\zeta)\simeq\frac{\zeta}{1-q}\,(\zeta\searrow 0)$ and $f(\zeta)\simeq\frac{1}{1-\zeta}\,(\zeta\nearrow 1)$. Thus the fugacity and density are related asymptotically as

[TABLE]

The difference $f(y)-f(\mu y)$ in (5.7) is similarly increasing monotonously from [math] to $\infty$ for $y\in(0,1)$. Thus $\rho\in(0,\infty)$ is in one-to-one correspondence with $y\in(0,1)$, which confirms the consistency of the assumption made before (5.4). It implies that there is no phase transition typically recognized as condensation (cf. [5]) in ZRPs. To summarize, (5.7) determines the relation $y=y(\rho)$ and $\rho=\rho(y)$.

The density $\rho(r)$ (4.7) is also independent of $r$. From (5.6) it is evaluated as

[TABLE]

in terms of the average density $\rho$ in (5.7) confirming the consistency.

Now we are ready to evaluate the currents of the second class particles (4.8). Since they are independent of the site $r$, we simply write it as

[TABLE]

For instance $J_{-}$ is computed from (2.17) and (5.6) as

[TABLE]

where the last equality is due to (3.7). A similar calculation for $J_{+}$ leads to current-density relation, a basic characteristic of the system, given as

[TABLE]

via $y=y(\rho)$ (5.7). Here the function $h(\zeta)$ is the derivative of the $q$-digamma function $f(\zeta)$ in (5.8):

[TABLE]

It obviously satisfies the relation

[TABLE]

The result (5.10) is restated as

[TABLE]

which is $\mu$-independent despite that $w_{-}((0,l)|(0,n))$ (2.17) and $P(n)$ (5.6) depend on $\mu$ individually. An expression similar to the total current $aJ_{+}-bJ_{-}$ (4.9) with $J_{\pm}$ given by (5.11) has also been obtained in [2, Sec.4.1].

The function $h(\zeta)$ is monotonously increasing and tends to $\infty$ as $\zeta$ approaches 1 from below. Consequently the both currents in (5.11) grow monotonously with the density $\rho$ via (5.7). Their leading asymptotic behavior is given by

[TABLE]

So $J_{+}$ converges to a finite value whereas $J_{-}$ diverges as the density $\rho$ gets large. Such large $\rho$ behavior is out of question in asymmetric simple exclusion processes where the well known relation $J=\mathrm{const}\,\rho(1-\rho)$ makes sense only for $0\leq\rho\leq 1$.

We close the section with the description on the limiting cases $q,\mu\rightarrow 0,1$. We shall only give the leading terms.

[TABLE]

6. Density and currents: Results in general case

Let us proceed to the general case in which $d_{1},\ldots,d_{s}$ defect (first class) particles are present at the sites $1,\ldots,s$. We present the final results on the conditional probability $P(r,n)$ (4.2), the local density $\rho(r)$ (4.7) and the currents $J(r)_{\pm}$ (4.8) in the regions I, II and III in Theorem 4, 7 and 9. Some of them look bit messy but the point is that they are always finite sums of appropriate building blocks, which allow accurate numerical evaluations. Another point is that they are expressed in the form that elucidates the difference from the defect-free case $s=0$. We continue to stay in the range $0<q,\mu,y<1$ and use the following quantities that have already appeared in the previous section:

$\rho$: average density of the second class particles in the entire system,

$y\in(0,1)$: fugacity of the second class particles determined from $\rho$ via (5.7),

$P(n)$: probability in the defect-free case (5.6),

$J_{\pm}$: currents in the defect-free case (5.11),

$h(\zeta)$: derivative of $q$-digamma function describing the currents (5.12),

$\eta_{j}$: quantity controlling the decay of correlations (5.3).

In addition to them the following functions, detailed in Section 7, will be the basic ingredients:

$\phi(l|m)$: constituent of $G_{m,l}(d_{1},\ldots,d_{s})$ (7.4),

$G_{m,l}(d_{1},\ldots,d_{s})$: building block incorporating the effect of defects (7.9).

The relation (5.7) between the fugacity $y$ and the average density $\rho$ of the second class particles was originally derived for the defect-free case. We will justify its use in the presence of defects in Proposition 10 and the comments following it.

6.1. Main results

First we consider the region inside the defect cluster.

Theorem 4**.**

In the region I $(1\leq r\leq s)$ the following formulas are valid:

[TABLE]

where the sums $\sum_{m}$ extend over $m\in{\mathbb{Z}}_{\geq 0}$ satisfying $d_{r}\leq m\leq d_{1}+\cdots+d_{r}$.

The proof will be given in Section 7.2. By means of (3.7) one finds that the total probability is expressed as $\sum_{n\geq 0}P_{\rm I}(r,n)=\sum_{m}G_{0,m}(d_{1},\ldots,d_{r})$. This indeed gives $1$ thanks to (7.14). In Section 7.5 we will show that (6.2) can also be expressed as

[TABLE]

where $K(0)=0$ and the sum over $m$ is taken in the same way as in Theorem 4. The difference structure $K(r)-K(r-1)$ in (6.5) matches the sum rule in Proposition 10.

Example 5**.**

At $r=1$ which is the left boundary of the defect cluster, Theorem 4 simplifies to

[TABLE]

due to (7.11). Thus we always have $\rho_{\rm I}(1)>\rho$ because of $d_{1}\geq 1$. The last result $J_{\rm I}(1)_{-}=h(y)$ may seem strange. It follows from the properties

[TABLE]

which are easily seen in (5.6) and (2.17) and the remark on (5.14).

Example 6**.**

At $r=2$ which is one step inside the defect cluster from the left, Theorem 4 takes the form

[TABLE]

where $G_{0,m}(d_{1},d_{2})=\phi(d_{1}+d_{2}-m|d_{1})$ has been used by (7.12).

Let us remark on the low and high density asymptotic behavior. They correspond to $y\searrow 0$ and $y\nearrow 1$, respectively. See (5.9). From the properties of the unperturbed density and currents in (5.9), (5.15) and (5.16), one can derive the following behavior by utilizing (7.17) and (7.19).

[TABLE]

The result (6.7) with $r=1$ agrees with the $y\rightarrow 0$ case of Example 5 implied by (5.9).

Next we consider the right region II.

Theorem 7**.**

In the region II $(r>s)$ the following formulas are valid:

[TABLE]

where the sum in (6.12) extends over $i,j,m\in{\mathbb{Z}}_{\geq 0}$ such that $0\leq i\leq j\leq m$ and $d_{s}\leq m\leq d_{1}+\cdots+d_{s}$. The sums in (6.13) and (6.14) are taken for $1\leq j\leq m$ and $d_{s}\leq m\leq d_{1}+\cdots+d_{s}$.

The proof will be given in Section 7.3 It is based on the choice $(d_{1},\ldots,d_{r})=(d_{1},\ldots,d_{s},\overset{r-s}{\overbrace{0,\ldots,0}})$ in Theorem 4. In this sense Theorem 4 covers Theorem 7 save the extraction of the dependence on the distance $r-s$ in (7.31).

Example 8**.**

Consider the case $s=1$. From (7.11) we know $G_{0,m}(d_{1})=\delta_{m,d_{1}}$. Thus (6.12), (6.13) and (6.14) for $d_{1}=1$ and $d_{1}=2$ read

[TABLE]

where $\eta_{1}=\frac{1-y}{1-\mu y}$, $\eta_{2}=\frac{(1-y)(1-qy)}{(1-\mu y)(1-q\mu y)}$ by (5.3).

In general the effect of defects disappears in the long distance, i.e.,

[TABLE]

These are derived from $\lim_{r\rightarrow\infty}\eta_{j}^{r-s}=\delta_{j,0}$ by (5.4) and (7.14). As seen in Example 8, deviation from the defect-free case comes in various mode proportional to $\eta_{j}^{r-s}$ having the correlation (or decay) length $-(\log\eta_{j})^{-1}\,(j\geq 1)$. Thus from (5.4) and (5.3) the influence of the defects reaches the distance of order \bigl{(}\log\frac{1-\mu y}{1-y}\bigr{)}^{-1}. Since $y\rightarrow 0$ as the density tends to [math], the correlation is longer for smaller density and $\mu$ closer to $1$.

By a calculation similar to the region I, one can show for $r>s$ the asymptotic behavior as follows.

[TABLE]

These results have the form that naturally extends (6.7) – (6.11) to $r>s$. In particular in the large $\rho$ limit, the region II feels the rightmost $d_{s}$ defects only.

Finally we consider the left region III.

Theorem 9**.**

In the region III $(r\leq 0)$ the following formulas are valid:

[TABLE]

The proof will be given in Section 7.4.

Let us introduce the quantity

[TABLE]

which represents the total excess of the number of the second class particles from the average value in the entire system. From (6.8) and (6.16) one can easily check $\lim_{\rho\rightarrow\infty}\Delta\rho_{\rm tot}=-(d_{1}+\cdots+d_{s})$. Our next result shows that this holds in general.

Proposition 10**.**

The total excess of the second class particles is given by

[TABLE]

The proof will be given in Section 7.5. Since $\Delta\rho_{\rm tot}$ is finite, the parameter $\rho$ indeed remains to be the average density of the second class particles in the infinite volume limit. It is curious that $\Delta\rho_{\rm tot}$ exactly cancels the total excess $+(d_{1}+\cdots+d_{s})$ coming from the first class particles. To find a simple explanation of this fact is an interesting open problem.

The results in Theorem 4 (resp. Theorem 9) are independent of $d_{r+1},\ldots,d_{s}$ (resp. $d_{1},\ldots,d_{s}$). They imply that the defects only influence their right in the large volume limit despite that the component $H_{-}$ in the Markov matrix (2.12) represents the left moving particles. We do not have an intuitive explanation of this fact. One should however remember that the both $H_{\pm}$ in (2.13) originate in the discrete time evolution $T(\lambda|\mu,\ldots,\mu)$ described by the right moving particles only as in (2.8). Moreover the derivative in (2.13) gives rise to the denominator $T(\lambda|\mu,\ldots,\mu)^{-1}$ which is $1$ for $H_{+}$ whereas it is the left cyclic shift for $H_{-}$ due to $\mathscr{S}(\mu,\mu)(|\alpha\rangle\otimes|\beta\rangle)=|\beta\rangle\otimes|\alpha\rangle$. Thus we could have had the right moving dynamics only, if not so neatly described, just by taking the derivative rather than the logarithmic derivative at $\lambda=\mu$. In other words, the curiosity is attributed to the commutativity $[H_{+},H_{-}]=0$, a basic consequence of the integrability of the model, which implies that the left and the right moving dynamics should possess the same stationary states on the ring of finite size.

6.2. Case of homogeneous defect

To grasp the results in Theorem 4 and 7 qualitatively, it is helpful to consider the homogeneous case $d_{1}=\cdots=d_{s}=d$ and the limits $\rho\rightarrow 0$ and $\rho\rightarrow\infty$. The behavior for general $\rho$ can then be inferred as something in between. In this setting we still have the integer parameters $d,s\in{\mathbb{Z}}_{\geq 1}$ and the continuous parameters $0<q,\mu<1$. Set

[TABLE]

They satisfy $D>d$ and $\nu>1$ obviously. The results on the local density (6.7), (6.8), (6.16) and (6.19) are summarized as

[TABLE]

Schematic plots of them look as Figure 4 below.

Actually in the limit $\rho\rightarrow 0$, either $\rho(s)/\rho\gtrless 1$ can happen according to $\nu q^{(s-1)d}\gtrless 1$. In the other limit $\rho\rightarrow\infty$, the density completely resumes the average value for $r\geq s+2$. For general $0<\rho<\infty$, the local density $\rho(r)$ gets larger than $\rho$ at the left boundary $r=1$ of the defect cluster forming a peak. It then decreases until the site $r=s+1$ forming a valley at the left boundary of the region II. Then in the region II, it recovers toward the average value $\rho$ exponentially as mentioned after (6.15). Such behavior and sensitivity to the defects are sharper in higher density case.

As for the local currents we are to compare $J(r)_{+}$ and $J(r+1)_{-}$ in view of (4.9). The results given in (6.9) – (6.11) and (6.17) –(6.19) are summarized as

[TABLE]

Thus we see that $J_{+}(r)$ and $J_{-}(r+1)$ have similar asymptotics aside from the overall normalization for large $\rho$. As mentioned under (5.16) about $J_{\pm}$, the former is convergent whereas the latter is divergent as $\rho$ tends to infinity. They do not exhibit a peak at $r=1$ like the density $\rho(r)$, but other behavior is more or less similar. In particular when $\rho$ is sufficiently large, the both currents get smaller inside the defect cluster than the defect-free case $J_{\pm}$.

6.3. Comparison with numerical

evaluation in canonical ensemble

One can evaluate the density numerically by the matrix product formula (3.8) in a fixed sector with system size $L$. It formally corresponds to the canonical ensemble average which will be denoted by $\rho^{\rm c}_{L}(r)$. In the actual calculation of $\rho^{\rm c}_{L}(r)$, we have truncated the operators (3.6) and (4.1) to the matrices acting on the finite dimensional subspace of $F$ of the form $\bigoplus_{m=0}^{d_{1}+\cdots+d_{s}+t}{\mathbb{R}}|m\rangle$ and checked the convergence has been achieved sufficiently already for $t=1,2$, etc. Figure 5 compares $\rho(r)$ and $\rho^{\rm c}_{L}(r)$ for $1\leq r\leq L$.

As remarked after (6.15), the correlation length of the system is smaller for larger average density $\rho$. Thus the agreement of $\rho(r)$ and $\rho^{\rm c}_{L}(r)$ is expected to be better for larger $\rho$, therefore is harder to observe in numerical calculations. Admittedly the agreement in Figure 5 is not quite excellent, but always exhibits the tendency to improve as $L$ gets large. In general fluctuations in the density is of order $L^{-\frac{1}{2}}$. So the grand canonical approach should coincide with the canonical one in the large volume limit for there is no symptom of phase transition in the range $0<q,\mu<1$. See the remarks following (5.9).

6.4. Density and current profiles

Here we present the result of numerical evaluation of the formulas in Theorem 4, 7 and 9 in a number of figures. We first consider the profile of local density $\rho(r)$ (6.2), (6.13) and (6.19) in the presence of various defects in Figure 6–8. In general the density can possibly break the monotonicity inside the defect cluster depending on the inhomogeneity of $d_{i}$’s.

Examples of the density profile for a finite range of $\rho$ are given in Figure 7. To display the defect region I in the center, the lattice coordinate $r$ has been shifted. A similar convention will be employed in the subsequent figures in this subsection except Figure 9.

For a fixed average density $\rho$, dependence of $\rho(r)$ on $q$ and $\mu$ are shown in Figure 8.

One sees that the inhomogeneous defects make the density profile in the region I irregular, but they still tend to produce a peak and a valley at their boundaries as observed in the homogeneous case in Figure 4.

Let us turn to the local currents $J(r)_{\pm}$ (4.8). Recall that they have been determined as $J(r)_{\pm}=J_{\rm I}(r)_{\pm}$ (6.3), (6.4) for $1\leq r\leq s$ and $J(r)_{\pm}=J_{\rm II}(r)_{\pm}$ (6.14) for $r>s$ and $J(r)_{\pm}=J_{\pm}$ (6.19) , (5.11) for $r\leq 0$. In view of (4.9) we plot $J_{+}(r)$ and $J_{-}(r+1)$.

In Figure 6 and 9, one observes that $\rho(r)$ reaches its bottom at $r=s+1$ whereas $J_{+}(r)$ and $J_{-}(r+1)$ do not.

7. Derivation of main results

In this section we derive the main results given in Theorem 4, 7, 9 and Proposition 10.

7.1. Preliminary

Our first step is to reduce the trace $\mathrm{Tr}(\cdots)^{\prime}$ over the Fock space to the “vacuum expectation value” $\langle 0|(\cdots)|0\rangle$ in the infinite volume limit.

Proposition 11**.**

The traces in the probabilities (4.4)–(4.6) are reduced to the following:

[TABLE]

This is a corollary of Lemma 18 in Appendix A.

Let us prepare the functions that will serve as building blocks to express (7.1) and (7.2). For $l,m\in{\mathbb{Z}}_{\geq 0}$ set

[TABLE]

where $\Phi_{q}$ was defined in (2.3). Although the dependence on $\mu$ and $y$ is suppressed in this notation, it deserves attention that the fugacity $y$ plays the role of a spectral parameter here. From (2.6) we know

[TABLE]

where the summand is nonzero only for $0\leq l\leq m$.

Lemma 12**.**

For any $m\in{\mathbb{Z}}_{\geq 0}$ the following equality is valid:

[TABLE]

Proof.

Explicitly it reads as

[TABLE]

where we have replaced $j$ by $m-j$. We prove (7.6) by induction on $m$. The case $m=0$ is obvious. In what follows we assume (7.6) is valid when $m$ is replaced by $0,1,\ldots,m-1$.

Upon multiplication by $(\mu y)_{m}$, the both sides of (7.6) become polynomials in $\mu$ of order $m-1$. Thus it suffices to check the equality at $m$ points, say $\mu=q^{-r}$ with $r=0,1,\ldots,m-1$. So we set $\mu=q^{-r}$ in (7.6) and further replace $y$ by $q^{r}y$ to simplify the formula slightly. The result reads

[TABLE]

where $A=\sum_{i=r}^{m-1}\frac{1}{1-yq^{i}}$ and the upper bound of the $j$ sum has been reduced from $m-1$ to $r$ owing to the factor $(q^{-r})_{j}$. The coefficient of $A$ is $1$ due to (7.5) with $(\mu,y)$ replaced by $(q^{-r},q^{r}y)$. Let us subtract $A$ from (7.7). The upper bound of the $k$ sum in the RHS becomes $r-1$. The second sum in the LHS restricts the upper bound of $j$ to $r-1$. Further applying

[TABLE]

in the process, we find that the result is equivalent to

[TABLE]

This coincides with (7.6) with $(m,\mu,y)$ replaced by $(r,q^{-m},q^{m}y)$. Since $r\leq m-1$, its validity is assured by the induction hypothesis. ∎

We note that the limit $y\rightarrow 1$ in Lemma 12 leads to the identity

[TABLE]

The following function plays the basic role in our working.

[TABLE]

where $l,m\in{\mathbb{Z}}_{\geq 0}$. The sum in (7.9) extends over $l_{1},\ldots,l_{s}\in{\mathbb{Z}}_{\geq 0}$ obeying the specified condition. Since $\phi(l|m)=0$ unless $0\leq l\leq m$, nonzero summands in (7.9) are only those satisfying $l_{1}+\cdots+l_{i}\leq m+d_{1}+\cdots d_{i-1}$ for $i=1,\ldots,s$. The definition (7.9) is depicted as

[TABLE]

if each vertex is interpreted as an element of the $n=1$ stochastic $R$ matrix according to (7.4), (2.2) and (2.4). According to this diagram, the $G_{m,l}(d_{1},\ldots,d_{s})$ may be regarded a sum of elements of a column monodromy matrix of the $U_{q}(A^{(1)}_{1})$ ZRP containing the fugacity $y$ as a spectral parameter via (7.4). Some special cases which will be of frequent use are

[TABLE]

where the sum in (7.12) extends over $l_{1},\ldots,l_{s-1}\in{\mathbb{Z}}_{\geq 0}$ under the specified condition. It is easy to see

[TABLE]

In our working there appear many sums involving $G_{0,m}(d_{1},\ldots,d_{s})$. They always range over those $m$’s that satisfy the non-vanishing condition implied by (7.13). The sum rule (7.14) is derived by successive use of (7.5). The relation (7.15) follows directly from the diagrammatic representation (7.10). The sum over $k$ in it is finite due to (7.13). In particular when $t=1$ it implies $G_{0,l}(d_{1},\ldots,d_{s})=G_{d_{1},l}(d_{2},\ldots,d_{s})$ due to (7.11). The relation (7.16) is derived by applying $\phi(\gamma|\beta)|_{y\rightarrow q^{i}y}=\phi(\gamma|\beta+i)\frac{(q)_{\beta-\gamma+i}(q)_{\beta}}{(q)_{\beta+i}(q)_{\beta-\gamma}}q^{i\gamma}\eta_{i}^{-1}$ to the representation (7.10). Combining (7.18) with (7.8), we get

[TABLE]

This relation will be utilized to extract the large $\rho$ limits (6.8), (6.11), (6.16) and (6.18). We note that neither $G_{m,l}(0,d_{1},\ldots,d_{s-1})$ nor $G_{m,l}(d_{1},\ldots,d_{s-1},0)$ are equal to $G_{m,l}(d_{1},\ldots,d_{s-1})$. The function $G_{m,l}(d_{1},\ldots,d_{s})$ originates in the following quantity representing the effect of defects.

Lemma 13**.**

[TABLE]

Proof.

Substitute the expansion $A_{d_{i}}=g_{d_{i}}\sum_{l_{i}\geq 0}g_{l_{i}}{\bf b}^{l_{i}}\frac{(q^{d_{i}}\mu y{\bf k})_{\infty}}{(y{\bf k})_{\infty}}{\bf c}^{d_{i}}$ of (4.1) into the LHS. By sending all the factors $\frac{(q^{d_{i}}\mu y{\bf k})_{\infty}}{(y{\bf k})_{\infty}}$ in the bracket to the left by (3.3), we find that $\frac{\langle m|A_{d_{1}}\cdots A_{d_{s}}|l\rangle}{g_{d_{1}}\cdots g_{d_{s}}}$ is expressed as

[TABLE]

where the sum is over $l_{1},\ldots,\l_{s}\in{\mathbb{Z}}_{\geq 0}$ and the constraint is imposed to pick the non-vanishing bracket. From $\frac{(q^{M+d_{i}}\mu y)_{\infty}}{(q^{M}y)_{\infty}}=\Lambda(y)\frac{(y)_{M}}{(\mu y)_{M+d_{i}}}$ and

[TABLE]

the claimed formula follows. ∎

Lemma 14**.**

[TABLE]

Proof.

First we show the $s=0$ case, i.e.,

[TABLE]

by induction on $m$. At $m=0$, it is valid since $\langle 0|A_{0}^{L}|0\rangle=\Lambda(y)^{L}$ for any finite $L$. Assume that it is valid for $0,1,\ldots,m-1$ and let $U_{m}$ denote the LHS of (7.21) to be determined. Then we have

[TABLE]

Solving this taking (7.5) into account we find $U_{m}=y^{-m}(\mu y)_{m}$ establishing (7.21). Now we prove (7.20) for general $s$. By inserting $1=\sum_{l\geq 0}\frac{|l\rangle\langle l|}{(q)_{l}}$ into the LHS, it becomes

[TABLE]

Due to Lemma 13 and (7.21), this is equal to $y^{-m-d_{1}-\cdots-d_{s}}g_{d_{1}}\cdots g_{d_{s}}(\mu y)_{m}\sum_{l\geq 0}G_{m,l}(d_{1},\ldots,d_{s})$. Thus the proof is completed by (7.14). ∎

Remark 15**.**

By Lemma 18, the result (7.20) with $m=0$ is equivalent to

[TABLE]

The factor $g_{d_{1}}\cdots g_{d_{s}}$ here is proportional to the stationary probability of the configuration $(d_{1},\ldots,d_{s})$ of the single species model [16] of the first class particles on the length $s$ ring. In this sense the effect of the second class particles in the grand canonical picture is “renormalized” into the extra factor $y^{-1}$ for the first class particles’ fugacity, reducing the $U_{q}(A^{(1)}_{2})$ ZRP effectively to the $U_{q}(A^{(1)}_{1})$ ZRP. We expect a similar recursive feature in the $U_{q}(A^{(1)}_{n})$ ZRP with general $n$ [11], which may be viewed as a reminiscent of the nested Bethe ansatz. The result (7.22) is valid including $d_{i}=0$. It implies that separated clusters of the first class particles do not attract nor repel depending on their distance under the grand canonical background of the second class particles. They feel each other only when they merge together at the same site or split from a common site.

The following result reduces the multiple sum in $G_{m,l}(d_{1},\ldots,d_{s})$ (7.9) into a single sum when $d_{1}=\cdots=d_{s}=0$ and elucidates the $s$-dependence explicitly. It will be utilized in the region II to extract the large distance behavior of the density and currents away from the defects.

Lemma 16**.**

For $0\leq i\leq m$, the following formula is valid:

[TABLE]

Proof.

Introduce the function

[TABLE]

where the sum is over $l_{1},\ldots,l_{r}\in{\mathbb{Z}}_{\geq 0}$ under the specified condition. It is related to the LHS by

[TABLE]

where $F_{m-i,r}(q^{i}y)$ is given by (7.23) with $\eta_{j}=\eta_{j}(y)$ (5.3) replaced by $\eta_{j}(q^{i}y)$. The relation (7.24) is easily checked by means of (7.16) and $\eta_{i+j}=\eta_{i+j}(y)=\eta_{i}\eta_{j}(q^{i}y)$.

Replacing $y,m$ and the summation variable $j$ with $q^{-i}y,m+i,j+i$ and using

[TABLE]

with $y=q$, the formula to show becomes

[TABLE]

We prove this by induction on $r$. When $r=1$, the RHS is expressed in terms of the $q$-hypergeometric [7, eq.(1.2.14)] as $y^{-m}\frac{(\mu y)_{m}}{(q)_{m}}{}_{2}\phi_{1}\left({q^{-m},y\atop\mu y};q,q\right)$. Due to the formula [7, p354(II.6)]

[TABLE]

we obtain $\frac{(\mu)_{m}}{(q)_{m}}$, which agrees with (7.23) with $r=1$. We now prove the $r+1$ case assuming the $r$ case holds. From (7.23) we have

[TABLE]

Using the induction hypothesis and (7.25), exchanging the orders of the summations with respect to $l$ and $j$, it can be written as

[TABLE]

Using (7.27) and $(q^{1-m}y^{-1})_{m}=(-y)^{-m}q^{-m(m-1)/2}(y)_{m}$, we arrive at the expression (7.26) with $r$ replaced by $r+1$, which finishes the proof. ∎

Lemma 16 and (7.24) with $i=0$ yield the formula

[TABLE]

We remark that by the definition (7.9) and (7.4) the LHS in Lemma 16 is regular at $q=0$ and manifestly positive in the range $0<q,\mu,y<1$ under consideration, whereas such features are highly nontrivial in the RHS. In particular, individual terms therein can become $O(q^{-\frac{1}{2}m^{2}})$ which needs care in numerical evaluation for small $q$.

7.2. Proof of Theorem 4

Theorem 4 is concerned with the region I (4.3). Let us derive (6.1) from (7.1). Substituting (3.6) into it and applying (7.20) to the denominator we get

[TABLE]

Insert $1=\sum_{i\geq 0}\frac{|i\rangle\langle i|}{(q)_{i}}$ at the two places so as to use

[TABLE]

By Lemma 13 and (7.20) the result reads

[TABLE]

where $\phi(l-m+d_{r}|l)$ is defined in (7.4). Thus (6.1) follows by taking the $l$ sum using (7.15) with $(s,t)=(r,r-1)$ and (7.11).

Next we show (6.2) for the local density. After substituting (6.1) into (4.7), the sum is taken by

[TABLE]

where $\rho$ is the average density (5.7). The resulting expression for $\rho_{\rm I}(r)$ agrees with (6.2) except that $\rho$ comes with the coefficient $\sum_{m}G_{0,m}(d_{1},\ldots,d_{r})$. But this equals 1 by (7.14).

Let us proceed to the currents (6.3) and (6.4). From (4.8) they are expressed as

[TABLE]

Substitution of (2.16), (2.17) and (6.1) leads to the calculation similar to (5.10). Its essential part is

[TABLE]

where $h(y)$ is the derivative of the $q$-digamma function defined in (5.12). Now we have

[TABLE]

Expand $h(q^{m}\mu y)$ and $h(q^{m-d_{r}}y)$ into the sum of $h(y)$ and the remainder terms by (5.13). From (5.11) the contribution from the $h(y)$ is expressed as $\sum_{m}G_{0,m}(d_{1},\ldots,d_{r})J_{\pm}=J_{\pm}$ by (7.14). The remainder terms give the RHS of (6.3) and (6.4).

7.3. Proof of Theorem 7

Theorem 7 is concerned with the region II (4.3). The conditional probability $P_{\rm II}(r,n)$ for $r>s$ with the defects $(d_{1},\ldots,d_{s})$ is deduced from the region I result $P_{\rm I}(r,n)$ (6.1) by setting $(d_{1},\ldots,d_{r})=(d_{1},\ldots,d_{s},\overset{r-s}{\overbrace{0,\ldots,0}})$. It yields the expression

[TABLE]

where $P(n)$ is the probability in the defect-free case (5.6) and the sum extends over $0\leq i\leq d_{1}+\cdots+d_{s}$. On the other hand the decomposition (7.15) and Lemma 16 lead to

[TABLE]

Substituting $\binom{m}{j}_{q}=(-1)^{j}q^{mj-j(j-1)/2}\frac{(q^{-m})_{j}}{(q)_{j}}$ into this, we get (6.12). The normalization $\sum_{n\geq 0}P_{\rm II}(r,n)=1$ is attributed to that of $P_{\rm I}(r,n)$. In fact it can also be directly verified from (7.30) by noting $\sum_{n\geq 0}P(n)q^{ni}=\eta_{i}$ and (7.14).

To derive the density (6.13), we again set $(d_{1},\ldots,d_{r})=(d_{1},\ldots,d_{s},\overset{r-s}{\overbrace{0,\ldots,0}})$ in the region I result (6.2) and use (7.31) to get

[TABLE]

where the sum extends over $i,j,m\in{\mathbb{Z}}_{\geq 0}$ such that $0\leq i\leq j\leq m$ and $d_{s}\leq m\leq d_{1}+\cdots+d_{s}$. We may restrict $j$ to $1\leq j\leq m$. Then the sums over $i$ and $k$ are taken by applying the identity

[TABLE]

where the $q$-binomial theorem and (3.7) are used. After some simplification the result becomes (6.13).

To derive the current (6.14), we set $(d_{1},\ldots,d_{r})=(d_{1},\ldots,d_{s},\overset{r-s}{\overbrace{0,\ldots,0}})$ in the region I result (6.3), (6.4) and use (7.31) to get

[TABLE]

where the sum $\sum_{i,j,m}$ is taken in the same way as (7.32). This time we apply the identity

[TABLE]

which follows from (7.33) by differentiation. The result yields (6.14).

7.4. Proof of Theorem 9

Theorem 9 is concerned with the region III (4.3). From the definition of $X_{0,n}$ in (3.6) and $\langle 0|\frac{(\mu{\bf b})_{\infty}}{({\bf b})_{\infty}}=\langle 0|{\bf k}^{n}=\langle 0|$ by (3.4), the conditional probability (7.3) in the region III becomes

[TABLE]

due to (7.20). This is independent of $r$ and coincides with the probability $P(n)$ (5.6) for the defect-free case. Consequently the local density and current also reduce to the average density $\rho$ and $J_{\pm}$ in (5.11).

7.5. Proof of Proposition 10

From (6.13), (5.3) and (5.4), it follows that

[TABLE]

where the sum $\sum_{j,m}$ ranges over $m$ such that $G_{0,m}\neq 0$ implied by (7.13) and $1\leq j\leq m$. On the other hand we have

[TABLE]

by replacing $(q,d,\mu)$ in (7.8) with $(q^{-1},m,q^{m-1}\mu y)$. Thus the above expression is simplified to

[TABLE]

Now we introduce

[TABLE]

From (6.2) it is easy to see

[TABLE]

Thus Proposition 10 follows from

[TABLE]

The first equality is obvious from (7.11) and the last one has already been derived in (7.34). So we are left to show $K(r)={\tilde{K}}(r+1)$ only. By using (7.15) with $(t,s)=(r,r+1)$ and (7.11) we find

[TABLE]

which tells that ${\tilde{K}}(r+1)$ is actually independent of $d_{r+1}$ as with $K(r)$. Now the proof is finished by Lemma 12.

8. Discussion

In this paper we have investigated the stationary properties of the $U_{q}(A^{(1)}_{2})$ ZRP [11] on a ring whose stationary probability is not factorized but expressed in a nontrivial matrix product form. The second class particles are treated in the grand canonical ensemble subject to the influence of the $d_{1},\ldots,d_{s}$ first class particles fixed as defects at the sites $1,\ldots,s$. The local density and current of the second class particles in the infinite volume limit are obtained in Theorem 4, 7 and 9. The density profile is shown to exhibit a peak and a valley at the boundaries of the defect cluster. The currents are locally suppressed by the defects. The detail is dependent on the inhomogeneity of the data $(d_{1},\ldots,d_{s})$. Outside the defect cluster, its influence reaches longer distance when the average density of the second class particles is lower. It was shown that the second class particles are decreased by the number of the defect particles in Proposition 10.

The effect of the defects is expressed through the function $G_{l,m}(d_{1},\ldots,d_{s})$, which is essentially a column monodromy matrix of the $U_{q}(A^{(1)}_{1})$ ZRP containing the fugacity $y$ as the spectral parameter. Compare the diagrams (7.10) and (2.8). We regard it as another reminiscent of the nested Bethe ansatz structure in addition to Remark 15.

It is natural to generalize the analysis in this paper to the grand canonical ensemble for the both classes of particles simultaneously. Denoting the fugacity of the first class particles by $x$, one is led to the analogue of the grand canonical partition function in two variables as

[TABLE]

where we have used the cyclicity of the trace and the commutativity $[\Gamma(x,y),\Gamma(\mu x,\mu y)]=0$ which follows from (3.6) by applying $(\mu)_{m+n}=(\mu)_{m}(q^{m}\mu)_{n}=(\mu)_{n}(q^{n}\mu)_{m}$. See also [14, Rem.3]. We leave the large $L$ asymptotics of $Z_{L}(x,y)$ and related issues for a future study.

Appendix A Proof of (5.5) and Proposition 11

We invoke the interchangeability of limits and infinite sums in elementary calculus adapted in the following form.

Lemma 17**.**

For the sequence $\{S_{L,l}\in{\mathbb{R}}_{\geq 0}\mid L,l\in{\mathbb{Z}}_{\geq 0}\}$, the interchangeability $\lim_{L\rightarrow\infty}\sum_{l=0}^{\infty}S_{L,l}=\sum_{l=0}^{\infty}\lim_{L\rightarrow\infty}S_{L,l}$ holds if the following (i) or (ii) is satisfied for $L$ sufficiently large:

1. (i)

$\lim_{L\rightarrow\infty}S_{L,l}<\infty$* and $S_{L,l}\leq T_{l}$ for some $T_{l}$ such that $\sum_{l=0}^{\infty}T_{l}<\infty$.* 2. (ii)

$S_{L+1,l}\leq S_{L,l}$* and $\sum_{l=0}^{\infty}S_{L,l}<B$ for some $B$ independent of $L$.*

Proof.

(i) Set $s_{l}=\lim_{L\rightarrow\infty}S_{L,l}$ By the definition $s_{l}\leq T_{l}$. For any $\varepsilon>0$ there is $M_{\varepsilon}$ such that $\sum_{l>M_{\varepsilon}}T_{l}<\varepsilon$, and there is also $L_{\varepsilon}$ such that $|S_{L,l}-s_{l}|<\frac{\varepsilon}{1+M_{\varepsilon}}\,(0\leq l\leq M_{\varepsilon})$ for $L>L_{\varepsilon}$. Then one has

[TABLE]

(ii) Without loss of generality we assume the condition holds for $L\geq 0$. From the assumption $\lim_{L\rightarrow\infty}S_{L,l}$ exists. Set $a_{i,l}=S_{i,l}-S_{i-1,l}\,(i\geq 1)$ and $a_{0,l}=S_{0,l}$. Then $|a_{i,l}|=(-1)^{\theta(i\geq 1)}a_{i,l}$. Since $\sum_{i=0}^{L}\sum_{l=0}^{\infty}|a_{i,l}|=2\sum_{l=0}^{\infty}S_{0,l}-\sum_{l=0}^{\infty}S_{L,l}<2B$, The sum of $a_{i,l}$ over $(i,l)\in{\mathbb{Z}}_{\geq 0}^{2}$ is absolutely convergent and can be taken in any order. Thus we get $\lim_{L\rightarrow\infty}\sum_{i=0}^{L}\sum_{l=0}^{\infty}a_{i,l}=\lim_{m\rightarrow\infty}\sum_{l=0}^{m}\sum_{i=0}^{\infty}a_{i,l}$. This means $\lim_{L\rightarrow\infty}\sum_{l=0}^{\infty}S_{L,l}=\sum_{l=0}^{\infty}\lim_{L\rightarrow\infty}S_{L,l}$. ∎

Proof of (5.5). For $n\geq 1$, set $S_{L,l}=q^{ln}\eta^{L}_{l}$ and $T_{l}=q^{ln}$. Then the condition (i) in Lemma 17 is satisfied. For $n=0$, set $S_{L,l}=\eta^{L}_{l}-\eta^{L}_{\infty}$ which satisfies $S_{L+1,l}\leq S_{L,l}$ for sufficiently large $L$. Set further $B_{0}=\sum_{l=0}^{\infty}S_{1,l}=\eta_{\infty}\sum_{j\geq 1}\frac{g_{j}y^{j}}{1-q^{j}}$ which is finite. From (5.4), we have $\sum_{l=0}^{\infty}S_{L,l}=1-\eta_{\infty}^{L}+\sum_{l=1}^{\infty}(\eta_{l}-\eta_{\infty})(\eta^{L-1}_{l}+\eta^{L-2}_{l}\eta_{\infty}+\cdots+\eta^{L-1}_{\infty})<1+B_{0}L\eta_{1}^{L-1}$. Thus the condition (ii) in Lemma 17 is satisfied by taking $B=1+B_{0}\max\{L\eta_{1}^{L-1}\mid L\in{\mathbb{Z}}_{\geq 1}\}<\infty$. ∎

Let $Q$ be the operators of the form

[TABLE]

for some $m_{1},\ldots,m_{\kappa},e,h\geq 0$ and $1\leq t\leq\kappa$. They are contained in the matrix product formula of the probabilities (4.4)–(4.6) in the form $\mathrm{Tr}(QA_{0}^{L})^{\prime}$ with $L$ large, where $(\cdots)^{\prime}$ is explained after (3.5).

Lemma 18**.**

For the operator $Q$ of the form (A.1), the following limits exist and are equal:

[TABLE]

Proof.

As for the RHS we have

[TABLE]

where $l_{1},\ldots,l_{L}$ are summed over ${\mathbb{Z}}_{\geq 0}$ with the specified condition. We have used the expansion (3.7), the quantity $\eta_{m}$ (5.3) and the definition of the function $F_{m,L}(y)$ (7.23). By the assumption on $Q$, $\langle 0|Q|m\rangle\neq 0$ holds only for finitely many $m\in{\mathbb{Z}}_{\geq 0}$. Thus we get

[TABLE]

where the second equality is due to (7.24) and Lemma 16 with $(i,r)=(0,L)$. This result shows that the RHS is finite.

As for the LHS we expand $Q$ by $\frac{(\mu{\bf b})_{\infty}}{({\bf b})_{\infty}}=\sum_{j\geq 0}g_{j}{\bf b}^{j}$ (3.7) and apply the commutation relation (3.3). The resulting terms contain ${\bf k}$ as an overall factor and are grouped into two types as

[TABLE]

Here $Q_{0}$ and $Q_{1}$ are finite linear combinations of the operators of the form

[TABLE]

with various $\alpha_{i},\beta_{i}\in{\mathbb{Z}},\,\lambda\in{\mathbb{Z}}_{\geq 0},\,n\in{\mathbb{Z}}_{\geq 1}$, and $\kappa$ is the one in (A.1). Let us pick one such term for each type:

[TABLE]

These are actually single sums over $l$ since $m$ is frozen to $m=l+\lambda$. Denote them by $W_{0}=\sum_{l\geq 0}U_{L,l}$ and $W_{1}=\sum_{l\geq 0}V_{L,l}$. A direct calculation gives

[TABLE]

We also introduce

[TABLE]

Since the $l$’s in $U_{L,l}$ and $V_{L,l}$ have the same meaning as those in (A.3), the proof is reduced to showing that the only $l=0$ term survives in the large $L$ limit.

[TABLE]

From (5.4) the equalities $\lim_{L\rightarrow\infty}u_{L,l}=\delta_{l,0}\lim_{L\rightarrow\infty}u_{L,0}$ and $\lim_{L\rightarrow\infty}v_{L,l}=\delta_{l,0}\lim_{L\rightarrow\infty}v_{L,0}=\delta_{l,0}\lim_{L\rightarrow\infty}\overline{v}_{L}$ are valid obviously. Therefore (A.5) follows from the next Lemma 19. ∎

Lemma 19**.**

The following interchangeability holds for the series involving $U_{L,l},V_{L,l}$ in (A.4):

[TABLE]

Proof.

(a). From (7.28) and (5.4) we have

[TABLE]

Thus by setting $S_{L,l}=U_{L,l}$ and $T_{l}=y^{-\lambda}q^{nl}\Gamma_{l}\frac{(\mu y)_{\lambda}(-q^{1-\lambda})_{\lambda}}{(q)_{\lambda}}$, the condition (i) in Lemma 17 is satisfied, hence (a) is valid. In particular, $\sum_{l=0}^{\infty}T_{l}<\infty$ because the series $\sum_{l=0}^{\infty}w^{l}\Gamma_{l}$ is convergent for $|w|<1$ due to $\lim_{l\rightarrow\infty}\Gamma_{l}=1$.

(b) We prove that Lemma 17 (ii) applies to $S_{L,l}=V_{L,l}$. Thus we are to verify (b)${}_{1}:\,v_{L+1,l}\leq v_{L,l}$ for sufficiently large $L$, and (b)${}_{2}:\,\sum_{l\geq 0}V_{L,l}<B$ for some $B$. To show (b)1, note the large $L$ asymptotic behavior

[TABLE]

where the first and the second equalities follow from (7.23) and (7.28) with $\eta_{l+m}=\eta_{l}\eta_{m}(q^{l}y)$, and the last estimation is due to (5.4). As for the second term in $v_{L,l}$ containing $-\eta^{L}_{\infty}$, one can estimate the coefficient of it as $\binom{L-1+\lambda}{\lambda}c_{1}\leq\sum_{l_{1}+\cdots+l_{L}=\lambda}g_{l_{1}}\cdots g_{l_{L}}\leq\binom{L}{\lambda}c_{2}$, where $c_{1},c_{2}$ are $L$-independent quantities $c_{1}=\min\{g_{l_{1}}\cdots g_{l_{\lambda}}\mid l_{i}\in{\mathbb{Z}}_{\geq 0},l_{1}+\cdots+l_{\lambda}=\lambda\}$ and $c_{2}=\sum_{l_{1}+\cdots+l_{\lambda}=\lambda}g_{l_{1}}\cdots g_{l_{\lambda}}$. From these results and the Stirling formula one finds that $v_{L,l}$ has the large $L$ asymptotic behavior $v_{L,l}\simeq c_{3}\eta_{l}^{L}-c_{4}L^{\lambda}\eta_{\infty}^{L}$ with $L$-independent $c_{3},c_{4}>0$. Thus the inequality (b)1 for sufficiently large $L$ follows from (5.4). To show (b)2, we use

[TABLE]

From (A.4), there exists $l_{0}\geq 0$ such that $\Gamma_{l}>0$ holds for all $l\geq l_{0}$333 Such $l_{0}$ is not unique but the non-uniqueness does not spoil the subsequent argument.. Then the above estimation leads to $\sum_{l\geq 0}V_{L,l}=C+\sum_{l\geq l_{0}}\Gamma_{l}v_{L,l}\leq C+L\eta_{\infty}\sum_{l_{1}+\cdots+l_{L}=\lambda}g_{l_{1}}\cdots g_{l_{L}}\sum_{j\geq 1}f_{j}g_{j}y^{j}$, where $C=\sum_{0\leq l<l_{0}}V_{L,l}$ is a finite constant and we have set $f_{j}=\sum_{l\geq l_{0}}\Gamma_{l}q^{jl}$. This series $f_{j}$ is convergent because of $\lim_{l\rightarrow\infty}\Gamma_{l}=1$ as noted in (a). The series $\sum_{j\geq 1}f_{j}g_{j}y^{j}$ in $y$ is also convergent due to $0<y<1$ and $\lim_{j\rightarrow\infty}\frac{f_{j}g_{j}}{f_{j+1}g_{j+1}}=q^{-l_{0}}\lim_{j\rightarrow\infty}\frac{1-q^{j+1}}{1-q^{j}\mu}=q^{-l_{0}}>1$. Thus $\sum_{l\geq 0}V_{L,l}$ is convergent for any fixed $L$. On the other hand (b)1 tells that $\sum_{l\geq 0}V_{L,l}$ is monotonously decreasing with respect to $L$ for sufficiently large $L$. This proves (b)2. ∎

Acknowledgments

The authors thank Masato Okado for discussion and Tomohiro Sasamoto for communication. This work is supported by Grants-in-Aid for Scientific Research No. 15K13429 from JSPS.