Many-server Gaussian limits for overloaded non-Markovian queues with customer abandonment

Abstract

Extending Ward Whitt’s pioneering work “Fluid Models for Multiserver Queues with Abandonments, Operations Research, 54(1) 37–54, 2006,” this paper establishes a many-server heavy-traffic functional central limit theorem for the overloaded $G{/}GI{/}n+GI$ queue with stationary arrivals, nonexponential service times, n identical servers, and nonexponential patience times. Process-level convergence to non-Markovian Gaussian limits is established as the number of servers goes to infinity for key performance processes such as the waiting times, queue length, abandonment and departure processes. Analytic formulas are developed to characterize the distributions of these Gaussian limits.

Appendices

Appendix

Additional Proofs

1.1 Proof of Proposition 3.1

To make sure (3.11) is well defined and to be able to characterize its distribution, we define a sequence of discrete versions of (3.11), that is, $\{L^{(m)}:m\ge 1\}$, where

$$\begin{aligned} L^{(m)}(t)\equiv \sum _{i=0}^{m-1} J(t,u_i) \left( Z(\omega ,u_{i+1}) - Z(\omega , u_{i}) \right) , \quad t \ge 0, \end{aligned}$$

(A.1)

for a given partition on the interval [0, t], $0=u_0< u_1< \cdots < u_m = t$. Suppose $\omega \in \varOmega$ is such that $(Z(\omega ,t):t \ge 0)$ has Hölder continuous sample paths. For simplicity, we suppress $\omega$ hereafter. For $m > 0$, consider the partition $0=u_0< u_1< \cdots < u_m=t$ and

$$\begin{aligned} L^{(m)}(t) =&\sum _{i=0}^{m-1} J(t,u_i) \left( Z(u_{i+1}) - Z(u_i) \right) \nonumber \\ =&\sum _{i=1}^{m} J(t,u_{i-1}) Z(u_i) - \sum _{i=0}^{m-1} J(t,u_i) Z(u_i) \nonumber \\ =&J(t,u_{m-1}) Z(u_m) - J(t,0) Z(0) \nonumber \\&- \sum _{i=1}^{m-1} \left[ \, J(t,u_i) - J(t,u_{i-1}) \, \right] Z(u_i). \end{aligned}$$

(A.2)

Because Z(t) is continuous, the summation converges to the Riemann–Stieltjes integral as the partition mesh goes to 0 if J(t, u) is monotone in the second component for each t. Moreover, if J(t, u) is differentiable for each t, we can replace the integrator $\mathrm{d}J(t,u)$ of the Riemann–Stieltjes integral with $J_u(t,u)\mathrm{d}u$, where the subscript denotes derivative with respect to the second component. The Riemann–Stieltjes integral is well defined if the derivative as a function of u for fixed t is continuous. (In general, finitely many jumps are allowed.) Therefore,

$$\begin{aligned} \int _0^t J(t,u)\,\mathrm{d}Z(u) \equiv&\lim _{m \rightarrow \infty } L^{(m)}(t) {\mathop {=}\limits ^\mathrm{a.s.}}J(t,t) Z(\omega ,t) - J(t,0) Z(\omega ,0) \\&- \int _0^t Z(\omega ,u)\,\mathrm{d}J(t,u). \end{aligned}$$

Moreover, with $\varDelta \equiv \max \{u_i-u_{i-1}: 1 \le i \le m\}$, we have

$$\begin{aligned}&\sum _{i=1}^{m-1} \left[ \, J(t,u_i)-J(t,u_{i-1})\, \right] Z(u_i) - \int _0^t J_u(t,u)Z(u)\,\mathrm{d}u \\&\quad = \sum _{i=1}^{m-1} \int _{u_{i-1}}^{u_i} \left[ \, \frac{J(t,u_i)-J(t,u_{i-1}) }{u_i-u_{i-1}}\, - J_u(t,u)\, \right] \left( Z(u)-Z(u_i) \right) \,\mathrm{d}u \quad \\&\quad \le \frac{1}{\varDelta } \cdot c_1 \varDelta \cdot c_2 \varDelta ^\alpha \rightarrow 0 \end{aligned}$$

as $\varDelta \rightarrow 0$, where the inequality holds because $\hat{Z}$ has Hölder continuous sample paths and J(t, u) is differentiable with respect to the second component.

We prove Proposition 3.1 in two steps. First, we show in Lemma A.1 that if the sequence of covariance functions associated with the processes $\{L^{(m)}:m \ge 1\}$ converges to some limit function, then the sequence $\{L^{(m)}:m \ge 1\}$ converges in distribution to a Gaussian process. Moreover, the covariance function of the limit Gaussian process coincides with the limit of the covariance function associated with $\{L^{(m)}:m \ge 1\}$. Then, in the second step, we show that the covariance functions associated with $\{L^{(m)}:m \ge 1\}$ indeed converge.

Lemma A.1

Let $X^{(m)} \equiv \left( X^{(m)}_1,\ldots ,X^{(m)}_l \right)$ be a sequence of centered Gaussian random vector in ${\mathbb {R}}^l$ and let $\varSigma ^{(m)}$ be the covariance matrix of $X^{(m)}$. If $\varSigma ^{(m)} \rightarrow \varSigma$ as $m \rightarrow \infty$, then $X^{(m)} \Rightarrow X$, where the limit X is Gaussian with mean zero and covariance $\varSigma$.

Proof

Consider the characteristic function $\phi _m(\theta ) \equiv \mathbb {E}\left[ e^{ i \theta ^T X^{(m)} } \right]$ of the vector $X^{(m)}$. The convergence $\varSigma ^{(m)} \rightarrow \varSigma$ implies the convergence of characteristic functions

$$\begin{aligned} \phi _m(\theta ) = e^{ -\frac{1}{2} \theta ^T \varSigma ^{(m)} \theta } \rightarrow \phi (\theta ) \equiv e^{ -\frac{1}{2} \theta ^T \varSigma \theta } \end{aligned}$$

due to continuity of $\phi _m$. Then the result follows from Lévy’s continuity theorem. $\square$

We next show that the covariance functions associated with the sequence $\{L^{(m)}:m \ge 1\}$ in (A.1) converge. We consider a partition of the interval $[0,t_2]$ such that there are a total of $m_2$ intervals partitioning $[0,t_2]$ and $m_1$ intervals partitioning $[0,t_1]$. We use the form in (A.2) to compute the covariance of $L^{(m)}(t)$. Let $C_Z(\cdot ,\cdot )$ be the covariance function associated with the process Z. Then, for $0 \le t_1 < t_2$ and the partition $0=s_0<s_1<\cdots<s_{m_1-1}<s_{m_1}=t_1<s_{m_1+1}<\cdots<s_{m_2-1}<s_{m_2}=t_2$,

$$\begin{aligned}&\mathbb {E}[L^{(m)}(t_1)L^{(m)}(t_2)] \\&\quad =\mathbb {E}\left[ \left( J(t_1,s_{m_1-1})Z(t_1) - J(t_1,0)Z(0) - \sum _{i=1}^{m_1-1} (J(t_1,s_i)-J(t_1,s_{i-1})) Z(s_i) \right) \right. \\&\qquad \times \left. \left( J(t_2,s_{m_2-1}) Z(t_2) - J(t_2,0) Z(0) - \sum _{i=1}^{m_2-1} (J(t_2,s_i)-J(t_2,s_{i-1})) Z(s_i) \right) \right] \\&\quad = J(t_1,s_{m_1-1})J(t_2,s_{m_2-1})C_Z(t_1,t_2) + J(t_1,0)J(t_2,0)C_Z(0,0) \\&\qquad - J(t_2,s_{m_2-1})J(t_1,0)C_Z(0,t_2) - J(t_1,s_{m_1-1})J(t_2,0)C_Z(0,t_1) \\&\qquad - \sum _{i=1}^{m_1-1} J(t_2,s_{m_2-1})(J(t_1,s_i)-J(t_1,s_{i-1})) C_Z(s_i,t_2) \\&\qquad + \sum _{i=1}^{m_1-1} J(t_2,0)(J(t_1,s_i)-J(t_1,s_{i-1})) C_Z(0,s_i) \\&\qquad - \sum _{i=1}^{m_2-1} J(t_1,s_{m_1-1})(J(t_2,s_i)-J(t_2,s_{i-1})) C_Z(t_1,s_i) \\&\qquad + \sum _{i=1}^{m_2-1} J(t_1,0)(J(t_2,s_i)-J(t_2,s_{i-1})) C_Z(0,s_i) \\&\qquad + \sum _{i=1}^{m_1-1} \sum _{j=1}^{m_2-1} (J(t_1,s_i)-J(t_1,s_{i-1})) (J(t_2,s_j)-J(t_2,s_{j-1})) C_Z(s_i,s_j). \end{aligned}$$

