∞ queue with batch arrivals and batch-dedicated servers

Operations Research Letters 32 (2004) 431 – 438

Operations Research Letters www.elsevier.com/locate/dsw

Analysis of an M=G=∞ queue with batch arrivals and batch-dedicated servers Re&k G)ull)u Industrial Engineering Department, Middle East Technical University, 06531 Ankara, Turkey Received 30 October 2003; accepted 3 December 2003

Abstract We analyze an M=G=∞ queue with batch arrivals, where jobs belonging to a batch have to be processed by the same server. The number of jobs in the system is characterized as a compound Poisson random variable through a scaling of the original arrival and batch size processes. c 2003 Elsevier B.V. All rights reserved.  Keywords: M=G=∞ queue; Batch arrivals

1. Introduction and motivation In this article we consider an M=G=∞ queueing system with batch arrivals. Our main di7erence from the mainstream research on M=G=∞ queues is that we restrict the arrivals from the same batch to be processed by the same server. Therefore, unlike classical M=G=∞ queues, the number of busy servers does not necessarily correspond to the number of jobs in the system. Fundamental characteristics of our model are as follows. Customers arrivals are governed by a homogeneous Poisson process. Each arrival brings a positive number of jobs to be processed. We assume that the number of jobs (the batch size) arriving with a customer is a random variable, independent of and identically distributed with the batch size associated with another customer. Service times of the jobs may depend on (1) the batch size, and (2) the sequence of

E-mail address: [email protected] (R. G)ull)u).

the jobs in the batch. Each batch is processed on a separate server, but jobs within the same batch have be served by the same server. The following points outline our main motivation for the study of M=G=∞ queues under batch arrivals, batch and sequence dependent processing times, and batch-dedicated servers: 1. M=G=∞ queues are commonly used in the analysis of repairable item inventory systems (see [7]). Consider a system where di7erent batch sizes correspond to di7erent types of repair orders to be performed. That is, the batch size represents the number of distinct activities to be accomplished to complete a repair. Each repair is handled by a separate server, but all the activities of a repair order has to be performed by the same server. Depending on the type of repair (the batch size) the processing time required to complete each activity may change. In this situation, a batch represents the number of work orders to be completed,

c 2003 Elsevier B.V. All rights reserved. 0167-6377/$ - see front matter  doi:10.1016/j.orl.2003.12.002

432

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

and the number of jobs in the system corresponds to the total number of outstanding work orders. 2. M=G=∞ queues are frequently used for obtaining approximate solutions to stochastic systems with &nitely many servers (telecommunication systems, transportation models, manufacturing systems, repairable item inventory systems, networks of queues). For instance, in repairable item inventory systems with batch arrivals the distribution of the number jobs in the system is used to approximate the distribution of outstanding orders under an (S − 1; S) inventory policy (see [4]). In these approximations however, each item arriving within a batch is processed on a separate server [7]. Therefore, these approximations tend to underestimate the workload in the system and may yield serious errors in system performance measures (see [12]). M=G=∞ queues are extensively studied in the literature. In the most elementary case, where arrivals occur according to a Poisson process, the result known as Palm’s theorem characterizes the distribution (either time-dependent or stationary) of the number of customers in the system as a Poisson distribution (whose rate function depends on the service time distribution) (see, for example [8]). When arrivals occur in batches, an extension of Palm’s theorem states that the distribution of the number in the system can be represented as a compound Poisson random variable [9,4,2]. Several authors have investigated more general arrival processes. Shanbhag [11] considers non-homogeneous Poisson arrival streams for the batches. Liu, Kashyap, and Templeton [5] investigate batch arrivals where the inter-arrival times follow a renewal process. In a recent article Masuyama and Takine [6] consider a system with multiple-type Markovian arrival streams and customer-type dependent service times. Under phase-type service time distribution Masuyama and Takine [6] derive formulas for the time-dependent and stationary moments of the number in the system. Best to our knowledge our study is the &rst attempt to analyze a batch arrival M=G=∞ queue where all the jobs within a batch have to be processed by the same server. In this study our major result is to show that under fairly general conditions the distribution of the number of jobs in the system (time dependent and steady state) is obtained as a compound Poisson

