- Email: [email protected]
Approximating nonstationary ph(t)⧸ph(t)⧸1⧸c queueing systems
Mathematics
and Computers
in Simulation
30 (1988)
441
441-452
North-Holland
APPROXIMATING
NONSTATIONARY
Ph( t ) / Ph( t ) / 1 /c
QUEUEING
SYSTEMS
*
Kim L. ONG BDM Corporation, McLean,
VA 22102, U.S.A.
Michael R. TAAFFE School of Industrial Engineering,
A state
space
time-dependent
partitioning behavior
with time-dependent tion requires in the arrival
phase
Kolmogorov-forward mean and standard
and arrival solution
and
distribution
systems
and service
equations
required
indicate
processes.
degardless
of entities
that the approximation
server
capacity,
where k,
process,
of solution.
in the system is extremely
approach single
of the system equations,
in the service
method
IN 47907, U.S.A
@DA)
for finite-capacity,
differential
of phases
for the classic
of the number
West Lafayette,
approximation
is described
of only k, + 3k, k,
is the number
k,
deviation
wide range of systems
surrogate
of queueing
the numerical process
Grissom Hall, Purdue University,
analyzing
queueing
Time-dependent Empirical
the
systems
c, the approxima-
is the number
compared
are obtained.
for
to the
of phases
k, + ck,k2
approximations test results
of
over a
accurate.
1. Introduction This paper considers finite-capacity single-server queueing systems in which the arrival and service processes are nonstationary phase processes. An approximation is developed based on a state-space partitioning and surrogate distribution approximation (SDA) methods. The nonstationary, phase process, and an approximation for the Ph(t)/M(t)/s/c model using another state-space partitioning and SDAs are described in [7] and [4] (see also [6] and [5] for a review of related nonstationary queueing approximations for simpler systems). The inter-arrival and service time random variables for the Ph( t)/Ph( t)/l/c system are the times till absorption in nonstationary Markov processes with exactly one absorbing state. Phase distributions, phase renewal processes and their use in modeling stationary queueing systems are described in [3]. Nonstationary phase processes are defined in [7]. Phase distributions can assume a wide variety of shapes. Nonstationary generalizations of Erlang and hyperexponential point processes are both special cases of Ph( t) processes. Nonstationary phase processes not only allow for the variety of shapes in the instantaneous distributions but also allow modeling of time-dependencies. Thus, the Ph( t)/Ph( t)/l/c is a very flexible system which is of practical value in modeling real time-dependent phenomena. Because phase distributions are dense over the set of distributions with support on the nonnegative reals, one can think of the Ph(t)/Ph(t)/l/c queueing model as a good approximation to the nonstationary or time-dependent generalization of the GI/G/l/c model. * This research
0378-4754/88/$3.50
was partially
0 1988,
supported
Elsevier
by the National
Science
Publishers
Science
Foundation
B.V. (North-Holland)
under Grant
No. ECS-8404409.
442
K. L. Ong M. R. Taaff
/ Ph( t)/Ph(
t )/l/c
queueing systems
Time-dependent state probabilities for the Ph( t)/Ph(t)/l/c model can be evaluated by numerical integration of nonhomogeneous Kolmogorov-forward differential difference equations. Unfortunately, the number of Kolmogorov-forward equations can be too large for convenient computation. The goal in approximating the time-dependent state probabilities and moments of the queue size distribution for the Ph( t)/Ph( t)/l/c system is to greatly reduce the number of differential equations requiring numerical integration to a small number with little loss in accuracy. Observing the qualitative as well as the quantitative time-dependent behavior of performance queues can be measures, such as E[N( t)] or P{ N(t) = 0) versus time, for nonstationary insightful. Observing only period averages of performance measures such as
‘IJ,k[
N(r)]
dr
or
t-‘]‘P{
N(T)
= O}
dr
0
can disguise much of the underlying variation of the process and lead to an incomplete if not distorted analysis of the system. Neuts [3] has discussed many aspects of the insight to be gained by studying queue dynamics.
2. The Ph(t)/Ph(t)/l/c
queueing system
Consider a single server queueing system with nonstationary phase arrival process (A(t), X(t)) Markov and nonstationary service process (B(t), p(t)), where A( t ) and B(t) are underlying chains of the arrival and service Markov processes, A(t)
= [ 2;;; ~~~~~‘__~~~~ 1 “:“‘],
B(t)
= [!!f!;!~-Bb(,!],
and a(t) and /3(t) are the initial state probability vectors after an absorption (an arrival or service completion) for the arrival and service processes at time t, respectively. The queue discipline is independent of actual customer service times. The elements of matrices A,(t), A2( t), B,(t), and B2( t), are the Markov routing probabilities for these nonstationary Markov processes at time t. The matrices have dimensions k, x k,, k, X 1, k, X k,, and k, X 1, respectively. The dimensions of a(t) and /3(t) are respectively k, and k,. The h(t) and p(t) vectors are time-dependent and state-dependent rate vectors for the phases of the arrival and service processes with dimensions k, and k,, respectively (see [7] for a description of the relationship of this parameterization of the nonstationary phase process and Neuts’ [3] parameterization of stationary phase processes). c model is described by the triple (N(t), I(t), J(t)), The random state of a Ph(t)/Ph(t)/l/ where N(t) is the number of customers in the system at time t, I(t) is the current phase of the next arrival at time t, and J(t) is the service phase of the customer in service at time t. The value J(t) is 0 if and only if the server is idle. The state space S is {(0, i, 0) 1i E I } U {(n, i, j) 1n = 1, 2,. . . , c, i E I, j E J}, where I = k,} is the set of all arrival phase indices and J = { 1, 2,. . . , k, } is the set of all service (1, 2,..., phase indices. The arrival phases are labeled 1 through k, and state k, + 1 is the instantaneous arrival state, i.e., the entry into state k, + 1 represents the customer arrival. Similarly, the service
K. L. Ong, M. R. Taaffe / Ph( t )/Ph( t )/l/c
phases are labeled 1 through k, and Because states (n, i, k, + l), ( n, k, + n = 0, 1,. . . , c are instantaneous states, system. The size of the state space, which cardinality of S:
443
queueing systems
state k, + 1 is the instantaneous end-of-service phase. 1, j), and (n, k, + 1, k, + 1) for i E I, j E J U {0}, and they are not included in the state space of the queueing is the number
of Kolmogorov-forward
equations,
1S I= k, + ck,k,.
is the
(1)
1S 1increases by k,k, for each additional queue space. Thus, for a system with a finite capacity of 100 and 5 phases respectively in the arrival and service processes, the number of differential equations to be numerically integrated is 2505, a significant amount of computing. The ) S 1 Kolmogorov-forward equations for the Ph( t)/Ph( t)/l/c model are
Ww&)/dt = -h(t)%&)
+
c
u,.;(t>h,(t)P,.,,,(t)
lEI,(O
(4 W,,,,,(t)/dt = -
[X,(t)+ P,(t)]fL,,(t) + 4,,,,(l)
c
~l.k,+lw,(fPl(4
,=f,,(f)
x [~nlp,m.,,o(t) + +
c
c
47n,~c ,,,,(t))]
b,,,(t)~.,(t)P,,,.,,(t)
m@,(t)
s,,= 1 -
6,,=
a,, ,
1, i 0,
i=A,
otherwise,
Ti,, = 1 - CVrA,
and ICI,(t)= Z,(t)=
{iliEI, (Zll~1,
ai(t)>O}, uI,i(t)>O},
ICI,(t) is the set of all possible
iEI,
initial arrival
and
I,(t)=
phase indices
(111~1,
a,,,,+,(t)>O).
at time
t. t. The sets Jo,<
K. L. Ong, hf. R. Taaffe
444
/ Ph( t)/Ph(
The p th partial moment of the distribution in phases i and j, respectively, is E,(y)(t)=
~r~~P,,;,~(t), n=O
iEI,
t)/l/c
of N(t)
jEJ,
queueing
systems
for the arrival and service processes
p=O,l,...
.
being
(4)
So the pth conditional moment of N(t) conditioned on i and j is E/~)(t)/E$‘( t). The unconditional pth moment of N(t) is ECP’(t) = X~&Z$&E~~‘( t). The follo&ng theorem is the foundation of the approximation to ECp’( t). Theorem. The PMDEs
for the Ph( t)/Ph(
dE:,;‘(t)/dt=
- [A,(t)
t)/l/c
system are
+/+)]E,rf(t)
+ 41&
c
al,k,+l(t)ai(t)Xl(t)
/Elo(t)
P,(t)Po,r,o(d +i
q=o
(;)E,‘;‘(t) - ~,,ii(t)~~~[~)~~ q=o
i j:)(-l)P-‘Ei%t) - q,,~,,~,,(t)
Lq=rJ
iEI,
jEJ,
p=O,
’
l,...
2 ’
.
1,
(5)
Proof dE,‘,;)(t)/dt
= c
np dPn,,,j(t)/dt,
i E I,
p=o,
j E J,
1 ,..a
3
n=l =
-
[ xi(t>
+
P,tt)]
E1lT’(t>
c-1
[
x P,wpo,t,o(d+
L
c nPPn+l,r,m (t )I
fl=l
1
c nPPn-l,,,,(t)+ cPPc,r,,(f)
n=2
(6)
K. L. Ong, M. R. Taaffe
proof is completed
/ Ph( t)/Ph(
t)/l/c
queueing
445
systems
by using the relationships
and c-1 c n=l
in (6).
nPPn+l
, im ,
Cl>
=
It
Cn
-
l)PPn,l,m(t)
=
fl
Cn
-
l)PPn,i,m(t)
-
8ppOPl,*,m(t>
fl=l
n=2
0
Note that the only probability ~o,,,o(t),
PI,,,,(t)
terms appearing and
PC,,,,(t),
in the right-hand
side of PMDEs
i ~1, jE J, lo Io(f), m E Jo(f).
(5) are:
(7)
The approximation exploits this structure. To numerically integrate PMDEs (5), values of (7) need to be available at every time t. The algorithm developed in the next section is based on an approximation to these probabilities.
3. The approximation Our approach in approximating multivariate queueing systems is to artfully partition the state space into subspaces of related states. Associated with each subspace is a distribution of the number of customers in the system, conditioned on being in that subspace. SDAs are then used to match two moments and numerically integrate the first two PMDEs for each subspace. Differential equations for the probabilities of being in each of the subspaces (the 0th partial moments) are also numerically integrated to obtain the joint state probabilities. Crucial to this kind of SDA approach is the choice of an efficient partitioning scheme. One needs to choose a partition such that the effect of largest dimension (in this case, the largest dimension is the number in the system which has size c + 1) is diminished and at the same time the resulting number of probability terms on the right-hand side of the PMDEs is small. Consider partitioning the Ph( t)/Ph( t)/l/c system state space into k,k, + 1 subspaces, one subspace for each pair of arrival and service phases, (I(t), J(t)), and one subspace consisting of the empty system states, {(N(t) = 0, I(t) = i, J(t) = 0);i E I}. The number of conditional
K. L. Ong M. R. Taaff
446
/ Ph( t )/Ph( t )/l/c
queue@
systems
distributions to describe a nonempty system is k,k,. To evaluate the three moments of each conditional distribution, only three sets of probabilities must be computed or approximated. The first probabilities necessary to evaluate the PMDEs at time t are the P,,,,O( t)‘s, where 1 E 1,(t) of equation (7). Although a surrogate distribution could be used to approximate these probabilities, numerically integrating the set of Kolmogorov-forward equations (2) is convenient since the number of equations k, is usually small. The second set of state probabilities at time t, PI,,,,(t), can be approximated by PE surrogate probabilities. The surrogate distributions for P{N(t)=nlI(t)=i,
J(t)=m},
n=l,2
,...,
c,
iEI,
mEJo(
are PE distributions, denoted by PE,( i, m ; t). The number of PE,( i, m ; t) surrogate distributions at time t is k, 1J,(t) 1, i.e., the product of the number of arrival phases times and the number of terminal service phases at time t. The number of PE,( i, m ; t) surrogate distributions at time t ranges over k,, . . . , k,k,, depending on 1.I,( t) 1. Further, only one probability from each of these surrogate distributions is needed at time t to evaluate PMDEs (5). The third set of probabilities needed at time t to evaluate PMDEs (5) is the set of P,,,,,(t) probabilities. The surrogate distributions for P{N(t)=nlI(t)=Z,
J(t)=j},
n=l,2
,...)
c,
/q)(t),
jEJ,
r,j;
are also PE, denoted by PE,.( t). The number of surrogate distributions at time t is k, ( I,(t) 1, so the number of surrogate distributions ranges over k,, . . . , k,k,, depending on 1I,(t) I. Again only one probability from each surrogate distribution PE,( I, j ; t) is needed at time t to evaluate the PMDEs (5). Some of the PE,(i, m ; t) and the PE,(I, j; t) distributions are the same. These distributions correspond to pairs (I, m) such that 1 E 1,(t) and m E J,(t)_ Both the probability of the high and low order states are approximated for these distributions. The total number of surrogate distributions is k, 1J,(t) 1+ k, ( I,(t) I - I I,,(t) I I.&,(t) 1, and the total number of surrogate probabilities at time t is k, I Jo(t) I + k, I lo(t) I. Three PMDEs (Oth, lst, and 2nd) for each conditional distribution of a nonempty system are equations for states needed in the approximation. In addition, the k, Kolmogorov-forward (N(t) = 0, I(t) = i, J(t) = 0) where i E I, are needed. Therefore, 3k,k, + k, differential equations yield the conditional moments for the conditional distributions associated with each subspace, as well as the unconditional partial moments. The surrogate probabilities, k, Kolmogorov-forward equations, and the PMDEs comprise a closed set of differential equations which can conveniently be numerically integrated. The approximation for the Ph( t)/Ph( t)/l/c system is summarized as follows: (1) At time t, match E~fl,(t)/E~~(t) and E$(t)/E/,z(t) to the first two moments of PE,( i, m ; t), for every pair (i, m) where i E I and m E J,(t). (2) At time t, match Ej:)( t)/El0,‘( t) and E$)( t)/E{;)( t) to the first two moments of PE,(I, j; t), for every pair (I, j) where I E lo(t) and j E J. (3) Approximate the P,,i,m( t) terms by the product of the probability of n = 1 in the PE,(i, m ; t) surrogate distribution and EyA( t), for each pair (i, m) where i E I and m E J,(t).
K. L. Ong, M. R. Taajje / Ph( t)/Ph(
t)/l/c
queueing systems
447
(4) Approximate the Pc,,,j (t) terms by the product of the probability of n = c in the PE,( 1, j; t) surrogate distribution and E{‘,)(t), for each pair (I, j) where I E 1,(t) and j E J. (5) Simultaneously numerically integrate the PMDEs (5) for p = 0, 1, 2 and the Kolmogorovforward equations (2) corresponding to states (0, i, 0) for i E I to get approximations for E:y’( t + At), E,‘j’( t + At), and E,‘:‘( t + At) for every pair (i, j) where i E I and j E J. (6) Sum the partial moments to get the first two moments of the distribution of the number of customers in the system at time t. (7) Set t = t + At and go to step (1). The above procedure is repeated for the time interval of interest.
4. Test cases The time-dependent mean and standard deviation of the number of customers in the system are compared to the values obtained by numerically integrating the full set of Kolmogorov-forward equations. In the set of twenty test cases, the following input parameters were varied: capacity (30, 20) number of arrival phases (3, 5), number of service phases (3, 5, 8, 15), average traffic intensity (light, medium, heavy), and structure of the arrival and service processes (series, parallel, mixed, and feedback). The initial condition for all cases reported here was an empty and idle system, although the initial arrival-process phase varied. The rates of the arrival and service nonstationary phase processes were sine and cosine functions of time plus constant terms and the routing matrices had fixed structures, i.e., (A, X(t)) and (B, p(t)). Four error measurements are shown in Tables 1 and 2: the maximum absolute and relative errors across time in the time-dependent expected number of customers in the system for each case and the maximum absolute and relative errors across time in the time-dependent standard deviation of the number of customers in the system for each case. Table 1 lists the maximum absolute error and the corresponding relative error as well as the maximum relative error and the corresponding absolute error in the expected number of customers in the system. Table 2 is a similar list of maximum relative and absolute errors for the time-dependent standard deviation of the number of customers in the system. Note that these tables report the sup norm errors; a conservative measure of accuracy. To indicate the error level at times other than the time at which the maximum error occurs two plots tracking the errors across time are included. The cases plotted are of the cases with the greatest maximum error in E [N( t)] and SD[ N( t)]. No pattern of error is evident in this set of test cases. The approximation is accurate in light, heavy, and medium traffic. Accuracy is slightly better for light and heavy traffic than for medium traffic. Across all twenty cases, the greatest absolute error in the time-dependent expected number of customers in the system, in case 3, is 0.748 customers, which corresponds to a 4.9% error. The greatest relative error for this performance measure across all twenty cases, in case 20, is 6.8%, which corresponds to an absolute error of 0.403 customers. The greatest absolute error in the time-dependent standard deviation of the number of customers across all twenty cases, in case 14, is 0.907 customers, which corresponds to a 29.2% error. And, across all twenty cases, the greatest relative error for this performance measure, again in case 14, is 48.9% which
K. L. Ong, M. R. Taajje / Ph( t )/Ph( t )/l/c
448 Table 1 Maximum
errors in Ph( t)/l/c
queueing systems
E[ N( t)] approximation
Case number
Traffic load
Maximum absolute error
Relative value
Maximum relative error
Absolute value
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
heavy light medium light heavy heavy heavy medium light heavy light medium medium heavy heavy light light medium medium light
0.364 0.071 0.748 0.036 0.256 0.312 0.279 0.252 0.075 0.407 0.050 0.174 0.208 0.435 0.209 0.207 0.282 0.434 0.427 0.358
0.014 0.020 0.049 0.017 0.009 0.011 0.011 0.028 0.034 0.022 0.015 0.026 0.035 0.016 0.007 0.039 0.019 0.045 0.029 0.037
0.014 0.030 0.052 0.018 0.009 0.011 0.011 0.040 0.036 0.022 0.016 0.048 0.036 0.016 0.008 0.043 0.019 0.048 0.029 0.068
0.364 0.060 0.572 0.035 0.259 0.312 0.279 0.224 0.064 0.382 0.044 0.106 0.172 0.435 0.193 0.188 0.281 0.413 0.426 0.403
Table 2 Maximum Case number 1
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
errors in Ph( t)/Ph(
t)/l/c
SD[ N( t)] approximation
Traffic load
Maximum absolute error
Relative value
Maximum relative error
Absolute value
heavy light medium light heavy heavy heavy medium light heavy light medium medium heavy heavy light light medium medium light
0.626 0.171 0.188 0.073 0.427 0.512 0.591 0.204 0.097 0.629 0.204 0.151 0.226 0.907 0.392 0.234 0.440 0.160 0.232 0.276
0.084 0.071 0.033 0.037 0.097 0.229 0.124 0.032 0.072 0.082 0.069 0.107 0.054 0.292 0.083 0.069 0.185 0.029 0.040 0.066
0.107 0.073 0.033 0.037 0.472 0.237 0.298 0.032 0.074 0.082 0.071 0.111 0.054 0.489 0.084 0.070 0.197 0.029 0.040 0.066
0.274 0.168 0.188 0.073 0.398 0.489 0.201 0.204 0.094 0.629 0.203 0.143 0.226 0.819 0.390 0.233 0.435 0.160 0.232 0.276
K. L. Ong, M.R. Taaff
/ Ph( t)/Ph(
t)/l/c
queueing systems
corresponds to an absolute error of 0.81 customers. So at no time in any case tested was the approximation for either performance measure as much as 1. customer off. No structure of arrival or service process emerged as being significantly more accurate than any other. Feedback in the arrival or service process also had no appreciable affect on the accuracy of the approximation. Likewise capacity did not appear to be important in predicting the accuracy of the approximation. Although capacity was varied over a large number of cases, only capacities of 20 and 30 are reported here since accuracy appeared to be insensitive to this parameter. Part of the test case selection strategy was to examine cases one would expect to be accurately approximated (call these ‘good cases’) and compare the results to cases one would expect to be least accurately approximated (call these ‘bad cases’). The heuristic strategy to determine the robustness of the approximation is to evaluate the difference in accuracy between these two sets of cases. If the difference is significant, then further exploration would be needed to determine the regions of acceptable accuracy. If the difference in accuracy is not significant, then the approximation is robust and the region of accuracy is large. Define the ‘bad cases’ as those that heavily rely on surrogate probabilities. ‘Bad cases’ were investigated to determine if they could serve as empirically determined worst case error bounds. Likewise, define ‘good cases’ as those that rely very little on surrogate probabilities. ‘Good cases’ were evaluated to investigate limits of accuracy of the approximation. Heavy dependence on the approximated probabilities can be due to either a large number of nonzero surrogate probabilities or to a large (significantly different from 0.0) total value of surrogate probabilities. Consider the ratio of the number of surrogate probabilities at time t to the number of PMDEs and Kolmogorov-forward equations numerically integrated. This ratio is (k, 1J,(t) 1+ k, 1I,(t) I)/( k, + 3k,k,). The set of surrogate probabilities is Pi,,,,(t) and P,,,,,(t), where i E I, m E .I,( t), 1 E I,(t), and j E J. Recall that I,,(t) and J,(t) are the sets of indices corresponding to the terminal phases of the arrival and service processes at time t. Systems that are very busy and have a small I I,(t) I do not heavily depend on approximations. Busy systems have small total value for the P I,r,m( t) terms, where i E I and m E J,(t), regardless of how many there are (e.g., even for large I J,(t) I) and if the system also has a small number of phases in I I,(t) 1, then there will be a small number of P,,,.,(t) terms, where 1 E I,(t) and j E J. These cases one would expect to be accurately approximated and will be called ‘good cases’. Likewise for systems that are not busy and have small I J,(t) I. Systems that are not busy have small total value of the P,,,,,(t), where 1 E I,,(t) and j E J regardless of the number of phase indices in 1,(t); and if the system also has a small number of phases in J,(t), then there will be a small number of P1,,,,( t) t erms where i E I and m E J,(t). These cases one would also expect to be accurately approximated and will also be called ‘good cases’. In fact, when traffic is such that P{ N(t) = c} + 0 and traffic is also heavy enough such that P{ N(t) = l} + 0, the approximation does not rely on surrogate probabilities and again one would expect extremely accurate approximations. Analogously, for systems that are not busy and have large I J,(t) lone would expect less accurate approximations than those for the systems described above; these systems will be called ‘bad cases’. Likewise, systems that are busy and have large I I,(t) ( should also be less accurately approximated; these will also be called ‘bad cases’. Another condition where one might reasonably expect to find ‘good cases’ are cases where the mean time spent in individual states (N(t) = 1, I(t) = i, J(t) = m) and (N(t) = C, I(t) = 1, J(t)
450
K. L. Ong, M. R. Taaffe / Ph( t)/Ph(
t)/l/c
queueing systems
=j), where i E I, m E .I,( t), 1 E I,(t), and j E J, is very small relative to other states and 1J,(t) 1 = ) Z,(t) I= 1. A system with the arrival and service processes both having infinite rates in their single final phases should be very accurately approximated. Test results for ‘good cases’ and ‘bas cases’ all were approximated extremely accurately. The approximation procedure is robust because both kinds of cases produced results with no significant difference in accuracy. In other words, cases expected to be most accurately approximated were, and cases expected to be least accurately approximated were almost as accurately approximated. A reasonable conclusion is that the approximation is uniformly accurate. One note about the last ‘good case’ described above. Consider cases with artificial final arrival and service phases having significantly faster rates than other phases (e.g., 100 times faster or approximately instantaneous phases). Our heuristic to find ‘good cases’ suggests that these cases would be close to exactly approximated. One might be tempted to include such artificial phases in an attempt to marginally increase the already excellent accuracy. In these cases, accurate and
13.07 Case 20
ll.W
___
approximation
-
actual
\ 9.806
8.171
E[N(t)l 6.537
3.269
I-
1.634
0.000 0
I l-
.ooo
1.250
2.500
3.750
5.000 Time
Fig. 1.
6.250
7.500
8.750
10 .oo
K. L. Ong, M. R. Taaffe / Ph( r)/Ph(
r )/l/c
451
queueing systems Case 14
8.000
-_-
approximation
-
actual
SD[N(t)l Lt.000
3.000
2.000
1.000
o.ooc
0.000
2.500
5.000
7.500
10.00
12.50
15.00
17.50
20 .oo
Time
Fig. 2.
stable numerical integration requires much smaller step sizes because of the large rates at the artificial final phases. Our test cases demonstrate that accuracy was not significantly improved. The negligible gain in accuracy is achieved at a great increase in computing cost. These results further support the conclusion that the approximation is uniformly accurate. The two plots in Figs. 1 and 2 also illustrate the accuracy of the approximation. Fig. 1 shows the actual and approximated values of E[ N( t)] versus time for case 20, which has the greatest maximum relative error in the time-dependent expected number of customers in the system. Fig. 2 shows the actual and approximated values of SD[ N( t)] versus time, for case 14, which has the greatest absolute and relative error in the approximated standard deviation of the time-dependent number of customers in the system. Both of these plots clearly indicate that even for these worst cases the approximation tracks the actual values quite well. The approximation for E [N( t)] is superior to the SD[ N( t)] approximation. Because of the wide variety of shapes that phase distributions can assume and the time-dependent input, this SDA-based approximation for the Ph( t)/Ph( t)l/c model appears to provide an efficient method of approximating a large class of single server queueing systems.
452
K. L. Ong, M. R. Taaffe / Ph( t)/Ph( t)/l/c
queueing systems
Acknowledgment
The authors
wish to thank Bruce W. Schmeiser
for his may helpful
suggestions.
References [l] G.M. Clark, Use of Polya distributions in approximate solutions to nonstationary M/M/s queues, Comm. ACM 24 (4) (1981) 206-218. [2] N.L. Johnson and S. Kotz, Urn Models and Their Applications (Wiley, New York, 1977). [3] M.F. Neuts, Matrix-Geometric Solutions in Stochastic ModelsAn Algorithmic Approach (The Johns Hopkins University Press, Baltimore, 1981). [4] K.L. Ong, Approximating Nonstationaty Multivariate Queueing Models, Ph.D. Dissertation, Purdue Univ., 1985. [5] M.R. Taaffe, Approximating Nonstationary Queueing Models, Ph.D. Dissertation, The Ohio State Univ., 1982. [6] M.R. Taaffe and G.M. Clark, Approximating nonstationary two priority non-preemptive queueing systems, Naval Research Logistics 35 (1) (1988) 125-145.. [7] M.R. Taaffe and K.L. Ong, Approximating nonstationary Ph(t)/M(t)/s/c queueing systems, Ann. Oper. Res. 9 (1987) 103-116.
and Computers
in Simulation
30 (1988)
441
441-452
North-Holland
APPROXIMATING
NONSTATIONARY
Ph( t ) / Ph( t ) / 1 /c
QUEUEING
SYSTEMS
*
Kim L. ONG BDM Corporation, McLean,
VA 22102, U.S.A.
Michael R. TAAFFE School of Industrial Engineering,
A state
space
time-dependent
partitioning behavior
with time-dependent tion requires in the arrival
phase
Kolmogorov-forward mean and standard
and arrival solution
and
distribution
systems
and service
equations
required
indicate
processes.
degardless
of entities
that the approximation
server
capacity,
where k,
process,
of solution.
in the system is extremely
approach single
of the system equations,
in the service
method
IN 47907, U.S.A
@DA)
for finite-capacity,
differential
of phases
for the classic
of the number
West Lafayette,
approximation
is described
of only k, + 3k, k,
is the number
k,
deviation
wide range of systems
surrogate
of queueing
the numerical process
Grissom Hall, Purdue University,
analyzing
queueing
Time-dependent Empirical
the
systems
c, the approxima-
is the number
compared
are obtained.
for
to the
of phases
k, + ck,k2
approximations test results
of
over a
accurate.
1. Introduction This paper considers finite-capacity single-server queueing systems in which the arrival and service processes are nonstationary phase processes. An approximation is developed based on a state-space partitioning and surrogate distribution approximation (SDA) methods. The nonstationary, phase process, and an approximation for the Ph(t)/M(t)/s/c model using another state-space partitioning and SDAs are described in [7] and [4] (see also [6] and [5] for a review of related nonstationary queueing approximations for simpler systems). The inter-arrival and service time random variables for the Ph( t)/Ph( t)/l/c system are the times till absorption in nonstationary Markov processes with exactly one absorbing state. Phase distributions, phase renewal processes and their use in modeling stationary queueing systems are described in [3]. Nonstationary phase processes are defined in [7]. Phase distributions can assume a wide variety of shapes. Nonstationary generalizations of Erlang and hyperexponential point processes are both special cases of Ph( t) processes. Nonstationary phase processes not only allow for the variety of shapes in the instantaneous distributions but also allow modeling of time-dependencies. Thus, the Ph( t)/Ph( t)/l/c is a very flexible system which is of practical value in modeling real time-dependent phenomena. Because phase distributions are dense over the set of distributions with support on the nonnegative reals, one can think of the Ph(t)/Ph(t)/l/c queueing model as a good approximation to the nonstationary or time-dependent generalization of the GI/G/l/c model. * This research
0378-4754/88/$3.50
was partially
0 1988,
supported
Elsevier
by the National
Science
Publishers
Science
Foundation
B.V. (North-Holland)
under Grant
No. ECS-8404409.
442
K. L. Ong M. R. Taaff
/ Ph( t)/Ph(
t )/l/c
queueing systems
Time-dependent state probabilities for the Ph( t)/Ph(t)/l/c model can be evaluated by numerical integration of nonhomogeneous Kolmogorov-forward differential difference equations. Unfortunately, the number of Kolmogorov-forward equations can be too large for convenient computation. The goal in approximating the time-dependent state probabilities and moments of the queue size distribution for the Ph( t)/Ph( t)/l/c system is to greatly reduce the number of differential equations requiring numerical integration to a small number with little loss in accuracy. Observing the qualitative as well as the quantitative time-dependent behavior of performance queues can be measures, such as E[N( t)] or P{ N(t) = 0) versus time, for nonstationary insightful. Observing only period averages of performance measures such as
‘IJ,k[
N(r)]
dr
or
t-‘]‘P{
N(T)
= O}
dr
0
can disguise much of the underlying variation of the process and lead to an incomplete if not distorted analysis of the system. Neuts [3] has discussed many aspects of the insight to be gained by studying queue dynamics.
2. The Ph(t)/Ph(t)/l/c
queueing system
Consider a single server queueing system with nonstationary phase arrival process (A(t), X(t)) Markov and nonstationary service process (B(t), p(t)), where A( t ) and B(t) are underlying chains of the arrival and service Markov processes, A(t)
= [ 2;;; ~~~~~‘__~~~~ 1 “:“‘],
B(t)
= [!!f!;!~-Bb(,!],
and a(t) and /3(t) are the initial state probability vectors after an absorption (an arrival or service completion) for the arrival and service processes at time t, respectively. The queue discipline is independent of actual customer service times. The elements of matrices A,(t), A2( t), B,(t), and B2( t), are the Markov routing probabilities for these nonstationary Markov processes at time t. The matrices have dimensions k, x k,, k, X 1, k, X k,, and k, X 1, respectively. The dimensions of a(t) and /3(t) are respectively k, and k,. The h(t) and p(t) vectors are time-dependent and state-dependent rate vectors for the phases of the arrival and service processes with dimensions k, and k,, respectively (see [7] for a description of the relationship of this parameterization of the nonstationary phase process and Neuts’ [3] parameterization of stationary phase processes). c model is described by the triple (N(t), I(t), J(t)), The random state of a Ph(t)/Ph(t)/l/ where N(t) is the number of customers in the system at time t, I(t) is the current phase of the next arrival at time t, and J(t) is the service phase of the customer in service at time t. The value J(t) is 0 if and only if the server is idle. The state space S is {(0, i, 0) 1i E I } U {(n, i, j) 1n = 1, 2,. . . , c, i E I, j E J}, where I = k,} is the set of all arrival phase indices and J = { 1, 2,. . . , k, } is the set of all service (1, 2,..., phase indices. The arrival phases are labeled 1 through k, and state k, + 1 is the instantaneous arrival state, i.e., the entry into state k, + 1 represents the customer arrival. Similarly, the service
K. L. Ong, M. R. Taaffe / Ph( t )/Ph( t )/l/c
phases are labeled 1 through k, and Because states (n, i, k, + l), ( n, k, + n = 0, 1,. . . , c are instantaneous states, system. The size of the state space, which cardinality of S:
443
queueing systems
state k, + 1 is the instantaneous end-of-service phase. 1, j), and (n, k, + 1, k, + 1) for i E I, j E J U {0}, and they are not included in the state space of the queueing is the number
of Kolmogorov-forward
equations,
1S I= k, + ck,k,.
is the
(1)
1S 1increases by k,k, for each additional queue space. Thus, for a system with a finite capacity of 100 and 5 phases respectively in the arrival and service processes, the number of differential equations to be numerically integrated is 2505, a significant amount of computing. The ) S 1 Kolmogorov-forward equations for the Ph( t)/Ph( t)/l/c model are
Ww&)/dt = -h(t)%&)
+
c
u,.;(t>h,(t)P,.,,,(t)
lEI,(O
(4 W,,,,,(t)/dt = -
[X,(t)+ P,(t)]fL,,(t) + 4,,,,(l)
c
~l.k,+lw,(fPl(4
,=f,,(f)
x [~nlp,m.,,o(t) + +
c
c
47n,~c ,,,,(t))]
b,,,(t)~.,(t)P,,,.,,(t)
m@,(t)
s,,= 1 -
6,,=
a,, ,
1, i 0,
i=A,
otherwise,
Ti,, = 1 - CVrA,
and ICI,(t)= Z,(t)=
{iliEI, (Zll~1,
ai(t)>O}, uI,i(t)>O},
ICI,(t) is the set of all possible
iEI,
initial arrival
and
I,(t)=
phase indices
(111~1,
a,,,,+,(t)>O).
at time
t. t. The sets Jo,<
K. L. Ong, hf. R. Taaffe
444
/ Ph( t)/Ph(
The p th partial moment of the distribution in phases i and j, respectively, is E,(y)(t)=
~r~~P,,;,~(t), n=O
iEI,
t)/l/c
of N(t)
jEJ,
queueing
systems
for the arrival and service processes
p=O,l,...
.
being
(4)
So the pth conditional moment of N(t) conditioned on i and j is E/~)(t)/E$‘( t). The unconditional pth moment of N(t) is ECP’(t) = X~&Z$&E~~‘( t). The follo&ng theorem is the foundation of the approximation to ECp’( t). Theorem. The PMDEs
for the Ph( t)/Ph(
dE:,;‘(t)/dt=
- [A,(t)
t)/l/c
system are
+/+)]E,rf(t)
+ 41&
c
al,k,+l(t)ai(t)Xl(t)
/Elo(t)
P,(t)Po,r,o(d +i
q=o
(;)E,‘;‘(t) - ~,,ii(t)~~~[~)~~ q=o
i j:)(-l)P-‘Ei%t) - q,,~,,~,,(t)
Lq=rJ
iEI,
jEJ,
p=O,
’
l,...
2 ’
.
1,
(5)
Proof dE,‘,;)(t)/dt
= c
np dPn,,,j(t)/dt,
i E I,
p=o,
j E J,
1 ,..a
3
n=l =
-
[ xi(t>
+
P,tt)]
E1lT’(t>
c-1
[
x P,wpo,t,o(d+
L
c nPPn+l,r,m (t )I
fl=l
1
c nPPn-l,,,,(t)+ cPPc,r,,(f)
n=2
(6)
K. L. Ong, M. R. Taaffe
proof is completed
/ Ph( t)/Ph(
t)/l/c
queueing
445
systems
by using the relationships
and c-1 c n=l
in (6).
nPPn+l
, im ,
Cl>
=
It
Cn
-
l)PPn,l,m(t)
=
fl
Cn
-
l)PPn,i,m(t)
-
8ppOPl,*,m(t>
fl=l
n=2
0
Note that the only probability ~o,,,o(t),
PI,,,,(t)
terms appearing and
PC,,,,(t),
in the right-hand
side of PMDEs
i ~1, jE J, lo Io(f), m E Jo(f).
(5) are:
(7)
The approximation exploits this structure. To numerically integrate PMDEs (5), values of (7) need to be available at every time t. The algorithm developed in the next section is based on an approximation to these probabilities.
3. The approximation Our approach in approximating multivariate queueing systems is to artfully partition the state space into subspaces of related states. Associated with each subspace is a distribution of the number of customers in the system, conditioned on being in that subspace. SDAs are then used to match two moments and numerically integrate the first two PMDEs for each subspace. Differential equations for the probabilities of being in each of the subspaces (the 0th partial moments) are also numerically integrated to obtain the joint state probabilities. Crucial to this kind of SDA approach is the choice of an efficient partitioning scheme. One needs to choose a partition such that the effect of largest dimension (in this case, the largest dimension is the number in the system which has size c + 1) is diminished and at the same time the resulting number of probability terms on the right-hand side of the PMDEs is small. Consider partitioning the Ph( t)/Ph( t)/l/c system state space into k,k, + 1 subspaces, one subspace for each pair of arrival and service phases, (I(t), J(t)), and one subspace consisting of the empty system states, {(N(t) = 0, I(t) = i, J(t) = 0);i E I}. The number of conditional
K. L. Ong M. R. Taaff
446
/ Ph( t )/Ph( t )/l/c
queue@
systems
distributions to describe a nonempty system is k,k,. To evaluate the three moments of each conditional distribution, only three sets of probabilities must be computed or approximated. The first probabilities necessary to evaluate the PMDEs at time t are the P,,,,O( t)‘s, where 1 E 1,(t) of equation (7). Although a surrogate distribution could be used to approximate these probabilities, numerically integrating the set of Kolmogorov-forward equations (2) is convenient since the number of equations k, is usually small. The second set of state probabilities at time t, PI,,,,(t), can be approximated by PE surrogate probabilities. The surrogate distributions for P{N(t)=nlI(t)=i,
J(t)=m},
n=l,2
,...,
c,
iEI,
mEJo(
are PE distributions, denoted by PE,( i, m ; t). The number of PE,( i, m ; t) surrogate distributions at time t is k, 1J,(t) 1, i.e., the product of the number of arrival phases times and the number of terminal service phases at time t. The number of PE,( i, m ; t) surrogate distributions at time t ranges over k,, . . . , k,k,, depending on 1.I,( t) 1. Further, only one probability from each of these surrogate distributions is needed at time t to evaluate PMDEs (5). The third set of probabilities needed at time t to evaluate PMDEs (5) is the set of P,,,,,(t) probabilities. The surrogate distributions for P{N(t)=nlI(t)=Z,
J(t)=j},
n=l,2
,...)
c,
/q)(t),
jEJ,
r,j;
are also PE, denoted by PE,.( t). The number of surrogate distributions at time t is k, ( I,(t) 1, so the number of surrogate distributions ranges over k,, . . . , k,k,, depending on 1I,(t) I. Again only one probability from each surrogate distribution PE,( I, j ; t) is needed at time t to evaluate the PMDEs (5). Some of the PE,(i, m ; t) and the PE,(I, j; t) distributions are the same. These distributions correspond to pairs (I, m) such that 1 E 1,(t) and m E J,(t)_ Both the probability of the high and low order states are approximated for these distributions. The total number of surrogate distributions is k, 1J,(t) 1+ k, ( I,(t) I - I I,,(t) I I.&,(t) 1, and the total number of surrogate probabilities at time t is k, I Jo(t) I + k, I lo(t) I. Three PMDEs (Oth, lst, and 2nd) for each conditional distribution of a nonempty system are equations for states needed in the approximation. In addition, the k, Kolmogorov-forward (N(t) = 0, I(t) = i, J(t) = 0) where i E I, are needed. Therefore, 3k,k, + k, differential equations yield the conditional moments for the conditional distributions associated with each subspace, as well as the unconditional partial moments. The surrogate probabilities, k, Kolmogorov-forward equations, and the PMDEs comprise a closed set of differential equations which can conveniently be numerically integrated. The approximation for the Ph( t)/Ph( t)/l/c system is summarized as follows: (1) At time t, match E~fl,(t)/E~~(t) and E$(t)/E/,z(t) to the first two moments of PE,( i, m ; t), for every pair (i, m) where i E I and m E J,(t). (2) At time t, match Ej:)( t)/El0,‘( t) and E$)( t)/E{;)( t) to the first two moments of PE,(I, j; t), for every pair (I, j) where I E lo(t) and j E J. (3) Approximate the P,,i,m( t) terms by the product of the probability of n = 1 in the PE,(i, m ; t) surrogate distribution and EyA( t), for each pair (i, m) where i E I and m E J,(t).
K. L. Ong, M. R. Taajje / Ph( t)/Ph(
t)/l/c
queueing systems
447
(4) Approximate the Pc,,,j (t) terms by the product of the probability of n = c in the PE,( 1, j; t) surrogate distribution and E{‘,)(t), for each pair (I, j) where I E 1,(t) and j E J. (5) Simultaneously numerically integrate the PMDEs (5) for p = 0, 1, 2 and the Kolmogorovforward equations (2) corresponding to states (0, i, 0) for i E I to get approximations for E:y’( t + At), E,‘j’( t + At), and E,‘:‘( t + At) for every pair (i, j) where i E I and j E J. (6) Sum the partial moments to get the first two moments of the distribution of the number of customers in the system at time t. (7) Set t = t + At and go to step (1). The above procedure is repeated for the time interval of interest.
4. Test cases The time-dependent mean and standard deviation of the number of customers in the system are compared to the values obtained by numerically integrating the full set of Kolmogorov-forward equations. In the set of twenty test cases, the following input parameters were varied: capacity (30, 20) number of arrival phases (3, 5), number of service phases (3, 5, 8, 15), average traffic intensity (light, medium, heavy), and structure of the arrival and service processes (series, parallel, mixed, and feedback). The initial condition for all cases reported here was an empty and idle system, although the initial arrival-process phase varied. The rates of the arrival and service nonstationary phase processes were sine and cosine functions of time plus constant terms and the routing matrices had fixed structures, i.e., (A, X(t)) and (B, p(t)). Four error measurements are shown in Tables 1 and 2: the maximum absolute and relative errors across time in the time-dependent expected number of customers in the system for each case and the maximum absolute and relative errors across time in the time-dependent standard deviation of the number of customers in the system for each case. Table 1 lists the maximum absolute error and the corresponding relative error as well as the maximum relative error and the corresponding absolute error in the expected number of customers in the system. Table 2 is a similar list of maximum relative and absolute errors for the time-dependent standard deviation of the number of customers in the system. Note that these tables report the sup norm errors; a conservative measure of accuracy. To indicate the error level at times other than the time at which the maximum error occurs two plots tracking the errors across time are included. The cases plotted are of the cases with the greatest maximum error in E [N( t)] and SD[ N( t)]. No pattern of error is evident in this set of test cases. The approximation is accurate in light, heavy, and medium traffic. Accuracy is slightly better for light and heavy traffic than for medium traffic. Across all twenty cases, the greatest absolute error in the time-dependent expected number of customers in the system, in case 3, is 0.748 customers, which corresponds to a 4.9% error. The greatest relative error for this performance measure across all twenty cases, in case 20, is 6.8%, which corresponds to an absolute error of 0.403 customers. The greatest absolute error in the time-dependent standard deviation of the number of customers across all twenty cases, in case 14, is 0.907 customers, which corresponds to a 29.2% error. And, across all twenty cases, the greatest relative error for this performance measure, again in case 14, is 48.9% which
K. L. Ong, M. R. Taajje / Ph( t )/Ph( t )/l/c
448 Table 1 Maximum
errors in Ph( t)/l/c
queueing systems
E[ N( t)] approximation
Case number
Traffic load
Maximum absolute error
Relative value
Maximum relative error
Absolute value
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
heavy light medium light heavy heavy heavy medium light heavy light medium medium heavy heavy light light medium medium light
0.364 0.071 0.748 0.036 0.256 0.312 0.279 0.252 0.075 0.407 0.050 0.174 0.208 0.435 0.209 0.207 0.282 0.434 0.427 0.358
0.014 0.020 0.049 0.017 0.009 0.011 0.011 0.028 0.034 0.022 0.015 0.026 0.035 0.016 0.007 0.039 0.019 0.045 0.029 0.037
0.014 0.030 0.052 0.018 0.009 0.011 0.011 0.040 0.036 0.022 0.016 0.048 0.036 0.016 0.008 0.043 0.019 0.048 0.029 0.068
0.364 0.060 0.572 0.035 0.259 0.312 0.279 0.224 0.064 0.382 0.044 0.106 0.172 0.435 0.193 0.188 0.281 0.413 0.426 0.403
Table 2 Maximum Case number 1
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
errors in Ph( t)/Ph(
t)/l/c
SD[ N( t)] approximation
Traffic load
Maximum absolute error
Relative value
Maximum relative error
Absolute value
heavy light medium light heavy heavy heavy medium light heavy light medium medium heavy heavy light light medium medium light
0.626 0.171 0.188 0.073 0.427 0.512 0.591 0.204 0.097 0.629 0.204 0.151 0.226 0.907 0.392 0.234 0.440 0.160 0.232 0.276
0.084 0.071 0.033 0.037 0.097 0.229 0.124 0.032 0.072 0.082 0.069 0.107 0.054 0.292 0.083 0.069 0.185 0.029 0.040 0.066
0.107 0.073 0.033 0.037 0.472 0.237 0.298 0.032 0.074 0.082 0.071 0.111 0.054 0.489 0.084 0.070 0.197 0.029 0.040 0.066
0.274 0.168 0.188 0.073 0.398 0.489 0.201 0.204 0.094 0.629 0.203 0.143 0.226 0.819 0.390 0.233 0.435 0.160 0.232 0.276
K. L. Ong, M.R. Taaff
/ Ph( t)/Ph(
t)/l/c
queueing systems
corresponds to an absolute error of 0.81 customers. So at no time in any case tested was the approximation for either performance measure as much as 1. customer off. No structure of arrival or service process emerged as being significantly more accurate than any other. Feedback in the arrival or service process also had no appreciable affect on the accuracy of the approximation. Likewise capacity did not appear to be important in predicting the accuracy of the approximation. Although capacity was varied over a large number of cases, only capacities of 20 and 30 are reported here since accuracy appeared to be insensitive to this parameter. Part of the test case selection strategy was to examine cases one would expect to be accurately approximated (call these ‘good cases’) and compare the results to cases one would expect to be least accurately approximated (call these ‘bad cases’). The heuristic strategy to determine the robustness of the approximation is to evaluate the difference in accuracy between these two sets of cases. If the difference is significant, then further exploration would be needed to determine the regions of acceptable accuracy. If the difference in accuracy is not significant, then the approximation is robust and the region of accuracy is large. Define the ‘bad cases’ as those that heavily rely on surrogate probabilities. ‘Bad cases’ were investigated to determine if they could serve as empirically determined worst case error bounds. Likewise, define ‘good cases’ as those that rely very little on surrogate probabilities. ‘Good cases’ were evaluated to investigate limits of accuracy of the approximation. Heavy dependence on the approximated probabilities can be due to either a large number of nonzero surrogate probabilities or to a large (significantly different from 0.0) total value of surrogate probabilities. Consider the ratio of the number of surrogate probabilities at time t to the number of PMDEs and Kolmogorov-forward equations numerically integrated. This ratio is (k, 1J,(t) 1+ k, 1I,(t) I)/( k, + 3k,k,). The set of surrogate probabilities is Pi,,,,(t) and P,,,,,(t), where i E I, m E .I,( t), 1 E I,(t), and j E J. Recall that I,,(t) and J,(t) are the sets of indices corresponding to the terminal phases of the arrival and service processes at time t. Systems that are very busy and have a small I I,(t) I do not heavily depend on approximations. Busy systems have small total value for the P I,r,m( t) terms, where i E I and m E J,(t), regardless of how many there are (e.g., even for large I J,(t) I) and if the system also has a small number of phases in I I,(t) 1, then there will be a small number of P,,,.,(t) terms, where 1 E I,(t) and j E J. These cases one would expect to be accurately approximated and will be called ‘good cases’. Likewise for systems that are not busy and have small I J,(t) I. Systems that are not busy have small total value of the P,,,,,(t), where 1 E I,,(t) and j E J regardless of the number of phase indices in 1,(t); and if the system also has a small number of phases in J,(t), then there will be a small number of P1,,,,( t) t erms where i E I and m E J,(t). These cases one would also expect to be accurately approximated and will also be called ‘good cases’. In fact, when traffic is such that P{ N(t) = c} + 0 and traffic is also heavy enough such that P{ N(t) = l} + 0, the approximation does not rely on surrogate probabilities and again one would expect extremely accurate approximations. Analogously, for systems that are not busy and have large I J,(t) lone would expect less accurate approximations than those for the systems described above; these systems will be called ‘bad cases’. Likewise, systems that are busy and have large I I,(t) ( should also be less accurately approximated; these will also be called ‘bad cases’. Another condition where one might reasonably expect to find ‘good cases’ are cases where the mean time spent in individual states (N(t) = 1, I(t) = i, J(t) = m) and (N(t) = C, I(t) = 1, J(t)
450
K. L. Ong, M. R. Taaffe / Ph( t)/Ph(
t)/l/c
queueing systems
=j), where i E I, m E .I,( t), 1 E I,(t), and j E J, is very small relative to other states and 1J,(t) 1 = ) Z,(t) I= 1. A system with the arrival and service processes both having infinite rates in their single final phases should be very accurately approximated. Test results for ‘good cases’ and ‘bas cases’ all were approximated extremely accurately. The approximation procedure is robust because both kinds of cases produced results with no significant difference in accuracy. In other words, cases expected to be most accurately approximated were, and cases expected to be least accurately approximated were almost as accurately approximated. A reasonable conclusion is that the approximation is uniformly accurate. One note about the last ‘good case’ described above. Consider cases with artificial final arrival and service phases having significantly faster rates than other phases (e.g., 100 times faster or approximately instantaneous phases). Our heuristic to find ‘good cases’ suggests that these cases would be close to exactly approximated. One might be tempted to include such artificial phases in an attempt to marginally increase the already excellent accuracy. In these cases, accurate and
13.07 Case 20
ll.W
___
approximation
-
actual
\ 9.806
8.171
E[N(t)l 6.537
3.269
I-
1.634
0.000 0
I l-
.ooo
1.250
2.500
3.750
5.000 Time
Fig. 1.
6.250
7.500
8.750
10 .oo
K. L. Ong, M. R. Taaffe / Ph( r)/Ph(
r )/l/c
451
queueing systems Case 14
8.000
-_-
approximation
-
actual
SD[N(t)l Lt.000
3.000
2.000
1.000
o.ooc
0.000
2.500
5.000
7.500
10.00
12.50
15.00
17.50
20 .oo
Time
Fig. 2.
stable numerical integration requires much smaller step sizes because of the large rates at the artificial final phases. Our test cases demonstrate that accuracy was not significantly improved. The negligible gain in accuracy is achieved at a great increase in computing cost. These results further support the conclusion that the approximation is uniformly accurate. The two plots in Figs. 1 and 2 also illustrate the accuracy of the approximation. Fig. 1 shows the actual and approximated values of E[ N( t)] versus time for case 20, which has the greatest maximum relative error in the time-dependent expected number of customers in the system. Fig. 2 shows the actual and approximated values of SD[ N( t)] versus time, for case 14, which has the greatest absolute and relative error in the approximated standard deviation of the time-dependent number of customers in the system. Both of these plots clearly indicate that even for these worst cases the approximation tracks the actual values quite well. The approximation for E [N( t)] is superior to the SD[ N( t)] approximation. Because of the wide variety of shapes that phase distributions can assume and the time-dependent input, this SDA-based approximation for the Ph( t)/Ph( t)l/c model appears to provide an efficient method of approximating a large class of single server queueing systems.
452
K. L. Ong, M. R. Taaffe / Ph( t)/Ph( t)/l/c
queueing systems
Acknowledgment
The authors
wish to thank Bruce W. Schmeiser
for his may helpful
suggestions.
References [l] G.M. Clark, Use of Polya distributions in approximate solutions to nonstationary M/M/s queues, Comm. ACM 24 (4) (1981) 206-218. [2] N.L. Johnson and S. Kotz, Urn Models and Their Applications (Wiley, New York, 1977). [3] M.F. Neuts, Matrix-Geometric Solutions in Stochastic ModelsAn Algorithmic Approach (The Johns Hopkins University Press, Baltimore, 1981). [4] K.L. Ong, Approximating Nonstationaty Multivariate Queueing Models, Ph.D. Dissertation, Purdue Univ., 1985. [5] M.R. Taaffe, Approximating Nonstationary Queueing Models, Ph.D. Dissertation, The Ohio State Univ., 1982. [6] M.R. Taaffe and G.M. Clark, Approximating nonstationary two priority non-preemptive queueing systems, Naval Research Logistics 35 (1) (1988) 125-145.. [7] M.R. Taaffe and K.L. Ong, Approximating nonstationary Ph(t)/M(t)/s/c queueing systems, Ann. Oper. Res. 9 (1987) 103-116.