A single-server queue with batch arrivals and semi-Markov services

Queueing Systems, May 2017

We investigate the transient and stationary queue length distributions of a class of service systems with correlated service times. The classical \(M^X/G/1\) queue with semi-Markov service times is the most prominent example in this class and serves as a vehicle to display our results. The sequence of service times is governed by a modulating process J(t). The state of \(J(\cdot )\) at a service initiation time determines the joint distribution of the subsequent service duration and the state of \(J(\cdot )\) at the next service initiation. Several earlier works have imposed technical conditions, on the zeros of a matrix determinant arising in the analysis, that are required in the computation of the stationary queue length probabilities. The imposed conditions in several of these articles are difficult or impossible to verify. Without such assumptions, we determine both the transient and the steady-state joint distribution of the number of customers immediately after a departure and the state of the process J(t) at the start of the next service. We numerically investigate how the mean queue length is affected by variability in the number of customers that arrive during a single service time. Our main observations here are that increasing variability may reduce the mean queue length, and that the Markovian dependence of service times can lead to large queue lengths, even if the system is not in heavy traffic.

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:

https://link.springer.com/content/pdf/10.1007%2Fs11134-017-9531-4.pdf

A single-server queue with batch arrivals and semi-Markov services

A single-server queue with batch arrivals and semi-Markov services Abhishek 0 1 Marko A. A. Boon 0 1 Onno J. Boxma 0 1 Rudesindo Núñez-Queija 0 1 B Abhishek 0 1 0 Department of Mathematics and Computer Science, Eindhoven University of Technology , Eindhoven , The Netherlands 1 Korteweg-de Vries Institute for Mathematics, University of Amsterdam , Amsterdam , The Netherlands We investigate the transient and stationary queue length distributions of a class of service systems with correlated service times. The classical M X /G/1 queue with semi-Markov service times is the most prominent example in this class and serves as a vehicle to display our results. The sequence of service times is governed by a modulating process J (t ). The state of J (·) at a service initiation time determines the joint distribution of the subsequent service duration and the state of J (·) at the next service initiation. Several earlier works have imposed technical conditions, on the zeros of a matrix determinant arising in the analysis, that are required in the computation of the stationary queue length probabilities. The imposed conditions in several of these articles are difficult or impossible to verify. Without such assumptions, we determine both the transient and the steady-state joint distribution of the number of customers immediately after a departure and the state of the process J (t ) at the start of the next service. We numerically investigate how the mean queue length is affected by variability in the number of customers that arrive during a single service time. Our main observations here are that increasing variability may reduce the mean - queue length, and that the Markovian dependence of service times can lead to large queue lengths, even if the system is not in heavy traffic. Mathematics Subject Classification 60K25 · 90B22 1 Introduction Service systems with correlated service durations have a long tradition in the queueing literature. Such systems enjoy a large variety of application domains, including logistics, production management and telecommunications [2,11,14,19]. Our main motivation stems from road traffic analysis, where traffic flows may interact at junctions or crossings [1,18]. Focus, for illustration, on a traffic flow that merges into a main flow (very similar considerations are valid for road intersections). If the traffic density on the main flow is high, vehicles in the secondary flow may queue up before merging into the main flow. The merging times required for two subsequent vehicles will be strongly correlated as they experience similar traffic conditions on the main flow. In this paper, we will capture this dependence in a queueing model in which the sequence of service times is governed by a modulating Markovian process. Although our analysis allows for a slightly larger class of models, we will use the classical M/G/1 queue with semi-Markov service times [14], and more specifically its extension to batch arrivals [15] to compare our results with existing literature. The first to have investigated this class of queueing models was Gaver [11], who derived the waiting time in a single-server queue with two types of customers arriving according to independent Poisson processes. In that model, service times are classspecific and when service switches from one type to the other, an additional switch-over time is required. This framework was generalized by Neuts [14], allowing for more than two customer types and the sequence of service times forming a semi-Markov process. Under technical assumptions (these will be discussed later in detail), Neuts obtained the transient and stationary distributions of queue lengths, waiting times and busy periods. Subsequently, Çinlar [4] obtained the transient and stationary queue length distributions under less restrictive assumptions, and Purdue [17] showed that the assumptions imposed by Neuts and Çinlar are not necessary for the analysis of the busy period, presenting an alternative approach. The literature on extensions of this model steadily expanded in the next two decades. In [16], Neuts studied the multitype M/G/1 queue with change-over times when switching service from one type of customer to another. A further generalization allowing for Poisson arrivals of groups (batches) of customers of arbitrary random size was investigated by Neuts in [15], obtaining the busy period, queue length and waiting time distributions. The departure process of a related model with single Poisson arrivals and exponential service times was determined by Magalhães and Disney [13]. In that model, the rate of the exponential service times depends on the type of the customer being served as well as that of its predecessor. Models with single arrivals, but with both the arrivals and the services depending on a common semi-Markov process have been investigated by De Smit [7] and Adan and Kulkarni [2]. Using the Wiener–Hopf factorization technique, De Smit [7] obtained the waiting time and queue length distributions. Adan and Kulkarni [2] considered a similar setting, but with the customer type being determined at arrival instants (independent of the service durations). In this paper, we investigate the transient and stationary queue length distributions in a single-server model with semi-Markov service times and with batch arrivals (our framework includes Poisson arrivals of batches as the most prominent example). In order to explain the technical contribution of our work, it is best to compare with the expositions of Neuts [14] and Çinlar [4]. In those papers only single Poisson arrivals were allowed, but the subsequent analysis is very similar. The earlier mentioned technical assumptions made by Neuts entail that the zeros of a particular matrix determinant appearing in the transient analysis are either strictly separated or completely coincide. This ensures that the zeros are analytic functions of the entries of the matrix and, consequently, that the stationary distribution can be obtained from the transient distribution. The assumptions were relaxed by Çinlar [4] while maintaining the analyticity of the zeros. Unfortunately, it remains hard, if not impossible, to verify the required conditions in practice, as they must hold for the zeros as functions of the matrix entries. As noted earlier, Purdue [17] showed that the assumptions imposed by Neuts and Çinlar are not necessary for the analysis of the busy period. Our work show that these assumptions are not needed for the analysis of the queue length distribution either. This comes at the expense of a separate analysis for the stationary distribution, which is more involved than that of the transient distribution. Specifically, we determine the generating function of the number of customers immediately after the departure of an arbitrary customer, considering both transient and steady-state behavior. For Poisson batch arrivals, in steady state we further obtain the queue length distribution at batch arrival instants and at arbitrary times, which are identical due to PASTA. Note that this distribution is in general not the same as that at departure times (for single arrivals, they would coincide). A further contribution is an extensive numerical investigation of the mean queue length in steady state. We show that due to the dependence between service times, the mean number of customers may be very large, even if the load on the system is not large. A noteworthy observation is that increasing the variability in the number of customers arriving during a service time may in fact decrease the mean queue length. The remainder of this paper is organized as follows. Section 2 gives the model description in two layers. First we describe the M X /G/1 queue with semi-Markov services and then present a somewhat more general framework. In Sect. 3, we derive the transient and stationary probability generating functions of the number of customers in the system immediately after a departure. In Sect. 4, we derive the generating functions of the stationary number of customers at an arbitrary epoch, at batch arrival epochs and at customer arrivals. The special case with only two customer types is specified in Sect. 5. Finally, in Sect. 6, we present numerical examples to demonstrate the impact of the correlated arrivals, and of the variability of the number of customers arriving during a service time, on the expected number of customers in the system. 2 Model description We start by describing the M X /G/1 queueing model with semi-Markov services, which is the most natural example in our framework. Our analysis extends directly to any model that satisfies the dynamics described in the recurrence relation (2.9) below. 2.1 The M X / G/1 queue with semi-Markov service times Customers arrive in batches at a single-server queue according to a Poisson process with rate λ; the batch size is denoted by the random variable B with generating function B(z), for |z| ≤ 1. Customers are served in order of arrival, with speed 1. Customers within a batch are assumed to be ordered arbitrarily. The service times are governed by a Markov process Jn, n = 0, 1, . . . , that can take values in {1, 2, . . . , N }, for some integer N . It will be convenient to refer to Jn as the type of the nth customer; thus, there are N customer types. The service time of the nth customer is denoted with G(n). An essential feature of our model is that the type of the (n + 1)th customer depends both on the type of the nth customer and on the service duration of the nth customer. This exactly matches the framework of semi-Markov service times introduced by Neuts [14]. We define For future use, we introduce the Laplace–Stieltjes transform (LST) where 1{.} denotes the indicator function. In particular, The type of a customer, and its service time, do not depend on the arrival process. It should be observed that { Jn, n = 1, 2, . . . } forms a finite-state Markov chain. We shall restrict ourselves to irreducible Markov chains. The stationary distribution P( J = j ) of the Markov chain Jn is given by the unique solution of the set of equations with normalizing condition Nj=1 P( J = j ) = 1. The mean service time of an arbitrary customer is given by The stability condition for this model is given by This can be formalized using Theorem 3 from Loynes [12], by describing the workload process in terms of “super customers” whose service times are the aggregate service times of customers in a single batch. Let G(m) be the service time of the super customer corresponding to the mth arriving batch, and Jm the type of the first customer in the mth batch. Starting from a stationary version of the sequence (G(n), Jn+1), one can readily construct a stationary sequence (G(m), Jm+1) for the super customers. Note that by construction G(m) is also stationary and, together with the arrival epochs of batches (which form an independent Poisson process), this sequence completely determines the workload process. This description of the workload process satisfies the criteria to use the characterization for stability in Loynes [12]. We will investigate the queue length process at departure times of customers. For that it will be convenient to define An as the number of customers arriving during the service time of the nth customer and Bn as the size of the batch in which the nth customer arrived. Note that for i, j = 1, 2, . . . , N , |z| ≤ 1, Ai j (z) := E[z An 1{Jn+1= j}| Jn = i ] = G˜ i j (λ(1 − B(z))). The queue length distribution at customer departure times is fully determined by the sequences An and Bn. For the analysis, it is not needed that the arrivals during service times occur in batches at Poisson instants. For that reason, we will now formulate our general model in terms of the An and Bn only; to specify our later results for the M X /G/1 queue with semi-Markov services, we will simply substitute the relation given in (2.7). 2.2 General model The inputs to our general model are probability generating functions of non-negative discrete random variables Ai j (z), i, j ∈ {1, 2, . . . , N }, and B(z). From the Ai j (z), we construct a Markov process ( An, Jn+1), n = 1, 2, . . . , satisfying E[z An 1{Jn+1= j}| Jn = i ] = Ai j (z). In this construction, it is implicit that ( An, Jn+1) conditional on Jn is independent of An−1. The sequence Bn is i.i.d. with generating function B(z) and independent of the sequence An. Next we define the recurrence relation Xn = Note: If the Ai j (z) are set equal to (2.7), then the sequence Xn follows the same law as the number of customers at departure times in the M X /G/1 queue with semiMarkov services. The role of the Bn is subtle in this representation: Bn is only included if the (n − 1)th customer leaves the system empty upon departure. The nth customer is therefore the first customer in a batch that arrives into an empty system. Only for that reason, the sequence Bn can be taken independent of the An in the M X /G/1 queue with semi-Markov services. In the sequel, we will study the transient and stationary distributions of Xn defined by (2.9). Again using Theorem 3 of Loynes [12], we may conclude that the stability condition in this case is Here E[ A] denotes the expectation of the An in stationarity: E[ A] = i=1 j=1 αi j = E[ An1{Jn+1= j}| Jn = i ] = Ai j (1). Note that at first sight (2.9) does not seem to fit the framework in Loynes [12], because of the special condition when the system is empty. For stability, however, the behavior of an empty system is irrelevant. 3 The queue length distribution at departure epochs We shall determine the transient and steady-state joint distribution of the number of customers immediately after a departure, and the type of the next customer to be served. From the recurrence relation (2.9), we find, for the probability generating functions, i=1 1 N i=1 + E z An+Bn−11{Jn+1= j}1{Xn−1=0} + E z An+Bn−11{Jn+1= j}1{Xn−1=0} E z Xn−1−1+An 1{Jn+1= j}| Jn = i P( Jn = i ) E z An 1{Jn+1= j}1{Xn−1=0}| Jn = i P( Jn = i ) E z An+Bn−11{Jn+1= j}1{Xn−1=0}| Jn = i P( Jn = i ), Now we exploit the fact that Xn−1 and ( An, Jn+1) are conditionally independent given Jn, and the Bn are also independent of all other random variables: E z Xn 1{Jn+1= j} = 1 N = z i=1 i=1 i=1 E z Xn−1 1{Jn=i} E z An 1{Jn+1= j}| Jn = i 3.1 Steady-state analysis In this subsection, we restrict ourselves to the steady-state queue length distribution, assuming that the stability condition (2.10) holds. In the next subsection, we will analyze the transient behavior of the queue length. It will be useful to introduce some further notation: for i = 1, 2, . . . , N , Ai (z) = j=1 j=1 f j (z) = limn→∞E z Xn 1{Jn+1= j} , f j (0) = limn→∞P(Xn = 0, Jn+1 = j ), f j (1) = limn→∞P( Jn+1 = j ) = P( J = j ). The probability generating function of the steady-state queue length distribution immediately after a departure is denoted by F (z) = j=1 M (z)T f (z) = b(z), In steady state, Eq. (3.1) leads to the following N equations: We can also write these N linear equations in matrix form as i=1,i = j i=1 (z − A j j (z)) f j (z) − Ai j (z) fi (z) = (B(z) − 1) Ai j (z) fi (0), j = 1, 2, . . . , N. (3.8) M (z) = ⎢⎢ ⎣ f (z) = ⎢⎢ ⎣ Therefore, solutions of the non-homogeneous linear system M (z)T f (z) = b(z) are of the form cof M (z)T T b(z), provided det M (z) = 0. Here cof M (z)T is the cofactor matrix of M (z)T . It remains to find the values of f1(0), f2(0), . . . , f N (0). We shall derive N linear equations for f1(0), f2(0), . . . , f N (0). 1 1 limz→1 z − 1 eˆ M (z)T f (z) = limz→1 z − 1 eˆb(z), where eˆ is a row vector with all entries one. After simplification, we can write this as z − 1 = limz→1 B(z) − 1 z − 1 j=1 i=1 Ai j (z) fi (0). i=1 (N-1) remaining equations: To find the remaining N − 1 equations, we first prove that det M (z) has exactly N − 1 zeros in |z| < 1 and the zero z = 1 on |z| = 1. Since fi (z) is an analytic function in |z| < 1, the numerator of fi (z) also has N − 1 zeros in the unit disk |z| < 1. As a consequence, these N − 1 zeros provide N − 1 linear equations for f1(0), f2(0), . . . , f N (0). To find the N − 1 zeros, we use a method that has also been applied in [2,6,9]. It is based on the concept of (strict) diagonal dominance in a matrix. The proof consists of four steps: Step 1: Prove that each element on the diagonal of M (z) has exactly one zero in |z| < 1. It follows from (3.9) that M (z) = D(z) + O(z), where D(z) is the diagonal matrix D(z) = ⎢⎢ ⎣ and O(z) is the off-diagonal matrix which corresponds to M (z). Proposition 1 det D(z) has exactly N zeros (counting multiplicities) in |z| < 1 and none satisfying |z| = 1. Step 2: Prove diagonal dominance for matrix M (t, z). Proof Consider an arbitrary i ∈ {1, 2, . . . , N }. |z − Aii (z)| ≥ |z| − | Aii (z)| ≥ 1 − Pii = for 0 ≤ t < 1, |z| = 1. j =i j =i We next turn to the case t = 1, |z| = 1, z = 1, again considering an arbitrary i ∈ {1, 2, . . . , N }. Now (3.13) is replaced by |z − Aii (z)| > j =i Pi j for |z| = 1, z = 1. On the other hand, j =i | Ai j (z)| < j =i Pi j . Therefore, |z − Aii (z)| > | j =i Ai j (z)| for |z| = 1, z = 1. This holds for i = 1, 2, . . . , N . In this way, we have proven the strict diagonal dominance, and hence the non-singularity, also for t = 1, |z| = 1, z = 1. Step 3: Prove that det M (t, z) has exactly N zeros in |z| < 1 and none on |z| = 1 for 0 ≤ t < 1. Proposition 3 The function det M (t, z) has exactly N zeros in |z| < 1 and none on |z| = 1 for 0 ≤ t < 1. Proof Let n(t ) be the number of zeros of det M (t, z) in |z| < 1. By the argument principle, see Evgrafov [8, p. 97], ∂∂z det M (t, z) where it should be noticed that det M (t, z) = 0 on |z| = 1 for 0 ≤ t < 1 according to Proposition 2. Here, n(t ) is a continuous integer-valued function of t for 0 ≤ t < 1 and n(0) = N according to Proposition 1. So n(t ) = n(0) = N . From the above, we may conclude that det M (1, z) = M (z) has at least N zeros in the closed unit disk, because the zeros of det M (t, z) are continuous functions for 0 ≤ t ≤ 1. Finally we need to prove that there are exactly N zeros in |z| ≤ 1, one of which (z = 1) lies on |z| = 1. Step 4: Use continuity of det M (t, z) in t for 0 ≤ t ≤ 1 to prove that det M (z) has N − 1 zeros in |z| < 1 and one zero z = 1 on |z| = 1. Proof Firstly, z = 1 is a zero of det M (z). Now we show that it is a simple zero. Use that limz→1 detz−M1(z) = ddz {det M (z)}|z=1 > 0, where the inequality is a consequence of the stability condition. Hence, z = 1 is a simple zero of det M (z). Proposition 5 det M (t, 1) > 0 for 0 ≤ t < 1. Proof We shall exploit the fact that det M (t, 1) is the product of all eigenvalues of M (t, 1). So we need to prove that the product of these eigenvalues is positive. Consider the matrix I − M (t, 1), where I is the identity matrix: Note that I − M (t, 1) is a substochastic matrix, so every eigenvalue of the matrix I − M (t, 1) lies in |z| < 1. Hence, every eigenvalue of the matrix M (t, 1) lies in |z − 1| < 1. M (t, 1) is a real matrix, so if M (t, 1) has a complex eigenvalue, then the conjugate of this complex eigenvalue is also one of the eigenvalues of M (t, 1). This implies that if M (t, 1) has complex eigenvalues, then the product of these complex eigenvalues is positive. The product of the real eigenvalues is also positive because every eigenvalue of the matrix M (t, 1) lies in |z − 1| < 1. This concludes the proof. Proposition 6 The function det M (z) has exactly N − 1 zeros in |z| < 1 and one zero on |z| = 1 (at z = 1). Proof We follow the argument of Gail et al. [9, p. 372]. By letting t → 1 in Proposition 3, it follows that det M (z) has at least N zeros in |z| ≤ 1. By Proposition 4, given > 0, there is a real z , 1 − < z < 1, such that det M (z ) is negative. By continuity, there is a real t , 1 − < t < 1, such that det M (t , z ) is negative. Since det M (t , 1) is positive according to Proposition 5, there is a real z , z < z < 1 with det M (t , z ) = 0. Thus, the zero of det M (z) at z = 1 is the limit of a zero of det M (t, z) from inside the unit disk. As t → 1, the limiting positions of the N zeros of det M (t, z) are: one at z = 1 and the other N − 1 in |z| < 1. 3.2 Transient analysis In this subsection, we shall determine the transient behavior of the probability generating function of the number of customers. The analysis proceeds largely analogously to the stationary case. In fact, for the transient analysis, it turns out to be less involved to demonstrate the location of the roots. We define n=0 f j (r, z) = r nE z Xn 1{Jn+1= j} for |r | < 1, j = 1, 2, ..., N , n=0 f j (r, 0) = r nP(Xn = 0, Jn+1 = j ). 1 N f j (r, z) =E z X0 1{J1= j} + z r nE z Xn−1 1{Jn=i} 1 N i=1 i=1 i=1,i = j n=1 n=1 n=0 r nP(Xn−1 = 0, Jn = i ) r n+1E z Xn 1{Jn+1=i} Ai j (z) fi (r, 0), provided the initial number of customers in the system is deterministic and equal to Using (3.15) and after simplification, we get the following N equations: (z − r A j j (z)) f j (r, z) − r Ai j (z) fi (r, z) = z X0+1P( J1 = j ) + r (B(z) − 1) Ai j (z) fi (r, 0), j = 1, 2, . . . , N . We can also write these N linear equations in matrix form as f (r, z) = ⎢⎢ ⎣ −r A12(z) . . . −r A1N (z) ⎤ z − r A22(z) . . . −r A2N (z) . . . . . . . . . ⎦⎥⎥ , −r AN 2(z) . . . z − r AN N (z) Therefore, solutions of the non-homogeneous linear system M (r, z)T f (r, z) = b(r, z) are of the form It remains to find the values of f1(r, 0), f2(r, 0), . . . , f N (r, 0). We shall derive N linear equations for f1(r, 0), f2(r, 0), . . . , f N (r, 0). To find N linear equations for f1(r, 0), f2(r, 0), . . . , f N (r, 0), we first prove that det M (r, z) has exactly N zeros for fixed r in |z| < 1. Since M (r, z) = z I − r A(z), det M (r, z) is a continuous function in r for 0 ≤ r ≤ 1, and therefore the zeros are continuous in 0 ≤ r ≤ 1. Remark It is worth emphasizing that it is at this point that our approach is different from the analysis by Neuts [14] and Çinlar [4]. We do not require for each pair of elementary roots that they either be strictly different for all values of 0 ≤ r ≤ 1 or coincide for all 0 ≤ r ≤ 1. The main price to pay is that we can not use that the roots are analytic in r and we can therefore not obtain the stationary distribution from the transient distribution as r → 1. Compared to the steady-state analysis, the proof is simpler and only consists of two steps: Step 1: Prove diagonal dominance of the matrix M (r, z). Proposition 7 det M (r, z) = 0 for 0 ≤ r < 1, |z| = 1. Proof Consider an arbitrary i ∈ {1, 2, . . . , N }. |z − r Aii (z)| ≥ |z| − r | Aii (z)| > 1 − Pii = Pi j for 0 ≤ r < 1, |z| = 1. (3.19) j =i j =i On the other hand, j =i |r Ai j (z)| ≤ r j =i Pi j for 0 ≤ r < 1, |z| = 1. Therefore, |z − r Aii (z)| > |r j =i Ai j (z)| for 0 ≤ r < 1, |z| = 1. This holds for i = 1, 2, . . . , N . Thus, M (r, z) is strictly diagonally dominant. This implies that M (r, z) is a nonsingular matrix, i.e., det M (r, z) = 0, for 0 ≤ r < 1, |z| = 1. This completes the proof. Proof Let n(r ) be the number of zeros of det M (r, z) in |z| < 1. As before, by the argument principle [8, p. 97], |z|=1 where it should be noticed that det M (r, z) = 0 on |z| = 1 for 0 ≤ r < 1 according to Proposition 7. Here, n(r ) is a continuous integer-valued function of r for 0 ≤ r < 1 and n(0) = N because det M (0, z) = zn . So n(r ) = n(0) = N . 4 Poisson batch arrivals: stationary queue length at arrival and arbitrary epochs In the previous section, we determined the stationary and the transient queue length distributions at departure times of customers. In the general framework, the exact arrival process of customers is not specified, but for the model with Poisson batch arrivals, we can obtain the stationary queue length distribution at arbitrary time, at batch arrival instants and at customer arrival instants. Because of PASTA, the distribution of the number of customers already in system just before a new batch arrives (let us denote this by a generic random variable X ba ) coincides with the distribution of the number of customers in the system at an arbitrary time ( X ar b). The number of customers at customer arrival instants (denoted with X ca ) needs to be further specified, because with batch arrivals all customers in the same batch have the same arrival time. As noted previously, customers within one batch are assumed to be (randomly) ordered. Although they arrive at the same time, they see different numbers of customers in front of them. In particular, the last customer in a batch sees all the customers that were already in the system plus all other customers (excluding him/her) arriving in the same batch. In the customer average distribution at arrival times, this must be taken into account. In Fig. 1 we depict three batch arrivals, two of which contain multiple customers and thus coincide with more than one customer arrival. Applying a simple level crossing argument with the aid of Fig. 1, it is readily seen that the distributions of X (at departure times) and X ca must coincide: indeed, for each level k = 1, 2, . . . , customer departures that decrease the queue length from k to k − 1 must be matched by customer arrivals increasing the level from k − 1 to k (since the arrival of each customer within a batch is counted separately, the difference can be at most 1, which is negligible in the long run). Fig. 1 Up- and down-crossing We can also link the distributions of X ba and X ca : A customer in an arriving batch sees in front of him the number of customers already in the system ( X ba ) and the number of customers in front of him in the same batch. For an arbitrary customer in the batch, the number of customers in front of him in the same batch has the forward recurrence distribution of B. Summarizing: where we use independence of the batch size and the number of customers already in system, and From these relations, we can obtain all the required distributions. It can be verified that these distributions agree with the results from Chaudhry[3] for the model without dependencies between successive service times. 5 The queueing model with two customer types : departure epochs In this section, we restrict ourselves to the case of two customer types, i.e., N = 2. In this case, we are able to give an explicit expression for the probability generating function of the number of customers in the system immediately after a departure. For the steady-state behavior, it follows from (3.8) that f1(z) = B(z) − 1 f1(0) (z A11(z) + A12(z) A21(z) − A11(z) A22(z)) + z f2(0) A21(z) (z − A11(z))(z − A22(z)) − A12(z) A21(z) (z − A11(z))(z − A22(z)) − A12(z) A21(z) f2(z) = F (z) = limn→∞E z Xn . F (z) = f1(z) + f2(z). After substituting the values of f1(z) and f2(z) from Eqs. (5.1) and (5.2), respectively, we obtain F (z) as F (z) = (z − A11(z))(z − A22(z)) − A12(z) A21(z) (B(z) − 1)( f1(0) + f2(0)) A12(z) A21(z) − A11(z) A22(z) (z − A11(z))(z − A22(z)) − A12(z) A21(z) Equation (3.2) states that Ai (z) = Ai1(z) + Ai2(z) for i = 1, 2. After substituting the values of fi (0) and Ai (z) for i = 1, 2, F (z) becomes F (z) = E[B] (z − A11(z))(z − A22(z)) − A12(z) A21(z) (B(z) − 1)(1 − ρ) A12(z) A21(z) − A11(z) A22(z) E[B] (z − A11(z))(z − A22(z)) − A12(z) A21(z) where c1 = A11(Azˆ)1+1(Azˆ)1−2(zˆzˆ)−zˆ , c2 = A21(Azˆ)2+2(Azˆ)2−2(zˆzˆ)−zˆ . After simplification, we can write F (z) as E[B] (z − A11(z))(z − A22(z)) − A12(z) A21(z) After differentiating F (z) w.r.t. z and taking the limit z → 1, we get ρ Var( A) E[B(B − 1)] E[X ] = 2 + 2(1 − ρ) + 2E[B] + −ρ + E[B]( f1(0)α1 + f(2P(012)α+2)P+21ρ)((1α−11 ρ+) α22) + α12α21 − α11α22 . For the transient distribution, it follows from (3.17) that F (z) = f1(r, z) = f2(r, z) = z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) r z(B(z) − 1) z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) r 2(B(z) − 1) A12(z) A21(z) − A11(z) A22(z) f1(r, 0) z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) r z(B(z) − 1) z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) r 2(B(z) − 1) A12(z) A21(z) − A11(z) A22(z) f2(r, 0) z − r A11(z) z − r A22(z) − r 2 A12(z) A21(z) (Bˆ (1) − 1)(Bˆ (2) − 1) Aˆ(221)(zˆ1 − r Aˆ(212)) − Aˆ(211)(zˆ2 − r Aˆ(222)) r zˆ2X0 (Bˆ (1) − 1) − zˆ1X0 (Bˆ (2) − 1) Aˆ(211) Aˆ(221)P(J1 = 2) (Bˆ (1) − 1)(Bˆ (2) − 1) Aˆ(221)(zˆ1 − r Aˆ(212)) − Aˆ(211)(zˆ2 − r Aˆ(222)) (Bˆ (1) − 1)(Bˆ (2) − 1) Aˆ(221)(zˆ1 − r Aˆ(212)) − Aˆ(211)(zˆ2 − r Aˆ(222)) zˆ2 − r Aˆ(222) P(J1 = 1) (Bˆ (1) − 1)(Bˆ (2) − 1) Aˆ(221)(zˆ1 − r Aˆ(212)) − Aˆ(211)(zˆ2 − r Aˆ(222)) f1(r, 0) = i, j = 1, 2. Remark 1 It can be observed that the first three terms on the right-hand-side of Eq. (5.5) are exactly equal to the mean queue length at departure epochs of the standard M X /G/1 queue without dependencies, cf. Gaver [10] and Cohen [5, Sect. III.2.3], and the remaining term appears due to the dependent service times. Remark 2 It can be shown, after some straightforward but tedious algebraic manipulations, that the queue length distribution in the system considered in the present paper also reduces to the distribution of the number of customers in an M X /G/1 queuing model if A1(z) = A2(z) = A(z), again cf. Gaver [10] and Cohen [5, Sect. III.2.3]. Similarly, we can also prove that the expected number of customers in the system considered in the present paper is equal to the expected number of customers in the corresponding M X /G/1 queuing model if α1 = α2 = E[ A]. 6 Numerical results In this section, we present four numerical examples in order to get more insight into the consequences of introducing dependencies between the service times of consecutive customers. For simplicity, we restrict ourselves to two customer types (N = 2). In all four examples, we assume that the overall batch arrival process is a Poisson process with rate λ and the load ρ equals 43 . 6.1 Example 1 In this example, we consider an almost symmetric system, with P( J = 1) = P( J = 2) = 21 and αi j = 83 for i, j = 1, 2. It follows that E[ A] = 43 , P11 = P22 and we shall vary P11. The batch sizes are geometrically distributed with P(B = k) = pk−1(1 − p), k = 1, 2, . . . We take p = 3/4, resulting in a mean batch size of E[B] = 4. The conditional service times are, respectively, exponential and Erlang distributed random variables, with m=0 for μi j > 0, i, j = 1, 2. In this example, we will take an Erlang distribution with four phases. If we define we can use Eq. (2.8) to obtain k j = Ai j (z) = Pi j for i = 1, 2 and j = 1, 2. The variance of the number of arrivals during one arbitrary service time, written as a function of P11, directly follows. For 0 < P11 < 1, 75 117 Var( A) = 16 + 512(1 − P11) P11 . We observe that α1 = α2, but A1(z) = A2(z). From Remark 2, we know that the mean queue length in our model is equal to the mean queue length of a standard M X /G/1 queue, but for higher moments of the queue length this equality is not true unless we can construct a case with A1(z) = A2(z). This is confirmed by Table 1, which depicts numerical values for the means and variances of the queue lengths in our model and in the corresponding M X /G/1 queue. Indeed, the mean queue lengths of both systems are equal, whereas the variances of the queue lengths are only equal in the case P11 = 21 , where A1(z) = A2(z). Since α1 = α2, we immediately conclude that the mean queue length and the variance of A are minimal when P11 = 1/2 (see Remark 2). Table 1 Means and variances of X and X M X /G/1 for various values of P11 in Example 1 E[X] = E X M X /G/1 Table 2 Mean queue length and variance of the number of arrivals during an arbitrary service time, for various values of P11 in Example 3. 6.2 Example 2 In this example, we take a similar setting as in the previous example, but we make two adjustments. First, for even more simplicity, we assume that all conditional service times are exponentially distributed, i.e., Gi j (x ) = (1 − e−μi j x ) Pi j , i, j = 1, 2. Secondly, we take α11 = α12 = 21 and α21 = α22 = 41 . As in the previous example, we let P( J = 1) = P( J = 2) = 21 . We observe that the difference with Example 1 is that all conditional service time distributions are exponential now, but with different parameters. Moreover, in this model α1 = α2. An interesting question is, how the mean queue length and the variance of the number of arrivals during an arbitrary service time are related. Since α1 = α2, the setting of Remark 2 does not apply. In Fig. 2, we show E[ X ] and Var( A) plotted versus P11. When studying the two plots carefully, one can see that the plots are not completely symmetric, which is obviously caused by the asymmetric service times. However, another observation that is not visible to the human eye is that the minima of both plots are not attained at the same value of P11. It can be shown analytically that the variance of A is minimal at exactly P11 = 1/2, and, numerically, that E[ X ] is minimal for P11 ≈ 0.500411. Although this is a small difference, it means that this system exhibits an interesting, rare feature: it is possible to obtain a smaller mean queue length by having a greater variance in the number of arrivals during one service time. In Example 3, we will create a setting in which this effect is even bigger. From Fig. 2a, b, we can observe that, except for the small region where 0.5 < P11 < 0.500411, the expected number of customers is increasing when the variance of the number of arrivals during a customer service time is increasing and conversely. This means that a bigger variance of the number of arrivals implies a larger expected number of customers. This also implies that the expected number of customers can grow beyond any bound in a stable system due to the very large variance of the number of arrivals during one service time. This scenario occurs when P11 tends to 0 or 1 in Fig. 2. Therefore, we can observe dependencies when P11 or (1 − P11) is small. Otherwise, E[ X ] and Var( A) appear to be rather insensitive to the value of P11. Fig. 2 Mean queue length E[X] and the variance of A in Example 2 Of course, the reason for the large variance in the number of arrivals during a customer service time lies in the dependence. When, for example, P11 = P22 is very small, services alternate for a long time between exp(μ12) and exp(μ21) services with small mean; rarely is there an exp(μ11) or exp(μ22) service which has a huge mean. 6.3 Example 3 Once again, we assume that the conditional service times are exponentially distributed, but in this example we choose less symmetric settings. Let P( J = 1) = 176 , P( J = Fig. 3 Variance of the number of arrivals versus the expected number of customers during an arbitrary customer service time. This implicit plot is obtained by varying P11. Figure (b) is a zoomed in version of Figure (a) Fig. 4 Numerical example 4: transient mean queue length analysis 7 P12, α1 = 0.3, and α2 = 1.1. The interesting phenomenon observed in P21 = 9 Example 2 is also taking place here. In fact, in this example there is a bigger difference between the value of P11 for which the mean queue length is minimal ( P11 ≈ 0.65), and the value resulting in a minimum variance of the number of arrivals during an arbitrary service time ( P11 ≈ 0.788) (in bold). More details can be found in Table 2. The interesting region is obviously 0.650 < P11 < 0.788, because in this region we know that an increase in Var( A) results in a decrease in E[ X ]. This is illustrated even better in Fig. 3, where Var( A) and E[ X ] are plotted against each other, for varying values of P11. 6.4 Example 4: transient-state analysis We return to the system in Example 2, but now we study the transient analysis. In this example, we start with an empty system, E[z X0 ] = 1, and set P11 = 1/10. Next, we repeatedly apply Eq. (3.1) to express E[z Xn ] in terms of E[z Xn−1 ]. We have taken four different distributions for the conditional service times, namely exponential, gamma with shape parameter 1/2, gamma with shape parameter 5, and deterministic. The results are shown in Fig. 4, where we depict the mean queue length after the departure of the nth customer, for n = 0, 1, 2, . . . , 200. In this example, it can clearly be seen that service time distributions with higher coefficients of variation result in longer queues. Also, it seems to take longer to reach steady state. For completeness, we give the steady-state mean queue lengths for the four systems below: Acknowledgements The research of Abhishek and Onno J. Boxma is partly funded by NWO Gravitation project Networks, Grant Number 024.002.003. The authors thank Michel Mandjes (University of Amsterdam) for helpful discussions. Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. 1. Abhishek , Boon, M.A.A., Mandjes , M.R.H. , Núñez-Queija , R. : Congestion analysis of unsignalized intersections . In: COMSNETS 2016 : Intelligent Transportation Systems Workshop ( 2016 ) 2. Adan , I.J.B.F., Kulkarni , V.G. : Single-server queue with Markov-dependent inter-arrival and service times . Queueing Syst . 45 , 113 - 134 ( 2003 ) 3. Chaudhry , M.L. : The queueing system M X /G/1 and its ramifications. Naval Res . Logist . Q. 26 , 667 - 674 ( 1979 ) 4. Çinlar , E. : Time dependence of queues with semi-Markovian services . J. Appl. Probab . 4 , 356 - 364 ( 1967 ) 5. Cohen , J.W.: The Single Server Queue . North-Holland, Amsterdam ( 1969 ) 6. de Smit, J.H.A. : The queue G I /M/s with customers of different types or the queue G I /Hm /s. Adv. Appl. Probab . 15 , 392 - 419 ( 1983 ) 7. de Smit, J.H.A. : The single server semi-Markov queue . Stoch. Process. Appl . 22 , 37 - 50 ( 1986 ) 8. Evgrafov , M.A.: Analytic Functions. Dover, New York ( 1978 ) 9. Gail , H.R. , Hantler , S.L. , Taylor , B.A.: On a preemptive Markovian queue with multiple servers and two priority classes . Math. Oper. Res . 17 , 365 - 391 ( 1992 ) 10. Gaver , D.P. : Imbedded Markov chain analysis of a waiting-line process in continuous time . Ann. Math. Stat . 30 , 698 - 720 ( 1959 ) 11. Gaver , D.P. : A comparison of queue disciplines when service orientation times occur . Naval Res . Logist . Q. 10 , 219 - 235 ( 1963 ) 12. Loynes , R.M. : The stability of a queue with non-independent inter-arrival and service times . Proc. Camb. Philos. Soc . 58 , 497 - 520 ( 1962 ) 13. Magalhães , M.N. , Disney , R.L. : Departures from queues with changeover times . Queueing Syst . 5 , 295 - 312 ( 1989 ) 14. Neuts , M.F. : The single server queue with Poisson input and semi-Markov service times . J. Appl. Probab . 3 , 202 - 230 ( 1966 ) 15. Neuts , M.F. : Some explicit formulas for the steady-state behavior of the queue with semi-Markovian service times . Adv. Appl. Probab . 9 , 141 - 157 ( 1977 ) 16. Neuts , M.F. : The M /G/1 queue with several types of customers and change-over times . Adv. Appl. Probab . 9 , 604 - 644 ( 1977 ) 17. Purdue , P. : A queue with Poisson input and semi-Markov service times: busy period analysis . J. Appl. Probab . 12 , 353 - 357 ( 1975 ) 18. Tachet , R. , Santi , P. , Sobolevsky , S. , Reyes-Castro , L.I. , Frazzoli , E. , Helbing , D. : Revisiting street intersections using slot-based systems . PLoS ONE 11 ( 3 ), e0149607 ( 2016 ). doi:10.1371/journal.pone. 0149607 19. Takács , L. : The transient behavior of a single server queueing process with a Poisson input . In: Proceedings of the Fourth Berkeley Symposium on Mathematical Statistics and Probability , vol. 2 , pp. 535 - 567 ( 1961 )


This is a preview of a remote PDF: https://link.springer.com/content/pdf/10.1007%2Fs11134-017-9531-4.pdf

Abhishek, Marko A. A. Boon, Onno J. Boxma, Rudesindo Núñez-Queija. A single-server queue with batch arrivals and semi-Markov services, Queueing Systems, 2017, 217-240, DOI: 10.1007/s11134-017-9531-4