# On Truncation of the Matrix-Geometric Stationary Distributions

^{1}

^{2}

^{3}

^{*}

Next Article in Journal / Special Issue

Previous Article in Journal

Previous Article in Special Issue

Previous Article in Special Issue

Service Innovation Research Institute, 00120 Helsinki, Finland

Department of Applied Informatics and Probability, Peoples’ Friendship University of Russia (RUDN University), Miklukho-Maklaya St. 6, Moscow 117198, Russia

Institute of Informatics Problems, Federal Research Center “Computer Science and Control” of the Russian Academy of Sciences, Vavilov St. 44-2, Moscow 119333, Russia

Author to whom correspondence should be addressed.

Received: 30 July 2019 / Revised: 25 August 2019 / Accepted: 27 August 2019 / Published: 1 September 2019

(This article belongs to the Special Issue Stochastic Processes: Theory and Applications)

In this paper, we study queueing systems with an infinite and finite number of waiting places that can be modeled by a Quasi-Birth-and-Death process. We derive the conditions under which the stationary distribution for a loss system is a truncation of the stationary distribution of the Quasi-Birth-and-Death process and obtain the stationary distributions of both processes. We apply the obtained results to the analysis of a semi-open network in which a customer from an external queue replaces a customer leaving the system at the node from which the latter departed.

Consider a queueing system with an infinite number of waiting places which may be modeled by a Quasi-Birth-and-Death (QBD) process $X(t)$ on a state space ${\cup}_{l=0}^{\infty}{X}_{l}},$ where ${X}_{0}=\{(i,0),\text{\hspace{0.17em}}1\le i\le n+m\}$ and ${X}_{l}=\{(i,l),\text{\hspace{0.17em}}1\le i\le m\}$, $l\ge 1$. Here l represents the number of customers in the queue and $i$ represents the internal state of the system. The generator of the process $X(t)$ has the block tri-diagonal form
where ${{\rm N}}_{0}$ is a square matrix of order n, the matrices $\Lambda $, ${\rm N}$, ${\rm M}$ and ${{\rm N}}_{1}$ are square of order m, and the matrices ${\Lambda}_{0}$ and ${{\rm M}}_{1}$ are rectangular of appropriate dimensions. We denote by $x$ the stationary probability vector associated with $Q$ and we partition it as $x=(p,{q}_{0},{q}_{1},\dots )$, where $p$ is an n-vector, and ${q}_{l}$, $l\ge 0$, are m-vectors. We also define the generator $A=\Lambda +{\rm N}+{\rm M}$.

$$Q=\left[\begin{array}{cccccc}{{\rm N}}_{0}& {\Lambda}_{0}& 0& 0& 0& \cdots \\ {{\rm M}}_{1}& {{\rm N}}_{1}& \Lambda & 0& 0& \cdots \\ 0& {\rm M}& {\rm N}& \Lambda & 0& \cdots \\ 0& 0& {\rm M}& {\rm N}& \Lambda & \cdots \\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \end{array}\right]$$

According to [1], the positive recurrent Markov chain $X(t)$ has a modified matrix-geometric stationary probability vector,
where $R$, which is called the rate matrix, has the spectral radius $\rho \text{\hspace{0.17em}}<\text{\hspace{0.17em}}1$ and is the minimal nonnegative solution of the matrix equation

$${q}_{l}={q}_{1}{R}^{l},\text{}l\ge 0,$$

$${R}^{2}{\rm M}+R{\rm N}+\Lambda =0.$$

Once $R$ is known, vectors $p$ and ${q}_{0}$ are obtained by solving the linear system
with the normalizing condition
where $u$ is the column vector of all ones of appropriate size.

$$p{{\rm N}}_{0}+{q}_{0}{{\rm M}}_{1}=0$$

$$p{\Lambda}_{0}+{q}_{0}({{\rm N}}_{1}+R{\rm M})=0$$

$$pu+{q}_{0}{(I-R)}^{-1}u=1,$$

The rate matrix $R$ is usually computed via iterative algorithms [1,2,3,4] among which the fastest is the logarithmic reduction algorithm [2]. Available methods for computation of the rate matrix are effective only for relatively small values of m. Exceptions are processes with off-diagonal blocks of a special type. If matrices $\Lambda $ or ${\rm M}$ have rank one [5,6], or matrix $A$ is triangular [7,8], it is possible to find an explicit solution of Equation (3).

