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

We investigate the transient and stationary queue length distributions of a class of service systems with correlated service times. The classical MX/G/1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$M^X/G/1$$\end{document} 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(·)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$J(\cdot )$$\end{document} at a service initiation time determines the joint distribution of the subsequent service duration and the state of J(·)\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$J(\cdot )$$\end{document} 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.


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,9,12,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 [12], and more specifically its extension to batch arrivals [13] to compare our results with existing literature. The first to have investigated this class of queueing models was Gaver [9], 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 class-specific and when service switches from one type to the other, an additional switch-over time is required. This framework was generalized by Neuts [12], 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 [15] 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 [14], Neuts studied the multi-type 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 [13], 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 [11]. 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 [17] and Adan and Kulkarni [2]. Using the Wiener-Hopf factorization technique, De Smit [17] 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 [12] 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 [15] 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 Section 3, we derive the transient and stationary probability generating functions of the number of customers in the system immediately after a departure. In Section 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 Section 5. Finally, in Section 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.

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 J n , n = 0, 1, . . . , that can take values in {1, 2, . . . , N }, for some integer N . It will be convenient to refer to J n 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 [12]. We define (2.1) For future use we introduce the Laplace-Stieltjes transform (LST) where 1 {.} denotes the indicator function. In particular, 3) The type of a customer, and its service time, do not depend on the arrival process.
It should be observed that {J n , 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 J n is given by the unique solution of the set of equations with normalizing condition N j=1 P(J = j) = 1. The mean service time of an arbitrary customer is given by (2.5) The stability condition for this model is given by This can be formalized using Theorem 3 from Loynes [10], 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 J m the type of the first customer in the mth batch. Starting from a stationary version of the sequence (G (n) , J n+1 ), one can readily construct a stationary sequence (G (m) , J m+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 [10].
We will investigate the queue length process at departure times of customers. For that it will be convenient to define A n as the number of customers arriving during the service time of the nth customer and B n as the size of the batch in which the nth customer arrived. Note that for i, j = 1, 2, . . . , N , |z| ≤ 1, The queue length distribution at customer departure times is fully determined by the sequences A n and B n . 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 A n and B n 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).

General model
The inputs to our general model are probability generating functions of non-negative discrete random variables A ij (z), i, j ∈ {1, 2, . . . , N }, and B(z). From the A ij (z), we construct a Markov process (A n , J n+1 ), n = 1, 2, . . . , satisfying In this construction it is implicit that (A n , J n+1 ) conditional on J n is independent of A n−1 . The sequence B n is i.i.d. with generating function B(z) and independent of the sequence A n .
Next we define the recurrence relation Note: If the A ij (z) are set equal to (2.7), then the sequence X n follows the same law as the number of customers at departure times in the M X /G/1 queue with semi-Markov services. The role of the B n is subtle in this representation: B n 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 B n can be taken independent of the A n in the M X /G/1 queue with semi-Markov services.
In the sequel we will study the transient and stationary distributions of X n defined by (2.9). Again using Theorem 3 of Loynes [10], we may conclude that the stability condition in this case is Here E[A] denotes the expectation of the A n in stationarity: Note that at first sight (2.9) does not seem to fit the framework in Loynes [10], because of the special condition when the system is empty. For stability, however, the behavior of an empty system is irrelevant.

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: Now we exploit the fact that X n−1 and (A n , J n+1 ) are conditionally independent given J n , and the B n are also independent of all other random variables: for n = 1, 2, 3, . . . , j = 1, 2, . . . , N. (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 , and, where the α ij are defined in (2.11). Furthermore, for j = 1, 2, . . . , N , |z| ≤ 1: and note that, The probability generating function of the steady-state queue length distribution immediately after a departure is denoted by In steady state, Equation (3.1) leads to the following N equations: We can also write these N linear equations in matrix form as Therefore, solutions of the non-homogeneous linear system M (z) T f (z) = b(z) are in the form: whereê is a row vector with all entries one. After simplification, we can write this as and after simplification, we get, . (3.11)

(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 f i (z) is an analytic function in |z| < 1, the numerator of f i (z) also has N − 1 zeros in the unit disc |z| < 1. As a consequence, these N − 1 zeros provide To find the N − 1 zeros, we use a method that has also been applied in [2,7,16]. It is based on the concept of (strict) diagonal dominance in a matrix. The proof consists of 4 steps: Step 1: Prove that each element on the diagonal of M (z) has exactly one zero in |z| < 1. Step , and prove strict diagonal dominance of M (t, z), i.e., each diagonal element of M (t, z) is in absolute value larger than the sum of the absolute values of the non-diagonal terms in the same row of the matrix.
Step 4: Use continuity of det M (t, z) in t for 0 ≤ t < 1 to prove that, indeed, det M (z) has N − 1 zeros in |z| < 1 and one zero z = 1 on |z| = 1.
Step 1: Prove that each element on the diagonal of M (z) has exactly one zero in |z| < 1.
and O(z) is the off-diagonal matrix which corresponds to M (z).
Rouché's theorem implies that the numbers of zeros of z and z − A ii (z) are the same in |z| < 1. z has exactly one zero in |z| < 1, and hence z − A ii (z) also has exactly one zero in |z| < 1, for i = 1, 2, . . . , N .
Hence det D(z) has N zeros in |z| < 1 and none on |z| = 1.
Step 2: Prove diagonal dominance for matrix M (t, z).

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 [6, p. 97], 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 disc, 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.
where the inequality is a consequence of the stability condition. Hence, z = 1 is a simple zero of det M (z).
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.
Proof. We follow the argument of Gail et al. [7, 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 disc. 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.

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 analogous 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 so that, provided the initial number of customers in the system is deterministic and equal to x 0 .
Using (3.15) and after simplification, we get the following N equations: We can also write these N linear equations in matrix form as . . .
To find N linear equations for f 1 (r, 0), f 2 (r, 0), . . . , f N (r, 0), we fist prove that det M (r, z) has exactly N zeros for fixed r in |z| < 1. Since M (r, z) = zI − rA(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 [12] 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).
Proof. Let n(r) be the number of zeros of det M (r, z) in |z| < 1. As before, by the argument principle [6, p. 97], where it should be noticed that det M (r, z) = 0 on |z| = 1 for 0 ≤ r < 1 according to Proposition 7.

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 arb ). 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 Figure 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 argu- ment with the aid of Figure 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). 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 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.

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: It is noted that the probability generating function of X n in steady state is From Equation (3.7), for N = 2, we can write F (z) as the sum of f 1 (z) and f 2 (z), i.e., After substituting the values of f 1 (z) and f 2 (z) from Equations (5.1) and (5.2) respectively, we obtain F (z) as .

Equation (3.2) states that
After substituting the values of f i (0) and A i (z) for i = 1, 2, F (z) becomes After simplification, we can write F (z) as .

(5.4)
Let us now determine the expected number of customers E[X] = F (1).
After differentiating F (z) w.r.t. z and taking the limit z → 1, we get For the transient distribution it follows from (3.17) that where , (5.8)
Remark 1 It can be observed that the first three terms in the right-hand-side of Equation (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 [8] and Cohen [5,Section 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 A 1 (z) = A 2 (z) = A(z), again cf. Gaver [8] and Cohen [5,Section 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].

Numerical results
In this section, we present four numerical examples in order to get more insight in 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 3 4 .

Example 1
In this example we consider an almost symmetric system, with P(J = 1) = P(J = 2) = 1 2 and α ij = 3 8 for ∀i, j = 1, 2. It follows that E[A] = 3 4 , P 11 = P 22 and we shall vary P 11 . The batch sizes are geometrically distributed with 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 for µ ij > 0, i, j = 1, 2. In this example we will take an Erlang distribution with four phases. If we define we can use Equation (2.8) to obtain 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 P 11 , directly follows. For 0 < P 11 < 1, Var(A) = 75 16 + 117 512(1 − P 11 )P 11 .
We observe that α 1 = α 2 , but A 1 (z) = A 2 (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 A 1 (z) = A 2 (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 P 11 = 1 2 , where A 1 (z) = A 2 (z). Since α 1 = α 2 , we immediately conclude that the mean queue length and the variance of A are minimal when P 11 = 1/2 (see Remark 2).  Table 1: Means and variances of X and X M X /G/1 for various values of P 11 in Example 1.

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., G ij (x) = (1 − e −µij x )P ij , ∀i, j = 1, 2.
Secondly, we take α 11 = α 12 = 1 2 and α 21 = α 22 = 1 4 . As in the previous example, we let P(J = 1) = P(J = 2) = 1 2 . We observe that the difference with Example 1 is that all conditional servicetime distributions are exponential now, but with different parameters. Moreover, in this model 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 Figure 2 we show E[X] and Var(A) plotted versus P 11 . 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 P 11 . It can be shown analytically, that the variance of A is minimal at exactly P 11 = 1/2, and numerically, that E[X] is minimal for P 11 ≈ 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 Figures 2(a) and (b), we can observe that, except for the small region where 0.5 < P 11 < 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 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 P 11 tends to 0 or 1 in Figure 2. Therefore, we can observe dependencies when P 11 or (1 − P 11 ) is small. Otherwise, E[X] and Var(A) appear to be rather insensitive to the value of P 11 . Of course, the reason for the large variance in the number of arrivals during a customer service time lies in the dependence. When, e.g., P 11 = P 22 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.

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) = 7 16 , P(J = 2) = 9 16 , α 11 = α 12 = α 21 = 3 20 and α 22 = 19 20 . From these settings we obtain P 21 = 7 9 P 12 , α 1 = 0.3, and α 2 = 1.1. The interesting phenomenon observed in Example 2, is also taking place here. In fact, in this example there is a bigger difference between the value of P 11 for which the mean queue length is minimal (P 11 ≈ 0.65), and the value resulting in a minimum variance of the number of arrivals during an arbitrary service time (P 11 ≈ 0.788). More details can be found in Table 2. The interesting region is obviously 0.650 < P 11 < 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 Figure 3, where Var(A) and E[X] are plotted against each other, for varying values of P 11 .

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 P 11 = 1/10. Next, we repeatedly apply Equation (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 Figure 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: