# Designing Low-Complexity Heavy-Traffic Delay-Optimal Load Balancing Schemes:
Theory to Algorithms

Xingyu Zhou
Department of ECE
The Ohio State University
<EMAIL_ADDRESS>Fei Wu
Department of CSE
The Ohio State University
<EMAIL_ADDRESS>Jian Tan
Department of ECE
The Ohio State University
<EMAIL_ADDRESS>
Yin Sun
Department of ECE
Auburn University
<EMAIL_ADDRESS>
Ness Shroff
Department of ECE and CSE
The Ohio State University
<EMAIL_ADDRESS>
###### Abstract

We establish a unified analytical framework for load balancing systems, which
allows us to construct a general class $\Pi$ of policies that are both
throughput optimal and heavy-traffic delay optimal. This general class $\Pi$
includes as special cases popular policies such as join-shortest-queue and
power-of-$d$, but not the join-idle-queue (JIQ) policy. In fact we show that
JIQ, which is not in $\Pi$, is actually not heavy-traffic delay optimal. Owing
to the significant flexibility offered by class $\Pi$, we are able to design a
new policy called join-below-threshold (JBT-d), which maintains the simplicity
of pull-based policies such as JIQ, but updates its threshold dynamically. We
prove that JBT-$d$ belongs to the class $\Pi$ when the threshold is picked
appropriately and thus it is heavy-traffic delay optimal. Extensive
simulations show that the new policy not only has a low complexity in message
rates, but also achieves excellent delay performance, comparable to the
optimal join-shortest-queue in various system settings.

##  1 Introduction

Load balancing, which is responsible for dispatching jobs on parallel servers,
is a key component in computer networks and distributed computing systems. For
a large number of practical applications, such as, Web service [6],
distributed caching systems (e.g., Memcached [13]), large data stores (e.g.,
HBase [5]), embarrassingly parallel computations [1] and grid computing [4],
the system performance critically depends on the load balancing algorithm it
employs.

In a load balancing system, there are two directions of message flows: push
messages (from the dispatcher to the servers) and pull messages (from the
servers to the dispatcher). In a push-based policy, the dispatcher actively
sends query messages to the servers and waits for their responses; In a pull-
based policy, the dispatcher passively listens to the report from the servers.
The job dispatching decision is conducted at the dispatcher based on the pull-
messages sent from the servers. Push-based policies (e.g., the join-shortest-
queue (JSQ) policy [19], [2] and the power-of-$d$ policy [10], [17]) have been
shown to be delay optimal in the heavy-traffic regime [2], [9]. Recently, the
pulled-based policies such as join-idle-queue (JIQ) [8] and the equivalent one
in [15], have been proposed. Compared with the push-based policies, these
pull-based policies not only achieve good delay performance, but also have
some nice features, such as, lower message overhead, lower computational
complexity, and zero dispatching delay. However, as shown in the simulations
of [8], the delay performance of existing pull-based polices will degrade
substantially as the load gets higher. In fact, as shown in Theorem 2 of this
paper, JIQ is not heavy-traffic delay optimal even for homogeneous servers.
Therefore, one key question is how to design load balancing policies that are
heavy-traffic delay optimal and meanwhile possess all the nice features of
pull-based policies such as zero dispatching delay, low message overhead and
low computational complexity.

In this paper, we take a systematic approach to address this question. To that
end, the main contributions of this paper are summarized as follows:

  * •

We derive inner-product based sufficient conditions for proving that a load-
balancing policy is throughput optimal and heavy-traffic delay optimal. Using
these sufficient conditions, we obtain a general class $\Pi$ of load balancing
policies that are throughput optimal and heavy-traffic delay optimal. This
class of load balancing policies contains the famous (push-based) JSQ and the
power-of-d policies as special cases, but not the (pull-based) JIQ policy.

  * •

On the other hand, we show that JIQ, which is not in $\Pi$, is not heavy-
traffic delay optimal even for homogeneous servers. While it has been
empirically shown in the past that the delay using JIQ is quite bad at high
loads, the question of whether it was heavy-traffic delay optimal in
homogeneous servers has been previously left unsolved. Furthermore, our novel
Lyapunov-drift approach offers a new avenue to show a policy is not heavy-
traffic delay optimal.

  * •

Owing to the significant flexibility offered by class $\Pi$, we are able to
design a new policy called Join-Below-Threshold (JBT-$d$). To the best of our
knowledge, this is the first load balancing policy that guarantees heavy-
traffic delay optimality while enjoying nice features of pull-based policy,
e.g., zero dispatching delay, low message overhead and low computational
complexity. Through extensive simulations, we demonstrate that JBT-$d$ has
excellent delay performance for different system sizes and various arrival and
service processes over a large range of traffic loads.

The rest of the paper is organized as follows. Section 1.1 reviews the related
work on load balancing schemes. Section 1.2 introduces the necessary notation
in the paper. Section 2 describes the system model and the related
definitions. Section 3 presents the main results of the paper. In particular,
a class $\Pi$ of flexible load balancing policies are introduced, containing
as special cases the popular existing ones and motivating new ones. Sufficient
conditions are derived to guarantee throughput and heavy-traffic delay
optimality. Section 4 contains the simulation results on comparing different
policies, demonstrating the performance and simplicity of our new policy.
Section 5 contains the proof of the main results.

###  1.1 Related work: push versus pull

This section reviews state-of-the-art load balancing policies with a focus on
the system performance in heavy traffic. We group these policies mainly into
two categories: push-based and pull-based as shown in Fig. 1.

Push-based policy: Under a push-based policy, the dispatcher tries to “push”
jobs to servers. More specifically, upon each job arrival, the dispatcher
sends probing messages to the servers, which feed back the required
information for dispatching decisions, e.g., queue lengths. After receiving
the feedback, the dispatcher sends the incoming jobs to servers based on a
dispatching distribution. A classical example in this category is the JSQ
policy, under which the dispatcher queries the queue length information of
each server upon new job arrivals, and sends the incoming jobs to the server
with the shortest queue, with ties broken randomly. It has been shown [19]
that for homogeneous servers this policy is delay optimal in a stochastic
ordering sense under the assumption of renewal arrival and non-decreasing
failure rate service. In the heavy-traffic regime, it has been proved that it
is heavy-traffic delay optimal for both heterogeneous and homogeneous servers
[2]. Nevertheless, the performance of this policy comes at the cost of
substantial overhead as it has to sample the queue lengths of all the servers,
which is undesirable in large-scale systems. To overcome this problem, an
alternative load balancing policy called power-of-$d$ has been introduced
[10], [17]; see also related works [20], [16]. Under this policy, the
dispatcher routes all the incoming jobs to the server that has the shortest
queue length, with ties broken randomly, out of the $d$ servers sampled
uniformly at random. This policy has also been shown to be heavy-traffic
optimal for homogeneous servers [9]. However, for heterogeneous servers, the
power-of-$d$ policy is neither throughput optimal, nor delay optimal in heavy
traffic.

Pull-based policy: Under a pull-based policy, the servers spontaneously send
messages to “pull” jobs from the dispatcher according to a fixed policy. One
illustrative example is the JIQ policy [8] and the equivalent one in [15].
Under the JIQ policy, each server sends a pull message to the dispatcher
whenever it becomes idle. Upon job arrivals, the dispatcher checks the
available pull messages in memory. If such messages exist, it removes one
uniformly at random, and sends the jobs to the corresponding server.
Otherwise, the new jobs will be dispatched uniformly at random to one of the
servers in the system. This policy has several favorable properties. The most
important property is that the required number of messages in steady-state is
at most one for each job arrival, which is smaller than the $2d$ of the power-
of-$d$-choices ($d$ for query and $d$ for response per job). However, as
already shown in [8], when the load becomes heavy, the performance of JIQ
keeps empirically degrades substantially, and in fact, in Theorem 2 we show
that it is not heavy traffic delay optimal even for homogeneous servers.

###  1.2 Notations

We use boldface letters to denote vectors in $\mathbb{R}^{N}$ and ordinary
letters for scalers. Denote by $\overline{{\mathbf{Q}}}$ the random vector
whose probability distribution is the same as the steady-state distribution of
$\\{{\mathbf{Q}}(t),t\geq 0\\}$. The dot product in $\mathbb{R}^{N}$ is
denoted by $\langle\mathbf{x},{\mathbf{y}}\rangle:=\sum_{i=1}^{N}x_{i}y_{i}$.
For any $\mathbf{x}\in\mathbb{R}^{N}$, the $l_{1}$ norm is denoted by
$\left\lVert\mathbf{x}\right\rVert_{1}:=\sum_{n=1}^{N}|x_{n}|$ and $l_{2}$
norm is denoted by
$\left\lVert\mathbf{x}\right\rVert:=\sqrt{\langle\mathbf{x},\mathbf{x}\rangle}$.
The parallel and perpendicular component of the queue length vector
${\mathbf{Q}}$ with respect to a vector $\mathbf{c}$ with unit norm is denoted
by ${\mathbf{Q}}_{\parallel}:=\langle\mathbf{c},{\mathbf{Q}}\rangle\mathbf{c}$
and ${\mathbf{Q}}_{\perp}:={\mathbf{Q}}-{\mathbf{Q}}_{\parallel}$,
respectively.

Figure 1: System model of general load balancing. (a) For push-based policy,
we have $m(t)=\emptyset$ for all $t$ since it does not require any memory. The
message exchange is bidirectional: probing from the dispatcher and feedback
messages from servers. (b) For pull-based policy, $m(t)$ stores the ID of the
servers that satisfy a certain condition at time $t$. The message exchange is
unidirectional, i.e., only the pull-message is sent from the servers to the
dispatcher.

##  2 Model and Definitions

This section describes a general model for the load balance systems as shown
in Fig. 1, and introduces necessary definitions.

###  2.1 Model Description

Consider a time-slotted load balancing system, with one central dispatcher and
$N$ parallel servers. These servers are indexed by its ID $n=1,2,\ldots,N$.
Each server $n$ is associated with a FIFO (first-in, first-out) queue of
length $Q_{n}(t)$ at the beginning of time slot $t$, $t=0,1,2,\ldots$. Thus,
we use index $n$ to represent both the server and the associated queue. Once a
job joins a queue, it will remain in that queue until its service is
completed.

######  Assumption 1 (Arrival Process).

Let $A_{\Sigma}(t)$ and $A_{n}(t)$ denote the number of exogenous job arrivals
and the number of arrivals routed to queue $n$ at time slot $t$, respectively.
We assume that all the exogenous arrivals at time $t$ are routed to one
selected queue $s$, using the standard model as in [2], [9], i.e.,
$A_{s}(t)=A_{\Sigma}(t),s\in\mathcal{N}=\\{1,2,\ldots,N\\}$ and $A_{i}(t)=0,$
for all $i\in\mathcal{N}\setminus\\{s\\}$. The job arrival process
$\\{A_{\Sigma}(t),t\geq 0\\}$ is a nonnegative integer valued stochastic
process that is _i.i.d_ across time $t$, with mean
$\mathbb{E}\left[A_{\Sigma}(t)\right]=\lambda_{\Sigma}$ and variance
$\text{Var}(A_{\Sigma}(t))=\sigma_{\Sigma}^{2}$. We further assume that the
number of exogenous arrivals at each time slot is bounded by a constant, i.e.,
$A_{\Sigma}(t)\leq A_{\max}<\infty$ for all $t\geq 0$.

######  Assumption 2 (Service Process).

Let $S_{n}(t)$ denote the potential service offered to queue $n$ at time $t$,
which represents the maximum number of jobs that can be served in time slot
$t$. Therefore, if the offered service $S_{n}(t)$ is larger than the number of
pending jobs in queue $n$ at time slot $t$, it will cause an unused service
$U_{n}(t)$, as defined in (1). For each $n$, the process $\\{S_{n}(t),t\geq
0\\}$ is a nonnegative integer valued _i.i.d._ stochastic process with mean
$\mathbb{E}\left[S_{n}(t)\right]=\mu_{n}$ and variance
$\text{Var}(S_{n}(t))=\nu_{n}^{2}$. Moreover,
$\lambda_{\Sigma}<\sum_{i=1}^{N}\mu_{n}$. Furthermore, the processes
$\\{S_{n}(t),t\geq 0\\},n\in\mathcal{N}$ are mutually independent across
different queues, which are also independent of the arrival processes. The
offered service $S_{n}(t)$ to each queue is uniformly bounded by a constant,
i.e., $S_{n}(t)\leq S_{\max}<\infty$ for all $t\geq 0$ and all
$n\in\mathcal{N}$.

Let ${\mathbf{Q}}(t)=\\{Q_{1}(t),\ldots,Q_{N}(t)\\}$ be the queue lengths
observed at the beginning of time $t$. Define $m(t)$ to be the set of server
IDs stored in the dispatcher at the beginning of time slot $t$. In general,
the dispatcher makes the decision of $A_{n}(t)$ based on
$({\mathbf{Q}}(t),m(t))$ for each time slot $t$. This includes the cases that
the dispatching decision depends only on ${\mathbf{Q}}(t)$ (e.g., JSQ),
partial information of ${\mathbf{Q}}(t)$ (e.g., power-of-$d$) or only on
$m(t)$ (e.g., JIQ). In each time slot, the queueing dynamics evolves according
to the following procedure. The job arrivals occur at the beginning of time
slot $t$. Then, the dispatching decision $A_{n}(t)$ is selected based on
$({\mathbf{Q}}(t),m(t))$. Further, the routed jobs are processed by the
allocated servers. Thus, the queueing dynamics is given by the following
equation,

$\begin{split}Q_{n}(t+1)&=\left[Q_{n}(t)+A_{n}(t)-S_{n}(t)\right]\\\
&=Q_{n}(t)+A_{n}(t)-S_{n}(t)+U_{n}(t),\end{split}$ (1)

where $[x]^{+}=\max(0,x)$, $U_{n}(t)=\max(S_{n}(t)-Q_{n}(t)-A_{n}(t),0)$
denotes the unused service of queue $n$.

###  2.2 Definitions

The load balancing system is modeled as a discrete-time Markov chain
$\\{Z(t)=({\mathbf{Q}}(t),m(t)),t\geq 0\\}$ with state space $\mathcal{Z}$,
using queue length vector ${\mathbf{Q}}(t)$ together with the memory state
$m(t)$. We consider a system $\\{Z^{(\epsilon)}(t),t\geq 0\\}$ parameterized
by $\epsilon$, i.e., the exogenous arrival process is
$\\{A_{\Sigma}^{(\epsilon)}(t),t\geq 0\\}$ with
$\lambda_{\Sigma}^{(\epsilon)}=\mu_{\Sigma}-\epsilon=\sum_{n}\mu_{n}-\epsilon$.
That is, we use $\epsilon$ to indicate the distance of arrival rate to the
capacity boundary, and it is also adopted as a superscript to represent the
corresponding random variables and processes.

######  Definition 1 (Stability).

$\\{Z^{(\epsilon)}(t),t\geq 0\\}$ is said to be stable if we have

$\limsup_{C\to\infty}\limsup_{t\to\infty}\mathbb{P}\left(\sum_{n}Q_{n}^{(\epsilon)}(t)>C\right)=0.$

A load balancing policy is said to be throughput optimal if it stabilizes the
system under any arrival rate in the capacity region. Since the capacity
region in our model is simply $\lambda_{\Sigma}<\mu_{\Sigma}$, the definition
of throughput optimality is given as follows.

######  Definition 2 (Throughput Optimality).

A load balancing policy is said to be throughput optimal if it stabilizes
$\\{Z^{(\epsilon)}(t),t\geq 0\\}$ for any $\epsilon>0$.

For the definition of heavy-traffic delay optimality, we need the following
definition and property.

######  Definition 3 (Resource-pooled System).

A single-server FCFS (first-come, first-serve) system
$\\{q^{(\epsilon)}(t),t\geq 0\\}$ is said to be the resource-pooled system
with respect to $\\{Z^{(\epsilon)}(t),t\geq 0\\}$, if its arrival and service
process satisfy $a^{(\epsilon)}(t)=A_{\Sigma}^{(\epsilon)}(t)$ and $s(t)=\sum
S_{n}(t)$ for all $t\geq 0$. Then, we have

$\mathbb{E}\left[q^{(\epsilon)}(t)\right]\leq\mathbb{E}\left[\sum
Q_{n}^{(\epsilon)}(t)\right],$ (2)

for all $t\geq 0$ and $\epsilon>0$.

In words, a resource-pooled system is a system that merges the total resource
of $N$ servers and queues to a single server with a single queue. Eq. (2)
holds due to the fact for any $t$, the overall arrivals to the resource-pooled
system and to load balancing system are the same, and the overall service in
the resource-pooled system is stochastically larger than the overall service
in the load balancing system. This is due to the fact that the jobs in load
balancing system cannot be moved from one queue to another, which often
results in a strict inequality in Eq. (2). However, in the heavy-traffic
regime, this lower bound can be achieved under some policy in an asymptotic
sense as defined in the next definition, and hence based on Little’s law this
policy achieves the minimum average delay of the system.

######  Definition 4 (Heavy-traffic Delay Optimality).

A load balancing policy is said to be heavy-traffic delay optimal if the
stationary workload of $\\{Z^{(\epsilon)}(t),t\geq 0\\}$ under all the arrival
and service processes in Assumptions 1 and 2, satisfies 111Assume
$(\sigma_{\Sigma}^{(\epsilon)})^{2}$ converges to a constant.

$\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\sum_{n}\overline{Q}_{n}^{(\epsilon)}\right]=\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\overline{q}^{(\epsilon)}\right],$ (3)

where $\overline{{\mathbf{Q}}}$ is the random vector whose probability
distribution is the same as the steady-state distribution of
$\\{{\mathbf{Q}}(t),t\geq 0\\}$.