distribution through a suitable scaling of the arrival and batch size distributions. We also provide examples and special cases where explicit characterization of the number in the system probabilities is possible. The rest of the article is organized as follows. In Section 2 we introduce the notation and essential variables that will be used throughout the paper. In Section 3 we derive the distribution of the number of jobs in the system for transient (Theorem 1) and steady-state (Corollary 1) cases. We re&ne our results by assuming stationary service times and present them in Section 4. Our conclusions and directions for future results are presented in Section 5.

2. Notation and the basic model We consider the arrival and service processes. Customer arrivals occur according to a homogeneous Poisson process with rate . Let N (t) be the number of customer arrivals by time t ¿ 0: For k = 0; 1; : : : Pr{N (t) = k} =

e−t (t)k : k!

Each arrival brings to the system a batch of random size. We assume that the batch sizes associated with arrivals are independent and identically distributed. Let M be the random variable denoting the size of a batch with k = Pr{M = k}, k = 1; 2; : : : . Let S1 ; S2 ; : : : be the arrival times for successive customers, and let M1 ; M2 ; : : : be their associated batch sizes (each having identical distribution with M ). Let {Yj; m ; m=1; 2; : : : ; j=1; 2; : : : ; m} denote the set of random variables with respective distribution functions {Gj; m (y); m = 1; 2; : : : ; j = 1; 2; : : : ; m}. Given that a batch of size m is realized in an arrival instance, the service times of the jobs are characterized by the distribution functions {Gj; m (y); j=1; 2; : : : ; m}. Therefore, Yj; m is the random variable denoting the processing time of jth job in a batch of size m. By convention, job j is the jth job to enter the service. We assume independent service times for the jobs within a batch and across batches. Let Lj; m = E[Yj; m ]. For the ease of exposition we assume that Yj; m ’s are continuous. Our derivations can be replicated for the discrete case. We also assume that there is a constant Lmax such that Lj; m 6 Lmax for all m = 1; 2; : : : and j = 1; 2; : : : ; m.

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

Let Gm(k) (y) be the distribution function of Y1; m + Y2; m + · · · + Yk; m (total processing times of the &rst k jobs within a batch of size m) for m = 1; 2; : : : , and (k) k = 1; 2; : : : ; m. Also let GL m (y) = 1 − Gm(k) (y). By (0) convention, GL m (y) ≡ 0. Let W (t; ; m) be the number of jobs in the system at time t for a batch whose arrival time is  ( 6 t) and with batch size m. Let {Wi (t; ; m); i = 1; 2; : : :} be independent and identically distributed copies of W (t; ; m). Assuming that the system is empty and idle at time t = 0, the total number of jobs in the system at time t ¿ 0, X (t), is de&ned as: X (t) =

N (t)


Wi (t; Si ; Mi ):

i=1

We also de&ne X := limt→∞ X (t). Similarly, let NB(t) be the number of busy servers at time t: NB(t) =

N (t)


1{Wi (t;Si ;Mi )¿0} ;

i=1

and NB = limt→∞ NB(t), where 1{A} is the indicator function for the event A.

433

Proof. Let PX (z; t) = E[z X (t) ]. Then, PX (z; t) = E[z X (t) ] = E[E[z X (t) |N (t)]] =

∞
e−t (t)n

n!

n=0

A(z; t; n);

(1)

where A(z; t; n) = E[z X (t) |N (t) = n] = E[z

n

i=1

Wi (t;Si ;Mi )

]:

Note that, conditioned on N (t) = n, (S1 ; S2 ; : : : ; Sn ) are distributed as ordered uniform random variables over (0; t). Let U1 ; U2 ; : : : ; Un be independent uniform random variables over (0; t). Then, A(z; t; n) = E[z

n

i=1

Wi (t;Ui ;Mi )

]:

Since batch sizes and service times associated with customers are identically distributed and independent, W1 ; : : : ; Wn is an independent and identically distributed collection of random variables. Let W (t; U; M ) be the generic random variable with the same distribution as any Wi , i = 1; 2; : : : ; n. Then, A(z; t; n) = (B(z; t))n ;

