text:
\begin{document}

\title{Scalable load balancing in networked systems:\\ A survey of recent advances}
\author[1]{Mark van der Boor}
\author[1,2]{Sem C.~Borst}
\author[1]{\\ Johan S.H.~van Leeuwaarden,}
\author[1]{Debankur Mukherjee}
\affil[1]{Eindhoven University of Technology, The Netherlands}
\affil[2]{Nokia Bell Labs, Murray Hill, NJ, USA}

\date{\today}

\maketitle

\begin{abstract}{
The basic load balancing scenario involves a single dispatcher where
tasks arrive that must immediately be forwarded to one
of $N$~single-server queues. We discuss recent advances on scalable load balancing schemes
which provide favorable delay performance when $N$ grows large,
and yet only require minimal implementation overhead.

Join-the-Shortest-Queue (JSQ) yields vanishing delays as $N$ grows large,
as in a centralized queueing arrangement, but involves a prohibitive communication burden. In contrast,
power-of-$d$ or JSQ($d$) schemes that assign an incoming task
to a server with the shortest queue among $d$~servers selected
uniformly at random require little communication, but lead to constant delays.
In order to examine this fundamental trade-off between delay
performance and implementation overhead, we consider
JSQ($d(N)$) schemes where the diversity parameter $d(N)$ depends on~$N$ and investigate
what growth rate of $d(N)$ is required to asymptotically
match the optimal JSQ performance on fluid and diffusion
scale.

Stochastic coupling techniques and stochastic-process limits play an instrumental role
in establishing the asymptotic optimality. We demonstrate how this methodology carries over
to infinite-server settings, finite buffers, multiple dispatchers, servers arranged on graph topologies, and token-based load balancing including the popular Join-the-Idle-Queue (JIQ) scheme. In this way we provide a broad overview of the many recent advances in the field. This survey extends the short review presented at ICM 2018~\cite{BBLM18}.}
\end{abstract}

\newpage
\setcounter{tocdepth}{2}
{\small \tableofcontents}

\newpage

\section{Introduction}\label{sec:intro}

In this survey we review scalable load balancing algorithms
(LBAs) which achieve excellent delay performance in large-scale
systems and yet only involve low implementation overhead.
LBAs play a critical role in distributing service requests or tasks
(e.g.~compute jobs, data base look-ups, file transfers) among servers
or distributed resources in parallel-processing systems.
The analysis and design of LBAs has attracted strong attention in recent
years, mainly spurred by crucial scalability challenges arising
in cloud networks and data centers with massive numbers of servers.
LBAs can be broadly categorized as static, dynamic, or some intermediate
blend, depending on the amount of feedback or state information
(e.g.~congestion levels) that is used in allocating tasks.
The use of state information naturally allows dynamic policies
to achieve better delay performance, but also involves higher
implementation complexity and a substantial communication burden.
The latter issue is particularly pertinent in cloud networks and data
centers with immense numbers of servers handling a huge influx
of service requests.
In order to capture the large-scale context, we examine scalability
properties through the prism of asymptotic scalings where the system
size grows large, and identify LBAs which strike an optimal balance
between delay performance and implementation overhead in that regime.

The most basic load balancing scenario consists of $N$~identical
parallel servers and a dispatcher where tasks arrive that must
immediately be forwarded to one of the servers.
Tasks are assumed to have unit-mean exponentially distributed service
requirements, and the service discipline at each server is supposed
to be oblivious to the actual service requirements.
In this canonical setup, the celebrated Join-the-Shortest-Queue (JSQ)
policy has several strong stochastic optimality properties.
In particular, the JSQ policy achieves the minimum mean overall
delay among all non-anticipating policies that do not have any
advance knowledge of the service requirements \cite{EVW80,Winston77}.
In order to implement the JSQ policy however, a dispatcher requires
instantaneous knowledge of all the queue lengths, which may involve
a prohibitive communication burden with a large number of servers~$N$.

This poor scalability has motivated consideration of JSQ($d$) policies,
where an incoming task is assigned to a server with the shortest queue
among $d \geq 2$ servers selected uniformly at random.
Note that this involves an exchange of $2 d$ messages per task,
irrespective of the number of servers~$N$.
Results in Mitzenmacher~\cite{Mitzenmacher01} and Vvedenskaya
{\em et al.}~\cite{VDK96} indicate that even sampling as few as $d = 2$
servers yields significant performance enhancements over purely
random assignment ($d = 1$) as $N$ grows large, which is commonly
referred to as the \emph{power-of-two} or \emph{power-of-choice} effect.
Specifically, when tasks arrive at rate $\lambda N$,
the queue length distribution at each individual server exhibits
super-exponential decay for any fixed $\lambda < 1$ as $N$ grows large, a considerable improvement
compared to exponential decay for purely random assignment.

The diversity parameter~$d$ thus induces
a fundamental trade-off between the amount of communication overhead
and the delay performance.
Specifically, a random assignment policy does not entail any
communication burden, but the mean waiting time remains \emph{constant}
as $N$ grows large for any fixed $\lambda > 0$.
In contrast, a nominal implementation of the JSQ policy (without
maintaining state information at the dispatcher) involves $2 N$
messages per task, but the mean waiting time \emph{vanishes}
as $N$ grows large for any fixed $\lambda < 1$.
Although JSQ($d$) policies with $d \geq 2$ yield major performance
improvements over purely random assignment while reducing the
communication burden by a factor O($N$) compared to the JSQ policy,
the mean waiting time \emph{does not vanish} in the limit.
Hence, no fixed value of~$d$ will provide asymptotically optimal
delay performance.
This is evidenced by results of Gamarnik {\em et al.}~\cite{GTZ16}
indicating that in the absence of any memory at the dispatcher the
communication overhead per task \emph{must increase} with~$N$ in order
for any scheme to achieve a zero mean waiting time in the limit.

We will explore the intrinsic trade-off between delay performance
and communication overhead as governed by the diversity parameter~$d$,
in conjunction with the relative load~$\lambda$.
The latter trade-off is examined in an asymptotic regime where not
only the overall task arrival rate is assumed to grow with~$N$,
but also the diversity parameter is allowed to depend on~$N$.
We write $\lambda(N)$ and $d(N)$, respectively, to explicitly reflect
that, and investigate what growth rate of $d(N)$ is required,
depending on the scaling behavior of $\lambda(N)$, in order to achieve
a zero mean waiting time in the limit.
The analysis covers both fluid-scaled and diffusion-scaled versions of
the queue length process in regimes where $\lambda(N) / N \to \lambda < 1$
and $(N - \lambda(N)) / \sqrt{N} \to \beta > 0$ as $N \to \infty$,
respectively.
We establish that the limiting processes are insensitive to the exact
growth rate of $d(N)$, as long as the latter is sufficiently fast,
and in particular coincide with the limiting processes for the JSQ policy.
This demonstrates
that the optimality of the JSQ policy can asymptotically be preserved
while dramatically lowering the communication overhead.

We will also consider network
scenarios where the $N$~servers are assumed to be inter-connected
by some underlying graph topology $G_N$.
Tasks arrive at the various servers as independent Poisson processes
of rate~$\lambda$, and each incoming task is assigned to whichever
server has the shortest queue among the one where it appears and its
neighbors in $G_N$.
In case $G_N$ is a clique (fully connected graph), each incoming task is assigned to the
server with the shortest queue across the entire system,
and the behavior is equivalent to that under the JSQ policy.
The stochastic optimality properties of the JSQ policy
thus imply that the queue length process in a clique will be `better'
than in an arbitrary graph~$G_N$.
We will establish sufficient conditions for the fluid-scaled
and diffusion-scaled versions of the queue length process in
an arbitrary graph to be equivalent to the limiting processes
in a clique as $N \to \infty$.
The conditions demonstrate that the optimality of a clique can
be asymptotically preserved while dramatically reducing the number
of connections, provided the graph $G_N$ is suitably random.

While a zero waiting time can be achieved in the limit by sampling
only $d(N) \ll N$ servers, the amount of communication overhead
in terms of $d(N)$ must still grow with~$N$.
This may be explained from the fact that a large number of servers
need to be sampled for each incoming task to ensure that at least one
of them is found idle with high probability.
This can be avoided by introducing memory
at the dispatcher, in particular maintaining a record of vacant servers,
and assigning tasks to idle servers, if there are any.
This so-called Join-the-Idle-Queue (JIQ) scheme \cite{BB08,LXKGLG11}
has gained huge popularity recently, and can be implemented through
a simple token-based mechanism generating at most one message per task.
As shown by Stolyar~\cite{Stolyar15}, the fluid-scaled queue length
process under the JIQ scheme is equivalent to that under the JSQ policy
as $N \to \infty$, and we will extend this result to the
diffusion-scaled queue length process.
Thus, the use of memory allows the JIQ scheme to achieve asymptotically
optimal delay performance with minimal communication overhead.
In particular, ensuring that tasks are assigned to idle servers
whenever available is sufficient to achieve asymptotic optimality,
and using any additional queue length information yields no meaningful
performance benefits on the fluid or diffusion levels.

Stochastic coupling techniques play an instrumental role in the proofs
of the above-described universality and asymptotic optimality properties.
A direct analysis of the queue length processes under a JSQ($d(N)$)
policy, in a load balancing graph $G_N$, or under the JIQ scheme
is confronted with formidable obstacles, and does not seem tractable.
As an alternative route, we leverage novel stochastic coupling
constructions to relate the relevant queue length processes to the
corresponding processes under a JSQ policy, and show that the
deviation between these processes is asymptotically negligible
under suitable assumptions on $d(N)$ or $G_N$.

While the stochastic coupling schemes provide an effective
and overarching approach, they defy a systematic recipe and involve
some degree of ingenuity and customization.
Indeed, the specific coupling arguments that we develop are not only
different from those that were originally used in establishing the
stochastic optimality properties of the JSQ policy, but also differ
in critical ways between a JSQ($d(N)$) policy, a load balancing graph
$G_N$, and the JIQ scheme.
Yet different coupling constructions are devised for model variants
with infinite-server dynamics that we will discuss in Section~\ref{bloc}.

The survey is organized as follows.
In Section~\ref{spec} we discuss various LBAs and evaluate
their scalability properties.
In Section~\ref{sec:powerofd} we introduce some useful preliminary concepts,
and then review fluid and diffusion limits for the JSQ policy as well
as JSQ($d$) policies with a fixed value of~$d$.
In Section~\ref{univ} we explore the trade-off between delay
performance and communication overhead as function of the diversity
parameter~$d$, in conjunction with the relative load.
In particular, we establish asymptotic universality properties
for JSQ($d$) policies, which are extended to systems with server pools
and network scenarios in Sections~\ref{bloc} and~\ref{networks},
respectively.
In Section~\ref{token} we establish asymptotic optimality properties
for the JIQ scheme.
We discuss somewhat related redundancy policies and alternative scaling
regimes and performance metrics in Section~\ref{miscellaneous}.
The survey is concluded in Section~\ref{sec:ext} with a discussion of yet
further extensions and several open problems and emerging research
directions.

\section{Scalability spectrum}
\label{spec}

In this section we review a wide spectrum of LBAs and examine their
scalability properties in terms of the delay performance vis-a-vis
the associated implementation overhead in large-scale systems.

\subsection{Basic model}

Throughout this section and most of the paper, we focus on a basic
scenario with $N$ parallel single-server infinite-buffer queues
and a single dispatcher where tasks arrive as a Poisson process
of rate~$\lambda(N)$, as depicted in Figure~\ref{figJSQ}.
Arriving tasks cannot be queued at the dispatcher,
and must immediately be forwarded to one of the servers.

Tasks are assumed to have unit-mean exponentially distributed service
requirements, and the service discipline at each server is supposed
to be oblivious to the actual service requirements.
This is the supermarket model described in Section~\ref{sec:intro}.

When tasks do not get served and never depart but simply accumulate,
the above setup corresponds to a so-called balls-and-bins model,
and we will further elaborate on the connections and differences
with work in that domain in Subsection~\ref{ballsbins}.

\subsection{Asymptotic scaling regimes}
\label{asym}

An exact analysis of the delay performance is quite involved,
if not intractable, for all but the simplest LBAs.
Numerical evaluation or simulation are not straightforward either,
especially for high load levels and large system sizes.
A common approach is therefore to consider various limit regimes,
which not only provide mathematical tractability and illuminate the
fundamental behavior, but are also natural in view of the typical
conditions in which cloud networks and data centers operate.
One can distinguish several asymptotic scalings that have been used
for these purposes:

\begin{enumerate}[{\normalfont (i)}]
\item In the classical heavy-traffic regime, $\lambda(N) = \lambda N$
with a fixed number of servers~$N$ and a relative load~$\lambda$
that tends to one in the limit.
\item In the conventional large-capacity or many-server regime,
the relative load $\lambda(N) / N$ approaches a constant $\lambda < 1$
as the number of servers~$N$ grows large.
\item The popular Halfin-Whitt regime~\cite{HW81} combines heavy traffic
with a large capacity, with
\begin{equation}
\label{eq:HW}
\frac{N - \lambda(N)}{\sqrt{N}} \to \beta > 0 \mbox{ as } N \to \infty,
\end{equation}
so the relative capacity slack behaves as $\beta / \sqrt{N}$
as the number of servers~$N$ grows large.
\item The so-called non-degenerate slow-down regime~\cite{Atar12}
involves $N - \lambda(N) \to \gamma > 0$, so the relative capacity
slack shrinks as $\gamma / N$ as the number of servers~$N$ grows large.
\end{enumerate}

The term non-degenerate slow-down refers to the fact that in the
context of a centralized multi-server queue, the mean waiting time
in regime (iv) tends to a strictly positive constant as $N \to \infty$,
and is thus of similar magnitude as the mean service requirement.
In contrast, in regimes (ii) and (iii), the mean waiting time in
a multi-server queue decays exponentially fast in~$N$ or is of the
order~$1 / \sqrt{N}$, respectively as $N \to \infty$,
while in regime (i) the mean waiting time grows arbitrarily large
relative to the mean service requirement.

In the context of a centralized M/M/$N$ queue, scalings (ii), (iii) and (iv) are commonly referred to as
Quality-Driven (QD), Quality-and-Efficiency-Driven (QED)
and Efficiency-Driven (ED) regimes.
These terms reflect that (ii) offers excellent service quality
(vanishing waiting time), (iv) provides high resource efficiency
(utilization approaching one), and (iii) achieves a combination
of these two, providing the best of both worlds.

In the present paper we will focus on scalings (ii) and (iii),
and occasionally also refer to these as fluid and diffusion scalings,
since it is natural to analyze the relevant queue length process
on fluid scale ($1 / N$) and diffusion scale ($1 / \sqrt{N}$)
in these regimes, respectively.
We will not provide a detailed account of scalings (i) and (iv),
which do not capture the large-scale perspective and do not allow
for low delays, respectively, but we will briefly mention some results
for these regimes in Subsections~\ref{classical}
and~\ref{nondegenerate}.

An important issue in the context of scaling limits is the rate
of convergence and the accuracy for finite-size systems.
Some interesting results for the accuracy of mean-field approximations
for interacting-particle systems including load balancing
models may be found in recent work of Gast~\cite{Gast17}, Gast \& Van Houdt~\cite{GH18}, and Ying \cite{Ying16,Ying17}.

\subsection{Basic load balancing algorithms}

\subsubsection{Random assignment: N independent M/M/1 queues}
\label{random}

One of the most basic LBAs is to assign each arriving task to a server
selected uniformly at random.
In that case, the various queues collectively behave as
$N$~independent M/M/1 queues, each with arrival rate $\lambda(N) / N$
and unit service rate.
In particular, at each of the queues, the total number of tasks in
stationarity has a geometric distribution with parameter $\lambda(N) / N$.
By virtue of the PASTA property, the probability that an arriving task
incurs a non-zero waiting time is $\lambda(N) / N$.
The mean number of waiting tasks (excluding the possible task in service)
at each of the queues is $\frac{\lambda(N)^2}{N (N - \lambda(N))}$,
so the total mean number of waiting tasks is
$\frac{\lambda(N)^2}{N - \lambda(N)}$, which by Little's law implies that
the mean waiting time of a task is $\frac{\lambda(N)}{N - \lambda(N)}$.
In particular, when $\lambda(N) = N \lambda$, the probability that
a task incurs a non-zero waiting time is $\lambda$,
and the mean waiting time of a task is $\frac{\lambda}{1 - \lambda}$,
independent of~$N$, reflecting the independence of the various queues.

As we will see later, a broad range of queue-aware LBAs can deliver
a probability of a non-zero waiting time and a mean waiting time
that vanish asymptotically.
While a random assignment policy is evidently not competitive
with such queue-aware LBAs, it still plays a relevant role due to
the strong degree of tractability inherited from its simplicity.
For example, the queue process under purely random assignment can
be shown to provide an upper bound (in a stochastic majorization sense)
for various more involved queue-aware LBAs for which even stability
may be difficult to establish directly, yielding conservative
performance bounds and stability guarantees.

A slightly better LBA is to assign tasks to the servers in
a Round-Robin manner, dispatching every $N$-th task to the same server.
In the fluid regime where $\lambda(N) = N \lambda$, the inter-arrival
time of tasks at each given queue will then converge to a constant
$1 / \lambda$ as $N \to \infty$.
Thus each of the queues will behave as a D/M/1 queue in the limit,
and the probability of a non-zero waiting time and the mean waiting
time will be somewhat lower than under purely random assignment.
However, both the probability of a non-zero waiting time and the mean
waiting time will still tend to strictly positive values and not vanish
as $N \to \infty$.

\subsubsection{Join-the-Shortest Queue (JSQ)}
\label{ssec:jsq}

Under the Join-the-Shortest-Queue (JSQ) policy, each arriving task is
assigned to the server with the currently shortest queue.
In the basic model described above, the JSQ policy has several strong
stochastic optimality properties, and yields the `most balanced
and smallest' queue process among all non-anticipating policies that
do not have any advance knowledge of the service requirements
\cite{EVW80,Winston77}.

\subsubsection{Join-the-Smallest-Workload (JSW): centralized M/M/N queue}
\label{ssec:jsw}

Under the Join-the-Smallest-Workload (JSW) policy, each arriving task
is assigned to the server with the currently smallest workload.
Note that this is an anticipating policy, since it requires advance
knowledge of the service requirements of all the tasks in the system.
Further observe that this policy (myopically) minimizes the waiting
time for each incoming task, and mimicks the operation of a centralized
$N$-server queue with a FCFS discipline.
The equivalence with a centralized $N$-server queue with a FCFS
discipline yields a strong optimality property of the JSW policy:
The vector of joint workloads at the various servers observed by each
incoming task is smaller in the Schur convex sense than under any
alternative admissible policy~\cite{FC01}.