######  Remark 1.

Based on the definition above, in order to show a policy, say
$\mathcal{P}_{1}$, is not heavy-traffic delay optimal, it is sufficient to
find a class of $\\{A_{\Sigma}^{(\epsilon)}(t)\\}$ and $\\{S_{n}(t)\\}$ such
that Eq. (3) does not hold. In other words, there exists a class of arrival
and service processes for which policy $\mathcal{P}_{1}$ cannot achieve the
lower bound (i.e., the resource-pooled system) while JSQ can (since it is
heavy-traffic delay optimal).

##  3 Main Results

In this section, we introduce a class $\Pi$ of load balancing policies which
are proven to be delay-optimal in the heavy-traffic regime. Popular load
balancing policies, such as JSQ and power-of-$d$, are special cases in $\Pi$;
but the JIQ policy does not belong to $\Pi$ as we will show in Theorem 2 that
it is not heavy-traffic delay optimal. In order to improve the delay
performance of JIQ while maintaining its low message overhead and simplicity,
we develop a new load balancing policy named join-below-threshold (JBT-$d$),
which is heavy-traffic delay-optimal as we can show JBT-$d$ is in $\Pi$ and
has a low message overhead similar to JIQ.

###  3.1 The Class of Load Balancing Policies $\Pi$

Let us denote $\mathbf{p}(t)=(p_{1}(t),\ldots,p_{N}(t))$, where $p_{n}(t)$ is
the probability that the new arrivals in time slot $t$ are dispatched to queue
$n$ such that $\sum_{n=1}^{N}p_{n}(t)=1$. We consider a class of load
balancing policies in which $\mathbf{p}(t)$ is a function of the system state
$Z(t)=\\{{\mathbf{Q}}(t),m(t)\\}$. Consider a permutation $\sigma_{t}(\cdot)$
of $(1,2,\ldots,N)$ that satisfies $Q_{\sigma_{t}(1)}(t)\leq
Q_{\sigma_{t}(2)}(t)\leq\ldots\leq Q_{\sigma_{t}(N)}(t)$ for all $t$, i.e.,
the queues are sorted according to an increasing order of the queue lengths in
time slot $t$ with ties broken randomly. Define
$\mathbf{P}(t)=(P_{1}(t),\ldots,P_{N}(t))$ such that $\mathbf{P}(t)$ is a
permutation of $\mathbf{p}(t)$ with $P_{n}(t)=p_{\sigma_{t}(n)}(t)$. Let

$\displaystyle\Delta_{n}(t)$
$\displaystyle=p_{\sigma_{t}(n)}(t)-{\mu_{\sigma_{t}(n)}}/{\mu_{\Sigma}}$
$\displaystyle=P_{n}(t)-{\mu_{\sigma_{t}(n)}}/{\mu_{\Sigma}}.$ (4)

######  Definition 5 (Equivalence in inner-product).

A dispatching distribution $\mathbf{\hat{P}}(t)$ is said to be _equivalent to
another dispatching distribution $\mathbf{P}(t)$ in inner product_, if

$\sum_{n}Q_{\sigma_{t}(n)}\Delta_{n}(t)=\sum_{n}Q_{\sigma_{t}(n)}\hat{\Delta}_{n}(t),$
(5)

or equivalently, if

$\sum_{n}Q_{\sigma_{t}(n)}P_{n}(t)=\sum_{n}Q_{\sigma_{t}(n)}\hat{P}_{n}(t).$
(6)

The equivalence between (5) and (6) follows immediately from (4). Intuitively
speaking, a load-balancing policy is ‘good’ if the inner product between
${\mathbf{Q}}_{\sigma_{t}}(t)$ and $\mathbf{P}(t)$ is as small as possible
such that more packets are dispatched to shorter queues. If $\mathbf{P}(t)$ is
equivalent to $\mathbf{\hat{P}}(t)$ in inner-product, we can replace
$\mathbf{\hat{P}}(t)$ by $\mathbf{P}(t)$ without affecting the property of the
policy in heavy-traffic regime, which will be explained in details later.

The following definitions enable us to distinguish different load balancing
policies based on $\mathbf{P}(t)$ or equivalently $\Delta(t)$:

######  Definition 6 (Tilted distribution).

A dispatching distribution $\mathbf{P}(t)$ is said to be _tilted_ , if there
exists $k\in\\{2,\ldots,N\\}$ such that $\Delta_{n}(t)\geq 0$ for all $n<k$
and $\Delta_{n}(t)\leq 0$ for all $n\geq k$.

######  Definition 7 ($\delta$-tilted distribution).

A dispatching distribution $\mathbf{P}(t)$ is said to be _$\delta$ -tilted_,
if (i) $\mathbf{P}(t)$ is tilted and (ii) these exists a constant $\delta>0$
such that $\Delta_{1}(t)\geq\delta$ and $\Delta_{N}(t)\leq-\delta$.

Some examples are presented in Fig. 2 to facilitate the understanding of
tilted distribution, $\delta$-tilted distribution, and equivalence in inner-
product. Fig. 2 (a)-(f) illustrate six dispatching distributions
$\mathbf{P}(t)$. The queue state ${\mathbf{Q}}(t)$ is given by (i) or (ii).
The service rates are $\mu_{A}=\mu_{B}=\mu_{C}=\mu_{D}=1$ such that
$\mu_{i}/{\mu_{\Sigma}}=1/4$ for $i=A,B,C,D$. By direct computation, one can
obtain that $P_{n}(t)$ is tilted in scenario (a), (b), (d), (e), and (f), and
is $\delta$-tilted in scenario (d), (e), and (f). If ${\mathbf{Q}}(t)$ is in
the State (i), there is no tie in the queue length and hence the permutation
$\sigma_{t}(\cdot)$ is unique, which means that $\mathbf{P}(t)$ is fully
determined by $\mathbf{p}(t)$. If ${\mathbf{Q}}(t)$ is in the State (ii), all
queue lengths are equal and hence the permutation $\sigma_{t}(\cdot)$ is non-
unique, which means that $\mathbf{P}(t)$ is determined by both $\mathbf{p}(t)$
and $\sigma_{t}(\cdot)$. In this case, however, the inner product between
${\mathbf{Q}}_{\sigma_{t}}(t)$ and $\mathbf{P}(t)$ is one in all (a)-(f), and
hence the dispatching distributions $\mathbf{P}(t)$ in (a)-(f) are mutually
equivalent in inner product. For example, in this case even though
$\mathbf{P}(t)$ in (c) is neither tilted nor $\delta$-tilted, it is equivalent
in inner product to $\mathbf{P}(t)$ in (d) which is both tilted and
$\delta$-tilted.

From the perspective of heavy-traffic delay performance, tilted distribution
is a dispatching distribution that is not worse than random routing and
$\delta$-tilted distribution is a dispatching distribution that is strictly
better than random routing. In addition, the equivalence in inner-product
allows us to transfer a tilted dispatching distribution to a $\delta$-tilted
dispatching distribution when there are ties in queue lengths, that is, it
allows to merges probability in $\mathbf{P}(t)$ from longer queues to shorter
queues without changing the inner product.

Figure 2:  Illustrations of tilted distribution, $\delta$-tilted distribution,
and equivalence in inner-product.

We now introduce a class of load balancing algorithms $\Pi$ based on the
property of $\mathbf{P}(t)$ or its equivalent distributions in inner product.

######  Definition 8.

A load balancing policy is said to belong to class $\Pi$ if it satisfies the
following two conditions:

  1. (i)

$\mathbf{P}(t)$ or one of its equivalent distributions in inner product is
tilted for all $Z(t)$ and $t\geq 0$.

  2. (ii)

For some finite positive constants $T$ and $\delta$ that both are independent
of $\epsilon$, there exists a time slot
$t_{k}\in\\{kT,kT+1,\ldots,(k+1)T-1\\}$ for each $k\in\mathbb{N}$ such that
$\mathbf{P}(t_{k})$ or one of its equivalent distributions in inner product is
$\delta$-tilted for all $Z(t_{k})$.

In the sequel, we will show that any policy in $\Pi$ satisfies the following
two sufficient conditions for throughput and heavy-traffic delay optimality,
which are obtained via the Lyapunov-drift based approach developed in [2].

######  Lemma 1.

If there exist $T_{1}>0$, $K_{1}\geq 0$, and $\gamma>0$ such that for all
$t_{0}=1,2,\ldots$, all $Z\in\mathcal{Z}$ and $\lambda_{\Sigma}<\mu_{\Sigma}$

$\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{1}-1}\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t_{0})=Z\right]\leq-\gamma\left\lVert{\mathbf{Q}}\right\rVert+K_{1},$ (7)

then the system is throughput-optimal. Moreover, the stationary distribution
of the queueing system has bounded moments, i.e., there exist finite $M_{r}$
such that for all $\epsilon>0$ and $r\in\mathbb{N}$

$\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}^{(\epsilon)}\right\rVert^{r}\right]\leq
M_{r}.$

###### Proof.

See Section A of the Appendix.∎

######  Lemma 2.

Under the assumptions of Lemma 1, if there further exist $T_{2}>0$, $K_{2}\geq
0$ and $\eta>0$ that are independent of $\epsilon$, such that for all
$t_{0}=1,2,\ldots$ and all $Z\in\mathcal{Z}$

$\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{2}-1}\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t_{0})=Z\right]\leq-\eta\left\lVert{\mathbf{Q}}_{\perp}\right\rVert+K_{2}$
(8)

holds for all $\epsilon\in(0,\epsilon_{0})$, $\epsilon_{0}>0$, where
${\mathbf{Q}}_{\perp}={\mathbf{Q}}-\langle{\mathbf{Q}},\mathbf{c}\rangle\mathbf{c}$
is the perpendicular component of ${\mathbf{Q}}$ with respect to the line
${\bf c}=\frac{1}{\sqrt{N}}{(1,1,\ldots,1)}$, then the system is heavy-traffic
delay optimal, i.e.,

$\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\sum_{n}\overline{Q}_{n}^{(\epsilon)}\right]=\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\overline{q}^{(\epsilon)}\right].$

###### Proof.

See Section B of the Appendix.∎

######  Remark 2.

Note that these two sufficient conditions distilled from the Lyapunov-drift
based approach not only provide a unified approach for throughput and heavy-
traffic optimality analysis, but also enable us to abstract a class of heavy-
traffic delay optimal policies. In particular, using Lemma 1 and Lemma 2, we
are able to prove the main result of this paper.

######  Theorem 1.

Any load balancing policy in $\Pi$ is throughput optimal and heavy-traffic
delay optimal.

###### Proof sketch of Theorem 1.

The insight for a policy in $\Pi$ to satisfy the sufficient condition in Eq.
(7) is that under tilted dispatching distribution the performance is no worse
than random dispatching. This follows from the following property of tilted
distribution

$\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\leq 0.$ (9)

The equality is obtained when all $\Delta_{n}(t)$ is zero, which is the case
of random dispatching as shown in (b) of Fig. 2. Note that for all other cases
of a tilted distribution, Eq. (9) is strictly less than zero. This is true
since $\sum_{n=1}^{N}\Delta_{n}(t)$ is always zero and the permutation is in
the non-decreasing order of the queue length.

The intuition for a policy in $\Pi$ to satisfy the sufficient condition in Eq.
(8) is that the performance under any $\delta$-tilted dispatching distribution
is strictly better than random dispatching, under which the term in Eq. (8) is
$0$ for homogeneous servers and of order $\epsilon$ for heterogeneous servers.
Note that under a $\delta$-tilted distribution, we have