Now, consider a queueing system with $L\ge 0$ waiting places, that can be modeled by a Markov chain ${X}_{L}(t)$ on a state space $\cup}_{l=0}^{L}{X}_{l$ with a generator of the block tri-diagonal form

$${Q}_{L}=\left[\begin{array}{cccccc}{{\rm N}}_{0}& {\Lambda}_{0}& 0& 0& \cdots & 0\\ {{\rm M}}_{1}& {{\rm N}}_{1}& \Lambda & 0& \cdots & 0\\ 0& {\rm M}& {\rm N}& \ddots & \ddots & \vdots \\ 0& 0& \ddots & \ddots & \Lambda & 0\\ \vdots & \vdots & \ddots & {\rm M}& {\rm N}& \Lambda \\ 0& 0& \cdots & 0& {\rm M}& {\rm N}+\Lambda \end{array}\right].$$

The change in the lower diagonal block is due to the fact that if an arriving customer finds all waiting places occupied it is lost, but the internal state of the system changes in the same manner as if the customer were not lost.

We partition the stationary probability vector ${x}_{L}$ associated with ${Q}_{L}$ as ${x}_{L}=({p}_{L},{q}_{L,0},{q}_{L,1},\dots ,{q}_{L,L})$, where ${p}_{L}$ is an n-vector, and ${q}_{L,l}$, $0\text{\hspace{0.17em}}\le \text{\hspace{0.17em}}l\text{\hspace{0.17em}}\le \text{\hspace{0.17em}}L,$ are m-vectors. Solving finite QBD processes is more complicated than solving infinite processes [9,10,11]. If for some number s matrices $\Lambda $, ${\rm N}$, ${\rm M}$ satisfy $\mathrm{det}({s}^{2}{\rm M}+s{\rm N}+\Lambda )\ne 0$ then the stationary probability vector ${x}_{L}$ can be expressed using two matrix-geometric terms as [12]
which require computation of the minimal nonnegative solutions of two matrix equations

$${q}_{L,l}={a}_{L}{R}^{l}+{b}_{L}{S}^{L-l},\text{}0\le l\le L,$$

$${R}^{2}{\rm M}+R{\rm N}+\Lambda =0,\text{}{\rm M}+S{\rm N}+{S}^{2}\Lambda =0.$$

The vectors ${a}_{L}$, ${b}_{L}$ and ${p}_{L}$ can be found by solving the linear system
and normalizing vector ${x}_{L}$. The condition $\mathrm{det}({s}^{2}{\rm M}+s{\rm N}+\Lambda )\ne 0$ holds, for example, if the generator $A=\Lambda +{\rm N}+{\rm M}$ is irreducible and $\mathsf{\pi}\Lambda u\ne \mathsf{\pi}{\rm M}u$, where $\mathsf{\pi}$ is the stationary probability vector of $A$ [13].

$${p}_{L}{{\rm N}}_{0}+({a}_{L}+{b}_{L}{S}^{L}){{\rm M}}_{1}=0$$

$${p}_{L}{\Lambda}_{0}+({a}_{L}+{b}_{L}{S}^{L}){{\rm N}}_{1}+({a}_{L}R+{b}_{L}{S}^{L-1}){\rm M}=0$$

$$({a}_{L}{R}^{L-1}+{b}_{L}S)\Lambda +({a}_{L}{R}^{L}+{b}_{L})({\rm N}+\Lambda )=0$$

The stationary probability vector $x=(p,{q}_{0},{q}_{1},\dots )$ of the generator (Equation (1)) satisfies the equation
for all $l\ge 0$. In the finite case the stationary probability vector ${x}_{L}=({p}_{L},{q}_{L,0},{q}_{L,1},\dots ,{q}_{L,L})$ also satisfies the equation
for all $0\le l<L$, but the boundary equation is different

$${q}_{l+1}{\rm M}+{q}_{l}{\rm N}+{q}_{l-1}\Lambda =0$$

$${q}_{L,l+1}{\rm M}+{q}_{L,l}{\rm N}+{q}_{L,l-1}\Lambda =0$$

$${q}_{L,L}({\rm N}+\Lambda )+{q}_{L,L-1}\Lambda =0.$$

It is clear that if the rate matrix $R$ satisfies ${R}^{L+1}{\rm M}={R}^{L}\Lambda $ then the stationary distribution of ${X}_{L}(t)$ can be obtained by the truncation of the matrix-geometric stationary distribution (Equation (2)) of $X(t)$ at the level of $l=L$. In turn, this is enough to fulfill the condition ${R}^{l+1}{\rm M}={R}^{l}\Lambda $ for some $l\le L.$ In this paper we study the properties of QBD processes with rate matrix $R$ satisfying $R{\rm M}=\Lambda $. If it is so, then for any L the stationary distribution of the process ${X}_{L}(t)$ can be obtained by truncation of the stationary distribution of $X(t)$.

In the next section we show when the rate matrix $R$ satisfies $R{\rm M}=\Lambda $. In Section 3, we analyze the stationary distributions of infinite and finite QBD processes. In Section 4, we use the obtained results for solving a semi-open network in which customers from an external queue may only arrive to the node from which a departed customer left the system. We assume the process $X(t)$ to be positive recurrent.

The following notation conventions are used throughout the article. Bold lowercase letters denote vectors, and bold capitalized letters denote matrices, $\delta (i,j)=1$ if $i=j$ and $\delta (i,j)=0$ otherwise; $I=[\delta (i,j)]$ is the identity matrix; the ith component of vector ${e}_{i}$ is equal to one and the others are zero.

Here we derive the explicit solution of matrix Equation (3) obtained for a rate matrix satisfying $R{\rm M}=\Lambda $.

Let$A$be an irreducible generator with stationary probability vector$\mathsf{\pi}$, and let$R$be the minimal nonnegative solution of Equation (3). Then the matrix$R$satisfies the relation$R{\rm M}=\Lambda $if and only if
where $\mathsf{\lambda}=\Lambda u$ and $\mathsf{\mu}={\rm M}u$. Moreover, the matrix $R$ and its spectral radius $\eta $ are given by

$$\Lambda =\frac{1}{\mathsf{\pi}\mathsf{\mu}}\mathsf{\lambda}\mathsf{\pi}{\rm M},$$

$$R=\frac{1}{\mathsf{\pi}\mathsf{\mu}}\mathsf{\lambda}\mathsf{\pi},$$

$$\eta =\frac{\mathsf{\pi}\mathsf{\lambda}}{\mathsf{\pi}\mathsf{\mu}}.$$

Let matrix $R$ satisfy the relation $R{\rm M}=\Lambda $. Then, Equation (3) can be rewritten as

$$RA=R\Lambda +R{\rm N}+\Lambda ={R}^{2}{\rm M}+R{\rm N}+\Lambda =0.$$

Since the generator $A$ is irreducible, each column of the matrix $R$ is proportional to the vector $\mathsf{\pi}$. Therefore $R$ has the form $R=w\mathsf{\pi}$, where $w$ is a row vector satisfying

$$w\mathsf{\pi}{\rm M}=\Lambda .$$

By postmultiplying Equation (15) by vector $u$ we obtain $w\mathsf{\pi}\mathsf{\mu}=\mathsf{\lambda},$ which implies Equations (11) and (12).

Let matrix $\Lambda $ be given by Equation (11). Since the process $X(t)$ is positive recurrent, the matrix $N+RM$ is nonsingular [1] and we have

$$R=-\Lambda {(N+RM)}^{-1}.$$

This relation implies that the rate matrix is of the form $R=\mathsf{\lambda}\mathsf{\xi}$, where $\mathsf{\xi}$ is the unique solution of the equation

$$\mathsf{\xi}=-\frac{1}{\mathsf{\pi}\mathsf{\mu}}\mathsf{\pi}{\rm M}{({\rm N}+\mathsf{\xi}\mathsf{\lambda}{\rm M})}^{-1}.$$

It is easy to check that the vector
satisfies Equation (16). Therefore, the rate matrix and its maximum eigenvalue are given by Equations (12) and (13). This completes the proof of the theorem. □

$$\mathsf{\xi}=\frac{1}{\mathsf{\pi}\mathsf{\mu}}\mathsf{\pi}$$

Let matrices$\Lambda $,${\rm M}$,${\rm N}$have the block-diagonal form
and suppose that ${A}_{i}={\Lambda}_{i}+{{\rm N}}_{i}+{{\rm M}}_{i}$ is an irreducible generator with a stationary probability vector ${\mathsf{\pi}}_{i}$. Then the minimal nonnegative solution
of Equation (3) satisfies the relation $R{\rm M}=\Lambda $, if and only if
where ${\mathsf{\lambda}}_{i}={\Lambda}_{i}u$ and ${\mathsf{\mu}}_{i}={{\rm M}}_{i}u$. The matrices ${R}_{i}$ and their spectral radii ${\eta}_{i}$ are given by

$$\Lambda =\left[\begin{array}{cccc}{\Lambda}_{1}& 0& \cdots & 0\\ 0& {\Lambda}_{2}& \ddots & \vdots \\ \vdots & \ddots & \ddots & 0\\ 0& \cdots & 0& {\Lambda}_{k}\end{array}\right],\text{\hspace{0.17em}}{\rm N}=\left[\begin{array}{cccc}{{\rm N}}_{1}& 0& \cdots & 0\\ 0& {{\rm N}}_{2}& \ddots & \vdots \\ \vdots & \ddots & \ddots & 0\\ 0& \cdots & 0& {{\rm N}}_{k}\end{array}\right],\text{\hspace{0.17em}}{\rm M}=\left[\begin{array}{cccc}{{\rm M}}_{1}& 0& \cdots & 0\\ 0& {{\rm M}}_{2}& \ddots & \vdots \\ \vdots & \ddots & \ddots & 0\\ 0& \cdots & 0& {{\rm M}}_{k}\end{array}\right],$$

$$R=\left[\begin{array}{cccc}{R}_{1}& 0& \cdots & 0\\ 0& {R}_{2}& \ddots & \vdots \\ \vdots & \ddots & \ddots & 0\\ 0& \cdots & 0& {R}_{k}\end{array}\right]$$

$${\Lambda}_{i}=\frac{1}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}{\mathsf{\lambda}}_{i}{\mathsf{\pi}}_{i}{{\rm M}}_{i},\text{\hspace{0.17em}}1\le i\le k,$$