3. Distribution of the number of jobs in the system In this section we derive the probability distribution of the number of jobs in the system, and in particular show that its distribution is compound Poisson under a suitable scaling of the original arrival process and batch size distribution. Theorem 1. X (t) is a compound Poisson random variable with rate v(t) and compounding random variable M˜ t , where  t
∞ v(t) = GL (m) m (y)m dy;

where B(z; t) = E[z W (t; U; M ) ]. Plugging A(z; t; n) back in (1), we obtain the form of PX (z; t) as PX (z; t) = e−t(1−B(z; t)) : For m = 1; 2; : : : , we de&ne another auxiliary variable B(z; t; m): B(z; t; m) = E[z W (t; U; m) ]  1 t = E[z W (t; ; m) ] d; t 0 and rewrite B(z; t) as

0 m=1

∞


B(z; t; m)m :

and the distribution of M˜ t is given by

B(z; t) =

k (t) := Pr{M˜ t = k}  t ∞ (m−k) L (m−k+1) (y)−GL m (y))m dy m=k (G m 0 = ; v(t)

Next we derive the distribution of W (t; ; m). For k = 0; 1; : : : ; m

for k = 1; 2; : : : .

m=1

{W (t; ; m) 6 k} if and only if  +

m−k
j=1

Yj; m ¡ t; (2)

434

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

with the convention that that

0 1

=0. Using (2) it follows

= −v(t) +

Pr{W (t; ; m) = 0} = Gm(m) (t − );

− ) −

Gm(m−k+1) (t

− ) where

Pr{W (t; ; m) = m} = 1 − Gm(1) (t − ): Consequently, we can obtain E[z W (t; ; m) ] as z k Pr{W (t; ; m) = k}

k=0

= Gm(m) (t − ) +

m−1


(m−k+1) z k (GL m (t − )

k=1 (m−k) (t − )) + z m GL (1) −GL m m (t − ): (3)

De&ne C(z; t) = −t(1 − B(z; t)), so that PX (z; t) = eC(z; t) . Using (3) we obtain  t
 t
∞ ∞ C(z; t) = − m d + Gm(m) (t − )m d 0 m=1

+

0 m=1

 t
∞ m−1
0 m=1 k=1

(m−k+1) z k (GL m (t − )

(m−k) (t − ))m d −GL m

+

=−

+

 t
∞ 0 m=1

 t
∞ 0 m=1

z m GL (1) m (t − )m d GL (m) m (y)m dy

 t
∞
m 0 m=1 k=1

(m−k+1) z k (GL m (y)

(m−k) (y))m dy −GL m

(m−k+1) (GL m (y)

z k k (t);

k=1

v(t) =

E[z W (t; ; m) ] =

∞
m=k

∞


= −v(t) + v(t)

for k = 1; : : : ; m − 1

m


0 k=1

zk

(m−k) (y))m dy −GL m