$\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\leq-\delta(Q_{\sigma_{t}(N)}(t)-Q_{\sigma_{t}(1)}(t).$
(10)

This inequality comes from the definition of the $\delta$-tilted distribution
and fact that the permutation is in the non-decreasing order of the queue
length. In order to have the term of
$\left\lVert{\mathbf{Q}}_{\perp}\right\rVert$, the following inequality would
be quite useful

$\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\leq\sqrt{N}(Q_{\sigma_{t}(N)}(t)-Q_{\sigma_{t}(1)}(t)).$
(11)

This is true since
${\mathbf{Q}}_{\perp}(t)={\mathbf{Q}}(t)-{\mathbf{Q}}_{\parallel}(t)={\mathbf{Q}}(t)-\frac{\sum
Q_{n}(t)}{N}\mathbf{1}={\mathbf{Q}}(t)-Q_{\text{avg}}(t)\mathbf{1}$, in which
$Q_{\text{avg}}(t)$ is the average queue length among the $N$ servers at time
slot $t$.

The details of the proof are presented in Section 5.1. ∎

From Eqs. (9) and (10), it can be seen that the important property of a given
policy is fully characterized by the inner product of
${\mathbf{Q}}_{\sigma_{t}}(t)$ and $\Delta(t)$ under the system state $Z(t)$,
which is actually the motivation to define equivalent distribution in inner
product. That is, even though the dispatching distribution $\mathbf{P}(t)$ is
not unique when there are ties in queue lengths, the inner product is actually
the same if two dispatching distributions are equivalent in inner product,
hence preserving the same property in heavy-traffic regime.

Note that class $\Pi$ is sufficient but not necessary for heavy-traffic delay
optimality. Nevertheless, in the next section, we will show that it not only
contains many well-known heavy-traffic delay optimal policies but also allows
us to design new heavy-traffic delay optimal policies which enjoy nice
features of pull-based policies.

###  3.2 Important Policies in $\Pi$

####  3.2.1 Join-shortest-queue (JSQ) policy

Under JSQ policy, all the incoming jobs are dispatched to the queue that has
the shortest queue length, ties are broken uniformly at random, out of all the
servers.

######  Proposition 1.

The JSQ policy belongs to $\Pi$, and hence is throughput optimal and heavy-
traffic delay optimal.

The result that JSQ is throughput and heavy traffic delay optimal has been
first proven via diffusion limits for two servers in [3] and via Lyapunov-
drift argument for $N$ servers in [2]. Here, we present another simple proof
based on our main result.

###### Proof.

Note that when there are no ties in queue lengths, the dispatching
distribution $\mathbf{P}(t)$ under JSQ satisfies that for all $t$

$P_{1}(t)=1\text{ and }P_{n}(t)=0,2\leq n\leq N.$ (12)

In other words, all the arrivals are dispatched to the shortest queue, which
is always the queue $\sigma_{t}(1)$ if there are no ties in queue lengths. If
there are ties in queue lengths, this $\mathbf{P}(t)$ is equivalent in inner
product to other dispatching distribution under the state $Z(t)$ in which ties
exist. In particular, if there are $m\leq N$ queues that all have the shortest
queue length, then in this case by random routing the dispatching distribution
under JSQ is given by $\hat{P}_{i}=\frac{1}{m}$ for all $1\leq i\leq m$, and
$\hat{P}_{i}=0$ for all $i>m$. It can be seen that $\mathbf{P}(t)$ in Eq. (12)
is equivalent in inner product to $\mathbf{\hat{P}}(t)$ according to the
definition because
$Q_{\sigma_{t}(1)}=Q_{\sigma_{t}(2)}=\ldots=Q_{\sigma_{t}(m)}$. Thus, for all
$Z(t)$, under JSQ the dispatching distribution or its equivalent distribution
in inner product is in the form of Eq. (12). Hence, we have
$\Delta_{1}(t)=1-\mu_{\sigma_{t}(1)}/\mu_{\Sigma}>0$, and
$\Delta_{n}(t)=-{\mu_{\sigma_{t}(n)}}/{\mu_{\Sigma}}<0$ for all $2\leq n\leq
N$, which implies that $P_{n}(t)$ is a $\delta$-tilted probability with
$\delta=\mu_{min}/{\mu_{\Sigma}}$ for all $Z(t),t\geq 0$, where
$\mu_{min}=\min_{n\in\mathcal{N}}\mu_{n}$. Therefore, the JSQ policy is
contained in the class $\Pi$ under both heterogeneous and homogeneous servers.
∎

####  3.2.2 The power-of-$d$ policy

Under the power-of-$d$ policy, all the incoming jobs are dispatched to the
queue that has the shortest queue length, ties are broken uniformly at random,
out of $d\geq 2$ servers, which are chosen uniformly at random.

######  Proposition 2.

The power-of-$d$ policy belongs to $\Pi$ under homogeneous servers, and hence
is throughput-optimal and heavy-traffic delay optimal.

The power-of-$d$ policy has been proven to be heavy-traffic delay optimal via
Lyapunov drift condition in [9]. Here, we will present another proof based on
our main result.

###### Proof.

Note that when there are no ties in queue lengths, the dispatching
distribution $\mathbf{P}(t)$ under the power-of-$d$ policy satisfies that for
all $t\geq 0$

$P_{n}(t)={\binom{N-n}{d-1}}\bigg{/}{\binom{N}{d}},1\leq n\leq N-d+1,$ (13)

and $P_{n}(t)=0$, for all $n>N-d+1$. This comes from the fact that all
arrivals are dispatched to the queue with shortest queue length among $d$
uniformly randomly sampled servers. Thus, if the queue $\sigma_{t}(n)$ is the
one with shortest queue length among $d$ samples, the remaining $d-1$ samples
must come from queues $\sigma_{t}(n+1)$, $\sigma_{t}(n+2)$,
$\ldots\sigma_{t}(N)$ if all the queue lengths are different in $Z(t)$. If
there are ties in queue lengths, it can be easily shown that this
$\mathbf{P}(t)$ is equivalent in inner product to other dispatching
distributions under any given $Z(t)$ in which there are ties in queue lengths.
Thus, for all $Z(t)$, the dispatching distribution or its equivalent
distribution in inner product under the power-of-$d$ policy can be fully
determined by Eq. (13). Since $P_{n}(t)$ is decreasing and
$\mu_{\sigma_{t}(n)}=\mu$ under homogeneous servers, $\mathbf{P}(t)$ is a
tilted distribution. Note that $\Delta_{1}(t)=\frac{d-1}{N}$ and
$\Delta_{N}(t)=-\frac{1}{N}$. As a result, $\mathbf{P}(t)$ is a
$\delta$-tilted distribution with $\delta=\frac{1}{N}$ for all $Z(t)$, which
implies that power-of-$d$ policy is included in the class $\Pi$ for
homogeneous servers. ∎

####  3.2.3 Join-idle-queue policy is not in $\Pi$

Now we will show that the JIQ policy is not contained in the class $\Pi$
because it is in fact not heavy-traffic delay optimal in homogeneous servers.
For the heterogeneous case, it is well-known that JIQ is not heavy-traffic
delay optimal since it is not even throughput optimal for a fixed number of
servers [15]. However, for the homogeneous case, it is still open whether it
is heavy-traffic optimal for a fixed number of servers, although it has been
shown to be heavy-traffic optimal when the number of servers goes to infinity
in the Halfin-Whitt regime [12]. It turns out that when the number of servers
is fixed, there exists a class of arrival process, under which the delay
performance of JSQ is strictly better than that of JIQ in the heavy-traffic
limit. More specifically, as shown in the proof of Theorem 2, for a class of
arrival process, the delay under JIQ cannot achieve the common lower bound
(i.e., the resource-pooled system), while JSQ can, which implies that JIQ is
not heavy-traffic delay optimal for homogeneous case.

In particular, we consider the two-server case with constant service process
with rate $1$. We are able to find a class of arrival process such that Eq.
(3) under JIQ does not hold. Let us first introduce the class of arrival
process $\mathcal{A}$.

######  Definition 9.

An arrival process $A_{\Sigma}(t)$ is said to belong to $\mathcal{A}$ if

  1. (i)

$\mathbb{P}(A_{\Sigma}^{(\epsilon)}(t)=0)=p_{0}$, which is a constant
independent of $\epsilon$.

  2. (ii)

$(\sigma_{\Sigma}^{(\epsilon)})^{2}$ approaches a constant
$\sigma_{\Sigma}^{2}$ which satisfies that
$\sigma_{\Sigma}^{2}>\frac{8}{p_{0}}-4$.

More concretely, we are able to show the following result.

######  Theorem 2.

Consider a load balancing system with two homogeneous servers, JIQ is not
heavy-traffic delay optimal in this case.

###### Proof.

See Section C of the Appendix. ∎

###  3.3 Designing New Policies in $\Pi$

It has been shown in the last section that the state-of-art push-based
policies, e.g., JSQ and power-of-$d$, are all included in $\Pi$. Recall that,
both of them need to sample the queue length information upon each new
arrival, which directly results in the following two problems.

  1. (a)

The message exchange rate between dispatcher and servers is high, especially
for join-shortest-queue.

  2. (b)

Each arrival has to wait for completion of the message exchange before being
dispatched, which increases the actual response time for each job.

To resolve the problem, the pull-based policies, join-idle-queue (JIQ) in [8]
and an equivalent algorithm called PULL in [15] are proposed, which have been
shown to enjoy low message rate (at most one message per job) and have a
better performance than the power-of-$d$ policy from light to moderate loads.
However, as shown via numerical results in [8] and the proof of Theorem 2 in
this paper, when the load becomes high, the performance of JIQ is much worse
than the power-of-$d$ policy, which motivates us to design policies that enjoy
low message rates, while still guaranteeing throughput and heavy-traffic delay
optimality.

Table 1: Summary of load balancing policies Policy | Message | Throughput-Optimal | Heavy-traffic Delay-Optimal
---|---|---|---
Homogeneous | Heterogeneous | Homogeneous | Heterogeneous
Random | 0 | $\surd$ | $\times$ | $\times$ | $\times$
JSQ [2] | 2$N$ | $\surd$ | $\surd$ | $\surd$ | $\surd$
Power-of-$d$ [9], [10] | 2$d$ | $\surd$ | $\times$ | $\surd$ | $\times$
JIQ[15], [8] |  $\leq$ 1 | $\surd$ | $\times$ | $\times$ | $\times$
JBT-$d$ | $\leq\frac{N+2d}{T}+1$ | $\surd$ | $\times$ | $\surd$ | $\times$
JBTG-$d$ | $\leq\frac{N+2d}{T}+1$ | $\surd$ | $\surd$ | $\surd$ | $\surd$
*The message rate for JBT-$d$ and JBTG-$d$ in this table is just a crude upper bound. When the new threshold is larger
than the old one, there is no need for the servers that are already recorded
in memory to resend pull-messages.

######  Definition 10.

Join-below-threshold-$d$ (JBT-$d$) policy is composed of three components:

  1. 1.

A threshold is updated every $T$ units of time by uniformly at random sampling
$d$ servers, and taking the shortest queue length among the $d$ servers as the
new threshold.

  2. 2.

Each server sends its ID to the dispatcher when its queue length is not larger
than the threshold for the first time.

  3. 3.

Upon a new arrival, the dispatcher checks the available IDs in the memory. If
they exist, it removes one uniformly at random, and sends all the new arrivals
to the corresponding server. Otherwise, all the new arrivals will be
dispatched uniformly at random to one of the servers in the system.

To be more specific, we explain the connections of the three components as
follows. At the beginning of each time slot, the dispatcher immediately routes
the new arrivals to a server only based on its memory state, i.e., no
sampling. If there are available IDs in memory, it removes one uniformly at
random and sends the newly arrived jobs to the corresponding server.
Otherwise, it sends the new jobs to a server selected uniformly at random
among all the servers. At the end of each time slot, if there is no update of
threshold, each server will immediately report its ID if its queue length is
not larger than the threshold for the first time, i.e., only reporting once
for each server before dispatched. Otherwise, the dispatcher updates the
threshold by uniformly at random sampling $d$ servers, and the new threshold
is set as the shortest queue length among $d$ samples. Then, each server
decides to whether or not to report based on its queue length and the new
threshold, using the same way as before.

######  Definition 11.

The JBT-$d$ policy can be easily generalized for heterogeneous servers, denote
by JBTG-$d$, as follows. The only difference is that the dispatching
probability distribution for the case of non-empty and empty memory is given
by

$\boldsymbol{\psi}_{i}(t):=\frac{\mu_{i}}{\sum_{j\in
m(t)}\mu_{j}}{\mathbbm{1}}_{\\{i\in m(t)\\}}\text{ and
}\boldsymbol{\phi}_{i}(t):=\frac{\mu_{i}}{\mu_{\Sigma}}\text{ for all }i$

That is, the probability to be selected for a server that has its ID in memory
is weighted by its service rate. This can be easily done by requiring the
server to report its service rate $\mu_{n}$ as well as its ID.

In the following, we will show that JBT-$d$ and JBTG-$d$ belong to $\Pi$, and
hence throughput and heavy-traffic delay optimal. More specifically, we have
the following result.

######  Proposition 3.

For any finite $T$ and $d\geq 1$, the following two assertions are true:

  1. 1.

JBT-$d$ is in $\Pi$ for homogeneous servers, and hence throughput and heavy-
traffic delay optimal.

  2. 2.

JBTG-$d$ is in $\Pi$ for both homogeneous and heterogeneous servers, and hence
throughput and heavy-traffic delay optimal.

###### Proof sketch of Proposition 3.

Let us look at JBT-$d$ for some key insights behind this proof. In order to
show it is in $\Pi$, we only need to show that it satisfies the two conditions
(i) and (ii). For the condition (i), we will show that at any time slot $t$,
the dispatching is no worse than the random routing. For the condition (ii),
we will show that at time slots $rT+1$, $r\in\\{0,1,2,\ldots\\}$, the
dispatching decision is strictly better than the random routing.

Note that under the JBT-$d$ policy, if the ID of the server $\sigma_{t}(n+1)$
is in $m(t)$, we must have that the ID of the server $\sigma_{t}(n)$ is also
in $m(t)$ as the permutation is in the non-decreasing order of the queue
length. Denote by $\tilde{p}_{k}(t)$ the probability that there are $k$ IDs in
the memory $m(t)$ for time $t$, i.e., $\tilde{p}_{k}(t)=\Pr(|m(t)|=k)$. Then,
the probability for the server $\sigma_{t}(n)$ to be selected at time $t$,
i.e., $P_{n}(t)$ is given by

$P_{n}(t)=\sum_{i=n}^{N}\tilde{p}_{i}(t)\frac{1}{i}.$ (14)

This is true since for the server $\sigma_{t}(n)$ to be selected, there should
be at least $n$ IDs in memory, i.e., $|m(t)|\geq n$ and in each case the
probability for the server $\sigma_{t}(n)$ to be chosen is $\frac{1}{|m(t)|}$.
Therefore, we can see that the probability of $P_{n}(t)$ satisfies

$P_{1}(t)\geq P_{2}(t)\geq\ldots\geq P_{N}(t),$ (15)

which directly implies that for all $t\geq 0$ there exists a $k$ between $2$
and $N$ such that $\Delta_{n}(t)=P_{n}(t)-\frac{1}{N}\geq 0$ for all $n<k$ and
$\Delta_{n}(t)\leq 0$ for all $n\geq k$. Therefore, condition (i) of $\Pi$ is
satisfied.

For condition (ii), we will show that there exists a lower bound for $\delta$
such that $\mathbf{P}(rT+1)$ (or an inner product equivalent distribution when
there are ties in queue lengths), $r\in\\{0,1,2,\ldots\\}$ is at least a
$\delta$-tilted distribution. In this case, we need only to show that
$P_{N}(rT+1)$ is strictly less than $\frac{1}{N}$ for all the system state
$Z(rT+1)$.

The full proof is presented in Section 5.2. ∎

###  3.4 Features of JBT-d

This section summarizes the main features of JBT-d policy and compares it with
existing policies in Table 1. In particular, we compare the number of messages
for each new arrival under different policies. For push-based polices, e.g.,
JSQ and power-of-$d$, there are $d$ query and $d$ response messages for each
new arrival ($d=N$ for JSQ policy). For JIQ policy, for each new arrival, it
requires at most one pull-message since when there are no pull-messages in
memory, the arrival is dispatched randomly without costing any pull-message.
Similarly, our JBT-$d$ policy requires $2d$ push-messages every $T$ time slots
to update the threshold. Due to the threshold update, the old pull-messages
may be discarded, which is upper bounded by $N$. Hence, the pull-message for
each new arrival under JBT-$d$ is at most $1+\frac{2d+N}{T}$.

In sum, the JBT-$d$ policy has the following nice features: a) It is
throughput and heavy-traffic delay optimal. b) It is able to guarantee heavy-
traffic delay optimal with very low complexity when $T$ is relatively large.
c) The computation overhead is small. It only needs to keep a list of the
available IDs and choose randomly. d) The arrival is immediately dispatched,
i.e., there is no dispatching delay as compared to push-based policies such as
JSQ and Power-of-$d$.

It is worth pointing that by just changing the way of updating the threshold
in JBT-d, we can design other new policies which also enjoy the nice features
above. For example, it can be easily shown via similar arguments that if the
threshold is updated by sampling all the servers and taking the average value
of the queue length as the new threshold, this corresponding new policy is
still in the class $\Pi$.

##  4 Numerical Results

In this section, we use simulations to compare our proposed policies JBT-$d$
and JBTG-$d$ with join-shortest-queue (JSQ), join-idle-queue (JIQ), power-
of-$d$ (SQ($d$)) and power-of-$d$ with memory (SQ($d$,$m$)). The power-of-$d$
with memory policy (SQ($d$,$m$)) improves power-of-$d$ by using extra memory
to store the $m$ shortest queues sampled at the previous time slot [11].

We compare the throughput performance, delay performance, heavy-traffic delay
performance and message overhead performance under various arrival and service
processes as well as different system sizes. Moreover, the $95\%$ confidence
intervals for all the simulation results can be found in Section E of the
Appendix, which justify the accuracy of the simulation results. Some of the
confidence intervals are also included in figures and similar accuracy goes
for other points as well. The exogenous arrival $A_{\Sigma}(t)$ and potential
service $S_{n}(t)$ are drawn from a Poisson distribution with rate
$\lambda_{\Sigma}$ and $\mu_{n}$ for each time slot unless otherwise
specified. The traffic load is equal to $\rho=\lambda_{\Sigma}/\mu_{\Sigma}$.
The parameter $T$ is the threshold update interval for JBT-$d$ and JBTG-$d$.

Before we dive into each case, let us first summarize some key observations in
these simulation results.

  1. (i)

Throughput performance:

    1. (a)

Our proposed policy JBT-$d$ continues to stabilize all the considered loads in
heterogeneous systems under all the different settings.

    2. (b)

JIQ and SQ($d$) cannot stabilize the system when the load is high in all the
cases.

    3. (c)

JIQ appears to have a larger capacity region as the number of servers
increases. This agrees with the theoretical result in [15].

  2. (ii)

Delay performance:

    1. (a)

Our proposed policy JBT-$d$ continues to have good performance from light to
heavy traffic among all the cases.

    2. (b)

As the system size increases, JBT-$d$ achieves the same performance as JSQ for
a larger range of loads. Meanwhile, the gains of JBT-$d$ over SQ($d$) and
SQ($d$,$m$) become larger as the system size increases.

    3. (c)

The gain of JBT-$d$ over JIQ decreases as the system size increase. This is
also intuitive since as $N$ goes to infinity, it is more likely to find an
idle server, which results in the fact that JIQ is heavy-traffic delay optimal
in the Halfin-Whitt regime [12].

    4. (d)

The gain of JBT-$d$ over JIQ increases when burstness is introduced in arrival
or service process. This agrees with the insight in the proof of Theorem 2
that larger variance of arrival or service process will degrade the
performance of JIQ.

  3. (iii)

Message overhead performance:

    1. (a)

Our proposed policy JBT-$d$ continues to have a low message overhead among all
the cases.

    2. (b)

Push-based policies such as SQ($d$) and SQ($d$,$m$) have to increase their
message overhead linearly with respect to $d$ to achieve good delay
performance as the system size increases. In contrast, our proposed JBT-$d$ is
able to achieve good performance with a message rate that is less than $1$ for
all the cases when $T$ is large.

  4. (iv)

Confidence interval:

    1. (a)

Our proposed policy JBT-$d$ continues to have good $95\%$ confidence intervals
in all the various settings.

###  4.1 Throughput Performance

We investigate the throughput region of different load balancing policies in
the case of heterogeneous servers. In particular, we consider the case that
the system consisting of two server pools each with five servers and the rates
are $1$ and $10$, respectively. A turning point in the curve indicates that
the load approaches the throughput region boundary of the corresponding
policy.

Figure 3: Throughput performance in $10$ heterogeneous servers.

Figure 3 shows that the system becomes unstable when $\rho>0.5$ under the
policy power-of-2 (SQ($2$)), and it becomes unstable under JIQ when
$\rho>0.9$. In contrast, our proposed JBTG-$d$ policy remains stable for all
the considered loads which agrees with the theoretical results. It can be seen
that JBT-$2$ is also able to stabilize the system for all the considered loads
in this case. Note that the system remains stable under the power-of-$2$ with
memory policy SQ($2$,$3$), which demonstrates the benefit of using memory to
obtain maximum throughput as first discussed in [14].

We further provide additional simulation results on throughput performance
under different arrival and service process as well as different system sizes
in Section D of the Appendix.

###  4.2 Delay Performance

We investigate the mean response time under different load balancing policies
in homogeneous servers with different system sizes and various arrival and
service processes. The time interval for threshold update of JBT-$d$ is set
$T=1000$.

Let us first look at the regime when $\rho$ is from $0.3$ to $0.99$, which
ranges from light traffic to heavy traffic. Figure 5 shows that our proposed
policy JBT-$d$ outperforms both power-of-$2$ and power-of-$2$ with memory
(SQ($2$,$3$), which uses the same amount of memory as in JBT-$d$) for nearly
the whole regime. Moreover, JBT-$d$ policy achieves nearly the same response
time of JIQ when the load is not too high. However, as the load becomes
heavier, the performance of JIQ gets worse and worse, and its mean response
time is as large as two times of the response time under JBT-$d$ policy when
the load is $0.99$.

Figure 4: Delay performance under $10$ homogeneous servers.

Figure 5: Heavy-traffic delay performance under $10$ homogeneous servers.

Figure 6: Heavy-traffic delay performance under $50$ homogeneous servers.

Figure 7: Heavy-traffic delay performance under $10$ homogeneous servers with
Poisson arrival and bursty service.

Now, let us get a closer look at the delay performance in heavy-traffic
regime, i.e., $\rho>0.9$, as shown in Figure 5. It can be seen that JBT-$10$
outperforms JIQ when $\rho>0.9$ and JBT-$2$ outperforms JIQ when $\rho>0.95$
in this case. More importantly, the gap between them keeps increasing as the
load gets higher. Note that power-of-$2$ with memory (SQ($2$,$3$)) also has
good performance in this case, which, however, uses a much higher message rate
compared to our JBT-$d$ policy, as discussed in the next section.

Last, we further provide some results on heavy-traffic delay performance for a
larger system size and a bursty service process, respectively. Due to space
limitation, the comprehensive results can be found in Section D of Appendix.
Figure 7 illustrates the heavy-traffic performance under Poisson arrival and
Poisson service when $N=50$. In this case, first thing to note is that even
though the power-of-$d$ with memory policy (SQ($2$,$9$)) uses the same amount
of memory as in JBT-$d$, it has a much poorer performance with a much higher
message overhead since the message overhead of JBT-$d$ is strictly less than
$1$ when $T=1000$ in this case. This means that to improve delay performance
in large system size, power-of-$d$ with memory has to increase its message
overhead linearly with respect $d$, while our JBT-$d$ policy is able to
achieve good performance with message rate less than $1$ even for $d=N$.
Moreover, as $\rho$ approaches to $1$, the performance of JIQ degrades
substantially while our proposed JBT-$d$ remains quite close to JSQ. In Figure
7, the potential number of jobs served in each time slot is either $0$ or
$10$. In this bursty service case, JIQ degrades much faster than that in the
Poisson service process. Moreover, in this setting we can easily observe the
difference between non-heavy-traffic policy (JIQ) and heavy-traffic optimal
policies (all the others). Note that the message overhead of SQ($2$,$3$) is
nearly $8$ times as large as that of JBT-$d$, as shown in next section, though
its delay is slightly better than JBT-$d$.