Convergence of the first four terms follows from continuity of $u \mapsto J(t,u)$ for each fixed t as $s_{m_1-1}\rightarrow t_1$ and $s_{m_2-1}\rightarrow t_2$ as $m \rightarrow \infty$. Convergence of the last four summations follows from the fact that $C_Z$ is bounded over compact intervals and J(t, u) is differentiable and, therefore, bounded for each t over compact intervals. Hence the limits of these terms are the Riemann–Stieltjes integrals given in (3.13). Finally, the last summation term converges to the two-dimensional Riemann–Stieltjes integral in (3.13) due to similar reasoning. $\square$

1.2 Proof of Theorem 4.1

We first establish a FWLLN for $W_n$ following the compactness approach, i.e., (i) the sequence $W_n$ is $\mathcal {C}$-tight, which implies that every subsequence has a convergent subsequence with a limit in $\mathcal {C}$; and (ii) every convergent subsequence converges to the same limit, which in our case uniquely solves the ODE in (4.2). Finally, we establish convergence for the other processes and characterize their limits. We remark that the tightness for $W_n$ is quite straightforward, but the tightness for the CLT-scaled processes (for example, $\widehat{W}_n$) is complicated (which is why we adopt a new approach to prove the FCLT).

The proof closely follows the arguments in [24] and Sect.  6.6 of [26]. We, hereby, redo the steps therein for the new representation of the enter-service process $E_n$; we use the decomposition in (3.3)–(3.6) that is different than the expressions for $E_n$ in [26].

Tightness of $\{ W_n \}$ To prove tightness, first we show that $W_n$ is stochastically bounded and then show that $W_n$ has controlled modulus of continuity, that is, for each $T>0$ and $\epsilon >0$,

$$\begin{aligned} \lim _{\delta \downarrow 0} \limsup _{n \rightarrow \infty } \mathbb {P}(w(W_n,\delta ,T)>\epsilon )=0, \end{aligned}$$

(A.3)

where $w(W_n,\delta ,T)$ is the modulus of continuity of $W_n$, i.e., $\sup \{ w(W_n,[t_1,t_2]): 0\le t_1<t_2 \le (t_1+\delta ) \wedge T \}$ with $w(W_n,A) \equiv \sup \{ W_n(s_1)- W_n(s_2): s_1,s_2 \in A \}$.

The stochastic boundedness is obvious, because, in any finite interval [0, T], we immediately see that HOL satisfies $0 \le W_n(t) \le T$ for all $n \ge 1$, $t \in [0,T]$.

To treat the modulus of continuity, we first see that $W_n(t+\delta )-W_n(t) \le \delta$ for $\delta >0$ and $0\le t\le T$, because the HWT can increase at most at rate 1. Therefore, it remains to find a bound on $W_n(t)-W_n(t+\delta )$. To this end, define

$$\begin{aligned} \bar{E}_{n,3}(t,\delta )&\equiv \bar{E}_{n,3}(t+\delta )-\bar{E}_{n,3}(t) = \int _{t-W_n(t)}^{t+\delta -W_n(t+\delta )}F^c(V_n(s)) \lambda \, \mathrm{d}s. \end{aligned}$$

(A.4)

Because the ccdf $F^c(x)>0$ for all $x \ge 0$, let $c\equiv \inf _{x \in [0,T]} \{ F^c(x) \} >0$. Hence, the integrand in (A.4) is bounded below by a constant $c\lambda >0$, which yields a lower bound on $\bar{E}_{n,3}(t,\delta )$:

$$\begin{aligned} W_n(t)-W_n(t+\delta ) + \delta \le \frac{\bar{E}_{n,3}(t,\delta )}{c\lambda }, \quad t \ge 0. \end{aligned}$$

From the FCLT in Theorem 2 of [35], we know that $\bar{D}_n(t) \Rightarrow D(t) = \mu t$ so that $\bar{E}_{n,3}(t) \rightarrow E_3(t) = D(t) = \mu t$. Therefore, we have $\limsup _{n \rightarrow \infty } \{ W_n(t)-W_n(t+\delta ) \} \le \mu \delta /c\lambda$ so that

$$\begin{aligned} \limsup _{n \rightarrow \infty } \left| W_n(t+\delta ) -W_n(t)\right| \le c^* \delta , \qquad c^* \equiv \max \left( \mu /c\lambda , 1\right) . \end{aligned}$$

(A.5)

Hence, $W_n$ is tight. In addition, (A.5) also implies that the limit of every convergent subsequence of $W_n$ is in $\mathcal {C}$ and is Lipschitz continuous.

Limit of Convergent Subsequence of $\{ W_n \}$ The C-tightness implies that every subsequence of $W_n$ has a convergent subsequence. Let $W_{n_k}$ be a convergent subsequence with the limit $w^{*}$, i.e., $W_{n_k} \Rightarrow w^{*}$. From (2.3) and (2.4), we deduce that the PWT on the subsequence also converges, that is, $V_{n_k} \Rightarrow v^{*}$, with the limit $v^*$ satisfying

$$\begin{aligned} v^*(t)= & {} w^*(t+v^*(t)) \qquad \text {and}\qquad v^*(t-w^*(t)) = w^*(t), \quad t \ge 0. \end{aligned}$$

(A.6)

We now show that $w^*$ solves the ODE (4.2). On the one hand, the FCLT in Theorem 2 of [35] implies that $(\bar{E}_n, \bar{D}_n) \Rightarrow (D, D)$ with $D(t)=\mu t$. On the other hand, (3.3) implies that $\bar{E}_{n}$ along the subsequence associated with $W_{n_k}$ and $V_{n_k}$ converges to a limit $E^*$. Specifically,

$$\begin{aligned} \bar{E}_{n_k}(t)\Rightarrow & {} E^*(t) = E^{*}_3(t) \equiv \int _{0}^{t-w^*(t)} F^c(v^{*}(s))\lambda \, \mathrm{d}s = D(t) = \mu t. \end{aligned}$$

(A.7)

Because the prelimit process is $\mathcal {C}$-tight, we know the derivative $\dot{w}^*(t)$ exists. Taking derivative in (A.7) yields

$$\begin{aligned} \mu&= (1 - \dot{w}^*(t))F^c(v^{*}(t-w^{*}(t))) \lambda = (1 - \dot{w}^*(t)) F^c(w^{*}(t))\lambda , \end{aligned}$$

(A.8)

which coincides with the ODE (4.2).

FWLLN for the other processes To prove full convergence of $V_n$, we write

$$\begin{aligned}&|V_n(t-W_n(t))-v(t-w(t))| \nonumber \\&\quad \le |V_n(t-W_n(t))-V_n(t-w(t))| +|V_n(t-w(t))-v(t-w(t))| \nonumber \\&\quad = | W_n(t)-w(t)+O(1/n) | + |w(t)+O(1/n)-w(t)| \nonumber \\&\quad \le | W_n(t)-w(t) | +O(1/n) \end{aligned}$$

(A.9)

Apply the change of variable to (A.9) with $u_n \equiv t-W_n(t)$ and $u \equiv t-w(t)$ to obtain

$$\begin{aligned} \Vert V_n -v \Vert \le \frac{ \Vert W_n - w \Vert }{\gamma } + O(1/n)=O(1/n) \end{aligned}$$

(A.10)

for a constant $\gamma >0$, where the equality holds because $u_n=u+o(1)$.

The limit of the sequences of processes (3.8)–(3.10) can be obtained the same way it is done in [26], which makes use of Theorem 3.1. of [31] and then applies the continuous mapping theorem given $W_n \Rightarrow w$. From (6.17) of [26], we immediately write

$$\begin{aligned} \bar{Q}_{n,i} \Rightarrow 0 \quad \text{ for }\quad i=1,2; \quad \bar{Q}_{n,3} \Rightarrow Q_{3}(t) \equiv \int _{t-w(t)}^{t}F^c(t-s)\lambda \,\mathrm{d}s,\quad \text {as}\quad n \rightarrow \infty . \end{aligned}$$

(A.11)

1.3 Proof of Theorem 4.3

The expressions for $C_{\widehat{W}_1}(t,t')$ and $C_{\widehat{W}_2}(t,t')$ are obtained by applying rules of the Itô integral. Derivation of these functions follows from standard arguments and therefore the details are omitted.

To compute $C_{\widehat{W}_3}(t,t')$ we make use of (3.13) with $J(t,u)\equiv H(t,u)/q(u,w(u))$, where H(t, u) and q(u, w(u)) are as in (4.17). In particular, for $0 \le t < t'$,

$$\begin{aligned} C_{\widehat{W}_3}(t,t') =&J(t,t) J(t',t') C_{{E}}(t,t') - \int _0^{t'} J(t,t)J_u(t',u)C_{{E}}(t,u)\mathrm{d}u \\&- \int _0^{t} J(t',t') J_u(t,u)C_{{E}}(t',u)\mathrm{d}u \\&+ \int _0^{t}\int _0^{t'} J_u(t,u) J_v(t',v) C_{{E}}(u,v)dv \mathrm{d}u \\ =&\frac{1}{\lambda ^2 F^c(w(t))F^c(w(t'))}\, C_{{E}}(t,t') - \frac{1}{\lambda F^c(w(t))} \int _0^{t'} J_u(t',u)C_{{E}}(t,u)\mathrm{d}u \\&- \frac{1}{\lambda F^c(w(t))} \int _0^{t} J_u(t,u)C_{{E}}(t',u)\mathrm{d}u \\&+ \int _0^{t}\int _0^{t'} J_u(t,u) J_v(t',v) C_{{E}}(u,v)\mathrm{d}v \mathrm{d}u, \end{aligned}$$

