The M/M/𝑁 Repairable Queueing System with Variable Breakdown Rates

Abstract

This paper considers the M/M/𝑁 repairable queuing system. The customers' arrival is a Poisson process. The servers are subject to breakdown according to Poisson processes with different rates in idle time and busy time, respectively. The breakdown servers are repaired by repairmen, and the repair time is an exponential distribution. Using probability generating function and transform method, we obtain the steady-state probabilities of the system states, the steady-state availability of the servers, and the mean queueing length of the model.

1. Introduction

In queuing researches, many researchers have studied the queuing system with repairable servers. Most of the works of the repairable queuing system deal with the single-server models [1–9]. The works about multiserver repairable systems is not sufficient. Mitrany and Avi-Itzhak [10] analyzed the model with 𝑁 units of servers and the same amount of repairmen, they obtained the steady-state mean queuing length of customers. Neuts and Lucantoni [11] studied the model with 𝑁 units of repairable servers and 𝑐 (0≤𝑐≤𝑁) repairmen by matrix analysis method and obtained the steady-state properties of the model.

In recent years, many flexible policies have been introduced to the repairable systems. Gray et al. [5] studied the model with a single server which may take a vacation in idle times and may breakdown in busy times; they obtained the mean queue length. Altman and Yechiali [12] presented a comprehensive analysis of the M/M/1 and M/G/1 queues, as well as of the M/M/𝑐 queue with server vacations; they obtained various closed-form results for the probability generating function (PGF) of the number of the customers. Zhang and Hou [6] analyzed an M/G/1 queue with working vacations and vacation interruptions; they obtained the queue length distribution and steady-state service status probability. Yang et al. [7] analyzed an M/G/1 queuing system with second optional service, server breakdowns, and general startup times under (𝑁,𝑝)-policy, they obtained the explicit closed-form expression of the joint optimum threshold values of (𝑁,𝑝) at the minimum cost. Chang et al. [8] studied the optimal management problem of a finite capacity M/H2/1 queuing system, where the unreliable server operates 𝐹-policy, a cost model is developed to determine the optimal capacity 𝐾, the optimal threshold 𝐹, the optimal setup rate, and the optimal repair rate at a minimum cost. Wang [9] used a quasi-birth-and-death (QBD) modeling approach to model queuing-inventory systems with a single removable server, performance measures are obtained by using both hybrid and standard procedures; an optimal control policy is proposed and verified through numerical studies.

The most works of repairable queuing system assumed that the server breakdown rate is constant, but the breakdown rate of a server may be variable in a real system. It is well known that many kinds of machine are easy to breakdown at their busy times, and some equipments may be easy to fail after a long idle period. For example, the tires of the truck prefer to breakdown when the truck is running on the road. On the other hand, the storage battery in an automobile may not work if the automobile is idle for long period. For the actual demands of the above cases, we study a multiserver repairable queuing system in this paper, and assume that the unreliable servers have different breakdown rates in their busy times and idle times, respectively.

The rest of this paper is organized as follows. Section 2 describes the model and gives the balance equations. Section 3 presents the equations of PGF. The steady-state availability of the system is derived in Section 4. The steady-state probabilities of the system states and mean queuing length are obtained in Section 5. Case analysis is given in Section 6. Section 7 presents the conclusions.

2. Model Description

The model characteristics are as follows. (1)There are 𝑁 units of identical servers in the system. The servers are subject to breakdown according to Poisson processes with different rates which are 𝜉₁ in idle times and 𝜉₂ in busy times, respectively. (2)Customers arrive according to a Poisson process with rate 𝜆. The service discipline is first come first served (FCFS). The service which is interrupted by a server breakdown will become the first one of the queue of customers. The service time distribution is an exponential distribution with parameter 𝜇. (3)There are 𝑐(1≤𝑐≤𝑁) reliable repairmen to maintain the unreliable servers. The repair discipline is first come first repaired (FCFR). The repair time distribution is an exponential distributions with parameter 𝜂. A server is as good as a new one after repair.

We define𝑋(𝑡)≡ the number of available servers in the system at the moment 𝑡(0⩽𝑋(t)⩽𝑁), 𝑌(𝑡)≡ the number of customers in the system at the moment 𝑡(0⩽𝑌(𝑡)).The stochastic process {𝑋(𝑡),𝑌(𝑡),𝑡≥0} is a two-dimensional Markov process which is called quasi-birth-and-death (QBD) process [11] with state space {(𝑖,𝑗),0≤𝑖≤𝑁,𝑗≥0}.