Figure 8: Message per new job arrival under $10$ homogeneous servers with
respect to $T$.

Figure 9: Delay performance under $10$ homogeneous servers with respect to
$T$.

###  4.3 Message Overhead

We use simulations to further show the low message rate of our proposed
JBT-$d$ policy, though a crude upper bound has been established. Here, we
consider the $10$ homogeneous servers with Poisson arrival and Poisson
service, and more results for different settings can be found in Section D of
Appendix. More specifically, we investigate the impact of different values of
$T$, i.e., the time interval for updating the threshold, on the message rate
and its corresponding delay performance at a fixed load $\rho=0.99$. In
particular, we calculate the average number of messages per new job arrival
under each policy. For push-based policies, e.g., JSQ, power-of-$2$ (SQ(2))
and power-of-$2$ with memory (SQ($2$,$3$)), the message only includes the
push-message and is easily calculated as $20$, $4$, and $4$, which is
independent of $T$. For JIQ, we know that the rate is at most one for each new
job arrival, which is also independent with $T$ and serves as the benchmark.

Figure 9 shows the message rate of JBT-$d$ with respect to $T$ for different
values of $d$, and the corresponding delay performance is shown in Figure 9.
The first thing to note is that the message rate of JIQ is much smaller than
one since the traffic is heavy and hence there are few idle servers in this
case, which directly results in the poor performance in the heavy-traffic
regime. Second, the message rate of JBT-$d$ is smaller than all push-based
policies and becomes less than one when $T>100$ in this case, which means that
it is able to achieve throughput and heavy-traffic delay optimality by
requiring a slightly more message than JIQ. Moreover, it can be seen that as
$T$ increases, there is no significant change of the delay performance, which
indicates that we are allowed to adopt a sufficiently large $T$ while not
incurring the loss of performance very much in this case. Last, it is worth
noting that a larger $d$ does not necessarily mean a larger message overhead
when $T$ is large. This is because when $T$ is large, the push-message in
JBT-$d$ will be dominated by the pull-message. For a small $d$, the number of
pull-message may be larger since the threshold may be higher than that under a
larger $d$. As shown in the additional results in Section D of Appendix, the
observations above hold almost for all the considered cases. The exact impact
and relationship of $T$ and $d$ would be one of our future research focus.

##  5 Proof of Main Results

The high-level insight for class $\Pi$ to be heavy-traffic delay optimal is
that it always has a preference to shorter queues in the way that is specified
by the $\delta$-tilted distribution. The key step behind the proof that
JBT-$d$ is heavy-traffic delay optimal is to show that the dispatching
distribution for the time slot that is immediately after the threshold update
is always a $\delta$-titled distribution.

###  5.1 Proof of Theorem 1

Before we adopt the sufficient conditions in Lemma 1 and Lemma 2 to prove
Theorem 1, we first present the following lemmas on the tilted distribution
and $\delta$-tilted distribution, respectively.

######  Lemma 3.

For a system with mean arrival rate $\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$
and a tilted distribution $\mathbf{P}(t)$ under $Z(t)$, we have

$\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\leq-\epsilon\frac{\mu_{min}}{\mu_{\Sigma}}\left\lVert{\mathbf{Q}}(t)\right\rVert$
(16)

and

$\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\leq\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert$
(17)

###### Proof.

Consider the left-hand-side (LHS) of Eq. (16)

$\begin{split}&\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\\\
=&\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\left[(\Delta_{n}(t)+\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\lambda_{\Sigma}-\mu_{\sigma_{t}(n)}\right]\\\
\mathop{=}\limits^{(a)}&\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\lambda_{\Sigma}+\sum_{n}^{N}Q_{\sigma_{t}(n)}(t)(-\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\\\
\mathop{\leq}\limits^{(b)}&\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)(-\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\\\
\mathop{\leq}\limits^{(c)}&-\epsilon\frac{\mu_{min}}{\mu_{\Sigma}}\left\lVert{\mathbf{Q}}(t)\right\rVert,\end{split}$

in which equation (a) holds since $\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$;
inequality (b) comes from the fact that
$\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\leq 0$ under a tilted
distribution. This fact is true since $Q_{\sigma_{t}(1)}(t)\leq
Q_{\sigma_{t}(2)}(t)\leq\ldots\leq Q_{\sigma_{t}(N)}(t)$ and
$\sum_{n=1}^{N}\Delta_{n}(t)=0$; inequality (c) follows from the fact that
$\left\lVert\mathbf{x}\right\rVert_{1}\geq\left\lVert\mathbf{x}\right\rVert$
for any $\mathbf{x}\in\mathbb{R}^{N}$.

Note that
${\mathbf{Q}}_{\perp}(t)={\mathbf{Q}}(t)-{\mathbf{Q}}_{\parallel}(t)={\mathbf{Q}}(t)-\frac{\sum
Q_{n}(t)}{N}\mathbf{1}={\mathbf{Q}}(t)-Q_{\text{avg}}(t)\mathbf{1}$, in which
$Q_{\text{avg}}(t)$ is the average queue length among the $N$ servers at time
slot $t$. Then, consider the left-hand-side (LHS) of Eq. (17)

$\begin{split}&\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\\\
=&\sum_{n=1}^{N}(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))\left[(\Delta_{n}(t)+\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\lambda_{\Sigma}-\mu_{\sigma_{t}(n)}\right]\\\
\mathop{=}\limits^{(a)}&\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\lambda_{\Sigma}+\sum_{n}^{N}(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))(-\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\\\
\mathop{\leq}\limits^{(b)}&\sum_{n=1}^{N}(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))(-\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\\\
\mathop{\leq}\limits^{(c)}&\epsilon\sum_{n=1}^{N}|(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))|\\\
\mathop{\leq}\limits^{(d)}&\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert,\end{split}$
(18)

where equation (a) comes from the facts that $\sum_{n=1}^{N}\Delta_{n}(t)=0$
and $\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$; inequality (b) holds since
$\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\leq 0$ under a tilted
distribution; inequality (c) is true since $x\leq|x|$ for any $x\in\mathbb{R}$
and $|\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}}|\leq\epsilon$ for all
$n\in\mathcal{N}$; inequality (d) is true since
$\left\lVert\mathbf{x}\right\rVert_{1}\leq\sqrt{N}\left\lVert\mathbf{x}\right\rVert$
for any $\mathbf{x}\in\mathbb{R}^{N}$. ∎

######  Lemma 4.

For a system with mean arrival rate $\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$
and a $\delta$-tilted distribution $\mathbf{P}(t)$ under $Z(t)$, we have

$\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\leq\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert(\epsilon-\frac{\delta\lambda_{\Sigma}}{N})$
(19)

###### Proof.

Consider the left-hand-side (LHS) of Eq. (19), we have

$\begin{split}&\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\\\
=&\sum_{n=1}^{N}(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))\left[(\Delta_{n}(t)+\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\lambda_{\Sigma}-\mu_{\sigma_{t}(n)}\right]\\\
\mathop{=}\limits^{(a)}&\sum_{n=1}^{N}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\lambda_{\Sigma}+\sum_{n=1}^{N}(Q_{\sigma_{t}(n)}(t)-Q_{avg}(t))(-\epsilon\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}})\\\
\mathop{\leq}\limits^{(b)}&\sum_{n=1}^{N}{}Q_{\sigma_{t}(n)}(t)\Delta_{n}(t)\lambda_{\Sigma}+\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\\\
\mathop{\leq}\limits^{(c)}&-\lambda_{\Sigma}\delta(Q_{\sigma_{t}(N)}(t)-Q_{\sigma_{t}(1)}(t))+\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\\\
\mathop{\leq}\limits^{(d)}&-\lambda_{\Sigma}\frac{\delta}{\sqrt{N}}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert+\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\\\
=&\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert(\epsilon-\frac{\delta\lambda_{\Sigma}}{N})\end{split}$

where equation (a) holds since $\sum_{n=1}^{N}\Delta_{n}(t)=0$ and
$\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$; inequality (b) follows from steps
(c) and (d) in Eq. (18); inequality (c) follows from the definition of
$\delta$-tilted probability and the fact that $Q_{\sigma_{t}(1)}(t)\leq
Q_{\sigma_{t}(2)}(t)\leq\ldots\leq Q_{\sigma_{t}(N)}(t)$; inequality (d) is
follows from the fact that
$\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\leq\sqrt{N}(Q_{\sigma_{t}(N)}(t)-Q_{\sigma_{t}(1)}(t))$.
∎

Now we are ready to present the proof of Theorem 1

Proof of Theorem 1: The proof is a direct application of the sufficient
conditions for throughput and heavy-traffic delay optimality, i.e., we need
only to show Eq. (7) and Eq. (8) hold.

Fix a load balancing policy $p$ in $\Pi$. Let us first consider the left-hand-
side (LHS) of Eq. (7) with $T_{1}=T$,

$\begin{split}LHS&\mathop{=}\limits^{(a)}{\sum_{t=t_{0}}^{t_{0}+T-1}\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t_{0})=Z\right]}\\\
&\mathop{=}\limits^{(b)}{\sum_{t=t_{0}}^{t_{0}+T-1}\mathbb{E}\left[\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]|Z(t_{0})=Z\right]}\\\
&\mathop{\leq}\limits^{(c)}{\sum_{t=t_{0}}^{t_{0}+T-1}\mathbb{E}\left[-\epsilon\frac{\mu_{min}}{\mu_{\Sigma}}\left\lVert{\mathbf{Q}}(t)\right\rVert\mid
Z(t_{0})=Z\right]}\\\
&\leq-\epsilon\frac{\mu_{min}}{{}\mu_{\Sigma}}\left\lVert{\mathbf{Q}}(t_{0})\right\rVert\end{split}$

where equation (a) comes from the linearity of condition expectation; equation
(b) follows from the tower property of conditional expectation and the fact
that ${\mathbf{Q}}(t)$, ${\mathbf{A}}(t)$ and ${\mathbf{S}}(t)$ are
conditionally independent of $Z(t_{0})$ when given $Z(t)$. inequality (c)
follows from Lemma 3 since the policy $p$ adopts a tilted distribution within
every time slot for all $Z(t)$. Hence, the condition of Lemma 1 is satisfied
and hence policy $p$ is throughput optimal.

Let us now turn to consider the left-hand-side (LHS) of Eq. (8) with $T_{2}=T$
and $\epsilon<\epsilon_{0}=\frac{\delta\mu_{\Sigma}}{2TN+2\delta}$.

$\begin{split}LHS\mathop{=}\limits^{(a)}&{\sum_{t=t_{0}}^{t_{0}+T-1}\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t_{0})=Z\right]}\\\
\mathop{=}\limits^{(b)}&{\sum_{t=t_{0}}^{t_{0}+T-1}\mathbb{E}\left[\mathbb{E}\left[\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)\right]\mid Z(t_{0})=Z\right]}\\\ \mathop{\leq}\limits^{(c)}&\sum_{t\neq
t^{*}}\mathbb{E}\left[\epsilon\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t)\right\rVert\mid
Z(t_{0})=Z\right]\\\
&+\mathbb{E}\left[\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t^{*})\right\rVert(\epsilon-\frac{\delta\lambda_{\Sigma}}{N})\mid
Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(d)}&(T-1)\epsilon\sqrt{N}(\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert+M)\\\
&+(\epsilon-\frac{\delta\lambda_{\Sigma}}{N})\sqrt{N}(\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert-M)\\\
=&(T\epsilon-\frac{\delta\lambda_{\Sigma}}{N})\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert+\sqrt{N}M(\frac{\delta\lambda_{\Sigma}}{N}+(T-2)\epsilon)\\\
\mathop{\leq}\limits^{(e)}&(T\epsilon-\frac{\delta\lambda_{\Sigma}}{N})\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert+K_{2}\\\
\mathop{\leq}\limits^{(f)}&-\frac{\delta\mu_{\Sigma}}{2N}\sqrt{N}\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert+K_{2}\end{split}$

where equation (a) comes from the linearity of condition expectation; equation
(b) follows from the tower property of conditional expectation and the fact
that ${\mathbf{Q}}_{\perp}(t)$, ${\mathbf{A}}(t)$ and ${\mathbf{S}}(t)$ are
conditionally independent of $Z(t_{0})$ when given $Z(t)$; inequality (c)
follows from Lemmas 3 and 4 since under policy $p$ there exists at least one
time slot $t^{*}$ within which at least a $\delta$-tilted distribution (or one
of its equivalent distribution in inner-product is $\delta$-tilted) is adopted
and all the other time slots a tilted distribution is used; inequality (d)
follow from the fact
$|\left\lVert{\mathbf{Q}}_{\perp}(t_{0}+T)\right\rVert-\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert|\leq
M=2T\sqrt{N}\max\\{A_{max},S_{max}\\}$ which is shown in Eq. (B) and the fact
$\epsilon-\frac{\delta\lambda_{\Sigma}}{N}<0$ for all $\epsilon<\epsilon_{0}$;
inequality (e) comes from the fact that
$\sqrt{N}M(\frac{\delta\lambda_{\Sigma}}{N^{2}}+(T-2)\epsilon)\leq
K_{2}=\sqrt{N}M(\frac{\delta\mu_{\Sigma}}{N}+T\mu_{\Sigma})$, which is
independent of $\epsilon$; inequality (f) holds since $\epsilon<\epsilon_{0}$
and $\lambda_{\Sigma}=\mu_{\Sigma}-\epsilon$. Therefore, since both
$-\frac{\delta\mu_{\Sigma}}{2N}\sqrt{N}$ and $K_{2}$ are independent of
$\epsilon$, the condition of Lemma 2 is satisfied, hence policy $p$ is heavy-
traffic delay optimal.∎

###  5.2 Proof of Proposition 3

Let us first look at assertion 1, i.e., JBT-$d$ is in $\Pi$ under homogeneous
servers. Based on Eq. (15), we can conclude that for any $t\geq 0$, the
dispatching distribution is a tilted distribution for all $Z(t)$. We are left
to show that at time slot $rT+1$, $r\in\\{0,1,2,\ldots\\}$, the dispatching
distribution is at least a $\delta$-distribution for some positive $\delta$.
This is equivalent to finding the maximum value for $P_{N}(rT+1)$ and the
minimum value of $P_{1}(rT+1)$ for all queue length state. In fact, they are
achieved at the same time when $\tilde{p}_{N}$ is in its largest value based
on Eq. (14), which is repeated as follows.

$P_{n}(t)=\sum_{i=n}^{N}\tilde{p}_{i}(t)\frac{1}{i}$

Then, there are two cases to consider.

(a) Note that $\tilde{p}_{N}(rT+1)=1$ if and only if all the servers have the
same queue length at the end of time slots $rT$ (sampling slots), which are
also the queue length state at the beginning of $rT+1$, i.e.,
${\mathbf{Q}}(rT+1)$. In this case, it can be easily seen that
$P_{n}(rT+1)=\frac{1}{N}$ for all $n$, which is not a $\delta$-tilted
distribution. However, it is the an equivalent distribution in inner-product
to $\hat{P}_{1}(rT+1)=1$ and $\hat{P}_{n}=0$ for $2\leq n\leq N$ as all the
queue lengths are equal, which is indeed a $\delta$-distribution.

(b) If the queue lengths are not all equal at the end of time slots $rT$, then
the maximum value for $\tilde{p}_{N}(rT+1)$ is strictly less than $1$ and it
is obtained when the queue length in the state that there are $N-1$ servers
that have the same queue length, which is strictly larger than the remaining
one. In this case, by sampling $d$ servers uniformly at random at the end of
times slots $rT$, the probability for the event that there are $N$ IDs in
memory, i.e., $\tilde{p}_{N}(rT+1)$ is given by

$\tilde{p}_{N}(rT+1)=1-\tilde{p}_{1}(rT+1)=1-\frac{d}{N}$

Therefore, we have $P_{1}(rT+1)=\frac{d}{N}+\frac{N-d}{N^{2}}$ and
$P_{n}(rT+1)=\frac{N-d}{N^{2}}$, which is inner product equivalent to
$\hat{P}_{1}(rT+1)=\frac{d}{N}+\frac{N-d}{N^{2}}$,
$\hat{P}_{2}(rT+1)=(N-1)\frac{N-d}{N^{2}}$ and $\hat{P}_{n}(rT+1)=0$ for all
$3\leq n\leq N$ as the $N-1$ queues have the same queue length. As a result,
we have for this state $Z(rT+1)$

$\hat{\Delta}_{1}(rT+1)=\frac{d}{N}(1-\frac{1}{N})\text{ and
}\hat{\Delta}_{N}(rT+1)=-\frac{1}{N}$

Thus, it is a $\delta$-distribution with
$\delta=\min\\{\frac{d}{N}(1-\frac{1}{N}),\frac{1}{N}\\}$, which is the lower
bound for $\delta$. That is, for any state $Z(rT+1)$, the dispatching
distribution is at least a $\delta$-distribution. Therefore, every $T+1$ time
slots, there exists one time slot in which the dispatching distribution (or an
inner product equivalent distribution) is at least a $\delta$-tilted
distribution with $\delta=\min\\{\frac{d}{N}(1-\frac{1}{N}),\frac{1}{N}\\}$.

