DTIC ADA442652: Interpolation Approximations for M/G/infinity Arrival Processes

17 Pages·0.33 MB·English

by  Defense Technical Information Center

T R R
ECHNICAL ESEARCH EPORT
Interpolation approximations for $M|G|infty$ arrival
processes
by Konstantinos P. Tsoukatos, Armand M. Makowski
CSHCN T.R. 99-35
(ISR T.R. 99-69)
The Center for Satellite and Hybrid Communication Networks is a NASA-sponsored Commercial Space
Center also supported by the Department of Defense (DOD), industry, the State of Maryland, the University
of Maryland and the Institute for Systems Research. This document is a technical report in the CSHCN
series originating at the University of Maryland.
Web site http://www.isr.umd.edu/CSHCN/
Report Documentation Page Form Approved
OMB No. 0704-0188
Public reporting burden for the collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and
maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information,
including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 1204, Arlington
VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to a penalty for failing to comply with a collection of information if it

does not display a currently valid OMB control number.
1. REPORT DATE 3. DATES COVERED
1999 2. REPORT TYPE -
4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER
Interpolation approximations for $M|G|infty$ arrival processes
5b. GRANT NUMBER
5c. PROGRAM ELEMENT NUMBER
6. AUTHOR(S) 5d. PROJECT NUMBER
5e. TASK NUMBER
5f. WORK UNIT NUMBER
7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION
Army Research Laboratory,2800 Powder Mill Road,Adelphi,MD,20783 REPORT NUMBER
9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S)
11. SPONSOR/MONITOR’S REPORT
NUMBER(S)
12. DISTRIBUTION/AVAILABILITY STATEMENT
Approved for public release; distribution unlimited
13. SUPPLEMENTARY NOTES
The original document contains color images.
14. ABSTRACT
see report
15. SUBJECT TERMS
16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 18. NUMBER 19a. NAME OF
ABSTRACT OF PAGES RESPONSIBLE PERSON