where $J_u(t,u)$ is as in (4.20).

We next derive the covariance function for the limit queue-length process. First, $C_{\widehat{Q}_1}(t,t')$ can be obtained from isometry property of the Itô integral. The function $C_{\widehat{Q}_2}(t,t')$ can be obtained by $\widehat{U}(\lambda s,y) = \mathcal {W}(\lambda s,y)- y \mathcal {W}(\lambda s,1)$, where $\mathcal {W}(\cdot .\cdot )$ is a two-dimensional Brownian motion. We refer interested readers to the long version of [31] and the references therein for a definition of the Kiefer process and of stochastic integrals with respect to two-parameter martingales. The last term easily follows by definition. $\square$

1.4 Proof of Lemma 4.1

First, we prove the existence of $\widetilde{E}(t)$. It suffices to show that for any $n\ge 1$ and $-t_1<-t_2<...<-t_n\le 0$, the matrix $M=(\widetilde{C}(-t_i,-t_j))^n_{i,j=1}$ is nonnegative definite. Let $r_1=t_1$ and $r_j=t_{j-1}-t_{j}$ for $j=2,...,n$, and define $N=(C_E(r_i,r_j))^n_{i,j=1}$. For any $z=(z_1,z_2,...,z_n)^T\in \mathbb {R}^n$, define $y=(y_1,y_2,..,y_n)^T$ such that $y_1=\sum _{i=1}^nz_i$ and $y_j=-\sum _{i=j}^n z_i$ for $j=2,...,n$. Given that $\widetilde{C}(t,s)=C_E(-t,-t)-C_E(-t,s-t)$, we can compute

$$\begin{aligned} z^TMz=\sum _{i=1}^n\widetilde{C}(-t_i,-t_i)z_i^2+2\sum _{1\le i<j\le n}\widetilde{C}(-t_i,-t_j)z_iz_j=y^TNy. \end{aligned}$$

We shall explain how to derive the above equation for $n=2$.

$$\begin{aligned}&z_1^2\widetilde{C}(-t_1,-t_1)+2z_1z_2\widetilde{C}(-t_1,-t_2)+z_2^2\widetilde{C}(-t_2,-t_2)\\&\quad =z_1^2C_E(r_1,r_1)+2z_1z_2(C_E(r_1,r_1)-C_E(r_1,r_2))+z_2^2C_E(r_1-r_2,r_1-r_2)\\&\quad =z_1^2C_E(r_1,r_1)+2z_1z_2(C_E(r_1,r_1)-C_E(r_1,r_2))\\&\qquad +z_2^2(C_E(r_1,r_1)-2C_E(r_1,r_2)+C_E(r_2,r_2))\\&\quad =(z_1+z_2)^2C_E(r_1,r_1)-2z_2(z_1+z_2)C_E(r_1,r_2)+z_2^2C_E(r_2,r_2) =y^T Ny. \end{aligned}$$

Since $C_E$ is the covariance function of a Gaussian process, the matrix N is nonnegative definite and hence $y^TNy\ge 0$. As the vector z is any vector in $\mathbb {R}^n$, we can conclude that M is also nonnegative definite and the existence of $\widetilde{E}$ follows. The argument is similar for $n\ge 3$, therefore, the details are omitted.

Next we show that (4.24) holds. Since a Gaussian process is fully characterized by its covariance function, it suffices to show that, for any fixed $t>0$ and $0\le r<s\le t$,

$$\begin{aligned} \mathrm{Cov}(\widetilde{E}(-t+r)-\widetilde{E}(-t),\widetilde{E}(-t+s)-\widetilde{E}(-t))=C_E(r,s). \end{aligned}$$

By our definition of $\widetilde{C}(t,s)$, we can compute

$$\begin{aligned}&\mathrm{Cov}(\widetilde{E}(-t+r)-\widetilde{E}(-t),\widetilde{E}(-t+s)-\widetilde{E}(-t)) \nonumber \\&\quad =C_E(t-r,t-r)-C_E(t-r,s-r)+C_E(t,s)+C_E(t,r)-C_E(t,t). \end{aligned}$$

(A.12)

By the stationary increments of $\widehat{E}$, we have

$$\begin{aligned} C_E(t-r,t-r)&=\mathrm{Var}(\widehat{E}(t-r)) = \mathrm{Var}(\widehat{E}(t)-\widehat{E}(r))\\&=C_E(t,t)-2C_E(t,r)+C_E(r,r),\\ C_E(t-r,s-r)&= \mathrm{Cov}(\widehat{E}(t-r),\widehat{E}(s-r))= \mathrm{Cov}(\widehat{E}(t)-\widehat{E}(r),\widehat{E}(s)-\widehat{E}(r))\\&=C_E(t,s)-C_E(t,r)-C_E(r,s)+C_E(r,r), \end{aligned}$$

which along with (A.12) implies that

$$\begin{aligned} \mathrm{Cov}(\widetilde{E}(-t+r)-\widetilde{E}(-t),\widetilde{E}(-t+s)-\widetilde{E}(-t))=C_E(r,s). \end{aligned}$$

This completes the proof.$\square$

1.5 Proof of Theorem 4.4

Steady state of $\widehat{W}$ Let $\mathcal {N}(0,\sigma ^2)$ denote the normal distribution with mean 0 and variance $\sigma ^2$. First, we treat $\widehat{W}_1(t)$ in (4.16) by applying a change of variable with $u=s+v(s)$. Let $\kappa (t)$ be the inverse of the function $\beta (t) = t+v(t)$. We write

$$\begin{aligned} \widehat{W}_1(t)&= \int _{\kappa (0)}^{\kappa (t)} \frac{F^c(w(s+v(s))) H(t,s+v(s))}{\lambda F^c(w(t))} c_{\lambda }\,\mathrm{d}{\mathcal {B}}_{\lambda }(\varLambda (s+v(s)-w(s+v(s)))) \\&= \int _{0}^{\kappa (t)} \frac{F^c(v(s)) H(t,s+v(s))}{\lambda F^c(w(t))} c_{\lambda }\,\mathrm{d}{\mathcal {B}}_{\lambda }(\varLambda (s))\\&{\mathop {=}\limits ^\mathrm{d}}\int _{0}^{\kappa (t)} \frac{c_{\lambda }}{\sqrt{\lambda }} \, e^{-h_F(w)(t-s-v(s))} \mathrm{d}{\mathcal {B}}_{\lambda }(s) \\&{\mathop {=}\limits ^\mathrm{d}}\widetilde{{\mathcal {B}}}_{\lambda }\left( \frac{c^{2}_{\lambda }}{\lambda }\int _0^{\kappa (t)} e^{-2h_{F}(w)(t-s-v(s)) }\,\mathrm{d}s \right) \\&{\mathop {=}\limits ^\mathrm{d}}\widetilde{{\mathcal {B}}}_{\lambda }\left( \frac{c^{2}_{\lambda }}{\lambda }\int _0^{t} e^{-2h_{F}(w)(t-s) }\,\mathrm{d}\kappa (s) \right) {\mathop {=}\limits ^\mathrm{d}}\widetilde{{\mathcal {B}}}_{\lambda }\left( \frac{c^{2}_{\lambda }}{2h_{F}(w)\lambda }\left( 1-e^{-2t\, h_{F}(w)}\right) \right) \\&\Rightarrow \widehat{W}_1(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}\left( 0,\frac{c_{\lambda }^2}{2h_{F}(w)\lambda }\right) ,\qquad \text {as}\quad t\rightarrow \infty , \end{aligned}$$

where the second equality follows from (4.3). Similarly, an application of Theorem 3.4.6 of [16] yields

$$\begin{aligned} \widehat{W}_2(t)&= \int _{0}^{t} \frac{\sqrt{F(w)}}{ \sqrt{\lambda F^c(w)}} \, e^{-h_F(w)(t-u)} \mathrm{d}{\mathcal {B}}_{a}(u) + o(1) \\&{\mathop {=}\limits ^\mathrm{d}}\widetilde{{\mathcal {B}}}_{a}\left( \frac{F(w)}{ 2 \lambda f(w)}\left( 1-e^{-2t\,h_{F}(w)}\right) \right) \\&\Rightarrow \widehat{W}_2(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}\left( 0,\frac{F(w)}{ 2 \lambda f(w) }\right) ,\quad \quad \text{ as }\quad t \rightarrow \infty . \end{aligned}$$

Next, (4.22) and (4.16) imply that