The proof for heterogeneous servers follows exact the same idea with
additional care on the service rate. The probability for the server
$\sigma_{t}(n)$ to be selected at time $t$, i.e., $P_{n}(t)$ is given by

$P_{n}(t)=\mu_{\sigma_{t}(n)}\sum_{i=n}^{N}\frac{\tilde{p}_{i}(t)}{\sum_{\\{j\in
m(t),|m(t)|=i\\}}\mu_{j}}$

From it we can easily see that if
$\Delta_{n}(t)=P_{n}(t)-\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}}$ is positive,
then we must have that $\Delta_{n-1}(t)$ is also positive as it has one more
term in the equation above. Therefore, we can find a $k$ between $2$ and $N$
such that $\Delta_{n}(t)=P_{n}(t)-\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}}\geq
0$ for all $n<k$ and $\Delta_{n}(t)\leq 0$ for all $n\geq k$. Therefore,
condition (i) of $\Pi$ is satisfied.

For the condition (ii), we need to find the maximum value of
$\tilde{p}_{N}(rT+1)$ to bound $\delta$. There are also two cases as before.

(a) If $\tilde{p}_{N}(rT+1)=1$, the we must have that the queue lengths are
all equal at the end of time slots $rT$, which is the same as the beginning of
time slot $rT+1$. In this case,
$P_{n}(rT+1)=\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}}$ for all $n$. Note that
this dispatching distribution an equivalent distribution in inner-product to
$\hat{P}_{1}(rT+1)=1$ and $\hat{P}_{n}=0$ for $2\leq n\leq N$ as all the queue
lengths are equal, which is a $\delta$-distribution.

(b) If $\tilde{p}_{N}(rT+1)\neq 1$, the maximum value of $\tilde{p}_{N}(rT+1)$
is obtained when there are $N-1$ servers that have the same queue length,
which is strictly larger than the remaining one. In this case, we have
$\tilde{p}_{N}(rT+1)=1-\tilde{p}_{1}(rT+1)=1-\frac{d}{N}$ as before. Thus, we
can obtain

$P_{1}(rT+1)=\frac{d}{N}+(1-\frac{d}{N})\frac{\mu_{\sigma_{t}(1)}}{\mu_{\Sigma}}$

and $P_{n}(rT+1)=(1-\frac{d}{N})\frac{\mu_{\sigma_{t}(n)}}{\mu_{\Sigma}}$,
which is inner product equivalent to $\hat{P}_{1}(rT+1)=P_{1}(rT+1)$,
$\hat{P}_{2}(rT+1)=\sum_{n=2}^{N}P_{n}(rT+1)$ and $\hat{P}_{n}(rT+1)=0$ for
all $3\leq n\leq N$ since the last $N-1$ servers have the same queue lengths.
As a result, we have for this $Z(rT+1)$

$\hat{\Delta}_{1}(rT+1)=\frac{d}{N}(1-\frac{\mu_{\sigma_{t}(1)}}{\mu_{\Sigma}})\text{
and }\hat{\Delta}_{N}(rT+1)=-\frac{\mu_{\sigma_{t}(N)}}{\mu_{\Sigma}}$

Thus, it is a $\delta$-distribution with
$\delta=\min\\{\frac{d}{N}(1-\frac{\mu_{max}}{\mu_{\Sigma}}),\frac{\mu_{min}}{\mu_{\Sigma}}\\}$,
in which $\mu_{max}=\max_{n\in\mathcal{N}}\mu_{n}$ and
$\mu_{min}=\min_{n\in\mathcal{N}}\mu_{n}$, which is the lower bound of
$\delta$. Hence, for any $Z(rT+1)$, the dispatching probability distribution
(or its inner product equivalent one) is at least a $\delta$-distribution. ∎

##  6 Conclusion

We introduce a class $\Pi$ of flexible load balancing policies, which are
shown to be throughput and heavy-traffic delay optimal. This class includes as
special cases JSQ, power-of-d, and also allows flexibility in designing other
new policies. The JIQ policy, albeit exhibiting a good performance when the
traffic load is not heavy, is not in $\Pi$ since it is not heavy-traffic delay
optimal even for homogeneous servers. A new policy called JBT-$d$ is proposed
in the class $\Pi$, which enjoys the simplicity of JIQ while guaranteeing
heavy-traffic delay optimal. A unified analytic framework is established to
characterize this class of policies by exploring their common characteristics
and provide sufficient conditions that guarantee the heavy-traffic delay
optimality. Extensive simulations are used to demonstrate the good performance
and low complexity of the proposed policy compared to other existing ones.

## References

  * [1] Jeffrey Dean and Sanjay Ghemawat.  Mapreduce: simplified data processing on large clusters.  Communications of the ACM, 51(1):107–113, 2008.
  * [2] Atilla Eryilmaz and R Srikant.  Asymptotically tight steady-state queue length bounds implied by drift conditions.  Queueing Systems, 72(3-4):311–359, 2012.
  * [3] G Foschini and JACK Salz.  A basic dynamic routing problem and diffusion.  IEEE Transactions on Communications, 26(3):320–327, 1978.
  * [4] Ian Foster, Yong Zhao, Ioan Raicu, and Shiyong Lu.  Cloud computing and grid computing 360-degree compared.  In 2008 Grid Computing Environments Workshop, pages 1–10. Ieee, 2008.
  * [5] Lars George.  HBase: the definitive guide.  ” O’Reilly Media, Inc.”, 2011.
  * [6] Varun Gupta, Mor Harchol Balter, Karl Sigman, and Ward Whitt.  Analysis of join-the-shortest-queue routing for web server farms.  Performance Evaluation, 64(9):1062–1081, 2007.
  * [7] Bruce Hajek.  Hitting-time and occupation-time bounds implied by drift analysis with applications.  Advances in Applied probability, pages 502–525, 1982.
  * [8] Yi Lu, Qiaomin Xie, Gabriel Kliot, Alan Geller, James R Larus, and Albert Greenberg.  Join-idle-queue: A novel load balancing algorithm for dynamically scalable web services.  Performance Evaluation, 68(11):1056–1071, 2011.
  * [9] Siva Theja Maguluri, R Srikant, and Lei Ying.  Heavy traffic optimal resource allocation algorithms for cloud computing clusters.  Performance Evaluation, 81:20–39, 2014.
  * [10] Michael Mitzenmacher.  The power of two choices in randomized load balancing.  IEEE Transactions on Parallel and Distributed Systems, 12(10):1094–1104, 2001.
  * [11] Michael Mitzenmacher, Balaji Prabhakar, and Devavrat Shah.  Load balancing with memory.  In Foundations of Computer Science, 2002. Proceedings. The 43rd Annual IEEE Symposium on, pages 799–808. IEEE, 2002.
  * [12] Debankur Mukherjee, Sem C Borst, Johan SH Van Leeuwaarden, Philip A Whiting, et al.  Universality of load balancing schemes on the diffusion scale.  Journal of Applied Probability, 53(4):1111–1124, 2016.
  * [13] Rajesh Nishtala, Hans Fugal, Steven Grimm, Marc Kwiatkowski, Herman Lee, Harry C Li, Ryan McElroy, Mike Paleczny, Daniel Peek, Paul Saab, et al.  Scaling memcache at facebook.  In Presented as part of the 10th USENIX Symposium on Networked Systems Design and Implementation (NSDI 13), pages 385–398, 2013.
  * [14] Devavrat Shah and Balaji Prabhakar.  The use of memory in randomized load balancing.  In Information Theory, 2002. Proceedings. 2002 IEEE International Symposium on, page 125. IEEE, 2002.
  * [15] Alexander L Stolyar.  Pull-based load distribution in large-scale heterogeneous service systems.  Queueing Systems, 80(4):341–361, 2015.
  * [16] John N Tsitsiklis and Kuang Xu.  Queueing system topologies with limited flexibility.  In ACM SIGMETRICS Performance Evaluation Review, volume 41, pages 167–178. ACM, 2013.
  * [17] Nikita Dmitrievna Vvedenskaya, Roland L’vovich Dobrushin, and Fridrikh Izrailevich Karpelevich.  Queueing system with selection of the shortest of two queues: An asymptotic approach.  Problemy Peredachi Informatsii, 32(1):20–34, 1996.
  * [18] Weina Wang, Kai Zhu, Lei Ying, Jian Tan, and Li Zhang.  Maptask scheduling in mapreduce with data locality: Throughput and heavy-traffic optimality.  IEEE/ACM Transactions on Networking, 24(1):190–203, 2016.
  * [19] Richard R Weber.  On the optimal assignment of customers to parallel servers.  Journal of Applied Probability, pages 406–413, 1978.
  * [20] Lei Ying, R Srikant, and Xiaohan Kang.  The power of slightly more than one sample in randomized load balancing.  In 2015 IEEE Conference on Computer Communications (INFOCOM), pages 1131–1139. IEEE, 2015.

##  Appendix A Proof of Lemma 1

Before we present the proof of Lemma 1, we first introduce two lemmas which
will be the key ingredients in the proof. The first lemma enables us to bound
the moments of a stationary distribution based on drift condition, which can
be simplified by the second lemma.

The following lemma is introduced in [18], which is an extension of Lemma 1 in
[2] and can be proved from the results in [7].

######  Lemma 5.

For an irreducible aperiodic and positive recurrent Markov chain
$\\{X(t),t\geq 0\\}$ over a countable state space $\mathcal{X}$, which
converges in distribution to $\overline{X}$, and suppose
$V:\mathcal{X}\rightarrow\mathbb{R}_{+}$ is a Lyapunov function. We define the
$T$ time slot drift of $V$ at $X$ as

$\Delta V(X):=[V(X(t_{0}+T))-V(X(t_{0}))]\mathcal{I}(X(t_{0})=X),$

where $\mathcal{I}(.)$ is the indicator function. Suppose for some positive
finite integer $T$, the $T$ time slot drift of $V$ satisfies the following
conditions:

  * •

(C1) There exists an $\gamma>0$ and a $\kappa<\infty$ such that for any
$t_{0}=1,2,\ldots$ and for all $X\in\mathcal{X}$ with $V(X)\geq\kappa$,

$\mathbb{E}\left[\Delta V(X)\mid X(t_{0})=X\right]\leq-\gamma.$

  * •

(C2) There exists a constant $D<\infty$ such that for all $X\in\mathcal{X}$,

$\mathbb{P}(|\Delta V(X)|\leq D)=1.$

Then $\\{V(X(t)),t\geq 0\\}$ converges in distribution to a random variable
$\overline{V}$, and there exists constants $\theta^{*}>0$ and $C^{*}<\infty$
such that $\mathbb{E}\left[e^{\theta^{*}\overline{V}}\right]\leq C^{*}$, which
directly implies that all moments of random variable $\overline{V}$ exist and
are finite. More specifically, there exist finite constants
$\\{M_{r},r\in\mathbb{N}\\}$ such that for each positive $r$,
$\mathbb{E}\left[V(\overline{X})^{r}\right]\leq M_{r}$, where $M_{r}$ are
fully determined by $\kappa$, $\gamma$ and $D$.

######  Lemma 6.

For any $t\geq 0$, we have

$\left\lVert{\mathbf{Q}}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert^{2}\leq
2\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle+K$ (20)

where $K$ is a finite constant.

###### Proof.

Consider the left-hand-side (LHS) of Eq. (20).

$\begin{split}LHS&=\left\lVert{\mathbf{Q}}(t)+{\mathbf{A}}(t)-{\mathbf{S}}(t)+{\mathbf{U}}(t)\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert^{2}\\\
&\mathop{\leq}\limits^{(a)}\left\lVert{\mathbf{Q}}(t)+{\mathbf{A}}(t)-{\mathbf{S}}(t)\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert^{2}\\\
&=2\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle+\left\lVert{\mathbf{A}}(t)-{\mathbf{S}}(t)\right\rVert^{2}\\\
&\mathop{\leq}\limits^{(b)}2\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle+K\end{split}$

where inequality (a) holds as $[\max(a,0)]^{2}\leq a^{2}$ for any
$a\in\mathbb{R}$; in inequality (b), $K\triangleq
N\max(A_{\max},S_{\max})^{2}$ holds due to the assumptions that
$A_{\Sigma}(t)\leq A_{\max}$ and $S_{n}(t)\leq S_{\max}$ for all $t\geq 0$ and
all $n\in\mathcal{N}$, and independent of the queue length. ∎

We are now ready to prove Lemma 1.

Proof of Lemma 1: The proof follows from the application of Lemma 5 to the
Markov chain $\\{Z^{(\epsilon)}(t),t\geq 0\\}$ with Lyapunov function
$V(Z^{(\epsilon)}):=\left\lVert{\mathbf{Q}}^{(\epsilon)}\right\rVert$ and
$T=T_{1}$ since $m^{(\epsilon)}(t)$ is always finite. In particular, this
proof is completed in two steps, where the superscript (ϵ) will be omitted for
ease of notations.

(i) First, in order to apply Lemma 5, we need to show that the Markov chain
$\\{Z(t),t\geq 0\\}$ is irreducible, aperiodic and positive recurrent under
the hypothesis of Lemma 1. It can be easily seen that $\\{Z(t),t\geq 0\\}$ is
irreducible and aperiodic. Thus, we are left with the task to prove that the
Markov chain is positive recurrent. By the extension of Foster-Lyapunov
theorem, it suffices to find a Lyapunov function and a positive constant $T$
such that the expected $T$ time slot Lyapunov drift is bounded within a finite
subset of the state space and negative outside this subset.

Consider the Lyapunov function
$W(Z):=\left\lVert{\mathbf{Q}}\right\rVert^{2}$, and the corresponding
expected $T_{1}$ time slot mean conditional Lyapunov drift under the
hypothesis of Lemma 1.

$\begin{split}&\mathbb{E}\left[W(Z(t_{0}+T_{1}))-W(Z(t_{0}))\mid
Z(t_{0})\right]\\\
=&\mathbb{E}\left[\left\lVert{\mathbf{Q}}(t_{0}+T_{1})\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t_{0})\right\rVert^{2}\mid
Z(t_{0})\right]\\\
=&\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{1}-1}(\left\lVert{\mathbf{Q}}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert^{2})\mid
Z(t_{0})\right]\\\
\mathop{\leq}\limits^{(a)}&2\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{1}-1}\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle+K\mid
Z(t_{0})\right]\\\
\mathop{\leq}\limits^{(b)}&-2\gamma\left\lVert{\mathbf{Q}}(t_{0})\right\rVert+2K_{1}+2KT_{1}\end{split}$
(21)

where inequality (a) follows from Lemma 6, and inequality (b) results directly
from the hypothesis in Eq. (7). Pick any $\beta>0$ and let
$\mathcal{B}=\\{Z\in\mathcal{S}:\left\lVert{\mathbf{Q}}\right\rVert\leq\frac{K_{1}+KT_{1}+\beta}{\gamma}\\}$.
Then $\mathcal{B}$ is a finite subset of $\mathcal{S}$ as $m(t)$ is finite.
Moreover, for any $Z\in\mathcal{B}$, the conditional mean drift is less or
equal to $2K_{1}+2KT_{1}$, and for any $Z\in\mathcal{B}^{c}$, it is less than
or equal to $-\beta$. This finishes the proof of positive recurrence for any
$\epsilon>0$, and hence throughput optimal.

(ii) Second, in order to show that the hypothesis in Lemma 1 also ensures the
bounded moments for the stationary distribution, we will resort to Lemma 5.
Thus, we need to check Conditions (C1) and (C2), respectively.

For Condition (C1), we have

$\begin{split}&\mathbb{E}\left[\Delta V({Z})\mid Z(t_{0})=Z\right]\\\
=&\mathbb{E}\left[\left\lVert{\mathbf{Q}}(t_{0}+T_{1})\right\rVert-\left\lVert{\mathbf{Q}}(t_{0})\right\rVert\mid
Z(t_{0})=Z\right]\\\
=&\mathbb{E}\left[\sqrt{\left\lVert{\mathbf{Q}}(t_{0}+T_{1})\right\rVert^{2}}-\sqrt{\left\lVert{\mathbf{Q}}(t_{0})\right\rVert^{2}}\mid
Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(a)}&\frac{1}{2\left\lVert{\mathbf{Q}}(t_{0})\right\rVert}\mathbb{E}\left[\left\lVert{\mathbf{Q}}(t_{0}+T_{1})\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t_{0})\right\rVert^{2}\mid
Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(b)}&-\gamma+\frac{K_{1}+KT_{1}}{\left\lVert{\mathbf{Q}}(t_{0})\right\rVert}\end{split}$

where inequality (a) follows from the fact that $f(x)=\sqrt{x}$ is concave;
(b) comes from the upper bound in Eq. (21). Hence, (C1) in Lemma 5 is
verified.

For Condition (C2), we have

$\begin{split}|\Delta
V(Z)|&=|\left\lVert{\mathbf{Q}}(t_{0}+T_{1})\right\rVert-\left\lVert{\mathbf{Q}}(t_{0})\right\rVert|\mathcal{I}(Z(t_{0})=Z)\\\
&\mathop{\leq}\limits^{(a)}\left\lVert{\mathbf{Q}}(t_{0}+T_{1})-{\mathbf{Q}}(t_{0})\right\rVert\mathcal{I}(Z(t_{0})=Z)\\\
&\mathop{\leq}\limits^{(b)}T_{1}\sqrt{N}\max(A_{\max},S_{\max})\end{split}$

where inequality (a) follows from the fact that $|\left\lVert{\bf
x}\right\rVert-\left\lVert{\bf y}\right\rVert|\leq\left\lVert{\bf x}-{\bf
y}\right\rVert$ holds for any ${\bf x}$, ${\bf y}\in\mathbb{R}^{N}$;
inequality (b) holds due to the assumptions that the $A_{\Sigma}(t)\leq
A_{\max}$ and $S_{n}(t)\leq S_{\max}$ for all $t\geq 0$ and all
$n\in\mathcal{N}$, and independent of the queue length. This verifies
Condition (C2) and hence complete the proof of Lemma 1.∎