a. REPORT b. ABSTRACT c. THIS PAGE 16
unclassified unclassified unclassified
Standard Form 298 (Rev. 8-98)
Prescribed by ANSI Std Z39-18
Sponsored by: ARL, NASA, NSF
Interpolation approximations
j j1
for M G arrival processes
(cid:3) y
Konstantinos P. Tsoukatos Armand M. Makowski
Electrical Engineering Department and Institute for Systems Research
University of Maryland, College Park, MD 20742
Abstract
We present an approximate analysis of a discrete{time queue with correlated ar-
rival processes of the so{called MjGj1 type. The proposed heuristic approximations
are developed around asymptotic results in the heavy and light tra(cid:14)c regimes. Inves-
tigation of the system behavior in light tra(cid:14)c quanti(cid:12)es the di(cid:11)erences between the
gradual MjGj1 inputs and the point arrivals of a classical GIjGIj1 queue. In heavy
tra(cid:14)c, salient features are e(cid:11)ectively captured by the exponential distribution and the
Mittag{Le(cid:15)er special function, under short{ and long{range dependence respectively.
By interpolating between the heavy and light tra(cid:14)c extremes we derive approxima-
tions to thequeuesize distribution, applicable to all tra(cid:14)c intensities. We examine the
accuracy of these expressions and discuss possible extensions of our results in several
numerical examples.
1 Introduction
The conclusions of a series of measurement studies demonstrating that network tra(cid:14)c ex-
hibits persistent long term correlations have spurred recent activity in the study of queueing
systems with correlated arrival processes. Analytical results reveal that, when strong de-
pendencies are present, diverse queueing patterns may arise, in contrast to the familiar
exponential decay encountered in traditional tra(cid:14)c models with bounded exponential mo-
ments.
(cid:3)TheworkofthisauthorwassupportedthroughNSFGrantNSFDCDR-88-03012andtheArmyResearch
Laboratory under Cooperative Agreement No. DAAL01-96-2-0002.
yTheworkofthisauthorwassupportedpartiallythroughNSFGrantNSFDCDR-88-03012,NASAGrant
NAGW277S and the Army Research Laboratory under Cooperative Agreement No. DAAL01-96-2-0002.
1
In this paper we deal with a discrete time queue, viewed as a surrogate for a network
multiplexer, driven by an MjGj1 arrival stream. Both the discrete time MjGj1 process
considered here and its continuous time variant are among the tra(cid:14)c models arising from
large aggregations of on/o(cid:11) sources, that have attracted a great deal of attention. Reasons
for this, such as flexibility in representing correlation functions of actual tra(cid:14)c traces, are
discussed in [9]. A fluid queue fed by the continuous time version of the process has been
studied at least as early as 1974 [3]. Later, Cox notices that the MjGj1 busy server process
with heavy tailed G is a second order asymptotically self{similar process [4].
Here, we are interested in the entire steady{state queue length distribution at a multi-
plexer fed by the MjGj1 arrival process. For arbitrary pmf G, the system lacks the desired
Markovian structure and a calculation using numerical inversion techniques, e.g. [13], is not
possible, since no z-transform expressions are available. An exception is provided by the ge-
ometric case, where a two dimensional Markov chain formulation and a functional equation
for the z-transform is given in [2].
To circumvent the di(cid:14)culties of an exact analysis, one may rely on information gleaned
from various asymptotic regimes. A promising approach consists of deriving approximations
from the analysis of large bu(cid:11)er asymptotics [7, 8, 10]; these estimates are exact in the limit
as the bu(cid:11)er level goes to in(cid:12)nity.
Our objective is to explore alternative approximations to the queue length probabilities,
developed arounda combination of lightandheavy tra(cid:14)casymptotics. Such approximations
become exact in the limit as the tra(cid:14)c intensity goes to zero and one respectively. In light
tra(cid:14)c we take advantage of the fact that the MjGj1 arrival process is obviously \Poisson
driven", sothattheReiman{Simontheory[12]applies, underaboundedexponentialmoment
assumption. The resulting light tra(cid:14)c limits of the queue with MjGj1 arrivals di(cid:11)er from
those of a classical GIjGIj1 queue. This is a manifestation of the fact that work that joins
the system gradually, as is the case with MjGj1 inputs, generates less queueing than work
that arrives instantaneously. From the heavy tra(cid:14)c regime, we collect the associated limit
distribution of the queue size: This is given through the exponential function in the standard
short{range dependent setup, and the Mittag{Le(cid:15)er special function in the case where the
MjGj1 process is long{range dependent. The approximation to the queue size distribution
is subsequently generated by interpolating between the heavy and light tra(cid:14)c extremes. For
some common pmfs G the approximant assumes a simple (cid:12)nal form. More interestingly, it
has the potential of capturing accurately the queue size distribution at small bu(cid:11)er sizes,
for which approximations based on large bu(cid:11)er asymptotics are usually ill (cid:12)tted. On the
other hand, when G has (cid:12)nite exponential moment, we do not expect the heavy{light tra(cid:14)c
interpolation to be accurate for bu(cid:11)er sizes much larger than the maximum burst length: It
simply does not possess the correct decay rate { it does so only as the tra(cid:14)c intensity tends
to one, i.e., in the heavy tra(cid:14)c limit. Surprisingly, this drawback is often absent under long{
range dependence, since there are cases where the queue size distribution has hyperbolic
asymptotics with the same exponent for all tra(cid:14)c intensities! Then an approximation is
more valuable, especially when considering that, in the presence of heavy tails, alternative
estimates by means of simulation take an unreasonably long time to obtain. Yet, this is
2
somewhat compromised by the unaivalability of rigorously established light tra(cid:14)c limits,
under long{range dependence. In Section 6 we rely on a postulated relationship, but the
problem is still unresolved.
The paper is organized asfollows: The description of thesystem model is given in Section
2. Section 3 contains the main conclusions of the light and heavy tra(cid:14)c analyses. These are
the ingredients for the approximation, which is presented in Section 4 and discussed through
numerical examples in Section 5. Further extensions of the results are suggested in Section
6.
2 The System Model
We introduce the queueing model of interest, together with the required notation. We start
by presenting the MjGj1 arrival processes and several of their properties; additional facts
can be found in [4, 10].
Consider a population of in(cid:12)nitely many information sources, operating in discrete{time.
Sources can be in one of two states, active or idle. During time slot [n,n+1), n = 0,1,...,
β new sources become active. Source j, j = 1,...,β begins generating information by
n+1 n+1
thestartofslot[n+1,n+2),itsactivityperiodhasdurationσ (innumberofslots). While
n+1,j
active, each source emits information at a constant rate of one information unit (packet) per
time slot. After its activity period expires, each source switches o(cid:11) permanently, never to
generate packets again. Let b denote the number of active sources, or equivalently, the
n
number of packets generated by the active sources at the beginning of time slot [n,n+1).
If initially (i.e., at time n = 0) there were already b active sources, we denote by σ the
0,j
residual activity duration (in time slots) for the jth active source, j = 1,...,b.
Throughout, the IN{valued rvs b, fβ , n = 0,1,...g, fσ , n = 1,2,...; j = 1,2,...g
n+1 n,j
and fσ , j = 1,2,...g satisfy the following assumptions: (i) These rvs are mutually inde-
0,j
pendent; (ii) The rvs fβ , n = 0,1,...g are i.i.d. Poisson rvs with parameter λ > 0; (iii)
n+1
The rvs fσ , n = 1,...; j = 1,2,...g are i.i.d. with common pmf G on f1,2,...g. Let
n,j
σ be a generic IN{valued rv distributed according to the pmf G, assume throughout that
E[σ] < 1; (iv) The rvs fσ , j = 1,2,...g are i.i.d. IN{valued rvs distributed accord-
0,j
ing to the equilibrium pmf G associated with G, i.e., if σ denotes a generic IN{valued rv
e e
distributed according to the pmf G , then
e
P[σ (cid:21) n]
P[σ = n] = , n = 1,2,... (1)
e
E[σ]
In summary, the process fb , n = 0,1,...g results from discrete{time Poisson(λ) arrivals
n
of information sessions, where the session duration is distributed according to the pmf G
and the packet generation rate of an ongoing session is one packet per time slot. Under the
3
enforced assumptions fb , n = 0,1,...g can be identi(cid:12)ed as the busy server process of a
n
discrete{time MjGj1 queue; for this reason the packet arrival process fb , n = 0,1,...g is
n
referred to as the MjGj1 arrival process. The following proposition shows that fb , n =
n
0,1,...g is a correlated process, with time dependencies controlled by the tail of σ [10].
Proposition 2.1 If b is taken to be a Poisson rv with parameter λE[σ], then the process
fb , n = 0,1,...g is a (strictly) stationary ergodic process with the properties:
n
(a) For each n = 0,1,..., the rv b is a Poisson rv with parameter λE[σ];
n
(b) Its covariance function is given by
h i
cov(b ,b ) = λE [σ −j]+ = λE[σ]P[σ > j], n,j = 0,1,...
n+j n e
(c) Its index of dispersion of counts (IDC) is given by
X1 X1 λ
IDC (cid:17) cov(b ,b ) = λE[σ] P[σ > j] = E[σ(σ +1)],
n+j n e
2
j=0 j=0
and the process is short{range dependent (i.e., IDC (cid:12)nite) if and only if E[σ2] is (cid:12)nite.
We now feed this MjGj1 arrival stream fb , n = 0,1,...g into a discrete{time single
n
server queue with in(cid:12)nite bu(cid:11)er capacity. Such a queueing system routinely serves as a
model for a network multiplexer: If q denotes the number of packets remaining in the
n
multiplexer bu(cid:11)er by the end of slot [n−1,n), and the multiplexer output link can transmit
c packets/slot, then the bu(cid:11)er content sequence fq , n = 0,1,...g evolves according to the
n
Lindley recursion
q = 0; q = [q +b −c]+, n = 0,1,... (2)
0 n+1 n n+1
From Part (a) of Proposition 2.1 the average input rate to the multiplexer is E[b ] =
n
λE[σ], and the system is stable if the tra(cid:14)c intensity ρ (cid:17) λE[σ]/c satis(cid:12)es ρ < 1. In that
case q =) q, where the IR{valued rv q is the stationary queue size in the multiplexer
n n
bu(cid:11)er. We are interested in evaluating
P(b,ρ) (cid:17) P [q > b], b (cid:21) 0 and Q (ρ) (cid:17) E [qm], 0 (cid:20) ρ < 1, m = 1,2,
ρ m ρ
i.e., the probability that the stationary queue size exceeds b, and the queue size (cid:12)rst and sec-
ond moments, when the tra(cid:14)c intensity is ρ. To that end we develop simple approximations
that all flow from asymptotic results under heavy and light tra(cid:14)c conditions.
3 Heavy and Light Tra(cid:14)c
The interpolation approximation we have in mind hinges on the availability of explicit ex-
pressions for limits of system quantities as ρ ! 1 (heavy tra(cid:14)c limits), and derivatives with
respect to ρ as ρ ! 0 (light tra(cid:14)c derivatives). It thus requires examination of the behavior
of the queue with MjGj1 arrivals under each one of these two asymptotic regimes.
4
Light Tra(cid:14)c We start with the light tra(cid:14)c regime. The right-hand derivatives at ρ = 0
of the various metrics of interest are evaluated using the Reiman{Simon technique [12]. For
the system to be in the domain of applicability of the Reiman{Simon results, an assumption
on (cid:12)niteness of the exponential moment of σ is needed.
h i
Assumption (A) There exists θ(cid:63) > 0 such that E eθσ < 1 for θ < θ(cid:63).
The detailed light tra(cid:14)c analysis of the queue with MjGj1 arrivals is provided in [15]. We
summarize the conclusions in the following.
Proposition 3.1 Consider the Lindley recursion (2) with integer release rate c = 1,2,...,
and let b = 0,1, .... Under Assumption (A) it holds that
(a) For each n = 0,1,...,c
∂n dn
P(b,0+) = 0 and Q (0+) = 0, m = 1,2,.... (3)
∂ρn dρn m
(b) In addition, for c = 1,
E[σ]2 ∂2 P(b,0+) = Eh[σ −b]+2iP[σ > b]+2 Eh[σ −b]+i2
∂ρ2
h i
−3 E [σ −b]+ P[σ > b]+P[σ > b]2 (4)
d2 E[σ2]
Q (0+) = (5)
dρ2 1 E[σ]
and
!
d2 1 E[σ2]2
Q (0+) = 1+ . (6)
dρ2 2 2 E[σ]2
Proposition (3.1) delineates a light tra(cid:14)c behavior for the queue with MjGj1 arrivals
that is certainly di(cid:11)erent from the one of a classical GIjGIj1 queue. As seen from (3), when
c = 1, in which case the multiplexer can serve no more thanone source per time slot, the (cid:12)rst
derivative of the tail probability is zero. Hence, in a Taylor expansion of P [q > b] around
ρ
ρ = 0 the linear term in ρ would o(cid:11)er no contribution. In contrast, the stationary workload
W in a single server MjGj1 queue is known to satisfy P [W > x] (cid:24) ρ(1−G (x)) (ρ ! 0),
ρ e
that is, in the classical queueing setup the corresponding expansion starts with a non-zero ρ
term. For MjGj1 arrivals it is the second derivative (4) which is the most informative. This
highlights the role of the activity duration rv σ, through both its distribution and its (cid:12)rst
two moments. Notice that even if Assumption (A) were to be relaxed, (4) shows that for
P [q > b] to decay like ρ2 for small ρ it is necessary that E[σ2] be (cid:12)nite. If E[σ2] = 1, as
ρ
is the case for long{range dependent MjGj1 arrivals, expression (4) yields in(cid:12)nity and ρ2 is
no longer the correct order of decay. A di(cid:11)erent, perhaps smaller exponent should be sought
5
in the long{range dependent case (see (21)). Finally, relations (3) reflect (though in a rough
manner) the statistical multiplexing gain: Since the (cid:12)rst non-zero contribution to a the tail
probability is no lower than ρc+1, (3) implies that increasing the multiplexer release rate c
while maintaining the same tra(cid:14)c intensity ρ would result in a decreasing tail probability
P [q > b], as could be expected.
ρ
Heavy tra(cid:14)c Considered next is the behavior of the queue with MjGj1 arrivals in heavy
tra(cid:14)c, that is, as the arrival rate λE[σ] tends to the multiplexer release rate c from below.
Clearly, as the tra(cid:14)c intensity ρ converges to one the system becomes unstable and the
queue length grows unbounded. It is thus necessary to seek a suitable normalizer for the
queue length process, so that its normalized version has a non trivial heavy tra(cid:14)c limit.
This problem is typically addressed in a setup where the system of interest is embedded into
a family of queueing systems, parametrized by an integer, say l = 1,2,..., ensuring that, as
l " 1, the appropriate trend to instability is established. Such an approach was pursued
in [16], providing a complete characterization of the arising heavy tra(cid:14)c limits. We tacitly
assume here that the heavy tra(cid:14)c limit of the stationary distribution coincides with the
stationary distribution of the heavy tra(cid:14)c limit, and gather the required results from [16] in
a convenient form:
Proposition 3.2 The heavy tra(cid:14)c limits of the stationary queue length distribution associ-
ated with (2) can be classi(cid:12)ed as follows:
(a) If E[σ2] < 1 then
!
2E[σ]
limP [(1−ρ) q > x] = exp − x , x (cid:21) 0. (7)
ρ!1 ρ E[σ2]
(b) If P[σ > n] = n−α, n = 1,2,..., with 1 < α < 2, then
!
h i (α−1)E[σ]
limP (1−ρ)1/(α−1)q > x = E − xα−1 , x (cid:21) 0, (8)
ρ!1 ρ α−1 Γ(2−α)
where
X1 xn
E (x) (cid:17) , ν > 0, x 2 IR, (9)
ν
Γ(νn+1)
n=0
is the Mittag{Le(cid:15)er special function [5].
Part (a) of Proposition 3.2 addresses the classical short{range dependent case, for which
theheavytra(cid:14)cnormalizeris(1−ρ)andthelimitingheavytra(cid:14)cdistributionisexponential.
Part (b) deals with a long{range dependent MjGj1 arrival process. Under long{range
dependence, the heavy tra(cid:14)cqueue length distribution is expressed througha Mittag{Le(cid:15)er
function with hyperbolic decay, while the power{law behavior of the heavy tra(cid:14)c normalizer
is (1 − ρ)1/(α−1). We mention that this result can be stated in a more general manner to
cover the situation where the tail of σ is regularly varying of order α, 1 < α < 2.
6
The results under the light and heavy tra(cid:14)c regimes are subsequently combined into
approximations for all values of the tra(cid:14)c intensity.
4 Interpolation Approximations
Whenever Assumption (A) is satis(cid:12)ed, P [q > b] is in(cid:12)nitely di(cid:11)erentiable with respect to
ρ
ρ at ρ = 0, hence it can be approximated by bringing together heavy tra(cid:14)c limits and light
tra(cid:14)c derivatives into a Taylor series{like expansion. To this end we follow the approach
proposed in [6]. In passing, we also discuss approximations for Q (ρ), m = 1,2. The details
m
are as follows:
The interpolation method Consider the normalized queue length rv (1−ρ) q and de(cid:12)ne
F(x,ρ) (cid:17) P [(1−ρ) q > x], 0 (cid:20) ρ < 1, x (cid:21) 0 (10)
ρ
and
F(x,1) (cid:17) limP [(1−ρ) q > x]. (11)
ρ
ρ!1
Assume that partial derivatives of F(x,ρ) with respect to ρ, up to order n, at ρ = 0+,
b
are available. Construct F (x,ρ), the nth order interpolation approximation to F(x,ρ), by
n
means of the polynomial
!
Xn ρi ∂i Xn 1 ∂i
Fb (x,ρ) (cid:17) F(x,0+)+ F(x,1)− F(x,0+) ρn+1. (12)
n i! ∂ρi i!∂ρi
i=0 i=0
Observe that
∂i ∂i
b b
F (x,1) = F(x,1) and F (x,0+) = F(x,0+), i = 0,1,...,n,
n ∂ρi n ∂ρi
b
that is, F (x,ρ) is precisely that unique n + 1 degree polynomial in ρ which matches the
n
n+1 partial derivatives of F(x,ρ) at ρ = 0+ and its heavy tra(cid:14)c limit. Now, by reversing
the (1−ρ) normalization in Fb (x,ρ) we generate the nth order interpolation approximation
n
to P [q > b] as
ρ
P [q > b] (cid:25) Fb ((1−ρ) b,ρ). (13)
ρ n
Note that, in principle, this may lie outside [0,1], in which case it is obviously a poor
approximation. To calculate the quantities associated with (13) it remains to express the
partial derivatives appearing in (12) in terms of the light tra(cid:14)c derivatives of P [q > b]. We
ρ
have
∂ ∂ ∂
F(x,0+) = P(x,0+)+x P(x,0+) (14)
∂ρ ∂ρ ∂x
7
and
∂2 ∂2 ∂2
F(x,0+) = P(x,0+)+2x P(x,0+)
∂ρ2 ∂ρ2 ∂ρ∂x
∂ ∂2
+2x P(x,0+)+x2 P(x,0+). (15)
∂x ∂x2
In case additional light tra(cid:14)c information is available, repeated application of the chain rule
will yield higher order derivatives, as needed.
Approximate expressions We are now ready to write approximate expressions anchored
on the heavy and light tra(cid:14)c results of Section 3. Proposition 3.2(a) provides the limit (11)
that should be inserted in (12). Proposition 3.1(a) can be used to substitute for the partials
in (14) and then in (12). Thus, if the multiplexer release rate is c = 1,2,..., the cth order
interpolation approximation to P [q > b] is simply
ρ
!
2E[σ]
P [q > b] (cid:25) Fb ((1−ρ) b,ρ) = ρc+1exp − (1−ρ) b . (16)
ρ c
E[σ2]
More can be accomplished in the case c = 1, since Proposition 3.1(b) a(cid:11)ords us a promising
2nd order interpolation approximation. Insertion of (15) in (12) yields
!
1 ∂2 2E[σ]
Fb (b,ρ) = ρ2(1−ρ) P(b,0+)+ρ3exp − b (17)
2 2 ∂ρ2 E[σ2]
and the latter leads to the 2nd order approximant
P [q > b] (cid:25) Fb ((1−ρ) b,ρ), c = 1, (18)
ρ 2
where Proposition 3.1(b) is used to supply the second partial derivative in (17).
Next, we briefly deal with moment approximations. We restrict attention to the case
c = 1 and consider only the queue length (cid:12)rst and second moment. The relevant light tra(cid:14)c
limits are given by (3), (5) and (6). In heavy tra(cid:14)c it can be inferred from (7) that
!
E[σ2] m
lim(1−ρ)m Q (ρ) = m! , m = 1,2,....
m
ρ!1 2E[σ]
Moment interpolations are then developed in very much the same manner as distribution
interpolations. We skip the details of the derivation and list the (cid:12)nal expressions
E[σ2] ρ2
Q (ρ) = , c = 1 (19)
1 2E[σ] 1−ρ
and ! !
ρ2 E[σ2]2 E[σ2]2
Q (ρ) (cid:25) ρ −1 + +1 , c = 1. (20)
2 4(1−ρ)2 E[σ]2 E[σ]2
8