#### Corrected asymptotics for a multi-server queue in the Halfin-Whitt regime

A.J.E.M. Janssen
0
1
J.S.H. van Leeuwaarden
0
1
B. Zwart
0
1
0
B. Zwart H. Milton Stewart School of Industrial and Systems Engineering, Georgia Institute of Technology
, 765 Ferst Drive,
30332 Atlanta, USA
1
J.S.H. van Leeuwaarden ( ) Mathematics and Computer Science Department, Eindhoven University of Technology
, 5612 AZ Eindhoven,
The Netherlands
To investigate the quality of heavy-traffic approximations for queues with many servers, we consider the steady-state number of waiting customers in an M/D/s queue as s . In the Halfin-Whitt regime, it is well known that this random variable converges to the supremum of a Gaussian random walk. This paper develops methods that yield more accurate results in terms of series expansions and inequalities for the probability of an empty queue, and the mean and variance of the queue length distribution. This quantifies the relationship between the limiting system and the queue with a small or moderate number of servers. The main idea is to view the M/D/s queue through the prism of the Gaussian random walk: as for the standard Gaussian random walk, we provide scalable series expansions involving terms that include the Riemann zeta function.
1 Introduction
1.1 Goals, results and insights
P(Q = 0) P(M = 0) =
2 r=0 r!(2r + 1)
(1/2 r) 2 r
2
(s) = (2s(1 + ln ))1/2 ,
1.2 Methodology
l=1
1 (P(Al l) P ( l)) ,
l
P (x) =
eu2/2du
Gk(a) =
l=1
e 12 lsx2 z(x)dx,
a, s R+, k Z,
1.3 Organization
2 The M/D/s queue and the Halfin-Whitt regime
Q,n+1 = (Q,n + A,n s)+,
n = 0, 1, . . . ,
f (x)
fn(x),
n=0
k1
n=0
f (x)
fn(x) = fk(x)(1 + o(1)).
Let d denote convergence in distribution.
2.1 Basic results
The following theorem can be proved using a similar approach as in Jelenkovic et
al. [18]. We give a separate proof, because our setting is slightly different and more
specific.
d
(i) Q M ;
(ii) P(Q = 0) P(M = 0);
(iii) E[Q k] E[Mk ] for any k > 0.
Proof Proof of (i): Note that
d
Q = (Q + A )+,
with A = (A )/. Since A converges in distribution to the standard normal
random variable as , (i) follows from Theorem X.6.1 in Asmussen [2], if the
family (A, 0) is uniformly integrable. But this follows simply from the fact that
E[A] = 1 for all .
2
Proof of (ii): The result lim sup P(Q = 0) P(M = 0) follows from (i). To
show the lim inf, note that from Spitzers identity (see (16))
1
1 x2ex , we see that any s is admissible that satisfies
Since ex 1 x 2
s22 es/ s 0.
P(Q > x) ex
1
l=1
P ( l) = ln P(M = 0),
now follows directly using for example the formula
ln P(M = 0) =
P(Sl > 0),
1
l=1
which gives for the M/D/s queue
ln P(Q = 0) =
1
l=1
1
P(Q = 0)
P(M = 0) exp
1
l=1
1
l=1
(P ( l) P(Al > ls))
exp
|P ( l) P(Al > ls)| .
4 30.6 1 4
5 , 1 + 3l3/2 l 5l
yields, upon substituting (21) into (19), the second inequality in (18). The first
inequality in (18) follows in a similar way.
2.2 Quasi-Gaussian form: motivation and outline
which leads to the approximation
1
l=1
l=1
l=1
Theorem 2
in which
l=1
ex2/2y (x/ ls)dx,
2 l
2s 1 s + ln s
p(n) = nnen2 n/n!,
and y is a function analytic in |x| < 2 (see (140)).
pk
k=0
n ,
and for y there is the power series representation
i=0
ln P(Q = 0)
pksk+1/2G(k+1)(/s),
l=1
Gk(a) =
e 21 lsx2 y (x)dx.
Similar expressions, though somewhat more complicated than the one in (29),
exist for EQ and VarQ (see Sect. 3.2) and these involve Gk with k = 0, 1, 2 . . .
and k = 1, 0, 1 . . . , respectively. We shall study Gk thoroughly, leading to series
expansions, asymptotics and bounds.
We close this section by giving an interpretation of the quasi-Gaussian form (24).
Using, see (134),
e 12 nx2 y (x)dx
we find from (24) that
elsx2/2y (x)dx
1
l=1
elx2/2y (x/s)dx
3 From Spitzer-type expressions to quasi-Gaussian forms In this section we show how to obtain the expression (24). In addition, we present similar results for the mean and variance of the queue length.
3.1 Proof of Theorem 2
For n = 0, 1, . . . we let
sn(z) =
k=0
z C.
j=n+1
= 1 el sn(l) =
nn+1en
with p(n) as defined in (27). We then consider the equation
ln P(Q = 0) = s1/2
l=1
1 2 1 3 1 4
f (y) = 2 y + 3 y + 4 y + ,
Then it holds that
Substituting (42) into (37) yields (24).
qls () = e 21 ls
e 12 lsx2 y (x)dx =
Table 1 Interrelations between
some parameters and functions
s = +
= /s
= 2(1 + ln )
= s
a = /s =
y(x) = x 31 x2 + 316 x3 + ; |x| < 2
= s y(/s)(1 y(/s))1/2
= + 61 s1/22 + 752 s13 + ; |/s| < 2
p(n) = nnen2n/n! 1 + 112 n1 + 2818 n2 +
/s = (1 yy((//ss)))1/2 .
= + 6 1s 2 + 752s 3 + .
3.2 Mean and variance of the queue length
This leads after considerable rewriting to
1
1
qls () ls(1 )
s
l=1
l=1
l=1
l=1
ex2/2y (x/ ls)dx ,
2 l
e 12 2l
lR() / ls p(ls)
p(ls)
(2lR2() + )
ex2/2y (x/ ls)dx
2 l
l=1
(ls(1 ) 1)qls ()
+ ((1 )2l2s2 + ls)
For the Gaussian random walk we have that (see [13])
l=1
l=1
l1/2e 12 2l .
4 Results for Gk
4.1 Principal series expansion
We let z : [0, ) C be a continuous function satisfying z(x) = O(exp(x2)) for
any > 0, and we assume that there is an r0 > 0 such that z(x) is represented by
its Taylor series j=0 bj xj for 0 x < r0. We consider for s > 0 and integer k the
function
Gk(a) =
l=1
e 12 lsx2 z(x)dx,
a > 0.
where Tk,i is defined as
Gk(a) = s i+21 Tk,i (as),
Tk,i (b) =
l=1
e 21 lx2 xi dx
Gk(a) =
(k + 3/2)
a z(x)
2k+1 bj aj2k2
j=0
2k + 2 j
j2k=+02 bj xj
b2k+2 ln a
+ Lk
r=0
x2r z(x)dx,
Proof We have
(ln 1/z)t1 + zv
The right-hand side of (60) can be expressed in terms of Lerchs transcendent ,
defined as the analytic continuation of the series
which converges for any real number v = 0, 1, 2, . . . if z and t are any complex
numbers with either |z| < 1, or |z| = 1 and Re(t ) > 1. Note that (t ) = ( 1, t, 1).
Thus,
Gk(a) = z(a)e 12 sa2
We then use the important result derived by Bateman [9], Sect. 1.11(8) (with
(t, v) := ( 1, t, v) the Hurwitz zeta function)
l=1
j=0
Lk =
e 21 lsx2 z(x)
j=0
(k + 3/2)
j=0
Gk(a) = z(a)
lk+1/2e 12 lsa2 .
(z, t, v) =
(v + n)t zn,
l=1
n=0
r=0
r=0
Gk(a) = z(a) (k + 3/2)
a2k3 +
(k r 1/2) ( 12 sa2)r .
r!
Gk(a) + (k + 3/2)
j=0
bj aj2k3
j=0
z(a)
bj aj a2k3
z(a)
r=0
Gk(a) + (k + 3/2)
= Lk (k + 3/2)
2 k+3/2 2k+1 bj aj2k2
j=0
z(x)
bj xj x2k3dx
j=0
x2r z(x)dx,
r=0
Gk(a) =
e 12 lsx2 z(x)
We shall determine Lk. It holds that, as a 0,
2 k+3/2 2k+1 bj aj2k2
l=1
l=1
i=0
+ o(1) +
l=1
i=0
j=0
i=0
i=0
e 21 lsx2 xidx
e 12 lsx2 z(x)
bi s i+21 Tk,i (as).
i2=k+0 2 bixi =
i=0
1
bixi dx = O (ls)k+2 .
Now from Janssen and Van Leeuwaarden [14], Sect. 5,
for i = 0, 1, . . . , 2k + 1 and
= (k + 3/2) s
ln(as)
+ s(k+3/2)Lk,2k+2 + O(a).
Lk,2k+2 = (k + 3/2)2k+3/2
Therefore, as a 0, we get
Gk(a) =
e 21 lsx2 z(x)
i=0
l=1
i=0
2 k+3/2 ai2k2
2k + 2 i
+ b2k+2 (k + 3/2) s
j=0
2 k+3/2 2k+1 bj aj2k2
j 2k 2 + b2k+2 ln a
Table 2 Some values of the
Riemann zeta function
0.50000000000000
1.46035450880959
2.61237534868549
1.64493406684823
1.34148725725092
1.20205690315959
1.12673386731706
1.08232323371114
1.05470751076145
1.03692775514337
1.02520457995469
at either side of (74) and letting a 0, we find that Lk has the required value (59).
Then (58) follows from (66).
Some values of the Riemann zeta function are given in Table 2.
We now give several special cases of Theorem 3. The next two corollaries focus
on negative values of k.
Gk (a) =
(k + 3/2)
r=0
+ Lk
x2r z(x)dx,
x2k3z(x)dx
( 21 s)r
( 21 s)r
G1(a) =
r=0
+ L1
x2r z(x)dx,
where L1 =
Theorem 3 is meant for the case that a and the convergence radius r0 of
j =0 bj xj are general. In the case that a < r0 the expressions can be simplified
considerably, as demonstrated below. If a < r0 we have
a z(x)
j2k=+02 bj xj
dx =
x2r z(x)dx =
j=2k+3
bj aj2k2
j 2k 2
j=0
bj aj+2r+1
simplifies to
2k+1 bj aj2k2
j=0
a z(x)
j2k=+02 bj xj
j=0,j=2k+2
Lemma 3 For the first line of (59)
1 2 (j+1)/2
e 12 lsx2 xj dx = 2 ls
euu(j1)/2du
l=1
j=0
e 12 lsx2 z(x)
j=0
there is the asymptotic expression
j=0,j=2k+2
s .
In case that bj ( j +21 ) = O(Bj ) for some B > 0, the asymptotic series in (83) is
convergent when s > 2B2, with sum equal to (82).
Proof Using
we find that
l=1
j=2k+3
j=2k+3
l=1
e 12 lsx2 z(x)
j=0
e 21 lsx2 xj dx
l=1
lj/2+k.
l=1 lj/2+k .
Remark 5 Chang and Peres [6], Theorem 1.1, proved that
2 r=0 r!(2r + 1)
(1/2 r) 2 r
2
l=1
lkP ( l) =
(k r 21 ) 2 r
r!(2r + 1) 2
where R1() = ln 2 and
k = 1.
4.2 Optimal truncation value
j=0,j=2k+2
We replace, furthermore, bJ +1 = (J + 2)aJ +2 by its asymptotic bound, see the
appendix, Lemma 13,
The right-hand side of (94) decreases in J until J /2 + 1 = 2 s; this J is (near to)
the optimal truncation point. At this point we estimate the right-hand side of (92) by
Stirlings formula as
(2 s)2s1/2e2s 2
is of the order
The factor (1/2(J + 2))1/2 is rather unimportant for determination of the optimal
truncation value J , and we focus on
J + 2 1/2
2(J + 2)
(J +2)/2
. (92)
|bJ +1|
DJ =
(J +2)/2
4.3 Accurate approximations for the M/D/s queue
We can apply Theorem 3 to obtain accurate approximations for the emptiness
probability and the mean and variance of the queue length. By way of illustration, we do
this in some detail for P(Q = 0) and briefly indicate at the end of this section how
one can proceed for the other cases.
We have from (27) and (29) that
ln P(Q = 0)
pksk+1/2G(k+1)(/s)
1
s1/2G1(/s) 12 s1/2G2(/s)
s3/2G3(/s) + .
Accurate approximations to ln P(Q = 0) are obtained by including 1, 2, 3, . . .
terms of the second line of (96) in which the Gs must be approximated. For the
number of terms of the asymptotic series in (96) to be included one could follow
a truncation strategy (based on (141), (154) and the bound Gk(a) ( 2s )1/2 (k),
k = 2, 3, . . .) pretty much as was done in Sect. 4.2. We shall not pursue this point
here.
We shall compute accurate approximations to Gk(a) for k = 1, 2, . . . . We have
from (76) and (77) for < 2
G1(/s) =
and for k = 2, 3, . . . ,
Gk(/s) =
+ L1
+ Lk
r=0
r=0
/s y (x) 1
x2r y (x)dx,
2 k+3/2
( k + 3/2)
x2k3y (x)dx
L1 =
l=1
l1/2
and for k = 2, 3, . . . ,
e 21 lsx2 (y (x) 1)dx
ln 2s,
Lk =
l=1
lk+1/2
e 21 lsx2 y (x)dx.
Below we specify the missing ingredients in (98)(101).
We have
/s y (x) 1
dx =
j=1
bj xj1dx =
bj
j=1
xny (x)dx =
bj xn+j dx =
j=0
j=0
n+j+1
x2r y (x)dx =
j=0
. (104)
Since, see [13], Sect. 6, (r + 21 )/r! = O(1/(2 )r ), the computation of the series
over r at the right-hand side of (98) is feasible when < 2 . A similar result
holds for the series over r at the right-hand side of (99).
We have by Lemma 3
l=1
l1/2
2
1
j=1
j=0
e 21 lsx2 (y (x) 1)dx
k = 2, 3, . . . , (106)
on the second line of (96) (so that the truncation error is O(s7/2)). We then
include in the right-hand side of (105) the terms with j = 1, 2, 3, 4, and in the
righthand side of (106) the terms with j = 1, 2.
5 Bounds and approximations for the emptiness probability
l=1
l=1
e 21 lsx2 y (x)dx
=: LB,
e 12 lsx2 y (x)dx
Lemma 5
P(Q = 0) exp
s1/2
P(Q = 0) exp
s1/2
=: U B.
Proof Follows directly from rewriting (24) as
and applying
s1/2
l=1
e 21 lsx2 y (x)dx
1 nn1/2 en
n! 2
= 0.01 (0.0141)
True (96)-1
= 0.1 (0.1334)
True (96)-1
Lemma 6 There is the inequality
l=1
= 0.2 (0.2518)
True (96)-1
= 0.5 (0.5293)
True (96)-1
Proof Letting = 122 we have
, (113)
where we have set x = . Now
(equality at t = 0) and thus
the proof is complete.
Lemma 8 The following inequalities hold:
(ii) For 0 < < 2,
ULBB exp 1442s .
l=1
l=1
(ii) Follows from rewriting the right-hand side of (112) in terms of Lerchs
transcendent and applying the Bateman series (63).
Lemma 9 There are the inequalities
LB exp
LB exp
l=1
l=1
1 1
1 1
1 2 2
e 2 x dx + 32 s l=1
l 2 l
2
Proof Follows from (108) and 1 3 x y (x) 1 (see Lemma 15).
Lemma 10
U B exp
l=1
1 1
18s2 s l=1
l=1
2
Proof Follows from (109) and 1 3 x y (x) (see Lemma 15).
The right-hand side of (122) can be written as
Proof (i) Follows from (112) and observing that
exp
l=1
1
l1P ( l) + 12
l=1
12s 1
18s2 s
2 exp 1424 + 23
12s 1
18s2 s
( 12 ) ( 32 )
+ 2 122
12s 1
18s2 s
which yields a sharp approximation for (123) for small values of .
Considering the leading component in the exponents of (121) and (122), it makes
sense to use the approximation
P(Q = 0) exp
1
l=1
1 1 2 2
l 2 l e 2 x dx + 32 s
l=1
= P(M = 0) exp
= P(M = 0) exp
l=1
. (126)
P(Q = 0) P(M = 0) 2
P(M = 0) = P(M = 0) exp 32 s
( 12 r)( 12 )r
r!(2r + 1)
r=0
( 32 r)( 12 )r 2r
5.1 Numerical experiments
= 1 (0.8005)
(120) True
= 0.1 (0.1334)
(120) True
6 Conclusions and outlook
The approach in this paper consists of three major steps:
= 0.5 (0.5293)
(120) True
form (24). The key facilitator is, see (42),
e 21 sx2 y (x)dx,
d
Q = (Q + A s)+.
Appendix: Analysis of the function y
We shall present an analysis to obtain some results on the function y(x) that
appeared in Sect. 2.2, especially in (24). As before, y(x) is, for |x| sufficiently small,
the solution y of the equation
(t + 1) =
euut du
(t + 1) = et t t+1
et (v+ln(1+v))dv
= et t t+1
etf (y)dy = et t t+1
e 21 tx2 y (x)dx.
n=1
y(x) =
n=1
converging for |x| sufficiently small. Then by Watsons lemma, see e.g. [23], pp. 112
116, the asymptotics of (t + 1) follows from (134) as
whence we can write (129) for small x and y as
with the principal value of the square root. Hence,
2 1 2
y 1 + 3 y + 2 y +
= x
1 2
y(x) = x 3 x + O(x3).
(t + 1) et t t+1
xn1e 12 tx2 dx
k=0
= et t t+12
(u)t eudu,
substitutions as in (134) we now get
etf (y)dy,
(t ) 2 et t t1
n=1
k=0
k = 0, 1, . . . .
xn1e 21 tx2 dx
e 21 tx2 y (ix)dx
nanin1
(1)k (2k + 1)(2k 1) 1 a2k+1,
nk
n ,
(t + 1) t t+1/2
k=0
ck (k + 3/2)
, t ,
with ck the coefficients of a function ( ) that is regular inside the circle | |1= 2 .
Apparently, since t (t ) = (t + 1) and (k + 3/2) = (k + 1/2)(k 1/2) 2 ,
we have that
a2k+1 = (1)k2k1/2ck,
k = 0, 1, . . . ,
y(x) = 1 + W (e(1+ 21 x2)),
x > 0.
It appears that Szegs analysis of the mapping in (144) has largely escaped the
attention of the Lambert W community. In the proof of the lemma below, we heavily rely
on this analysis and omit some of the details that are contained implicitly or explicitly
in [28], Sect. 2.
sin
cos 1 ln = tan 2 .
Fig. 2 Curves C in the
y-plane. = 0.1
lim
:= |y|,yC
Also, |x| = |2f (y)|1/2 in these cases. We thus conclude that y(x) can be
continued analytically along all rays x = rei , r 0, when [, ], except for
= /4. In the latter cases the analytic continuation can only be carried out until
r = 2 .
Lemma 12 We have
y(x) =
n=1
(n + 1)an+1ak+2n ,
k = 0, 1, . . . .
y (x)y(x) = x xy(x).
k1
n=0
(n + 1)an+1akn = ak1,
k = 2, 3, . . . ,
and this gives (150). Also see [23], Ex. 8.3 on p. 88 and [31], p. 70.
The first five coefficients an are given by
a1 = 1,
Lemma 13 There is the asymptotic form (rapid decay to 0)
an =
sin 4 (n 1) 6n29 cos 4 (n 1)
1 + O(1/n) .
and similarly,
Fig. 3 The ratio of an and the
approximation in (154),
minus 1, for the first 40
coefficients
Note that the leading term at the right-hand side of (158) vanishes when n = 4k + 1,
k = 0, 1, . . . , and so more precise information is required. This can be obtained when
we write
n=1
an(ei/4(4 v2)1/2)n,
Lemma 14 It holds that
y(x) = 1
mm1
zez = u,
so that by the Brmann-Lagrange inversion formula
z(u) =
1
m=1 m!
m1
z=0
um =
m=1
mm1 um.
An interesting consequence of (161) is that
e 21 lsx2 y (x)dx =
m=1
m=1 m!em
mme 21 (m+ls)a2
m=1 m!em(m + ls)
e 21 (m+ls)x2 xdx
This series expansion can be then be inserted into (55) to yield
Gk(a) =
m=1 m!
e 21 lsa2 .
l=1 m + ls
Proof Solve f(y) = y ln(1 y) = 21x2 for the unique solution y(x) [0,1).
From the inequalities, to be proved below,
and monotonicity of f we infer that
From (151) we get
Inserting the inequalities (170) into (171), we find
1
y (x) = x y(x) 1 , x 0.
as required. We now prove the two inequalities in (169). As to the second one, we
need to show that
Lemma 15 There is the inequality
2
1 3x y (x) 1, x 0.
1 3
g2(x) = 1 + x 1 + 2x = h2(x), x 0,
we get (173). Note that the second inequality in (169) is the sharpest in its kind: when
> 0, the inequality y(x) x/(1+(1)x) fails to hold for large x 0, since this
would require
x 0, (177)