##  Appendix B Proof of Lemma 2

We now proceed to prove Lemma 2. Before we present the proof, the following
lemmas which serve as useful preliminary steps are first introduced. Denote by
${\mathbf{Q}}_{\parallel}$ and ${\mathbf{Q}}_{\perp}$ the parallel and
perpendicular components of the queue length vector ${\mathbf{Q}}$ with
respect to the line $\mathbf{c}=\frac{1}{\sqrt{N}}\mathbf{1}$, i.e.,

${\mathbf{Q}}_{\parallel}:=\langle\mathbf{c},{\mathbf{Q}}\rangle\mathbf{c}\quad\quad{\mathbf{Q}}_{\perp}:={\mathbf{Q}}-{\mathbf{Q}}_{\parallel}$
(22)

The following lemma is a natural extension of Lemma 7 in [2] to $T$ time
slots.

######  Lemma 7.

Define the following Lyapunov functions

$V_{\perp}(Z):=\left\lVert{\mathbf{Q}}_{\perp}\right\rVert,W(Z):=\left\lVert{\mathbf{Q}}\right\rVert^{2}\text{
and }W_{\parallel}(Z):=\left\lVert{\mathbf{Q}}_{\parallel}\right\rVert^{2}$

with the corresponding $T$ time-slot drift given by

$\begin{split}&\Delta
V_{\perp}(Z):=[V_{\perp}(Z(t_{0}+T))-V_{\perp}(Z(t_{0}))]\mathcal{I}(Z(t_{0})=Z)\\\
&\Delta W(Z):=[W(Z(t_{0}+T))-W(Z(t_{0}))]\mathcal{I}(Z(t_{0})=Z)\\\ &\Delta
W_{\parallel}(Z):=[W_{\parallel}(Z(t_{0}+T))-W_{\parallel}(Z(t_{0}))]\mathcal{I}(Z(t_{0})=Z)\\\
\end{split}$

Then, the drift of $V_{\perp}(.)$ can be bounded in terms of $W(.)$ and
$W_{\parallel}(.)$ as follows.

$\Delta
V_{\perp}(Z)\leq\frac{1}{2\left\lVert{\mathbf{Q}}_{\perp}\right\rVert}(\Delta
W(Z)-\Delta W_{\parallel}(Z))$

for all $Z\in\mathcal{S}$.

######  Lemma 8.

For any $t\geq 0$, we have

$\left\lVert{\mathbf{Q}}_{\parallel}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}_{\parallel}(t)\right\rVert^{2}\geq
2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle.$

###### Proof.

$\begin{split}&\left\lVert{\mathbf{Q}}_{\parallel}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}_{\parallel}(t)\right\rVert^{2}\\\
=&2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}_{\parallel}(t+1)-{\mathbf{Q}}_{\parallel}(t)\rangle+\left\lVert{\mathbf{Q}}_{\parallel}(t+1)-{\mathbf{Q}}_{\parallel}(t)\right\rVert^{2}\\\
\geq&2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}_{\parallel}(t+1)-{\mathbf{Q}}_{\parallel}(t)\rangle\\\
=&2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}(t+1)-{\mathbf{Q}}(t)\rangle-2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}_{\perp}(t+1)-{\mathbf{Q}}_{\perp}(t)\rangle\\\
\mathop{\geq}\limits^{(a)}&2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}(t+1)-{\mathbf{Q}}(t)\rangle\\\
\mathop{\geq}\limits^{(b)}&2\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\end{split}$

where the inequality (a) is true because
$\langle{\mathbf{Q}}_{\parallel}(t),{\mathbf{Q}}_{\perp}(t)\rangle=0$ and
$\langle{\mathbf{Q}}_{\perp}(t+1),{\mathbf{Q}}_{\parallel}(t)\rangle=0$; (b)
follows from the fact that all the components of ${\mathbf{Q}}_{\parallel}(t)$
and ${\mathbf{U}}(t)$ are nonnegative. ∎

We are now ready to prove the following result, which is often called state
space collapse and is the key ingredient for establishing heavy traffic delay
optimality. It shows that under the hypothesis of Lemma 2, the multi-dimension
space for the queue length vector will reduce to one dimension in the sense
that the deviation from the line $\mathbf{c}$ is bounded by a constant, which
is independent with the heavy-traffic parameter $\epsilon$.

######  Lemma 9.

If the hypothesis of Lemma 2 holds, then we have that ${\mathbf{Q}}_{\perp}$
is bounded in the sense that in steady state there exists finite constants
$\\{L_{r},r\in\mathbb{N}\\}$ such that

$\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}_{\perp}^{(\epsilon)}\right\rVert^{r}\right]\leq
L_{r}$

for all $\epsilon\in(0,\epsilon_{0})$ and $r\in\mathbb{N}$.

###### Proof.

It suffices to show that $V_{\perp}(Z)$ satisfies the Conditions (C1) and (C2)
in Lemma 5. Fix $\epsilon\in(0,\epsilon_{0})$, and the superscript will be
omitted for simplicity in the following arguments.

(i) For the Condition (C1), let
$\Lambda(t):=\left\lVert{\mathbf{Q}}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert^{2}$
and
$\Lambda_{\parallel}(t):=\left\lVert{\mathbf{Q}}_{\parallel}(t+1)\right\rVert^{2}-\left\lVert{\mathbf{Q}}_{\parallel}(t)\right\rVert^{2}$.
Then, we have

$\begin{split}&\mathbb{E}\left[\Delta V_{\perp}(Z)\mid Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(a)}&\frac{1}{2\left\lVert{\mathbf{Q}}_{\perp}\right\rVert}\mathbb{E}\left[\Delta
W(Z)-\Delta W_{\parallel}(Z)\mid Z(t_{0})=Z\right]\\\
=&\frac{1}{2\left\lVert{\mathbf{Q}}_{\perp}\right\rVert}\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{2}-1}\Lambda(t)-\Lambda_{\parallel}(t)\mid
Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(b)}&\frac{1}{2\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert}\mathbb{E}\left[\sum_{t=t_{0}}^{t_{0}+T_{2}-1}2\langle{\mathbf{Q}}_{\perp}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle+K\mid
Z(t_{0})=Z\right]\\\
\mathop{\leq}\limits^{(c)}&-\eta+\frac{2K_{2}+KT_{2}}{2\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert}\end{split}$

where the inequality (a) follows from Lemma 7; the inequality (b) holds as a
result of Lemmas 6 and 8; the inequality (c) follows directly from the
hypothesis Eq. (8). Hence, the Condition (C1) is verified.

(ii) For the Condition (C2), we have

$\displaystyle|\Delta V_{\perp}(Z)|$ $\displaystyle=$
$\displaystyle|\left\lVert{\mathbf{Q}}_{\perp}(t_{0}+T_{2})\right\rVert-\left\lVert{\mathbf{Q}}_{\perp}(t_{0})\right\rVert|\mathcal{I}(Z(t_{0})=Z)$
$\displaystyle\mathop{\leq}\limits^{(a)}$
$\displaystyle\left\lVert{\mathbf{Q}}_{\perp}(t_{0}+T_{2})-{\mathbf{Q}}_{\perp}(t_{0})\right\rVert\mathcal{I}(Z(t_{0})=Z)$
$\displaystyle=$
$\displaystyle\left\lVert{\mathbf{Q}}(t_{0}+T_{2})-{\mathbf{Q}}_{\parallel}(t_{0}+T_{2})-{\mathbf{Q}}(t_{0})+{\mathbf{Q}}_{\parallel}(t_{0})\right\rVert\mathcal{I}(Z(t_{0})=Z)$
$\displaystyle\mathop{\leq}\limits^{(b)}$
$\displaystyle\left\lVert{\mathbf{Q}}(t_{0}+T_{2})-{\mathbf{Q}}(t_{0})\right\rVert+\left\lVert{\mathbf{Q}}_{\parallel}(t_{0}+T_{2})-{\mathbf{Q}}_{\parallel}(t_{0})\right\rVert\mathcal{I}(Z(t_{0})=Z)$
$\displaystyle\mathop{\leq}\limits^{(c)}$ $\displaystyle
2\left\lVert{\mathbf{Q}}(t_{0}+T_{2})-{\mathbf{Q}}(t_{0})\right\rVert\mathcal{I}(Z(t_{0})=Z)$
$\displaystyle\mathop{\leq}\limits^{(d)}$ $\displaystyle
2T_{2}\sqrt{N}\max(A_{\max},S_{\max})$ (23)

where the inequality (a) follows from the fact that $|\left\lVert{\bf
x}\right\rVert-\left\lVert{\bf y}\right\rVert|\leq\left\lVert{\bf x}-{\bf
y}\right\rVert$ holds for any ${\bf x}$, ${\bf y}\in\mathbb{R}^{N}$;
inequality (b) follows from triangle inequality; (c) holds due to the non-
expansive property of projection to a convex set. (d) holds due to the
assumptions that the $A_{\Sigma}(t)\leq A_{\max}$ and $S_{n}(t)\leq S_{\max}$
for all $t\geq 0$ and all $n\in\mathcal{N}$, and independent of the queue
length. This verifies Condition (C2) and hence complete the proof of Lemma 9.
∎

The following result on the unused service is another key ingredient for
establishing heavy-traffic delay optimal.

######  Lemma 10.

For any $\epsilon>0$ and $t\geq 0$, we have

$Q_{n}^{(\epsilon)}(t+1)U_{n}^{(\epsilon)}(t)=0\text{ and
}q^{(\epsilon)}(t+1)u^{(\epsilon)}(t)=0.$

If the system has a finite first moment, then we have for some constants
$c_{1}$ and $c_{2}$

$\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert^{2}\right]\leq
c_{1}\epsilon\text{ and
}\mathbb{E}\left[({\overline{u}^{(\epsilon)}})^{2}\right]\leq c_{2}\epsilon$

###### Proof.

According to the queue dynamic in Eq. (1), we can see when $U_{n}(t)$ is
positive, $Q_{n}(t+1)$ must be zero, which gives the results
$Q_{n}^{(\epsilon)}(t+1)U_{n}^{(\epsilon)}(t)=0$ for all $n\in\mathcal{N}$ and
all $t\geq 0$, and the corresponding result for the resource-pooled system
$q^{(\epsilon)}(t+1)u^{(\epsilon)}(t)=0$.

Then, let us consider the Lyapunov function
$W_{1}(Z(t))=\left\lVert{\mathbf{Q}}(t)\right\rVert_{1}$. In the steady state
with a finite first moment, the mean drift of $W_{1}(Z(t))$ is zero. Then, we
have

$\begin{split}0=\mathbb{E}\left[\left\lVert{\mathbf{A}}^{(\epsilon)}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}+\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert_{1}\right]\end{split}$

which directly implies

$\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert_{1}\right]=\epsilon$
(24)

Moreover, due to the fact that $U_{n}(t)\leq S_{\max}$ for all
$n\in\mathcal{N}$ and $t\geq 0$, we have
$\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert^{2}\leq
S_{\max}\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert_{1}$.
Therefore, we can conclude that
$\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert^{2}\right]\leq
S_{\max}\epsilon$ and
$\mathbb{E}\left[({\overline{u}^{(\epsilon)})}^{2}\right]\leq
NS_{\max}\epsilon$. ∎

Now, we are well prepared to prove Lemma 2

Proof of Lemma 2: First, let us consider the Lyapunov function
$V_{1}(Z):=\left\lVert{\mathbf{Q}}\right\rVert_{1}^{2}$ and the corresponding
conditional mean drift, defined as
$D_{1}(Z(t)):=\mathbb{E}\left[V_{1}(Z(t+1))-V_{1}(Z(t))\mid Z(t)=Z\right].$

Then, we have the following equation, in which the time reference $(t)$ will
be omitted after the second step for brevity and
${\mathbf{Q}}^{+}:={\mathbf{Q}}(t+1)$.

$\begin{split}&D_{1}(Z(t))\\\
=&\mathbb{E}\left[\left\lVert{\mathbf{Q}}(t+1)\right\rVert_{1}^{2}-\left\lVert{\mathbf{Q}}(t)\right\rVert_{1}^{2}\mid
Z(t)=Z\right]\\\
=&\mathbb{E}\left[\left(\left\lVert{\mathbf{Q}}(t)\right\rVert_{1}+\left\lVert{\mathbf{A}}(t)\right\rVert_{1}-\left\lVert{\mathbf{S}}(t)\right\rVert_{1}+\left\lVert{\mathbf{U}}(t)\right\rVert_{1}\right)^{2}\mid
Z(t)=Z\right]\\\
&-\mathbb{E}\left[\left\lVert{\mathbf{Q}}(t)\right\rVert_{1}^{2}\mid
Z(t)=Z\right]\\\
=&\mathbb{E}\left[2\left\lVert{\mathbf{Q}}{}\right\rVert_{1}\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)+\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}\right.\\\
&+\left.{}2\left(\left\lVert{\mathbf{Q}}\right\rVert_{1}+\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)\left\lVert{\mathbf{U}}\right\rVert_{1}+\left\lVert{\mathbf{U}}\right\rVert_{1}^{2}\vphantom{\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}}\mid
Z\right]\\\
=&\mathbb{E}\left[2\left\lVert{\mathbf{Q}}\right\rVert_{1}\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)+\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}\right.\\\
&+\left.{}2\left\lVert{\mathbf{Q}}^{+}\right\rVert_{1}\left\lVert{\mathbf{U}}\right\rVert_{1}-\left\lVert{\mathbf{U}}\right\rVert_{1}^{2}\vphantom{\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}}\mid
Z\right]\\\
\leq&\mathbb{E}\left[2\left\lVert{\mathbf{Q}}\right\rVert_{1}\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)+\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}\right.\\\
&+\left.{}2\left\lVert{\mathbf{Q}}^{+}\right\rVert_{1}\left\lVert{\mathbf{U}}\right\rVert_{1}\vphantom{\left(\left\lVert{\mathbf{A}}\right\rVert_{1}-\left\lVert{\mathbf{S}}\right\rVert_{1}\right)^{2}}\mid
Z\right]\end{split}$ (25)

Under the hypothesis of Lemma 2, there exists a steady-state distribution with
finite moments for any $\epsilon>0$. Therefore, the mean drift in steady-state
is zero, i.e., $\mathbb{E}\left[{D_{1}}(\overline{Z}^{(\epsilon)})\right]=0$.
Therefore, taking the expectation of both sides of Eq. (25) with respect to
the steady-state distribution $\overline{Z}^{(\epsilon)}$, yields

$\begin{split}\epsilon\mathbb{E}\left[\sum_{n=1}^{N}\overline{Q}_{n}^{(\epsilon)}\right]\leq\frac{\zeta^{(\epsilon)}}{2}+\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}^{(\epsilon)}(t+1)\right\rVert_{1}\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}(t)\right\rVert_{1}\right]\end{split}$

where
$\zeta^{(\epsilon)}=(\sigma_{\Sigma}^{(\epsilon)})^{2}+\nu_{\Sigma}^{2}+\epsilon^{2}$.
For the resource-pooled system, by letting $N=1$ in Eq. (25) and taking the
expectation with respect to $\overline{q}^{(\epsilon)}$, we have

$\begin{split}\epsilon\mathbb{E}\left[\overline{q}^{(\epsilon)}\right]&=\frac{\zeta^{(\epsilon)}}{2}+\mathbb{E}\left[\overline{q}^{(\epsilon)}(t+1)u^{(\epsilon)}(t)\right]-\frac{1}{2}\mathbb{E}\left[(u^{(\epsilon)})^{2}\right].\\\
\end{split}$

Then, based on the property on the unused service in Lemma 10, we have

$\frac{\zeta^{(\epsilon)}}{2}-\frac{c_{2}}{2}\epsilon\leq\epsilon\mathbb{E}\left[\overline{q}^{(\epsilon)}\right]\leq\epsilon\mathbb{E}\left[\sum_{n=1}^{N}\overline{Q}_{n}^{(\epsilon)}\right]\leq\frac{\zeta^{(\epsilon)}}{2}+\overline{B}^{(\epsilon)}$
(26)

where
$\overline{B}^{(\epsilon)}:=\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}^{(\epsilon)}(t+1)\right\rVert_{1}\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}(t)\right\rVert_{1}\right]$.

Therefore, in order to show heavy-traffic delay optimality, all we need to
show is that $\lim_{\epsilon\downarrow 0}\overline{B}^{(\epsilon)}=0$. Note
that $\overline{B}^{(\epsilon)}$ can be bounded as follows.

$\begin{split}\overline{B}^{(\epsilon)}&\mathop{=}\limits^{(a)}N\mathbb{E}\left[\langle\overline{{\mathbf{U}}}^{(\epsilon)}(t),-\overline{{\mathbf{Q}}}^{(\epsilon)}_{\perp}(t+1)\rangle\right]\\\
&\mathop{\leq}\limits^{(b)}N\sqrt{\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}_{\perp}(t)\right\rVert^{2}\right]\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}^{(\epsilon)}_{\perp}(t+1)\right\rVert^{2}\right]}\\\
&\mathop{=}\limits^{(c)}N\sqrt{\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}_{\perp}(t)\right\rVert^{2}\right]\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}^{(\epsilon)}_{\perp}(t)\right\rVert^{2}\right]},\\\
\end{split}$