$$\begin{aligned} \mathrm{Var}(\widehat{W}_3(t))&= \frac{1}{\lambda ^2 F^c(w)^2} \mathrm{Var}\left( \int _{0}^{t} e^{-h_{F}(w)(t-u)}\, \mathrm{d}\widehat{E}(u) \right) \nonumber \\&=\frac{1}{\lambda ^2 F^c(w)^2}\mathrm{Var}\left( -e^{- h_{F}(w)t }\widetilde{E}(-t)-\int _0^t h_{F}(w) e^{-h_{F}(w)s}\widetilde{E}(-s)\mathrm{d}s \right) \nonumber \\&=\frac{1}{\lambda ^2 F^c(w)^2}\bigg [e^{-2 h_F(w)t}\widetilde{C}(-t,-t)\nonumber \\&\quad + 2h_F(w)e^{-h_F(w)t}\int _0^t e^{-h_F(w)s}\widetilde{C}(-s,-t)\mathrm{d}s \nonumber \\&\quad + 2h_F(w)^2 \int _0^t\int _0^x e^{-h_F(w)(x+y)}\widetilde{C}(-x,-y)dy \mathrm{d}x\bigg ], \end{aligned}$$

(A.13)

where the second equality follows from (4.24). Note that $\widetilde{C}(-t,-s)= \mathrm{Cov}(\widetilde{E}(-s),\widetilde{E}(-t))\le \sqrt{\mathrm{Var}(\widetilde{E}(t))\mathrm{Var}(\widetilde{E}(s))}=\sqrt{\mathrm{Var}(\widehat{E}(t))\mathrm{Var}(\widehat{E}(s))}$. As $\mathrm{Var}(\widehat{E}(t))=O(t^2)$ as $t\rightarrow \infty$, we can conclude that $\widetilde{C}(-t,-s)=O(st)$ as $s,t\rightarrow \infty$. As a result, the first term in (A.13) is $O(e^{-2h_F(w)t } t^2) \rightarrow 0$ and the second term is $O(e^{-h_F(w)t}t^2)\rightarrow 0$ as $t\rightarrow \infty$. Hence, we conclude that

$$\begin{aligned} \widehat{W}_3(t) \Rightarrow \widehat{W}_3(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}(0, \sigma ^{2}_{W_3}), \quad \text{ as }\quad t\rightarrow \infty , \end{aligned}$$

where

$$\begin{aligned} \sigma ^{2}_{W_3} \equiv 2 \frac{h_F(w)^2}{\lambda ^2 F^c(w)^2}\int _{0}^\infty \int _0^x e^{-h_F(w)(x+y)}\widetilde{C}(-x,-y)\mathrm{d}y\mathrm{d}x, \end{aligned}$$

(A.14)

and $\widetilde{C}(\cdot ,\cdot )$ is as defined in Lemma 4.1.

Finally, by independence, we conclude that

$$\begin{aligned} \widehat{W}(t)&\Rightarrow \widehat{W}(\infty )\equiv \widehat{W}_1(\infty )+\widehat{W}_2(\infty )+\widehat{W}_3(\infty ){\mathop {=}\limits ^\mathrm{d}}\mathcal {N}\left( 0,\sigma ^2_{W}\right) ,\quad \text {as}\quad t\rightarrow \infty ,\\&\qquad \text {where}\quad \sigma ^2_{W} \equiv \frac{c_{\lambda }^2}{2h_{F}(w)\lambda }+\frac{F(w)}{2\lambda f(w)}+\sigma ^2_{W_3}. \end{aligned}$$

Steady state of $\widehat{Q}$ We next characterize the steady state for the queue-length process.

$$\begin{aligned} \widehat{Q}_1(t)&= c_{\lambda }\int _{t-w}^{t}F^c(t-s)\mathrm{d}{\mathcal {B}}_{\lambda }(\varLambda (s)) {\mathop {=}\limits ^\mathrm{d}}c_{\lambda } \sqrt{\lambda } \int _{0}^{w}F^c(w-s)\mathrm{d}{\mathcal {B}}_{\lambda }( t-w+s)\nonumber \\&{\mathop {=}\limits ^\mathrm{d}}c_{\lambda } \sqrt{\lambda } \int _{0}^{w}F^c(w-s)\mathrm{d}{\mathcal {B}}_\lambda (s) \Rightarrow \widehat{Q}_1(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}(0, \sigma ^{2}_{Q_1}), \quad \text{ as }\quad t \rightarrow \infty , \quad \nonumber \\&\qquad \text {where } \quad \sigma ^2_{Q_1} \equiv \lambda \, c_{\lambda }^2 \int _{0}^{w} F^c(u)^2 \mathrm{d}u. \end{aligned}$$

(A.15)

Next, the expression in (5.18) implies that, as $t\rightarrow \infty$,

$$\begin{aligned} \widehat{Q}_2(t) {\mathop {=}\limits ^\mathrm{d}}{\mathcal {B}}_a \left( \int ^{w(t)}_0 F^c(u) F(u)\lambda \,\mathrm{d}u \right) \Rightarrow \widehat{Q}_2(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}(0, \sigma ^{2}_{Q_2}), \end{aligned}$$

where $\sigma ^2_{Q_2} \equiv \lambda \int _{0}^{w} F(u) F^c(u)\, \mathrm{d}u$. Finally, (4.10) yields that

$$\begin{aligned} \widehat{Q}_{3}(t)&= \lambda \,F^c(w) \,\widehat{W}(t) \Rightarrow \widehat{Q}_3(\infty )\equiv \lambda \, F^{c}(w) \,\widehat{W}(\infty ) \\&{\mathop {=}\limits ^\mathrm{d}}\mathcal {N}\left( 0, \sigma ^2_{Q_3}\right) ,\quad \text {as}\quad t\rightarrow \infty , \end{aligned}$$

where $\sigma ^2_{Q_3} \equiv \lambda ^2 F^c(w)^2 \sigma ^2_{W}.$

The independence of $\widehat{Q}_1$,$\widehat{Q}_2$ and $\widehat{Q}_3$ yields

$$\begin{aligned} \widehat{Q}(t) \Rightarrow \widehat{Q}(\infty ) {\mathop {=}\limits ^\mathrm{d}}\mathcal {N}\left( 0, \sigma ^2_{Q}\right) , \quad \text {as} \quad t\rightarrow \infty , \qquad \text {where}\quad \sigma _{Q}^2 \equiv \sigma _{Q_1}^2+\sigma _{Q_2}^2+\sigma _{Q_3}^2. \end{aligned}$$

1.6 Proof of Corollary 4.2

Remaining service times are exponentially distributed due to lack of memory if the service-time distribution is exponential. Consequently, service completions at each server are a Poisson process with constant rate $\mu >0$, which implies by [35] that the sequence $\widehat{E}_n$ converges to a centered Gaussian process with covariance function $C_{{E}}(s,t) = \mu (s \wedge t)$ for $s,t \ge 0$. Then (4.23) becomes $\widetilde{C}(-x,-y)=\mu (x \vee y) - \mu | x-y |$ for $x \ge 0$, $y \ge 0$. Consequently, (A.14) becomes

$$\begin{aligned} \sigma ^2_{W_3}&= 2\frac{h_{F}(w)^2}{\lambda ^2 F^c(w)^2}\int _{0}^{\infty }\int _{0}^{x} e^{-h_{F}(w)(x+y)} \mu y\,\mathrm{d}y\mathrm{d}x \\&= 2\frac{h_{F}(w)^2}{\lambda ^2 F^c(w)^2}\int _{0}^{\infty } \mu e^{-h_{F}(w)x} \int _{0}^{x} y e^{-h_{F}(w) y} \,\mathrm{d}y\mathrm{d}x \\&= 2\frac{h_{F}(w)^2}{\lambda ^2 F^c(w)^2}\int _{0}^{\infty } \mu e^{-h_{F}(w)x} \left( -\frac{x}{h}e^{-h_{F}(w) x} + \frac{1}{h^2}\left( 1-e^{-h_{F}(w) x}\right) \right) \mathrm{d}x \\&= 2\frac{h_{F}(w)^2}{\lambda ^2 F^c(w)^2} \left( \frac{-\mu }{h_F(w)}\int _{0}^{\infty } x e^{-2h_{F}(w)x} \mathrm{d}x \right. \\&\quad \left. + \frac{\mu }{h_F(w)^2}\int _{0}^{\infty } e^{-h_{F}(w)x} \left( 1- e^{-h_{F}(w)x} \right) \mathrm{d}x \right) \\&= 2\frac{h_{F}(w)^2}{\lambda ^2 F^c(w)^2} \left( \frac{-\mu }{2h_F(w)^2}\frac{1}{2h_F(w)} + \frac{\mu }{h_F(w)^3} - \frac{\mu }{2h_F(w)^3} \right) = \frac{1}{2\lambda f(w)}. \end{aligned}$$

Summing up $\sigma ^2_{W_3}$ with $\sigma ^2_{W_i}(\infty )$ for $i=1,2$ yields (4.27).