Let 𝑃_𝑖,𝑗(𝑡) denote the probability that the system is in a state of (𝑖,𝑗) at the moment 𝑡, and 𝑃_𝑖,𝑗 denote the steady-state probability of 𝑃_𝑖,𝑗(𝑡), then we have 𝑃_𝑖,𝑗=lim_𝑡→∞𝑃_𝑖,𝑗(𝑡),𝑖=0,1,2,…,𝑁,𝑗=0,1,2,…,0,other. (1) Assuming that the system is positive recurrent, the balance equations are as follows: 𝜆+𝑐𝜂𝑃_0,0=𝜉₁𝑃_1,0,𝜆+𝑐𝜂𝑃_0,𝑗=𝜆𝑃_0,𝑗−1+𝜉₂𝑃_1,𝑗,𝑗>0,𝜆+𝑐𝜂+𝑖𝜉₁𝑃_𝑖,0=𝑐𝜂𝑃_𝑖−1,0+𝜇𝑃_𝑖,1+(𝑖+1)𝜉₁𝑃_𝑖+1,0,𝑃0<𝑖≤𝑁−𝑐,𝑗=0,_𝑖,𝑗𝜆+𝑐𝜂+𝑗𝜇+𝑖−𝑗𝜉₁+𝑗𝜉₂=𝜆𝑃_{𝑖,𝑗−1}+𝑐𝜂𝑃_{𝑖−1,𝑗}+𝑗+1𝜇𝑃_𝑖,𝑗+1+𝑖+1−𝑗𝜉₁+𝑗𝜉₂𝑃_𝑖+1,𝑗,𝑃0<𝑖≤𝑁−𝑐,0<𝑗<𝑖,_𝑖,𝑗𝜆+𝑐𝜂+𝑖𝜇+𝑖𝜉₂=𝜆𝑃_{𝑖,𝑗−1}+𝑐𝜂𝑃_{𝑖−1,𝑗}+𝑖𝜇𝑃_𝑖,𝑗+1+𝜉₁+𝑗𝜉₂𝑃_𝑖+1,𝑗,𝑃0<𝑖≤𝑁−𝑐,𝑗=𝑖,_𝑖,𝑗𝜆+𝑐𝜂+𝑖𝜇+𝑖𝜉₂=𝜆𝑃_{𝑖,𝑗−1}+𝑐𝜂𝑃_{𝑖−1,𝑗}+𝑖𝜇𝑃_𝑖,𝑗+1+(𝑖+1)𝜉₂𝑃_𝑖+1,𝑗,𝑃0<𝑖≤𝑁−𝑐,𝑗>𝑖,_𝑖,0𝜆+(𝑁−𝑖)𝜂+𝑖𝜉₁=(𝑁−𝑖+1)𝜂𝑃_𝑖−1,0+𝜇𝑃_𝑖,1+(𝑖+1)𝜉₁𝑃_𝑖+1,0,𝑃𝑁−𝑐<𝑖<𝑁,𝑗=0,_𝑖,𝑗𝜆+(𝑁−𝑖)𝜂+𝑗𝜇+𝑖−𝑗𝜉₁+𝑗𝜉₂=(𝑁−𝑖+1)𝜂𝑃_{𝑖−1,𝑗}+𝜆𝑃_{𝑖,𝑗−1}+𝑗+1𝜇𝑃_𝑖,𝑗+1+𝑖+1−𝑗𝜉₁+𝑗𝜉₂𝑃_𝑖+1,𝑗,𝑃𝑁−𝑐<𝑖<𝑁,0<𝑗<𝑖,_𝑖,𝑗𝜆+(𝑁−𝑖)𝜂+𝑖𝜇+𝑖𝜉₂=𝜆𝑃_{𝑖,𝑗−1}+(𝑁−𝑖+1)𝜂𝑃_{𝑖−1,𝑗}+𝑖𝜇𝑃_𝑖,𝑗+1+𝜉₁+𝑗𝜉₂𝑃_𝑖+1,𝑗𝑃,𝑁−𝑐<𝑖<𝑁,𝑗=𝑖,_𝑖,𝑗𝜆+(𝑁−𝑖)𝜂+𝑖𝜇+𝑖𝜉₂=𝜆𝑃_{𝑖,𝑗−1}+(𝑁−𝑖+1)𝜂𝑃_{𝑖−1,𝑗}+𝑖𝜇𝑃_𝑖,𝑗+1+(𝑖+1)𝜉₂𝑃_𝑖+1,𝑗𝑃,𝑁−𝑐<𝑖<𝑁,𝑗>𝑖,_𝑁,0𝜆+𝑁𝜉₁=𝜂𝑃_𝑁−1,0+𝜇𝑃_𝑁,1𝑃,𝑖=𝑁,𝑗=0,_𝑁,𝑗𝜆+𝑗𝜇+𝑁−𝑗𝜉₁+𝑗𝜉₂=𝜆𝑃_{𝑁,𝑗−1}+𝜂𝑃_{𝑁−1,𝑗}+𝑗+1𝜇𝑃_𝑁,𝑗+1,𝑃𝑖=𝑁,0<𝑗<𝑁,_𝑁,𝑗𝜆+𝑁𝜇+𝑁𝜉₂=𝜆𝑃_{𝑁,𝑗−1}+𝜂𝑃_{𝑁−1,𝑗}+𝑁𝜇𝑃_𝑁,𝑗+1,𝑖=𝑁,𝑗≥𝑁. (2)

Here, we give the derivation of the second equation of (2). Since the process {𝑋(𝑡),𝑌(𝑡),𝑡≥0} is a vector Markov process of continuous time, we write the equations of the state of (0,𝑗) by considering the transitions occurring between the moments 𝑡 and 𝑡+Δ𝑡(Δ𝑡>0) as follows: 𝑃_0,𝑗(𝑡+Δ𝑡)=𝑃_0,𝑗−1(𝑡)𝜆Δ𝑡+𝑃_1,𝑗(𝑡)𝜉₂Δ𝑡+𝑃_0,𝑗(𝑡)1−𝜆+𝑐𝜂Δ𝑡+𝑜(Δ𝑡), (3) then we have 𝑃_0,𝑗(𝑡+Δ𝑡)−𝑃_0,𝑗(𝑡)=𝑃_0,𝑗−1(𝑡)𝜆Δ𝑡+𝑃_1,𝑗(𝑡)𝜉₂Δ𝑡−𝑃_0,𝑗(𝑡)𝜆+𝑐𝜂Δ𝑡+𝑜(Δ𝑡),𝑃_0,𝑗(𝑡+Δ𝑡)−𝑃_0,𝑗(𝑡)Δ𝑡=𝑃_0,𝑗−1(𝑡)𝜆+𝑃_1,𝑗(𝑡)𝜉₂−𝑃_0,𝑗(𝑡+𝑜())𝜆+𝑐𝜂Δ𝑡.Δ𝑡 (4)(5) Letting Δ𝑡→0 in (5), we have 𝑃_0,𝑗(𝑡)^′=𝑃_0,𝑗−1(𝑡)𝜆+𝑃_1,𝑗(𝑡)𝜉₂−𝑃_0,𝑗(𝑡)𝜆+𝑐𝜂. (6) If the system is positive recurrent, we have the formulas lim_𝑡→∞𝑃_0,𝑗(𝑡)^′=0 [13]. Letting 𝑡→0 in (6), we obtain the second equation of (2). The derivations of other formulas in (2) are similar.

3. Equations of Probability Generating Functions

The PGFs of the number of customers are defined as follows: 𝐺_𝑖(𝑧)≡^∞_𝑗=0𝑧^𝑗𝑃_𝑖,𝑗,𝐺(𝑧)≡^𝑁_𝑖=0𝐺_𝑖(𝑧),0≤𝑖≤𝑁,|𝑧|≤1. (7) Then 𝐺_𝑖(1)=^∞_𝑗=0𝑃_𝑖,𝑗(𝑖=0,1,2,…,𝑁), (8) where 𝐺_𝑖(1) is the steady-state probability that the number of the available servers of the system is 𝑖. Hence, ^𝑁_𝑖=0𝐺_𝑖(1)=1. (9)