where the equality (a) comes from the property
$Q_{n}^{(\epsilon)}(t+1)U_{n}^{(\epsilon)}(t)=0$ for all $n\in\mathcal{N}$ and
all $t\geq 0$ in Lemma 10 and the definition of ${\mathbf{Q}}_{\perp}$; the
inequality (b) holds due to Cauchy-Schwartz inequality; the last equality (c)
is true since the distributions of
$\overline{{\mathbf{Q}}}^{(\epsilon)}_{\perp}(t+1)$ and
$\overline{{\mathbf{Q}}}^{(\epsilon)}_{\perp}(t)$ are the same in steady
state.

As shown in Lemma 9,
$\mathbb{E}\left[\left\lVert\overline{{\mathbf{Q}}}_{\perp}^{(\epsilon)}\right\rVert^{2}\right]\leq
L_{2}$ holds for all $\epsilon\in(0,\epsilon_{0})$ and some constant $L_{2}$
which is independent of $\epsilon$. Meanwhile, note that
$\mathbb{E}\left[\left\lVert\overline{{\mathbf{U}}}^{(\epsilon)}\right\rVert^{2}\right]\leq
c_{1}\epsilon$ for some $c_{1}$ independent of $\epsilon$ based on Lemma 10.
Then, we have for all $\epsilon\in(0,\epsilon_{0})$

$\overline{B}^{(\epsilon)}\leq N\sqrt{c_{1}\epsilon L_{2}}$ (27)

Therefore, it can be seen from Eq. (27) that $\lim_{\epsilon\downarrow
0}\overline{B}^{(\epsilon)}=0$, which directly implies
$\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\sum_{n}\overline{Q}_{n}^{(\epsilon)}\right]=\lim_{\epsilon\downarrow
0}\epsilon\mathbb{E}\left[\overline{q}^{(\epsilon)}\right]$, and thus the
proof of Lemma 2 is completed.∎

##  Appendix C Proof of theorem 2

Before we present the proof of Theorem 2, we first show that both the
following Lyapunov functions have finite expectations in steady state under
JIQ, i.e., throughput optimal in a strong sense.

$W(Z):=\left\lVert{\mathbf{Q}}\right\rVert^{2},V_{1}(Z):=\left\lVert{\mathbf{Q}}\right\rVert_{1}^{2}$

######  Lemma 11.

Consider a load balancing system with homogeneous servers under JIQ policy,
the steady state means $\mathbb{E}\left[W(\overline{Z}^{(\epsilon)})\right]$
and $\mathbb{E}\left[V_{1}(\overline{Z}^{(\epsilon)})\right]$ are both finite
for any $\epsilon>0$.

###### Proof.

This proof is a direct application of Lemma 1. Let us consider $T_{1}=1$ in
Lemma 1, we have

$\begin{split}&\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)-{\mathbf{S}}(t)\rangle\mid
Z(t)=Z\right]\\\
=&\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)\rangle\mid
Z(t)=Z\right]-\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{S}}(t)\rangle\mid
Z(t)=Z\right]\\\
\mathop{=}\limits^{(a)}&\mathbb{E}\left[\langle{\mathbf{Q}}(t),{\mathbf{A}}(t)\rangle\mid
Z(t)=Z\right]-\mu\left\lVert{\mathbf{Q}}\right\rVert_{1}\\\
\mathop{\leq}\limits^{(b)}&\left(\mu-\frac{\epsilon}{N}\right)\left\lVert{\mathbf{Q}}\right\rVert_{1}-\mu\left\lVert{\mathbf{Q}}\right\rVert_{1}\\\
\mathop{\leq}\limits^{(c)}&-\frac{\epsilon}{N}\left\lVert{\mathbf{Q}}\right\rVert\end{split}$

where step (a) follows from the fact that the servers are homogeneous and each
service process is independent with the system state; inequality (b) holds due
to the property of JIQ; inequality (c) holds since $l_{1}$ norm
$\left\lVert\mathbf{x}\right\rVert_{1}$ of any vector
$\mathbf{x}\in\mathbb{R}^{N}$ is no smaller than its $l_{2}$ norm
$\left\lVert\mathbf{x}\right\rVert$. Thus, according to Lemma 1, we have
$\mathbb{E}\left[W(\overline{Z}^{(\epsilon)})\right]$ is finite for any
$\epsilon>0$.

Note that $\left\lVert{\mathbf{Q}}\right\rVert_{1}^{2}\leq
N\left\lVert{\mathbf{Q}}\right\rVert^{2}$. Therefore, it follows that
$\mathbb{E}\left[V_{1}(\overline{Z}^{(\epsilon)})\right]$ is also finite. ∎

Now, we are well prepared to prove Theorem 2

Proof of Theorem 2: Recall Remark 1, to prove JIQ is not heavy-traffic delay
optimal, it suffices to find a certain combination of arrival and service
processes under which the heavy-traffic limit under JIQ cannot achieve the
heavy-traffic limit of resource-pooled system, while JSQ can. Mathematically,
it is sufficient to show that there exists a subsequence of $\epsilon$ such
that Eq. (3) (definition of heavy-traffic delay optimality) does not hold. In
particular, in this proof we consider the case of two homogeneous servers with
arrival process in $\mathcal{A}$ and constant service with rate one.

We will prove the result by contradiction. Suppose the result in Theorem 2
does not hold, i.e., JIQ is heavy-traffic delay optimal. Then it means that
$\lim_{\epsilon\downarrow 0}\overline{B}^{(\epsilon)}=0$, i.e.,
$\limsup_{\epsilon\downarrow 0}\overline{B}^{(\epsilon)}=0$ must hold in our
considered case. This is true since

$\frac{\zeta^{(\epsilon)}}{2}+\overline{B}^{(\epsilon)}-c_{1}\epsilon\leq\epsilon\mathbb{E}\left[\sum_{n=1}^{2}\overline{Q}_{n}^{(\epsilon)}\right]\leq\frac{\zeta^{(\epsilon)}}{2}+\overline{B}^{(\epsilon)},$

which comes from the fact that
$\mathbb{E}\left[V_{1}(\overline{Z}^{(\epsilon)})\right]$ is finite under JIQ
by Lemma 11 and the property of unused service in Lemma 10.

Note that ${\mathbf{Q}}^{+}:={\mathbf{Q}}(t+1)$, we have

$\begin{split}\overline{B}^{(\epsilon)}&\mathop{=}\limits^{(a)}\mathbb{E}\left[(\overline{Q}_{1}^{+})^{(\epsilon)}\overline{U}_{2}^{(\epsilon)}+(\overline{Q}_{2}^{+})^{(\epsilon)}\overline{U}_{1}^{(\epsilon)}\right]\\\
&\mathop{=}\limits^{(b)}2\mathbb{E}\left[(\overline{Q}_{1}^{+})^{(\epsilon)}\overline{U}_{2}^{(\epsilon)}\right]\\\
&\mathop{=}\limits^{(c)}2\sum_{k=1}^{\infty}k\mathbb{P}\left((\overline{Q}_{1}^{+})^{(\epsilon)}=k,\overline{U}_{2}^{(\epsilon)}>0\right),\\\
\end{split}$

where inequality (a) follows from the property shown in Lemma 10 that
$Q_{n}(t+1)U_{n}(t)=0$ for all $t\geq 0$; equality (b) holds due to the
symmetry property of JIQ policy; inequality (c) comes from the fact that the
unused service is one when it is positive since the service rate is constant
$1$. Since $\limsup_{\epsilon\downarrow 0}\overline{B}^{(\epsilon)}=0$, we
have

$\lim_{\epsilon\downarrow
0}\sum_{k=1}^{\infty}k\mathbb{P}\left((\overline{Q}_{1}^{+})^{(\epsilon)}=k,\overline{U}_{2}^{(\epsilon)}>0\right)=0.$
(28)

Now, let us consider another Lyapunov function
$V_{2}(Z(t))=(Q_{1}(t)-Q_{2}(t))^{2}$. By the fact that the service process is
constant, we have

$\begin{split}&\mathbb{E}\left[V_{2}(Z(t+1))-V_{2}(Z(t))\mid Z(t)=Z\right]\\\
=&2\mathbb{E}\left[(Q_{1}-Q_{2})(A_{1}(t)-A_{2}(t))\mid
Z\right]+\mathbb{E}\left[(A_{1}(t)-A_{2}(t))^{2}\mid Z\right]\\\
&-2\mathbb{E}\left[Q_{1}^{+}U_{2}(t)+Q_{2}^{+}U_{1}(t)\mid
Z\right]-\mathbb{E}\left[(U_{1}(t)-U_{2}(t))^{2}\mid Z\right]\\\
=&2\mathbb{E}\left[(Q_{1}-Q_{2})(A_{1}(t)-A_{2}(t))\mid
Z\right]+\mathbb{E}\left[(A_{1}(t)-A_{2}(t))^{2}\mid Z\right]\\\
&-4\mathbb{E}\left[Q_{1}^{+}U_{2}(t)\mid
Z\right]-\mathbb{E}\left[(U_{1}(t)-U_{2}(t))^{2}\mid Z\right].\end{split}$

Since
$\mathbb{E}\left[V_{2}(Z(t))\right]\leq\mathbb{E}\left[W(Z(t))\right]<\infty$
by Lemma 11, the mean drift of $V_{2}(.)$ in steady state is zero, which
implies the following equation in steady state for any $\epsilon>0$.

$\begin{split}&4\mathbb{E}\left[(\overline{Q}_{1}^{+})^{(\epsilon)}\overline{U}_{2}^{(\epsilon)}\right]\\\
=&2\mathbb{E}\left[\left(\overline{Q}_{1}^{(\epsilon)}-\overline{Q}_{2}^{(\epsilon)}\right)\left(\overline{A}_{1}^{(\epsilon)}-\overline{A}_{2}^{(\epsilon)}\right)\right]+\mathbb{E}\left[\left(\overline{A}_{1}^{(\epsilon)}-\overline{A}_{2}^{(\epsilon)}\right)^{2}\right]\\\
&-\mathbb{E}\left[\left(\overline{U}_{1}^{(\epsilon)}-\overline{U}_{2}^{(\epsilon)}\right)^{2}\right].\end{split}$
(29)

In the following, we will analyze each term on the right-hand side of Eq.
(29).

To start with, let us look at the first term on the right-hand side of Eq.
(29), denoted as $\mathcal{T}_{1}$. For ease of exposition, we shall omit the
superscript (ϵ) after the first step, and reintroduce it when necessary. In
the following term $\mathcal{T}_{1}$, the equality (a) follows from the
symmetry property of JIQ policy; (b) holds since when both queues are non-idle
random routing with equal probability is adopted in JIQ policy; (c) comes from
the fact that when one queue is idle, the arrival is always routed to the idle
queue under JIQ; equality (d) holds since $\overline{{\mathbf{Q}}}(t+1)$ has
the same distribution as $\overline{{\mathbf{Q}}}(t)$ in steady state.

$\begin{split}\mathcal{T}_{1}&=2\mathbb{E}\left[\left(\overline{Q}_{1}^{(\epsilon)}-\overline{Q}_{2}^{(\epsilon)}\right)\left(\overline{A}_{1}^{(\epsilon)}-\overline{A}_{2}^{(\epsilon)}\right)\right]\\\
&\mathop{=}\limits^{(a)}2\mathbb{E}\left[\left(\overline{Q}_{1}-\overline{Q}_{2}\right)\left(\overline{A}_{1}-\overline{A}_{2}\right)\mathbbm{1}\left(\overline{Q}_{1}>0,\overline{Q}_{2}>0\right)\right]\\\
&\quad+4\mathbb{E}\left[\left(\overline{Q}_{1}-\overline{Q}_{2}\right)\left(\overline{A}_{1}-\overline{A}_{2}\right)\mathbbm{1}\left(\overline{Q}_{1}\geq
0,\overline{Q}_{2}=0\right)\right]\\\
&\mathop{=}\limits^{(b)}4\mathbb{E}\left[\left(\overline{Q}_{1}-\overline{Q}_{2}\right)\left(\overline{A}_{1}-\overline{A}_{2}\right)\mathbbm{1}\left(\overline{Q}_{1}\geq
0,\overline{Q}_{2}=0\right)\right]\\\
&\mathop{=}\limits^{(c)}-4\lambda_{\Sigma}\sum_{k=1}^{\infty}k\mathbb{P}\left(\overline{Q}_{1}=k,\overline{Q}_{2}=0\right)\\\
&\mathop{=}\limits^{(d)}-4\lambda_{\Sigma}\sum_{k=1}^{\infty}k\mathbb{P}\left(\overline{Q}_{1}^{+}=k,\overline{Q}_{2}^{+}=0,\overline{U}_{2}=0\right)\\\
&\quad-4\lambda_{\Sigma}\sum_{k=1}^{\infty}k\mathbb{P}\left(\overline{Q}_{1}^{+}=k,\overline{Q}_{2}^{+}=0,\overline{U}_{2}>0\right).\\\
\end{split}$

Now, let us define three events as follows where $k\geq 1$.

$\begin{split}&E_{1,k}:=\\{\overline{Q}_{1}^{+}=k,\overline{Q}_{2}^{+}=0,\overline{U}_{2}=0\\},\\\
&E_{2,k}:=\\{\overline{Q}_{1}(t+2)=k,\overline{Q}_{2}(t+2)=0,\overline{U}_{2}^{+}>0\\},\\\
&E_{3,k}:=\\{\overline{Q}_{1}^{+}=k,\overline{Q}_{2}^{+}=0,\overline{U}_{2}>0\\}.\\\
\end{split}$

Then, it can be easily seen that for each occurrence of event $E_{1,k}$, due
to the fact that arrival process is in $\mathcal{A}$, there exists a positive
probability $p_{0}$, i.e., the probability of no arrival, such that
$E_{2,k-1}$ will happen. This is because when there is no arrival, the unused
service for $Q_{2}$ is $1$ and the queue length of $Q_{1}$ will decrease by
$1$ due to the constant service rate $1$. Since

$\mathbb{P}\left(E_{2,k-1}\right)\geq\mathbb{P}\left(E_{1,k}\right)\mathbb{P}\left(E_{2,k-1}\mid
E_{1,k}\right),$

which directly implies that $\mathbb{P}\left(E_{2,k-1}\right)\geq
p_{0}\mathbb{P}\left(E_{1,k}\right)$. Thus, we have

$\mathbb{P}\left(E_{3,k-1}\right)=\mathbb{P}\left(E_{2,k-1}\right)\geq
p_{0}\mathbb{P}\left(E_{1,k}\right).$ (30)

which holds since these events are defined on the steady-state distribution.

Now, let us first rewrite $\mathcal{T}_{1}$ with respect to events $E_{1,k}$
and $E_{3,k}$ as follows.

$\begin{split}\mathcal{T}_{1}&=-4\lambda_{\Sigma}\left(\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{1,k}\right)+\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}\right)\right)\\\
&\mathop{\geq}\limits^{(a)}-4\lambda_{\Sigma}\left(\sum_{k=1}^{\infty}k\frac{1}{p_{0}}\mathbb{P}\left(E_{3,k-1}\right)+\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}\right)\right)\\\
&=-4\lambda_{\Sigma}\left(\sum_{k=1}^{\infty}(k-1+1)\frac{1}{p_{0}}\mathbb{P}\left(E_{3,k-1}\right)+\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}\right)\right)\\\
&=-4\lambda_{\Sigma}\left(\left(\frac{1}{p_{0}}+1\right)\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}\right)+\frac{1}{p_{0}}\sum_{k=1}^{\infty}\mathbb{P}\left(E_{3,k-1}\right)\right)\\\
&\mathop{\geq}\limits^{(b)}-8\left(\left(\frac{1}{p_{0}}+1\right)\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}\right)+\frac{1}{p_{0}}\right),\end{split}$
(31)

where inequality (a) follows from Eq. (30); inequality (b) holds since
$\sum_{k=1}^{\infty}\mathbb{P}\left(E_{3,k-1}\right)\leq 1$ and
$\lambda_{\Sigma}\leq 2$.

Now, let us reintroduce the superscript (ϵ). Note that

$\lim_{\epsilon\downarrow
0}\sum_{k=1}^{\infty}k\mathbb{P}\left(E_{3,k}^{(\epsilon)}\right)\mathop{=}\limits^{(a)}\lim_{\epsilon\downarrow
0}\sum_{k=1}^{\infty}k\mathbb{P}\left((\overline{Q}_{1}^{+})^{(\epsilon)}=k,\overline{U}_{2}^{(\epsilon)}>0\right)\mathop{=}\limits^{(b)}0,$

where equality (a) follows from the fact that $U_{2}(t)>0$ implies that
$Q_{2}(t+1)=0$ for all $t\geq 0$ and $\epsilon>0$; equality (b) follows
directly from Eq. (28). Taking liminf on both sides of Eq. (31), yields

$\liminf_{\epsilon\downarrow 0}\mathcal{T}_{1}\geq-\frac{8}{p_{0}}.$ (32)

Next, let us turn to consider the second term on the right-hand side of Eq.
(29)

$\begin{split}\mathcal{T}_{2}&=\mathbb{E}\left[\left(\overline{A}_{1}^{(\epsilon)}-\overline{A}_{2}^{(\epsilon)}\right)^{2}\right]\\\
&\mathop{=}\limits^{(a)}\mathbb{E}\left[\left(\overline{A}_{1}^{(\epsilon)}+\overline{A}_{2}^{(\epsilon)}\right)^{2}\right]\\\
&=\mathbb{E}\left[\left(A_{\Sigma}^{(\epsilon)}\right)^{2}\right]\\\
&=\left(\sigma_{\Sigma}^{(\epsilon)}\right)^{2}+\left(\lambda_{\Sigma}^{(\epsilon)}\right)^{2},\end{split}$

where (a) holds since $A_{1}(t)A_{2}(t)=0$ for all $t\geq 0$ and $\epsilon>0$.
Taking liminf of both sides, we obtain

$\liminf_{\epsilon\downarrow 0}\mathcal{T}_{2}=\sigma_{\Sigma}^{2}+4,$ (33)