The variance $\sigma ^2_{W_3}$ for the $M/M/n+M$ queue can be immediately obtained by letting $c_\lambda ^2=1$ and $h_F(w)=\theta$ in (4.27). Finally, we obtain $\sigma ^2_{Q}$ in (4.28) as follows:

$$\begin{aligned} \sigma ^2_{Q}(\infty )&= \lambda \int _0^w F^c(u)^2\,\mathrm{d}u + \lambda \int _0^w F(u) F^c(u)^2\,\mathrm{d}u + \lambda ^2 F^c(w)^2 \sigma _{W}^2 \\&= \lambda \int _0^w F^c(u)\,\mathrm{d}u + \lambda ^2 F^c(w)^2 \sigma _{W}^2 = \frac{\lambda }{\theta } (1-e^{-\theta w}) + \frac{\lambda }{\theta }\cdot \frac{1}{\rho } = \frac{\lambda }{\theta }, \end{aligned}$$

where the last equality holds because $w = F^{-1}\left( 1-1/\rho \right)$, so that $1-1/\rho = F(w) = 1-e^{-\theta w}$. $\square$

1.7 Proof of Lemma 5.1

1.7.1 Proof of the convergence in (5.2)

We consider the modified processes $\widehat{E}^{'}_{n,1}(t)$ given below. We first prove convergence for the sequence $\widehat{E}^{'}_{n,1}$ and then show that the difference between the modified sequence $\widehat{E}^{'}_{n,1}$ and the desired sequence $\widehat{E}_{n,1}$ is asymptotically negligible (see (A.18)), which proves the desired convergence in (5.2).

Now define, for $t \ge 0$,

$$\begin{aligned} \widehat{E}^{'}_{n,1}(t)&\equiv \frac{1}{\sqrt{n}} E^{'}_{n,1}(t) = \int _{0}^{t-w(t)} F^c(v(s))\,\mathrm{d}\widehat{N}_n(s) \nonumber \\&= F^c(v(t-w(t))) \widehat{N}_n(t-w(t)) - F^c(v(0)) \widehat{N}_n(0) \nonumber \\&\quad - \int _{0}^{t-w(t)} \widehat{N}_n(s-) \,\mathrm{d}F^c(v(s)) \nonumber \\&= F^c(w(t)) \widehat{N}_n(t-w(t)) - \widehat{N}_n(0) - \int _{0}^{t} \widehat{N}_n(s-w(s)) \,\mathrm{d}F^c(w(s)) .\nonumber \\ \end{aligned}$$

(A.16)

The second equality holds by integration by parts. The last equality follows from (4.3). Next we define a mapping $\psi :\mathcal {D}\rightarrow \mathcal {D}$ such that, for $z \in \mathcal {D}$,

$$\begin{aligned} \psi (z)(t) \equiv F^c(w(t)) z(t) - z(0) - \int _{0}^{t} z(s) \,\mathrm{d}F^c(w(s)), \quad 0 \le t \le T. \end{aligned}$$

We now prove that the mapping $\psi$ is continuous in $\mathcal {D}$. Let $\{x_n\}$ be a sequence in $\mathcal {D}$ such that $\Vert x_n - x \Vert _T \rightarrow 0$. Then

$$\begin{aligned}&| \psi (x_n)(t) - \psi (x)(t) | \\&\quad = \bigg | F^c(w(t)) x_n(t) - x_n(0) - \int _{0}^{t} x_n(s) \,\mathrm{d}F^c(w(s)) \\&\qquad - F^c(w(t)) x(t) + x(0) + \int _{0}^{t} x(s) \,\mathrm{d}F^c(w(s)) \bigg | \\&\quad \le F^c(w(t)) | x_n(t)-x(t) | + | x_n(0)-x(0) | \\&\qquad + \Vert x_n-x\Vert _T \bigg | \int _{0}^{t} \,\mathrm{d}F^c(w(s)) \bigg | \le 4\, \Vert x_n-x \Vert . \end{aligned}$$

Hence the mapping $\psi$ is continuous. In general, proving convergence with respect to the uniform topology does not necessarily imply $J_1$ convergence because there may be measurability issues (see, for example, [6, 36]). However, we will be interested in the case where the limit x is continuous, i.e., $x \in \mathcal {C}$. Therefore, we will not have any measurability issues and obtain the desired convergence in $\mathcal {D}$ with respect to Skorokhod’s $J_1$ metric.

Convergence of the modified process in (A.16) follows by the continuous mapping theorem with composition. In particular, let $Z_n(t) \equiv \widehat{N}_n(t -W_n(t))$. Then $Z_n:[0,T] \rightarrow {\mathbb {R}}$ and $Z_n \Rightarrow Z$, where $Z(t) \equiv \widehat{N}(t -w(t))$. Convergence of $\{Z_n\}$ follows from the continuous mapping theorem with composition. Then we have $n^{-1/2} \widetilde{E}_{n,1}(t ) = \psi (Z_n)(t) \Rightarrow \psi (Z)(t)$ in $\mathcal {D}$ with

$$\begin{aligned} \psi (Z) \equiv&\int _{0}^{t} F^c(w(s))\,\mathrm{d}\widehat{N}(s-w(s)) \equiv F^c(w(t)) \widehat{N}(t-w(t)) - \widehat{N}(0) \\&- \int _{0}^{t} \widehat{N}(s-w(s)) \,\mathrm{d}F^c(w(s))\\ =&F^c(w(t))c_{\lambda } {\mathcal {B}}_{\lambda }(\varLambda (t-w(t))) - c_{\lambda }{\mathcal {B}}_{\lambda }(0) \\&- \int _{0}^{t} c_{\lambda }{\mathcal {B}}_{\lambda }(\varLambda (s-w(s))) \,\mathrm{d}F^c(w(s)). \end{aligned}$$

Finally, to establish (5.2), we show that the difference between the processes $n^{-1/2} E_{n,1}(t)$ and $n^{-1/2} E^{'}_{n,1}(t)$ is asymptotically negligible. In particular,

$$\begin{aligned}&\left| \widehat{E}_{n,1}(t) - \widehat{E}^{'}_{n,1}(t)\right| \end{aligned}$$

(A.17)

$$\begin{aligned}&\quad = \frac{1}{\sqrt{n}} \left| \int _{0}^{t-W_n(t)} F^c(V_n(s-))\,\mathrm{d}\widehat{N}_n(s) - \int _{0}^{t-w(t)} F^c(v(s))\,\mathrm{d}\widehat{N}_n(s)\right| \nonumber \\&\quad \le \frac{1}{\sqrt{n}} \left| \int _{t-w(t)}^{t-W_n(t)} F^c(V_n(s-))\,\mathrm{d}\widehat{N}_n(s) \right| \nonumber \\&\qquad + \frac{1}{\sqrt{n}} \int _{0}^{t-w(t)} \left| F^c(V_n(s-)) -F^c(v(s))\right| \,\mathrm{d}\widehat{N}_n(s) \nonumber \\&\quad \le \frac{1}{\sqrt{n}} \left| \widehat{N}_n(t-W_n(t)) - \widehat{N}_n(t-w(t)) \right| +\frac{1}{\sqrt{n}} \left| \widehat{N}_n(t-w(t)) - \widehat{N}_n(0) \right| \Rightarrow 0. \end{aligned}$$

(A.18)

1.7.2 Proof of the convergence in (5.3)

To prove convergence in (5.3), we apply the martingale FCLT in [32] (see also [9, 14] for applications of the martingale FCLT). First we define a sequence of discrete-time processes (see (A.19)) and argue that it is a sequence of martingales adapted to a specific filtration ${\mathscr {H}}_k^n$ as defined below. Next, we define continuous-time martingales using the discrete-time martingales in (A.19). Then we invoke Theorem 7.1.4. on p.339 in [10] to establish convergence and characterize the limit.

Consider the discrete-time processes

$$\begin{aligned} \widehat{H}^n_k \equiv \frac{1}{\sqrt{n}} \sum _{i=1}^{k} \left( {\mathbf {1}}(\gamma ^n_{i} > w^n_{i}) - F^c(w^n_{i}) \right) \quad \text{ for }\quad k=1,2,\dots . \end{aligned}$$

(A.19)

Also, consider the filtration ${\mathscr {H}}_k^n \equiv \sigma \{ \tau ^n_{i+1}, \nu _i^n, \gamma _{i}^n : 1 \le i \le k \}$. Then $\mathbb {E}[ | \widehat{H}^n_k |] \le k/ \sqrt{n}$ and

$$\begin{aligned} \mathbb {E}[H_k^n - H_{k-1}^n | {\mathscr {H}}_{k-1}^n] = \frac{1}{\sqrt{n}} \left( \mathbb {E}[ {\mathbf {1}}(\gamma ^n_{k} > w^n_{k}) | {\mathscr {H}}_{k-1}^n ] - F^c(w^n_{k}) \right) = 0, \end{aligned}$$

which implies that the process $\{(\widehat{H}^n_k,{\mathscr {H}}_k^n): k\ge 1 \}$ is a discrete-time martingale for each $n\ge 1$.

Our next step is to replace k with $\lfloor nt \rfloor$ for $t \ge 0$ to obtain a continuous-time martingale. By a direct application of Lemma 4.2 of [9], we deduce that the continuous-time process $(\widehat{H}^n(t),{\mathscr {H}}^n(t): t \ge 0) \equiv (\widehat{H}_{\lfloor nt \rfloor }^n, {\mathscr {H}}_{\lfloor nt \rfloor }^n : t \ge 0)$ is a martingale with quadratic variation

$$\begin{aligned} \langle \widehat{H}^n\rangle (t)= \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } \left( {\mathbf {1}}(\gamma ^n_{i} > w^n_{i}) - F^c(w^n_{i}) \right) ^2. \end{aligned}$$

(A.20)

We next show that the sequence of martingales $(\widehat{H}^n(t),{\mathscr {H}}^n(t): t \ge 0)$ satisfies the conditions of Theorem 7.1.4. of [10]. In particular, it is required that (i) jumps of the processes $\widehat{H}^n(y)$ are asymptotically negligible and (ii) the quadratic quadratic variation of the processes converges in probability to a limit characterized in Theorem 7.1.1. of [10].

(i) Negligibility of jumps. We now show that condition (a) of Theorem 7.1.4. holds. Let $\widehat{H}^n(t-)\equiv \lim _{s \uparrow t} \widehat{H}^n(s)$. Then, for each $T>0$, we have $\sup _{0 \le t \le T} | \widehat{H}^n(t) - \widehat{H}^n(t-) |\le 1/\sqrt{n}$ and hence

$$\begin{aligned} \lim _{n \rightarrow \infty } \mathbb {E}\left[ \sup _{0 \le t \le T} | \widehat{H}^n(t) - \widehat{H}^n(t-) | \right] = 0, \end{aligned}$$

which is the desired condition.

(ii) Convergence of quadratic variations. We now prove that the quadratic variation processes given in (A.20) converges in the $L^2$ sense as $n \rightarrow \infty$. In particular,

$$\begin{aligned}&\mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } \left( {\mathbf {1}}(\gamma ^n_{i}> w^n_{i}) - F^c(w^n_{i}) \right) ^2 - \int _0^{\varLambda ^{-1}(t)} F^c(v(u)) F^c(v(u)) \,\mathrm{d}\varLambda (u) \right) ^2 \right] \nonumber \\&\quad \le 2 \, \mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } \left[ \left( {\mathbf {1}}(\gamma ^n_{i}> w^n_{i}) - F^c(w^n_{i}) \right) ^2 - F^c(w^n_{i})F(w^n_{i}) \right] \right) ^2 \right] \nonumber \\&\qquad + 4 \, \mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } \left[ F^c(w^n_{i}) F(w^n_{i})- F^c(v(\tau ^n_{i}-)) F(v(\tau ^n_{i}-))\right] \right) ^2 \right] \nonumber \\&\qquad + 4 \, \mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } F^c(v(\tau ^n_{i}-)) F(v(\tau ^n_{i}-))\right. \right. \nonumber \\&\qquad \left. \left. - \int _0^{\varLambda ^{-1}(t)} F^c(v(u-)) F(v(u-))\,\mathrm{d}\varLambda (u) \right) ^2 \right] \nonumber \\&\quad \le \frac{2}{n^2}\, \sum _{i=1}^{\lfloor nt \rfloor } \mathbb {E}\left[ \left( {\mathbf {1}}(\gamma ^n_{i}> w^n_{i}) - F^c(w^n_{i}) \right) ^2 \Big ( F(w^n_{i}) - F^c(w^n_{i}) \Big )^2 \right] \nonumber \\&\qquad + \frac{2}{n^2} \mathbb {E}\sum _{i \ne j} \left[ \left( {\mathbf {1}}(\gamma ^n_{i}> w^n_{i}) - F^c(w^n_{i}) \right) \left( {\mathbf {1}}(\gamma ^n_{j} > w^n_{j}) - F^c(w^n_{j}) \right) \right. \nonumber \\&\qquad \left. \times \Big (F(w^n_{i}) - F^c(w^n_{i})\Big ) \Big ( F(w^n_{j}) - F^c(w^n_{j}) \Big ) \right] \nonumber \\&\qquad + 4 \, \mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } \left[ F^c(w^n_{i}) F(w^n_{i})- F^c(v(\tau ^n_{i}-)) F(v(\tau ^n_{i}-))\right] \right) ^2 \right] \end{aligned}$$