$${R}_{i}=\frac{1}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}{\mathsf{\lambda}}_{i}{\mathsf{\pi}}_{i},$$

$${\eta}_{i}=\frac{{\mathsf{\pi}}_{i}{\mathsf{\lambda}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}.$$

This Corollary 1 follows directly from Theorem 1 applied to each diagonal block of the matrices (17) and (18).

Consider a queueing system with $L\ge 0$ waiting places that can be modeled by a Markov chain ${X}_{L}(t)$ with a generator of the form (7) with the block-diagonal matrices $\Lambda $, ${\rm M}$, ${\rm N}$ as in Equation (17). We partition the stationary probability vector ${x}_{L}$ corresponding to ${Q}_{L}$ as ${x}_{L}=(p,{q}_{0},\dots ,{q}_{L})$, where $p$ is an n-vector, and ${q}_{l}$ is partitioned as ${q}_{l}=({q}_{l,1},{q}_{l,2},\dots ,{q}_{l,k})$, with ${m}_{i}$-vectors ${q}_{l,i}$, $1\le i\le k,$$0\le l\le L$.

Let a generator${Q}_{L}$be irreducible, matrices$\Lambda $,${\rm M}$,${\rm N}$have the block-diagonal form (Equation (17)),${A}_{i}={\Lambda}_{i}+{{\rm N}}_{i}+{{\rm M}}_{i}$be an irreducible generator with the stationary probability vector${\mathsf{\pi}}_{i},$the matrices${\Lambda}_{i},$${R}_{i}$be given by Equations (19) and (20) respectively, and vectors$a$and$b=({b}_{1},{b}_{2},\dots ,{b}_{k})$satisfy the linear system

$$a{{\rm N}}_{0}+b{{\rm M}}_{1}=0,$$

$$a{\Lambda}_{0}+b({{\rm N}}_{1}+\Lambda )=0,$$

$$au+bu=1.$$

Then the stationary probability vector of ${Q}_{L}$ is given by
where ${c}_{L}$ is the normalization constant,

$$p={c}_{L}a,\text{\hspace{0.17em}}{q}_{0,i}={c}_{L}{b}_{i},\text{\hspace{0.17em}}{q}_{l,i}={c}_{L}{\rho}_{i}^{l-1}\frac{{b}_{i}{\mathsf{\lambda}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}{\mathsf{\pi}}_{i},\text{\hspace{0.17em}}1\le l\le L,\text{\hspace{0.17em}}1\le i\le k,$$

$${c}_{L}={(1+{\displaystyle \sum _{i=1}^{k}{s}_{L,i}})}^{-1},\text{\hspace{0.17em}}{s}_{L,i}=\{\begin{array}{cc}\frac{{b}_{i}{\mathsf{\lambda}}_{i}(1-{\rho}_{i}^{L})}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}(1-{\rho}_{i})},& {\rho}_{i}\ne 1,\\ L\frac{{b}_{i}{\mathsf{\lambda}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}},& {\rho}_{i}=1.\end{array}$$

Since the matrices ${\Lambda}_{i}$ and ${R}_{i}$ are given by Equations (19) and (20), under the Corollary 1 of Theorem 1, matrix (18) is the solution of Equation (3) and satisfies $R{\rm M}=\Lambda $. It implies that vectors $p$ and ${q}_{l}=({q}_{l,1},{q}_{l,2},\dots ,{q}_{l,k})$ given by Equation (25) satisfy the following relations:

$$p{{\rm N}}_{0}+{q}_{0}{{\rm M}}_{1}={c}_{L}(a{{\rm N}}_{0}+b{{\rm M}}_{1})=0,$$

$$p{\Lambda}_{0}+{q}_{0}{{\rm N}}_{1}+{q}_{1}{\rm M}={c}_{L}(a{\Lambda}_{0}+b{{\rm N}}_{1}+bR{\rm M})={c}_{L}(a{\Lambda}_{0}+b{{\rm N}}_{1}+b\Lambda )=0,$$

$${q}_{l+1}{\rm M}+{q}_{l}{\rm N}+{q}_{l-1}\Lambda ={q}_{l-1}({R}^{2}{\rm M}+R{\rm N}+\Lambda )=0,\text{\hspace{0.17em}}0<l<L,$$

$${q}_{L}(\Lambda +{\rm N})+{q}_{L-1}\Lambda ={q}_{L-1}(R\Lambda +R{\rm N}+\Lambda )={q}_{L-1}({R}^{2}{\rm M}+R{\rm N}+\Lambda )=0.$$

It follows that ${x}_{L}=(p,{q}_{0},\dots ,{q}_{L})$ is the stationary probability vector of the generator ${Q}_{L}.$ Power ${R}_{i}^{l}$ of the matrix ${R}_{i}$ is given by
Using this formula and equality (Equation (24)) we obtain the normalization constant (Equation (26)). This completes the proof of the theorem. □

$${R}_{i}^{l}={\eta}_{i}^{l-1}\frac{{\mathsf{\lambda}}_{i}{\mathsf{\pi}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}.$$

Note that the vector $(a,b)$ in Equations (22)–(24) gives the stationary distribution of the Markov chain ${X}_{0}(t)$ which models a queueing system without waiting places, and the stationary distribution ${x}_{L}$ of the Markov chain ${X}_{L}(t)$ depends on the number of waiting places L only via the normalizing constant ${c}_{L}$. Similar results are valid for the infinite QBD process $X(t)$ with generator $Q$ of the form (1).

Let generator$Q$be irreducible, matrices$\Lambda $,${\rm M}$,${\rm N}$have the block-diagonal form (17),${A}_{i}={\Lambda}_{i}+{{\rm N}}_{i}+{{\rm M}}_{i}$be an irreducible generator with stationary probability vector${\mathsf{\pi}}_{i},$matrices${\Lambda}_{i},$${R}_{i}$be given by (19) and (20) respectively, and vectors$a$and$b=({b}_{1},{b}_{2},\dots ,{b}_{k})$satisfy the linear system (22)–(24). Then process$X(t)$is positive recurrent if and only if${\mathsf{\pi}}_{i}{\mathsf{\lambda}}_{i}<{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}$for all$1\le i\le k.$In this case the stationary probability vector of$Q$is given by
where$c$is the normalization constant,

$$p=ca,\text{}{q}_{0,i}=c{b}_{i},\text{}{q}_{l,i}=c{\eta}_{i}^{l-1}\frac{{b}_{i}{\mathsf{\lambda}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}}{\mathsf{\pi}}_{i},\text{}l\ge 1,1\le i\le k,$$

$$c={\left(1+{\displaystyle \sum _{i=1}^{k}\frac{{b}_{i}{\mathsf{\lambda}}_{i}}{{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}-{\mathsf{\pi}}_{i}{\mathsf{\lambda}}_{i}}}\right)}^{-1}.$$

The proof of Equations (27) and (28) is similar to the proof of Equations (25) and (26). The spectral radius of the rate matrix $R$ is the maximum of the spectral radii if its diagonal blocks are given by Equation (21). Therefore, the process $X(t)$ is positive recurrent if and only if ${\mathsf{\pi}}_{i}{\mathsf{\lambda}}_{i}<{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}$ for all $1\le i\le k$. □

The following result is a direct consequence of Equation (25) for the stationary distribution of the processes ${X}_{L}(t)$ and of Equation (27) for the stationary distribution of the process $X(t)$.

If the conditions of Theorem 2 hold and${\mathsf{\pi}}_{i}{\mathsf{\lambda}}_{i}<{\mathsf{\pi}}_{i}{\mathsf{\mu}}_{i}$for all$1\le i\le k$then the stationary distributions of the processes${X}_{L}(t)$can be obtained by truncation of the stationary distribution of the process$X(t)$.

Since the first papers [14,15] on semi-open networks appeared, the use of these models has become increasingly widespread [16,17,18]. However, the exact analysis even of the simplest network, with two exponential single-server nodes, Poisson arrival process, and an infinite queue, remains a difficult task [15]. The semi-open networks with external queue are mainly analyzed by approximate methods and the matrix-geometric method [19,20,21,22,23].

Consider a semi-open network consisting of an internal network with v nodes and an external queue. The number of customers in the internal network cannot exceed N. If an arriving customer finds the internal network busy it waits in the external queue or is lost if the number of customers in the queue has reached its maximum value $L$. Such a network can be considered as a queuing system with N resource units. If all N resource units are occupied, the arriving customer waits in the external queue until a resource unit becomes available. Upon completion of its service in the internal network, the customer releases the occupied resource unit.

We consider semi-open networks in which the customers waiting in the external queue select a node for entering the internal network according to the following rule, which we call the replacement. If a customer leaves the system upon completing service at node i, then a customer from the external queue will start service also at node i, thus “replacing” the leaving customer. As a result, the number of customers in the external queue will decrease by one, but the number of customers at each node of the internal network remains unchanged. An example of a system in which only the replacement of customers leaving the system is possible is shown in Figure 1. It is a semi-open network with only one node where customers can complete service and through which waiting customers can enter the internal network.

Let the arrival process be Poisson with rate $\lambda $. The service rate ${\mu}_{{k}_{i},i}$ at node i depends on the number ${k}_{i}$ of customers at this node: ${\mu}_{0,i}=0$ and ${\mu}_{{k}_{i},i}>0$ for $1\le {k}_{i}\le N.$ A customer entering the internal network is sent to node i with probability ${a}_{i}$, ${\sum}_{i=1}^{n}{a}_{i}}=1$. A customer completing service at node i will either require service at node j with a probability ${q}_{i,j}$ or leave the system with probability ${q}_{i,0}=1-{\displaystyle {\sum}_{j=1}^{n}{q}_{i,j}}$. The visit ratios ${q}_{i},$ $i=1,2,\dots ,v,$ can be found as a solution of the linear system
and the nodes’ throughputs are given by ${\lambda}_{i}=\lambda {q}_{i}$, $i=1,2,\dots ,v$. Summing both sides of Equation (1) over j, we obtain a useful relation

$${q}_{j}={a}_{j}+{\displaystyle \sum _{i=1}^{n}{q}_{i}}{q}_{i,j},\text{}j=1,2,\dots ,v,$$

$$\sum _{i=1}^{n}{q}_{i}}{q}_{i,0}=1.$$