since $\left(\sigma_{\Sigma}^{(\epsilon)}\right)^{2}$ approaches to
$\sigma_{\Sigma}^{2}$.

We are left with the third term on the right-hand side of Eq. (29).

$\begin{split}\mathcal{T}_{3}&=-\mathbb{E}\left[\left(\overline{U}_{1}^{(\epsilon)}-\overline{U}_{2}^{(\epsilon)}\right)^{2}\right]\\\
&\geq-\mathbb{E}\left[\left(\overline{U}_{1}^{(\epsilon)}\right)^{2}+\left(\overline{U}_{2}^{(\epsilon)}\right)^{2}\right]\\\
&\mathop{=}\limits^{(a)}-\mathbb{E}\left[\left(\overline{U}_{1}^{(\epsilon)}\right)+\left(\overline{U}_{2}^{(\epsilon)}\right)\right]\\\
&\mathop{=}\limits^{(b)}-\epsilon,\end{split}$

where (a) follows from the fact the unused service is at most one; (b) follows
directly from the property of unused service shown in Eq. (24) of Lemma 10.
Taking liminf of both sides, yields,

$\liminf_{\epsilon\downarrow 0}\mathcal{T}_{3}\geq 0.$ (34)

Now, taking liminf on both sides of Eq. (C) and using Eq. (29), yields,

$\begin{split}\liminf_{\epsilon\downarrow
0}B^{(\epsilon)}&\geq\frac{1}{2}\liminf_{\epsilon\downarrow
0}\left(\mathcal{T}_{1}+\mathcal{T}_{2}+\mathcal{T}_{3}\right)\\\
&\mathop{\geq}\limits^{(a)}\frac{1}{2}\left(\sigma_{\Sigma}^{2}+4-\frac{8}{p_{0}}\right).\\\
&\mathop{>}\limits^{(b)}0.\end{split}$

where (a) follows from the super-additivity of liminf and Eqs. (32), (33) and
(34); inequality (b) comes from the fact that the arrival process is in
$\mathcal{A}$. Therefore, we arrive that $\limsup_{\epsilon\downarrow
0}\overline{B}^{(\epsilon)}>0$, which contradicts with our assumption that
$\limsup_{\epsilon\downarrow 0}\overline{B}^{(\epsilon)}=0$. Therefore, the
result in Theorem 2 must hold, i.e., JIQ is not heavy-traffic delay optimal
for a two-server homogeneous settings.∎

##  Appendix D Additional Simulation Results

In this section, we provide additional simulation results with various
parameters and system sizes. More specifically, we provide the same set of
simulation results for larger system sizes, e.g., $50$ and $100$ servers, as
well as different combinations of arrival and service distributions for each
time slot, including Poisson arrival-constant service, Poisson arrival-bursty
service and bursty arrival-Poisson service.

###  D.1 Throughput Performance

In this subsection, we further explore the throughput performance for
different system sizes and different combinations of arrival and service
process. In particular, we present the results for $50$ and $100$
heterogeneous servers under Poisson arrival and Poisson service. Moreover,
simulation results for $10$ heterogeneous servers under Poisson arrival-
constant service, Poisson arrival-bursty service and bursty arrival-Poisson
arrival are all presented. Similar trend can be observed in all these results.

Figure 10: Throughput performance in $50$ heterogeneous servers.

Figure 11: Throughput performance in $100$ heterogeneous servers.

Let us first look at the throughput performance under larger system sizes.
Figures 11 and 11 demonstrate the throughput performance under $50$ and $100$
heterogeneous servers, respectively. As before, the servers are equally
divided into two server pools with rate $1$ and $10$. A turning point in each
curve indicates that the load approaches the throughput region boundary of the
corresponding policy. It can be seen that power-of-$d$ policy has nearly the
same throughput region as in the case of $10$ servers. For JIQ policy, it
remains unstable for heavy loads, yet it tends to have a larger throughput
region than that in the case of $10$ servers, which is intuitive since it is
more likely for an arrival to find an idle queue in a larger system. As
expected, our proposed JBTG-$d$ policy is able to stabilize the system for all
the considered loads which again agrees with our theoretical results. The
JBT-$d$ policy can also ensure stability in these two cases as in the case of
$10$ servers. The use of memory for improving stability of power-of-$d$ is
verified by the power-of-$d$ with memory policies for both cases
(SQ($2$,$9$),SQ($2$,$15$)).

Figure 12: Throughput performance in $10$ heterogeneous servers with Poisson
arrival and constant service.

Figure 13: Throughput performance in $10$ heterogeneous servers with Poisson
arrival and bursty service.

Figure 14: Throughput performance in $10$ heterogeneous servers with bursty
arrival and Poisson service.

Figure 15: Delay performance under $50$ homogeneous servers.

Figure 16: Delay performance under $100$ homogeneous servers.

Figure 17: Heavy-traffic Delay performance under $50$ homogeneous servers.

Figure 18: Heavy-traffic Delay performance under $100$ homogeneous servers.

Figure 19: Delay performance under $10$ homogeneous servers with Poisson
arrival and constant service.

Figure 20: Heavy-traffic delay performance under $10$ homogeneous servers with
Poisson arrival and constant service.

Figure 21: Delay performance under $10$ homogeneous servers with Poisson
arrival and bursty service.

Figure 22: Heavy-traffic delay performance under $10$ homogeneous servers with
Poisson arrival and bursty service.

Figure 23: Delay performance under $10$ homogeneous servers with bursty
arrival and Poisson service.

Figure 24: Heavy-traffic delay performance under $10$ homogeneous servers with
bursty arrival and Poisson service.

Figure 25: Message per new job arrival under $50$ homogeneous servers with
respect to $T$.

Figure 26: Delay performance under $50$ homogeneous servers with respect to
$T$.

Figure 27: Message per new job arrival under $100$ homogeneous servers with
respect to $T$.

Figure 28: Delay performance under $100$ homogeneous servers with respect to
$T$.

Figure 29: Message per new job arrival under $10$ homogeneous servers with
Poisson arrival and constant service versus $T$.

Figure 30: Delay performance under $10$ homogeneous servers with Poisson
arrival and constant service versus $T$.

Figure 31: Message per new job arrival under $10$ homogeneous servers with
Poisson arrival and bursty service versus $T$.

Figure 32: Delay performance under $10$ homogeneous servers with Poisson
arrival and bursty service versus $T$.

Figure 33: Message per new job arrival under $10$ homogeneous servers with
bursty arrival and Poisson service versus $T$.

Figure 34: Delay performance under $10$ homogeneous servers with bursty
arrival and Poisson service versus $T$.

Let us turn to look at the throughput performance for different combinations
of arrival and service distribution in $10$ heterogeneous servers. As before,
half of them have rate $1$ and half of them have rate $10$. Figure 13 shows
the case of Poisson arrival and constant service. Similar results as before
can be observed in this case. Figure 13 shows the results for Poisson arrival
and bursty service. In particular, the exogenous arrival at each time slot is
drawn from Poisson distribution. The number of potential service at each time
slot for servers with rate $1$ can be either $0$ or $10$, and the number of
potential service at each time for servers with rate $10$ is either $0$ or
$15$ in this setting. Again, similar results are observed in this case. Figure
15 illustrates the case of bursty arrival and Poisson service. In particular,
the exogenous arrival at each time slot is either $0$ or $60$, and the
potential service at each time slot of each server is drawn from a Poisson
distribution with the corresponding rate. It can be easily seen that similar
results hold in this case too.

###  D.2 Delay Performance

In this subsection, we provide additional delay performance for larger system
sizes and different combinations of arrival and service process in homogeneous
servers. In particular, we present the results for $50$ and $100$ homogeneous
servers under Poisson arrival and Poisson service. Moreover, simulation
results for $10$ homogeneous servers under Poisson arrival-constant service,
Poisson arrival-bursty service and bursty arrival-Poisson arrival are all
presented. As can be seen, our proposed policy JBT-$d$ achieves good
performance in all these cases. All the servers have rate $1$.

Let us first look at the delay performance under a load which ranges from
light to heavy traffic regime for $50$ and $100$ homogeneous servers, as shown
in Figures 15 and respectively. It can be easily seen that our proposed policy
JBT-$d$ achieves a much smaller response time over all the considered loads
than that of power-of-$2$ and power-of-$2$ with memory (SQ($2$,$9$), which
utilizes almost the same amount of memory as the JBT-$d$). Moreover, the
performance under JBT-$d$ is as good as JSQ for an even larger range of loads
compared to the $10$ servers case, and the same trend holds for JIQ, which can
be explained that as $N$ becomes larger it is more likely to find an idle
queue for a given load. However, as load approaches to $1$, the performance of
JIQ degrades badly as discussed in the following.

Now, let us get a closer look for the heavy-traffic regime, i.e., $\rho>0.9$,
under $50$ and $100$ servers, as shown in Figures 17 and 19, respectively. We
can easily observe that, as $\rho$ approaches to $1$, the performance of JIQ
degrades substantially due to the lack of idle servers while our proposed
policy JBT-$d$ remains quite close to JSQ in both cases.

Next, it is time for us to take a look at the delay performance under
different combinations of arrival and service distribution in $10$ servers
case. Figures 19 and 21 shows the same set of results for the case of Poisson
arrival and constant service. Same trend can be observed as before. Figures 21
and 23 demonstrates the same set of results for the case of Poisson arrival
and bursty service. In particular, the potential number of jobs served at each
time slot is either $0$ or $10$ in this bursty setting. It can be easily seen
that JIQ degrades much faster in this case as the load increases than that in
the Poisson case. Figures 23 and 25 illustrates the same set of results for
the case of bursty arrival and Poisson service. In particular, the exogenous
number of arrivals for each time slot is taken either $0$ or $12$. From the
figures, we can observe that burst in arrival process also aggravates the
performance of response time, and yet our proposed policy JBT-$d$ still has
superior performance over other polices from light traffic to heavy traffic.

###  D.3 Message Overhead

In this subsection, further results on message overhead for different settings
are provided. As before, we present results for larger system sizes, e.g.,
$50$ and $100$ server, with Poisson arrival-Poisson service. Meanwhile, the
message overhead for Poisson arrival-constant service, Poisson arrival-bursty
service and bursty arrival-Poisson service in $10$ servers are also given. All
these cases demonstrate the low message rate property of our proposed JBT-$d$
policy.

Figure 25 illustrates the message overhead with respect to $T$ under $50$
homogeneous servers at load $\rho=0.99$, and the corresponding delay
performance is shown in Figure 27. It can be easily seen that the JBT-$d$
policy is able to achieve a much better performance with a much lower message
rate in a broad range of $T$ when compared to power-of-$d$ and power-of-$d$
with memory since both of them have message rate of $4$. Moreover, the message
overhead of JBT-$d$ for $T>100$ is quite close to that of JIQ, which indicates
that JBT-$d$ is capable of guaranteeing heavy-traffic delay optimality with
just a slightly more message. In addition, we can observe that the delay
performance is almost invariant with respect to $T$, especially when $T>100$,
which implies that we are allowed to use a large $T$ without sacrificing the
delay performance. Figures 27 and 29 shows the case of $100$ homogeneous
servers. We can see that the same results hold in this case as before.

Now, let us look at the message overhead under different arrival and service
distributions. Figures 29 and 31 show the case for Poisson arrival and
constant service with rate $1$. Nearly the same results hold in this case as
well. Figures 31 and 33 demonstrates the case of Poisson arrival and bursty
service. In particular, the burst in service is simulated by allowing the
number of potential service at each time slot to take either $0$ or $10$.
First thing to note is that burst in service degrades the performance of JIQ
significantly, which agrees with our theoretical results. From the perspective
of message rate, we can see that in this case message rate of JIQ is quite
close to $0$, which accounts for the poor performance. In contrast, our
proposed JBT-$d$ policy still achieve good performance with a much lower
performance in this case. It is worth noting that although the delay
performance of power-of-$d$ with memory (SQ($2$,$3$)) is slightly better than
that of our JBT-$d$ policy, its message rate is four fold as much as the rate
of JBT-$d$ when $T>50$. As before, a larger $d$ does not necessarily mean a
larger message rate when $T$ is large. Figures 33 and 34 illustrate the case
of bursty arrival and Poisson service. More specifically, the exogenous
arrival for each time slot is either $0$ or $12$, and the potential service
for each time slot at each server is drawn from Poisson distribution with mean
$1$. The behavior of message overhead is nearly the same as before. That is,
for a large range of $T$, the message rate of JBT-$d$ is strictly less than
$1$, which implies that our proposed JBT-$d$ is able to achieve heavy-traffic
delay optimality by using a message rate that is less than $1$. From the delay
performance in Figure 34, we can observe that JBT-$d$ still has good
performance in bursty arrival case and it also suggests that a large $d$ would
be better for this bursty arrival setting.

##  Appendix E Confidence Intervals

In this section, we provide the $95\%$ confidence interval for all average-
based simulation results. In particular, it is obtained by running the
simulation for $9,000,000$ time slots, which are divided into $30$ batches.
Then, the standard batch means method is adopted to compute the confidence
interval for all the results.

###  E.1 Throughput performance

In this subsection, we present the $95\%$ confidence intervals for all the
simulation results of throughput performance, which is summarized in Table 2.

###  E.2 Delay performance

In this subsection, we present the $95\%$ confidence intervals for all the
simulation results of delay performance, which is summarized in Table 3.

###  E.3 Message overhead

In this subsection, we present the $95\%$ confidence intervals for all the
simulation results of message overhead and the corresponding delay
performance, which are summarized in Table 4.

Table 2: The simulation result of throughput performance and its corresponding confidence interval. Figures | Confidence interval table
---|---
Figure 3 | Table 5
Figure 11 | Table 6
Figure 11 | Table 7
Figure 13 | Table 8
Figure 13 | Table 9
Figure 15 | Table 10
Table 3: The simulation result of delay performance and its corresponding confidence interval. Figures | Confidence interval table
---|---
Figures 5 and 5 | Table 11
Figures 15 and 7(17) | Table 12
Figures 17 and 19 | Table 13
Figures 19 and 21 | Table 14
Figures 21 and 7(23) | Table 15
Figures 23 and 25 | Table 16
Table 4: The simulation result of message overhead (delay performance) and its corresponding confidence interval. Figures | Confidence interval table
---|---
Figures 9 and 9 | Table 17 and Table 18
Figures 25 and 27 | Table 19 and Table 20
Figures 27 and 29 | Table 21 and Table 22
Figures 29 and 31 | Table 23 and Table 24
Figures 31 and 33 | Table 25 and Table 26
Figures 33 and 34 | Table 27 and Table 28
Table 5: The $95\%$ confidence interval of response time under $10$ heterogeneous servers with Poisson arrival and Poisson service. Policy | $\rho=0.5$ | $\rho=0.7$ | $\rho=0.9$ | $\rho=0.95$ | $\rho=0.99$
---|---|---|---|---|---
JSQ | $3.363\pm 0.003$ | $3.894\pm 0.003$ | $4.998\pm 0.007$ | $5.804\pm 0.010$ | $8.093\pm 0.083$
JBTG-10 | $2.156\pm 0.005$ | $3.246\pm 0.005$ | $5.732\pm 0.027$ | $8.151\pm 0.100$ | $13.562\pm 0.380$
JBTG-2 | $2.192\pm 0.008$ | $3.414\pm 0.034$ | $6.159\pm 0.079$ | $8.294\pm 0.142$ | $13.902\pm 0.412$
JBT-2 | $3.715\pm 0.084$ | $5.288\pm 0.113$ | $16.882\pm 1.463$ | $28.904\pm 4.272$ | $53.830\pm 16.384$
SQ(2,3) | $3.968\pm 0.008$ | $6.398\pm 0.151$ | $12.942\pm 0.052$ | $18.200\pm 0.103$ | $31.200\pm 0.291$
SQ(2) | $1.9\mathrm{e}4\pm 0.4\mathrm{e}4$ | $4.2\mathrm{e}4\pm 0.9\mathrm{e}4$ | $5.5\mathrm{e}4\pm 1.1\mathrm{e}4$ | $5.7\mathrm{e}4\pm 1.2\mathrm{e}4$ | $5.9\mathrm{e}4\pm 1.2\mathrm{e}4$
JIQ | $3.112\pm 0.003$ | $3.896\pm 0.006$ | $1.1\mathrm{e}4\pm 2.3\mathrm{e}3$ | $2.7\mathrm{e}4\pm 5.6\mathrm{e}3$ | $3.9\mathrm{e}4\pm 8.1\mathrm{e}3$
Table 6: The $95\%$ confidence interval of response time under $50$ heterogeneous servers with Poisson arrival and Poisson service.  Policy | $\rho=0.5$ | $\rho=0.7$ | $\rho=0.9$ | $\rho=0.95$ | $\rho=0.99$
---|---|---|---|---|---
JSQ | $17.841\pm 0.007$ | $20.564\pm 0.008$ | $23.369\pm 0.008$ | $24.510\pm 0.016$ | $27.290\pm 0.051$
JBTG-50 | $12.115\pm 0.019$ | $17.134\pm 0.021$ | $22.449\pm 0.015$ | $26.091\pm 0.069$ | $42.996\pm 0.685$
JBTG-10 | $12.119\pm 0.019$ | $17.147\pm 0.022$ | $22.837\pm 0.085$ | $26.036\pm 0.146$ | $35.153\pm 0.441$
JBT-10 | $16.659\pm 0.011$ | $20.184\pm 0.044$ | $25.544\pm 0.159$ | $34.323\pm 0.496$ | $64.107\pm 5.690$
SQ(2,9) | $22.554\pm 0.077$ | $39.086\pm 0.071$ | $77.287\pm 0.223$ | $103.829\pm 0.511$ | $170.932\pm 1.552$