(A.21)

$$\begin{aligned}&\qquad + 4 \, \mathbb {E}\left[ \left( \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } F^c(v(\tau ^n_{i}-)) F(v(\tau ^n_{i}-))\right. \right. \nonumber \\&\qquad \left. \left. - \int _0^{\varLambda ^{-1}(t)} F^c(v(u)-) F(v(u)-)\,\mathrm{d}\varLambda (u) \right) ^2 \right] . \end{aligned}$$

(A.22)

The first sum vanishes as $n \rightarrow \infty$ because the summands are bounded by 1 and, therefore, the first term is bounded by $2\lfloor nt \rfloor /n^2 \rightarrow 0$ as $n \rightarrow \infty$. The summands of the second term are independent. Therefore, the second term is equal to 0.

To prove convergence of (A.21), we first rewrite the summands of (A.21) as

$$\begin{aligned}&F^c(w^n_{i})F(w^n_{i})- F^c(v(\tau ^n_{i}-))F(v(\tau ^n_{i}-)) = F^c(w^n_{i}) - F^c(v(\tau ^n_{i}-)) \nonumber \\&\quad - \left( F^c(w^n_{i})^2 -F^c(v(\tau ^n_{i}-))^2 \right) , \end{aligned}$$

(A.23)

Next we make use of the FWLLN for PWT $V_n(t)$, i.e., $V_n \Rightarrow v$ in $\mathcal {D}$, and continuity of the function F to show that (A.21) converges to 0. In particular, for all $i \ge 1$,

$$\begin{aligned} F^c(w^n_{i}) = F^c(V_n(\tau ^n_{i}-)) = F^c(v(\tau ^n_{i}-)+o(1)). \end{aligned}$$

Combining with (A.23), this implies that the summands in (A.21) can be bounded above by

$$\begin{aligned}&| F^c(v(\tau ^n_{i}-)+o(1)) - F^c(v(\tau ^n_{i}-)) | + | F^c(v(\tau ^n_{i}-)+o(1))^2 \\&\quad - F^c(v(\tau ^n_{i}-))^2 | \le |o(1)|, \end{aligned}$$

where the inequality holds by continuity of cdf F. This implies that the squared sum inside the expectation in (A.21) is bounded above by $( |o(1)| \lfloor nt \rfloor /n )^2 \le t^2 |o(1)| = o(1)$ for all $t \ge 0$. Convergence of (A.21) to 0 then follows from the dominated convergence theorem.

The summation in (A.22) can be alternatively represented as

$$\begin{aligned} \frac{1}{n} \sum _{i=1}^{\lfloor nt \rfloor } F^c(v(\tau ^n_{i}-)) F(v(\tau ^n_{i}-))&= \int _0^{\varLambda ^{-1}_n(t)} F^c(v(u-)) F(v(u-))\,d\bar{N}_n(u) \nonumber \\&\Rightarrow \int _0^{\varLambda ^{-1}(t)} F^c(v(u-)) F(v(u-))\,\mathrm{d}\varLambda (u), \end{aligned}$$

(A.24)

where the convergence (A.24) follows from the continuous mapping theorem. Having established the convergence in (A.24), convergence in mean square is obtained by first applying the continuous mapping theorem with the function $f(x) = x^2$ and then applying the dominated convergence theorem by using the fact that both the summation and the limit integral in (A.24) are bounded by t. Hence (A.22) converges to 0. That completes the proof of convergence of the quadratic variation (A.20).

Having proved conditions (i) and (ii) are indeed satisfied, by Theorem 7.1.4 of [10], we deduce that $\widehat{H}^n \Rightarrow \widehat{H}$ in $\mathcal {D}$, where $\widehat{H}$ is a Gaussian process with independent increments and continuous sample paths. Moreover, as implied by the proof of Theorem 7.1.1. of [10], the limit $\widehat{H}$ is indeed a time-changed Brownian motion, where the time change is the limit of the quadratic variation, i.e.,

$$\begin{aligned} \widehat{H}(t) = {\mathcal {B}}_a(\langle \widehat{H} \rangle (t)) = {\mathcal {B}}_a \left( \int _0^{\varLambda ^{-1}(t)} F^c(v(u-))F(v(u-))\, \mathrm{d}\varLambda (u) \right) , \end{aligned}$$

where ${\mathcal {B}}_a$ is the standard Brownian motion.

Finally, to complete the proof, we note that $\widehat{E}_{n,2}(t) = \widehat{H}^n(\bar{N}_n(t-W_n(t)))$. Then, by the convergence-together theorem, we have $\widehat{H}^n(\bar{N}_n) \Rightarrow \widehat{H}(\varLambda )$ in $\mathcal {D}$. Consequently, as $n\rightarrow \infty$,