Pr{W (t; ; m) = k} =Gm(m−k) (t

 t
∞

 t
∞ 0 m=1

and

GL (m) m (y)m dy;

 t ∞

k (t) :=

0

L (m−k+1) (y) m=k (G m

(m−k) − GL m (y))m dy : v(t)

∞ Note that k (t) ¿ 0 and k=1 k (t)=1. We let M˜ t be the random variable whose probability mass function is given by {k (t); k = 1; 2; : : :}. Also let H (z; t) as the z-transform for M˜ t . We obtain PX (z; t) as PX (z; t) = e−v(t)(1−H (z; t)) ;

(4)

which is the z-transform of a compound Poisson distribution, where the associated Poisson random variable has mean v(t), and M˜ t is the compounding random variable. This observation completes the proof. We can easily obtain E[X (t)] as  dPX (z; t)  E[X (t)] = dz  z=1

 dC(z; t)  = dz  =

z=1

 t
∞
m 0 m=1 k=1

eC(z; 1)

(k) (y)m dy: GL m

(5)

The probability mass function of X (t) can be obtained by utilizing (4). Let M˜ t(k) be the k-fold convolution of M˜ t with itself. Notice that M˜ t(k) ¿ k for all k =1; 2; : : : . Then, for i = 0; 1; : : : Pr{X (t) = i} = e−v(t)

i
(v(t))k k=0

k!

Pr{M˜ t(k) = i}; (6)

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

where M˜ (0) ≡ 0. A numerically more attractive ext pression for (6) can be obtained using properties of compound Poisson random variables ([10], Corollary 2.5.4). For i = 1; 2; : : : Pr{X (t) = i} =

i v(t)


i

k Pr{X (t) = i − k}

k=1

for n = 1; 2; : : : . By following similar steps as in the proof of Theorem 1 we can easily show that the z-transform of NB(t) can be obtained as :

Therefore, NB(t) has a non-homogeneous Poisson distribution with rate v(t). It is interesting to note that, using (5)  t
∞ m−1
(k) E[X (t)] =  (y)m dy + E[NB(t)] GL m 0 m=1 k=1

and E[X (t)] = E[NB(t)] whenever M ≡ 1. The z-transform for the stationary number of jobs in the system, X = limt→∞ X (t) can be obtained as PX (z)=limt→∞ PX (z; t)=limt→∞ eC(z; t) . We outline the respective result in the following corollary. Corollary 1. X is a compound Poisson random variable with rate v and compounding random variable M˜ , where ∞
m
v= Lk; m m ; (7) m=1 k=1

and the distribution of M˜ is given by ∞ Lm−k+1; m m ˜ k = Pr{M = k} = m=k ; v for k = 1; 2; : : : .

(8)

Proof. The result is obtain by taking limits of v(t) and k (t): v := lim v(t) t→∞

∞

0

∞
m


GL (m) m (y) dym

Lk; m m ;

m=1 k=1

where Pr{X (t) = 0} = e−v(t) . Similarly, moments of X (t) can be obtained as ([10], Corollary 2.5.3):   n−1
n−1 n E[X (t)] = v(t) E[X k (t)]E[M˜ tn−k ] k k=0

PNB (z; t) = e

∞ 
m=1

=

×Pr{M˜ t = k}

−v(t)(1−z)

=

435

and similarly for k = 1; 2; : : : ∞ Lm−k+1; m m k := lim k (t) = m=k : t→∞ v The limits above are justi&ed by the Dominated Convergence Theorem (see, for example, [3]) since expected times are bounded. Also note that  t (m) processing  L m (y) dy m Lk; m as t → ∞. We let M˜ be G k=1 0 the random variable whose probability mass function is given by {k ; k = 1; 2; : : :}. Also let H (z) be the z-transform for M˜ . Therefore PX (z) is obtained as PX (z) = e−v(1−H (z)) ; which is the z-transform of a compound Poisson distribution, where the associated Poisson random variable has mean v, and M˜ is the compounding random variable. This observation completes the proof. We note that v is the expected total processing time of a batch. In the following example we demonstrate the e7ect of scaling the original process when the batch size is constant. Example. Let M = K ¿ 1 be constant (so that K = 1 and m = 0 for all m = K). Batches arrive according to a Poisson process with rate , and each batch has exactly K jobs. Let L1 ; L2 ; : : : ; LK be the associated mean processing times for the jobs in a batch. Then, K by using (7) and (8), v = m=1 Lm and LK−k+1 k = L1 + · · · + L K for k = 1; 2; : : : ; K. Even though the original batch size is constant, in the scaled process the resulting compounding random variable M˜ takes values in {1; 2; : : : ; K} with probabilities proportional to mean processing times. 4. Identically distributed processing times Considerable simpli&cation is obtained when the processing times are not a7ected by the batch size and

436

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

the sequence of a job in the batch (that is, for Yj; m = Y for all m ¿ 1 and j = 1; 2; : : : ; m). In this case, Lj; m = L and Gm(k) (y) = G (k) (y) for all m ¿ 1, k = 0; 1; : : : , and j = 1; 2; : : : ; m. Some of the identities of the previous section become: v=L

∞


mm = LE[M ];

m=1

k =

Pr{M ¿ k} ; E[M ]

H (z) =

∞


k = 1; 2; : : : ;

z k Pr{M ¿ k}=E[M ];

k=1

PX (z) = e−LE[M ](1−H (z)) : In particular, the &rst two moments of X can be obtained in terms of the moments of M as

and Var(X ) = E[X 2 ] − (E[X ])2 = LE[M ]E[M˜ 2 ] L {2E[M 3 ] + 3E[M 2 ] + E[M ]}: 6 One can deduce that, E[X n ], for n ¿ 1 can be obtained in terms of E[M k ], for k = 1; 2; : : : ; n + 1. We should note that in an ordinary M=G=∞ system, where each job is processed on a separate server, the nth moment of X can be obtained in terms of the &rst n moments of the batch size random variable (for instance, the mean number of jobs in the system is LE[M ]), whereas here one also needs to use E[M n+1 ]. =

4.1. Geometrically distributed batch size

E[X 2 ] = LE[M ]{E[M˜ 2 ] + E[X ]E[M˜ ]};

Let k = p(1 − p)k−1 for k = 1; 2; : : : be the probability mass function for the size of a batch. We can obtain k (t) as

where

k (t) =

E[X ] = LE[M ]E[M˜ ];

∞

1
E[M˜ ] = k Pr{M ¿ k} E[M ]

=

k=1 ∞

=

j=1



1 1 E[M 2 ] + E[M ] ; 2 2

=

=

and E[M˜ 2 ] =

1 E[M ]

∞


j(j + 1)(2j + 1) 1
Pr{M = j} E[M ] 6 j=1

=

1 E[M ]

1 1 1 E[M 3 ] + E[M 2 ] + E[M ] : 3 2 6

Therefore, E[X ] =

0

0

(GL (m−k+1) (y) − GL (m−k) (y))p(1 − p)m−1 dy

m=k

j=0

v(t) (GL ( j+1) (y) − GL ( j) (y))p(1 − p)k+j−1 dy v(t)

p(1 − p)


∞ k−1 t

p(1 − p)


∞ k−1 t

0

j=0

(GL ( j+1) (y) − GL ( j) (y))(1 − p) j dy v(t)

0

j=1

GL ( j) (y)p(1 − p) j−1 dy

v(t) k−1

;

since v(t) = 0 m=1 GL (m) (y)p(1 − p)m−1 dy. Therefore, M˜ t also has geometric distribution (independent of t), and therefore M˜ t(k) has negative binomial distribution:   i−1 (k) ˜ Pr{M t = i} = pk (1 − p)i−k k −1 for i = k; k + 1; : : : . As a result PX (z; t) is obtained as

L {E[M 2 ] + E[M ]} 2

L = {Var(M ) + E[M ](E[M ] + 1)} 2

v(t)


t ∞

 t ∞

k=1



(GL (m−k+1) (y) − GL (m−k) (y))m dy

m=k

= p(1 − p)

k 2 Pr{M ¿ k}

∞

=

0


t ∞

j(j + 1) 1
Pr{M = j} = E[M ] 2 1 = E[M ]


t ∞

PX (z; t) = e−v(t)(1−(pz=1−(1−p)z)) : (9)

Hence X (t) is a non-homogeneous compound Poisson process with arrival function v(t) and geometric compounding distribution. We can easily obtain

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

the time-dependent probabilities for the number in the system: Pr{X (t) = 0} = e−v(t) ; Pr{X (t) = i} −v(t)

=e

i
(v(t))k k=1



k!

k −1

pk (1 − p)i−k

for i = 1; 2; : : : . On the other hand, as t → ∞, v(t) → LE[M ] = L=p; and consequently, Pr{X = 0} = e Pr{X = i} = e

−L=p

−L=p

i
(L)k

k!

k=1



For the special case outlined above, if one uses a standard M=G=∞ system with batch arrivals (where each job immediately enters service without waiting), the expected number of jobs in the system, E[Xs (t)], can be obtained as



i−1

i−1 k −1

 (1 − p)i−k

437

E[Xs (t)] =

L (1 − e−t=L ): p

Therefore, the standard model underestimates the mean number of jobs in the system by a percentage: E[X (t)] − E[Xs (t)] × 100 E[X (t)] 

1 − e−t=L × 100: = 1−p 1 − e−tp=L

(10)

for i = 1; 2; : : : . An explicit characterization of v(t) is generally dif&cult. In the following example we derive v(t) for the special case of exponential processing times.

As t → ∞, (10) converges to (1 − p)100. In general, by using (9) we can show that for any service time or batch size distribution chosen the stationary percent underestimation is found as

Example (exponential processing times): If we assume that the processing times are exponentially distributed with mean L, then for m = 1; 2; : : : m−1
e−y=L (y=L)k : GL (m) (y) = k!

E[M 2 ] − E[M ] × 100: E[M 2 ] + E[M ]

k=0

Therefore,  t
∞ v(t) = GL (m) (y)p(1 − p)m−1 dy 0 m=1

=

=

=

 t
∞ m−1
e−y=L (y=L)k k!

0 m=1 k=0

 t
∞
∞ 0 k=0 m=k+1

 t
∞ 0 k=0



t

p(1 − p)m−1 dy

p(1 − p)m−1

(1 − p)k

e−y=L (y=L)k dy k!

e−y=L (y=L)k dy k!

L (1 − e−tp=L ): p 0 By following similar steps one can also show that L E[X (t)] = 2 (1 − e−tp=L ): p =

e−yp=L dy =

5. Conclusions This work can be extended in several directions. One possible future research is to consider more general arrival processes, such as Markov modulated arrival streams with di7erent customer classes as in [6]. We can easily impose correlations among processing times of the jobs within a batch. Another possible research direction is to analyze systems with &nitely many batch-dedicated servers. M=G=∞ queues are commonly used in the analysis of inventory systems. Consider an inventory system operating with (S − 1; S) policy, where upon each demand arrival a new order is placed whose amount equals that demand size. Suppose a customer arrives at time t with a demand size of M . Let I (t) be the on-hand inventory at the system at time t. Then the waiting time of the arriving customer, W (t) is W (t)=0 if I (t) ¿ M and W (t) ¿ 0 if I (t) ¡ M (see [1] for a similar analysis under a standard M=G=∞ queue). One interesting extension would be the analysis of the

438

R. G(ull(u / Operations Research Letters 32 (2004) 431 – 438

customer waiting time under batch-dedicated server restriction of our model. References [1] M. Berg, M.J.M. Posner, Customer delay in M=G=∞ repair systems with spares, Oper. Res. 38 (1990) 344–348. [2] M.L. Chaudry, J.G.C. Templeton, A First Course in Bulk Queues, Wiley, New York, 1983. [3] R. Durrett, Probability: Theory and Examples, Wadsworth & Brooks/Cole, Paci&c Grove, 1991. [4] G.J. Feeney, C.C. Sherbrooke, The (s − 1; s) inventory policy under compound Poisson demand, Manage. Sci. 12 (1966) 391–411. [5] L. Liu, B.R.K. Kashyap, J.G.C. Templeton, On the GI X =G=∞ system, J. Appl. Probab. 27 (1990) 671–683.

[6] H. Masuyama, T. Takine, Analysis of an in&nite-server queue with batch Markovian arrival streams, Queueing Systems 42 (2002) 269–296. [7] S. Nahmias, Managing repairable item inventory systems: a review, TIMS Stud. Manage. Sci. 16 (1981) 253–277. [8] N.U. Prabhu, Queues and Inventories, Wiley, New York, 1965. [9] J.F. Reynolds, Some results for the bulk-arrival in&nite-server Poisson queue, Oper. Res. 16 (1968) 186–189. [10] S.M. Ross, Stochastic Processes, Wiley, New York, 1996. [11] D.N. Shanbhag, On in&nite server queues with batch arrivals, J. Appl. Probab. 3 (1966) 274–279. [12] A. Sleptchenko, M.C. van der Heijden, A. van Harten, E7ects of &nite repair capacity in multi-echelon, multi-indenture service part supply systems, Internat. J. Product. Econom. 79 (2002) 209–230.