On: 26 April 2014, At: 00:15 Publisher: Taylor & Francis
Informa Ltd Registered in England and Wales Registered Number: 1072954 Registered office: Mortimer House, 37-41 Mortimer Street, London W1T 3JH, UK
Optimization: A Journal of
Mathematical Programming and
Operations Research
Publication details, including instructions for authors and subscription information:
http://www.tandfonline.com/loi/gopt20
A heuristic algorithm for the
optimization of a retrial system with
Bernoulli vacation
Jau-Chuan Ke a , Chia-Huang Wu b & Wen Lea Pearn b a
Department of Applied Statistics , National Taichung Institute of Technology , Taichung , Taiwan , ROC
b
Department of Industrial Engineering and Management , National Chiao Tung University , Taiwan , ROC
Published online: 25 May 2011.
To cite this article: Jau-Chuan Ke , Chia-Huang Wu & Wen Lea Pearn (2013) A heuristic algorithm for the optimization of a retrial system with Bernoulli vacation, Optimization: A Journal of Mathematical Programming and Operations Research, 62:3, 299-321, DOI: 10.1080/02331934.2011.579966
To link to this article: http://dx.doi.org/10.1080/02331934.2011.579966
PLEASE SCROLL DOWN FOR ARTICLE
Taylor & Francis makes every effort to ensure the accuracy of all the information (the “Content”) contained in the publications on our platform. However, Taylor & Francis, our agents, and our licensors make no representations or warranties whatsoever as to the accuracy, completeness, or suitability for any purpose of the Content. Any opinions and views expressed in this publication are the opinions and views of the authors, and are not the views of or endorsed by Taylor & Francis. The accuracy of the Content should not be relied upon and should be independently verified with primary sources of information. Taylor and Francis shall not be liable for any losses, actions, claims, proceedings, demands, costs, expenses, damages, and other liabilities whatsoever or howsoever caused arising directly or indirectly in connection with, in relation to or arising out of the use of the Content.
This article may be used for research, teaching, and private study purposes. Any substantial or systematic reproduction, redistribution, reselling, loan, sub-licensing, systematic supply, or distribution in any form to anyone is expressly forbidden. Terms &
and-conditions
A heuristic algorithm for the optimization of a retrial system with
Bernoulli vacation
Jau-Chuan Kea*, Chia-Huang Wub
and Wen Lea Pearnb
a
Department of Applied Statistics, National Taichung Institute of Technology, Taichung, Taiwan, ROC;bDepartment of Industrial Engineering and Management, National Chiao
Tung University, Taiwan, ROC
(Received 12 November 2010; final version received 6 April 2011)
In this study, we consider an M/M/c retrial queue with Bernoulli vacation under a single vacation policy. When an arrived customer finds a free server, the customer receives the service immediately; otherwise the customer would enter into an orbit. After the server completes the service, the server may go on a vacation or become idle (waiting for the next arriving, retrying customer). The retrial system is analysed as a quasi-birth-and-death process. The sufficient and necessary condition of system equilibrium is obtained. The formulae for computing the rate matrix and stationary probabilities are derived. The explicit close forms for system performance measures are developed. A cost model is constructed to determine the optimal values of the number of servers, service rate, and vacation rate for minimizing the total expected cost per unit time. Numerical examples are given to demonstrate this optimization approach. The effects of various parameters in the cost model on system performance are investigated.
Keywords: Bernoulli vacation schedule; matrix-geometric method; quasi-Newton method; retrial; single vacation policy
1. Introduction
Retrial queueing system is characterized by the feature that the arriving customers who, on encountering the busy server, join a retrial queue called orbit. An arbitrary customer in the orbit generates a stream of repeated requests that is independent of the rest of customers in the orbit. This situation arises in telephony, where an arriving call is not allowed to await the termination of a busy signal. Such queueing systems play important roles in the analysis of many telephone switching systems, telecommunication networks and computer systems. Review of retrial queue literature could be found in Yang and Templeton [44], Falin and Templeton [23] and Artalejo [2]. A number of applications of retrial queues in science and engineering can be found in Kulkarni and Liang [28].
*Corresponding author. Email: jauchuan@ntit.edu.tw
ß 2013 Taylor & Francis
Vol. 6 , No. , 299–2 3 321, http://dx.doi.org/10.1080/02331934.2011.57
Apart from its practical interest, due to its more accurate representation of several congestion phenomena, the multi-server retrial queue raises interesting mathematical and computational remarks. The investigation of the multi-server retrial queues is essentially more difficult than those with single server. Explicit formulae for the stationary distribution of M/M/c retrial queue are known only when the number of servers c is no more than two. Most multi-server retrial queues can be modelled by a level-dependent quasi-birth-and-death (QBD) process. The main feature of its infinitesimal generator is the spatial heterogeneity caused by transitions due to repeated attempts. This lack of homogeneity supports the analytical complexity of retrial models. Many interesting studies were devoted to an approximate approach of the stationary probabilities for system states (Falin [22], Bright and Taylor [9], Neuts and Rao [36], Stepanov [38], Artalejo and Pozo [7], Breuer et al. [8], Chakravarthy and Dudin [10]). Recently, Gomez-Corral [24] gave a detailed bibliographical guide to the analysis of retrial queues through matrix analytic techniques.
It is worth noting that the truncation models seem to be the most convenient method for obtaining reliable numerical solutions for the M/M/c retrial queue. For example, Falin [22] assumed that the retrial rate becomes infinite when the number of customers in orbit exceeds a level M. It means that, when the number of customers in the system is greater than M, the system performs as an ordinary M/M/1 queue with arrival rate and service rate c, so that 5 c is a sufficient and necessary condition for system ergodicity. Neuts and Rao [36] and Artalejo and Pozo [7] proposed several models in this direction and provided efficient approximate solutions to the stationary distribution of the M/M/c retrial queue. As related works, a number of studies investigated the computation of the other system characteristics, such as the distributions of busy period, successful and blocked (unsuccessful) retrials, for the multi-server retrial queue of type M/M/c. The readers can refer to Artalejo et al. [4], Amador and Artalejo [1] and others. Artalejo et al. [5,6] presented an algorithmic analysis of the maximum number of customers in orbit (and in the system) during a busy period for the M/M/c retrial queue. The multi-server retrial queueing problems are extensively studied as mentioned earlier. However, in the literature, there are no detailed studies on multi-server retrial queue with a vacation at each service completion instantly.
Alternatively, queueing models with server vacations are practical models for performance analysis of manufacturing systems, local area networks and data communication systems. Past works on vacation queueing models include those with single-server and multiple-server systems. Surveys on the single-server vacation models have been reported by Doshi [21] and Takagi [41]. The variations and extensions of these vacation models were developed by several researchers such as Lee et al. [31,32], Krishna Reddy et al. [27], Choudhury [13,14], Shomrony and Yechiali [37], Ke and Chu [26] and many others. For the multiple-server vacation models, there are only a limited number of studies due to the complexity of the systems. The M/M/c queue with exponential vacations was first studied by Levy and Yechiali [33]. Chao and Zhao [11] investigated a GI/M/c vacation system and provided an algorithm to compute the performance measures. Tian et al. [42] gave a detailed study of the M/M/c vacation systems in which all servers take
multiple vacation policy when the system is empty. Zhang and Tian [45,46] and Xu and Zhang [43] analysed the M/M/c vacation systems with a ‘partial server multiple vacation policy’ in which some servers (only the idle ones) take single or multiple vacations.
Studies on various queueing models in the past are characterized by common feature; servers always serve the waiting customers in the queue until all customers are served exhaustively or the number of the waiting customers is dropped to predetermined level. In reality, however, it may occur that the service process requires to be temporarily stopped for overhauling at the end of a service. This overhauling can be utilized as a vacation in the presented model. For example, consider a production process with a number of machines (or c machines). A number of investigations (Madan et al. [34], Choudhury and Madan [18,19], Tadj et al. [39,40], and Choudhury et al. [20]) have recently appeared in queueing literature in which the single server provides to each service with Bernoulli schedule vacation (BSV). The so-called BSV means that when the service of a unit is completed, the server may leave for a vacation of random length with probability p to serve the next unit with probability 1 p (Choudhury and Madan [18,19]). Analytic steady-state solutions of a multi-sever retrial queue with Bernoulli schedules under a single vacation policy (BSV) have not been found. Multi-server vacation models are more flexible and applicable in practice than the queueing models with single server. Existing research works, including those mentioned above, have not addressed the analytical study and optimization issue in the multiple-server retrial queues in which the server may take a vacation upon his each service completion. This motivates us to discuss an M/M/c retrial queue with BSV by applying matrix analytic approach.
Recently, Choudhury [15,16] investigated the M/G/1 and M[x]/G/1 queue with two phases of heterogeneous service and Bernoulli vacation schedule which operate under various retrial policies. Some extensive stationary analyses of the queueing system were carried out including the system size distribution and orbit size distribution. In the following year, Choudhury and Deka [17] dealt with the steady-state behaviour of M[x]/G/1 retrial queue with second optional service, unreliable server and Bernoulli admission mechanism. The above-mentioned model generalizes both M[x]/G/1 retrial queue with server breakdown and Bernoulli admission mechanism as well as M[x]/G/1 queue with second optional service and unreliable server. Furthermore, Ke and Chang [25] derived the mathematical model of M[x]/(G1, G2)/1 retrial queue under Bernoulli vacation schedules with general
repeated attempts and starting failures. A practical mail system example was presented. Later, Langaris and Dimitriou [29] investigated a single-server queueing with n-phases of service and (n 1) types of retrial customers. Some numeri-cal results under exponentially distributed service time were provided. Artalejo [3] presented a bibliography on retrial queues made during the past decade 2000–2009.
This study considers an M/M/c retrial queue where primary customers arrive as a Poisson process with parameter . An arriving primary customer finding one or more servers available (free) gets service immediately. On the other hand, if the primary customer finds all servers busy, he joins the orbit and tries to get the service later.
There are c channels (servers) that provide service for the arrivals and the service times are assumed to be exponentially distributed with mean 1=. Each server can serve only one customer at a time. At each service completion instant of a server, the server may take a vacation of random length with probability p or wait to serve the next arrival with probability q(¼1 p). The vacation times follow an exponen-tially distributed with a parameter . Furthermore, each customer staying in the orbit makes repeated attempts independently and the inter-retrial time is assumed to be exponentially distributed with parameter . Upon requesting service from the orbit, customer who finds all c servers busy always rejoins the orbit; this manner continues until he is eventually served. It is assumed that the number of customers in the orbit that is allowed to conduct retrials have an upper bound N (Neuts and Rao [36] and Artalejo and Pozo [7]). This implies that the probability of a repeated attempt during ðt, t þ dtÞ, given that j customers in the orbit at time t, is jdt þ oðdtÞ, where
j¼minf j, Ng Moreover, the process of primary arrivals, service times, and
inter-retrial times are assumed to be mutually independent. Conveniently, we represent this multi-server system with Bernoulli vacation as M/M/c/BSV retrial queue.
This article is organized as follows. In Section 2, the QBD model of the M/M/c/ BSV retrial queue is set up. The computable form of the rate matrix is derived and the stable condition is obtained using the matrix-geometric property. In Section 3, an efficient algorithm is developed to find the stationary probabilities by matrix-geometric method. In Section 4, some system performance measures are derived. In Section 5, a cost model is developed to determine the optimal number of servers, service rate and vacation rate, simultaneously, in order to minimize the total expected cost per unit time. The quasi-Newton method and direct search method are implemented to deal with the optimization tasks. Some numerical examples are provided to illustrate the optimization procedures. In Section 6, conclusions are made with some remarks.
2. M/M/c/BSV retrial queue
For M/M/c/BSV retrial queue system, the state of the system can be described by the pair ði, j, kÞ, i ¼ 0, 1, 2, . . . , c, j ¼ 0, 1, 2, . . . and k ¼ 0, 1, 2, . . . , c i, where i denotes the number of busy servers, j the number of customers in orbit (sources of repeated demands) and k the number of vacation servers. According to system assumptions, the number of customers in orbit allowed to conduct retrials is restricted to an appropriate number N (N4c), so the retrial rate is j¼min{ j, N},
j 0 and one server would go on vacation with probability p ( p40) or resumes service with probability q ¼ 1 p at a service completion instant. The cus-tomers upon the server get services immediately as i þ k5c. The new arriving customer who finds all c servers busy (i þ k c) always rejoins the retrial group (orbit).
In steady state, the steady-state probability is defined as
Pki,j probability, that is, there are i busy servers, j customers in orbit and k vacation servers, where 0 i þ k c and j ¼ 0, 1, 2, . . ..
2.1. Matrix representation of M/M/c/BSV retrial queue
The infinitesimal generator Q of the QBD describing the M/M/c/BSV retrial queueing system is Q ¼ A0 B C1 A1 B C2 A2 B . . . . . . . . . CN AN1 B CN1 AN B CN AN B . . . . . . . . . 2 6 6 6 6 6 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 7 7 7 7 7 5 : ð1Þ
The entries B, Ajðj4 0Þ and Cjðj4 1Þ are block-diagonal matrices of order
ðc þ1Þðc þ 2Þ=2 defined by B ¼ b0 b1 . . . bc1 bc 2 6 6 6 6 4 3 7 7 7 7 5 and Cj¼ c0j c1 j . . . cc1 j cc j 2 6 6 6 6 6 6 4 3 7 7 7 7 7 7 5 , j ¼ 1, 2, . . .
where sub-matrices bi and ci
j are ðc þ 1 iÞ ðc þ 1 iÞ square matrices with
elements bi½c þ1 i, c þ 1 i ¼ 0 e:w ( and c i j½k, k þ 1 ¼ j, 1 k c i 0 e:w Aj¼ Y0j X0 Z1 Y1j X1 Z2 Y2j X2 . . . . . . . . . . . . . . . . . . Zc1 Yc1j Xc1 Zc Ycj 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 , j ¼0, 1, 2, . . .
where Xi is a ðc þ 1 iÞ ðc iÞ matrix with Xi½k þ1, k ¼ kp, 1 k c i, Zi a ðc iÞ ðc þ 1 iÞ matrix with Zi½k, k ¼ i, 1 k c i and Yij a square matrix of order ðc þ 1 iÞ with elements
Yij½k, k þ 1 ¼ , 1 k c i Yij½k þ1, k ¼ kð1 pÞ, 1 k c i Yij½1, 1 ¼ ½ þ ði 1Þ þ j Yij½k, k ¼ ½ þ ðk þ 1Þ þ i þ j, 2 k c i Yij½c þ1 i, c þ 1 i ¼ ½ þ i þ ðc iÞ 8 > > > > > < > > > > > : :
The detailed descriptions of the above matrices (for c ¼ 3) are given in the Appendix.
Let & ¼ ½&0, &1, &2, . . . with &i¼ ½P00,i, P01,i, . . . , P0c,i, P10,i, P11,i, . . . , P1c1,i,
P00,i, . . . , Pc1
0,i , Pc11,i , Pc0,i, i ¼ 0, 1, 2, . . . be the unique solution to &Q ¼ 0 and
&e ¼ 1, where e is a column vector with all elements equal to 1. It is noted that the vector & ¼ ½&0, &1, &2, &3, . . . with the following properties
&Nþk¼&NRk, for k 1: ð2Þ
The matrix R is the unique non-negative solution with spectral radius less than 1 of the equation
B þ RANþR2CN¼0: ð3Þ
From Neuts [35] and Latouche and Ramaswami [30], it is known that R is given by limn!1Rn, where the sequence {Rn} is defined by
R0¼0, and Rnþ1¼ BA1N R2nCNA1N , for n 0: ð4Þ
The sequence {Rn} is monotone so that R could be evaluated from (4) by
successive substitutions.
2.2. Stability condition
It is also known (Theorem 3.1.1 of Neuts [35]) that the steady-state probability vector exists if and only if
xBe 5 xCNe, ð5Þ
where x is the invariant probability of the matrix F ¼ CNþANþB. Here, x satisfies
xF ¼ 0 and xe ¼ 1. First we solve xF ¼ 0, where x ¼ ½x0
0, x01. . . , x0c, . . . ,
xc10 , xc1
1 , xc0. We can get the following ðc þ 1Þðc þ 2Þ=2 equations:
For k ¼ 0,
ð þ NÞx00þqx01þx10¼0, ð6-1aÞ
ð þ NÞx0i1 ð þ i þ NÞx0i þ ði þ1Þqx0iþ1þx1i ¼0, 1 i c 1, ð6-1bÞ
ð þ NÞx0c1cx0c ¼0: ð6-1cÞ For 1 k c 1,
pxk11 ð þ N þ kÞxk0þqxk1þ ðk þ1Þxkþ10 ¼0, ð6-2aÞ
ð þ NÞxki2þipxk1i ½ þ N þ ði 1Þ þ kxki1
þiqxki þ ðk þ1Þxkþ1i1 ¼0, 2 i c k, ð6-2bÞ
ðc þ1 kÞ pxk1cþ1kþ ð þ NÞxkck1 ½ðc kÞ þ kxkck¼0: ð6-2cÞ
For k ¼ c,
pxc11 cxc0¼0: ð6-3aÞ Using a effective Maple software to solve Equations (6-1a)–(6-2c), the following results are derived
xki ¼ c!
ck
i!k!ð þ NÞcikipckx c
0, 0 i þ k c: ð7Þ
Then, using the normalization condition xe ¼ 1, xc
0 can be determined as xc0¼ X c k¼0 Xck i¼0 c!ck i!k!ð þ NÞcikipck " #1 ð8Þ
Substituting B and CNinto Equation (5) and doing some routine manipulations,
then we have Nð1 PFÞ4 PFull, ð9Þ where PFull¼ Xc i¼0 xcii ¼X c i¼0 c!i i!ðc iÞ!ipix c 0 ¼ 1 þ p c Xc k¼0 Xck i¼0 c!ck i!k!ð þ NÞcikipck " #1 , ð10Þ
which is referred to the probability that all normal working (non-vacation) servers are busy (i.e. i þ k ¼ c). That is, the system will be stable if the expected successful retrial rate is greater then the expected arrival rate of ‘orbit’.
3. Steady-state solution
Under the stability condition, the stationary probability vector & exists. We deal with the steady-state equations using matrix technique. The steady-state equations are given by
&0A0þ&1C1¼0, ð11aÞ
&i1B þ&iAiþ&iþ1Ciþ1¼0, 1 i N 1, ð11bÞ
&N1B þ&NANþ&NRCN¼0, ð11cÞ
&NRi1NB þ&NRiNANþ&NRiþ1NCN¼0, N þ1 i, ð11dÞ
X1 i¼0
&ie ¼ 1: ð12Þ
After doing some routine manipulations to Equations (11a)–(11c) recursively, we have
&0¼&1C1ðA0Þ1¼&11,
&i1¼&iCi½ði1B þ Ai1Þ1¼&ii, 2 i N,
ð13Þ
and
&NNB þ&NANþ&NRCN¼0: ð14Þ
Consequently, &ið0 i N 1Þ in Equation (13) can be written in terms of
&N as &0¼&N1i¼Ni, &1 ¼&N2i¼Ni, . . . , &N1¼&NNi¼Ni and the rest
steady-state vector ½&N, &Nþ1, &Nþ2, . . . can be determined recursively as
&i¼&NRiN, for i N. Therefore, once the steady-state probability &N is
obtained, the steady-state solutions ½&0, &1, &2, . . . , &N1, &N, &Nþ1, . . . are
deter-mined. The steady-state probability &N can be solved by Equation (14) with the
following normalization equation X1
i¼0
&ie ¼ ½&0þ&1þ þ&N1þ&Nþ&Nþ1þ&Nþ2þ e
¼ ½&N 1 i¼Niþ&N 2 i¼Niþ þ&N N
i¼Niþ&Nþ&NR þ&NR
2þ e ¼&N XN k¼1 k i¼Niþ ðI RÞ 1 " # e ¼ 1: ð15Þ
where I denotes the identity matrix with suitable size. Solving Equations (14) and (15) in accordance with Cramer’s rule, &N can be obtained. Then, the prior state
probabilities ½&0, &1, &2, . . . , &N1 are computed from (13) and
½&Nþ1, &Nþ2, &Nþ3, . . . are gained by the formula &i¼&NRiN, i N þ 1. The
solution procedure of steady-state probabilities is summarized as follows: Algorithm Recursive Solver
Step 1 Set 1¼C1ðA0Þ1
Step 2 For i from 2 to N, set i¼Ci½ði1B þ Ai1Þ1.
Step 3 For k from 1 to N, set (k¼ k i¼Ni.
Step 4 Solving &NNB þ&NANþ&NRCN¼0, &N½PNk¼1(kþ ðI RÞ1e ¼ 1
and obtain steady-state probability &N.
Step 5 Construct steady-state probability &i as follows:
(a) if 0 i N, assign &i¼&N(iþ1,
(b) if N i, assign &iþ1¼&iR,
4. System performance measures
There are several system descriptors (system performance measures) of the M/M/c/BSV retrial queue, such as the expected number of busy servers (denoted by E½B), the expected number of vacation servers (denoted by E ½V) and the expected number of customers in orbit (denoted by E ½Orbit), which can be
evaluated from the steady-state probabilities. The explicit expressions for E ½B, E ½V, and E ½Orbit are given by
E ½B ¼X 1 j¼0 &jv ¼ X N1 j¼0
&jv þ&Nv þ&NRv þ&NR2v þ
¼X
N1
j¼0
&N(jþ1v þ&Nv þ&NRv þ&NR2v þ
¼&N XN j¼1 (jþ ðI RÞ1 " # v ð16Þ E ½V ¼X 1 j¼0 &ju ¼ X N1 j¼0
&ju þ&Nu þ&NRu þ&NR2u þ
¼X
N1
j¼0
&N(jþ1u þ&NðI RÞ1u
¼&N XN j¼1 (jþ ðI RÞ1 " # u ð17Þ E ½Orbit ¼X 1 j¼1 j&je ¼ X N1 j¼1
j&N(jþ1e þ N&Ne þ ðN þ 1Þ&NRe þ ðN þ 2Þ&NR2e þ
¼X
N
j¼2
ðj 1Þ&N(je þ&N½NðI RÞ1þRðI RÞ2e
¼N XN j¼2 ðj 1Þ(jþNðI RÞ1þRðI RÞ2 " # e ð18Þ where v ¼ 0, 1, . . . , c |fflfflfflfflfflffl{zfflfflfflfflfflffl} #¼cþ1 , 0, 1, . . . , c 1 |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} #¼c , . . . , 0, 1 |{z} #¼2 , 0 2 4 3 5 and u ¼ 0, 0, . . . , 0 |fflfflfflfflfflffl{zfflfflfflfflfflffl} #¼cþ1 , 1, 1, . . . , 1 |fflfflfflfflfflffl{zfflfflfflfflfflffl} #¼c , . . . , c 1, c 1 |fflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflffl} #¼2 , c 2 4 3 5
are column vectors with dimension ðc þ 1Þðc þ 2Þ=2.
4.1. System performance versus system parameters
For an M/M/c/BSV retrial queue, the numerical results of E ½Orbit are obtained by considering the following four cases with different values of c.
Case 1 N ¼30, ¼ 5, ¼ 10, p ¼ 0.5, ¼ 5, varying from 10 to 15.
Case 2 N ¼30, ¼ 5, ¼ 10, p ¼ 0.5, ¼ 10, varying from 10 to 15. Case 3 N ¼30, ¼ 15, ¼ 15, p ¼ 0.5, ¼ 10, varying from 5 to 10. Case 4 N ¼30, ¼ 5, ¼ 15, ¼ 15, p ¼ 0.5, varying from 10 to 15.
Results of E ½Orbit are depicted in Figures 1–4 for Cases 1–4, respectively. One sees from Figures 1 and 2 that E ½Orbit drastically decreases as or increases for c ¼1, while E ½Orbit is not sensitive to or for c 2. It reveals from Figure 3 that E ½Orbit increases violently as increases for c ¼ 1, while E ½Orbit slightly increases as increases for c 2. Figure 4 reports that E ½Orbit decreases as increases for c ¼1, while E ½Orbit is not sensitive to for c 2.
Figure 2. The expected number of customers in orbit E[Orbit] versus . Figure 1. The expected number of customers in orbit E[Orbit] versus .
There are several general descriptors of retrial queues, some of which are listed below:
(1) The overall rate of retrials 1 ¼ XN j¼1 jX c k¼0 Xck i¼0 Pki,jþ X 1 j¼Nþ1 NX c k¼0 Xck i¼0 Pki,j¼X N j¼1 jje þ X1 j¼Nþ1 NNRjNe ¼X N j¼1
j&je þ N&NRðI RÞ1e ¼
XN
j¼1
j&jþN&NRðI RÞ1
" # e: ¼&N X N1 j¼1 j(jþ1þNðI RÞ1 " # e ð19Þ
Figure 3. The expected number of customers in orbit E[Orbit] versus .
Figure 4. The expected number of customers in orbit E[Orbit] versus .
(2) The rate of retrials that are successful 2 ¼ XN j¼1 jX c k¼0 X ck1 i¼0 Pki,jþ X 1 j¼Nþ1 NX c k¼0 X ck1 i¼0 Pki,j: ð20Þ
(3) The fraction of retrials that are successful F ¼ 2 1 : ð21Þ
(4) The marginal distribution of the number of busy servers X1
j¼0
Pki,j, 0 i þ k c: ð22Þ
(5) Busy period: The busy period T of a retrial queue is defined as the period that starts at the epoch when an arriving customer finds an empty system (all servers are idle and no customer in the orbit) and ends at the departure epoch at which the system is empty again.
The mean busy period
EðT Þ ¼1 1 P0 0,0 1 ! ¼1 1 &N(1½1 1 ð23Þ
where the symbol ‘&N(1½1’ denotes the first element of the column
vector &N(1.
(6) Vain retrials: A vain retrial is an unsuccessful retrial when all servers are busy.
The steady-state probability of vain retrial PV
PV¼ P1 j¼1 P iþk¼cPki,j P1 j¼1 Pc k¼0 Pck i¼0 Pki,j ¼ P1 j¼1 P iþk¼cPki,j 1 &0e : ð24Þ
4.2. System performance versus truncated parameters
To understand how system performance measures listed above vary with N, we also perform a numerical investigation to the measures based on changing the value of N from 5 to 25, which is based on ¼ 5, ¼ 15, p ¼ 0:5, ¼ 10, and ¼ 10. The numerical illustration is graphically presented in Figures 5–8.
From Figures 5–8, it is clear that increasing the retrial rate beyond a certain point does not result in a commensurate improvement in the system performance, which is according with the result of Neuts and Rao [36].
5. Optimization analysis
In this section, we construct the total expected cost function per unit time based on the system performance measures for the M/M/c/BSV retrial queue, in which the
number of servers (c) is a discrete decision variable, and the service rate () and the vacation rate () are continuous decision variables. Let us define the following cost elements:
Chholding cost per unit time per customer present in orbit;
Cscost per unit time of providing a service rate ;
Cvcost per unit time when one server is on vacation;
Crcost per unit time of providing a vacation rate ; and
Cpfixed cost for purchasing one server: Figure 6. The fraction of successful retrials F versus N.
Figure 5. The expected number of customers in orbit E[Orbit] versus N.
Based on the definition of the cost parameters, the total expected cost function per unit time can be expressed as
Fðc, , Þ ¼ ChE ½Orbit þ Cs þ CvE ½V þ Cr þ Cpc ð25Þ
where E ½Orbit and E ½V are defined previously.
The main objective is to find the optimal number of servers c , and the optimal
values of service rate and vacation rate ð , Þsimultaneously which minimize the
cost function Fðc, , Þ. The analytical study of the optimization behaviour of the
Figure 7. The mean busy period E[T] versus N.
Figure 8. The steady-state probability of vain retrial PVversus N.
expected cost function would have been an arduous task to undertake since the decision variables appear in an expression which is a highly non-linear and complex and non-linear in terms of ðc, , Þ. Next, two methods are provided to deal with this problem heuristically.
In the next section, we first use the quasi-Newton method to find the approximate optimal value of continuous variable ð, Þ, say ð , Þ, and then use
direct search method to search the optimal value of discrete variable c, say c .
5.1. Quasi-Newton method for optimal ðl, gÞ
For practice situation of purchase budget, the number of servers is bounded by a positive integer cU1. We want to find the joint optimal value ( , ) for each given
cin the feasible set {1, 2, . . . , cU}. The cost minimization problem can be illustrated
mathematically as
Fðc, , Þ ¼ min
and s:t: ð9Þ Fðc, , Þ cj
, c ¼1, 2, . . . , cU ð26Þ
For the problem of (26), we should show the convexity of Fðc, , Þ in ð, Þ. However, this study is difficult to implement. It is noted that the derivative of the cost function F with respect to ð, Þ indicates the direction at which the cost function increases. It means that, the optimal value ð , Þcan be found along this
opposite direction of the gradient (Chong and Zak [12]). That is, for a fixed c, quasi-Newton method is employed to search ð, Þ until the approximate minimum value of Fðc, , Þ is achieved, say Fðc, , Þ. An effective procedure that makes it possible
to calculate the optimal value ðc, , Þis presented as follows:
Algorithm Quasi-Newton Method
Step 1 Set the initial trial solution for ~hð0Þ, and compute Fðc, ð0Þ, ð0ÞÞ.
Step 2 Compute the cost gradient ~rFð~hÞ ¼ ½@F=@, @F=@T and the cost Hessian matrix Hð~hÞ ¼ @ 2F=@2 @2F=@@ @2F=@@ @2F=@2 at point ~hðiÞ:
Step 3 While j@F=@j 4 " or j@F=@j 4 ", set the new trial solution ~hðiþ1Þ¼ ~hðiÞ ½Hð~hðiÞÞ1rFð~~ hðiÞÞand return to Step 2.
To demonstrate the validness and the approximate optimization solution, we perform some computation and analysis on the examples given in Table 1 by considering the following cost parameters as
Ch¼ $25=customer=unit time, Cs¼ $45=unit time,
Cv¼ $120=server=unit time, Cr¼ $90=unit time, Cp¼ $120=server
From Table 1, it can be seen that the minimum expected cost per unit time of 1474.377 is achieved at ð , Þ ¼(11.54626, 6.305710) by using six iterations, which
is based on Case (i) with initial value ðc, , Þ ¼ (1, 15, 5). Based on Case (ii) with
Table 1. The illustration of the implement process of quasi-Newton method. Iterations 0 1 2 3 4 5 6 Case (i): ð ,p , Þ¼ (5, 0.5, 10) with initial value ðc , , Þ¼ (1, 15, 5) F ðc , , Þ 1544.435 1517.015 1482.721 1474.921 1474.380 1474.377 1474.377 15 10.74763 11.11560 11.41594 11.53441 11.54617 11.54626 5 5.932174 6.131345 6.263916 6.303111 6.305700 6.305710 @F @ 15.31879 78.2392 25.8695 5.64068 0.43039 0.00300 7 :8 10 8 @F @ 73.2424 133.720 43.6031 9.22994 0.66640 0.00424 1 :5 10 7 E [Orbit] 7.177405 10.75622 8.070767 6.782341 6.418249 6.388411 6.388210 E [V ] 0.500000 0.421422 0.407740 0.399111 0.396630 0.396467 0.396466 Case (ii): ð ,p , Þ¼ (10, 0.8, 15) with initial value ðc , , Þ¼ (2, 10, 10) F ðc , , Þ 2037.910 1988.860 1971.630 1968.793 1968.692 1968.692 1968.692 10 11.05421 11.93856 12.42039 12.52661 12.53093 12.53093 10 9.256253 8.869115 8.722289 8.697166 8.696282 8.696281 @F @ 98.0608 41.9620 13.3913 2.29042 0.09060 0.00016 7 :7 10 9 @F @ 35.0235 22.3227 9.22534 1.86890 0.08050 0.00014 1 :6 10 9 E [Orbit] 9.276428 7.785777 6.717369 6.192268 6.074761 6.069724 6.069715 E [V ] 0.799990 0.862781 0.902006 0.917190 0.919840 0.919933 0.919933
initial value ðc, , Þ ¼ (2, 10, 10), the minimum expected cost per unit time of 1968.692 is achieved at ð , Þ ¼(12.53093, 8.696281) using six iterations.
5.2. Direct search method for optimal c
After obtaining the joint approximate optimal value ð , Þ of the continuous
variable ð, Þ, we use direct search method to obtain the optimal c such that the expected cost function Fðc, , Þ attains a minimum, say Fðc , , Þ. Therefore,
the cost minimization problem can be illustrated mathematically as Fðc , , Þ ¼ min
1ccU
Fðc, , Þ
ð27Þ
The procedure to find the optimal solution is described in the following. A numerical example given in Table 2 is based on (i) ð, p, Þ ¼ (10, 0.8, 15) and (ii) ð, p, Þ ¼ (15, 0.5, 20).
Algorithm Direct Search Method
Step 1 Set F ¼Mwhich M is a sufficiently large number.
Step 2 For each i from 1 to cU, set a initial trial solution ð, Þ and use
Quasi-Newton method to find the optimal value ð , Þ and the cost function
Fðc, , Þ.
Step 3 If the quasi-Newton method diverges, try another initial trial solution and back to Step 1.
Step 4 If Fðc, , Þ5 F , set F ¼Fðc, , Þand S ¼ ðc, , Þ.
It is noted that the optimal value is ðc , , Þ ¼(4, 5.999552, 5.046493) and the
corresponding minimum cost is F ¼1708.284 for Case (i). For Case (ii),
ðc , , Þ ¼(4, 8.099802, 5.265980) and F ¼1819.241 are optimal.
Table 2. The optimal value ð , Þand the corresponding minimum expected cost.
c Initial value Coverage value ð , Þ Iteration Cost*
Case (i) ð, p, Þ ¼ (10, 0.8, 15) 1 [25, 15] [25.13488, 16.43305] 6 3118.635 2 [10, 10] [12.53093, 8.696281] 6 1968.692 3 [10, 5] [8.214208, 6.210196] 6 1725.728 4 [5, 5] [5.999552, 5.046493] 7 1708.284 5 [5, 5] [4.652035, 4.414643] 7 1779.094 Case (ii) ð, p, Þ ¼ (15, 0.5, 20) 1 [30, 20] [33.17698, 17.35916] 6 3601.021 2 [15, 10] [16.60255, 9.183037] 5 2210.467 3 [10, 5] [10.97471, 6.530226] 10 1882.075 4 [6, 6] [8.099802, 5.265980] 8 1819.241 5 [5, 5] [6.347280, 4.561196] 7 1861.652
Table 3. The optimal value ðc , , Þ and the minimum expected cost for various values of and p . ð ,p , Þ (5, 0.2, 10) (10, 0.2, 10) (20, 0.2, 10) (5, 0.8, 10) (10, 0.8, 10) (20, 0.8, 10) c* 2 34 44 5 ð , Þ [4.965695, 2.123714] [6.427349, 2.781059] [9.416220, 3.974561] [2.997995, 2.998664] [6.062298, 5.075460] [9.609657, 7.689420] F ðc , , Þ 901.7296 1245.806 1727.201 1325.523 1716.873 2386.602 E [Orbit] 2.825372 3.199280 4.199710 1.626472 3.125312 4.497047 E [V ] 0.470873 0.719505 1.006400 1.333927 1.576212 2.080781 ð ,p , Þ (5, 0.2, 10) (5, 0.5, 10) (5, 0.8, 10) (10, 0.2, 15) (10, 0.5, 15) (10, 0.8, 15) c* 2 34 33 4 ð , Þ [4.965695, 2.123714] [3.774111, 2.689427] [2.997995, 2.998664] [6.347744,2.767427] [7.295827, 4.645567] [5.999552, 5.046493] F ðc , , Þ 901.7296 1116.483 1325.523 1237.045 1511.634 1708.284 E [Orbit] 2.825372 2.122060 1.626472 3.024207 3.662626 2.955528 E [V ] 0.470873 0.929566 1.333927 0.722693 1.076295 1.585259 ð ,p , Þ (10, 0.2, 5) (10, 0.2, 10) (10, 0.2, 15) (10, 0.8, 5) (10, 0.8, 10) (10, 0.8, 15) c* 2 33 44 4 ð , Þ [10.00245, 3.820378] [6.427349, 2.781059] [6.347744,2.767427] [6.232824, 5.154912] [6.062298, 5.075460] [5.999552, 5.046493] F ðc , , Þ 1361.503 1245.806 1237.045 1739.966 1716.873 1708.284 E [Orbit] 5.789514 3.199280 3.024207 3.572681 3.125312 2.955528 E [V ] 0.5235084 0.719505 0.722693 1.551918 1.576212 1.585259
Finally, we perform a sensitivity investigation on the optimal values ðc , , Þ.
For various values of and p, the minimum expected cost Fðc , , Þ, the system
performance measures Ls and E ½V at the optimum values ðc , , Þare given in
Table 3.
From Table 3, it can be seen that (1) c is insensitive to or p; (2) increases as
increases; and (3) increases as or p increases. Moreover, the minimum expected
cost increases Fðc , , Þas or p increases.
6. Conclusions
An M/M/c retrial queue with Bernoulli vacation (M/M/c/BSV retrial queue) was investigated using the matrix-geometric method. The queueing system was formu-lated as a QBD process. The sufficient and necessary condition for the stability of the system was discussed. The stationary probability vectors were obtained. We also obtained some system performance in matrix forms. A cost model was constructed to calculate the optimal number of servers, the optimal service rate and vacation rate, so that the cost function is minimized. Two methods were provided to deal with the optimization problem heuristically. We performed a sensitivity analysis of the joint optimal values ðc , , Þwith respect to specific values of , p and .
Acknowledgements
The authors gratefully acknowledge the constructive comments of editors and the anonymous reviewers.
References
[1] J. Amador and J.R. Artalejo, On the distribution of the successful and blocked events in the M/M/c retrial queue: A computational approach, Appl. Math. Comput. 190 (2007), pp. 1612–1626.
[2] J.R. Artalejo, A classified bibliography of research on retrial queues: Progress in 1990–1999, Top 7 (1999), pp. 187–211.
[3] J.R. Artalejo, Accessible bibliography on retrial queues: Progress in 2000–2009, Math. Comput. Model. 51 (2010), pp. 1071–1081.
[4] J.R. Artalejo, S.R. Chakravarthy, and M.J. Lopez-Herrero, The busy period and the waiting time analysis of a MAP/M/c queue with finite retrial group, Stoch. Anal. Appl. 25 (2007), pp. 445–469.
[5] J.R. Artalejo, A. Economou, and M.J. Lopez-Herrero, Algorithmic approximations for the busy period distribution of the M/M/c retrial queue, Eur. J. Oper. Res. 176 (2007), pp. 1687–1702.
[6] J.R. Artalejo, A. Economou, and M.J. Lopez-Herrero, Algorithmic analysis of the maximum queue length in a busy period for the M/M/c retrial queue, INFORMS J. Comput. 19 (2007), pp. 121–126.
[7] J.R. Artalejo and M. Pozo, Numerical calculation of the stationary distribution of the main multiserver retrial queue, Ann. Oper. Res. 116 (2002), pp. 41–56.
[8] L. Breuer, A.N. Dudin, and V.I. Klimenok, A retrial BMAP/PH/N system, Queueing Syst. 40 (2002), pp. 433–457.
[9] L. Bright and P.G. Taylor, Calculating the equilibrium distribution in level dependent quasi-birth-and-death processes, Stoch. Models 11 (1995), pp. 497–525.
[10] S.R. Chakravarthy and A.N. Dudin, Multiserver retrial queue with BMAP arriving and group services, Queueing Syst. Theory Appl. 42 (2002), pp. 5–31.
[11] X. Chao and Y. Zhao, Analysis of multiple-server queues with station and server vacations, Eur. J. Oper. Res. 110 (1998), pp. 392–406.
[12] E.K.P. Chong and S.H. Zak, An Introduction to Optimization, 2nd ed., John Wiley and Sons, Inc., New York, NY, 2001.
[13] G. Choudhury, An M[x]/G/1 queueing system with a setup period and a vacation period, Queueing Syst. Theory Appl. 36 (2000), pp. 23–38.
[14] G. Choudhury, A batch arrival queue with a vacation time under single vacation policy, Comput. Oper. Res. 29 (2002), pp. 1941–1955.
[15] G. Choudhury, A two phase batch arrival retrial queueing system with Bernoulli vacation schedule, Appl. Math. Comput. 188 (2007), pp. 1455–1466.
[16] G. Choudhury, Steady state analysis of an M/G/1 queue with linear retrial policy and two phase service under Bernoulli vacation schedule, Appl. Math. Model. 32 (2008), pp. 2480–2489.
[17] G. Choudhury and K. Deka, An M[x]/G/1 unreliable retrial queue with two phases
of service and Bernoulli admission mechanism, Appl. Math. Comput. 215 (2009), pp. 936–949.
[18] G. Choudhury and K.C. Madan, A two phase batch arrival queueing system with a vacation time under Bernoulli schedule, Appl. Math. Comput. 149 (2004), pp. 337–349.
[19] G. Choudhury and K.C. Madan, A two-stage batch arrival queueing system with a modified Bernoulli schedule vacation under N-policy, Math. Comput. Model. 42 (2005), pp. 71–85.
[20] G. Choudhury, L. Tadj, and M. Paul, Steady state analysis of an M[x]
/G/1 queue with two phase service and Bernoulli vacation schedule under multiple vacation policy, Appl. Math. Model. 31 (2007), pp. 1079–1091.
[21] B.T. Doshi, Queueing systems with vacations-a survey, Queueing Syst. 1 (1986), pp. 29–66. [22] G.I. Falin, Calculations of probability characteristics of a multiline system with repeated
calls, Moscow Univ. Comput. Math. Cybernet. 1 (1983), pp. 43–49.
[23] G.I. Falin and J.G.C. Templeton, Retrial Queues, Chapman and Hall, London, 1997. [24] A. Gomez-Corral, A bibliographical guide to the analysis of retrial queues through matrix
analytic techniques, Ann. Oper. Res. 141 (2006), pp. 177–207.
[25] J.C. Ke and F.M. Chang, M[x]/(G1, G2)/1 retrial queue under Bernoulli vacation schedules
with general repeated attempts and starting failures, Appl. Math. Model., 33 (2009), pp. 3186–3196.
[26] J.C. Ke and Y.K. Chu, A modified vacation model M[x]/G/1 system, Appl. Stoch. Bus. Ind. 22 (2006), pp. 1–16.
[27] G.V. Krishna Reddy, R. Nadarajan, and R. Arumuganathan, Analysis of a bulk queue with N-policy multiple vacations and setup times, Comput. Oper. Res. 25 (1998), pp. 957–967.
[28] V.G. Kulkarni and H.M. Liang, Retrial queues revisited, in Frontiers in Queueing: Models and Applications in Science and Engineering, J.H. Dshalalow, ed., CRC Press, Inc., Boca Raton, FL, 1997, pp. 19–34.
[29] C. Langaris and I. Dimitriou, A queueing system with n-phases of service and (n-1)-types of retrial customers, Eur. J. Oper. Res. 205 (2010), pp. 638–649.
[30] G. Latouche and V. Ramaswami, Introduction to Matrix Analytic Methods in Stochastic Modeling (ASA-SIAM Series on Statistics and Applied Probability), Elsevier Science, Amsterdam, Netherlands, 1999.
[31] H.W. Lee, S.S. Lee, J.O. Park, and K.C. Chae, Analysis of M[x]/G/1 queue with N policy and multiple vacations, J. Appl. Prob. 31 (1994), pp. 467–496.
[32] S.S. Lee, H.W. Lee, S.H. Yoon, and K.C. Chae, Batch arrival queue with N policy and single vacation, Comput. Oper. Res., 22 (1995), pp. 173–189.
[33] Y. Levy and Y.U. Yechiali, An M/M/c queue with servers’ vacations, INFOR 14 (1976), pp. 153–163.
[34] K.C. Madan, W. Abu-Dayyeh, and F. Taiyyan, A two server queue with Bernoulli schedules and a single vacation policy, Appl. Math. Comput. 145 (2003), pp. 59–71.
[35] M.F. Neuts, Matrix Geometric Solutions in Stochastic Models: An Algorithmic Approach, The John Hopkins University Press, Baltimore, MD, 1981.
[36] M.F. Neuts and B.M. Rao, Numerical investigation of a multiserver retrial model, Queueing Syst. Theory Appl. 7 (1990), pp. 169–190.
[37] M. Shomrony and U. Yechiali, Burst arrival queues with server vacations and random timers, Math. Meth. Oper. Res. 53 (2001), pp. 117–146.
[38] S.N. Stepanov, Markov models with retrials: The calculation of stationary performance measures based on the concept of truncation, Math. Comput. Model. 30 (1999), pp. 207–228.
[39] L. Tadj, G. Choudhury, and C. Tadj, A quorum queueing system with a random setup time under N-policy and with Bernoulli vacation schedule, Qual. Technol. Quant. Manag. 3 (2006), pp. 145–160.
[40] L. Tadj, G. Choudhury, and C. Tadj, A bulk quorum queueing system with a random setup time under N-policy and with Bernoulli vacation schedule, Stoch.: Int. J. Probab. Stoch. Process. 78 (2006), pp. 1–11.
[41] H. Takagi, Queueing Analysis: A Foundation of Performance Evaluation, Vol. I, Vacation and Priority Systems, Part I, North-Holland, Amsterdam, 1991.
[42] N. Tian, Q. Li, and J. Cao, Conditional stochastic decompositions in the M/M/c queue with server vacations, Stoch. Models, 14 (1999), pp. 367–377.
[43] X. Xu and Z.G. Zhang, Analysis of multiple-server queue with a single vacation (e,d)-policy, Perform. Eval. 63 (2006), pp. 825–838.
[44] T. Yang and J.G.C. Templeton, A survey on retrial queues, Queueing Syst. Theory Appl. 2 (1987), pp. 201–233.
[45] Z.G. Zhang and N. Tian, Analysis of queueing systems with synchronous single vacation for some servers, Queueing Syst. Theory Appl., 45 (2003), pp. 161–175.
[46] Z.G. Zhang and N. Tian, Analysis on queueing systems with synchronous vacations of partial servers, Perform. Eval. 52 (2003), pp. 269–282.
Appendix
For instance, for c ¼ 3, the sub-matrices of B are
b0¼ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 6 6 6 6 4 3 7 7 7 7 5 , b1¼ 0 0 0 0 0 0 0 0 2 6 4 3 7 5, b2¼ 0 0 0 " # , b3¼:
The sub-matrices of C1, C2and C3are
c01¼ 0 0 0 0 0 0 0 0 0 0 0 0 0 2 6 6 6 6 4 3 7 7 7 7 5 , c11¼ 0 0 0 0 0 0 0 2 6 4 3 7 5, c21¼ 0 0 0 " # , c31¼0: c02¼ 0 2 0 0 0 0 2 0 0 0 0 2 0 0 0 0 2 6 6 6 4 3 7 7 7 5, c 1 2¼ 0 2 0 0 0 2 0 0 0 2 6 4 3 7 5, c22¼ 0 2 0 0 , c32¼0: c03¼ 0 3 0 0 0 0 3 0 0 0 0 3 0 0 0 0 2 6 6 6 4 3 7 7 7 5, c 1 3¼ 0 3 0 0 0 3 0 0 0 2 6 4 3 7 5, c23¼ 0 3 0 0 , c33¼0:
The diagonal sub-matrices of Aj, where j ¼ 0, 1, 2, 3 are described as follows. For A0: Y00¼ ð1 pÞ ð þ Þ 2ð1 pÞ ð þ2Þ 3ð1 pÞ ð þ3Þ 2 6 6 6 4 3 7 7 7 5, Y10¼ ð þ Þ ð1 pÞ ð þ þ Þ 2ð1 pÞ ð þ2 þ Þ 2 6 4 3 7 5, Y20¼ ð þ2Þ ð1 pÞ ð þ þ2Þ , Y30¼ ð þ3Þ: For A1: Y01¼ ð þ Þ ð1 pÞ ð þ þ Þ 2ð1 pÞ ð þ2 þ Þ 3ð1 pÞ ð þ3Þ 2 6 6 6 4 3 7 7 7 5, Y11¼ ð þ þ Þ ð1 pÞ ð þ þ þ Þ 2ð1 pÞ ð þ2 þ Þ 2 6 4 3 7 5, Y21¼ ð þ2 þ Þ ð1 pÞ ð þ þ2Þ , Y31¼ ð þ3Þ For A2:
Y02¼ ð þ2Þ ð1 pÞ ð þ þ2Þ 2ð1 pÞ ð þ2 þ 2Þ 3ð1 pÞ ð þ3Þ 2 6 6 6 6 4 3 7 7 7 7 5 , Y11¼ ð þ þ Þ ð1 pÞ ð þ þ þ Þ 2ð1 pÞ ð þ2 þ Þ 2 6 4 3 7 5, Y22¼ ð þ2 þ 2Þ ð1 pÞ ð þ þ2Þ " # , Y32¼ ð þ3Þ: For A3: Y03¼ ð þ3Þ ð1 pÞ ð þ þ3Þ 2ð1 pÞ ð þ2 þ 3Þ 3ð1 pÞ ð þ3Þ 2 6 6 6 4 3 7 7 7 5, Y13¼ ð þ þ3Þ ð1 pÞ ð þ þ þ3Þ 2ð1 pÞ ð þ2 þ Þ 2 6 4 3 7 5, Y23¼ ð þ2 þ 3Þ ð1 pÞ ð þ þ2Þ , Y33¼ ð þ3Þ:
For A0, A1, A2 and A3, the first super-diagonal sub-matrices and the first sub-diagonal sub-matrices are given by
X0¼ 0 p 0 0 0 0 2p 0 0 0 0 3p 2 6 6 4 3 7 7 5, X1¼ 0 p 0 0 0 2p 2 4 3 5, X2¼ 0 p , and Z1¼ 0 0 0 0 0 0 0 0 0 2 6 4 3 7 5, Z2¼ 2 0 0 2 0 0 , Z3¼3 0 , respectively.
![Figure 2. The expected number of customers in orbit E[Orbit] versus . Figure 1. The expected number of customers in orbit E[Orbit] versus .](https://thumb-ap.123doks.com/thumbv2/9libinfo/7635328.136262/12.739.218.523.91.356/figure-expected-customers-figure-expected-number-customers-versus.webp)
![Figure 3. The expected number of customers in orbit E[Orbit] versus .](https://thumb-ap.123doks.com/thumbv2/9libinfo/7635328.136262/13.739.220.521.89.356/figure-the-expected-number-customers-orbit-orbit-versus.webp)
![Figure 5. The expected number of customers in orbit E[Orbit] versus N.](https://thumb-ap.123doks.com/thumbv2/9libinfo/7635328.136262/15.739.220.522.81.354/figure-the-expected-number-customers-orbit-orbit-versus.webp)
![Figure 7. The mean busy period E[T] versus N.](https://thumb-ap.123doks.com/thumbv2/9libinfo/7635328.136262/16.739.218.522.89.362/figure-the-mean-busy-period-e-t-versus.webp)