$$\begin{aligned}&\frac{1}{\sqrt{n}} \sum _{i=1}^{N_n(t-W_n(t))} \left( {\mathbf {1}}(\gamma ^n_{i}>w^n_{i}) - F^c(w^n_{i}) \right) \nonumber \\&\quad \Rightarrow {\mathcal {B}}_a \left( \int _{0}^{t-w(t)} F^c(v(u-))F(v(u-))\,\mathrm{d}\varLambda (u) \right) . \end{aligned}$$

(A.25)

We next verify the other two expressions in (5.3). The last expression is obtained by a change of variable with $u = s-w(s)$. (Note that, according to (4.3), we have $v(s-w(s)) = w(s)$.) The Kiefer integral expression holds because it is a Gaussian process with zero mean and the same covariance function as the Brownian expression. Specifically, for $t,t'>0$, the first Brownian expression has the covariance

$$\begin{aligned} \int _0^{(t-w(t))\wedge (t'-w(t'))} F^c(v(u))F(v(u)) \,\mathrm{d}\varLambda (u). \end{aligned}$$

(A.26)

On the other hand, the Kiefer integral in (5.3) has the covariance

$$\begin{aligned}&\mathbb {E}\left[ \int _0^{t-w(t)} \int _{0}^{1} {\mathbf {1}}(y>F(v(s))) \mathrm{d}\widehat{U}(s,y) \right. \\&\qquad \left. \times \int _0^{t'-w(t')} \int _{0}^{1} {\mathbf {1}}(y>F(v(s))) \mathrm{d}\widehat{U}(s,y) \right] \\&\quad =\mathbb {E}\left[ \int _0^{t-w(t)} \int _{0}^{\infty } {\mathbf {1}}(x>v(s)) \mathrm{d}\widehat{U}(s,F(x)) \right. \\&\qquad \left. \times \int _0^{t'-w(t')} \int _{0}^{\infty } {\mathbf {1}}(x>v(s)) \mathrm{d}\widehat{U}(s,F(x)) \right] \\&\quad = \int ^{(t-w(t))\wedge (t'-w(t'))}_0 \int _{0}^{\infty } {\mathbf {1}}(x>v(s)) \mathrm{d}F(x)\mathrm{d}\varLambda (s)\\&\qquad + \int ^{(t-w(t))\wedge (t'-w(t'))}_0 F^c(v(s))F^c(v(s))\mathrm{d}\varLambda (s) \\&\qquad - 2 \int ^{(t-w(t))\wedge (t'-w(t'))}_0 \int _{0}^{\infty } {\mathbf {1}}(x>v(s)) F^c(v(s)) \mathrm{d}F(x)\mathrm{d}\varLambda (s) \\&\quad = \int ^{(t-w(t))\wedge (t'-w(t'))}_0 F^c(v(s))\mathrm{d}\varLambda (s) + \int ^{(t-w(t))\wedge (t'-w(t'))}_0 \left( F^c(v(s))\right) ^2 \mathrm{d}\varLambda (s) \\&\qquad - 2 \int ^{(t-w(t))\wedge (t'-w(t'))}_0 F^c(v(s))F^c(v(s))\mathrm{d}\varLambda (s) \\&\quad = \int ^{(t-w(t))\wedge (t'-w(t'))}_0 F^c(v(s))\mathrm{d}\varLambda (s) - \int ^{(t-w(t))\wedge (t'-w(t'))}_0 \left( F^c(v(s))\right) ^2 \mathrm{d}\varLambda (s) \\&\quad = \int ^{(t-w(t))\wedge (t'-w(t'))}_0 F^c(v(s))\left( 1-F^c(v(s))\right) \mathrm{d}\varLambda (s), \end{aligned}$$

which coincides with (A.26).

Refined staffing levels

In this section, we consider a refined staffing function given by

$$\begin{aligned} s_n \equiv \lceil ns_1 + \sqrt{n} s_2 \rceil , \quad \text {where} \quad s_1, \, s_2 >0. \end{aligned}$$

(B.1)

The general form of (B.1) enables us to recover two of the staffing functions considered in [30] that, respectively, lead to the ED and ED+QED operating regimes. More specifically, the two staffing functions in [30] are given by

$$\begin{aligned} n_{\text {ED}}&= \lceil (1- \gamma ) R_n \rceil , \end{aligned}$$

(B.2)

$$\begin{aligned} n_{\text {ED+QED}}&= \lceil (1- \gamma ) R_n + \delta \sqrt{R_n} \rceil , \end{aligned}$$

(B.3)

where $0<\gamma <1$ and $R_n$ is the offered load, defined as $R_n=n \lambda /\mu$. Letting $s_1 = (1-\gamma ) \lambda / \mu$ and $s_2 = 0$ yields (B.2), whereas letting $s_1 = (1-\gamma ) \lambda / \mu$ and $s_2 = \delta \sqrt{\lambda / \mu }$ yields (B.3). See also [21] and Sect. 10 in [26] for time-varying versions of the refined staffing (B.1).

We next briefly discuss the changes resulting from considering the staffing function $s_n$ instead of n. In the previous sections, the staffing function happens to coincide with our scaling factor n, i.e., $s_n=n$. In this section, we let n and $\sqrt{n}$ be the scaling factors for FWLLN and FCLT, respectively, and let the staffing function have a more general form $s_n = \lceil ns_1 + \sqrt{n} s_2 \rceil$, where $s_1,s_2>0$. To indicate the processes associated with the new staffing function, we use a superscript r, whereas to indicate the processes associated the case where $s_n=n$, we use notation without a superscript. Because the arrival process is independent of the staffing level, it holds that $\widehat{N}_n^r(t) = \widehat{N}_n(t)$ for all $n\ge 1$ and $t\ge 0$, and hence, $\widehat{N}^r(t) = \widehat{N}(t)$. The FWLLN limit and the FCLT limit for the service-completion process, on the other hand, becomes $D^r(t) = s_1 D(t)$, and $\widehat{D}^r(t) = \sqrt{s_1} \widehat{D}(t)+ s_2 D(t)$, respectively, where $\widehat{D}(t)$ is a centered Gaussian process with covariance function $C_E$ in Theorem 4.2, and $D(t) = \mu t$ as in (4.4). Hence, $\widehat{D}^r(t)$ is a Gaussian process with covariance function $C^r(\cdot ,\cdot ) = s_1 C_E(\cdot ,\cdot )$ and mean $s_2 \mu t$. Consequently, the enter-service process satisfies $\widehat{E}^r(t) = \sqrt{s_1} \widehat{D}(t)+ s_2 \mu t$.

The following theorem is an analog of Theorems 4.1 and 4.2 for the $G/GI/n+GI$ model having the refined staffing function $s_n$ in (B.1).

Theorem B.1

(FWLLN and FCLT with refined staffing) Consider the $G/GI/n+GI$ with staffing level $s_n$ given by (B.1) and $\rho ^r = \lambda /\mu s_1>1$.

$\mathrm{(i)}$ Under the conditions of Theorem 4.1, an analog of joint convergence in (4.1) holds as $n \rightarrow \infty$, where $\varLambda ^r(t)=\lambda t$,

$$\begin{aligned} D^r(t)&= E^r(t) = s_1 \mu t, \quad w^r(t) = \int _0^t \left( 1- \frac{s_1 \mu }{\lambda F^c(w^r(u))} \right) \,\mathrm{d}u, \nonumber \\ v^r(t)&= w^r(t+v^r(t)). \end{aligned}$$

(B.4)

The limits $Q^r(t)$, $X^r(t)$ and $A^r(t)$ have the same mathematical form as their counterparts in Theorem 4.1 with modified components.

$\mathrm{(ii)}$ Under the conditions of Theorem 4.1, an analog of joint convergence in (4.5) holds as $n \rightarrow \infty$, where

$$\begin{aligned} \widehat{W}^r(t) =&\int _0^t \frac{F^c(w^r(u)) H^r(t,u)}{q(t,w^r(t))}c_{\lambda }\, \mathrm{d}{\mathcal {B}}_{\lambda }(\varLambda (u-w^r(u))) \nonumber \\&+ \int _0^t \frac{\sqrt{\lambda F^c(v^r(u))F(v^r(u)) } H^r(t,u)}{q(t,w^r(t))}\, \mathrm{d}{\mathcal {B}}_a(u) \nonumber \\&- \sqrt{s_1} \int _0^t \frac{H^r(t,u)}{q(t,w^r(t))}\, d\hat{E}(u) - s_2 \mu \int _0^t \frac{H^r(t,u)}{q(t,w^r(t))}\,\mathrm{d}u, \end{aligned}$$

(B.5)