It is worth observing that the above optimality properties in fact do
not rely on Poisson arrival processes or exponential service
requirement distributions.
At the same time, these optimality properties do not imply that the
JSW policy minimizes the mean stationary waiting time.
In our setting with Poisson arrivals and exponential service
requirements, however, it can be shown through direct means that the
total number of tasks under the JSW policy is stochastically smaller
than under the JSQ policy.
Even though the JSW policy requires a similar excessive communication
overhead, aside from its anticipating nature, the equivalence
with a centralized FCFS queue means that there cannot be any idle
servers while tasks are waiting and that the total number of tasks behaves
as a birth-death process, which renders it far more tractable than
the JSQ policy.
Specifically, given that all the servers are busy, the total number
of waiting tasks is geometrically distributed with parameter
$\lambda(N) / N$.
The total mean number of waiting tasks is then
$\Pi_W(N, \lambda(N)) \frac{\lambda(N)}{N - \lambda(N)}$,
and the mean waiting time is
$\Pi_W(N, \lambda(N)) \frac{1}{N - \lambda(N)}$,
with $\Pi_W(N, \lambda(N)$ denoting the probability of all servers
being occupied and a task incurring a non-zero waiting time.
This immediately shows that the mean waiting time is smaller
by at least a factor $\lambda(N)$ than for the random assignment
policy considered in Subsection~\ref{random}.

In the large-capacity regime $\lambda(N) = N \lambda$, it can be shown
that the probability $\Pi_W(N, \lambda(N))$ of a non-zero waiting time
decays exponentially fast in~$N$, and hence so does the mean waiting time.
In the Halfin-Whitt heavy-traffic regime~\eqref{eq:HW}, the probability
$\Pi_W(N, \lambda(N))$ of a non-zero waiting time converges to a finite
constant $\Pi_W^*(\beta)$, implying that the mean waiting time of a task
is of the order $1 / \sqrt{N}$, and thus vanishes as $N \to \infty$.

\subsubsection{Power-of-d load balancing (JSQ(d))}
\label{ssec:powerd}

We have seen that the achilles heel of the JSQ policy is its
excessive communication overhead in large-scale systems.
This poor scalability has motivated consideration of so-called
JSQ($d$) policies, where an incoming task is assigned to a server
with the shortest queue among $d$~servers selected uniformly at random.
Results in Mitzenmacher~\cite{Mitzenmacher01} and Vvedenskaya
{\em et al.}~\cite{VDK96} indicate that

in the fluid regime where $\lambda(N) = \lambda N$,
the probability that there are $i$ or more tasks at a given queue is
proportional to $\lambda^{\frac{d^i - 1}{d - 1}}$ as $N \to \infty$,
and thus exhibits super-exponential decay as opposed to exponential decay
for the random assignment policy considered in Subsection~\ref{random}.

As alluded to in Section~\ref{sec:intro}, the diversity parameter~$d$ thus induces
a fundamental trade-off between the amount of communication overhead
and the performance in terms of queue lengths and delays.
A rudimentary implementation of the JSQ policy ($d = N$, without
replacement) involves $O(N)$ communication overhead per task,
but it can be shown that the probability of a non-zero waiting time
and the mean waiting \emph{vanish} as $N \to \infty$, just like
in a centralized queue.
Although JSQ($d$) policies with a fixed parameter $d \geq 2$ yield
major performance improvements over purely random assignment,

the probability of a non-zero waiting time
and the mean waiting time \emph{do not vanish} as $N \to \infty$.

\subsubsection{Token-based mechanisms: Join-the-Idle-Queue (JIQ)}
\label{ssec:jiq}

While a zero waiting time can be achieved in the limit by sampling
only $d(N) \ll N$ servers, the amount of communication overhead
in terms of $d(N)$ must still grow with~$N$.
This can be countered by introducing memory at the dispatcher,
in particular maintaining a record of vacant servers,
and assigning tasks to idle servers as long as there are any,
or to a uniformly at random selected server otherwise.
This so-called Join-the-Idle-Queue (JIQ) scheme \cite{BB08,LXKGLG11}
has received keen interest recently, and can be implemented through
a simple token-based mechanism.
Specifically, idle servers send tokens to the dispatcher to advertize
their availability, and when a task arrives and the dispatcher has
tokens available, it assigns the task to one of the corresponding
servers (and disposes of the token).
Note that a server only issues a token when a task completion leaves
its queue empty, thus generating at most one message per task.
Surprisingly, the mean waiting time and the probability of a non-zero
waiting time vanish under the JIQ scheme in both the fluid and diffusion
regimes, as we will further discuss in Section~\ref{token}.
Thus, the use of memory allows the JIQ scheme to achieve asymptotically
optimal delay performance with minimal communication overhead.

\subsection{Performance comparison}
\label{ssec:perfcomp}

We now present some simulation experiments
to compare the above-described LBAs in terms of delay performance.
\begin{figure}
\includegraphics[width=\linewidth]{DS.pdf}
\caption{Simulation results for mean waiting time $\mathbb{E}[W^N]$
and probability of a non-zero waiting time $p_{\textup{wait}}^N$,
for both a fluid regime and a diffusion regime.}
\label{differentschemes}
\end{figure}
Specifically, we evaluate the mean waiting time and the probability
of a non-zero waiting time in both a fluid regime ($\lambda(N) = 0.9 N$)
and a diffusion regime ($\lambda(N) = N - \sqrt{N}$).
The results are shown in Figure~\ref{differentschemes}.
An overview of the delay performance and overhead associated with
various LBAs is given in Table~\ref{table}.

We are specifically interested in distinguishing two classes of LBAs --
the ones delivering a mean waiting time and probability of a non-zero
waiting time that vanish asymptotically, and the ones that fail to do so
-- and relating that dichotomy to the associated communication
overhead and memory requirement at the dispatcher.
We give these classifications for both the fluid regime and the
diffusion regime.

\paragraph{JSQ, JIQ, and JSW.}
Three schemes that clearly have vanishing wait are JSQ, JIQ and JSW.
The optimality of JSW is observed in the figures;
JSW has the smallest mean waiting time, and all three schemes have
vanishing wait in both the fluid and diffusion regime.
There is a significant difference, however, between JSW and JSQ/JIQ.
We observe that the probability of positive wait does not vanish for JSW,
while it does vanish for JSQ/JIQ.
This implies that the mean of all positive waiting times is an order
larger in JSQ/JIQ compared to JSW.
Intuitively, this is clear since in JSQ/JIQ, when a task is placed
in a queue, it waits for at least one specific other task.
In JSW, which is equivalent to the M/M/N queue, a task that cannot
start service immediately, can start service when one of the $N$~servers
becomes idle.

\paragraph{Random and Round-Robin.}
The mean waiting time does not vanish for Random and Round-Robin
in the fluid regime, as already mentioned in Subsection \ref{random}.
Moreover, the mean waiting time grows without bound in the diffusion regime
for these two schemes.
This is because the system can still be decomposed into single-server queues, and the loads
of the individual M/M/1 and D/M/1 queues tend to~1.

\paragraph{JSQ($d$) policies.}
Three versions of JSQ($d$) are included in the Figure~\ref{differentschemes};
$d(N)=2$, $d(N)=\lfloor\log(N)\rfloor\to \infty$
and $d(N)=N^{2/3}$ for which $\frac{d(N)}{\sqrt{N}\log(N)}\to \infty$.
Note that the graph for $d(N)=\lfloor \log(N) \rfloor$ shows sudden
jumps when $d(N)$ increases by~1.
As can be seen in Figure~\ref{differentschemes}, the choices
for which $d(N)\to\infty$ have vanishing wait in the fluid regime,
while $d = 2$ has not.

Overall, we see that JSQ($d$) policies clearly outperform
Random and Round-Robin dispatching, while JSQ, JIQ, and JSW are better
in terms of mean wait.

\begin{table}\centering
\def\arraystretch{1.5}
\begin{tabular}{|C{3cm}|C{3cm}|C{2.5cm}|C{3cm}|C{1.75cm}|}
\hline
Scheme & Queue length & Waiting time (fixed $\lambda < 1$) &
Waiting time ($1 - \lambda \sim 1 / \sqrt{N}$) & Overhead per task \\
\hline\hline
Random & $q_i^\star = \lambda^i$ &  $\frac{\lambda}{1 - \lambda}$ &
$\Theta(\sqrt{N})$ & 0 \\
\hline
JSQ($d$) & $q_i^\star = \lambda^{\frac{d^i - 1}{d - 1}}$ &
$\Theta$(1) & $\Omega(\log{N})$ & $2 d$ \\
\hline
\mbox{$d(N)$} \mbox{$\to \infty$} & same as JSQ & same as JSQ & ?? &
$2d(N)$ \\
\hline
$\frac{d(N)}{\sqrt{N} \log(N)}\to \infty$ &
same as JSQ & same as JSQ & same as JSQ & $2d(N)$ \\
\hline
JSQ & $q_1^\star = \lambda$, $q_2^\star =$ o(1) & o(1) &
$\Theta(1 / \sqrt{N})$ & $2 N$ \\
\hline\hline
JIQ & same as JSQ & same as JSQ & same as JSQ & $\leq 1$ \\
\hline
\end{tabular}
\caption{Queue length distribution, waiting times, and communication
overhead for various~LBAs.}
\label{table}
\end{table}

\section{Preliminaries, JSQ policy, and power-of-d algorithms}
\label{sec:powerofd}

In this section we first introduce some useful notation and preliminary
concepts, and then review fluid and diffusion limits for the JSQ policy
as well as JSQ($d$) policies with a fixed value of~$d$.

We keep focusing on the basic scenario
where all the servers are homogeneous, the service requirements are
exponentially distributed, and the service discipline at each server
is oblivious of the actual service requirements.
In order to obtain a Markovian state description, it therefore
suffices to only track the number of tasks, and in fact we do not need
to keep record of the number of tasks at each individual server,
but only count the number of servers with a given number of tasks.
Specifically, we represent the state of the system by a vector
\begin{equation}
\label{sos}
\QQ(t) := \left(Q_1(t), Q_2(t), \dots\right)
\end{equation}
with $Q_i(t)$ denoting the number of servers with $i$ or more tasks
at time~$t$, including the possible task in service, $i = 1, 2 \dots$.
Note that if we represent the queues at the various servers as (vertical)
stacks, and arrange these from left to right in non-descending order,
then the value of $Q_i$ corresponds to the width of the $i$-th (horizontal)
row, as depicted in the schematic diagram in Figure~\ref{figB}.

\begin{figure}
\begin{minipage}{\textwidth}
\begin{center}
\hfill\parbox{0.48\textwidth}{\centering
\begin{tikzpicture}[scale=0.5]
\node[above] at (3.750,5) {\small$\lambda(N)$};
\draw[thick,black,->] (2.75,5)--(4.75,5);
\draw[thick,black,fill=black](5,5) circle [radius=0.1];
\draw[thick,black,->,dash pattern=on 2 off 1] (5,5)--(7,5);
\draw[thick,black,->,dash pattern=on 2 off 1] (5,5)--(7,6.5);
\draw[thick,black,->,dash pattern=on 2 off 1] (5,5)--(7,8);
\draw[thick,black,->,dash pattern=on 2 off 1] (5,5)--(7,3.5);
\draw[thick,black,->,dash pattern=on 2 off 1] (5,5)--(7,2);
\foreach \i in {1,...,4}
{
\draw[thick,black,fill=mygray](13-\i,8) circle [radius=0.4];
}
\foreach \i in {1,...,2}
{
\draw[thick,black,fill=mygray](13-\i,6.5) circle [radius=0.4];
}
\foreach \i in {1,...,3}
{
\draw[thick,black,fill=mygray](13-\i,5) circle [radius=0.4];
}
\foreach \i in {1,...,5}
{
\draw[thick,black,fill=mygray](13-\i,2) circle [radius=0.4];
}
\draw[thick,black](13.2,8) circle [radius=0.5];
\draw[thick,black](13.2,6.5) circle [radius=0.5];
\draw[thick,black](13.2,5) circle [radius=0.5];
\draw[thick,black](13.2,2) circle [radius=0.5];
\draw[thick,black,->] (13.7,8)--(14.5,8);
\draw[thick,black,->] (13.7,6.5)--(14.5,6.5);
\draw[thick,black,->] (13.7,5)--(14.5,5);
\draw[thick,black,->] (13.7,2)--(14.5,2);
\node at (13.2,8) {\small 1};
\node at (13.2,6.5) {\small 2};
\node at (13.2,5) {\small 3};
\node at (13.2,3.5) {$\vdots$};
\node at (13.2,2) {\small $N$};

\end{tikzpicture}
}
\parbox{0.48\textwidth}{\centering
\begin{tikzpicture}[scale=.60]
\foreach \x in {10, 9,...,1}
    \draw (\x,6)--(\x,0)--(\x+.7,0)--(\x+.7,6);
\foreach \x in {10, 9,...,1}
    \draw (\x+.35,-.45) node[circle,inner sep=0pt, minimum size=10pt,draw,thick] {{{\tiny $\mathsmaller{\x}$}}} ;
\foreach \y in {0, .5}
    \draw[fill=mygray,mygray] (1.15,.1+\y) rectangle (1.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5}
    \draw[fill=mygray,mygray] (2.15,.1+\y) rectangle (2.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5}
    \draw[fill=mygray,mygray] (3.15,.1+\y) rectangle (3.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5}
    \draw[fill=mygray,mygray] (4.15,.1+\y) rectangle (4.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3}
    \draw[fill=mygray,mygray] (5.15,.1+\y) rectangle (5.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3, 3.5, 4}
    \draw[fill=mygray,mygray] (6.15,.1+\y) rectangle (6.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3, 3.5, 4}
    \draw[fill=mygray,mygray] (7.15,.1+\y) rectangle (7.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5}
    \draw[fill=mygray,mygray] (8.15,.1+\y) rectangle (8.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5}
    \draw[fill=mygray,mygray] (9.15,.1+\y) rectangle (9.55,.5+\y);
\foreach \y in {0, .5, 1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5}
    \draw[fill=mygray,mygray] (10.15,.1+\y) rectangle (10.55,.5+\y);

\draw[thick] (.9,0) rectangle (10.8,.5);
\draw[thick] (.9,.5) rectangle (10.8,1);
\draw[thick] (3.9,2.5) rectangle (10.8,3);

\draw  (12, .2) node {{\scriptsize $\leftarrow Q_1=10$}};
\draw  (12, .9) node {{\scriptsize $\leftarrow Q_2=10$}};

\draw  (11.85, 1.3) node {{\tiny $\cdot$}};
\draw  (11.85, 1.8) node {{\tiny $\cdot$}};
\draw  (11.85, 2.3) node {{\tiny $\cdot$}};
\draw  (12, 2.7) node {{\scriptsize $\leftarrow Q_i=7$}};
\draw  (11.85, 3.3) node {{\tiny $\cdot$}};
\draw  (11.85, 3.8) node {{\tiny $\cdot$}};
\draw  (11.85, 4.3) node {{\tiny $\cdot$}};
\end{tikzpicture}
}\\[1ex]
\parbox{0.48\textwidth}{\centering
\captionof{figure}{Tasks arrive at the dispatcher as a Poisson process
of rate $\lambda(N)$, and are forwarded to one of the $N$ servers
according to some specific load balancing algorithm.
\label{figJSQ}}
}
\hfill\parbox{0.48\textwidth}{\centering
\captionof{figure}{The value of $Q_i$ represents the width of the
$i$-th row, when the servers are arranged in non-descending order
of their queue lengths.
\label{figB}}
}
\end{center}
\end{minipage}
\end{figure}

In order to examine the fluid and diffusion limits in regimes
where the number of servers~$N$ grows large, we consider a sequence
of systems indexed by~$N$, and attach a superscript~$N$ to the
associated state variables.
The fluid-scaled occupancy state is denoted by
$\qq^N(t) := (q_1^N(t), q_2^N(t), \dots)$, with $q_i^N(t) = Q_i^N(t) / N$
representing the fraction of servers in the $N$-th system with $i$
or more tasks as time~$t$, $i = 1, 2, \dots$.
Let
$\cS = \{\qq \in [0, 1]^\infty: q_i \leq q_{i-1} \forall i = 2, 3,\dots\}$
be the set of all possible fluid-scaled states.
Whenever we consider fluid limits, we assume the sequence of initial
states is such that $\qq^N(0) \to \qq^\infty \in \cS$ as $N \to \infty$.

The diffusion-scaled occupancy state is defined as
$\bar{\QQ}^N(t) = (\bar{Q}_1^N(t), \bar{Q}_2^N(t), \dots)$, with
\begin{equation}\label{eq:diffscale}
\bar{Q}_1^N(t) = - \frac{N - Q_1^N(t)}{\sqrt{{N}}}, \qquad
\bar{Q}_i^N(t) = \frac{Q_i^N(t)}{\sqrt{{N}}}, \quad i = 2,3, \dots.
\end{equation}
Note that $-\bar{Q}_1^N(t)$ corresponds to the number of vacant servers,
normalized by $\sqrt{N}$.
The reason why $Q_1^N(t)$ is centered around~$N$ while $Q_i^N(t)$,
$i = 2,3, \dots$, are not, is that for the scalable LBAs we consider the fraction of servers with exactly one task tends to one,
whereas the fraction of servers with two or more tasks tends to zero
as $N\to\infty$.
For convenience, we will assume that each server has an infinite-capacity
buffer, but all the results extend to the finite-buffer case.

\subsection{Fluid limit for JSQ(d) policies}

We first consider the fluid limit for JSQ($d$) policies
with an arbitrary but fixed value of~$d$ as characterized
by Mitzenmacher~\cite{Mitzenmacher01} and Vvedenskaya
{\em et al.}~\cite{VDK96}:

{\em The sequence of processes $\{\qq^N(t)\}_{t \geq 0}$ has a weak
limit $\{\qq(t)\}_{t \geq 0}$ that satisfies the system
of differential equations}
\begin{equation}
\label{fluid:standard}
\frac{\dif q_i(t)}{\dif t} =
\lambda (q_{i-1}^d(t) - q_i^d(t)) - (q_i(t) - q_{i+1}(t)).
\quad i = 1, 2, \dots.
\end{equation}
The fluid-limit equations may be interpreted as follows.
The first term represents the rate of increase in the fraction
of servers with $i$ or more tasks due to arriving tasks that are
assigned to a server with exactly $i - 1$ tasks.
Note that the latter occurs in fluid state $\qq \in \cS$ with probability
$q_{i-1}^d - q_i^d$, i.e., the probability that all $d$~sampled servers
have $i - 1$ or more tasks, but not all of them have $i$ or more tasks.
The second term corresponds to the rate of decrease in the fraction
of servers with $i$ or more tasks due to service completions from servers
with exactly $i$ tasks, and the latter rate is given by $q_i - q_{i+1}$.
While the system in~\eqref{fluid:standard} characterizes the functional law of large numbers (FLLN) behavior of systems under the JSQ($d$) scheme, weak convergence to a certain Ornstein-Ulenbeck process (both in the transient regime and in steady state) was shown in~\cite{Graham05}, establishing a functional central limit theorem (FCLT) result.
Strong approximations for systems under the JSQ($d$) scheme on any finite time interval by the deterministic system in~\eqref{fluid:standard}, a certain infinite-dimensional jump process, and a diffusion approximation were established in~\cite{LN05}.

When the derivatives in~\eqref{fluid:standard} are set equal to zero
for all~$i$, the unique fixed point for any $d \geq 2$ is obtained as
\begin{equation}
\label{eq:fixedpoint1}
q_i^* = \lambda^{\frac{d^i-1}{d-1}}.
\quad i = 1, 2, \dots.
\end{equation}
It can be shown that the fixed point is asymptotically stable in the
sense that $\qq(t) \to \qq^*$ as $t \to \infty$ for any initial fluid
state $\qq^\infty$ with $\sum_{i = 1}^{\infty} q_i^\infty < \infty$.
As mentioned earlier, the fixed point reveals that the stationary queue length distribution at
each individual server exhibits super-exponential decay as $N \to \infty$,
as opposed to exponential decay for a random assignment policy.
It is worth observing that this involves an interchange of the
many-server ($N \to \infty$) and stationary ($t \to \infty$) limits.
The justification is provided by the asymptotic stability of the fixed
point along with a few further technical conditions.

\subsection{Fluid limit for JSQ policy}
\label{ssec:jsqfluid}

We now turn to the fluid limit for the ordinary JSQ policy,
which rather surprisingly was not rigorously established until fairly
recently in~\cite{MBLW16-3}, leveraging martingale functional limit
theorems and time-scale separation arguments~\cite{HK94}.

In order to state the fluid limit starting from an arbitrary
fluid-scaled occupancy state, we first introduce some additional notation.
For any fluid state $\qq \in \cS$,
denote by $m(\qq) = \min\{i: q_{i + 1} < 1\}$ the minimum queue length
among all servers.
Now if $m(\qq)=0$, then define $p_0(m(\qq))=1$ and $p_i(m(\qq))=0$
for all $i=1,2,\ldots$.
Otherwise, in case $m(\qq)>0$, define
\begin{equation}
\label{eq:fluid-gen}
p_i(\qq) =
\begin{cases}
\min\big\{(1 - q_{m(\qq) + 1})/\lambda,1\big\} & \quad\mbox{ for }\quad i=m(\qq)-1, \\
1 - p_{m(\qq) - 1}(\qq) & \quad\mbox{ for }\quad i=m(\qq), \\
0&\quad \mbox{ otherwise.}
\end{cases}
\end{equation}

{\em Any weak limit of the sequence of processes $\{\qq^N(t)\}_{t \geq 0}$
is given by the deterministic system $\{\qq(t)\}_{t \geq 0}$ that satisfies the system of differential equations
\begin{equation}
\label{eq:fluid}
\frac{\dif^+ q_i(t)}{\dif t} =
\lambda p_{i-1}(\qq(t)) - (q_i(t) - q_{i+1}(t)),
\quad i = 1, 2, \dots,
\end{equation}
where $\dif^+/\dif t$ denotes the right-derivative.}
The reason we have used derivative in~\eqref{fluid:standard},
and right-derivative in~\eqref{eq:fluid} is that the limiting
trajectory for the JSQ policy may not be differentiable at all time points.
In fact, one of the major technical challenges in proving the fluid limit
for the JSQ policy is that the drift of the process is not continuous,
which leads to non-smooth limiting trajectories.

The fluid-limit trajectory in~\eqref{eq:fluid}
can be interpreted as follows.
The coefficient $p_i(\qq)$ represents the instantaneous fraction
of incoming tasks assigned to servers with a queue length of exactly~$i$
in the fluid state $\qq \in \mathcal{S}$.
Note that a strictly positive fraction $1 - q_{m(\qq) + 1}$
of the servers have a queue length of exactly~$m(\qq)$.
Clearly the fraction of incoming tasks that get assigned to servers
with a queue length of $m(\qq) + 1$ or larger is zero:
$p_i(\qq) = 0$ for all $i = m(\qq) + 1, \dots$.
Also, tasks at servers with a queue length of exactly~$i$ are completed
at (normalized) rate $q_i - q_{i + 1}$, which is zero for all
$i = 0, \dots, m(\qq) - 1$, and hence the fraction of incoming tasks
that get assigned to servers with a queue length of $m(\qq) - 2$ or less
is zero as well: $p_i(\qq) = 0$ for all $i = 0, \dots, m(\qq) - 2$.
This only leaves the fractions $p_{m(\qq) - 1}(\qq)$
and $p_{m(\qq)}(\qq)$ to be determined.
Now observe that the fraction of servers with a queue length of exactly
$m(\qq) - 1$ is zero.
If $m(\qq)=0$, then clearly the incoming tasks will join an empty queue,
and thus, $p_{m(\qq)}=1$, and $p_i(\qq) = 0$ for all $i \neq m(\qq)$.
Furthermore, if $m(\qq) \geq 1$, since tasks at servers with a queue
length of exactly $m(\qq)$ are completed at (normalized) rate
$1 - q_{m(\qq) + 1} > 0$, incoming tasks can be assigned to servers
with a queue length of exactly $m(\qq) - 1$ at that rate.
We thus need to distinguish between two cases, depending on whether the
normalized arrival rate $\lambda$ is larger than $1 - q_{m(\qq) + 1}$
or not.
If $\lambda < 1 - q_{m(\qq) + 1}$, then all the incoming tasks can be
assigned to a server with a queue length of exactly $m(\qq) - 1$,
so that $p_{m(\qq) - 1}(\qq) = 1$ and $p_{m(\qq)}(\qq) = 0$.
On the other hand, if $\lambda > 1 - q_{m(\qq) + 1}$, then not all
incoming tasks can be assigned to servers with a queue length of
exactly $m(\qq) - 1$ active tasks, and a positive fraction will be
assigned to servers with a queue length of exactly $m(\qq)$:
$p_{m(\qq) - 1}(\qq) = (1 - q_{m(\qq) + 1}) / \lambda$
and $p_{m(\qq)}(\qq) = 1 - p_{m(\qq) - 1}(\qq)$.

The unique fixed point
$\qq^\star = (q_1^\star,q_2^\star,\ldots)$ of the dynamical system
in~\eqref{eq:fluid} is given by
\begin{equation}
\label{eq:fpjsq}
q_i^* = \left\{\begin{array}{ll} \lambda, & i = 1, \\
0, & i = 2, 3,\dots. \end{array} \right.
\end{equation}
Note that the fixed point naturally emerges when $d \to \infty$ in the
fixed point expression~\eqref{eq:fixedpoint1} for fixed~$d$.
However, the process-level results in \cite{Mitzenmacher01,VDK96}
for fixed~$d$ cannot be readily used to handle joint scalings,
and do not yield the entire fluid-scaled sample path for arbitrary
initial states as given by~\eqref{eq:fluid}.

The fixed point in~\eqref{eq:fpjsq}, in conjunction with an interchange
of limits argument, indicates that in stationarity the fraction
of servers with a queue length of two or larger under the JSQ policy
is negligible as $N \to \infty$.

\subsection{Diffusion limit for JSQ policy}\label{ssec:diffjsq}

We next describe the diffusion limit for the JSQ policy in the
Halfin-Whitt heavy-traffic regime \eqref{eq:HW}, as derived
by Eschenfeldt \& Gamarnik~\cite{EG15}.

Recall the centered and diffusion-scaled processes in~\eqref{eq:diffscale}.

{\em For suitable initial conditions, the sequence of processes
$\big\{\bar{\QQ}^N(t)\big\}_{t \geq 0}$ converges weakly to the limit
$\big\{\bar{\QQ}(t)\big\}_{t \geq 0}$,
where
$(\bar{Q}_1(t), \bar{Q}_2(t),\ldots)$
is the unique solution
to the following system of SDEs
\begin{equation}
\label{eq:diffusionjsq}
\begin{split}
\dif\bar{Q}_1(t) &= \sqrt{2}\dif W(t) - \beta\dif t - \bar{Q}_1(t) + \bar{Q}_2(t)-\dif U_1(t), \\
\dif\bar{Q}_2(t) &= \dif U_1(t) -  \bar{Q}_2(t),
\end{split}
\end{equation}
and $\bQ_i(t) = 0$,
for $t \geq 0$, where $W$ is the standard Brownian motion and $U_1$ is
the unique non-decreasing non-negative process
satisfying
$\int_0^\infty \mathbbm{1}_{[\bar{Q}_1(t) < 0]} \dif U_1(t) = 0$.}

The above convergence of the scaled occupancy measure was established in~\cite{EG15} only in the transient regime on any finite time interval.
The tightness of the sequence of diffusion-scaled steady-state occupancy measures $\{(\bQ_1^N(\infty), \bQ_2^N(\infty))\}_{N\geq 1}$, the ergodicity of the limiting diffusion process~\eqref{eq:diffusionjsq}, and hence the interchange of limits were open until Braverman~\cite{Braverman18} recently further established that the weak-convergence result extends to the steady state as well, i.e., $\bar{\QQ}^N(\infty)$ converges weakly to $(\bQ_1(\infty), \bQ_2(\infty), 0, 0,\ldots)$ as $N\to\infty$, where $(\bQ_1(\infty), \bQ_2(\infty))$ has the stationary distribution of the process $(\bQ_1, \bQ_2)$.
 Thus, the steady state of the diffusion process in~\eqref{eq:diffusionjsq} is proved to capture the asymptotic behavior of large-scale systems under the JSQ policy.

In~\cite{Braverman18} a Lyapunov function is obtained via a generator expansion framework using Stein's method, which establishes exponential ergodicity of $(\bQ_1, \bQ_2)$.
Although this approach gives a good handle on the rate of convergence to stationarity, it sheds little light on the form of the stationary distribution of the limiting diffusion process~\eqref{eq:diffusionjsq} itself.
In two companion papers Banerjee \& Mukherjee~\cite{BM18, BM18b} perform a thorough analysis of the steady state of this diffusion process.

Using a classical regenerative process construction of the diffusion process in~\eqref{eq:diffusionjsq}, \cite[Theorem 2.1]{BM18} establishes that $\bQ_1(\infty)$ has a Gaussian tail, and the tail exponent is uniformly bounded by constants which do not depend on $\beta$, whereas $\bQ_2(\infty)$ has an exponentially decaying tail, and the coefficient in the exponent is linear in $\beta$.
More precisely, for any $\beta>0$ there exist positive constants $C_1, C_2, D_1, D_2$ not depending on $\beta$ and positive constants $C^{l}(\beta)$, $C^{u}(\beta)$, $D^{l}(\beta)$, $D^{u}(\beta)$, $C_R(\beta)$, $D_R(\beta)$ depending only on $\beta$ such that
\begin{equation}\label{eq:statail}
\begin{split}
C^{l}(\beta)e^{-C_1x^2} \le \P(\bQ_1(\infty) < -x) \le C^{u}(\beta)e^{-C_2x^2}, \ \ x \ge C_R(\beta)\\
D^{l}(\beta)e^{-D_1\beta y} \le \P(\bQ_2(\infty) > y) \le D^{u}(\beta)e^{-D_2\beta y}, \ \ y \ge D_R(\beta).
\end{split}
\end{equation}
It was further shown in \cite[Theorem 2.2]{BM18} that there exists a positive constant $\mathcal{C^*}$ not depending on $\beta$ such that almost surely along any sample path:
\begin{equation}\label{eq:fluc-width}
\begin{split}
-2\sqrt{2} &\le \liminf_{t \rightarrow \infty} \frac{\bQ_1(t)}{\sqrt{\log t}} \le -1,\\
\frac{1}{\beta} &\le \limsup_{t \rightarrow \infty} \frac{\bQ_2(t)}{\log t} \le \frac{2}{\mathcal{C^*} \beta}.
\end{split}
\end{equation}

Notice that the width of fluctuation of $\bQ_1$ does not depend on the value of $\beta$, whereas that of $\bQ_2$
is linear in $\beta^{-1}$.

Since the $N$-th system is ergodic and its arrival rate is $N-\beta\sqrt{N}$, it is straightforward to see that $\E(\bQ_1^N) = -\beta$ for all $N$, and hence as expected, it can also be derived from the evolution of the limiting diffusion process that $\E(\bQ_1) = -\beta$.
Thus, intuitively, for large enough $\beta$, the system has mostly many idle servers, and thus the number of servers with queue length at least two diminishes.
But the way $\bQ_2(\infty)$ scales as $\beta$ becomes large, is highly non-trivial.
Specifically, it was shown in~\cite{BM18b} that there exists $\beta_0 \ge 1$
and positive constants $C_1, C_2, D_1, D_2$
such that for all $\beta \ge \beta_0$,
\begin{eq}\label{eq:large-beta}
\e^{-C_1\beta^2} \le \mathbb{E}_{\pi}\left(Q_2(\infty)\right) \le \e^{-C_2\beta^2},\\
\P\Big(Q_2(\infty) \ge \e^{-\e^{D_1\beta^2}}\Big) \le \e^{-D_2\beta^2},
\end{eq}
i.e., the steady-state mean is of order $\e^{-C\beta^2}$, but most of the steady-state mass concentrates at a much smaller scale $\e^{-\e^{D_1\beta^2}}$.
This suggests \emph{intermittency} in the behavior of the $\bQ_2$ process, namely, $\bQ_2$ is typically of order $\e^{-\e^{D_1\beta^2}}$, but during rare events when it achieves higher values, it takes a long time to decay.

However, for small enough $\beta$, the behavior is qualitatively different.
Since $\E(\bQ_1) = -\beta$, the system is expected to become more congested as $\beta$ becomes smaller.
As a result, intuitively, $\bQ_2$ should increase.
In this regime as well, $\bQ_2$ exhibits some surprising behavior.
Specifically, it was shown in~\cite{BM18b} that
there exist positive constants $\beta^*, M_1,$ and $M_2$ such that for all $\beta \le \beta^*$
\begin{eq}\label{eq:small-beta}
\frac{M_1}{\beta} \le \mathbb{E}(\bQ_2(\infty)) \le \frac{M_2}{\beta}.
\end{eq}

\paragraph{Comparison with M/M/$N$ queue.}
It is worth mentioning that the M/M/N queue in the Halfin-Whitt heavy-traffic regime has been studied extensively (see~\cite{GG13a, GG13b, HW81, LK11, LK12, FKL14, LMZ17}, and the references therein).
In this case, the centered and scaled total number of tasks in the system $(\bar{S}^N(t)-N)/\sqrt{N}$ converges weakly to a diffusion process $\{\bar{S}(t)\}_{t\geq 0}$~\cite[Theorem 2]{HW81}

with

\begin{equation}\label{eq:diffusionmmn}
\dif \bar{S}(t) = \sqrt{2}\dif W(t) - \beta \dif t - \dif \bar{S}(t)\ind{\bar{S}(t)\leq 0},
\end{equation}
where $W$ is the standard Brownian motion.

As reflected in~\eqref{eq:diffusionjsq} and~\eqref{eq:diffusionmmn}, the JSQ policy and the M/M/$N$ system share some surprising similarities in terms of the qualitative behavior of the total number of tasks in the system.
In particular, both the number of idle servers and the number of waiting tasks are of the order
$\Theta(\sqrt{N})$.
This shows that despite the distributed queueing operation
a suitable load balancing policy can deliver a similar combination
of excellent service quality and high resource utilization efficiency
in the QED regime as in a centralized queueing arrangement.
Moreover,
the interchange of limits result in~\cite{Braverman18} implies that for systems under the JSQ policy, $\bQ_{tot}^N(\infty):= \sum_{i = 1}^{\infty} Q_i^N(\infty)$ converges weakly to $\bQ_1(\infty)+\bQ_2(\infty)$, which has an Exponential upper tail (large positive deviation) and a Gaussian lower tail (large negative deviation), see~\eqref{eq:statail}.
This is again reminiscent of the corresponding tail asymptotics for the M/M/$N$ queue.
Note that since $\bar{S}(\cdot)$ is a simple combination of a Brownian motion with a negative drift (when all servers are fully occupied) and an Ornstein Uhlenbeck (OU) process (when there are idle servers), the steady-state distribution $\bar{S}(\infty)$ can be computed explicitly, and is indeed a combination of an exponential distribution (from the Brownian motion with a negative drift) and a Gaussian distribution (from the OU process).

There are, however, some clear differences between the diffusion in~\eqref{eq:diffusionjsq} and~\eqref{eq:diffusionmmn}:
\begin{enumerate}[{\normalfont (i)}]
\item Observe that in case of M/M/N systems, whenever there are some waiting tasks (equivalent to $Q_2$ being positive in our case), the queue length has a constant negative drift towards zero. This leads to the exponential upper tail of $\bar{S}(\infty)$, by comparing with the stationary distribution of a reflected Brownian motion with constant negative drift. In the JSQ case, however, the rate of decrease of $Q_2$ is always proportional to itself, which makes it somewhat counter-intuitive that its stationary distribution has an exponential tail.
\item In the M/M/N system, the number of idle servers can be non-zero only when the number of waiting tasks is zero.
Thus, the dynamics of both the number of idle servers and the number of waiting tasks are completely captured by the one-dimensional process $\bar{S}^N$ and by the one-dimensional diffusion  $\bar{S}$ in the limit.
But in the JSQ case, $\bQ_2$ is never zero, and the dynamics of $(\bQ_1, \bQ_2)$ are truly two-dimensional (although the diffusion is non-elliptic) with $\bQ_1$ and $\bQ_2$ interacting with each other in an intricate manner.
\item From \eqref{eq:diffusionjsq} we see that $\bQ_2$ {\em never} hits zero.
Thus, in steady state, there is no mass at $\bQ_2=0$, and the system always has waiting tasks.
This is in sharp contrast with the M/M/N case, where the system has no waiting tasks in steady state with positive probability.
\item In the M/M/$N$ system, a positive fraction of the
tasks incur a non-zero waiting time as $N \to \infty$, but a non-zero
waiting time is only of length $1 / (\beta \sqrt{N})$ in expectation.
Moreover, in the JSQ case, it is easy to see that $\bQ_1$ (the limit of the scaled number of idle servers) spends zero time at the origin, i.e., in steady state  the fraction of arriving tasks that find all servers busy vanishes in the large-N limit (in fact, this is of order $1/\sqrt{N}$, see~\cite{Braverman18}).
However, such tasks will have to wait for the duration of a residual
service time (the time till the service of the task ahead of it in its queue finishes), yielding a waiting time of the order $O(1)$.
\item As $\beta\to 0$, \cite[Proposition 2]{HW81} implies that $\beta \bar{S}(\infty)$ for the M/M/$N$ queue converges weakly to a unit-mean exponential distribution.
In contrast, results in~\cite{BM18b} show that $\beta (\bQ_1(\infty)+\bQ_2(\infty))$ converges weakly to a Gamma$(2)$ random variable.
This indicates that despite similar order of performance, due to the distributed operation, in terms of the number of waiting tasks JSQ is a factor 2 worse in expectation than the corresponding centralized system.
\end{enumerate}

\subsection{JSQ(d) policies in heavy-traffic regime}

Finally, we briefly discuss the behavior of JSQ($d$) policies
with a fixed value of~$d$ in the Halfin-Whitt heavy-traffic
regime~\eqref{eq:HW}.
While a complete characterization of the occupancy process
for fixed~$d$ has remained elusive so far, significant partial results
were recently obtained by Eschenfeldt \& Gamarnik~\cite{EG16}.
In order to describe the transient asymptotics,
introduce the following rescaled processes
\begin{equation}
\label{eq:HTtransient}
\bQ_i^N(t) := \frac{N-Q_i^N(t)}{\sqrt{N}}, \quad i = 1, 2,\ldots.
\end{equation}
Note that in contrast to~\eqref{eq:diffscale}, in~\eqref{eq:HTtransient}
\emph{all} components are centered by $N$.
\emph{For suitable initial states, \cite[Theorem~2]{EG16}
establishes that on any finite time interval, $\bar{\QQ}^N(\cdot)$
converges weakly to a deterministic system $\bar{\QQ}(\cdot)$ that
satisfies the system of ODEs
\begin{equation}
\dif \bQ_i(t) = - d (\bQ_i(t) - \bQ_{i-1}(t)) + \bQ_{i+1}(t) - \bQ_i(t),
\quad i = 1, 2, \ldots
\end{equation}
with the convention that $\bQ_0(t) \equiv 0$.}
It is noteworthy that the scaled occupancy process loses its
diffusive behavior for fixed~$d$.
It is further shown in~\cite{EG16} that with high probability the
steady-state fraction of queues with length at least
$\log_d(\sqrt{N} / \beta) - \omega(1)$ tasks approaches unity,
which in turn implies that with high probability the steady-state delay
is {\em at least} $\log_d(\sqrt{N} / \beta) - O(1)$ as $N \to \infty$.
The diffusion approximation of the JSQ($d$) policy in the Halfin-Whitt regime~\eqref{eq:HW}, starting from a different initial scaling, has been studied by Budhiraja \& Friedlander~\cite{BF17}.

In the work of Ying~\cite{Ying17} a broad framework involving Stein's
method was introduced to analyze the rate of convergence of the
stationary distribution in a heavy-traffic regime where
$(N - \lambda(N))/\eta(N) \to \beta > 0$ as $N \to \infty$,
with $\eta(N)$ a positive function diverging to infinity as $N \to \infty$.
Note that the case $\eta(N)= \sqrt{N}$ corresponds to the Halfin-Whitt
heavy-traffic regime~\eqref{eq:HW}.
Using this framework, it was proved that when $\eta(N) = N^\alpha$
with some $\alpha>0.8$,
\begin{equation}
\label{eq:po2ht}
\mathbb{E}\Big(\sum_{i=1}^\infty \Big|q_i^N(\infty) - q_i^\star\Big|\Big)\leq\frac{1}{N^{2\alpha -1-\xi}},\qquad\mbox{where}\qquad q_i^\star = \Big(\frac{\lambda(N)}{N}\Big)^{2^k-1},
\end{equation}
and $\xi>0$ is an arbitrarily small constant.
Equation~\eqref{eq:po2ht} not only shows that asymptotically the
stationary occupancy measure concentrates at $\qq^\star$,
but also provides the rate of convergence.

\section{Universality of JSQ(d) policies}
\label{univ}

In this section we will further explore the trade-off between delay
performance and communication overhead as a function of the diversity
parameter~$d$, in conjunction with the relative load.
The latter trade-off will be examined in an asymptotic regime where
not only the total task arrival rate $\lambda(N)$ grows with~$N$,
but also the diversity parameter depends on~$N$,
and we write $d(N)$ to explicitly reflect that.
We will specifically investigate what growth rate of $d(N)$ is required,
depending on the scaling behavior of $\lambda(N)$,
in order to asymptotically match the optimal performance of the JSQ
policy and achieve a zero mean waiting time in the limit.
The results presented in the remainder of the section are based
on~\cite{MBLW16-3} unless specified otherwise.

\begin{theorem}[Universality fluid limit for JSQ($d(N)$)]
\label{fluidjsqd}
If $d(N)\to\infty$ as $N\to\infty$, then the fluid limit of the
JSQ$(d(N))$ scheme coincides with that of the ordinary JSQ policy,
and in particular, is given by the dynamical system in~\eqref{eq:fluid}.
Consequently, the stationary occupancy states converge to the unique
fixed point as in~\eqref{eq:fpjsq}.
\end{theorem}

\begin{theorem}[Universality diffusion limit for JSQ($d(N)$)]
\label{diffusionjsqd}
If $d(N) /( \sqrt{N} \log N)\to\infty$, then for suitable initial
conditions the weak limit of the sequence of processes
$\big\{\bar{\QQ}^{ d(N)}(t)\big\}_{t \geq 0}$ coincides with that
of the ordinary JSQ policy, and in particular, is given by the system
of SDEs in~\eqref{eq:diffusionjsq}.
\end{theorem}

The above universality properties indicate that the JSQ overhead can
be lowered by almost a factor O($N$) and O($\sqrt{N} / \log N$)
while retaining fluid- and diffusion-level optimality, respectively.
In other words, Theorems~\ref{fluidjsqd} and~\ref{diffusionjsqd} reveal that it is sufficient for $d(N)$ to grow at any rate,
and faster than $\sqrt{N} \log N$, in order to observe similar scaling
benefits as in a pooled system with $N$~parallel single-server queues
on fluid scale and diffusion scale, respectively.
The stated conditions are in fact close to necessary, in the sense
that if $d(N)$ is uniformly bounded and $d(N) /( \sqrt{N} \log N) \to 0$
as $N \to \infty$, then the fluid-limit and diffusion-limit paths
of the system occupancy process under the JSQ($d(N)$) scheme differ
from those under the ordinary JSQ policy.
In particular, if $d(N)$ is uniformly bounded, the mean steady-state
delay does not vanish asymptotically as $N \to \infty$.

\begin{remark}\normalfont
One implication of Theorem~\ref{fluidjsqd} is that in the subcritical regime any growth rate of $d(N)$ is enough to achieve asymptotically vanishing steady-state probability of wait.
This result is complimented by recent results of Liu and Ying~\cite{LY18} and Brightwell {\em et al.}~\cite{BFL18}, where steady-state analysis is extended in the heavy-traffic regime.
Specifically, it is established in~\cite{LY18} that when the system load of the $N$-th system scales as $N-N^{\alpha}$ with $\alpha\in (0,1/2)$ (i.e., the system is in heavy traffic, but the load is lighter than that in the Halfin-Whitt regime), the steady-state probability of wait for the JSQ($d(N)$) policy with $d(N)\geq N^{1-\alpha}\log N$ vanishes as $N\to\infty$.
The results of~\cite{BFL18} imply that when $\lambda(N) = N - N^{\alpha}$ and $d(N) = \lfloor N^{\beta}\rfloor$ with $\alpha,\beta\in (0,1]$, $k = \lceil (1-\alpha)/\beta\rceil$, and $2\alpha+\beta (k-1)>1$, with probability tending to 1 as $N\to\infty$, the proportion of queues with queue length equal to $k$ is at least $1-2N^{-1+\alpha + (k-1)\beta}$ and there are no longer queues.
It is important to note that a crucial difference between the result stated in Theorem~\ref{diffusionjsqd} and the results in~\cite{BFL18, LY18} is that the former analyzes the system on diffusion scale (and describes its behavior in terms of a limiting diffusion process), whereas~\cite{BFL18, LY18} analyze the system on fluid-scale (and characterize its behavior in terms of limiting fluid-scaled occupancy state).
\end{remark}

\paragraph{High-level proof idea.}
The proofs of both Theorems~\ref{fluidjsqd} and~\ref{diffusionjsqd}
rely on a stochastic coupling construction to bound the difference
in the queue length processes between the JSQ policy and a scheme
with an arbitrary value of $d(N)$.
This coupling is then exploited to obtain the fluid and diffusion
limits of the JSQ($d(N)$) policy, along with the associated fixed point,
under the conditions stated in Theorems~\ref{fluidjsqd}
and~\ref{diffusionjsqd}.

A direct comparison between the JSQ$(d(N))$ scheme and the ordinary JSQ policy is not straightforward, which is why the $\CJSQ(n(N))$ class of schemes is introduced as an intermediate scenario to establish the universality result.
Just like the JSQ$(d(N))$ scheme, the schemes in the class $\CJSQ(n(N))$ may be thought of as ``sloppy'' versions of the JSQ policy, in the sense that tasks are not necessarily assigned to a server with the shortest queue length but to one of the $n(N)+1$
lowest ordered servers, as graphically illustrated in Figure~\ref{fig:sfigCJSQ}.
In particular, for $n(N)=0$, the class only includes the ordinary JSQ policy.
Note that the JSQ$(d(N))$ scheme is guaranteed to identify the lowest ordered server, but only among a randomly sampled subset of $d(N)$ servers.
In contrast, a scheme in the $\CJSQ(n(N))$ class only guarantees that
one of the $n(N)+1$ lowest ordered servers is selected, but
across the entire pool of $N$ servers.
We will show that for sufficiently small $n(N)$, any scheme from the class $\CJSQ(n(N))$ is still `close' to the ordinary JSQ policy.
We will further prove that for sufficiently large $d(N)$ relative to $n(N)$ we can construct a scheme
called JSQ$(n(N),d(N))$, belonging to the $\CJSQ(n(N))$ class, which differs `negligibly' from the JSQ$(d(N))$ scheme.
Therefore,  for a `suitable' choice of $d(N)$ the idea is to produce a `suitable' $n(N)$.
This proof strategy is schematically represented in Figure~\ref{fig:sfigRelation}.

\begin{figure}
\begin{center}
\begin{subfigure}{.4\textwidth}
  \centering
  \includegraphics[scale=.55]{sfigCJSQ.pdf}
  \caption{$\CJSQ(n(N))$ scheme\vspace{16pt}}
  \label{fig:sfigCJSQ}
\end{subfigure}
\begin{subfigure}{.6\textwidth}
  \centering
  \includegraphics[scale=1.1]{sfigRelation.pdf}
  \caption{Asymptotic equivalence relations}
  \label{fig:sfigRelation}
\end{subfigure}
\caption{(a) High-level view of the $\CJSQ(n(N))$ class of schemes, where as in Figure~\ref{figB}, the servers are arranged in nondecreasing order of their queue lengths, and the arrival must be assigned through the left tunnel. (b) The equivalence structure is depicted for various intermediate load balancing schemes to facilitate the comparison between the JSQ$(d(N))$ scheme and the ordinary JSQ policy.
}
\label{fig:strategy}
\end{center}
\end{figure}

In order to prove the stochastic comparisons among the various schemes,
the many-server system is described as an ensemble of stacks,
in a way that two different ensembles can be ordered.
This stack formulation has also been considered in the literature
for establishing the stochastic optimality properties of the JSQ
policy \cite{towsley,Towsley95,Towsley1992}.
In Remark~\ref{rem:literatureJSQ} we will compare and contrast the
various stochastic comparison techniques.
In this formulation, at each step, items are added or removed
(corresponding to an arrival or departure) according to some rule.
From a high level, it is then shown that if two systems follow some
specific rules, then at any step, the two ensembles maintain some kind
of deterministic ordering.
This deterministic ordering turns into an almost sure ordering in the
probability space constructed by a specific coupling.
In what follows, each server along with its queue is thought of as
a stack of items, and the stacks are always considered to be arranged
in non-decreasing order of their heights.
The ensemble of stacks then represents the empirical CDF of the queue
length distribution, and the $i$-th horizontal bar
corresponds to $Q_i^\Pi$ (for some task assignment scheme~$\Pi$),
as depicted in Figure~\ref{figB}.
For the sake of full exposure, we will describe the coupling
construction in the scenario when the buffer capacity~$B$ at each
stack can possibly be finite.
If $B < \infty$ and an arriving item happens to land on a stack which
already contains $B$ items, then the item is discarded, and is added
to a special stack $L^\Pi$ of discarded items, where it stays forever.

Any two ensembles $\cA$ and $\cB$, each having $N$~stacks
and a maximum height $B$ per stack, are said to follow
Rule($n_{\cA},n_{\cB},k$) at some step, if either an item is removed
from the $k$-th stack in both ensembles (if nonempty),
or an item is added to the $n_{\cA}$-th stack in ensemble
$\cA$ and to the $n_{\cB}$-th stack in ensemble $\cB$.

\begin{proposition}
\label{prop:det ord}
For any two ensembles of stacks $\cA$ and $\cB$, if at any step \emph{Rule}$(n_{\cA},n_{\cB},k)$ is followed for some value of $n_{\cA}$, $n_{\cB}$, and $k$, with $n_{\cA}\leq n_{\cB}$, then the following ordering is always preserved: for all $m\leq B$,
\begin{equation}\label{eq:det ord}
\sum_{i=m}^B Q_i^{\cA} +L^{\cA}\leq \sum_{i=m}^B Q_i^{\cB} +L^{\cB}.
\end{equation}
\end{proposition}

This proposition says that, while adding the items to the ordered stacks,
if we ensure that in ensemble $\cA$ the item is always placed to the
left of that in ensemble $\cB$, and if the items are removed from the
same ordered stack in both ensembles, then the aggregate size of the
$B-m+1$ highest horizontal bars as depicted in Figure~\ref{figB}
plus the cumulative number of discarded items is no larger in $\cA$
than in $\cB$ throughout.

\paragraph{Another type of sloppiness.}
Recall that the class $\CJSQ(n(N))$ contains all schemes that assign
incoming tasks by some rule to any of the $n(N)+1$ lowest ordered servers.
Let MJSQ$(n(N))$ be a particular scheme that always assigns incoming
tasks to precisely the $(n(N)+1)$-th ordered server.
Notice that this scheme is effectively the JSQ policy when the system
always maintains $n(N)$ idle servers, or equivalently, uses only
$N-n(N)$ servers, and $\MJSQ(n(N)) \in \CJSQ(n(N))$.
For brevity, we will often suppress $n(N)$ in the notation where it is
clear from the context.
We call any two systems \emph{S-coupled}, if they have synchronized
arrival clocks and departure clocks of the $k$-th longest queue,
for $1 \leq k \leq N$ (`S' in the name of the coupling stands for `Server').
Consider three S-coupled systems following respectively the JSQ policy,
any scheme from the class $\CJSQ$, and the $\MJSQ$ scheme.
Recall that $Q_i^\Pi(t)$ is the number of servers with at least $i$~tasks
at time~$t$ and $L^\Pi(t)$ is the total number of lost tasks up to time~$t$,
for the schemes $\Pi=$ JSQ, $\CJSQ$, $\MJSQ$.
The following proposition provides a stochastic ordering for any scheme
in the class CJSQ with respect to the ordinary JSQ policy and the MJSQ scheme.

\begin{proposition}
\label{prop:stoch-ord}
For any fixed $m \geq 1$,
\begin{enumerate}[{\normalfont(i)}]
\item\label{item:jsq-cjsq} $
\left\{\sum_{i=m}^{B} Q_i^{\JSQ}(t) + L^{\JSQ}(t)\right\}_{t \geq 0} \leq_{\mathrm{st}}
\left\{\sum_{i=m}^{B} Q_i^{\CJSQ}(t) + L^{\CJSQ}(t)\right\}_{t \geq 0},
$
\item\label{item:cjsq-mjsq} $\left\{\sum_{i=m}^{B} Q_i^{\CJSQ}(t) +
L^{\CJSQ}(t)\right\}_{t \geq 0}\leq_{\mathrm{st}}
\left\{\sum_{i=m}^{B} Q_i^{\MJSQ}(t) + L^{\MJSQ}(t)\right\}_{t \geq 0},$
\end{enumerate}
provided the inequalities hold at time $t=0$.
\end{proposition}

The above proposition has the following immediate corollary,
which will be used to prove bounds on the fluid and the diffusion scale.

\begin{corollary}
\label{cor:bound}
In the joint probability space constructed by the S-coupling of the three systems under respectively \emph{JSQ}, \emph{MJSQ}, and any scheme from the class \emph{CJSQ}, the following ordering is preserved almost surely throughout the sample path: for any fixed $m \geq 1$
\begin{enumerate}[{\normalfont(i)}]
\item $Q_m^{\CJSQ}(t) \geq \sum_{i=m}^{B} Q_i^{\JSQ}(t) -
\sum_{i=m+1}^{B} Q_i^{\MJSQ}(t) + L^{\JSQ}(t) - L^{\MJSQ}(t)$,
\item $Q_m^{\CJSQ}(t) \leq \sum_{i=m}^{B} Q_i^{\MJSQ}(t) -
\sum_{i=m+1}^{B} Q_i^{\JSQ}(t) + L^{\MJSQ}(t)-L^{\JSQ}(t),$
\end{enumerate}
provided the inequalities hold at time $t=0$.
\end{corollary}

\begin{remark}
\label{rem:literatureJSQ}
\normalfont
Note that $\sum_{i=1}^{B} \min\big\{Q_i, k\big\}$ represents the
aggregate size of the rightmost $k$~stacks, i.e., the $k$~longest queues.
Using this observation, the stochastic majorization property of the
JSQ policy as stated in \cite{towsley,Towsley95,Towsley1992} can be
shown following similar arguments as in the proof
of Proposition~\ref{prop:stoch-ord}.
Conversely, the stochastic ordering between the JSQ policy and the
MJSQ scheme presented in Proposition~\ref{prop:stoch-ord} can also be
derived from the weak majorization arguments developed
in \cite{towsley,Towsley95,Towsley1992}.
But it is only through the stack arguments developed in~\cite{MBLW16-3}
as described above, that the results could be extended to compare
any scheme from the class CJSQ with the scheme MJSQ as well
as in Proposition~\ref{prop:stoch-ord}~\eqref{item:cjsq-mjsq}.
\end{remark}

\paragraph{Comparing two arbitrary schemes.}
To analyze the JSQ$(d(N))$ scheme, we need a further stochastic
comparison argument.
Consider two S-coupled systems following schemes $\Pi_1$ and $\Pi_2$.
Fix a specific arrival epoch, and let the arriving task join the
$n_{\Pi_i}$-th ordered server in the $i$-th
system following scheme $\Pi_i$, $i = 1, 2$ (ties can be broken
arbitrarily in both systems).
We say that at a specific arrival epoch the two systems
\emph{differ in decision}, if $n_{\Pi_1} \neq n_{\Pi_2}$,
and denote by $\Delta_{\Pi_1,\Pi_2}(t)$ the cumulative number of times
the two systems differ in decision up to time~$t$.

\begin{proposition}
\label{prop:stoch-ord2}
For two S-coupled systems under schemes $\Pi_1$ and $\Pi_2$ the
following inequality is preserved almost surely
\begin{equation}
\label{eq:stoch-ord2}
\sum_{i=1}^{B} |Q_i^{\Pi_1}(t) - Q_i^{\Pi_2}(t)| \leq
2 \Delta_{\Pi_1,\Pi_2}(t) \qquad \forall\ t \geq 0,
\end{equation}
provided the two systems start from the same occupancy state at $t=0$,
i.e., $Q_i^{\Pi_1}(0) = Q_i^{\Pi_2}(0)$ for all $i = 1, 2, \ldots, B$.
\end{proposition}

\paragraph{A bridge between two types of sloppiness.}
We will now introduce the JSQ$(n,d)$ scheme with $n, d \leq N$,
which is an intermediate blend between the CJSQ$(n)$ schemes
and the JSQ$(d)$ scheme.

At its first step, just as in the JSQ$(d)$ scheme, the JSQ$(d,n)$ scheme first chooses
the shortest of $d$~random candidates but only sends the arriving task
to that server's queue if it is one of the $n+1$ shortest queues.
If it is not, then at the second step it picks any of the $n+1$ shortest
queues uniformly at random and then sends the task to that server's queue.
Note that by construction, JSQ$(d,n)$ is a scheme
in CJSQ$(n)$.
Consider two S-coupled systems with a JSQ$(d)$ and a JSQ$(n,d)$ scheme.
Assume that at some specific arrival epoch, the incoming task is
dispatched to the $k$-th ordered server in the system
under the JSQ($d$) scheme.
If $k \in \{1, 2, \ldots, n+1\}$, then the system under the JSQ$(n,d)$
scheme also assigns the arriving task to the $k$-th ordered
server.
Otherwise, it dispatches the arriving task uniformly at random among
the first $(n+1)$ ordered servers.

The next proposition provides a bound on the number of times these two
systems differ in decision on any finite time interval.
For any $T \geq 0$, let $A(T)$ and $\Delta(T)$ be the total number
of arrivals to the system and the cumulative number of times that the
JSQ($d$) scheme and JSQ$(n,d)$ scheme differ in decision up to time~$T$.

\begin{proposition}
\label{prop:differ}
For any $T\geq 0$, and $M>0,$
\[
\Pro{\Delta(T) \geq
M \given A(T)} \leq \frac{A(T)}{M} \left(1-\frac{n}{N}\right)^d.
\]
\end{proposition}

\noindent
\emph{Proof-sketch of  Theorem~\ref{fluidjsqd}.}
The proof of Theorem~\ref{fluidjsqd} uses the S-coupling results,
and consists of three main steps:
\begin{enumerate}[{\normalfont (i)}]
\item First it is shown that if $n(N) / N \to 0$ as $N \to \infty$,
then the MJSQ$(n(N))$ scheme has the same fluid limit as the ordinary
JSQ policy.
\item Then the application of Corollary~\ref{cor:bound} proves that
as long as $n(N) / N \to 0$, \emph{any} scheme from the class
$\CJSQ(n(N))$ has the same fluid limit as the ordinary JSQ policy.
\item Next, Propositions~\ref{prop:stoch-ord2} and~\ref{prop:differ}
are used to establish that if $d(N) \to \infty$, then for \emph{some}
$n(N)$ with $n(N) / N \to 0$, the JSQ$(d(N))$ scheme and the
JSQ$(n(N),d(N))$ scheme have the same fluid limit.
The proposition then follows by observing that the JSQ$(n(N),d(N))$
scheme belongs to the class $\CJSQ(n(N))$. \hfill\qed
\end{enumerate}

The proof of Theorem~\ref{diffusionjsqd} follows the same arguments, but
uses the candidate $n(N)/\sqrt{N}\to 0$ (instead of $n(N)/N\to 0$) in Steps (i) and (ii),
and the candidate $d(N)/(\sqrt{N}\log(N))\to\infty$ (instead of $d(N)\to\infty$) in Step (iii).

\paragraph*{Extension to batch arrivals.}
Consider an extension of the model in which tasks arrive in batches.
We assume that the batches arrive as a Poisson process of rate
$\lambda(N) / \ell(N)$, and have fixed size $\ell(N) > 0$,
so that the effective total task arrival rate remains $\lambda(N)$.
We will show that even for arbitrarily slowly growing batch size,
fluid-level optimality can be achieved with $O(1)$ communication
overhead per task.
For that, we define the JSQ($d(N)$) scheme adapted for batch arrivals:
When a batch of size $\ell(N)$ arrives, the dispatcher samples
$d(N) \geq \ell(N)$ servers without replacement, and assigns the
$\ell(N)$ tasks to the $\ell(N)$ servers with the smallest queue
length among the sampled servers.

\begin{theorem}[Batch arrivals]
\label{th:batch}
Consider the batch arrival scenario with growing batch size $\ell(N)\to\infty$ and $\lambda(N)/N\to\lambda<1$ as $N\to\infty$. For the \jsq$(d(N))$ scheme with $d(N)\geq \ell(N)/(1-\lambda-\varepsilon)$ for any fixed $\varepsilon>0$, if $q^{\sss d(N)}_1(0)\to q_1(0)\leq \lambda$, and $q_i^{\sss d(N)}(0)\to 0$ for all $i\geq 2$, then the weak limit of the sequence of processes $\big\{\qq^{\sss d(N)}(t)\big\}_{t\geq 0}$ coincides with that of the ordinary JSQ policy, and in particular, is given by the system in~\eqref{eq:fluid}.
\end{theorem}

Observe that for a fixed $\varepsilon>0$, the communication overhead
per task is on average given by $(1-\lambda-\varepsilon)^{-1}$
which is $O(1)$.
Thus Theorem~\ref{th:batch} ensures that in case of batch arrivals
with growing batch size, fluid-level optimality can be achieved
with $O(1)$ communication overhead per task.
The result for the fluid-level optimality in stationarity can also be
obtained indirectly by exploiting the fluid-limit result in~\cite{YSK15}.
Specifically, it can be deduced from the result in~\cite{YSK15} that
for batch arrivals with growing batch size, the JSQ$(d(N))$ scheme
with suitably growing $d(N)$ yields the same fixed point of the fluid
limit as described in~\eqref{eq:fpjsq}.

\section{Blocking and infinite-server dynamics}
\label{bloc}

The basic scenario that we have focused on so far involved
single-server queues.
In this section we turn attention to a system with parallel server
pools, each with $B$~servers, where $B$ can possibly be infinite.
As before, tasks arrive at a single dispatcher and must immediately
be forwarded to one of the server pools, but also directly start
execution or be discarded otherwise.
The execution times are assumed to be exponentially distributed,
and do not depend on the number of other tasks receiving service
simultaneously, but the experienced performance (e.g.~in terms
of received throughput or packet-level delay) does degrade
in a convex manner with an increasing number of concurrent tasks.
In order to distinguish it from the single-server queueing dynamics
as considered earlier, the current scenario will henceforth be referred
to as the `infinite-server dynamics'.
These characteristics pertain for instance to video streaming sessions
and various interactive applications.
In contrast to elastic data transfers or computing-intensive jobs,
the duration of such sessions is hardly affected by the number
of contending service requests.
The perceived performance in terms of video quality or packet-level
latency however strongly varies with the number of concurrent tasks,
creating an incentive to distribute the incoming tasks across the
various server pools as evenly as possible.

As it turns out, the JSQ policy has similar stochastic optimality
properties as in the case of single-server queues, and in particular
stochastically minimizes the cumulative number of discarded tasks
\cite{STC93,J89,M87,MS91}.
However, the JSQ policy also suffers from a similar scalability issue
due to the excessive communication overhead in large-scale systems,
which can be mitigated through JSQ($d$) policies.
Results of Turner~\cite{T98} and recent papers by Mukhopadhyay
{\em et al.}~\cite{MKMG15,MMG15}, Karthik {\em et al.}~\cite{KMM17},
and Xie {\em et al.}~\cite{XDLS15} indicate that JSQ($d$) policies
provide similar \emph{power-of-choice} gains for loss probabilities.
It may be shown though that the optimal performance of the JSQ policy
cannot be matched for any fixed value of~$d$.

Motivated by these observations, we explore the trade-off between
performance and communication overhead for infinite-server dynamics.
We will demonstrate that the optimal performance of the JSQ policy
can be asymptotically retained while drastically reducing the
communication burden, mirroring the universality properties described
in Section~\ref{univ} for single-server queues.
The results presented in the remainder of the section are extracted
from~\cite{MBLW16-4}, unless indicated otherwise.

\subsection{Fluid limit for JSQ policy}
\label{ssec:jsqfluid-infinite}

As in Subsection~\ref{ssec:jsqfluid}, for any fluid state $\qq \in \cS$,
denote by $m(\qq) = \min\{i: q_{i + 1} < 1\}$ the minimum queue length
among all servers.
Now if $m(\qq)=0$, then define $p_0(m(\qq))=1$ and $p_i(m(\qq))=0$
for all $i=1,2,\ldots$.
Otherwise, in case $m(\qq)>0$, define
\begin{equation}
\label{eq:fluid-prob-infinite}
p_{i}(\qq) =
\begin{cases}
\min\big\{m(\qq)(1 - q_{m(\qq) + 1})/\lambda,1\big\} & \quad \mbox{ for }
\quad i=m(\qq)-1, \\
1 - p_{ m(\qq) - 1}(\qq) & \quad \mbox{ for } \quad i=m(\qq), \\
0 & \quad \mbox{ otherwise.}
\end{cases}
\end{equation}

{\em Any weak limit of the sequence of processes $\{\qq^N(t)\}_{t \geq 0}$
is given by the deterministic system $\{\qq(t)\}_{t \geq 0}$ satisfying
the following of differential equations
\begin{equation}
\label{eq:fluid-infinite}
\frac{\dif^+ q_i(t)}{\dif t} =
\lambda p_{i-1}(\qq(t)) - i (q_i(t) - q_{i+1}(t)),
\quad i = 1, 2, \dots,
\end{equation}
where $\dif^+/\dif t$ denotes the right-derivative.}

Equations~\eqref{eq:fluid-prob-infinite} and \eqref{eq:fluid-infinite}
are to be contrasted with Equations~\eqref{eq:fluid-gen}
and~\eqref{eq:fluid}.
While the form of the evolution equations~\eqref{eq:fluid-infinite}
of the limiting dynamical system remains similar to~\eqref{eq:fluid},
the rate of decrease of $q_i$ is now $i (q_i - q_{i+1})$,
reflecting the infinite-server dynamics.

Let $K := \lfloor \lambda \rfloor$ and $f := \lambda - K$ denote the
integral and fractional parts of~$\lambda$, respectively.
It is easily verified that, assuming $\lambda<B$, the unique fixed point
of the dynamical system in~\eqref{eq:fluid-infinite} is given by
\begin{equation}
\label{eq:fixed-point-infinite}
q_i^\star = \left\{\begin{array}{ll} 1 & i = 1, \dots, K \\
f & i = K + 1 \\
0 & i = K + 2, \dots, B, \end{array} \right.
\end{equation}
and thus $\sum_{i=1}^{B} q_i^\star = \lambda$.
This is consistent with the results in Mukhopadhyay
{\em et al.}~\cite{MKMG15,MMG15} and Xie {\em et al.}~\cite{XDLS15}
for fixed~$d$, where taking $d \to \infty$ yields the same fixed point.
However, the results in \cite{MKMG15,MMG15,XDLS15} for fixed~$d$
cannot be directly used to handle joint scalings, and do not yield the
universality of the entire fluid-scaled sample path for arbitrary
initial states as stated in~\eqref{eq:fluid-infinite}.

The fixed point in~\eqref{eq:fixed-point-infinite}, in conjunction
with an interchange of limits argument, indicates that in stationarity
the fraction of server pools with at least $K+2$ and at most $K-1$
active tasks is negligible as $N \to \infty$.

\subsection{Diffusion limit for JSQ policy}
\label{ssec:jsq-diffusion-infinite}

As it turns out, the diffusion-limit results may be qualitatively
different, depending on whether $f = 0$ or $f > 0$,
and we will distinguish between these two cases accordingly.
Observe that for any assignment scheme, in the absence of overflow events,
the total number of active tasks evolves as the number of jobs
in an M/M/$\infty$ system with arrival rate $\lambda(N)$ and unit
service rate, for which the diffusion limit is well-known~\cite{Robert03}.
For the JSQ policy we can establish, for suitable initial conditions,
that the total number of server pools with $K - 2$ or less and $K + 2$
or more tasks is negligible on the diffusion scale.
If $f > 0$, the number of server pools with $K - 1$ tasks is negligible
as well, and the dynamics of the number of server pools with $K$
or $K + 1$ tasks can then be derived from the known diffusion limit
of the total number of tasks mentioned above.
In contrast, if $f = 0$, the number of server pools with $K - 1$ tasks
is not negligible on the diffusion scale, and the limiting behavior is
qualitatively different, but can still be characterized.

\subsubsection{Diffusion-limit results for non-integral
\texorpdfstring{$\boldsymbol{\lambda}$}{lambda}}

We first consider the case $f > 0$, and define $f(N) := \lambda(N) - K N$.
Based on the above observations,
we define the following centered and scaled processes:
\begin{equation}
\begin{split}
\bar{Q}^N_i(t)&=N-Q^N_i(t)\geq 0\quad \mathrm{for}\quad i\leq K-1,\\
\bar{Q}_{K}^N(t)&:=\frac{N-Q_K^N(t)}{\log (N)}\geq 0,\\
\bar{Q}_{K+1}^N(t)&:=\frac{Q^N_{K+1}(t)-f(N)}{\sqrt{N}}\in\R,\\
\bar{Q}^N_i(t)&:=Q^N_i(t)\geq 0\quad \mathrm{for}\quad i\geq K+2.
\end{split}
\end{equation}

\begin{theorem}[{Diffusion limit for JSQ policy; $f>0$}]
\label{th:diffusion-infinite}
Assume $\bar{Q}^N_i(0)$ converges to $\bar{Q}_i(0)$ in $\R$, $i \geq 1$,
and $\lambda(N)/N \to \lambda > 0$ as $N \to \infty$. Then
\begin{enumerate}[{\normalfont(i)}]
\item $\lim_{N\to\infty}\Pro{\sup_{t\in[0,T]}\bar{Q}_{K-1}^N(t)\leq 1}=1$, and $\big\{\bar{Q}^N_i(t)\big\}_{t\geq 0}$ converges weakly to $\big\{\bar{Q}_i(t)\big\}_{t\geq 0}$, where $\bar{Q}_i(t)\equiv 0$, provided
$\lim_{N\to\infty}\Pro{\bar{Q}_{K-1}^N(0)\leq 1}=1$, and
$\bar{Q}_i^N(0)\pto 0$ for $i\leq K-2$.
\item $\big\{\bar{Q}^N_K(t)\big\}_{t\geq 0}$ is a stochastically bounded sequence of processes.
\item $\big\{\bar{Q}^N_{K+1}(t)\big\}_{t\geq 0}$ converges weakly to $\big\{\bar{Q}_{K+1}(t)\big\}_{t\geq 0}$, where $\bar{Q}_{K+1}(t)$ is given by the Ornstein-Uhlenbeck process satisfying the following stochastic differential equation:
$$d\bar{Q}_{K+1}(t)=-\bar{Q}_{K+1}(t)dt+\sqrt{2\lambda}dW(t),$$
where $W(t)$ is the standard Brownian motion,
provided $\bar{Q}_{K+1}^N(0)$ converges to $\bar{Q}_{K+1}(0)$ in $\mathbb{R}$.
\item For $i\geq K+2$, $\big\{\bar{Q}^N_i(t)\big\}_{t\geq 0}$ converges weakly to $\big\{\bar{Q}_i(t)\big\}_{t\geq 0}$, where $\bar{Q}_i(t)\equiv 0$, provided $\bar{Q}_i^N(0)$ converges to 0 in $\R$.\\
\end{enumerate}
\end{theorem}

Theorem~\ref{th:diffusion-infinite} implies that for suitable initial
states, for large $N$, there will be almost no server pool with $K-2$
or less tasks and $K+2$ or more tasks on any finite time interval.
Also, the number of server pools having fewer than $K$ tasks is
of order $\log(N)$, and there are $fN+O_P(\sqrt{N})$ server pools
with precisely $K + 1$ active tasks.

\paragraph{High-level proof idea.}

Observe that $\sum_{i=1}^{K} (N-Q_i^N(\cdot))$ increases by one at rate
\[
\sum_{i=1}^{K} i (Q_i(t)-Q_{i+1}(t)) =
\sum_{i=1}^{K} (Q_i(t)-Q_{K+1}(t)) \approx K (1 - f) N,
\]
which is when there is a departure from some server pool with at most
$K$~active tasks, and if positive, decreases by one at constant rate
$\lambda(N) = (K + f) N + o(N)$, which is whenever there is an arrival.
Thus, $\sum_{i=1}^{K} (N - Q_i^N(\cdot))$ roughly behaves as
a birth-and-death process with birth rate $K (1 - f) N$ and death rate
$(K + f) N$.
Since $f > 0$, we have $K + f > K (1 - f)$, and on any finite time
interval the maximum of such a birth-and-death process scales as $\log(N)$.

Similar to the argument above, the process
$\sum_{i=1}^{K-1} \bQ_i^N(\cdot)$ increases by one at rate
\begin{align*}
\sum_{i=1}^{K-1} i (Q^N_i(t)-Q^N_{i+1}(t)) &=
\sum_{i=1}^{K-1} Q^N_i(t) - (K-1) Q_K^N(t) \\
&\leq (K-1) (N-Q_K^N(t)) = O(\log(N)),
\end{align*}
which is when there is a departure from some server pool with at most
$K - 1$ active tasks, and if positive, decreases by one at rate
$\lambda(N)$, which is whenever there is an arrival.
Thus, $\sum_{i=1}^{K-1} \bQ_i^N(\cdot)$ roughly behaves as
a birth-and-death process with birth rate $O(\log(N))$ and death rate $O(N)$.
This leads to the asymptotic result for $\sum_{i=1}^{K-1} \bQ_i^N(\cdot)$,
and in particular for $\bQ_{K-1}^N(\cdot)$.
This completes the proof of Parts~(i) and~(ii)
of Theorem~\ref{th:diffusion-infinite}.

Furthermore, since $\lambda < K + 1$, the number of tasks that are
assigned to server pools with at least $K + 1$ tasks converges to zero
in probability and this completes the proof of Part (iv)
of Theorem~\ref{th:diffusion-infinite}.

Finally, all the above combined also means that on any finite time
interval the total number of tasks in the system behaves with high
probability as the total number of jobs in an M/M/$\infty$ system.
Therefore with the help of the following diffusion limit result
for the M/M/$\infty$ system in~\cite[Theorem 6.14]{Robert03},
we conclude the proof of Part (iii) of Theorem~\ref{th:diffusion-infinite}.

\begin{theorem}[{\cite[Theorem 6.14]{Robert03}}]
\label{th:robert-book-mmn}
Let $\big\{Y_\infty^N(t)\big\}_{t\geq 0}$ be the total number of jobs
in an M/M/$\infty$ system with arrival rate $\lambda (N)$
and unit-mean service time.
If $(Y_\infty^N(0) - \lambda(N))/\sqrt{N} \to v \in \R$,
then the process $\big\{\bar{Y}_\infty^N(t)\big\}_{t\geq 0}$, with
$$\bar{Y}^N_{\infty}(t)=\frac{Y^N_\infty(t)-\lambda( N)}{\sqrt{N}},$$
converges weakly to an Ornstein-Uhlenbeck process
$\big\{X(t)\big\}_{t\geq 0}$ described by the SDE
\begin{align*}
X(0) = v, \qquad \dif X(t) &= - X(t) \dif t + \sqrt{2 \lambda} \dif W(t).
\end{align*}
\end{theorem}

\subsubsection{Diffusion-limit results for integral
\texorpdfstring{$\boldsymbol{\lambda}$}{lambda}}

We now turn to the case $f = 0$, and assume that
\begin{equation}
\label{eq:f=0}
\frac{K N - \lambda(N)}{\sqrt{N}} \to \beta \in \R \quad \mbox{ as }
\quad N \to \infty,
\end{equation}
which can be thought of as an analog of the Halfin-Whitt regime in~\eqref{eq:HW}.

We now consider the following scaled quantities:
\begin{equation}
\begin{split}
\zeta_1^N(t) := \frac{1}{\sqrt{N}}\sum_{i=1}^{K} (N-Q_i^N(t)), \qquad
\zeta_2^N(t) := \frac{Q_{K+1}^N(t)}{\sqrt{N}}.
\end{split}
\end{equation}

\begin{theorem}
\label{th: f=0 diffusion}
Assuming the convergence of initial states, the process
$\big\{(\zeta_1^N(t),\zeta_2^N(t))\big\}_{t \geq 0}$ converges weakly
to the process $\big\{(\zeta_1(t),\zeta_2(t))\big\}_{t \geq 0}$
governed by the system of SDEs
\begin{align*}
\dif\zeta_1(t) &= \sqrt{2K} \dif W(t) - (\zeta_1(t) + K \zeta_2(t)) +
\beta \dif t + \dif V_1(t) \\
\dif\zeta_2(t) &= \dif V_1(t) - (K + 1) \zeta_2(t),
\end{align*}
where $W$ is the standard Brownian motion, and $V_1(t)$ is the unique
non-decreasing process satisfying
\begin{align*}
\int_0^t \ind{\zeta_1(s) \geq 0}\dif V_1(s)=0.
\end{align*}
\end{theorem}

Unlike the $f > 0$ case, the above theorem says that if $f = 0$,
then over any finite time horizon, there will be $O_P(\sqrt{N})$
server pools with fewer than $K$ or more than $K$~active tasks,
and hence most of the server pools have precisely $K$~active tasks.
The proof of Theorem~\ref{th: f=0 diffusion} uses the reflection
argument developed in~\cite{EG15}.

\begin{remark}
\normalfont
Let $Y^N(t)$ denote the total number of tasks in the system at time~$t$.
Note that $Y^N(t) - K N = Z_2^N(t) - Z_1^N(t)$.
Thus, under the assumption in~\eqref{eq:f=0}, the diffusion limit
in Theorem~\ref{th: f=0 diffusion} implies that
\begin{align*}
\frac{Y^N(\cdot) - \lambda(N)}{\sqrt{N}} =
\frac{Y^N(\cdot) - K N}{\sqrt{N}} + \frac{K N - \lambda(N)}{\sqrt{N}}
\end{align*}
converges weakly to $\zeta_2(\cdot) - \zeta_1(\cdot) + \beta$.
Writing $X(t) = \zeta_2(t) - \zeta_1(t) - \beta$,
from Theorem~\ref{th: f=0 diffusion}, one can note that the process
$\big\{X(t)\big\}_{t \geq 0}$ satisfies
$$\dif X(t) = - X(t) \dif t - \sqrt{2 K} \dif W(t),$$
which is consistent with the diffusion-level behavior of $Y^N(\cdot)$
stated in Theorem~\ref{th:robert-book-mmn}.
\end{remark}

\subsection{Universality of JSQ(d) policies in infinite-server dynamics}
\label{ssec:univ-infinite}

As in Section~\ref{univ}, we will now further explore the trade-off
between performance and communication overhead as a function
of the diversity parameter~$d(N)$, in conjunction with the relative load.
We will specifically investigate what growth rate of $d(N)$ is required,
depending on the scaling behavior of $\lambda(N)$, in order to
asymptotically match the optimal performance of the JSQ policy.

\begin{theorem}[Universality fluid limit for JSQ($d(N)$) and infinite-server dynamics]
\label{fluidjsqd-infinite}
If $d(N)\to\infty$ as $N\to\infty$, then the fluid limit of the
JSQ$(d(N))$ scheme coincides with that of the ordinary JSQ policy,
and in particular, is given by the dynamical system
in~\eqref{eq:fluid-infinite}.
Consequently, the stationary occupancy states converge to the unique
fixed point as in~\eqref{eq:fixed-point-infinite}.
\end{theorem}

In order to state the universality result on diffusion scale,
define in case $f > 0$,
\begin{equation}
\begin{split}
\bar{Q}_i^{d(N)}(t) & := \dfrac{N - Q_i^{d(N)}(t)}{\sqrt{N}} \geq 0,
\quad i \leq K, \\
\bar{Q}_{K+1}^{d(N)}(t) & := \dfrac{Q_{K+1}^{d(N)}(t) - f(N)}{\sqrt{N}}
\in \R, \\
\bar{Q}_i^{d(N)}(t) & := \frac{Q_i^{d(N)}(t)}{\sqrt{N}}\geq 0,
\quad \text{ for } \quad i \geq K + 2,
\end{split}
\end{equation}
and otherwise, if $f = 0$,
\begin{equation}

\begin{split}
\hQ_{K-1}^{d(N)}(t) & :=
\sum_{i=1}^{K-1} \dfrac{N - Q_i^{d(N)}(t)}{\sqrt{N}} \geq 0, \\
\hQ_K^{d(N)}(t) & := \dfrac{N - Q_K^{d(N)}(t)}{\sqrt{N}} \geq 0, \\
\hQ_i^{d(N)}(t) & := \dfrac{Q_i^{d(N)}(t)}{\sqrt{N}} \geq 0, \quad
\text{ for } \quad i \geq K + 1.
\end{split}
\end{equation}

\begin{theorem}[Universality diffusion limit for JSQ($d(N)$) and infinite-server dynamics]
\label{diffusionjsqd-infinite}
Assume $d(N) / (\sqrt{N} \log N) \to \infty$.
Under suitable initial conditions
\begin{enumerate}[{\normalfont (i)}]
\item If $f>0$, then $\bQ_i^{d(N)}(\cdot)$ converges to the zero process for $i\neq K+1$, and $\bQ^{d(N)}_{K+1}(\cdot)$ converges weakly to the Ornstein-Uhlenbeck process satisfying the SDE
$$d\bar{Q}_{K+1}(t)=-\bar{Q}_{K+1}(t)dt+\sqrt{2\lambda}dW(t),$$
where $W(t)$ is the standard Brownian motion.
\item If $f=0$, then $\hQ_{K-1}^{d(N)}(\cdot)$ converges weakly to the zero process, and $(\hQ_{K}^{d(N)}(\cdot), \hQ_{K+1}^{d(N)}(\cdot))$ converges weakly to $(\hQ_{K}(\cdot), \hQ_{K+1}(\cdot))$, described by the unique solution of the system of SDEs
\begin{align*}
\dif\hQ_{K}(t) &= \sqrt{2 K} \dif W(t) - (\hQ_K(t) + K \hQ_{K+1}(t)) +
\beta \dif t + \dif V_1(t) \\
\dif\hQ_{K+1}(t) &= \dif V_1(t) - (K + 1) \hQ_{K+1}(t),
\end{align*}
where $W$ is the standard Brownian motion, and $V_1(t)$ is the unique
non-decreasing process satisfying
\begin{align*}
\int_0^t \ind{\hQ_K(s)\geq 0} \dif V_1(s) = 0.
\end{align*}
\end{enumerate}
\end{theorem}

Having established the asymptotic results for the JSQ policy
in Subsections~\ref{ssec:jsqfluid-infinite}
and~\ref{ssec:jsq-diffusion-infinite}, the proofs of the asymptotic
results for the JSQ$(d(N))$ scheme in Theorems~\ref{fluidjsqd-infinite}
and~\ref{diffusionjsqd-infinite} involve establishing a universality
result which shows that the limiting processes for the JSQ$(d(N))$
scheme are `asymptotically equivalent' to those for the ordinary JSQ
policy for suitably large~$d(N)$.
The notion of asymptotic equivalence between different schemes is now
formalized in the next definition.
\begin{definition}
Let $\Pi_1$ and $\Pi_2$ be two schemes parameterized by the number of server pools $N$. For any positive function $g:\N\to\R_+$, we say that $\Pi_1$ and $\Pi_2$ are `$g(N)$-alike' if there exists a common probability space, such that for any fixed $T\geq 0$, for all $i\geq 1$,
$$\sup_{t\in[0,T]}(g(N))^{-1}|Q_i^{\Pi_1}(t)-Q_i^{\Pi_2}(t)|\pto 0\quad \mathrm{as}\quad N\to\infty.$$
\end{definition}

\noindent
Intuitively speaking, if two schemes are $g(N)$-alike, then in some sense,
the associated system occupancy states are indistinguishable
on $g(N)$-scale.
For brevity, for two schemes $\Pi_1$ and $\Pi_2$ that are $g(N)$-alike,
we will often say that $\Pi_1$ and $\Pi_2$ have the same process-level
limits on $g(N)$-scale.
The next theorem states a sufficient criterion for the JSQ$(d(N))$
scheme and the ordinary JSQ policy to be $g(N)$-alike, and thus,
provides the key vehicle in establishing the universality result.

\begin{theorem}
\label{th:pwr of d}
Let $g: \N \to \R_+$ be a function diverging to infinity.
Then the JSQ policy and the JSQ$(d(N))$ scheme are $g(N)$-alike,
with $g(N) \leq N$, if
\begin{align}
\label{eq:fNalike cond1}
\mathrm{(i)}&\quad d(N) \to \infty, \quad \text{for} \quad g(N) = O(N), \\
\mathrm{(ii)}&\quad d(N) \left(\frac{N}{g(N)}\log\left(\frac{N}{g(N)}\right)\right)^{-1} \to \infty, \quad \text{for} \quad g(N)=o(N).
\label{eq:fNalike cond2}
\end{align}
\end{theorem}

Theorem~\ref{th:pwr of d} yields the next two immediate corollaries.
\begin{corollary}
\label{cor:fluid}
If $d(N) \to \infty$ as $N \to \infty$, then the JSQ$(d(N))$ scheme
and the ordinary JSQ policy are $N$-alike.
\end{corollary}

\begin{corollary}
\label{cor-diff}
If $d(N)/(\sqrt{N} \log(N)) \to \infty$ as $N \to \infty$, then the
JSQ$(d(N))$ scheme and the ordinary JSQ policy are $\sqrt{N}$-alike.
\end{corollary}

Observe that Corollaries~\ref{cor:fluid} and~\ref{cor-diff} together
with the asymptotic results for the JSQ policy
in Subsections~\ref{ssec:jsqfluid-infinite}
and~\ref{ssec:jsq-diffusion-infinite} imply
Theorems~\ref{fluidjsqd-infinite} and~\ref{diffusionjsqd-infinite}.
The rest of the section will be devoted to the proof
of Theorem~\ref{th:pwr of d}.
The proof crucially relies on a novel coupling construction,
which will be used to (lower and upper) bound the difference
of occupancy states of two arbitrary schemes.

\paragraph{The coupling construction.}
Throughout the description of the coupling, we fix $N$, and suppress
the superscript $N$ in the notation.
Let $Q_i^{\Pi_1}(t)$ and $Q_i^{\Pi_2}(t)$ denote the number of server
pools with at least $i$~active tasks at time~$t$ in two systems
following schemes $\Pi_1$ and $\Pi_2$, respectively.
With a slight abuse of terminology, we occasionally use $\Pi_1$
and $\Pi_2$ to refer to systems following schemes $\Pi_1$ and $\Pi_2$,
respectively.
To couple the two systems, we synchronize the arrival epochs
and maintain a single exponential departure clock with instantaneous rate
at time $t$ given by $M(t) := \max\left\{\sum_{i=1}^{B} Q_i^{\Pi_1}(t),
\sum_{i=1}^{B} Q_i^{\Pi_2}(t)\right\}$.
We couple the arrivals and departures in the various server pools
as follows:

\noindent
\textit{Arrival:}
At each arrival epoch, assign the incoming task in each system to one
of the server pools according to the respective schemes.

\noindent
\textit{Departure:}
Define $$H(t) :=
\sum_{i=1}^{B} \min\left\{Q_i^{\Pi_1}(t), Q_i^{\Pi_2}(t)\right\}$$
and $$p(t):=
\begin{cases}
\dfrac{H(t)}{M(t)}, & \quad \text{if} \quad M(t) > 0, \\
0, & \quad \text{otherwise.}
\end{cases}
$$
At each departure epoch $t_k$ (say), draw a uniform$[0,1]$ random
variable $U(t_k)$.
The departures occur in a coupled way based upon the value of $U(t_k)$.
In either of the systems, assign a task index $(i,j)$, if that task is
at the $j$-th position of the $i$-th ordered
server pool.
Let $\mathcal{A}_1(t)$ and $\mathcal{A}_2(t)$ denote the set of all
task indices present at time~$t$ in systems $\Pi_1$ and $\Pi_2$,
respectively.
Color the indices (or tasks) in $\mathcal{A}_1 \cap \mathcal{A}_2$,
$\mathcal{A}_1 \setminus \mathcal{A}_2$
and $\mathcal{A}_2 \setminus \mathcal{A}_1$, green, blue and red,
respectively, and note that $|\mathcal{A}_1 \cap \mathcal{A}_2| = H(t)$.
Define a total order on the set of indices as follows:
$(i_1,j_1)<(i_2,j_2)$ if $i_1<i_2$, or $i_1=i_2$ and $j_1<j_2$.
Now, if $U(t_k)\leq p(t_k-)$, then select one green index uniformly
at random and remove the corresponding tasks from both systems.
Otherwise, if $U(t_k)> p(t_k-)$, then choose one integer~$m$,
uniformly at random from all the integers between~$1$
and $M(t) - H(t) = M(t) (1 - p(t))$, and remove the tasks corresponding
to the $m$-th smallest (according to the order defined above)
red and blue indices in the corresponding systems.
If the number of red (or blue) tasks is less than~$m$, then do nothing.

\begin{figure}
\begin{subfigure}{.4\textwidth}
\centering
\includegraphics[scale=.8]{superposition.pdf}
\caption{T-coupling}
\label{fig:tcoupling}
\end{subfigure}
\begin{subfigure}{.6\textwidth}
\begin{center}
  \begin{tikzpicture}[scale=.6]

  \draw(14,1) node[mybox, text width = 2.5cm,text centered]  {
   {\scriptsize JSQ$(n(N),d(N))$}};

  \draw(22,1) node[mybox, text width = 2.5cm,text centered]  {
   {\scriptsize  $\CJSQ(n(N))$}};

   \draw(14,6) node[mybox, text width = 2.5cm,text centered]  {
   {\scriptsize  JSQ$(d(N))$}};

   \draw(22,6) node[mybox, text width = 2.5cm,text centered]  {
   {\scriptsize  JSQ}};

\draw[doublearr] (16.5,6) to (19.5,6);
\draw[doublearr2] (14,2) to (14,5);
\draw[doublearr2] (16.5,1) to (19.5,1);
\draw[doublearr2] (22,2) to (22,5);

\node  at (18,6.5) {\tiny Theorem~\ref{th:pwr of d}};

\node [rotate=270] at (21.5,3.5) {\tiny Proposition~\ref{prop: modified JSL}};
\node [rotate=270] at (22.5,3.5) {\tiny Suitable $n(N)$};

\node [rotate=90] at (14.5,3.5) {\tiny Proposition~\ref{prop: power of d}};
\node [rotate=90] at (13.5,3.5) {\tiny Suitable $d(N)$};

\draw[text width=1.2cm, align=center] (18,1.75) node {\tiny Belongs to};
\draw[text width=1.1cm, align=center] (18,1.35) node {\tiny the class};
\end{tikzpicture}
\caption{Asymptotic equivalence relations}
\label{fig:sfigrelation-inf}
\end{center}
\end{subfigure}
\caption{(a) Superposition of the occupancy states at some particular time instant, of  schemes $\Pi_1$ and $\Pi_2$ when the server pools in both systems are arranged in nondecreasing order of the number of active tasks. The $\Pi_1$ system is the union of the green and blue tasks, and the $\Pi_2$ system is the union of the green and red tasks.
(b) The equivalence structure is depicted for various intermediate load balancing schemes to facilitate the comparison between the JSQ$(d(N))$ scheme and the ordinary JSQ policy.}
\end{figure}

The above coupling has been schematically represented
in Figure~\ref{fig:tcoupling}, and will henceforth be referred to as
T-coupling, where T stands for `task-based'.
Now we need to show that, under the T-coupling, the two systems,
considered independently, evolve according to their own statistical laws.
This can be seen in several steps.
Indeed, the T-coupling basically uniformizes the departure rate
by the maximum number of tasks present in either of the two systems.
Then informally speaking, the green region signifies the common portion of
tasks, and the red and blue regions represent the separate contributions.
Observe that
\begin{enumerate}[{\normalfont (i)}]
\item The total departure rate from $\Pi_i$ is
\begin{align*}
& M(t) \left[p(t)+(1-p(t)) \frac{|\mathcal{A}_i\setminus\mathcal{A}_{3-i}|}{M(t)-H(t)}\right]

= |\mathcal{A}_1\cap\mathcal{A}_2|+|\mathcal{A}_i\setminus\mathcal{A}_{3-i}|

=|\mathcal{A}_i|, \quad i= 1,2.
\end{align*}
\item Assuming without loss of generality
$|\mathcal{A}_1| \geq |\mathcal{A}_2|$, each task in $\Pi_1$ is equally
likely to depart.
\item Each task in $\Pi_2$ within $\mathcal{A}_1\cap \mathcal{A}_2$
and each task within $\mathcal{A}_2 \setminus \mathcal{A}_1$ is equally
likely to depart, and the probabilities of departures are proportional
to $|\mathcal{A}_1 \cap \mathcal{A}_2|$
and $|\mathcal{A}_2 \setminus \mathcal{A}_1|$, respectively.
\end{enumerate}

The T-coupling can be used to derive several stochastic inequality results
that will play an instrumental role in proving Theorem~\ref{th:pwr of d}.
In order to compare the JSQ policy with the CJSQ($n(N)$) schemes,
denote by $Q_i^{\Pi_1}(t)$ and $Q_i^{\Pi_2}(t)$ the number of server
pools with at least $i$~tasks under the JSQ policy and $\CJSQ(n(N))$
scheme, respectively.

\begin{lemma}
\label{lem:majorization}
For any $k \in \big\{1, 2, \ldots, B\big\}$,
\begin{equation}
\label{eq:majorization}
\left\{\sum_{i=1}^{k} Q_i^{\Pi_1}(t) - k n(N)\right\}_{t \geq 0} \leq_{st}
\left\{\sum_{i=1}^{k} Q_i^{\Pi_2}(t)\right\}_{t \geq 0} \leq_{st}
\left\{\sum_{i=1}^{k} Q_i^{\Pi_1}(t)\right\}_{t \geq 0},
\end{equation}
provided at $t=0$ the two systems start from the same occupancy states.
\end{lemma}

\begin{remark}
\label{rem:novelty}
{\normalfont
The stochastic ordering in Lemma~\ref{lem:majorization} is to be
contrasted with the weak majorization results
in~\cite{Winston77,towsley,Towsley95,Towsley1992,W78} in the context
of the ordinary JSQ policy in the single-server queueing scenario,
and in~\cite{STC93,J89,M87,MS91} in the scenario of state-dependent
service rates, non-decreasing with the number of active tasks.
In the current infinite-server scenario, the results
in~\cite{STC93,J89,M87,MS91} imply that for any non-anticipating
scheme~$\Pi$ taking assignment decisions based on the number of active
tasks only, for all $t \geq 0$,
\begin{align}
\label{eq: towsley}
\sum_{m=1}^{\ell} X_{(m)}^{\rm JSQ}(t) & \leq_{st}
\sum_{m=1}^\ell X_{(m)}^{\Pi}(t), \mbox{ for } \ell = 1, 2, \ldots, N, \\
\left\{L^{\rm JSQ}(t)\right\}_{t \geq 0} &\leq_{st}
\left\{L^{\Pi}(t)\right\}_{t \geq 0},
\end{align}
where $X_{(m)}^\Pi(t)$ is the number of tasks in the $m$-th
ordered server pool at time~$t$ in the system following scheme~$\Pi$
and $L^{\Pi}(t)$ is the total number of overflow events under policy~$\Pi$
up to time~$t$.
Observe that $X_{(m)}^\Pi$ can be visualized as the $m$-th
largest (rightmost) vertical bar (or stack) in Figure~\ref{figB}.
Thus~\eqref{eq: towsley} says that the sum of the lengths of the
$\ell$~largest \emph{vertical} stacks in a system following any
scheme~$\Pi$ is stochastically larger than or equal to that following
the ordinary JSQ policy for any $\ell = 1, 2, \ldots, N$.
Mathematically, this ordering can be equivalently written as
\begin{equation}
\label{eq:equiv-ord}
\sum_{i = 1}^{B} \min\big\{\ell, Q_i^{\rm JSQ}(t)\big\} \leq_{st}
\sum_{i = 1}^{B} \min\big\{\ell, Q_i^{\Pi}(t)\big\},
\end{equation}
for all $\ell = 1, \dots, N$.
In contrast, in order to show asymptotic equivalence on various scales,
we need to both upper and lower bound the occupancy states of the
$\CJSQ(n(N))$ schemes in terms of the JSQ policy, and therefore need
a much stronger hold on the departure process.
The T-coupling provides us just that, and has several useful properties
that are crucial for our proof technique.
For example, Proposition~\ref{prop:stoch-ord2} uses the fact that
if two systems are T-coupled, then departures cannot increase the sum
of the absolute differences of the $Q_i$-values, which is not true
for the coupling considered in the above-mentioned literature.
The left stochastic ordering in~\eqref{eq:majorization} also does not
remain valid in those cases.
Furthermore, observe that the right inequality in~\eqref{eq:majorization}
(i.e., $Q_i$'s) implies the stochastic inequality is \emph{reversed}
in~\eqref{eq:equiv-ord}, which is counter-intuitive in view of the
well-established optimality properties of the ordinary JSQ policy.
The fundamental distinction between the two coupling techniques is
also reflected by the fact that the T-coupling does not allow
for arbitrary nondecreasing state-dependent departure rate functions,
unlike the couplings in~\cite{STC93,J89,M87,MS91}.}
\end{remark}

\begin{proposition}
\label{prop: modified JSL}
For any function $g:\N\to\R_+$ diverging to infinity,
if $n(N)/ g(N)\to 0$ as $N \to \infty$, then the JSQ policy and the
$\CJSQ(n(N))$ schemes are $g(N)$-alike.
\end{proposition}

\begin{proof}[Proof of Proposition~\ref{prop: modified JSL}]
Using Lemma~\ref{lem:majorization}, there exists a common probability
space such that for any $k \geq 1$ we can write
\begin{equation}\label{eq:upperocc}
\begin{split}
Q_k^{\Pi_2}(t) &= \sum_{i=1}^{k} Q_i^{\Pi_2}(t) - \sum_{i=1}^{k-1} Q_i^{\Pi_2}(t) \\
&\leq \sum_{i=1}^{k} Q_i^{\Pi_1}(t) - \sum_{i=1}^{k-1} Q_i^{\Pi_1}(t) + k n(N) \\
&= Q_k^{\Pi_1}(t) + k n(N).
\end{split}
\end{equation}
Similarly, we can write
\begin{equation}
\label{eq:lowerocc}
\begin{split}
Q_k^{\Pi_2}(t) &= \sum_{i=1}^{k} Q_i^{\Pi_2}(t) - \sum_{i=1}^{k-1} Q_i^{\Pi_2}(t) \\
&\geq \sum_{i=1}^{k} Q_i^{\Pi_1}(t) - k n(N) - \sum_{i=1}^{k-1} Q_i^{\Pi_1}(t) \\
&= Q_k^{\Pi_1}(t) - k n(N).
\end{split}
\end{equation}
Therefore, for all $k \geq 1$, we have
$\sup_t |Q_k^{\Pi_2}(t) - Q_k^{\Pi_1}(t)| \leq k n(N)$.
Since $n(N)/ g(N) \to \infty$ as $N \to \infty$, the proof is complete.
\end{proof}

\begin{remark}[Comparison of T-coupling and S-coupling]
\label{rem:contrast}
\normalfont
As briefly mentioned earlier, in the current infinite-server scenario,
the departures of the ordered server pools cannot be coupled,
mainly since the departure rate at the $m^{\rm th}$ ordered server pool,
for some $m = 1, 2, \ldots, N$, depends on its number of active tasks.
It is worthwhile to mention that the T-coupling in the current section is
stronger than the S-coupling used in Section~\ref{univ} in the single-server queueing scenario.
Observe that due to Lemma~\ref{lem:majorization}, the absolute
difference of the occupancy states of the JSQ policy and any scheme
from the CJSQ class at any time point can be bounded deterministically
(without any terms involving the cumulative number of lost tasks).
It is worth emphasizing that the universality result on some specific
scale, stated in Theorem~\ref{th:pwr of d}, does not depend on the
behavior of the JSQ policy on that scale, whereas in the single-server queueing scenario
it does, mainly because the upper and lower bounds
in Corollary~\ref{cor:bound} involve tail sums of two different policies.
More specifically, in the single-server queueing scenario the fluid and diffusion limit results of CJSQ($n(N)$) class crucially use those of the MJSQ($n(N)$) scheme, while in the current scenario it does not -- the results for the MJSQ($n(N)$) scheme comes as a consequence of those for the CJSQ($n(N)$) class of schemes.
Also, the bounds in Lemma~\ref{lem:majorization}
do not depend on~$t$, and hence, apply in the steady state as well.
Moreover, the S-coupling compares the $k$ \emph{highest} horizontal bars, whereas the T-coupling
in the current section compares the $k$ \emph{lowest}
horizontal bars.
As a result, the bounds on the occupancy states established
in Corollary~\ref{cor:bound}
involve tail sums of the occupancy states of the ordinary JSQ policy,
which necessitates proving the convergence of tail sums of the
occupancy states of the ordinary JSQ policy.
In contrast, as we will see in the proof of Proposition~\ref{prop: modified JSL}, the bound in the infinite-server scenario involves only a single component
(see Equations \eqref{eq:upperocc} and \eqref{eq:lowerocc}),
and thus, proving convergence of each component suffices.
\end{remark}

Next we compare the CJSQ($n(N)$) schemes with the JSQ($d(N)$) scheme.
The comparison follows a somewhat similar line of argument as
in Section~\ref{univ}, and involves a JSQ$(n(N),d(N))$ scheme
which is an intermediate blend between the $\CJSQ(n(N))$ schemes
and the JSQ$(d(N))$ scheme.
Specifically, the JSQ$(n(N),d(N))$ scheme selects a candidate server
pool in the exact same way as the JSQ$(d(N))$ scheme.
However, it only assigns the task to that server pool if it belongs
to the $n(N)+1$ lowest ordered ones,
and to a randomly selected server pool among these otherwise.
By construction, the JSQ$(n(N),d(N))$ scheme belongs to the class
$\CJSQ(n(N))$.

We now consider two T-coupled systems with a JSQ$(d(N))$ scheme
and a JSQ$(n(N),d(N))$ scheme.
Assume that at some specific arrival epoch, the incoming task is
assigned to the $k$-th ordered server pool in the system
under the JSQ($d(N)$) scheme.
If $k \in \big\{1, 2, \ldots, n(N)+1\big\}$, then the scheme
JSQ$(n(N),d(N))$ also assigns the arriving task to the $k$-th
ordered server pool.
Otherwise it dispatches the arriving task uniformly at random
among the first $n(N)+1$ ordered server pools.

We will establish a sufficient criterion on $d(N)$ in order for the
JSQ$(d(N))$ scheme and JSQ$(n(N),d(N))$ scheme to be close in terms
of $g(N)$-alikeness, as stated in the next proposition.

\begin{proposition}
\label{prop: power of d}
Assume, $n(N) / g(N)\to 0$ as $N \to \infty$ for some function
$g: \N\to \R_+$ diverging to infinity.
The JSQ$(n(N),d(N))$ scheme and the JSQ($d(N)$) scheme are $g(N)$-alike
if the following condition holds:
\begin{equation}
\label{eq:condition-same}
\frac{n(N)}{N}d(N) - \log\frac{N}{g(N)} \to \infty, \quad \text{as}
\quad N \to \infty.
\end{equation}
\end{proposition}

Finally, Proposition~\ref{prop: power of d} in conjunction
with Proposition~\ref{prop: modified JSL} yields Theorem~\ref{th:pwr of d}.
The overall proof strategy as described above, is schematically
represented in Figure~\ref{fig:sfigrelation-inf}.

\begin{proof}[Proof of Theorem~\ref{th:pwr of d}]
Fix any $d(N)$ satisfying either~\eqref{eq:fNalike cond1}
or~\eqref{eq:fNalike cond2}.
From Propositions~\ref{prop: modified JSL} and~\ref{prop: power of d}
observe that it is enough to show that there exists an $n(N)$
with $n(N) \to \infty$ and $n(N) / g(N) \to 0$, as $N \to \infty$,
such that $$\frac{n(N)}{N}d(N)-\log\left(\frac{N}{g(N)}\right)\to\infty.$$

\noindent
(i) If $g(N) = O(N)$, then observe that $\log(N / g(N))$ is $O(1)$.
Since $d(N) \to \infty$, choosing $n(N) = N / \log (d(N))$ satisfies
the above criteria, and hence part (i) of the theorem is proved.

\noindent
(ii) Next we obtain a choice of $n(N)$ if $g(N)=o(N)$.
Note that, if
$$h(N) := \frac{d(N) \frac{g(N)}{N}}{\log\left(\frac{N}{g(N)}\right)} \to
\infty, \quad \text{as} \quad N \to \infty,$$
then choosing $n(N) = g(N) / \log(h(N))$, it can be seen that
as $N \to \infty$, $n(N) / g(N)\to 0$, and
\begin{equation}
\label{eq:equiv-nN-dN}
\begin{split}
&\frac{d(N)\frac{n(N)}{N}}{\log\left(\frac{N}{g(N)}\right)} =
\frac{h(N)}{\log(h(N))} \to \infty \\
\implies& \frac{n(N)}{N}d(N) - \log\left(\frac{N}{g(N)}\right) \to \infty,
\end{split}
\end{equation}
which completes the proof.
\end{proof}

\section{Load balancing in networks}
\label{networks}

In this section we return to the single-server queueing dynamics,
and extend the universality properties to network scenarios,
where the $N$~servers are assumed to be inter-connected by some
underlying graph topology~$G_N$.
Tasks arrive at the various servers as independent Poisson processes
of rate~$\lambda$, and each incoming task is assigned to whichever
server has the smallest number of tasks among the one where it arrives
and its neighbors in~$G_N$.
Thus, in case $G_N$ is a clique, each incoming task is assigned to the
server with the shortest queue across the entire system,
and the behavior is equivalent to that under the JSQ policy.
The stochastic optimality properties of the JSQ policy thus imply that
the queue length process in a clique will be better balanced
and smaller (in a majorization sense) than in an arbitrary graph~$G_N$.\\

\noindent
{\bf Sparse graph topologies.}
Besides the prohibitive communication overhead discussed earlier,
a further scalability issue of the JSQ policy arises when executing
a task involves the use of some data.
Storing such data for all possible tasks on all servers will typically
require an excessive amount of storage capacity.
These two burdens can be effectively mitigated in sparser graph
topologies where tasks that arrive at a specific server~$i$ are only
allowed to be forwarded to a subset of the servers ${\mathcal N}_i$.
For the tasks that arrive at server~$i$, queue length information
then only needs to be obtained from servers in ${\mathcal N}_i$,
and it suffices to store replicas of the required data on the
servers in ${\mathcal N}_i$.
The subset ${\mathcal N}_i$ containing the peers of server~$i$ can
be naturally viewed as its neighbors in some graph topology~$G_N$.
Here we consider the case of undirected graphs,
but most of the analysis can be extended to directed graphs.

While sparser graph topologies relieve the scalability issues
associated with a clique, the queue length process will be worse
(in the majorization sense) because of the limited connectivity.
Surprisingly, however, even quite sparse graphs can asymptotically
match the optimal performance of a clique, provided they are
suitably random, as we will further describe below.

The above model has been studied in~\cite{G15,T98}, focusing
on certain fixed-degree graphs and in particular ring topologies.
The results demonstrate that the flexibility to forward tasks
to a few neighbors, or even just one, with possibly shorter queues
significantly improves the performance in terms of the waiting time
and tail distribution of the queue length.
This resembles the \emph{power-of-choice} gains observed for JSQ($d$)
policies in complete graphs.

However, the results in \cite{G15,T98} also establish that
the performance sensitively depends on the underlying graph topology,
and that selecting from a fixed set of $d - 1$ neighbors typically
does not match the performance of re-sampling $d - 1$ alternate
servers for each incoming task from the entire population,
as in the power-of-$d$ scheme in a complete graph.
Further related problems have been investigated
in \cite{M96focs,ACM01,KM00,MPS02}.\\

\noindent
{\bf Balls and bins on graphs.}
If tasks do not get served and never depart but simply accumulate,
then the scenario described above amounts to a so-called balls-and-bins
problem on a graph.
Viewed from that angle, a close counterpart of our setup is studied
in Kenthapadi \& Panigrahy~\cite{KP06}, where in our terminology each
arriving task is routed to the shortest of $d \geq 2$ randomly
selected neighboring queues.
In this setup \cite{KP06} show that if the underlying graph is almost
regular with degree $N^{\varepsilon}$, where $\varepsilon$ is not too
small, the maximum number of balls in a bin scales as
$\log(\log(N)) / \log(d) + O(1)$, just like when the underlying graph
is a clique.
There are fundamental differences between the ball-and-bins and the
queueing scenarios, however, and an inherently different approach is
required in the current setup than what was developed in~\cite{KP06}.
Moreover, \cite{KP06} considers only the scaling of the maximum
queue length, whereas we analyze a more detailed time-varying
evolution of the entire system along with its stationary behavior.

We will further elaborate on the connections and differences
with balls-and-bins problems in Subsection~\ref{ballsbins}.\\

\noindent
{\bf Supermarket model on graphs.}
When each arriving task is routed to the shortest of $d \geq 2$ randomly
selected neighboring queues, the process-level convergence in the transient regime has been established recently in~\cite{BMW17}.
In this work, the authors analyze the evolution of the queue length process at an arbitrary tagged server as the system size becomes large.
The main ingredient is a careful analysis of local occupancy measures associated with the neighborhood of each server and to argue that under suitable conditions their asymptotic behavior is the same for all servers.
Under mild conditions on the graph topology $G_N$ (diverging minimum degree and the ratio between minimum degree and maximum degree in each connected component converges to 1), for a suitable initial occupancy measure, \cite[Theorem~2.1]{BMW17} establishes that for any fixed $d\geq 2$, the global occupancy state process for the JSQ($d$) scheme on $G_N$ has the same weak limit in~\eqref{fluid:standard} as that on a clique, as the number of vertices $N$ becomes large.
Also, the propagation of chaos property was shown to hold for this system, in the sense that the  queue lengths at any finite collection of tagged servers are asymptotically independent, and the
 queue length process for each server converges in distribution (in the path space) to the corresponding McKean-Vlasov process~\cite[Theorem~2.2]{BMW17}.
 Furthermore, when the graph sequence is random, with the $N$-th graph given as an Erd\H{o}s-R\'enyi random graph (ERRG) on $N$ vertices with average degree $c(N)$, annealed convergence of the occupancy process to the same deterministic limit as above, is established under the condition $c(N)\to\infty$, and under a stronger condition $c(N)/ \ln N\to\infty$, convergence (in probability) is shown for almost every realization of the random graph.\\

\noindent
{\bf Asymptotic optimality on graphs.}
We return to the case when each incoming task is assigned to whichever
server has the smallest number of tasks among the one where it arrives and its neighbors in~$G_N$.
The results presented in the remainder of the section are extracted
from~\cite{MBL17}, unless indicated otherwise.
As mentioned earlier, the queue length process in a clique will be
better balanced and smaller (in a majorization sense) than
in an arbitrary graph~$G_N$.
Accordingly, a graph $G_N$ is said to be $N$-optimal or
$\sqrt{N}$-optimal when the queue length process on $G_N$ is equivalent
to that on a clique on an $N$-scale or $\sqrt{N}$-scale, respectively.
Roughly speaking, a graph is $N$-optimal if the \emph{fraction}
of nodes with $i$~tasks, for $i = 0, 1, \ldots$, behaves as in a clique
as $N \to \infty$.
The fluid-limit results for the JSQ policy discussed
in Subsection~\ref{ssec:jsqfluid} imply that the latter fraction is
zero in the limit for all $i \geq 2$ in a clique in stationarity,
i.e., the fraction of servers with two or more tasks vanishes in any
graph that is $N$-optimal, and consequently the mean waiting time
vanishes as well as $N \to \infty$.
Furthermore, the diffusion-limit results of~\cite{EG15} for the JSQ
policy discussed in Subsection~\ref{ssec:diffjsq} imply that the
number of nodes with zero tasks and that with two tasks both scale as
$\sqrt{N}$ as $N \to \infty$.
Again loosely speaking, a graph is $\sqrt{N}$-optimal if in the
heavy-traffic regime the number of nodes with zero tasks and that
with two tasks when scaled by $\sqrt{N}$ both evolve as in a clique
as $N \to \infty$.
Formal definitions of asymptotic optimality on an $N$-scale or
$\sqrt{N}$-scale will be introduced in Definition~\ref{def:opt} below.

As one of the main results, we will demonstrate that, remarkably,
asymptotic optimality can be achieved in quite sparse ERRGs.

We prove that a sequence of ERRGs indexed by the number of vertices~$N$
with $d(N) \to \infty$ as $N \to \infty$, is $N$-optimal.
We further establish that the latter growth condition for the average
degree is in fact necessary in the sense that any graph sequence that
contains $\Theta(N)$ bounded-degree vertices cannot be $N$-optimal.
This implies that a sequence of ERRGs with finite average degree
cannot be $N$-optimal.
The growth rate condition is more stringent for optimality
on $\sqrt{N}$-scale in the heavy-traffic regime.
Specifically, we prove that a sequence of ERRGs indexed by the number
of vertices $N$ with $d(N) / (\sqrt{N} \log(N)) \to \infty$
as $N \to \infty$, is $\sqrt{N}$-optimal.

The above results demonstrate that the asymptotic optimality of cliques
on an $N$-scale and $\sqrt{N}$-scale can be achieved in far sparser graphs,
where the number of connections is reduced by nearly a factor $N$
and $\sqrt{N}/\log(N)$, respectively, provided the topologies are
suitably random in the ERRG sense.
This translates into equally significant reductions in communication
overhead and storage capacity, since both are roughly proportional
to the number of connections.\\

\noindent
{\bf Arbitrary graph topologies.}
While quite sparse graphs can achieve asymptotic optimality in the
presence of randomness, the worst-case graph instance may even
in very dense regimes (high average degree) not be optimal.
In particular, we prove that any graph sequence with minimum degree
$N-o(N)$ is $N$-optimal, but that for any $0 < c < 1/2$ one can
construct graphs with minimum degree $c N + o(N)$ which are not
$N$-optimal for some $\lambda < 1$.

The key challenge in the analysis of load balancing on arbitrary graph
topologies is that one needs to keep track
of the evolution of number of tasks at each vertex along with their
corresponding neighborhood relationship.
This creates a major problem in constructing a tractable Markovian
state descriptor, and renders a direct analysis of such processes
highly intractable.
Consequently, even asymptotic results for load balancing processes
on an arbitrary graph have remained scarce so far.
We take a radically different approach and aim to compare the load
balancing process on an arbitrary graph with that on a clique.
Specifically, rather than analyze the behavior for a given class
of graphs or degree value, we explore for what types of topologies
and degree properties the performance is asymptotically similar
to that in a clique.\\

\noindent
{\bf Stochastic coupling for graphs.}
Our proof arguments build on the stochastic coupling constructions
developed in Section~\ref{univ} for JSQ($d$) policies.
Specifically, we view the load balancing process on an arbitrary
graph as a `sloppy' version of that on a clique, and thus construct
several other intermediate sloppy versions.
By constructing novel couplings, we develop a method of comparing
the load balancing process on an arbitrary graph and that on a clique.
In particular, we bound the difference between the fraction of vertices
with $i$ or more tasks in the two systems for $i = 1, 2, \dots$,
to obtain asymptotic optimality results.
From a high level, conceptually related graph conditions for
asymptotic optimality were examined using quite different
techniques by Tsitsiklis and Xu \cite{TX11,TX17} in a dynamic
scheduling framework (as opposed to the load balancing context).\\

\noindent
{\bf Notation.}
For $k = 1, \ldots, N$, denote by $X_k(G_N, t)$ the queue length
at the $k$-th server at time~$t$ (including the task possibly in service),
and by $X_{(k)}(G_N, t)$ the queue length at the $k$-th ordered server
at time~$t$ when the servers are arranged in  non-decreasing order
of their queue lengths (ties can be broken in some way that will be
evident from the context).
Let $Q_i(G_N, t)$ denote the number of servers with queue length
at least~$i$ at time~$t$, $i = 1, 2, \ldots, B$.

It is important to note that $\{(q_i(G_N, t))_{i \geq 1}\}_{t \geq 0}$
is itself \emph{not} a Markov process, but the joint process
$\{(q_i(G_N, t))_{i \geq 1}, (X_{k}(G_N, t))_{k=1}^{N}\}_{t \geq 0}$
is Markov.
Also, in the Halfin-Whitt heavy-traffic regime~\eqref{eq:HW},
define the centered and scaled processes
\begin{equation}
\label{eq:HWOcc}
\bar{Q}_1(G_N,t) = - \frac{N - Q_1(G_N,t)}{\sqrt{N}}, \qquad
\bar{Q}_i(G_N,t) = \frac{ Q_i(G_N,t)}{\sqrt{N}},
\end{equation}
analogous to~\eqref{eq:diffscale}. \\

\noindent
{\bf Asymptotic optimality.}
As stated before, a clique is an optimal load balancing topology,
as the occupancy process is better balanced and smaller
(in a majorization sense) than in any other graph topology.
In general the optimality is strict, but it turns out that
near-optimality can be achieved asymptotically in a broad class
of other graph topologies.
Therefore, we now introduce two notions of \emph{asymptotic optimality},
which will be useful to characterize the performance in large-scale
systems.

\begin{definition}[{Asymptotic optimality}]
\label{def:opt}
A graph sequence $\GG = \{G_N\}_{N \geq 1}$ is called `asymptotically optimal on $N$-scale' or `$N$-optimal', if for any $\lambda < 1$,
the process $(q_1(G_N, \cdot), q_2(G_N, \cdot), \ldots)$
converges weakly, on any finite time interval, to the process
$(q_1(\cdot), q_2(\cdot),\ldots)$ given by~\eqref{eq:fluid}.

Moreover, a graph sequence $\GG = \{G_N\}_{N \geq 1}$ is called
`asymptotically optimal on $\sqrt{N}$-scale' or `$\sqrt{N}$-optimal',
if in the Halfin-Whitt heavy-traffic regime~\eqref{eq:HW},
on any finite time interval, the process
$(\bQ_1(G_N, \cdot), \bQ_2(G_N, \cdot), \ldots)$ as in~\eqref{eq:HWOcc}
converges weakly to the process $(\bQ_1(\cdot), \bQ_2(\cdot), \ldots)$
given by~\eqref{eq:diffusionjsq}.
\end{definition}

Intuitively speaking, if a graph sequence is $N$-optimal
or $\sqrt{N}$-optimal, then in some sense, the associated occupancy
processes are indistinguishable from those of the sequence of cliques
on $N$-scale or $\sqrt{N}$-scale.
In other words, on any finite time interval their occupancy processes
can differ from those in cliques by at most $o(N)$ or $o(\sqrt{N})$,
respectively.
For brevity, $N$-scale and $\sqrt{N}$-scale will henceforth be referred
to as \emph{fluid scale} and \emph{diffusion scale}, respectively.

In particular, exploiting interchange of the stationary ($t \to \infty$)
and many-server ($N \to \infty$) limits, we obtain that for any
$N$-optimal graph sequence $\{G_N\}_{N \geq 1}$, as $N \to \infty$
\begin{equation}
q_1(G_N, \infty) \to \lambda \quad \mbox{ and} \quad
q_i(G_N, \infty) \to 0 \quad \mbox{ for all } i = 2, \dots, B,
\end{equation}
implying that the stationary fraction of servers with queue length two
or larger and the mean waiting time vanish.

\subsection{Asymptotic optimality criteria for deterministic graph sequences}

We now proceed to develop a criterion for asymptotic optimality
of an arbitrary deterministic graph sequence on different scales.
Next this criterion will be leveraged to establish optimality
of a sequence of random graphs.
We start by introducing some useful notation, and two measures
of \emph{well-connectedness}.
Let $G = (V, E)$ be any graph.
For a subset $U \subseteq V$, define $\com(U) := |V\setminus N[U]|$
to be the cardinality of the set of all vertices that are disjoint from $U$, where
$N[U] := U\cup \{v \in V:\ \exists\ u \in U \mbox{ with } (u, v) \in E\}$.
For any fixed $\varepsilon > 0$ define
\begin{equation}
\label{def:dis}
\dis_1(G,\varepsilon) := \sup_{U\subseteq V, |U|\geq \varepsilon |V|}\com(U),
\qquad
\dis_2(G,\varepsilon) := \sup_{U\subseteq V, |U|\geq \varepsilon \sqrt{|V|}}\com(U).
\end{equation}

The next theorem provides sufficient conditions for asymptotic
optimality on $N$-scale and $\sqrt{N}$-scale in terms of the above
two well-connectedness measures.

\begin{theorem}
\label{th:det-seq}
For any graph sequence $\GG = \{G_N\}_{N \geq 1}$,
\begin{enumerate}[{\normalfont (i)}]
\item $\GG$ is $N$-optimal if for any $\varepsilon > 0$,
$\dis_1(G_N, \varepsilon) / N \to 0$ as $N \to \infty$.
\item $\GG$ is $\sqrt{N}$-optimal if for any $\varepsilon > 0$,
$\dis_2(G_N, \varepsilon) / \sqrt{N} \to 0$ as $N \to \infty$.
\end{enumerate}
\end{theorem}

From a high level, the conditions in Theorem~\ref{th:det-seq} (i) and (ii) require that neighborhoods of any $\Theta(N)$ and $\Theta(\sqrt{N})$ vertices contain at least $N-o(N)$ and $N-o(\sqrt{N})$ vertices, respectively.
As we will see below in Theorems~\ref{th:inhom} and~\ref{th:reg}, the conditions in Theorem~\ref{th:det-seq} impose right levels of connectivity in the graph topology in order for it to be asymptotically optimal on fluid and diffusion scales, while significantly reducing the total number of connections.
The next corollary is an immediate consequence
of Theorem~\ref{th:det-seq}.

\begin{corollary}
Let $\GG = \{G_N\}_{N \geq 1}$ be any graph sequence.
Then {\rm(i)} If $d_{\min}(G_N) = N - o(N)$, then $\GG$ is $N$-optimal,
and
{\rm(ii)} If $d_{\min}(G_N) = N - o(\sqrt{N})$, then $\GG$ is
$\sqrt{N}$-optimal.
\end{corollary}
The rest of the subsection is devoted to a discussion of the main
proof arguments for Theorem~\ref{th:det-seq}, focusing on the proof
of $N$-optimality.
The proof of $\sqrt{N}$-optimality follows along similar lines.
We establish in Proposition~\ref{prop:cjsq} that if a system is able
to assign each task to a server in the set $\cS^N(n(N))$ of the $n(N)$
nodes with shortest queues, where $n(N)$ is $o(N)$,
then it is $N$-optimal.
Since the underlying graph is not a clique however (otherwise there
is nothing to prove), for any $n(N)$ not every arriving task can be
assigned to a server in $\cS^N(n(N))$.
Hence we further prove in Proposition~\ref{prop:stoch-ord-new}
a stochastic comparison property implying that if on any finite time
interval of length~$t$, the number of tasks $\Delta^N(t)$ that are not
assigned to a server in $\cS^N(n(N))$ is $o_P(N)$, then the system is
$N$-optimal as well.
The $N$-optimality can then be concluded when $\Delta^N(t)$ is
$o_P(N)$, which we establish in Proposition~\ref{prop:dis-new}
under the condition that $\dis_1(G_N, \varepsilon) / N \to 0$ as
$N \to \infty$ as stated in Theorem~\ref{th:det-seq}.

To further explain the idea described in the above proof outline,
it is useful to adopt a slightly different point of view towards
load balancing processes on graphs.
From a high level, a load balancing process can be thought of as follows:
there are $N$~servers, which are assigned incoming tasks by some scheme.

The assignment scheme can arise from some topological structure,
in which case we will call it \emph{topological load balancing},
or it can arise from some other property of the occupancy process,
in which case we will call it \emph{non-topological load balancing}.
As mentioned earlier, the JSQ policy or the clique is optimal
among the set of all non-anticipating schemes, irrespective of being
topological or non-topological.
Also, load balancing on graph topologies other than a clique can be
thought of as a `sloppy' version of that on a clique, when each server
only has access to partial information on the occupancy state.
Below we first introduce a different type of sloppiness in the task
assignment scheme, and show that under a limited amount of sloppiness
optimality is retained on a suitable scale.
Next we will construct a scheme which is a hybrid of topological
and non-topological schemes, whose behavior is simultaneously close
to both the load balancing process on a suitable graph and that
on a clique.

\paragraph{A class of sloppy load balancing schemes.}

Fix some function $n: \N \to \N$, and recall the set $\cS^N(n(N))$
as before as well as the class $\CJSQ(n(N))$ where each arriving task
is assigned to one of the servers in $\cS^N(n(N))$.
It should be emphasized that for any scheme in $\CJSQ(n(N))$,
we are not imposing any restrictions on how the incoming task should
be assigned to a server in $\cS^N(n(N))$.
The scheme only needs to ensure that the arriving task is assigned
to some server in $\cS^N(n(N))$ with respect to \emph{some} tie
breaking mechanism.
Observe that using Corollary~\ref{cor:bound} and following the
arguments as in the proof of Theorems~\ref{fluidjsqd}
and~\ref{diffusionjsqd}, we obtain the next proposition,
which provides a sufficient criterion for asymptotic optimality
of any scheme in $\CJSQ(n(N))$.

\begin{proposition}
\label{prop:cjsq}
For $0 \leq n(N) < N$, let $\Pi \in \CJSQ(n(N))$ be any scheme.
{\rm(i)} If $n(N) / N \to 0$ as $N \to \infty$, then $\Pi$ is $N$-optimal,
and
{\rm(ii)} If $n(N) / \sqrt{N} \to 0$ as $N \to \infty$, then $\Pi$ is
$\sqrt{N}$-optimal.
\end{proposition}

\paragraph{A bridge between topological and non-topological load balancing.}
For any graph $G_N$ and $n \leq N$, we first construct a scheme
called $I(G_N, n)$, which is an intermediate blend between the
topological load balancing process on $G_N$ and some kind
of non-topological load balancing on $N$~servers.
The choice of $n = n(N)$ will be clear from the context.

To describe the scheme $I(G_N, n)$, first synchronize the arrival
epochs at server~$v$ in both systems, $v = 1, 2, \ldots, N$.
Further, synchronize the departure epochs at the $k$-th ordered
server with the $k$-th smallest number of tasks in the two systems,
$k = 1, 2, \ldots, N$.
When a task arrives at server~$v$ at time~$t$ say, it is assigned
in the graph $G_N$ to a server $v' \in N[v]$ according to its own
statistical law.
For the assignment under the scheme $I(G_N,n)$, first observe that if
\begin{equation}
\label{eq:criteria}
\min_{u \in N[v]} X_u(G_N,t) \leq \max_{u \in \cS(n)} X_u(G_N,t),
\end{equation}
then there exists \emph{some} tie-breaking mechanism for which
$v' \in N[v]$ belongs to $\cS(n)$ under $G_N$.
Pick such an ordering of the servers, and assume that $v'$ is the
$k$-th ordered server in that ordering, for some $k \leq n+1$.
Under $I(G_N,n)$ assign the arriving task to the $k$-th ordered
server (breaking ties arbitrarily in this case).
Otherwise, if \eqref{eq:criteria} does not hold, then the task is
assigned to one of the $n+1$ servers with minimum queue lengths
under $G_N$ uniformly at random.

Denote by $\Delta^N(I(G_N, n), T)$ the cumulative number of arriving
tasks up to time $T \geq 0$ for which Equation~\eqref{eq:criteria} is
violated under the above coupling.
The next proposition shows that the load balancing process under the
scheme $I(G_N,n)$ is close to that on the graph $G_N$ in terms of the
random variable $\Delta^N(I(G_N, n), T)$.

\begin{proposition}
\label{prop:stoch-ord-new}
The following inequality is preserved almost surely
\begin{equation}
\label{eq:stoch-ord-new}
\sum_{i=1}^{B} |Q_i(G_N, t) - Q_i(I(G_N, n), t)| \leq
2 \Delta^N(I(G_N, n), t) \qquad \forall\ t \geq 0,
\end{equation}
provided the two systems start from the same occupancy state at $t=0$.
\end{proposition}

In order to conclude optimality on $N$-scale or $\sqrt{N}$-scale,
it remains to be shown that the term $\Delta^N(I(G_N, n), T)$ is sufficiently small.
The next proposition provides suitable asymptotic bounds for
$\Delta^N(I(G_N, n), T)$ under the conditions on $\dis_1(G_N, \varepsilon)$
and $\dis_2(G_N, \varepsilon)$ stated in Theorem~\ref{th:det-seq}.

\begin{proposition}
\label{prop:dis-new}
\begin{enumerate}[{\normalfont (i)}]
\item For any $\varepsilon>0$, there exists $\varepsilon'>0$ and $n_{\varepsilon'}(N)$ with $n_{\varepsilon'}(N)/N\to 0$ as $N\to\infty$, such that if $\dis_1(G_N,\varepsilon')/N\to 0$ as $N\to\infty$, then for all $T>0$, $$\Pro{\Delta^N(I(G_N,n_{\varepsilon'}),T)/N>\varepsilon}\to 0.$$
\item For any $\varepsilon>0$, there exists $\varepsilon'>0$ and $m_{\varepsilon'}(N)$ with $m_{\varepsilon'}(N)/\sqrt{N}\to 0$ as $N\to\infty$, such that if $\dis_2(G_N,\varepsilon')/\sqrt{N}\to 0$ as $N\to\infty$, then for all $T>0$, $$\Pro{\Delta^N(I(G_N,m_{\varepsilon'}),T)/\sqrt{N}>\varepsilon}\to 0.$$
\end{enumerate}
\end{proposition}

The proof of Theorem~\ref{th:det-seq} then readily follows
by combining Propositions~\ref{prop:cjsq}-\ref{prop:dis-new}
and observing that the scheme $I(G_N,n)$ belongs to the class
$\CJSQ(n)$ by construction.

\begin{proof}[Proof of Proposition~\ref{prop:dis-new}]
Fix any $\varepsilon, T>0$ and choose
$\varepsilon' = \varepsilon / (2 \lambda T)$.
With the coupling described above, when a task arrives at some
vertex~$v$ say, Equation~\eqref{eq:criteria} is violated only if none
of the vertices in $\cS(n_{\varepsilon'}(N))$ is a neighbor of~$v$.
Thus, the total instantaneous rate at which this happens is
\[
\lambda \com(\cS(n_{\varepsilon'}(N),t)) \leq
\lambda \sup_{U \subseteq V_N, |U| \geq n_{\varepsilon'}(N)} \com(U),
\]
irrespective of what this set $\cS^N(n(N))$ actually is.
Therefore, for any fixed $T \geq 0$,
\[
\Delta^N(I(G_N,n_{\varepsilon'}),T) \leq
A\Big(\lambda \sup_{U \subseteq V_N, |U| \geq n_{\varepsilon'}(N)} \com(U)\Big),
\]
where $A(\cdot)$ represents a unit-rate Poisson process.
This can then be leveraged to show that
$\Delta^N(I(G_N, n_{\varepsilon'}), T)$ is small on an $N$-scale
and $\sqrt{N}$-scale, respectively, under the conditions stated in the
proposition, by choosing a suitable $n_{\varepsilon'}$.

Specifically, if $\dis_1(G_N,\varepsilon')/N \to 0$, then there exists
$n_{\varepsilon'}(N)$ with $n_{\varepsilon'}(N) / N \to 0$ such that
$\dis_1(G_N,\varepsilon') \leq n_{\varepsilon'}(N)$ for all $N \geq 1$,
and hence
$\sup_{U \subseteq V_N, |U| \geq n_{\varepsilon'}(N)} \com(U) \leq
\varepsilon' N$.
It then follows that with high probability,
\[
\limsup_{N \to \infty} \frac{1}{N} \Delta^N(I(G_N, n_{\varepsilon'}), T) \leq
\limsup_{N \to \infty} \frac{1}{N} A \Big(\lambda T \varepsilon' N\Big) \leq
2 \lambda T \varepsilon' = \varepsilon.
\]

Likewise, if $\dis_2(G_N,\varepsilon')/\sqrt{N} \to 0$, then there exists
$m_{\varepsilon'}(N)$ with $m_{\varepsilon'}(N)/\sqrt{N}\to 0$ such that
$\dis_2(G_N,\varepsilon') \leq m_{\varepsilon'}(N)$ for all $N \geq 1$,
and hence
$\sup_{U \subseteq V_N, |U| \geq m_{\varepsilon'}(N)} \com(U) \leq
\varepsilon' \sqrt{N}$.
It then follows that with high probability,
\[
\limsup_{N \to \infty} \frac{1}{\sqrt{N}} \Delta^N(I(G_N, m_{\varepsilon'}),T) \leq
\limsup_{N \to \infty} \frac{1}{\sqrt{N}} A\Big(\lambda T \varepsilon' \sqrt{N}\Big) \leq
2 \lambda T \varepsilon' = \varepsilon
\]
and the proof is completed.
\end{proof}

From the conditions of Theorem~\ref{th:det-seq} it follows that
if for all $\varepsilon > 0$, $\dis_1(G_N, \varepsilon)$
and $\dis_2(G_N, \varepsilon)$ are $o(N)$ and $o(\sqrt{N})$,
respectively, then the total number of edges in $G_N$ must be
$\omega(N)$ and $\omega(N\sqrt{N})$, respectively.
Theorem~\ref{th:bdd-deg} below states that the \emph{super-linear}
growth rate of the total number of edges is not only sufficient,
but also necessary in the sense that any graph with $O(N)$ edges is
asymptotically sub-optimal on $N$-scale.

\begin{theorem}
\label{th:bdd-deg}
Let $\GG = \{G_N\}_{N \geq 1}$ be any graph sequence, such that there
exists a fixed integer $M < \infty$ with
\begin{equation}
\label{eq:bdd-deg}
\limsup_{N \to \infty} \dfrac{\#\big\{v \in V_N: d_v \leq M\big\}}{N} > 0,
\end{equation}
where $d_v$ is the degree of the vertex~$v$.
Then $\GG$ is sub-optimal on $N$-scale.
\end{theorem}

To prove Theorem~\ref{th:bdd-deg}, we show that starting
from an all-empty state, in finite time, a positive fraction
of servers in $G_N$ will have at least two tasks.
This will prove that the occupancy processes when scaled by~$N$ cannot
agree with those in the sequence of cliques, and hence $\{G_N\}_{N \geq 1}$
cannot be $N$-optimal.
The idea of the proof can be explained as follows: If a system contains
$\Theta(N)$ bounded-degree vertices, then starting from an all-empty state,
in any finite time interval there will be $\Theta(N)$ servers~$u$ say,
for which all the servers in $N[u]$ have at least one task.
For all such servers an arrival at~$u$ must produce a server
with queue length two.
Thus, it shows that the instantaneous rate at which servers of queue
length two are formed is bounded away from zero, and hence $\Theta(N)$
servers of queue length two are produced in finite time.

\paragraph{Worst-case scenario.}
Next we consider the worst-case scenario.
Theorem~\ref{th:min-deg-negative} below asserts that a graph sequence
can be sub-optimal for some $\lambda < 1$ even when the minimum degree
$d_{\min}(G_N)$ is~$\Theta(N)$.

\begin{theorem}
\label{th:min-deg-negative}
For any $\big\{d(N)\big\}_{N \geq 1}$, such that $d(N) / N\to c$
with $0 < c < 1/6$, there exists $\lambda < 1$, and a graph sequence
$\GG = \{G_N\}_{N \geq 1}$ with $d_{\min}(G_N) = d(N)$, such that
$\GG$ is sub-optimal on $N$-scale.
\end{theorem}

To construct such a sub-optimal graph sequence, consider a sequence
of complete bipartite graphs $G_N = (V_N, E_N)$,
with $V_N = A_N \sqcup B_N$ and $|A_N| / N \to c \in (0, 1/2)$
as $N \to \infty$.
If this sequence were $N$-optimal, then starting from an all-empty state,
asymptotically the fraction of servers with queue length one would
converge to~$\lambda$, and the fraction of servers with queue length
two or larger should remain zero throughout.
Now note that for large $N$ the rate at which tasks join the empty
servers in $A_N$ is given by $(1-c) \lambda$, whereas the rate
of empty server generation in $A_N$ is at most~$c$.
Choosing $\lambda > c / (1-c)$, one can see that in finite time each
server in $A_N$ will have at least one task.
From that time onward with at least instantaneous rate
$\lambda (\lambda - c) - c$, servers with queue length two start forming.
The range for $c$ stated in Theorem~\ref{th:min-deg-negative} is only to
ensure that there exists $\lambda<1$ with $\lambda (\lambda - c) - c > 0$.

\subsection{Asymptotic optimality of random graph sequences}

Next we investigate how the load balancing process behaves on random
graph topologies.
Specifically, we aim to understand what types of graphs are
asymptotically optimal in the presence of randomness
(i.e., in an average-case sense).
Theorem~\ref{th:inhom} below establishes sufficient conditions
for asymptotic optimality of a sequence of inhomogeneous random graphs.
Recall that a graph $G' = (V', E')$ is called a supergraph
of $G = (V, E)$ if $V = V'$ and $E \subseteq E'$.

\begin{theorem}
\label{th:inhom}
Let $\GG= \{G_N\}_{N \geq 1}$ be a graph sequence such that for each~$N$,
$G_N = (V_N, E_N)$ is a supergraph of the inhomogeneous random graph
$G_N'$ where any two vertices $u, v \in V_N$ share an edge
with probability $p_{uv}^N$.
\begin{enumerate}[{\normalfont (i)}]
\item If $\inf\ \{p^N_{uv}: u, v\in V_N\}$ is $\omega(1/N)$, then $\GG$ is $N$-optimal.
\item If $\inf\ \{p^N_{uv}: u, v\in V_N\}$ is $\omega(\log(N)/\sqrt{N})$, then $\GG$ is $\sqrt{N}$-optimal.
\end{enumerate}
\end{theorem}

The proof of Theorem~\ref{th:inhom} relies on Theorem~\ref{th:det-seq}.
Specifically, if $G_N$ satisfies conditions~(i) and~(ii) in
Theorem~\ref{th:inhom}, then the corresponding conditions~(i) and~(ii)
in Theorem~\ref{th:det-seq} hold.

As an immediate corollary of Theorem~\ref{th:inhom} we obtain
an optimality result for the sequence of ERRGs.
Let $\ER_N(p(N))$ denote a graph on $N$ vertices, such that any pair of vertices share an edge with probability $p(N)$.

\begin{corollary}
\label{cor:errg}
Let $\GG = \{G_N\}_{N \geq 1}$ be a graph sequence such that for each~$N$,
$G_N$ is a super-graph of $\ER_N(p(N))$, and $d(N) = (N-1) p(N)$.
Then
{\normalfont (i)}
If $d(N)\to\infty$ as $N\to\infty$, then $\GG$ is $N$-optimal.
{\normalfont (ii)}
If $d(N)/(\sqrt{N}\log(N))\to\infty$ as $N\to\infty$, then $\GG$ is
$\sqrt{N}$-optimal.
\end{corollary}

Theorem~\ref{th:det-seq} can be further leveraged to establish the
optimality of the following sequence of random graphs.
For any $N \geq 1$ and $d(N) \leq N-1$ such that $N d(N)$ is even,
construct the \emph{erased random regular} graph on $N$~vertices
as follows:
Initially, attach $d(N)$ \emph{half-edges} to each vertex.
Call all such half-edges \emph{unpaired}.
At each step, pick one half-edge arbitrarily, and pair it to another
half-edge uniformly at random among all unpaired half-edges
to form an edge, until all the half-edges have been paired.
This results in a uniform random regular multi-graph with degree
$d(N)$~\cite[Proposition 7.7]{remco-book-1}.

Now the erased random regular graph is formed by erasing all the
self-loops and multiple edges, which then produces a simple graph.

\begin{theorem}
\label{th:reg}
Let $\GG = \{G_N\}_{N \geq 1}$ be a sequence of erased random regular
graphs with degree $d(N)$.
Then
{\normalfont (i)}
If $d(N)\to\infty$ as $N\to\infty$, then $\GG$ is $N$-optimal.
{\normalfont (ii)}
If $d(N)/(\sqrt{N}\log(N))\to\infty$ as $N\to\infty$, then $\GG$ is $\sqrt{N}$-optimal.
\end{theorem}

Note that due to Theorem~\ref{th:bdd-deg}, we can conclude that the
growth rate condition for $N$-optimality in Corollary~\ref{cor:errg}~(i)
and Theorem~\ref{th:reg}~(i) is not only sufficient, but necessary as well.
Thus informally speaking, $N$-optimality is achieved under the minimum
condition required as long as the underlying topology is suitably random.

\section{Token-based load balancing}
\label{token}

While a zero waiting time can be achieved in the limit by sampling
only $d(N) =o(N)$ servers as Sections~\ref{univ} and~\ref{networks} showed,
even in network scenarios, the amount of communication overhead
in terms of $d(N)$ must still grow with~$N$.
As mentioned earlier, this can be avoided by introducing memory at the
dispatcher, in particular maintaining a record of only vacant servers,
and assigning tasks to idle servers, if there are any,
or to a uniformly at random selected server otherwise.
This so-called Join-the-Idle-Queue (JIQ) scheme \cite{BB08,LXKGLG11}
can be implemented through a simple token-based mechanism generating
at most one message per task.
Remarkably enough, even with such low communication overhead,
the mean waiting time and the probability of a non-zero waiting time
vanish under the JIQ scheme in both the fluid and diffusion regimes,
as we will discuss in the next two subsections.

\subsection{Fluid-level optimality of JIQ scheme}
\label{ssec:fluidjiq}

We first consider the fluid limit of the JIQ policy.
Recall that $q_i^N(\infty)$ denotes a random variable denoting the
process $q_i^N(\cdot)$ in steady state.
Under significantly more general conditions (in the presence
of finitely many heterogeneous server pools and for general service
time distributions with decreasing hazard rate) it was proved
in~\cite{Stolyar15} that under the JIQ scheme
\begin{equation}
\label{eq:fpjiq}
q_1^N(\infty) \to \lambda, \qquad q_i^N(\infty) \to 0 \quad
\mbox{ for all } i \geq 2, \qquad \mbox{ as } \quad N \to \infty.
\end{equation}
The above equation in conjunction with the PASTA property yields that
the steady-state probability of a non-zero wait vanishes as $N \to \infty$,
thus exhibiting asymptotic optimality of the JIQ policy on fluid scale.

\paragraph{High-level proof idea.}
Loosely speaking, the proof of~\eqref{eq:fpjiq} consists of three
principal components:
\begin{enumerate}[{\normalfont (i)}]
\item Starting from an all-empty state, observe that the asymptotic
rate of increase of~$q_1$ is given by the arrival rate~$\lambda$.
Also, the rate of decrease is~$q_1$.
Thus, on a small time interval $\dif t$, the rate of change of~$q_1$
is given by
\begin{equation}\label{eq:fluidjiq}
\frac{\dif q_1(t)}{\dif t} = \lambda - q_1(t).
\end{equation}
Under the above dynamics, the system occupancy states converge to the
unique fixed point of the above ODE, given by $(\lambda, 0, 0, \ldots)$.
\item The occupancy process is monotone, in the sense that
(a) Starting from an all-empty state, the occupancy process is
componentwise stochastically nondecreasing in time, and
(b) The occupancy process at any fixed time~$t$ starting
from an arbitrary state is componentwise stochastically dominated
by the occupancy process at time~$t$ starting from an all-empty state.
\item Under the JIQ scheme, the system is stable, and hence the
occupancy process is ergodic.
Since $q_1(t)$ is the instantaneous rate of departure from the system,
ergodicity implies that in steady state there can be {\em at most}
$\lambda$ fraction of busy servers (containing at least one task).
In fact, it further establishes that the steady-state fraction
of servers with more than one tasks vanishes asymptotically.
\end{enumerate}

Points~(i) and~(ii) above imply that starting from any state the
system must have at least $\lambda$ fraction of busy servers,
and finally this along with Point (iii) establishes that the
steady-state occupancy process must converge to $(\lambda, 0, 0, \ldots)$.

\subsection{Diffusion-level optimality of JIQ scheme}

We now turn to the diffusion limit of the JIQ scheme.
Recall the centered and scaled occupancy process
as in~\eqref{eq:diffscale}, and the Halfin-Whitt heavy-traffic regime
in~\eqref{eq:HW}.

\begin{theorem}[Diffusion limit for JIQ]
\label{diffusionjiq}
Assume that $\lambda(N)$ satisfies~\eqref{eq:HW}.
Under suitable initial conditions the weak limit of the sequence of centered and diffusion-scaled occupancy process in~\eqref{eq:diffscale} coincides with that of the ordinary JSQ policy, and in particular, is given by the system of SDEs in~\eqref{eq:diffusionjsq}.
\end{theorem}

The above theorem implies that for suitable states,
on any finite time interval, the occupancy process of a system
under the JIQ policy is indistinguishable from that under the JSQ policy.

\paragraph{High-level proof idea.}
The proof of Theorem~\ref{diffusionjiq} relies on a novel coupling
construction introduced in~\cite{MBLW16-1} as described below in detail.
The idea is to compare the occupancy processes of two systems
following JIQ and JSQ policies, respectively.
Comparing the JIQ and JSQ policies is facilitated when viewed as follows:
(i) If there is an idle server in the system, both JIQ and JSQ perform
 similarly,
(ii)~Also, when there is no idle server and only $O(\sqrt{N})$ servers
with queue length two, JSQ assigns the arriving task to a server
with queue length one.
In that case, since JIQ assigns at random, the probability that the
task will land on a server with queue length two and thus JIQ acts
differently than JSQ is $O(1/\sqrt{N})$.
Since on any finite time interval the number of times an arrival finds
all servers busy is at most $O(\sqrt{N})$, all the arrivals except
an $O(1)$ of them are assigned in exactly the same manner in both JIQ
and JSQ, which then leads to the same scaling limit for both policies. \\

The diffusion limit result in Theorem~\ref{diffusionjiq} is in fact
true for an even broader class of load balancing schemes.
Recall that $B$ denotes the buffer capacity (possibly infinite)
of each server, and in case $B < \infty$, if a task is assigned
to a server with $B$~outstanding tasks, it is instantly discarded.
Define the class of schemes
$$\Pi^{(N)} := \{\Pi(d_0, d_1, \ldots, d_{B-1}): d_0 = N, 1 \leq d_i \leq N,
1 \leq i \leq B-1, B \geq 2\},$$ where
in the scheme $\Pi(d_0, d_1, \ldots, d_{B-1}),$ the dispatcher assigns
an incoming task to the server with the minimum queue length among $d_k$
(possibly function of~$N$) servers selected uniformly at random when the
minimum queue length across the system is~$k$, $k = 0, 1, \ldots, B - 1$.
The system analyzed in~\cite{EG15} (JSQ with $B = 2$) can be written
as $\Pi(N,N)$, JIQ can be expressed as $\Pi(N, 1, 1,\ldots)$,
and JIQ with a buffer capacity $B = 2$ is $\Pi(N, 1)$.

The crux of the argument in proving diffusion-level optimality for any
scheme in $\Pi^{(N)}$ goes as follows:
First the scheme $\Pi(N, d_1, \ldots, d_{B-1})$ is sandwiched
between $\Pi(N, 1)$ and $\Pi(N, d_1)$.
More specifically, the gap between $\Pi(N, d_1, \ldots, d_{B-1})$
and $\Pi(N, 1)$ is bounded by the number of items lost due to full buffers.
Next, this loss is bounded using the number of servers with queue length
$2$ in $\Pi(N,N)$.
This allows the use of the results in~\cite{EG15}, and yields that
on any finite time interval with high probability an $O(1)$ number of items
are lost due to full buffers, which is negligible on $\sqrt{N}$ scale.
Specifically, this shows that for suitable initial states,
the schemes $\Pi(N, 1)$ and $\Pi(N, d_1)$, along with any scheme
in the class $\Pi^{(N)}$ has the same diffusion limit in Halfin-Whitt
heavy-traffic regime.
We conclude this subsection by describing the coupling construction
stating the stochastic inequalities, and a brief proof sketch
for Theorem~\ref{diffusionjiq}.

\paragraph{The coupling construction.}
We now construct a coupling between two systems following any two schemes,
say $\Pi_1 = \Pi(l_0, l_1, \ldots, l_{B-1})$
and $\Pi_2 = \Pi(d_0, d_1, \ldots, d_{B'-1})$ in $\Pi^{(N)}$, respectively,
to establish the desired stochastic ordering results.
With slight abuse of notation we will denote by $\Pi_i$ the system
following scheme $\Pi_i$, $i = 1, 2$.

For the arrival process we couple the two systems as follows. First we synchronize the  arrival epochs of the two systems. Now assume that in the systems $A$ and $B$, the minimum queue lengths are $k$ and $m$, respectively, $k\leq B-1$, $m\leq B'-1$. Therefore, when a task arrives, the dispatchers in $\Pi_1$ and $\Pi_2$ have to select $l_k$ and $d_m$ servers, respectively, and then have to send the task to the one having the minimum queue length among the respectively selected servers. Since the servers are being selected uniformly at random we can assume without loss of generality, as in the stack construction, that the servers are arranged in non-decreasing order of their queue lengths and are indexed in increasing order.
Hence, observe that when a few server indices are selected, the server having the minimum of those indices will be the server with the minimum queue length among these. Hence, in this case the dispatchers in $\Pi_1$ and $\Pi_2$ select $l_k$ and $d_m$ random numbers (without replacement) from $\{1,2,\ldots,N\}$ and then send the incoming task to the servers having indices to be the minimum of those selected numbers. To couple the decisions of the two systems, at each arrival epoch a single random permutation of $\{1,2,\ldots,N\}$ is drawn, denoted by $\boldsymbol{\Sigma}^{(N)}:=(\sigma_1, \sigma_2,\ldots,\sigma_N)$. Define $\sigma_{(i)}:= \min_{j\leq i}\sigma_j$. Then observe that system $\Pi_1$ sends the task to the server with the index $\sigma_{(l_k)}$ and system $\Pi_2$ sends the task to the server with the index $\sigma_{(d_m)}$. Since at each arrival epoch both systems use a common random permutation, they take decisions in a coupled manner.

For the potential departure process, couple the service completion times of the $k^{th}$ queue in both scenarios, $k= 1,2,\ldots,N$. More precisely, for the potential departure process assume that we have a single synchronized exp($N$) clock independent of arrival epochs for both systems. Now when this clock rings, a number $k$ is uniformly selected from $\{1,2,\ldots,N\}$ and a potential departure occurs from the $k^{th}$ queue in both systems. If at a potential departure epoch an empty queue is selected, then we do nothing. In this way the two schemes, considered independently, still evolve according to their appropriate statistical laws.

\begin{proposition}
\label{prop: stoch_ord}
For two schemes $\Pi_1=\Pi(l_0,l_1,\ldots,l_{B-1})$ and $\Pi_2=\Pi(d_0, d_1,\ldots, d_{B'-1})$ with $B\leq B'$ assume $l_0=\ldots=l_{B-2}=d_0=\ldots=d_{B-2}=d$, $l_{B-1}\leq d_{B-1}$ and either $d=N$ or $d\leq d_{B-1}$. Then the following holds:
\begin{enumerate}[{\normalfont (i)}]
\item\label{component_ordering} $\{Q^{ \Pi_1}_i(t)\}_{t\geq 0}\leq_{st}\{Q^{ \Pi_2}_i(t)\}_{t\geq 0}$ for $i=1,2,\ldots,B$,
\item\label{upper bound} $\{\sum_{i=1}^B Q^{ \Pi_1}_i(t)+L^{ \Pi_1}(t)\}_{t\geq 0}\geq_{st} \{\sum_{i=1}^{B'} Q^{ \Pi_2}_i(t)+L^{ \Pi_2}(t)\}_{t\geq 0}$,
\item\label{delta_ineq} $\{\Delta(t)\}_{t\geq 0}\geq \{\sum_{i=B+1}^{B'}Q_i^{ \Pi_2}(t)\}_{t\geq 0}$ almost surely under the coupling defined above,
\end{enumerate}
for any fixed $N\in\mathbbm{N}$ where $\Delta(t):=L^{ \Pi_1}(t)-L^{ \Pi_2}(t)$, provided that at time $t=0$ the above ordering holds.
\end{proposition}

\begin{proof}[Proof of Theorem~\ref{diffusionjiq}]
Let $\Pi=\Pi(N,d_1,\ldots,d_{B-1})$ be a load balancing scheme in the class $\Pi^{(N)}$. Denote by $\Pi_1$ the scheme $\Pi(N,d_1)$ with buffer size $B=2$ and let $\Pi_2$ denote the JIQ policy $\Pi(N,1)$ with buffer size $B=2$.

Observe that from Proposition~\ref{prop: stoch_ord} we have under the coupling defined above,
\begin{equation}
\label{eq: bound}
\begin{split}
|Q_i^\Pi(t) - Q_i^{\Pi_2}(t)| & \leq |Q_i^{\Pi}(t) - Q_i^{\Pi_1}(t)| +
|Q_i^{\Pi_1}(t) - Q_i^{\Pi_2}(t)| \\ & \leq
|L^{\Pi_1}(t) - L^{\Pi}(t)| + |L^{\Pi_2}(t) - L^{\Pi_1}(t)| \\
& \leq 2 L^{\Pi_2}(t),
\end{split}
\end{equation}
for all $i \geq 1$ and $t \geq 0$ with the understanding that $Q_j(t)=0$
for all $j > B$, for a scheme with buffer capacity~$B$.
The third inequality above is
due to Proposition~\ref{prop: stoch_ord}\eqref{delta_ineq},
which in particular says that
$\{L^{\Pi_2}(t)\}_{t \geq 0} \geq \{L^{\Pi_1}(t)\}_{t \geq 0} \geq
\{L^{\Pi}(t)\}_{t \geq 0}$ almost surely under the coupling.
Now we have the following lemma which we will prove below.

\begin{lemma}
\label{lem: tight}
For all $t\geq 0$, under the assumptions of Theorem~\ref{diffusionjiq},
$\{L^{\Pi_2}(t)\}_{N \geq 1}$ forms a tight sequence.
\end{lemma}

Since $L^{\Pi_2}(t)$ is non-decreasing in~$t$, the above lemma
in particular implies that
\begin{equation}
\label{eq: conv 0}
\sup_{t\in[0,T]}\frac{L^{\Pi_2}(t)}{\sqrt{N}}\pto 0.
\end{equation}
For any scheme $\Pi\in\Pi^{(N)}$, from \eqref{eq: bound} we know that
$$\{Q_i^{\Pi_2}(t) - 2 L^{\Pi_2}(t)\}_{t \geq 0} \leq
\{Q_i^{\Pi}(t)\}_{t \geq 0} \leq \{Q_i^{\Pi_2}(t) +
2 L^{\Pi_2}(t)\}_{t \geq 0}.$$
Combining~\eqref{eq: bound} and~\eqref{eq: conv 0} shows that if the
weak limits under the $\sqrt{N}$ scaling exist, they must be the same
for all the schemes in the class $\Pi^{(N)}$.
Also, as described in Section~\ref{sec:powerofd}, the weak limit
for $\Pi(N, N)$ exists and the common weak limit can be described
by the unique solution of the SDEs in~\eqref{eq:diffusionjsq}.
Hence, the proof of Theorem~\ref{diffusionjiq} is complete.
\end{proof}

\subsection{Multiple dispatchers}
\label{multiple}

\newcommand{\Nb}{\bar{x}_0}

\begin{figure}\centering
\includegraphics[width=.5\linewidth]{modelfigure.pdf}
\caption{Schematic view of the model with $R$ dispatchers and $N$ servers.}
\label{modelfigure}
\end{figure}
So far we have focused on a basic scenario with a single dispatcher,
but it is not uncommon for LBAs to operate across multiple dispatchers.
While the presence of multiple dispatchers does not affect the
queueing dynamics of JSQ($d$) policies, it does matter for the JIQ
scheme which uses memory at the dispatcher.
In order to examine the impact, we consider in this subsection
a scenario with $N$ parallel identical servers as before
and $R \geq 1$ dispatchers, as depicted in Figure~\ref{modelfigure}.
Tasks arrive at dispatcher~$r$ as a Poisson process of rate
$\alpha_r \lambda N$, with $\alpha_r > 0$, $r = 1, \dots, R$,
$\sum_{r = 1}^{R} \alpha_r = 1$, and $\lambda$ denoting the task
arrival rate per server.
For conciseness, we denote $\alpha = (\alpha_1, \dots, \alpha_R)$,
and without loss of generality we assume that the dispatchers are
indexed such that $\alpha_1 \geq \alpha_2 \geq \dots \geq \alpha_R$.

When a server becomes idle, it sends a token to one of the
dispatchers selected uniformly at random, advertising its availability.
When a task arrives at a dispatcher which has tokens available,
one of the tokens is selected, and the task is immediately forwarded
to the corresponding server.

We distinguish two scenarios when a task arrives at a dispatcher
which has no tokens available, referred to as the {\em blocking\/}
and {\em queueing\/} scenario respectively.
In the blocking scenario, the incoming task is blocked
and instantly discarded.
In the queueing scenario, the arriving task is forwarded to one
of the servers selected uniformly at random.
If the selected server happens to be idle, then the outstanding
token at one of the other dispatchers is revoked.

In the queueing scenario we assume $\lambda < 1$, which is not only
necessary but also sufficient for stability.
It is not difficult to show that the joint queue length process is
stochastically majorized by a scheme that assigns each task to a server chosen uniformly at random.
In the latter case, the system decomposes into $N$~independent
M/M/1 queues, each of which has load $\lambda<1$ and is stable.

Scenarios with multiple dispatchers have received limited attention
in the literature, and the scant papers that exist
\cite{LXKGLG11,Mitzenmacher16,Stolyar17} almost exclusively assume
that the loads at the various dispatchers are strictly equal, i.e.,
$\alpha_1 = \dots = \alpha_R = 1 / R$.
In these cases the fluid limit, for suitable initial states, is the
same as in Equation~\eqref{eq:fluidjiq} for a single dispatcher,
and in particular the fixed point is the same,
hence, the JIQ scheme continues to achieve asymptotically optimal
delay performance with minimal communication overhead.
The results in~\cite{Stolyar17} in fact show that the JIQ scheme
remains asymptotically optimal even when the servers are heterogeneous,
while it is readily seen that JSQ($d$) policies cannot even be
maximally stable in that case for any fixed value of~$d$.
As one of the few exceptions, van der Boor {\em et al.} \cite{BBL17a}
allows the loads at the various dispatchers to be different.
It is not uncommon for such skewed load patterns to arise for example
when the various dispatchers receive tasks from external sources
making it difficult to perfectly balance the task arrival rates.

\paragraph{Results for blocking scenario.}
For the blocking scenario, Denote by $B(R, N, \lambda, \alpha)$ the steady-state blocking
probability of an arbitrary task.
It is established in \cite{BBL17a} that,
\begin{equation*}
B(R,N,\lambda,\alpha) \to \max\{1-R\alpha_R,1-1/\lambda\}
\mbox{ as } N \to \infty
\end{equation*}
This result shows that in the many-server limit the system performance
in terms of blocking is either determined by the relative load of the
least-loaded dispatcher, or by the aggregate load.
This may be informally explained as follows.
Let $\Nb$ be the expected fraction of busy servers in steady state,
so that each dispatcher receives tokens on average at a rate $\Nb N/R$.
We distinguish two cases, depending on whether a positive
fraction of the tokens reside at the least-loaded dispatcher~$R$
in the limit or not.
If that is the case, then the task arrival rate $\alpha_R \lambda N$
at dispatcher~$R$ must equal the rate $\Nb N / R$
at which it receives tokens, i.e., $\Nb / R = \alpha_R \lambda$.
Otherwise, the task arrival rate $\alpha_R \lambda N$ at
dispatcher~$R$ must be no less the rate $\Nb N / R$ at which it
receives tokens, i.e., $\Nb / R \leq \alpha_R \lambda$.
Since dispatcher~$R$ is the least-loaded, it then follows that
$\Nb / R \leq \alpha_r \lambda$ for all $r = 1, \dots, R$,
which means that the task arrival rate at all the dispatchers is
higher that the rate at which tokens are received.
Thus the fraction of tokens at each dispatcher is zero in the limit,
i.e., the fraction of idle servers is zero, implying $\Nb = 1$.
Combining the two cases, and observing that $\Nb \leq 1$,
we conclude $\Nb = \min\{R \alpha_R \lambda, 1\}$.
Because of Little's law, $\Nb$ is related to the blocking
probability~$B$ as $\Nb = \lambda (1 - B)$.
This yields $1 - B = \min\{R \alpha_R \lambda, 1 / \lambda\}$,
or equivalently, $B = \max\{1 - R \alpha_R, 1 - 1 / \lambda\}$.

The above explanation also reveals that, somewhat counter-intuitively,
it is the least-loaded dispatcher that throttles tokens and leaves
idle servers stranded, thus acting as bottleneck.
Specifically, in the limit dispatcher~$R$ (or the set of least-loaded
dispatchers in case of ties) inevitably ends up with all the available
tokens, if any.
The accumulation of tokens hampers the visibility of idle servers
to the heavier-loaded dispatchers, and leaves idle servers stranded
while tasks queue up at other servers.

\paragraph{Results for queueing scenario.}
For the queueing scenario, denote by $W(R, N, \lambda, \alpha)$ a random variable
with the steady-state waiting-time distribution of an arbitrary task.
It is shown in~\cite{BBL17a} that, for a fixed
$\lambda < 1$ and $N \to \infty$,
\[
\mathbb{E}[W(R, N, \lambda, \alpha)] \to
\frac{\lambda_2(R, \lambda, \alpha)}{1 - \lambda_2(R, \lambda, \alpha)},
\]
where
\[
\lambda_2(R,\lambda,\alpha)=1-\frac{1-\lambda\sum_{i=1}^{r^*}\alpha_i}{1-\lambda r^* / R}
\]
with
\[
r^*=\sup\big\{r \big| \alpha_r > \frac{1}{R}\frac{1-\lambda\sum_{i=1}^{r}\alpha_i}{1-\lambda r/R}\big\}
\]
and the convention that $r^*=0$ if $\alpha_1=\hdots=\alpha_R=1/R$.
In particular,
\[
\lambda_2(2, \lambda, (1 - \alpha_2, \alpha_2)) =
\lambda \frac{1 - 2 \alpha_2}{2 - \lambda},
\]
so that
\[
\mathbb{E}[W(2, N, \lambda, (1 - \alpha_2, \alpha_2))] \rightarrow
\frac{\lambda (1-2\alpha_2)}{2-2\lambda (1-\alpha_2)}.
\]
Here $\lambda_2$ can be interpreted as the rate at which tasks are forwarded to randomly selected servers. Furthermore, dispatchers $1, \hdots, r^*$ receive tokens at a lower rate than the incoming tasks, and in particular $\lambda_2^*=0$ if and only if $r^*=0$.

When the arrival rates at all dispatchers are strictly equal, i.e.,
$\alpha_1 = \dots = \alpha_R = 1 / R$, the above results indicate that
the stationary blocking probability and the mean waiting time
asymptotically vanish as $N \to \infty$, which is in agreement with
the observations in~\cite{Stolyar17} mentioned above.
However, when the arrival rates at the various dispatchers are not
perfectly equal, so that $\alpha_R < 1 / R$, the blocking probability
and mean waiting time are strictly positive in the limit, even for arbitrarily
low overall load and an arbitrarily small degree of skewness in the
arrival rates.
Thus, the ordinary JIQ scheme fails to achieve asymptotically optimal
performance for heterogeneous dispatcher loads.

\paragraph{Enhancements.}

In order to counter the above-described performance degradation for
asymmetric dispatcher loads,  \cite{BBL17a} proposes two enhancements.

\newtheorem{myenhancement}{Enhancement}

\begin{myenhancement}[Non-uniform token allotment]\label{alg1}
When a server becomes idle, it sends a token to dispatcher~$r$
with probability~$\beta_r$.
\end{myenhancement}

\begin{myenhancement}[Token exchange mechanism]\label{alg2}
Any token is transferred to a uniformly randomly selected
dispatcher at rate~$\nu$.
\end{myenhancement}

Note that the token exchange mechanism only creates a constant
communication overhead per task as long as the rate~$\nu$ does not
depend on the number of servers~$N$, and thus preserves the
scalability of the basic JIQ scheme.
The above enhancements can achieve asymptotically optimal performance for
suitable values of the $\beta_r$ parameters and the exchange rate~$\nu$.
Specifically, the stationary blocking probability in the blocking
scenario and the mean waiting time in the queueing scenario
asymptotically vanish as $N \to \infty$, upon using Enhancement~\ref{alg1}
with $\beta_r = \alpha_r$ or Enhancement~\ref{alg2}
with $\nu \geq \frac{\lambda}{1 - \lambda}(\alpha_1 R - 1)$.
Alternative options to deal with heterogeneous loads could be to have
idle servers issue copies of their availability tokens to multiple
dispatchers \cite{LXKGLG11}.

\paragraph{Large number of dispatchers.}

In the above set-up we assumed the number of dispatchers to remain fixed
as the number of servers grows large, but a further natural scenario
would be for the number of dispatchers $R(N)$ to scale with the number
of servers as considered by Mitzenmacher~\cite{Mitzenmacher16}.
He analyzes the case $R(N) = r N$ for some constant~$r$, so that the
relative load of each dispatcher is $\lambda r$.
The term `I-queue' is used for the queue of (idle) servers that is known by one of the dispatchers. A server is added to an I-queue when it becomes idle. With fluid limits and fixed-point analyses, Mitzenmacher determines the fraction of I-queues with $i$ queued servers and the fraction of servers with $i$ tasks in queue that are in the $j$th position in one of the I-queues. The fixed point can be computed numerically and is accurate.\\

\noindent
{\bf Anticipation}. Mitzenmacher also proposes to have servers issue
their availability tokens to the dispatchers already before they are idle,
e.g. when they have just one task remaining.
This appears beneficial at very high load when there are (on average)
fewer idle servers than dispatchers, and tasks would frequently be
assigned to uniformly at random selected servers otherwise.
Mitzenmacher introduces two variants. First, an LCFS-scheme in which the server that is in the I-queue the least amount of time is chosen for the incoming task. Second, a server that became idle, may probe $d$ I-queues after which it chooses the least loaded one of the $d$. Both variants lead to small performance improvements.

\subsection{Joint load balancing and auto-scaling}

Besides delay performance and implementation overhead, a further key
attribute in the context of large-scale cloud networks and data centers
is energy consumption.
So-called auto-scaling algorithms have emerged as a popular mechanism
for adjusting service capacity in response to varying demand levels
so as to minimize energy consumption while meeting performance targets,
but have mostly been investigated in settings with a centralized queue,
and queue-driven auto-scaling techniques have been widely investigated
in the literature \cite{ALW10,GDHS13,LCBWGWMH12,LLWLA11a,LLWLA11b,LLWA12,LWAT13,PP16,UKIN10,WLT12}.
In systems with a centralized queue it is very common to put servers
to `sleep' while the demand is low, since servers in sleep mode
consume much less energy than active servers.
Under Markovian assumptions, the behavior of these mechanisms can
be described in terms of various incarnations of M/M/$N$ queues
with setup times.
There are several further recent papers which examine on-demand server
addition/removal in a somewhat different vein~\cite{PS16,NS16}.

Unfortunately, data centers and cloud networks with massive numbers
of servers are too complex to maintain any centralized queue,
as it involves a prohibitively high communication burden to obtain
instantaneous state information.

Motivated by these observations, Mukherjee {\em et al.}~\cite{MDBL17}
propose a joint load balancing and auto-scaling strategy,
which retains the excellent delay performance and low implementation
overhead of the ordinary JIQ scheme,
and at the same time minimizes the energy consumption.
The strategy is referred to as TABS (Token-Based Auto-Balance Scaling)
and operates as follows:
\begin{itemize}
\item When a server becomes idle, it sends a `green' message to the dispatcher, waits for an $\exp(\mu)$ time (standby period), and turns itself off by sending a `red' message to the dispatcher (the corresponding green message is destroyed).
\item When a task arrives, the dispatcher selects a green message at random if there are any, and assigns the task to the corresponding server (the corresponding green message is replaced by a `yellow' message).
Otherwise, the task is assigned to an arbitrary busy  server, and if at that arrival epoch there is a red message at the dispatcher, then it selects one at random, and the setup procedure of the corresponding server is initiated, replacing its red message by an `orange' message.
Setup procedure takes $\exp(\nu)$ time after which the server becomes active.
\item Any server which activates due to the latter event, sends a green message to the dispatcher (the corresponding orange message is replaced), waits for an $\exp(\mu)$ time for a possible assignment of a task, and again turns itself off by sending a red message to the dispatcher.
\end{itemize}

\begin{figure}
\begin{center}
\includegraphics[scale=1]{flowchart.pdf}
\end{center}
\caption{Illustration of server on-off decision rules in the TABS scheme,
along with message colors and state variables.}
\label{fig:scheme}
\end{figure}

The TABS scheme gives rise to a distributed operation in which servers
are in one of four states (busy, idle-on, idle-off or standby),
and advertize their state to the dispatcher via exchange of tokens.
Figure~\ref{fig:scheme} illustrates this token-based exchange protocol.
Note that setup procedures are never aborted and continued even when
idle-on servers do become available.
Very recently dynamic scaling and load balancing with variable service
capacity and on-demand agents has been further examined in~\cite{GFPJ17}.

To describe systems under the TABS scheme,
we use $\mathbf{Q}^N(t) := (Q_1^N(t), Q_2^N(t), \dots, Q_B^N(t))$
to denote the system occupancy state at time~$t$ as before

Also, let $\Delta_0^N(t)$ and $\Delta_1^N(t)$ denote the number
of idle-off servers and servers  in setup mode at time $t$, respectively.
The fluid-scaled quantities are denoted by the respective small letters,
\emph{viz.} $q_i^{N}(t) := Q_i^{N}(t) / N$,
$\delta_0^N(t) = \Delta_0^N(t) / N$,
and $\delta_1^N(t) = \Delta_1^N(t)/N$.
For brevity in notation, we will write
$\mathbf{q}^N(t) = (q_1^N(t),\dots,q_B^N(t))$
and $\boldsymbol{\delta}^N(t) = (\delta_0^N(t), \delta_1^N(t))$.

\paragraph{Fluid limit.}
Under suitable initial conditions, on any finite time interval,
with probability~$1$, any sequence $\{N\}$ has a further subsequence
along which the sequence of processes $(\qq^N(\cdot), \dd^N(\cdot))$
converges to a deterministic limit $(\qq(\cdot), \dd(\cdot))$
that satisfies the following system of ODEs
\begin{equation}
\label{eq:tabsfluid}
\begin{split}
\frac{\dif^+ q_i(t)}{\dif t} &= \lambda(t) p_{i-1}(\qq(t),\dd(t),\lambda(t))
 - (q_i(t) - q_{i+1}(t)),\ i = 1, \ldots, B, \\

\frac{\dif^+\delta_0(t)}{\dif t} &= u(t) - \frac{\dif^+\xi(t)}{\dif t}, \qquad
\frac{\dif^+\delta_1(t)}{\dif t} = \frac{\dif^+\xi(t)}{\dif t} - \nu \delta_1(t),
\end{split}
\end{equation}
where by convention $q_{B+1}(\cdot) \equiv 0$, and
\begin{align*}
u(t) &= 1 - q_1(t) - \delta_0(t) - \delta_1(t), \quad
\frac{\dif^+\xi(t)}{\dif t} =
\lambda(t) (1 - p_0(\qq(t), \dd(t), \lambda(t))) \ind{\delta_0(t) > 0}.
\end{align*}
For any $(\qq, \dd)$ and $\lambda > 0$,
$(p_i(\qq, \dd, \lambda))_{i \geq 0}$ are given by
\begin{align*}
p_0(\qq,\dd,\lambda) &=
\begin{cases}
& 1 \qquad \text{ if } \qquad u = 1 - q_1 - \delta_0 - \delta_1 > 0, \\
& \min\{\lambda^{-1}(\delta_1 \nu + q_1 - q_2), 1\},\quad \text{otherwise,}
\end{cases} \\
\quad p_i(\qq, \dd, \lambda) &=
(1 - p_0(\qq, \dd, \lambda)) (q_i - q_{i+1})q_1^{-1},\ i = 1, \ldots, B.
\end{align*}

We now provide an intuitive explanation of the fluid limit stated above.
The term $u(t)$ corresponds to the asymptotic fraction of idle-on
servers in the system at time~$t$, and $\xi(t)$ represents the
asymptotic cumulative number of server setups (scaled by~$N$)
that have been initiated during $[0,t]$.
The coefficient $p_i(\qq, \dd, \lambda)$ can be interpreted as the
instantaneous fraction of incoming tasks that are assigned to some
server with queue length~$i$, when the fluid-scaled occupancy state is
$(\qq, \dd)$ and the scaled instantaneous arrival rate is~$\lambda$.

Observe that as long as $u > 0$, there are idle-on servers,
and hence all the arriving tasks will join idle servers.
This explains that if $u > 0$, $p_0(\qq, \dd, \lambda) = 1$
and $p_i(\qq, \dd, \lambda)=0$ for $i = 1, \ldots, B-1$.
If $u = 0$, then observe that servers become idle at rate $q_1 - q_2$,
and servers in setup mode turn on at rate $\delta_1 \nu$.
Thus the  idle-on servers are created at a total rate
$\delta_1 \nu + q_1 - q_2$.
If this rate is larger than the arrival rate~$\lambda$, then almost
all the arriving tasks can be assigned to idle servers.
Otherwise, only a fraction $(\delta_1 \nu + q_1 - q_2) / \lambda$
of arriving tasks join idle servers.
The rest of the tasks are distributed uniformly among busy servers,
so a proportion $(q_i - q_{i+1}) q_1^{-1}$ are assigned to servers
having queue length~$i$.
For any $i = 1, \ldots, B$, $q_i$ increases when there is an arrival
to some server with queue length $i - 1$, which occurs at rate
$\lambda p_{i-1}(\qq, \dd, \lambda)$, and it decreases when there is
a departure from some server with queue length~$i$, which occurs
at rate $q_i - q_{i-1}$.
Since each idle-on server turns off at rate~$\mu$, the fraction
of servers in the off mode increases at rate $\mu u$.
Observe that if $\delta_0 > 0$, for each task that cannot be assigned
to an idle server, a setup procedure is initiated  at one idle-off server.
As noted above, $\xi(t)$ captures the (scaled) cumulative number
of setup procedures initiated up to time~$t$.
Therefore the fraction of idle-off servers and the fraction of servers
in setup mode decreases and increases by $\xi(t)$, respectively,
during $[0,t]$.
Finally, since each server in setup mode becomes idle-on at rate~$\nu$,
the fraction of servers in setup mode decreases at rate $\nu \delta_1$.

\paragraph{Fixed point and global stability.}
In case of a constant arrival rate $\lambda(t) \equiv \lambda < 1$,
any fluid sample path in~\eqref{eq:tabsfluid} has a unique fixed point:
\begin{equation}
\label{eq:fixed point}
\delta_0^* = 1 - \lambda, \qquad \delta_1^* = 0, \qquad
q_1^* = \lambda \quad \mbox{ and } \quad q_i^*=0,
\end{equation}
for $i = 2, \ldots, B$.
Indeed, it can be verified that $p_0(\qq^*, \dd^*, \lambda) = 1$
and $u^* = 0$ for $(\qq^*, \dd^*)$ given by~\eqref{eq:fixed point}
so that the derivatives of $q_i$, $i = 1, \dots, B$, $\delta_0$,
and $\delta_1$ become zero, and that these cannot be zero at any other
fluid-scaled occupancy state.
Note that, at the fixed point, a fraction $\lambda$ of the servers
have exactly one task while the remaining fraction have zero tasks,
independently of the values of the parameters~$\mu$ and~$\nu$.

In order to establish the convergence of the sequence of steady states,
we need the global stability of the fluid limit, i.e., starting
from any fluid-scaled occupancy state, any fluid sample path described
by~\eqref{eq:tabsfluid} converges to the unique fixed
point~\eqref{eq:fixed point} as $t \to \infty$.
More specifically, irrespective of the starting state,
\begin{equation}
\label{eq:globalstab}
(\qq(t), \dd(t)) \to (\qq^*, \dd^*), \quad \mbox{ as } \quad t \to \infty,
\end{equation}
where $(\qq^*, \dd^*)$ is as defined in~\eqref{eq:fixed point}.

\paragraph{Interchange of limits.}
The global stability can be leveraged to show that the steady-state
distribution of the $N$-th system, for large $N$, can be
well approximated by the fixed point of the fluid limit
in~\eqref{eq:fixed point}.
Specifically, it justifies the interchange of the many-server
($N \to \infty$) and stationary ($t \to \infty$) limits.
Since the buffer capacity~$B$ at each server is supposed to be finite,
for every~$N$, the Markov process $(\QQ^N(t), \Delta_0^N(t), \Delta_1^N(t))$ is irreducible, has a finite state space, and thus has a unique
steady-state distribution.
Let $\pi^N$ denote the steady-state distribution of the
$N$-th system, i.e.,
$$\pi^{N}(\cdot) = \lim_{t \to \infty}
\mathbb{P}\ \big(\qq^{N}(t) = \cdot, \dd^N(t) = \cdot\big).$$
The fluid limit result and the global stability thus yield that
$\pi^N$ converges weakly to~$\pi$ as $N \to \infty$, where $\pi$ is
given by the Dirac mass concentrated upon $(\qq^*,\dd^*)$ defined
in~\eqref{eq:fixed point}.
\begin{remark}\normalfont
Note that the above interchange of limits result was obtained under the assumption that the queues have finite buffers, and analysis of the infinite-buffer scenario was left open.
The key challenge in the latter case stems from the fact that the system stability under the usual subcritical load assumption is not automatic.
In fact as explained in~\cite{MS18}, when the number of servers $N$ is fixed, the stability may {\em not} hold even under a subcritical load assumption.
In~\cite{MS18} the stability issue of the TABS scheme has been addressed and the convergence of the sequence of steady states was shown for the infinite-buffer scenario.
In particular, it was established that for a fixed choice of parameters $\lambda<1$, $\mu>0$, and $\nu>0$, the system with $N$ servers under the TABS scheme is stable for large enough $N$.
There the authors introduce an induction-based approach that uses both the conventional fluid limit (in the sense of a large starting state) and the mean-field fluid limit (when $N\to\infty$) in an intricate fashion to prove the large-$N$ stability of the system.
\end{remark}

\newcommand{\expt}{\mathbb{E}}

\paragraph{Performance metrics.}

As mentioned earlier, two key performance metrics are the expected
waiting time of tasks $\expt[W^N]$ and energy consumption $\expt[P^N]$
for the $N$-th system in steady state.
In order to quantify the energy consumption, we assume that the energy
usage of a server is $P_{\full}$ when busy or in set-up mode,
$P_{\idle}$ when idle-on, and zero when turned off.
Evidently, for any value of~$N$, at least a fraction $\lambda$ of the
servers must be busy in order for the system to be stable,
and hence $\lambda P_{\full}$ is the minimum mean energy usage per
server needed for stability.
We will define $\expt[Z^N] = \expt[P^N] - \lambda P_{\full}$ as the
relative energy wastage accordingly.
The interchange of limits result can be leveraged to obtain that
asymptotically the expected waiting time and energy consumption
for the TABS scheme vanish in the limit, for any strictly positive
values of~$\mu$ and~$\nu$.
More specifically, for a constant arrival rate
$\lambda(t) \equiv \lambda < 1$, for any $\mu > 0$, $\nu > 0$,
as $N \to \infty$,
\begin{enumerate}[{\normalfont (a)}]
\item Zero mean waiting time: $\expt[W^N] \to 0$,
\item Zero energy wastage: $\expt[Z^N] \to 0$.
\end{enumerate}
The key implication is that the TABS scheme, while only involving
constant communication overhead per task, provides performance
in a distributed setting that is as good at the fluid level as can
possibly be achieved, even in a centralized queue, or with unlimited
information exchange.

\paragraph{Comparison to ordinary JIQ policy.}
Consider again a constant arrival rate $\lambda(t) \equiv \lambda$.
It is worthwhile to observe that the component $\qq$ of the fluid limit
as in~\eqref{eq:tabsfluid} coincides with that for the ordinary JIQ
policy where servers always remain on, when the system following the
TABS scheme starts with all the servers being idle-on,
and $\lambda + \mu < 1$.
To see this, observe that the component $\qq$ depends on~$\dd$ only
through $(p_{i-1}(\qq, \dd))_{i \geq 1}$.
Now, $p_0 = 1$, $p_i = 0$, for all $i \geq 1$,
whenever $q_1 + \delta_0 + \delta_1 < 1$, irrespective of the precise
values of $(\qq, \dd)$.
Moreover, starting from the above initial state, $\delta_1$ can
increase only when $q_1 + \delta_0 = 1$.
Therefore, the fluid limit of $\qq$ in~\eqref{eq:tabsfluid} and the
ordinary JIQ scheme are identical if the system parameters
$(\lambda, \mu, \nu)$ are such that $q_1(t) + \delta_0(t) < 1$,
for all $t \geq 0$.
Let $y(t) = 1 - q_1(t) - \delta_0(t)$.
The solutions to the differential equations
\[
\frac{\dif q_1(t)}{\dif t} = \lambda - q_1(t), \quad
\frac{\dif y(t)}{\dif t} = q_1(t) - \lambda - \mu y(t),
\]
$y(0) = 1$, $q_1(0) = 0$ are given by
\[
q_1(t) = \lambda (1 - \e^{- t}), \quad
y(t) = \frac{\e^{- (1+\mu) t}}{\mu-1} \big(\e^t(\lambda+\mu-1) -
\lambda \e^{\mu t}\big).
\]
Notice that if $\lambda + \mu<1$, then $y(t) > 0$ for all $t \geq 0$
and thus, $q_1(t) + \delta_0(t) < 1$, for all $t \geq 0$.
The fluid-level optimality of the JIQ scheme was described
in Subsection~\ref{ssec:fluidjiq}.
This observation thus establishes the optimality of the fluid-limit
trajectory under the TABS scheme for suitable parameter values
in terms of response time performance.
From the energy usage perspective, under the ordinary JIQ policy,
since the asymptotic steady-state fraction of busy servers ($q_1^*$)
and idle-on servers are given by $\lambda$ and $1-\lambda$, respectively,
the asymptotic steady-state (scaled) energy usage is given by
\begin{align*}
\expt[P^{\mathrm{JIQ}}] = \lambda P_{\full} + (1-\lambda) P_{\idle} =
\lambda P_{\full}(1 + (\lambda^{-1}-1)f),
\end{align*}
where $f = P_{\idle}/P_{\full}$ is the relative energy consumption
of an idle server.
As described earlier, the asymptotic steady-state (scaled) energy
usage under the TABS scheme is $\lambda P_{\full}$.
Thus the TABS scheme reduces the asymptotic steady-state energy usage
by $\lambda P_{\full}(\lambda^{-1}-1)f = (1-\lambda)P_{\idle}$, which
amounts to a relative saving of $(\lambda^{-1}-1)f/(1+(\lambda^{-1}-1)f)$.
In summary, the TABS scheme performs as good as the ordinary JIQ policy
in terms of the waiting time and communication overhead while providing
a significant energy saving.

\section{Redundancy policies and alternative scaling}
\label{miscellaneous}

In this section we discuss somewhat related redundancy policies,
 alternative scaling regimes, and some additional performance metrics of interest.

\subsection{Redundancy-d policies}

So-called redundancy-$d$ policies involve a somewhat similar operation
as JSQ($d$) policies, and also share the primary objective of ensuring
low delays \cite{AGSS13,VGMSRS13}.
In a redundancy-$d$ policy, $d \geq 2$ candidate servers are selected
uniformly at random (with or without replacement) for each arriving task,
just like in a JSQ($d$) policy.
Rather than forwarding the task to the server with the shortest queue
however, replicas are dispatched to all sampled servers.
Note that the initial replication to $d$ servers selected uniformly
at random does not entail any communication burden, but the abortion
of redundant copies at a later stage does involve a significant amount
of information exchange and complexity.

Two common options can be distinguished for abortion of redundant clones.
In the first variant, as soon as the first replica starts service,
the other clones are abandoned.
In this case, a task gets executed by the server which had the smallest
workload at the time of arrival (and which may or may not have had the
shortest queue length) among the sampled servers.
This may be interpreted as a power-of-$d$ version of the Join-the-Smallest
Workload (JSW) policy discussed in Subsection~\ref{ssec:jsw}.
The optimality properties of the JSW policy mentioned in that subsection
suggest that redundancy-$d$ policies should outperform JSQ($d$) policies,
which appears to be supported by simulation experiments,
but has not been established by analytical comparisons so far.

In the second option the other clones of the task are not aborted
until the first replica has completed service (which may or may not
have been the first replica to start service).
While a task is only handled by one of the servers in the former case,
it may be processed by several servers in the latter case.
When the service times are exponentially distributed and independent for
the various clones, the aggregate amount of time spent by all the servers
until completion remains exponentially distributed with the same mean.
An exact analysis of the delay distribution in systems with $N = 2$
or $N = 3$ servers is provided in \cite{GZDHHS15,GZDHHS16},
and exact expressions for the mean delay with an arbitrary number
of servers are established in~\cite{GZHS16}.
The limiting delay distribution in a fluid regime with $N \to \infty$
is derived in \cite{GZVHS16,GHSVZ17}.
Based on an asymptotic independence assumption among the servers.
In general, the mean aggregate amount of time devoted to a task
and the resulting delay may be larger or smaller for less or more
variable service time distributions, also depending on the number
of replicas per task \cite{PC16,SLR16,WJW14,WJW15}.
In particular, for heavy-tailed service time distributions, the mean
aggregate time spent on a task may be considerably reduced by virtue
of the redundancy.
Indeed, even if the first replica to start service has an extremely long
service time, that is not likely to be case for the other clones as well.
In spite of the extremely long service time of the first replica,
it is therefore unlikely for the aggregate amount of time spent on the
task or its waiting time to be large.
This provides a significant performance benefit to redundancy-$d$
policies over JSQ($d$) policies, and has also motivated a strong
interest in adaptive replication schemes \cite{APS17a,Joshi17a, Joshi17b}.

A further closely related model is where $k$ of the replicas need to
complete service, $1 \leq k \leq d$, in order for the task to finish
which is relevant in the context of storage systems with coding
and MapReduce tasks \cite{JSW15a,JSW17}.
The special case where $k = d = N$ corresponds to a classical
fork-join system.

\subsection{Conventional heavy traffic}
\label{classical}

In this subsection we briefly discuss a few asymptotic results for LBAs
in the classical heavy-traffic regime as described in Subsection~\ref{asym}
where the number of servers~$N$ is fixed and the relative load tends
to one in the limit.

The papers \cite{Foschini77,FS78,Reiman84,ZHW95} establish diffusion
limits for the JSQ policy in a sequence of systems with Markovian
characteristics as in our basic model set-up, but where in the $K$-th
system the arrival rate is $K \lambda + \hat\lambda \sqrt{K}$, while
the service rate of the $i$-th server is $K \mu_i + \hat\mu_i \sqrt{K}$,
$i = 1, \dots, N$, with $\lambda = \sum_{i = 1}^{N} \mu_i$,
inducing critical load as $K \to \infty$.
It is proved that for suitable initial conditions the queue lengths are
of the order O($\sqrt{K}$) over any finite time interval and exhibit
a state-space collapse property.
In particular, a properly scaled version of the joint queue length
process lives in a one-dimensional rather than $N$-dimensional space,
reflecting that the various queue lengths evolve in lock-step,
with the relative proportions remaining virtually identical in the limit,
while the aggregate queue length varies.

Atar {\em et al.}~\cite{AKM17} investigate a similar scenario,
and establish diffusion limits for three policies: the JSQ($d$) policy,
the redundancy-$d$ policy (where the redundant clones are abandoned
as soon as the first replica starts service), and a combined policy
called Replicate-to-Shortest-Queues (RSQ) where $d$ replicas are
dispatched to the $d$-shortest queues.
Note that the latter policy requires instantaneous knowledge of all
the queue lengths, and hence involves a similar excessive communication
overhead as the ordinary JSQ policy, besides the substantial information
exchange associated with the abortion of redundant copies.
Conditions are derived for the values of the relative service rates
$\mu_i$, $i = 1, \dots, N$, in conjunction with the diversity
parameter~$d$, in order for the queue lengths under the JSQ($d$)
and redundancy-$d$ policies to be of the order O($\sqrt{K}$) over any
finite time interval and exhibit state-space collapse.
The conditions for the two policies are distinct, but in both cases
they are weaker for larger values of~$d$, as intuitively expected.
While the conditions for the values of~$\mu_i$ depend on~$d$,
whenever they are met, the actual diffusion-scaled queue length
processes do not depend on the exact value of~$d$ in the limit,
showing a certain resemblance with the universality property
as identified in Subsection~\ref{ssec:powerd} for the large-capacity
and Halfin-Whitt regimes.

Zhou {\em et al.} \cite{ZWTSS17} consider a slightly different model
set-up with a time-slotted operation, and identify a class~$\Pi$ of LBAs
that not only provide throughput-optimality (or maximum stability,
i.e., keep the queues stable in a suitable sense whenever feasible to
do so at all), but also achieve heavy-traffic delay optimality, in
the sense that the properly scaled aggregate queue length is the same
as that in a centralized queue where all the resources are pooled
as the load tends to one.
As it turns out, the class~$\Pi$ includes JSQ($d$) policies with $d \geq 2$,
but does \emph{not} include the JIQ scheme, which tends to degenerate
into a random assignment policy when idle servers are rarely available.
The authors further propose a threshold-based policy which has low
implementation complexity like the JIQ scheme, but \emph{does} belong
to the class~$\Pi$, and hence achieves heavy-traffic delay optimality.

\subsection{Non-degenerate slowdown}
\label{nondegenerate}

In this subsection we briefly discuss a few of the scarce asymptotic
results for LBAs in the so-called non-degenerate slow-down regime
described in Subsection~\ref{asym} where $N - \lambda(N) \to \gamma > 0$,
as the number of servers~$N$ grows large.
We note that in a centralized queue the process tracking the evolution
of the number of waiting tasks, suitably accelerated and normalized
by~$N$, converges in this regime to a Brownian motion with drift
$- \gamma$ reflected at zero as $N \to \infty$.

In stationarity, the number of waiting tasks, normalized by~$N$,
converges in this regime to an exponentially distributed random
variable with parameter~$\gamma$ as $N \to \infty$.
Hence, the mean number of waiting tasks must be at least of the order
$N / \gamma$, and the waiting time cannot vanish as $N \to \infty$
under any policy.

Gupta \& Walton~\cite{GW17} characterize the diffusion-scaled queue
length process under the JSQ policy in this asymptotic regime.
They further compare the diffusion limit for the JSQ policy with that
for a centralized queue as described above as well as several LBAs
such as the JIQ scheme and a refined version called Idle-One-First (I1F),
where a task is assigned to a server with exactly one task if no idle
server is available and to a randomly selected server otherwise.

It is proved that the diffusion limit for the JIQ scheme is no longer
asymptotically equivalent to that for the JSQ policy in this asymptotic
regime, and the JIQ scheme fails to achieve asymptotic optimality
in that respect, as opposed to the behavior in the large-capacity
and Halfin-Whitt regimes discussed in Subsection~\ref{ssec:jiq}.
In contrast, the I1F scheme does preserve the asymptotic equivalence
with the JSQ policy in terms of the diffusion-scaled queue length process,
and thus retains asymptotic optimality in that sense.

These results provide further indication that the amount and accuracy
of queue length information needed to achieve asymptotic equivalence
with the JSQ policy depend not only on the scale dimension (e.g.~fluid
or diffusion), but also on the load regime.
Put differently, the finer the scale and the higher the load,
the more strictly one can distinguish various LBAs in terms
of the relative performance compared to the JSQ policy.

\subsection{Sparse-feedback regime}
\label{sparse}

\iffalse

	Load balancing algorithms play a vital role in enhancing performance in data centers and cloud networks.
	Due to the massive size of these systems, scalability challenges,
	and especially the communication overhead associated with load
	balancing mechanisms, have emerged as major concerns.
	Motivated by these issues, we introduce and analyze a novel class
	of load balancing schemes where the various servers provide occasional
	queue updates to guide the load assignment.

	We show that the proposed schemes strongly outperform JSQ($d$)
	strategies with comparable overhead per task, and can achieve
	a vanishing waiting time in the many-server limit with just one
	message per task, just like the popular JIQ scheme.
	The proposed schemes are particularly geared however towards the
	sparse feedback regime with $d < 1$ message per task,
	where they outperform corresponding sparsified JIQ versions.

	We investigate fluid limits for synchronous updates as well as
	asynchronous exponential update intervals.
	The fixed point of the fluid limit is identified in the latter case,
	and used to derive the queue length distribution.
	We also demonstrate a remarkable dichotomy in the ultra-low feedback
	regime where the mean waiting time tends to a constant in the
	synchronous case, but grows without bound in the asynchronous case.

\fi

As described in Section~\ref{ssec:jiq}, the JIQ scheme involves
a communication overhead of at most one message per task, and yet
achieves optimal delay performance in the fluid and diffusion regimes.
However, even just one message per task may still be prohibitive,
especially when tasks do not involve big computational tasks,
but small data packets which require little processing.
In such situations the sheer message exchange in providing queue
length information may be disproportionate to the actual amount
of processing required.

Motivated by the above issues, \cite{BBL17b} proposes and examines
a novel class of LBAs which also leverage memory at the dispatcher,
but allow the communication overhead to be seamlessly adapted
and reduced below that of the JIQ scheme.
Specifically, in the proposed schemes, the various servers provide
occasional queue status notifications to the dispatcher,
either in a synchronous or asynchronous fashion.
The dispatcher uses these reports to maintain queue estimates,
and forwards incoming tasks to the server with the lowest queue estimate.
The queue estimate for a server is incremented for every task
assigned, and set to the true queue length at update moments,
but never lowered in between updates.
Note that when the update frequency per server is $\delta$,
the number of messages per task is $d = \delta / \lambda$,
with $\lambda < 1$ denoting the arrival rate per server.

The results in~\cite{BBL17b} demonstrate that the proposed schemes
markedly outperform JSQ($d$) policies with the same number
of $d \geq 1$ messages per task and they can achieve
a vanishing waiting time in the many-server limit when the
update frequency $\delta$ exceeds $\lambda / (1 - \lambda)$.
In case servers only report zero queue lengths and suppress
updates for non-zero queues, the update frequency required for
a vanishing waiting time can in fact be lowered to just~$\lambda$,
matching the one message per task involved in the JIQ scheme.

From a scalability viewpoint, the most pertinent regime is $d < 1$
where only very sparse server feedback is required.
It is shown in~\cite{BBL17b} that the proposed schemes then
outperform the corresponding sparsified versions of the JIQ scheme
where idle servers only provide notifications to the dispatcher
with probability~$d$.
In order to further explore the performance for $d < 1$ in the
many-server limit, \cite{BBL17b} investigates fluid limits for the
synchronous case as well as the asynchronous case with exponential
update intervals.
The fixed point of the fluid limit are leveraged to derive the
stationary queue length distribution as function of the update
frequency.

Additionally, \cite{BBL17b} examines the performance in the
ultra-low feedback regime where the update frequency~$\delta$ goes
to zero, and in particular establishes a somewhat counter-intuitive
dichotomy.
In the synchronous case, the behavior of each of the individual
queues approaches that of a single-server queue with
a near-deterministic arrival process and exponential service times,
with the mean waiting time tending to a finite constant.
In contrast, in the asynchronous case, the individual queues
experience saw-tooth behavior with oscillations and waiting times
that grow without bound.

\iffalse

\begin{enumerate}

\item Power-of-$d$ is a under-performer among these policies
by a wide margin.
To place Power-of-$d$ in the most favorable position,
the communication overhead is assumed to be $d$ in these results.
The true overhead could well be argued to be $2 d$ per task
(query and response), which would shift the Power-of-$d$ curve
to the right by a factor two, and show even worse performance.

\item $\sujsqdi$, where only idle servers send reports, further
decreases the communication overhead.
This scheme achieves a near-zero waiting time with just one message
per task, just like the JIQ scheme, and additionally outperforms
the sparsified version JIQ($p$) across most of the relevant domain.

\item Both $\sujsqd$ and $\sujsqe$ resemble Round-Robin as the
update frequency~$\delta$ approaches zero, and are the clear
winners in the ultra-low feedback regime, which is crucial
from a scalability perspective.

\item While $\aujsqd$ outperforms the synchronous variants for
large values of the update frequency~$\delta$, it produces a mean
waiting time that grows without bound as $\delta$ approaches zero.
The latter issue also occurs for $\sujsqdi$, and render the
asynchronous versions far inferior in the ultra-low feedback regime
compared to both synchronous variants.

\end{enumerate}

In order to interpret the resemblance in 3), we observe that the
initial queue lengths after an update will be small compared
to the number of arrivals until the next update.
Thus shortly after the update the dispatcher will essentially start
forwarding tasks in a (probabilistic) Round-Robin manner.
Specifically, the arrival time of the $i$-th task during an update
interval at a particular server~$n$ will be $T_n + (i + U_i) / \lambda$
for large~$N$, where $T_n$ is an offset determined by the initial
queue length at server~$n$ at the start of the update interval,
and $U_1, U_2, \dots$ are independent standard uniform random variables.
Even though the assignment is not strictly deterministic, the time
between the arrival of the $i$-th and $(i+j)$-th task is $j + U_j - U_i$,
and has a coefficient of variation that vanishes for large~$j$.

In order to explain the remarkable dichotomy in 4), we observe that
the queue estimates must inevitably grow to increasingly large values
of the order $\lambda / \delta$ and significantly diverge from the
true queue lengths as the update frequency becomes small,
both in the synchronous and asynchronous versions.
However, in the synchronous variants the queue estimates will all
be lowered and updated to the true queue lengths simultaneously,
prompting the dispatcher to evenly distribute incoming tasks over time.
In contrast, in the asynchronous versions, a server will be the
only one with a low queue estimate right after an update,
and almost immediately be assigned a huge pile of tasks to bring
its queue estimate at par, resulting in oscillatory effects.
This somewhat counter-intuitive dichotomy reveals that the
synchronous variants behave benignly in the presence of outdated
information, while the asynchronous versions are adversely impacted.

\fi

In order to achieve a vanishing waiting time, the dispatcher must assign each incoming task to an idle server with high probability, and thus be able to identify on average at least one idle server for every incoming task.  When the amount of memory at the dispatcher is limited, the dispatcher may in fact have to identify more idle servers on average to ensure that at least one is available with high probability for each incoming task, as also reflected in the results of \cite{GTZ16}.  These conditions, in conjunction with the fact that the fraction of idle servers in equilibrium is $ 1 - \lambda$, translate into a minimum required communication overhead for various families of algorithms.  For example, if the dispatcher samples a server at random, it will find that server idle with probability $1 - \lambda$, so in the absence of any memory it will need to sample a number of servers that grows with $N$ for each incoming task, while with unlimited memory, it will need to sample on average $1 / (1 - \lambda)$ servers per incoming task.  Likewise, if servers report their queue status to the dispatcher, then an arbitrary server will report to be idle with probability $1 - \lambda$, so they all need to do that every $\lambda / (1 - \lambda)$ time units on average, i.e.,~$1 / (1 - \lambda)$ times on average per incoming task.  When only idle servers report their status to the dispatcher, as in the JIQ algorithm, they only need to do so at most once per incoming task.  When servers report their status asynchronously rather than all simultaneously, or idle servers only after some delay, the associated memory requirement at the dispatcher can be reduced.

\subsection{Scaling of maximum queue length}
\label{ballsbins}

So far we have focused on the asymptotic behavior of LBAs in terms
of the number of servers with a certain queue length, either on fluid
scale or diffusion scale, in various regimes as $N \to \infty$.
A related but different performance metric is the maximum queue length
$M(N)$ among all servers as $N \to \infty$.
Luczak \& McDiarmid~\cite{LM06} showed that
for a fixed $d\geq 2$ the steady-state maximum queue length
$M(N)$ in a system under JSQ($d$) policy is given by $\ln(\ln(N)) / \ln(d) + O(1)$ and is concentrated on at most two adjacent values, whereas
for purely random assignment ($d=1$), it scales as $\ln(N) / \ln(1/\lambda)$ and does not concentrate on a bounded range of values.
This is yet a further manifestation of the power-of-choice effect.

The maximum queue length $M(N)$ is the central performance
metric in balls-and-bins models where arriving items (balls) do not
get served and never depart but simply accumulate in bins,
and (stationary) queue lengths are not meaningful.
In fact, the very notion of randomized load balancing and power-of-$d$
strategies was introduced in a balls-and-bins setting in the seminal
paper by Azar {\em et al.}~\cite{ABKU94}.
Several further variations and extensions in that context have been considered in \cite{V99,ACMR95,BCSV00,BCSV06, CMS95,DM93,FPSS05,PR04, P05, G81}.

As alluded to earlier, there are natural parallels between the
balls-and-bins setup and the queueing scenario that we have focused
on so far.
These commonalities are for example reflected in the fact that
power-of-$d$ strategies yield similar dramatic performance improvements
over purely random assignment in both settings.

However, there are also quite fundamental differences between the
balls-and-bins setup and the queueing scenario, besides the obvious
contrasts in the performance metrics.
The distinction is for example evidenced by the fact that a simple
Round-Robin strategy produces a perfectly balanced allocation in
a balls-and-bins setup but is far from optimal in a queueing scenario
as observed in Subsection~\ref{random}.
In particular, the stationary fraction of servers with two or more
tasks under a Round-Robin strategy remains positive in the limit
as $N \to \infty$, whereas it vanishes under the JSQ policy.
Furthermore, it should also be noted~\cite{LM05} that the maximum number of balls in a bin under the purely random assignment policy scales as $\ln(N)/\ln(\ln(N))$ and is concentrated on two adjacent values, which is again in contrast with the queueing scenario.
On a related account, since tasks get served and eventually depart
in a queueing scenario, less balanced allocations with a large
portion of vacant servers will generate fewer service completions
and result in a larger total number of tasks.
Thus different schemes yield not only various degrees of balance,
but also variations in the aggregate number of tasks in the system,
which is not the case in a balls-and-bins set-up.

\section{Extensions}\label{sec:ext}

Throughout most of the paper we have focused on the supermarket model
as a canonical setup and adopted several common assumptions in that context:
(i) all servers are identical;
(ii) the service requirements are exponentially distributed;
(iii) no advance knowledge of the service requirements is available;
(iv) in particular, the service discipline at each server is oblivious
to the actual service requirements.
As mentioned earlier, the stochastic optimality of the JSQ policy,
and hence its central role as an ideal performance benchmark,
critically rely on these assumptions.
The latter also broadly applies to the stochastic coupling techniques
and asymptotic universality properties that we have considered
in the previous sections.
In this section however we review some results for scenarios
where these assumptions are relaxed, in particular allowing for general
service requirement distributions and possibly heterogeneous servers,
along with some broader methodological issues.
In Subsection~\ref{ssec:powerofd-gen} we focus on the behavior of JSQ($d$)
policies in such scenarios, mainly in the large-$N$ limit,
while also briefly commenting on the JIQ policy.
In Subsection~\ref{anticipating} we discuss strategies which specifically exploit
knowledge of server speeds or service requirements of arriving tasks
in making task assignment decisions, and may not necessarily use
queue length information, mostly in a finite-$N$ regime.

\subsection{JSQ(d) policies with general service requirement
distributions}
\label{ssec:powerofd-gen}

Foss \& Chernova~\cite{FC91,FC98} use direct probabilistic methods
and fluid limits to obtain stability conditions
for finite-size systems with a renewal arrival process, a FCFS discipline
at each server, various state-dependent routing policies, including JSQ,
and general service requirement distributions, which may depend on the
task type, the server or both.
Using fluid limits as well as Lyapunov functions,
Bramson~\cite{Bramson98,Bramson11} shows that JSQ($d$) policies achieve
stability for any subcritical load in finite-size systems
with a renewal arrival process, identical servers, non-idling local
service disciplines and general service requirement distributions.
In addition, he derives uniform bounds on the tails of the marginal
queue length distributions, and uses these to prove relative
compactness of these distributions.

Bramson {\em et al.}~\cite{BLP10,BLP12} examine mean-field limits
for JSQ($d$) policies with generally distributed service requirements,
leveraging the above-mentioned tail bounds and relative compactness.
They establish that similar \emph{power-of-choice} benefits occur
as originally demonstrated for exponentially distributed service
requirements in the work of Mitzenmacher~\cite{Mitzenmacher96}
and Vvedenskaya {\em et al.}~\cite{VDK96}, provided a certain `ansatz'
holds asserting that finite subsets of queues become independent in the
large-$N$ limit.
The latter `propagation of chaos' property is shown to hold in several
settings, e.g. when the service requirement distribution has
a decreasing hazard rate and the discipline at each server is FCFS
or when the service requirement distribution has a finite second moment
and the load is sufficiently low.
The ansatz also always holds for the power-of-$d$ version of the JSW
rather than JSQ policy.

It is further shown in~\cite{BLP10,BLP12} that the arrival process
at any given server tends to a state-dependent Poisson process in the
large-$N$ limit, and that the queue length distribution becomes insensitive
with respect to the service requirement distribution when the service
discipline is either Processor Sharing or LCFS with preemptive resume.
This may be explained from the insensitivity property of queues
with state-dependent Poisson arrivals and symmetric service disciplines.

There are strong plausibility arguments that a similar asymptotic
insensitivity property should hold for the JIQ policy in a queueing
scenario, even if the discipline at each server is not symmetric
but FCFS for example.
So far, however, this has only been rigorously established for service
requirement distributions with decreasing hazard rate in~\cite{Stolyar15}.
This result was in fact proved for systems with heterogeneous server
pools, and was further extended in~\cite{Stolyar17} to systems
with multiple symmetric dispatchers.
As it turns out, general service requirement distributions
with an increasing hazard rate give rise to major technical challenges
due to a lack of certain monotonicity properties.
This has only allowed a proof of the asymptotic zero-wait property
for the JIQ policy for load values strictly below $1/2$ so far~\cite{FS17}.

A fundamental technical issue associated with any general service
requirement distribution is that the joint queue length no longer
provides a suitable state description, and that the state space
required for a Markovian description is no longer countable.
Aghajani \& Ramanan~\cite{AR17} and Aghajani {\em et al.}~\cite{ALR18}
introduce a particle representation for the state of the system
and describe the state dynamics for a JSQ($d$) policy via a sequence
of interacting measure-valued processes.
They prove that as $N$ grows large, a suitably scaled sequence
of state processes converges to a hydrodynamic limit which is
characterized as the unique solution of a countable system of coupled
deterministic measure-valued equations, i.e., a system of PDE rather
than the usual ODE equations.
They also establish a `propagation of chaos' result, meaning that
finite collections of queues are asymptotically independent.

Mukhopadhyay \& Mazumdar~\cite{MM14,MM16}
and Mukhopadhyay {\em et al.}~\cite{MKM16} analyzed the performance
and stability of static probabilistic routing strategies and power-of-$d$
policies in the large-$N$ limit in systems with exponential service
requirement distributions, but heterogeneous server pools
and a Processor-Sharing discipline at each server.
They also considered variants of the JSQ($d$) policy which account
for the server speed in the selection criterion as well as hybrid
combinations of the JSQ($d$) policy with static probabilistic routing.
Related results for heterogeneous loss systems rather than queueing
scenarios are presented in \cite{KMM17,MMG15,MKMG15}.
As the results in~\cite{MM14,MM16} reflect, ordinary JSQ($d$) policies
may fail to sample the faster servers sufficiently often in such scenarios,
and therefore fail to achieve maximum stability, let alone
asymptotic optimality.
In~\cite{MKM16} a weighted version of JSQ($d$) policies is presented that
does provide maximum stability, without requiring any specific knowledge
of the underlying system parameters and server speeds in particular.

Vasantam {\em et al.} \cite{VMM17a,VMM17b} examine mean-field limits
for power-of-$d$ policies in many-server loss systems with phase-type
service requirement distributions.
They observe that the fixed point suggests a similar insensitivity
property of the stationary occupancy distribution as mentioned above.
In view of the insensitivity of loss systems with possibly
state-dependent Poisson arrivals, this may be interpreted as an indirect
indication that the arrival process at any given server pool tends
to a state-dependent Poisson arrival process in the large-$N$ limit.
In a somewhat different strand of work, Jonckheere \& Prabhu~\cite{JP16}
investigate the behavior of blocking probabilities in various load
regimes in systems with many single-server finite-buffer queues,
a Processor-Sharing discipline at each server, and an \emph{insensitive}
routing policy.

\subsection{Heterogeneous servers and knowledge of service requirements}
\label{anticipating}

The bulk of the literature has focused on systems
with identical servers, and scenarios with non-identical server speeds
have received relatively limited attention.
A natural extension of the JSQ policy is to assign jobs to the server
with the normalized shortest queue length, or equivalently,
assuming exponentially distributed service requirements, the shortest
expected delay.
While such a Generalized JSQ (GJSQ) or Shortest Expected Delay (SED)
strategy tends to perform well~\cite{BZ89}, it is not strictly
optimal in general~\cite{EVW80}, and the true optimal strategy may in
fact have a highly complicated structure.
Selen {\em et al.}~\cite{SAK16} present approximations for the
performance of GJSQ policies in a finite-$N$ regime with generally
distributed service requirements and a Processor-Sharing discipline
at each server, extending the analysis in Gupta {\em et al.}~\cite{GHSW07}
for the ordinary JSQ policy with homogeneous servers.

In a separate line of work, Feng {\em et al.}~\cite{FMR05} consider
static dispatching policies in a finite-$N$ regime with heterogeneous
servers and a FCFS or Processor-Sharing discipline at each server.
The assignment decision may depend on the service requirement of the
arriving task, but not on the actual queue lengths or any other
state information.
In case of FCFS the optimal routing policy is shown to have a nested
size interval structure, generalizing the strict size interval
structure of the task assignment strategies in Harchol-Balter
{\em et al.}~\cite{HCM99} which are optimal for homogeneous servers.
In case of Processor Sharing, the knowledge of the service requirements
of arriving tasks is irrelevant, in the absence of any state information.

Altman {\em et al.}~\cite{AAP11} consider static probabilistic routing
policies in a somewhat similar setup of a finite-$N$ regime
with multiple task types, servers with heterogeneous speeds,
and a Processor-Sharing discipline at each server.
The routing probabilities are selected so as to either minimize the
global weighted holding cost or the expected holding cost for
an individual task, and may depend on the type of the task and its
service requirement, but not on any other state information.

When knowledge of the service requirements of arriving tasks is available,
it is natural to exploit that for the purpose of local scheduling
at the various servers, and for example use size-based disciplines.
The impact of the local scheduling discipline and server heterogeneity
on the performance and degree of efficiency of load balancing strategies
is examined in~\cite{CMW09}.
An interesting broader issue concerns the relative benefits provided
by exploiting knowledge of service requirements of arriving tasks versus
using information on queue lengths or workloads at the various servers,
which strongly depend on the service requirement distribution~\cite{HSY09}.

\bibliographystyle{abbrv}
\bibliography{bibliography}

\end{document}

config: arxiv