A semi-open network can be modeled by a Markov chain ${Y}_{L}(t)=({y}_{0}(t),{y}_{1}(t),\dots ,{y}_{n}(t)),$ where ${y}_{0}(t)$ represents the number of customers in the system at time t, and ${y}_{j}(t)$, $1\le j\le n$, represent the number of customers at node j of the internal network. The process ${Y}_{L}(t)$ has a state space $\cup}_{l=0}^{L}{Y}_{l$ with
where

$${Y}_{0}=\{(l,k)|k\in {Z}_{l},0\le l\le N\},\text{}{Y}_{l}=\{(l+N,k)|k\in {Z}_{N}\},\text{}1\le l\le L,$$

$${Z}_{l}=\{({k}_{1},\dots ,{k}_{v})|{k}_{j}\ge 0,\text{}j=1,\dots ,v,\text{}{k}_{1}+\dots +{k}_{v}=l\},\text{}0\le l\le N.$$

We assume that all states of the process $Y(t)$ communicate and denote its stationary distribution as ${\theta}_{l}=({\theta}_{l}(k))$, $(l,k)\in {Y}_{l}$, $0\le l\le N+L$.

The generator of the process $Y(t)$ has the block tri-diagonal form (Equation (18)) with the following nonzero transition rates
where

$$\gamma ((l,k),(l+1,k+{e}_{i}))=\lambda {a}_{i},\text{}k\in {Z}_{l},\text{}0\le lN,$$

$$\gamma ((l,k),(l+1,k))=\lambda ,\text{}k\in {Z}_{N},\text{}N\le lN+L,$$

$$\gamma ((l,k),(l-1,k-{e}_{i}))={\mu}_{{k}_{i}}{q}_{i,0},\text{}k\in {Z}_{l},\text{}1\le l\le N,$$

$$\gamma ((l,k),(l-1,k))=\mu (k),\text{}k\in {Z}_{N},\text{}Nl\le N+L,$$

$$\gamma ((l,k),(l,k+{e}_{j}-{e}_{i}))={\mu}_{{k}_{i}}{q}_{i,j},\text{}k\in {Z}_{l},\text{}1\le l\le N,$$

$$\mu (k)={\displaystyle \sum _{i=1}^{v}{\mu}_{{k}_{i},i}}{q}_{i,0}.$$

In Equation (8) the matrices $\Lambda $, ${\rm M}$, ${\rm N}$ are diagonal with elements $\gamma ((l,k),(l+1,n))=\lambda \delta (k,n)$, $\gamma ((l,k),(l-1,n))=\mu (k)\delta (k,n)$, $\gamma ((l,k),(l,n))=-(\lambda +\mu (k))\delta (k,n)$ respectively. Therefore, the stationary distributions for the finite and infinite semi-open networks with replacement can be found using Theorems 2 and 3.

For a semi-open network without an external queue the closed form solution can be easily derived because there exists an equivalent closed network of capacity N. The equivalent network can be obtained by adding a single-server node 0 with an infinite queue and exponentially distributed service times with parameter $\lambda $ to the original network. Node 0 plays the role of the external environment of the original network. When a customer arrives to or departs from the original network, it departs from or arrives to node 0 in the equivalent network. The stationary distribution ${x}_{l}=({x}_{l}(k))$, $k\in {Z}_{l}$, $l=0,1,\dots ,N$ of closed queueing networks with state-dependent service rates have a product form [24] and the solution $a=({x}_{0},{x}_{1},\dots ,{x}_{N-1})$ and $b={x}_{N}$ of the linear system (Equations (19)–(21)) is given by

$${x}_{l}(k)=\frac{1}{G(N)}{\displaystyle \prod _{i=1}^{n}\frac{{\lambda}_{i}^{{k}_{i}}}{{\mu}_{1,i}{\mu}_{2,i}\cdots {\mu}_{{k}_{i},i}}},k\in {Z}_{l},\text{}l=0,1,\dots ,N,$$

$$G(N)=1+{\displaystyle \sum _{l=1}^{N}{\displaystyle \sum _{k\in {X}_{l}}{\displaystyle \prod _{i=1}^{v}\frac{{\lambda}_{i}^{{k}_{i}}}{{\mu}_{1,i}{\mu}_{2,i}\cdots {\mu}_{{k}_{i},i}}}}}.$$

For closed queueing networks, there are several effective methods for computation of the normalization constant $G(N)$ [25].

By virtue of Theorem 2 for semi-open network with replacement and a finite external queue, the stationary distribution of the process ${Y}_{L}(t)$ is given by
where

$${\theta}_{l}(k)={c}_{L}{x}_{l}(k),k\in {Z}_{l},\text{}0\le l\le N,$$

$${\theta}_{l}(k)={c}_{L}{x}_{N}(k)\rho {(k)}^{l-N},\text{}k\in {Z}_{N},\text{}Nl\le N+L,$$

$${c}_{L}={(1+{\displaystyle \sum _{k\in {Z}_{N}}{s}_{L}(k)})}^{-1},\text{}{s}_{L}(k)=\{\begin{array}{cc}{x}_{N}(k)\rho (k)\frac{1-\rho {(k)}^{L}}{1-\rho (k)},& \rho (k)\ne 1,\\ {x}_{N}(k)L,& \rho (k)=1,\end{array}$$

$$\rho (k)=\frac{\lambda}{{\displaystyle \sum _{i=1}^{v}{\mu}_{{k}_{i},i}}{q}_{i,0}}.$$

It follows from Theorem 3 that a semi-open network with replacement and an infinite external queue is stable if and only if $\rho (k)<1$ for all states $k\in {Z}_{N}$. In this case, the stationary distribution of $Y(t)$ is given by

$${\theta}_{l}(k)={c}_{L}{x}_{l}(k),k\in {Z}_{l},\text{}0\le l\le N,$$

$${\theta}_{l}(k)={c}_{L}{x}_{N}(k)\rho {(k)}^{l-N},\text{}k\in {Z}_{N},\text{}lN,$$

$${c}_{L}={\left(1+{\displaystyle \sum _{k\in {Z}_{N}}\frac{{x}_{N}(k)\rho (k)}{1-\rho (k)}}\right)}^{-1}.$$

Similarly, it is possible to analyze semi-open queuing networks partitioned into subnetworks. Each such subnetwork consists of a subset of nodes of the internal network. If a customer leaves the system, then a customer from the external queue will proceed with a certain probability to a node in the same subnetwork in which the customer leaving the system was served. The generator of the related QBD process has the matrices $\Lambda $ of the form of Equation (17). The results of this paper can be used if the diagonal blocks of $\Lambda $ satisfy relation (19).

We have studied the properties of QBD processes with the rate matrix $R$ satisfying $R{\rm M}=\Lambda $, which guarantees that the stationary distribution of the finite QBD process is a truncation of the stationary distribution of the infinite QBD process. We have obtained the necessary and sufficient condition for the fulfillment of the equality $R{\rm M}=\Lambda $ and have derived the formulae for the rate matrix and for the superdiagonal blocks of the process generator. These matrices have rank one, and therefore the stationary distribution and the performance parameters of QBD processes can be easily found. The results obtained have been applied to a semi-open network in which customers from the external queue can enter the internal network only through the node from which the departing customer left the system.

Conceptualization, V.A.N. and K.E.S.; methodology and validation, V.A.N.; review and editing, original draft preparation, examples, Y.V.G.; supervision and project administration, K.E.S.

The publication has been supported by the Ministry of Education and Science of the Russian Federation (project No. 2.882.2017/4.6).

The authors declare no conflict of interest.

- Neuts, M.F. Matrix-Geometric Solutions an Algorithmic Approach; The Johns Hopkins University Press: Baltimore, MD, USA, 1981. [Google Scholar]
- Latouche, G.; Ramaswami, V. A logarithmic reduction algorithm for Quasi-Birth-Death processes. J. Appl. Probab.
**1993**, 30, 650–674. [Google Scholar] [CrossRef] - Meini, B. Solving QBD problems: The cyclic reduction algorithm versus the invariant subspace method. Adv. Perform. Anal.
**1998**, 1, 215–225. [Google Scholar] - Qiang, Y. On Latouche-Ramaswami’s logarithmic reduction algorithm for Quasi-Birth-and-Death processes. Stoch. Models
**2001**, 18, 449–467. [Google Scholar] - Neuts, M.F. Explicit steady-state solutions to some elementary queueing models. Oper. Res.
**1982**, 30, 480–489. [Google Scholar] [CrossRef] - Latouche, G.; Ramaswami, V. A general class of Markov processes with explicit matrix-geometric solutions. Oper. Res.
**1986**, 8, 209–218. [Google Scholar] - Leeuwaarden, J.S.H.; Squillante, M.; Winands, E. Quasi-birth-and-death processes, lattice path counting, and hypergeometric functions. J. Appl. Probab.
**2009**, 46, 507–520. [Google Scholar] [CrossRef] - Doroudi, S.; Fralix, B.; Harchol-Balter, M. Clearing analysis on phases: Exact limiting probabilities for skip-free, unidirectional, quasi-birth-death processes. Stoch. Syst.
**2016**, 6, 420–458. [Google Scholar] [CrossRef] - Hajek, B. Birth-and-Death Processes on the Integers with Phases and General Boundaries. J. Appl. Probab.
**1982**, 19, 488–499. [Google Scholar] [CrossRef] - Naumov, V. Modified Matrix-Geometric Solution for Finite QBD processes. In Advances in Algorithmic Methods for Stochastic Models; Latouche, G., Taylor, P., Eds.; Notable Publications: Neshanic Station, NJ, USA, 2000; pp. 257–264. [Google Scholar]
- Akar, N.; Sohraby, K. Finite and infinite QBD chains: A simple and unifying algorithmic approach. In Proceedings of the INFOCOM ‘97, Kobe, Japan, 7–11 April 1997; Volume 3, pp. 1105–1113. [Google Scholar]
- Naoumov, V. Matrix-multiplicative approach to quasi-birth-and-death processes analysis. In Matrix-Anal. Methods Stoch. Models; Alfa, A.S., Chakravarthy, S.R., Eds.; Marcel Dekker: New York, NY, USA, 1996; pp. 87–106. [Google Scholar]
- Krieger, U.R.; Naumov, V. Analysis of a delay-loss system with a superimposed markovian arrival process and state dependent service times. In Numer. Solut. Markov Chain; Plateau, B., Stewart, W.J., Silva, M., Eds.; Prensas Universitarias de Zaragoza: Zaragoza, Spain, 1999. [Google Scholar]
- Avi-Itzhak, B.; Heyman, D.P. Approximate Queuing Models for Multiprogramming Computer Systems. Oper. Res.
**1973**, 21, 1212–1230. [Google Scholar] [CrossRef] - Konheim, A.G.; Reiser, M. Finite Capacity Queuing Systems with Applications in Computer Modeling. SIAM J. Comput.
**1978**, 7, 210–229. [Google Scholar] [CrossRef] - Busacott, J.A.; Shanthikumar, J.G. Design of Manufacturing Systems Using Queueing Models. Queueing Syst.
**1992**, 12, 135–214. [Google Scholar] [CrossRef] - Cai, X.; Heragu, S.S.; Yang Liu, Y. Modeling Automated Warehouses Using Semi-Open Queueing Networks. In Handbook of Stochastic Models and Analysis of Manufacturing System Operations; Smith, J., Tan, B., Eds.; Springer: New York, NY, USA, 2003. [Google Scholar]
- Zou, B.P.; de Koster, R.; Xu, X.H. Operating Policies in Robotic Compact Storage and Retrieval Systems. Transp. Sci.
**2008**, 52, 788–811. [Google Scholar] [CrossRef] - Jia, J.; Heragu, S.S. Solving Semi-Open Queuing Networks. Oper. Res.
**2009**, 57, 391–401. [Google Scholar] [CrossRef] - Roy, D. Semi-Open Queuing Networks: A Review of Stochastic Models, Solution Methods and New Research Areas. Int. J. Prod. Res.
**2016**, 54, 1735–1752. [Google Scholar] [CrossRef] - Ekren, B.Y.; Heragu, S.S.; Krishnamurthy, A.; Malmborg, C.J. Matrix-Geometric Solution for Semi-Open Queuing Network Model of Autonomous Vehicle Storage and Retrieval System. Comput. Ind. Eng.
**2014**, 68, 78–86. [Google Scholar] [CrossRef] - Dhingraa, V.; Kumawat, G.L.; Roy, D.; de Koster, R. Solving Semi-open Queuing Networks with Time-varying Arrivals: An Application in Container Terminal Landside Operations. Eur. J. Oper. Res.
**2018**, 267, 855–876. [Google Scholar] [CrossRef] - Kim, J.; Dudin, A.; Dudin, S.; Kim, C. Analysis of a semi-open queueing network with Markovian arrival process. Perform. Eval.
**2018**, 120, 1–19. [Google Scholar] [CrossRef] - Baskett, F.; Chandy, K.M.; Muntz, A.A.; Palacios, F.G. Open, closed, and mixed networks of queues with different classes of customers. J. ACM
**1975**, 22, 248–260. [Google Scholar] [CrossRef] - Bolch, G.; Greiner, S.; de Meer, H.; Trivedi, K.S. Queueing Networks and Markov Chains: Modeling and Performance Evaluation with Computer Science Applications, 2nd ed.; John Wiley & Sons: Hoboken, NJ, USA, 2006. [Google Scholar]

© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).