$w^r(t)$ and $v^r(t)$ are as in (B.4), $H^r(\cdot ,\cdot )$ and $q(\cdot ,w^r(\cdot ))$ are as in (4.17) with w(t) replaced by $w^r(t)$. The virtual waiting time $\widehat{V}^r(t)$ and the queue-length process $\widehat{Q}^r(t)$ have the same mathematical forms as in (4.9) and (4.10), with w(t), v(t) and $\widehat{W}(t)$ replaced by their counterparts $w^r(t)$, $v^r(t)$ and $\widehat{W}^r(t)$. The FCLT limit for the abandonment process is $\widehat{A}^r(t) = \widehat{N}(t)- \widehat{Q}^r(t)- \widehat{E}^r(t)$.

Proof of Theorem B.1

The proof closely follows the arguments in the proofs of Theorem 4.1, Theorem 4.2, Corollary 4.1 and Theorem 4.4. Therefore, we mostly refer to proofs of those results in the proofs below and argue in what way the new staffing function $s_n = \lceil ns_1+\sqrt{n}s_2 \rceil$ changes the arguments. We skip lengthy details. Throughout this subsection, the processes with a superscript r correspond to those associated with staffing level $s_n$, whereas the processes without a superscript r correspond to those associated with staffing level n.

The LLN- and CLT-scaled departure process

$$\begin{aligned} \bar{D}^r_n(t)\equiv & {} \frac{\sum _{j=1}^{s_n} D_j(t)}{n} = \frac{s_n}{n} \cdot \frac{\sum _{j=1}^{s_n} D_j(t)}{s_n} \Rightarrow D^r(t) \equiv s_1 D(t) = s_1 \mu t,\end{aligned}$$

(B.6)

$$\begin{aligned} \widehat{D}^r_n(t)\equiv & {} \frac{\sum _{j=1}^{s_n} D_j(t)-n D^r(t)}{\sqrt{n}} = \frac{\sqrt{s_n}}{\sqrt{n}} \cdot \frac{\sum _{j=1}^{s_n} D_j(t) - s_n \mu t}{\sqrt{s_n}} + \frac{s_n \mu t- n D^r(t)}{\sqrt{n}} \nonumber \\= & {} \sqrt{s_1+ \frac{s_2}{\sqrt{n}}} \cdot \frac{\sum _{j=1}^{s_n} D_j(t) - s_n \mu t}{\sqrt{s_n}} + s_2 \mu t + O(1/\sqrt{n}) \Rightarrow \widehat{D}^r(t) \nonumber \\\equiv & {} \sqrt{s_1} \widehat{D}(t) - s_2\mu t, \end{aligned}$$

(B.7)

where $O(1/\sqrt{n})$ in the second equality accounts for the error caused by dropping $\lceil \cdot \rceil$ in $s_n$, and $\widehat{D}(t)$ is the Gaussian process in Theorem 4.2. Hence, we deduce from (B.7) that $\hat{D}^r(t)$ is a Gaussian process with negative drift $- s_2 \mu t$ and covariance function $C^r(\cdot ,\cdot ) = s_1 C_E(\cdot ,\cdot )$, with $C_E$ being the covariance function in Theorem 4.2.

Having obtained the modified fluid limits in (B.6) and established the joint convergence, we deduce that the proof in Sect. A.2 continues to hold with minor modifications. But the limit in (A.7) changes because the fluid limit of the departure process is now given by $D^r(t)=s_1 \mu t$. Consequently, the ODE in Theorem 4.1 has $s_1 \mu t$ in the numerator instead of $\mu t$.

Similarly, given the joint convergence $(\widehat{N}^r_n,\widehat{D}^r_n, \widehat{E}^r_n) \Rightarrow (\widehat{N}^r,\widehat{D}^r,\widehat{E}^r)$, we can prove the FCLT with a slightly modified proof. The arguments in Sect. 5 continue to hold for modified fluid limits and cause only minor changes in the final expressions. In particular, (5.4)–(5.7) have the same mathematical form with fluid limits and prelimit stochastic process replaced with their counterparts with a superscript r. Hence the steps of the proof in Sect. 5.1.2 can be replicated with counterpart processes. The only step that requires careful treatment is that the limit of the enter-service process is now $\widehat{E}^r(t) \equiv \sqrt{s_1} \widehat{E}(t) - s_2\mu t$. Since the additional term $- s_2\mu t$ is deterministic and $\sqrt{s_1} \widehat{E}(t)$ is a centered Gaussian process, we can use similar arguments to proof of Corollary 4.1 to deduce that (B.5) is indeed the desired solution. $\square$

Note that (B.5) is different than (4.16) in that the third term on the right-hand side is scaled by $\sqrt{s_1}$ and that there is an additional deterministic term. This implies that both the variance and mean of HWT change and so do those of the PWT and queue length. The corresponding steady-state formulas are given in the following corollary.

Corollary B.1

(Steady state of limits with refined staffing) Under the assumptions of Theorem 4.4, the steady-state random variables $\hat{W}^r(\infty )$, $\hat{V}^r(\infty )$ and $\hat{Q}^r(\infty )$ have Gaussian distributions with means and variances given below:

$$\begin{aligned} \mu _{W^r} \equiv \mathbb {E}\left[ \widehat{W}^r(\infty ) \right]&= \mathbb {E}\left[ \widehat{V}^r(\infty ) \right] = - \frac{s_2 \mu }{\lambda f(w^r) },\\ \mathbb {E}\left[ \widehat{Q}^r(\infty ) \right]&= \lambda F^c(w^r)\mathbb {E}\left[ \widehat{W}^r(\infty ) \right] = - \frac{s_2 \mu }{h_F(w^r)},\\ \mathrm{Var}(\widehat{W}^r(\infty )) =\mathrm{Var}(\widehat{V}^r(\infty ))&\equiv \sigma ^2_{W^r} \equiv \frac{c_{\lambda }^2}{2h_{F}(w^r)\lambda } + \frac{F(w^r)}{2 \lambda f(w^r)} + s_1 \sigma ^2_{W_3^r},\\ \text {where}&\quad \sigma ^2_{W_3^r} \equiv 2\frac{h_{F}(w^r)^2}{\lambda ^2 F^c(w^r)^2}\\&\qquad \times \int _{0}^{\infty }\int _{0}^{x} e^{-h_{F}(w^r)(x+y)}\widetilde{C}(-x, -y)\,\mathrm{d}y\mathrm{d}x, \end{aligned}$$

the covariance function $\widetilde{C}(-x, -y)$ is as in (4.23), and $w^r = F^{-1}(1-1/\rho ^r)$. The variance of the steady-state queue length $\widehat{Q}(\infty )$ has the same mathematical form with w and $\widehat{W}$ replaced by $w^r$ and $\widehat{W}^r$.

Remark B.1

(Optimal staffing problems) Heavy-traffic FWLLN and FCLT limits have proven useful in solving optimal staffing problems with respect to service-level constraints in large scale service systems [4, 30]. A general framework for this type of approach has two steps: First, a corresponding optimal staffing problem is formulated and solved using analytic FWLLN or FCLT limits (which are often more convenient than their corresponding stochastic versions); Next, an asymptotic optimality result is established by showing that the FWLLN- or FCLT-based optimal staffing problem is asymptotically equivalent to its desired stochastic version as the scale $n\rightarrow \infty$. We advocate that our new FCLT limit with refined staffing functions provides a basis for solving optimal staffing problems in $G/GI/n+GI$ queueing systems (note that two control factors $s_1$ and $s_2$ for the staffing function are preserved in the limit). For example, in the FCLT-based optimal staffing problem, we may choose the optimal $s_1^*$ and $s_2^*$ in order to minimize certain performance functions, for example, the mean waiting time, queue length, or abandonment probability; see the formulation in [4]. We leave this to future research.

Proof of Corollary B.1

We first derive the mean of $\widehat{W}(\infty )$ from (B.5). Since the first three terms in (B.5) have zero means, $\mathbb {E}[\widehat{W}(\infty )]$ is the limit of the last term in (B.5) as $t \rightarrow \infty$:

$$\begin{aligned} \mathbb {E}[\widehat{W}^r(\infty )]&= \lim _{t \rightarrow \infty } -s_2 \mu \int _0^t \frac{H^r(t,u)}{q(t,w^r(t))}\,\mathrm{d}u = \lim _{t \rightarrow \infty } -s_2 \mu \int _0^t \frac{e^{-h_F(w^r)(t-u)}}{\lambda F^c(w^r)}\,\mathrm{d}u\\&= \lim _{t \rightarrow \infty } \frac{-s_2 \mu }{\lambda f(w^r)} \left( 1-e^{-h_F(w^r)t} \right) \!. \end{aligned}$$

Having established $\mathbb {E}[\widehat{W}(\infty )]$, it is easy to establish $\mathbb {E}[\widehat{V}(\infty )]$ and $\mathbb {E}[\widehat{Q}(\infty )]$ by letting $t\rightarrow \infty$ in

$$\begin{aligned} \widehat{V}^r(t) = \frac{\widehat{W}^r(t)}{1-\dot{w}^r(t+v^r(t))} \quad \text{ and }\quad \widehat{Q}^r(t) = \lambda F^c(w^r(t)) \widehat{W}^r(t). \end{aligned}$$

Computation of variance is standard and as given in Sect. A.4. $\square$