Multiplying the two sides of every equation of (2) by 𝑧^𝑗+1, and summing over 𝑗(𝑗=0,1,2,…) for every 𝑖, we obtain 𝑧𝜆−𝜆𝑧+𝑐𝜂𝐺₀(𝑧)−𝑧𝜉₂𝐺₁(𝑧)=𝑧𝜉₁−𝜉₂𝑃_1,0,−𝑐𝜂𝑧𝐺_𝑖−1(𝑧)+𝑧𝑖𝜇+𝑖𝜉₂+𝜆+𝑐𝜂−𝜆𝑧²−𝑖𝜇𝐺_𝑖(𝑧)−(𝑖+1)𝑧𝜉₂𝐺_𝑖+1(𝑧)=^𝑖−1_𝑚=0𝜉₂−𝜉₁(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚+1+^𝑖_𝑚=0𝜉₁−𝜉₂(𝑖+1−𝑚)𝑃_𝑖+1,𝑚𝑧^𝑚+1+(𝑧−1)^𝑖−1_𝑚=0𝜇(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚−(),0<𝑖≤𝑁−𝑐,𝑁−𝑖+1𝜂𝑧𝐺_𝑖−1(𝑧)+𝑧𝑖𝜇+𝑖𝜉₂+𝜆+(𝑁−𝑖)𝜂−𝜆𝑧²−𝑖𝜇𝐺_𝑖(𝑧)−(𝑖+1)𝑧𝜉₂𝐺_𝑖+1(𝑧)=^𝑖−1_𝑚=0𝜉₂−𝜉₁(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚+1+^𝑖_𝑚=0𝜉₁−𝜉₂(𝑖+1−𝑚)𝑃_𝑖+1,𝑚𝑧^𝑚+1+(𝑧−1)^𝑖−1_𝑚=0𝜇(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚,𝑁−𝑐<𝑖<𝑁,−𝜂𝑧𝐺_𝑁−1(𝑧)+𝑧𝑁𝜇+𝑁𝜉₂+𝜆−𝜆𝑧²−𝑁𝜇𝐺_𝑁(𝑧)=^𝑁−1_𝑚=0𝜉₂−𝜉₁(𝑁−𝑚)𝑃_𝑁,𝑚𝑧^𝑚+1+(𝑧−1)^𝑁−1_𝑚=0𝜇(𝑁−𝑚)𝑃_𝑁,𝑚𝑧^𝑚. (10)

We give some explanations of (10), the first equation of (2) multiplied by 𝑧, we get 𝜆+𝑐𝜂𝑃_0,0𝑧=𝜉₁𝑃_1,0𝑧. (11) The second equation of (2) multiplied by 𝑧^𝑗+1, we get 𝜆+𝑐𝜂𝑃_0,𝑗𝑧^𝑗+1=𝜆𝑃_0,𝑗−1𝑧^𝑗+1+𝜉₂𝑃_1,𝑗𝑧^𝑗+1,𝑗>0. (12) Summing (11) and (12) over 𝑗 and using (7), we obtain the first equation of (10). The other equations of (10) are obtained in the same way.

In order to simplify (10), we defined the following notations: 𝑓_𝑖(𝑧)≡𝑖𝜇+𝑖𝜉₂+𝜆+𝑐𝜂𝑧−𝜆𝑧²𝑓−𝑖𝜇,𝑖=0,1,2,…,𝑁−𝑐,_𝑖(𝑧)≡𝑖𝜇+𝑖𝜉₂+𝜆+(𝑁−𝑖)𝜂𝑧−𝜆𝑧²−𝑖𝜇,𝑖=𝑁−𝑐+1,𝑁−𝑐+2,…,𝑁, (13)𝐴(𝑧)≡^{𝑓0(𝑧)−𝜉2𝑧−𝑐𝜂𝑧𝑓1(𝑧)−2𝜉2𝑧⋱⋱⋱⋱⋱⋱−𝑐𝜂𝑧𝑓𝑁−𝑐+1(𝑧)−()𝜉𝑁−𝑐+22𝑧}_{()𝑐−1𝜂𝑧𝑓𝑁−𝑐+2(𝑧)−()𝜉𝑁−𝑐+32𝑧⋱⋱⋱⋱⋱⋱−2𝜂𝑧𝑓𝑁−1(𝑧)−𝑁𝜉2𝑧−𝜂𝑧𝑓𝑁(𝑧)}, (14)𝑏₀(𝑧)≡𝑧𝜉₁−𝜉₂𝑃_1,0,𝑏_𝑖(𝑧)≡^𝑖−1_𝑚=0𝜉₂−𝜉₁(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚+1+^𝑖_𝑚=0𝜉₁−𝜉₂(𝑖+1−𝑚)𝑃_𝑖+1,𝑚𝑧^𝑚+1+(𝑧−1)^𝑖−1_𝑚=0𝜇(𝑖−𝑚)𝑃_𝑖,𝑚𝑧^𝑚𝑏,0<𝑖<𝑁,_𝑁(𝑧)≡^𝑁−1_𝑚=0𝜉₂−𝜉₁(𝑁−𝑚)𝑃_𝑁,𝑚𝑧^𝑚+1+(𝑧−1)^𝑁−1_𝑚=0𝜇(𝑁−𝑚)𝑃_𝑁,𝑚𝑧^𝑚,𝑏(𝑧)≡𝑏₀(𝑧)𝑏₁(𝑧)⋮𝑏_𝑁(𝑧),𝑔(𝑧)≡𝐺₀(𝑧)𝐺₁(𝑧)⋮𝐺_𝑁(𝑧). (15) Using the above notations, (10) is rewritten as follows: 𝐴(𝑧)𝑔(𝑧)=𝑏(𝑧). (16) Using Cramer’s rule we obtain |𝐴(𝑧)|𝐺_𝑖(𝑧)=𝐴_𝑖(𝑧),𝑖=0,1,2,…,𝑁, (17) where |𝐴(𝑧)| denotes the determinant of 𝐴(𝑧), and 𝐴_𝑖(𝑧) is a matrix obtained by replacing the (𝑖+1)th column of 𝐴(𝑧) with 𝑏(𝑧). In (17), the functions of 𝑧 are continuous and bounded in the interval [0,1], so the equations in (17) are valid in the interval [0,1] no matter |𝐴(𝑧)|=0 or not.

4. Steady-State Availability

In this section, we discuss the steady-state availabilities 𝐺_𝑖(1)(𝑖=0,1,2,…,𝑁).

Letting 𝑧=1 in (10) we obtain 𝑐𝜂𝐺₀(1)−𝜉₂𝐺₁(1)=𝜉₁−𝜉₂𝑃_1,0,𝑐𝜂𝐺_𝑖−1(1)−𝑖𝜉₂𝐺_𝑖(1)−𝑐𝜂𝐺_𝑖(1)−(𝑖+1)𝜉₂𝐺_𝑖+1(1=𝜉₁)−𝜉₂^𝑖_𝑚=1𝑚𝑃_{𝑖,𝑖−𝑚}−𝜉₁−𝜉₂^𝑖+1_𝑚=1𝑚𝑃_{𝑖+1,𝑖+1−𝑚},[(0<𝑖≤𝑁−𝑐,𝑁−𝑖−1)]𝜂𝐺_𝑖−1(1)−𝑖𝜉₂𝐺_𝑖(1)−(𝑁−𝑖)𝜂𝐺_𝑖(1)−(𝑖+1)𝜉₂𝐺_𝑖+1(1=𝜉₁)−𝜉₂^𝑖_𝑚=1𝑚𝑃_{𝑖,𝑖−𝑚}−𝜉₁−𝜉₂^𝑖+1_𝑚=1𝑚𝑃_{𝑖+1,𝑖+1−𝑚},𝑁−𝑐<𝑖<𝑁,𝜂𝐺_𝑁−1(1)−𝑁𝜉₂𝐺_𝑁(1)=𝜉₁−𝜉₂^𝑁_𝑚=1𝑚𝑃_{𝑁,𝑁−𝑚}. (18) The (𝑁+1) equations of (18) are simplified to 𝑁 independent equations, joined with (9), we have 𝑐𝜂𝐺_𝑖(1)−(𝑖+1)𝜉₂𝐺_𝑖+1(1)=𝜉₁−𝜉₂^𝑖+1_𝑚=1𝑚𝑃_{𝑖+1,𝑖+1−𝑚}(),0≤𝑖≤𝑁−𝑐,𝑁−𝑖𝜂𝐺_𝑖(1)−(𝑖+1)𝜉₂𝐺_𝑖+1(1)=𝜉₁−𝜉₂^𝑖+1_𝑚=1𝑚𝑃_{𝑖+1,𝑖+1−𝑚},𝑁−𝑐≤𝑖≤𝑁−1,_𝑁𝑖=0𝐺_𝑖(1)=1. (19) Using (2), 𝑃_𝑖,𝑗(0≤𝑗≤𝑖−1,1≤𝑖≤𝑁) are reduced to 𝑃_𝑖,0(1≤𝑖≤𝑁) which will be solved in Section 5. Given 𝑃_𝑖,0(1≤𝑖≤𝑁), we obtain 𝐺_𝑖(1)(𝑖=0,1,2,…,𝑁) by solving (19), then the steady-state availability of the system is as follows: 𝐴=1−𝐺₀(1). (20)

5. Steady-State Probabilities of the System States and Mean Queuing Length

5.1. The Roots of |𝐴(𝑧)| in the Interval (0,1)

In order to get the steady-state probabilities 𝑃_𝑖,0(1≤𝑖≤𝑁), we need all 𝑁 independent linear equations. Further, for getting the linear equations of 𝑃_𝑖,0(1≤𝑖≤𝑁) we need the roots of |𝐴(𝑧)| in the interval (0,1), so we discuss the roots of |𝐴(𝑧)| as follows.

We define the following notations: 𝑄₀(𝑧)≡1,𝑄₁(𝑧)≡𝑓_𝑁(𝑧), (21)(𝑧)≡𝑓_𝑁−1(𝑧)−𝑁𝜉₂𝑧−𝜂𝑧𝑓_𝑁(𝑧),𝑄⋮⋮_𝑁(𝑧)≡𝑓₁(𝑧)−2𝜉₂𝑧−𝑐𝜂𝑧𝑓₂(𝑧)−3𝜉₂𝑧⋱⋱⋱−𝑐𝜂𝑧𝑓_{𝑁−𝑐+1}(𝑧)−(𝑁−𝑐+2)𝜉₂𝑧(𝑐−1)𝜂𝑧𝑓_{𝑁−𝑐+2}(𝑧)−(𝑁−𝑐+3)𝜉₂𝑧⋱⋱⋱−2𝜂𝑧𝑓_𝑁−1(𝑧)−𝑁𝜉₂𝑧−𝜂𝑧𝑓_𝑁(𝑧), (22)𝑄_𝑁+1(𝑧)≡|𝐴(𝑧)|,(0≤𝑧≤∞), (23) where 𝑄_𝑖(𝑧)(𝑖=1,2,…,𝑁) is formed by the last 𝑖 rows and last 𝑖 columns of |𝐴(𝑧)|. According to (14) we have 𝑄_𝑘+1(𝑧)=𝑓_𝑁−𝑘(𝑧)𝑄_𝑘(𝑧)−(𝑁−𝑘+1)𝜉₂𝑘𝜂𝑧²𝑄_𝑘−1(𝑧),𝑄1≤𝑘≤𝑐−1,_𝑘+1(𝑧)=𝑓_𝑁−𝑘(𝑧)𝑄_𝑘(𝑧)−(𝑁−𝑘+1)𝜉₂𝑐𝜂𝑧²𝑄_𝑘−1(𝑧),𝑐≤𝑘≤𝑁. (24) The properties of 𝑄_𝑖(𝑧)(𝑖=1,2,…,𝑁) and the necessary proofs are as follows. (a)𝑄₀(𝑧) has no roots. (b)𝑄_𝑘(𝑧) and 𝑄_𝑘+1(𝑧) have no common roots in the interval (0,∞)(𝑘=1,2,…,𝑁).

Proof. Suppose that (b) is not true and 𝑧₀(>0) is a common roots of 𝑄_𝑘(𝑧) and 𝑄_𝑘+1(𝑧), then 𝑄_𝑘−1(𝑧₀)=0 due to (24). Similarly, 𝑄_𝑘−2(𝑧₀)=0, so we get 𝑄₀(𝑧₀)=0 which contradicts the statement (a).

(c) If 𝑧₀ is a positive root of 𝑄_𝑘(𝑧)(𝑘=1,2,…,𝑁), then 𝑄_𝑘−1(𝑧₀) and 𝑄_𝑘+1(𝑧₀) are opposite in sign. This property readily follows from (24). (d) 𝑄_𝑘(1)>0, 𝑘=0,1,2,…,𝑁.

Proof. Substituting 𝑧=1 into 𝑄_𝑘(𝑧), we get 𝑄_𝑘(1)=(𝑁−𝑘+1)(𝑁−𝑘+2)⋯(𝑁−1)𝑁𝜉^𝑘₂>0, 𝑘=1,2,…,𝑁, and 𝑄₀(1)=1.

(e) 𝑄_𝑁+1(1)=|𝐴(1)|=0, and 𝑄_𝑁+1(0)=|𝐴(0)|=0.

Proof. The first row of |𝐴(𝑧)| has a common factor 𝑧, so 𝑄_𝑁+1(0)=|𝐴(0)|=0. The sum of every column of |𝐴(𝑧)| has a common factor (𝑧−1). Replacing every element of the last row of |𝐴(𝑧)| with the sum of the corresponding column and extracting the common factor (𝑧−1), |𝐴(𝑧)| is written as follows: |𝐴(𝑧)|=𝑧(𝑧−1)×𝜆+𝑐𝜂−𝜆𝑧−𝜉₂0⋯00−𝑐𝜂𝑧𝑓₁(𝑧)−2𝜉₂𝑧⋯000−𝑐𝜂𝑓₂(𝑧)⋮⋯00000⋯𝑓_𝑁−1(𝑧)−𝑁𝜉₂𝑧().−𝜆𝑧−𝜆𝑧+𝜇−𝜆𝑧+2𝜇⋯−𝜆𝑧+𝑁−1𝜇−𝜆𝑧+𝑁𝜇 (25) If we define 𝐷(𝑧)≡𝜆+𝑐𝜂−𝜆𝑧−𝜉₂0⋯00−𝑐𝜂𝑧𝑓₁(𝑧)−2𝜉₂𝑧⋯000−𝑐𝜂𝑓₂(𝑧)⋮⋯00000⋯𝑓_𝑁−1(𝑧)−𝑁𝜉₂𝑧(),−𝜆𝑧−𝜆𝑧+𝜇−𝜆𝑧+2𝜇⋯−𝜆𝑧+𝑁−1𝜇−𝜆𝑧+𝑁𝜇 (26)
then |𝐴(𝑧)|=𝑧(𝑧−1)𝐷(𝑧). (27)

(f) Sign[𝑄_𝑘(0)]=(−1)^𝑘(𝑘=0,1,2,…,𝑁).

Proof. From the definitions of 𝑓_𝑖(𝑧)(𝑖=0,1,2,…,𝑁), we got 𝑓₀(0)=0 and 𝑓_𝑘(0)<0(𝑘=1,2,…,𝑁), so we get this property from (24).

(g) Sign[𝑄_𝑘(+∞)]=(−1)^𝑘(𝑘=0,1,2,…,𝑁+1).

Proof. It is since the highest power term of 𝑄_𝑘(𝑧) is (−𝜆𝑧²)^𝑘(𝑘=0,1,2,…,𝑁+1) and the sign of 𝑄_𝑘(+∞) is determined by its highest power term.

Theorem 1. If 𝐷(1)>0, the polynomial |𝐴(𝑧)| has exactly (𝑁−1) distinct roots in the interval (0,1).

Proof. Since 𝑄₁(𝑧)=𝑓_𝑁(𝑧)=[𝑁(𝜇+𝜉₂)+𝜆]𝑧−𝜆𝑧²−𝑁𝜇, 𝑄₁(𝑧) is a 2-power polynomial of 𝑧. Further, we find that 𝑄₁(1)=𝜉₂>0 and 𝑄₁(0)=−𝜇<0, so 𝑄₁(𝑧) has two distinct roots which are denoted by 𝑧_1,1(0<𝑧_1,1<1) and 𝑧_1,2(>1).
With the fact that 𝑧_1,1 and 𝑧_1,2 are roots of 𝑄₁(𝑧), and 𝑄₀(𝑧)=1>0, according to the property (c) or (24), we get 𝑄₂(𝑧_1,1)<0 and 𝑄₂(𝑧_1,2)<0.
𝑄₂(𝑧) is a 4-power polynomial of 𝑧. From the properties (c), (d), (f), and (g), we find that 𝑄₂(𝑧) has one and only one root in each interval of (0,𝑧_1,1), (𝑧_1,1,1), (1,𝑧_1,2), and (𝑧_1,2,+∞).
So on, we find that 𝑄_𝑁(𝑧) is a 2𝑁-power polynomial of 𝑧, it has 𝑁 distinct roots in the interval (0,1) and 𝑁 distinct roots in the interval (1,∞). We denote the 2𝑁 roots of 𝑄_𝑁(𝑧) by 𝑧_𝑁,𝑖(𝑖=1,2,…,2𝑁) orderly.
From the properties (c), (d), (e), and (f), we find that |𝐴(𝑧)| has one and only one root in each interval (𝑧_𝑁,𝑖,𝑧_𝑁,𝑖+1)(𝑖=1,2,…,𝑁−1,𝑁+1,…,2𝑁−1), all of them are 2(𝑁−1) distinct roots of |𝐴(𝑧)|.
Since |𝐴(1)|=0 and 𝐷(1)=(|𝐴(𝑧)|)^′_𝑧=1>0, it has a real number 𝜀 (>0) satisfies 1+𝜀<𝑧_𝑁,𝑁+1 and |𝐴(1+𝜀)|>0. On the other hand, from (c) and (d) we get Sign𝑄_𝑁+1𝑧_𝑁,𝑖=(−1)^𝑁+𝑖,𝑖=𝑁+1,𝑁+2,…,2𝑁, (28) then Sign[𝑄_𝑁+1(𝑧_𝑁,𝑁+1)]=(−1)^𝑁+𝑁+1 or |𝐴(𝑧_𝑁,𝑁+1)|=𝑄_𝑁+1(𝑧_𝑁,𝑁+1)<0, so |𝐴(𝑧)| has at least one root in the interval (1,𝑧_𝑁,𝑁+1).
From (28), we get Sign[|𝐴(𝑧_𝑁,2𝑁)|]=Sign[𝑄_𝑁+1(𝑧_𝑁,2𝑁)]=(−1)^𝑁+2𝑁=(−1)^3𝑁=(−1)^𝑁. From the property (g), we get Sign[|𝐴(+∞)|]=Sign[𝑄_𝑁+1(+∞)]=(−1)^𝑁+1. So we know that |𝐴(𝑧)| has at least one root in the interval (𝑧_𝑁,2𝑁,∞).
From the properties (e) and (f), we know that 0 and 1 are roots of |𝐴(𝑧)|.
In conclusion, |𝐴(𝑧)| is a 2(𝑁+1)-power polynomial of 𝑧, it has 2(𝑁+1) distinct roots at most. Now we make certain all roots of |𝐴(𝑧)| and find that it has 𝑁−1 distinct roots in the interval (0,1).

From the proof, we find that the (𝑁−1) distinct roots in the interval (0,1) of |𝐴(𝑧)| are also the roots of 𝐷(𝑧).

5.2. Steady-State Probabilities

Assuming that the system parameters meet 𝐷(1)>0. Letting 𝑧_𝑘(𝑘=1,2,…,𝑁−1) denote the roots of |𝐴(𝑧)| in the interval (0,1). Substituting 𝑧₁ in (17), we obtain a set of linear equations about the steady-state probabilities of 𝑃_𝑖,0(𝑖=1,2,…,𝑁), but these equations are similar to each other. However, the equations belong to different 𝑧_𝑘(𝑘=1,2,…,𝑁−1) are independent mutually, so we obtain (𝑁−1) independent equations by the (𝑁−1) different roots of 𝑧_𝑘, respectively.

In the following, we discuss about the 𝑁th-independent linear equation of 𝑃_𝑖,0(𝑖=1,2,…,𝑁). Similar to (27), |𝐴_𝑖(𝑧)| is written as follows: 𝐴_𝑖(𝑧)=𝑧(𝑧−1)𝐷_𝑖(𝑧),𝑖=0,1,2,…,𝑁, (29) where 𝐷_𝑖(𝑧)=𝑐𝜂+𝜆(1−𝑧)−𝜉₂⋯𝜉₁−𝜉₂𝑃_1,0⋯0−𝑐𝜂𝑧𝑓₁(𝑧)⋯𝑏₁(𝑧)⋯00−𝑐𝜂𝑧⋯𝑏₂(𝑧)⋯0⋮⋮00⋯𝑏_𝑁−1(𝑧)⋯−𝑁𝜉₂𝑧−𝜆𝑧−𝜆𝑧+𝜇⋯∑^𝑁𝑖_𝑖=1∑_𝑚=1𝑚𝜇𝑝_{𝑖,𝑖−𝑚}𝑧^𝑖−𝑚_{(𝑖+1)thcolumn}.⋯−𝜆𝑧+𝑁𝜇 (30) Substituting (27) and (29) into (17), we obtain 𝐷(𝑧)𝐺_𝑖(𝑧)=𝐷_𝑖(𝑧),𝑖=0,1,2,…,𝑁. (31) Substituting 𝑧=1 in (31), we obtain 𝐷(1)𝐺_𝑖(1)=𝐷_𝑖(1),𝑖=0,1,2,…,𝑁. (32) From (19), we know that 𝐺_𝑖(1)(𝑖=0,1,2,…,𝑁) can be expressed by 𝑃_𝑖,0(𝑖=1,2,…,𝑁), so (32) are linear equations of 𝑃_𝑖,0(𝑖=1,2,…,𝑁), but they are similar to each other. However, every equation of (32) is independent with the (𝑁−1) linear equations obtained by the roots of |𝐴_𝑖(𝑧)| in the interval (0,1). Then all independent linear equations of 𝑃_𝑖,0(𝑖=1,2,…,𝑁) are as follows: 𝐴₀𝑧_𝑘𝐷(1)𝐺=0,𝑘=1,2,…,𝑁−1,₀(1)=𝐷₀(1). (33) Further, from (29), (33) is equivalent to the follows: 𝐷₀𝑧_𝑘𝐷(1)𝐺=0,𝑘=1,2,…,𝑁−1,₀(1)=𝐷₀(1). (34)

The steady-state probabilities of 𝑃_𝑖,0(𝑖=0,1,2,…,𝑁) are obtained by solving (34). Using 𝑃_𝑖,0(𝑖=0,1,2,…,𝑁) and (2), we obtain the other steady-state probabilities of 𝑃_𝑖,𝑗(𝑖=0,1,2,…,𝑁,𝑗=1,2,…).

5.3. Mean Queuing Length

After getting the probabilities 𝑃_𝑖,𝑗(0≤𝑗≤𝑖−1,1≤𝑖≤𝑁) in (19) we obtain 𝐺_𝑖(1)(𝑖=0,1,2,…,𝑁) by solving (19) and obtain the steady-state availability (𝐴) of the model by (20).

From (31), we obtain 𝐺_𝑖(𝑧)=𝐷_𝑖(𝑧)𝐷(𝑧),𝐷(𝑧)≠0,|𝑧|≤1,𝑖=0,1,2,…,𝑁. (35) The PGF of 𝐺(𝑧) is obtained by (7). Using the property of PGF [13] we obtain the steady-state mean queuing length is as follows: (𝑧)𝐿=𝑑𝐺𝑑𝑧_𝑧=1. (36)

6. Case Analysis

We analyze the case of 𝑁=2 and 𝑐=1 in this section. According to the above discussion, the determinant |𝐴(𝑧)| (or 𝐷(𝑧)) of this case has only one root 𝑧₁ in the interval (0,1). The notations of this case are as follows: 𝑓₀(𝑧)=𝜆+𝜂𝑧−𝜆𝑧²,𝑓₁(𝑧)=𝜇+𝜉₂+𝜆+𝜂𝑧−𝜆𝑧²𝑓−𝜇,₂(𝑧)=2𝜇+2𝜉₂+𝜆𝑧−𝜆𝑧²𝑏−2𝜇,₀(𝑧)=𝜉₁−𝜉₂𝑃_1,0𝑏𝑧,₁(𝑧)=𝜇+𝜉₂−𝜉₁𝑃_1,0𝑧+²_𝑚=1𝑚𝜉₁−𝜉₂𝑃_2,2−𝑚𝑧^3−𝑚−𝜇𝑃_1,0=𝜇𝑃_1,0(𝑧−1)+𝑧𝜉₁−𝜉₂𝑃_2,1𝑧+2𝑃_2,0−𝑃_1,0,𝑏₂(𝑧)=²_𝑚=1𝑚𝜇+𝜉₂−𝜉₁𝑃_2,2−𝑚𝑧^3−𝑚−²_𝑚=1𝑚𝜇𝑃_2,2−𝑚𝑧^2−𝑚,𝐴(𝑧)=𝑓₀(𝑧)−𝜉₂𝑧0−𝜂𝑧𝑓₁(𝑧)−2𝜉₂𝑧0−𝜂𝑧𝑓₂(𝑧),|𝐴(𝑧()×)|=𝑧𝑧−1𝜆+𝜂−𝜆𝑧−𝜉₂0−𝜂𝑧𝜇+𝜉₂+𝜆+𝜂𝑧−𝜆𝑧²−𝜇−2𝜉₂𝑧,𝐷(𝑧)=−𝜆𝑧−𝜆𝑧+𝜇−𝜆𝑧+2𝜇𝜆+𝜂−𝜆𝑧−𝜉₂0−𝜂𝑧𝜇+𝜉₂+𝜆+𝜂𝑧−𝜆𝑧²−𝜇−2𝜉₂𝑧−𝜆𝑧−𝜆𝑧+𝜇−𝜆𝑧+2𝜇=−𝜂+𝜆−𝑧𝜆𝑧𝜆−2𝜇−𝜇+𝑧𝜂+𝜆−𝑧𝜆+𝜇+𝑧𝜉₂𝑧𝜆−2𝜂−3𝜆+3𝑧𝜆+2𝜂−2𝑧𝜆+2𝜆𝜇−2𝑧𝜆𝜉₂, (37) then 𝐷(1)=𝜂−𝜉₂0−𝜂𝜂+𝜉₂−2𝜉₂𝜉−𝜆−𝜆+𝜇−𝜆+2𝜇=2𝜇₂𝜂+𝜂²−𝜆2𝜉²₂+2𝜉₂𝜂+𝜂², (38) and 𝐷(1)>0 is equivalent to 2𝜇𝜉₂𝜂+𝜂²−𝜆2𝜉²₂+2𝜉₂𝜂+𝜂²>0, (39) or 𝜆𝜇<2𝜉₂𝜂+2𝜂²2𝜉²₂+2𝜉₂𝜂+𝜂²=2𝜂/𝜉₂+𝜂1+𝜉₂/𝜉₂+𝜂². (40) The left of (40) is the mean service quantities that all customers need per unit time. The right of (40) is the mean service quantities that the two servers provide per unit time. So (40) is the necessary and sufficient condition of recurrence of the system.

Equations (19) in this case are as follows: 𝜂𝐺₀(1)−𝜉₂𝐺₁(1)=𝜉₁−𝜉₂𝑃_1,0,𝜂𝐺₁(1)−2𝜉₂𝐺₂(1)=²_𝑚=1𝑚𝜉₁−𝜉₂𝑃_2,2−𝑚,_2𝑖=0𝐺_𝑖(1)=1. (41)

Equations (2) in this case are as follows: 𝜆+𝜂𝑃_0,0=𝜉₁𝑃_1,0𝑃,𝑖=0,𝑗=0,𝜆+𝜂_0,𝑗=𝜆𝑃_0,𝑗−1+𝜉₂𝑃_1,𝑗,𝑖=0,𝑗>0,𝜆+𝜂+𝜉₁𝑃_1,0=𝜂𝑃_0,0+𝜇𝑃_1,1+2𝜉₁𝑃_2,0,𝑖=1,𝑗=0,𝜆+𝜂+𝜉₂+𝜇𝑃_1,1=𝜆𝑃_1,0+𝜂𝑃_0,1+𝜇𝑃_1,2+𝜉₁+𝜉₂𝑃_2,1,𝑖=1,𝑗=1,𝜆+𝜂+𝜉₂+𝜇𝑃_1,𝑗=𝜆𝑃_1,𝑗−1+𝜂𝑃_0,𝑗+𝜇𝑃_1,𝑗+1+2𝜉₂𝑃_2,𝑗,𝑖=1,𝑗>1,𝜆+2𝜉₁𝑃_2,0=𝜂𝑃_1,0+𝜇𝑃_2,1,𝑖=2,𝑗=0,𝜆+𝜉₁+𝜉₂+𝜇𝑃_2,1=𝜆𝑃_2,0+𝜂𝑃_1,1+2𝜇𝑃_2,2,𝑖=2,𝑗=1,𝜆+2𝜇+2𝜉₂𝑃_2,𝑗=𝜆𝑃_2,𝑗−1+𝜂𝑃_1,𝑗+2𝜇𝑃_2,𝑗+1,𝑖=2,𝑗≥2. (42) Using (42), we obtain 𝑃_2,1=2𝜉₁+𝜆𝑃_2,0−𝜂𝑃_1,0𝜇. (43) Using (41) and (43), 𝐺₀(1), 𝐺₁(1), and 𝐺₂(1) are expressed in an algebraic expressions of 𝑃_1,0 and 𝑃_2,0.

If (40) is satisfied, we obtain 𝑧₁ by solving 𝐷(𝑧)=0. (44) Equations (34) in this case are as follows: 𝐷₀𝑧₁𝐷(1)𝐺=0,₀(1)=𝐷₀(1). (45) Solving (45), we obtain 𝑃_1,0 and 𝑃_2,0.

Using (42), we obtain 𝑃_𝑖,𝑗 (𝑖=0,1,2, 𝑗=1,2,…).

Using (20) and (41), we obtain the steady-state availability 𝐴.

For the mean queuing lengths, we have 𝐷₀(𝑧)=𝜉₁−𝜉₂𝑃_1,0−𝜉₂0𝑏₁(𝑧)𝜇+𝜉₂+𝜆+𝜂𝑧−𝜆𝑧²−𝜇−2𝜉₂𝜇𝑃_2,1𝑧+2𝑃_2,0+𝑃_1,0,𝐷−𝜆𝑧+𝜇𝜆𝑧+2𝜇₁(𝑧)=𝜆+𝜂−𝜆𝑧𝜉₁−𝜉₂𝑃_1,00−𝜂𝑧𝑏₁(𝑧)−2𝜉₂−𝜆𝑧𝜇𝑃_2,1𝑧+2𝑃_2,0+𝑃_1,0,𝐷𝜆𝑧+2𝜇₂(𝑧)=𝜆+𝜂−𝜆𝑧−𝜉₂𝜉₁−𝜉₂𝑃_1,0−𝜂𝑧𝜇+𝜉₂+𝜆+𝜂𝑧−𝜆𝑧²−𝜇𝑏₁(𝑧)−𝜆𝑧−𝜆𝑧+𝜇𝜇𝑃_2,1𝑧+2𝑃_2,0+𝑃_1,0, (46)𝐺(𝑧)=𝐷₀(𝑧)+𝐷₁(𝑧)+𝐷₂(𝑧)𝐷(𝑧). (47) Using (36), we obtain the mean queuing lengths 𝐿 is as follows: 𝐿=𝑑𝐺(𝑧)𝑑𝑧_𝑧=1=1𝜆𝜂²−2𝜇𝜂²+2𝜉₂𝜂𝜆−𝜂𝜇+𝜆𝜉₂²×−𝜂−2𝜂𝜆+2𝜆²+2𝜂𝜇−5𝜆𝜇+2𝜇²+𝜉₂−4𝜂𝜆+3𝜆²+2𝜂𝜇−4𝜆𝜇−4𝜆𝜉₂×𝜇𝜂²+2𝜉₁𝜂+2𝜉₁𝜉₂𝑃_1,0+𝜇𝜂²+𝜂𝜉₂+𝜂𝜉₁+2𝜉₁𝜉₂2𝑃_2,0+𝑃_2,1+𝜇−𝜂²𝜆+2𝜇𝜂²−2𝜉₂𝜂𝜆−𝜂𝜇+𝜆𝜉₂×𝜂²−2𝜂𝜆+2𝜂𝜇−2𝜆𝜉₂+2𝜉₁𝜂−𝜆𝜉₁+2𝜇𝜉₁+2𝜉₁𝜉₂𝑃_1,0+𝜂²+𝜂𝜉₂+𝜂𝜉₁+2𝜉₁𝜉₂𝑃_2,1+𝜂²−2𝜂𝜆+𝜂𝜇+𝜂𝜉₂−2𝜆𝜉₂+𝜂𝜉₁−𝜆𝜉₁+2𝜉₁𝜉₂2𝑃_2,0+𝑃_2,1. (48)

Numerical Example. Letting 𝑁=2, 𝑐=1, 𝜆=1, 𝜉₁=0.5, 𝜉₂=1, 𝜂=1, and 𝜇=2, we have 𝜆𝜇=12<45=2𝜂/𝜉₂+𝜂1+𝜉₂/𝜉₂+𝜂². (49) The roots of 𝐷(𝑧)=0 are as follows: 𝑧₁=0.349123,𝑧₂𝑧=1.84513,₃=4.46896,𝑧₄=8.33679, (50) and only 𝑧₁ in the interval (0,1).

Solving (45), we obtain 𝑃_1,0=0.151375,𝑃_2,0=0.297974. (51) Using (42) we obtain 𝑃_0,0=0.037844,𝑃_2,1=0.146599. (52) Solving (41), we obtain 𝐺₀(1)=0.280333,𝐺₁(1)𝐺=0.35602,₂(1)=0.363647. (53) Finally, the availability and mean queuing lengths of this example are as follows: 𝐴=0.719667,𝐿=6.04039. (54)

The other numerical results are shown in Table 1. All the system parameters in Table 1 satisfy (40).

Table 1 

The availability 𝐴 and mean queuing length 𝐿 (𝑁=2, 𝑐=1, and 𝜆=1).

We find that the mean queuing length (𝐿) decreases with the increasing of the parameter 𝜇 in Table 1, it is because of the greater service rate the less customers in the system. Furthermore, we find that the availability (𝐴) increases with the increasing of the parameter 𝜇, where 𝜉₁<𝜉₂ (the cases: 𝜉₁=0,𝜉₂=0.5;𝜉₁=0.3,𝜉₂=0.5;𝜉₁=0.5,𝜉₂=1); on the contrary, the availability decreases with the increasing of the parameter 𝜇, where 𝜉₁>𝜉₂ (the case: 𝜉₁=1,𝜉₂=0.8); otherwise, the availability is constant, where 𝜉₁=𝜉₂ (the case: 𝜉₁=0.5,𝜉₂=0.5).

7. Conclusions

In Section  5.1, the inequality 𝐷(1)>0 of Theorem 1 is the necessary and sufficient condition for the system to be positive recurrent, and a probability explanation of this condition is given by (40).

We find that the idle time breakdown rate 𝜉₁ does not appear in (40). This is because the busy time breakdown rate 𝜉₂ is at work when the number of the customers is greater than or equal to the number of the available servers, and the criteria of positive recurrence depends on the busy time breakdown rate.

A case analysis is given to illustrate the analysis of this paper, and the numerical results indicate that the variation of breakdown rates has a significant effect on the steady-state availability and steady-state queue length of the system.

Acknowledgments

This paper is supported by the National Natural Science Foundation of China no. 71071133 and the Natural Science Foundation of Hebei Province no. G2012203136.