Uniform stability of some large-scale parallel server networks

Uniform stabilit y of some large-scale parallel serv er net w orks HASSAN HMEDI ∗ , ARI ARAPOST A THIS ∗ , AND GUODONG P ANG † Abstract. In this pap er w e study the uniform stability prop erties of tw o classes of parallel serv er net works with m ultiple classes of jobs and m ultiple serv er po ols of a tree top ology . These include a class of netw orks with a single non-leaf serv er po ol, such as the ‘N’ and ‘M’ mo dels, and netw orks of any tree top ology with class-dep endent service rates. W e sho w that with √ n safety staffing, and no abandonment, in the Halfin–Whitt regime, the diffusion-scaled con trolled queueing pro cesses are exp onen tially ergo dic and their inv ariant probability distributions are tight, uniformly ov er all stationary Mark o v con trols. W e use a unified approac h in which the same Ly apunov function is used in the study of the prelimit and diffusion limit. A parameter called the sp ar e c ap acity (safety staffing) of the network pla ys a central role in c haracterizing the stability results: the parameter b eing p ositive is necessary and sufficient that the limiting diffusion is uniformly exp onentially ergo dic o ver all stationary Marko v controls. W e in tro duce the concept of “ system-wide work c onserving p olicies ”, which are defined as p olicies that minimize the num b er of idle serv ers at all times. This is stronger than the so-called joint work conserv ation. W e show that, pro vided the spare capacity parameter is p ositiv e, the diffusion-scaled pro cesses are geometrically ergo dic and the inv ariant distributions are tight, uniformly ov er all “system-wide work conserving p olicies”. In addition, when the spare capacity is negativ e we show that the diffusion-scaled pro cesses are transient under any stationary Marko v control, and when it is zero, they cannot b e p ositive recurrent. 1. Introduction Large-scale parallel server netw orks hav e b een the sub ject of intense study , due to their use in mo deling a v ariety of systems including telecommunications, data centers, customer services and man ufacturing systems; see, e.g., [ 1 , 9 , 14 , 19 , 20 , 25 , 30 , 31 ]. In suc h netw orks, there are multiple classes of jobs and m ultiple serv er p o ols where each job class can b e serv ed b y a subset of server p o ols while each server p o ol can serve a subset of job classes, thus requiring optimal routing and sc heduling decisions. Man y of these systems op erate in the so-called the Halfin–Whitt regime (or Qualit y-and-Efficiency-Driven (QED) regime [ 12 , 23 , 32 ]), where the arriv al rates and the n umbers of serv ers grow large as the scale of the system gro ws, while the service rates remain fixed in such a w ay that the system b ecomes critically loaded. Ensuring stability of these systems through allo cating av ailable resources b y means of adjusting con troller parameters is of great imp ortance. Existing w ork in the literature has addressed the follo wing imp ortant questions: (i) Uniform stability of the multiclass single-p o ol “V” netw ork. The study in [ 17 ] fo cused on the prelimit diffusion-scaled pro cess and sho wed that, with square-ro ot safet y staffing in the single-p ool of servers, the inv arian t probabilit y distributions under all w ork-conserving sc heduling p olicies are tight, and hav e a uniform exp onential tail when the mo del has no ∗ Dep ar tment of Electrical and Computer Engineering, The University of Texas a t Austin, Austin, TX 78712 † Dep ar tment of Comput a tional and Applied Ma thema tics, George R. Brown College of Engineer- ing, Rice University, Houston, TX 77005 E-mail addresses : { hmedi,ari } @utexas.edu, gdpang@rice.edu . 2000 Mathematics Subje ct Classific ation. Primary: 90B22. Secondary: 60K25; 49L20; 90B36. Key wor ds and phr ases. uniform exp onential ergo dicit y , parallel serv er (multiclass multi-po ol) netw orks, Halfin– Whitt regime, spare capacity , system-wide work conserv ation. 1 2 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG abandonmen t (or a sub-Gaussian tail with abandonment). In [ 4 ], a unified approach with a common Lyapuno v function is dev elop ed to establish a F oster-Lyapuno v equation for both the diffusion limit and the diffusion-scaled pro cesses, which shows that the asso ciated in- v ariant probabilit y measures ha v e exp onen tial tails, uniformly o ver the scale of the net work, and o ver all stationary (work-conserving) Mark ov controls. (ii) Stability of the ‘N’ netw ork under a static priority scheduling p olicy . With safet y staffing in one serv er p o ol and no abandonment, Stoly ar [ 27 ] employ ed a in tegral type of Lyapuno v function and established the tightness of stationary distributions of the diffusion-scaled pro cess (there is no analysis of the rate of conv ergence though). (iii) Counterexamples for stabilit y of m ulti-class multi-po ol netw orks. Stoly ar and Y udo vina [ 29 ] sho wed that the stationary distributions of the diffusion-scaled pro cesses may not b e tigh t in these regimes under a natural load balancing scheduling policy , “Longest-queue freest-serv er” (LQFS-LB) (also true in the underloaded regime). (iv) Stability of multi-class m ulti-p o ol net works with p o ol-dep endent service rates under the LQFS-LB p olicy [ 29 ]. W e also refer the reader to [ 26 , 28 ], even though these concern the underloaded case. (v) Stability of multiclass m ulti-p o ol net works under a family of Mark ov p olicies. In [ 6 , 7 ], it is shown that a class of state-dep endent p olicies, referred to as b alanc e d satur ation p olicies (BSP) are stabilizing for the prelimit diffusion-scaled queueing pro cess, when at least one abandonmen t rate is strictly p ositiv e. (vi) Stability of the limiting controlled diffusions for multiclass m ulti-p o ol netw orks under a constan t con trol. Arap ostathis and P ang [ 5 ] developed a leaf elimination algorithm to deriv e an explicit expression of the drift, and, consequen tly , by using the structural prop erties of the drift, a static priority scheduling and routing control is identified which stabilizes the limiting diffusion, when at least one of the classes has a p ositive abandonmen t rate. (vii) Stabilizability of multiclass m ulti-p o ol netw orks of any tree top ology without abandonment in the Halfin-Whitt regime. Hmedi, Arap ostahis and Pang [ 24 ] iden tified a system-wide safet y staffing parameter and sho wed that that parameter b eing p ositive is a necessary and sufficien t condition for the net work to b e stabilizable, that is, there exists a scheduling p olicy under which the stationary distributions of the controlled diffusion-scaled queueing pro cesses are tight ov er the size of the netw ork. The stabilit y results in (v) and (vi) are used in the aforementioned pap ers for the study of ergodic con trol problems for multiclass m ulti-p o ol net works. In [ 2 , 5 – 7 ], due to the lack of the “uniform stabilit y” (also called “blanket stability”) prop erty , ergo dic control problems were studied using a rather elab orate metho dology . The uniform stabilit y prop erties established in this pap er render the ergo dic control problem muc h simpler, and it can b e studied by applying the metho dology in [ 3 , Chapter 3.7]. Despite all the imp ortan t results in (i)–(vi), the ergo dic prop erties of m ulticlass multi-po ol net works in the Halfin–Whitt regime are far from b eing w ell understo o d. The stability analysis of m ulticlass multi-po ol net works in the Halfin–Whitt regime is considerably more c hallenging than the corresp onding one for the ‘V’ netw ork. The problem is particularly difficult when the system do es not hav e abandonmen t. Giv en the coun terexamples in [ 29 ], uniform stability , that is, tigh tness of the inv ariant probability distributions, do es not hold for multiclass multi-po ol net works of an y tree top ology . In this pap er w e iden tify a large class of such net works that are indeed uniformly stable: (a) netw orks with one dominant server p o ol, that is, a single non-leaf serv er p o ol, whic h include the ‘N’, ‘M’ and generalized ‘N’, ‘M’ net works with diameters equal to three or four, and (b) netw orks with class- dep enden t service rates. It migh t app ear to the reader that the top ology in (a) is restrictive. One should note though that even for simple net works with t wo non-leaf serv er p o ols [ 29 , Figure 2, p. 21] the parameters can b e chosen so that uniform stability fails. UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 3 The classes of netw orks in (a)–(b) share an imp ortan t structural prop ert y in their drift, that is, the matrix B 1 in (2.30) is diagonal. W e establish a necessary and sufficient condition for uniform stabilit y , via the so-called spare capacit y (safet y staffing) parameter defined in (2.35) . F or the net works under consideration, we show that if the spare capacity is negativ e, then the limiting diffusion is transient under any stationary Mark ov con trol, if it is zero, the diffusion cannot b e p ositiv e recurrent, and if it is p ositive, the diffusion limit is uniformly stable ov er all stationary Mark ov controls (see Prop osition 3.1 and Theorem 2.1 ). The analogous results for the diffusion- scaled pro cesses are also established. Lastly , we provide a characterization of the spare capacity parameter for the limiting diffusion when the latter is p ositiv e recurrent. W e show in Theorem 3.1 that the spare capacity is equal to an av erage ‘idleness’ weigh ted by the critical quan tity in (3.1) . T o prov e the uniform exp onential ergo dicity for the limiting con trolled diffusion, we use a com- mon Lyapuno v function giv en in (4.7) . This Ly apunov function consists of t wo comp onen ts that treat the p ositiv e and negative half spaces of the state space in a delicate manner. An imp or- tan t ‘tilting’ parameter must b e carefully c hosen to account for not only the differen t effects of queueing and idleness (p ositiv e and negativ e half state space), but also the second order deriv a- tiv es of the extended generator of the diffusion. Note that these Lyapuno v functions differ from the quadratic Lyapuno v functions used in [ 5 – 8 , 16 ] for the study of stability under either constan t con trols or assuming abandonmen t, and also differ from that used in [ 4 ] for the uniform stabilit y of the ‘V’ netw ork. In [ 16 ] for example, the stability analysis requires the existence of a common quadratic Lyapuno v function which cannot be sho wn for the mo dels under consideration. As it will b e clear to the reader later in the pap er, the uniform stability analysis without abandonment requires the use of the sum of tw o functions where eac h of them ‘dominates’ the other ov er a part of the state space. See for example the pro ofs of Lemmas 4.1 and 4.2 . The same Lyapuno v function is used to pro ve the uniform exp onen tial ergo dicity for the prelimit diffusion-scaled pro cesses. How ever, unlik e the ‘V’ netw ork studied in [ 4 ], the F oster-Lyapuno v equations for the limiting diffusion do not carry ov er to the analogous equations for the diffusion- scaled queueing pro cesses ov er the entire state space. The reason lies in the join tly work conserving (JW C) condition (that is, all the queues hav e to b e empty when there are idle servers) which is essen tial in establishing the w eak con vergence to the controlled limiting diffusion (see [ 10 , 11 ]). T o tackle this difficulty , we first provide an explicit ‘drift’ represen tation of the diffusion-scaled pro cesses which differs from the drift of the diffusion by an extra term that accoun ts for the deviation from the JWC condition in the n th system, and which v anishes in the limit. A natural extension of the concept of work conserv ation for m ulticlass m ulti-p o ol netw orks is minimization of the idle serv ers at all times. This defines an action space whic h w e call system-wide work c onserving (SW C). Establishing the “uniform” geometric ergo dicity ov er all SWC Marko v p olicies when the spare capacit y is positive, is accomplished by first pro ving a useful upper b ound for the minimum of idle servers and cumulativ e queue size for the n t h system, and then using this to derive the F oster– Ly apunov drift inequalities in the region of the state space where the drifts of the diffusion limit and the n t h system do not match. This facilitates establishing the drift inequalities for the diffusion- scaled processes. As a consequence of the F oster-Lyapuno v equations, the inv arian t probability measures of the diffusion-scaled queueing pro cesses hav e uniform exp onential tails. The prop ert y of interc hange of limits attests to the v alidity of the diffusion appro ximation for the queueing netw ork. F or sto chastic netw orks in the conv entional heavy traffic regime, we refer the readers to the pap ers [ 13 , 15 , 18 , 21 , 34 , 35 ] and references therein. F or the ‘V’ netw ork in the Halfin–Whitt regime, interc hange of limits is established in [ 4 , 17 ]. F or the ‘N’ net work, Stolyar [ 27 ] has shown the interc hange of limits under a sp ecific static priority p olicy . This prop ert y also holds for net works with p o ol-dep enden t service rates under the LQFS-LB scheduling p olicy , as sho wn in [ 29 , Section 7.2]. Stolyar and Y udo vina [ 28 ] and Stolyar [ 27 ] then pro ved tightness of the stationary distributions and in terchange of limits of a leaf-activit y priority p olicy in the sub- diffusion and diffusion scales, resp ectiv ely , in the underloaded regime. This pap er con tributes to this 4 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG literature b y establishing that the limit of the diffusion-scaled inv ariant distributions is equal to the in v arian t distribution of the limiting diffusion pro cess for the large classes of net works considered under an y stationary Marko v p olicy (see Remark 5.3 ). 1.1. Organization of the pap er. In the next subsection, we summarize the notation used in the pap er. In Subsection 2.1 , we describ e the mo del and state informally the assumptions used. W e define the diffusion scaled pro cesses, and c haracterize the corresp onding controlled generator in Subsection 2.2 . In Subsection 2.3 , the notion of system-wide work c onserving p olicies is introduced, and this is used in Subsection 2.4 to take limits and establish the diffusion appro ximation. In Section 3 , w e define the parameter of sp ar e c ap acity ( % ) for m ulticlass multi-po ol netw orks and sho w that whenever % < 0, the pro cess is transien t under any stationary Marko v control b oth for the diffusion limit and the n th system for the mo dels under consideration. In the same subsection, w e establish the relation b etw een the spare capacity and a verage idleness. In Section 4 we first pro vide equiv alent c haracterizations of uniform exp onen tial ergo dicit y of controlled diffusions, and then pro ceed to establish that the diffusion limits of the aforementioned classes of netw orks are uniformly exp onentially ergo dic and their inv ariant probability measures hav e uniform exp onen tial tails. Finally , Section 5 is devoted to the study of uniform exp onential ergo dicit y of the n th system of net works under consideration. 1.2. Notation. W e use R m (and R m + ), m ≥ 1, to denote real-v alued m -dimensional (nonnegative) v ectors, and write R for the real line. W e use z T to denote the transp ose of a v ector z ∈ R m . Throughout the pap er e ∈ R m stands for the vector whose elemen ts are equal to 1, that is, e = (1 , . . . , 1) T , and e i ∈ R m denotes the vector whose elemen ts are all 0 except for the i th elemen t whic h is equal to 1. F or x, y ∈ R , x ∨ y = max { x, y } , x ∧ y = min { x, y } , x + = max { x, 0 } and x − = max {− x, 0 } . F or a set A ⊆ R m , w e use A c , ∂ A , and 1 A to denote the complemen t, the b oundary , and the indicator function of A , resp ectively . A ball of radius r > 0 in R m around a p oint x is denoted by B r ( x ), or simply as B r if x = 0. W e also let B ≡ B 1 . The Euclidean norm on R m is denoted b y | · | , and h· , ·i stands for the inner pro duct. F or x ∈ R m , w e let k x k 1 : = P i | x i | , and by K r , or K ( r ), for r > 0, we denote the closed cub e K r : = { x ∈ R m : k x k 1 ≤ r } . (1.1) Also, w e define x max : = max i x i , and x min : = min i x i , and x ± : =  x ± 1 , . . . , x ± m  . F or a finite signed measure ν on R m , and a Borel measurable f : R m → [1 , ∞ ), the f -norm of ν is defined b y k ν k f : = sup g ∈B ( R m ) , | g |≤ f     Z R m g ( x ) ν (d x )     , (1.2) where B ( R m ) denotes the class of Borel measurable functions on R m . 2. The queueing network model and the diffusion limit In this section, we consider a sequence of parallel serv er net works whose pro cesses, parameters, and v ariables are indexed by n . W e recall some of the definitions and notations used in [ 5 , 7 ]. 2.1. Mo del and assumptions. Consider a general Marko vian parallel server (m ulticlass multi- p o ol) net work with m classes of customers and J serv er p o ols. Customer classes tak e v alues in I = { 1 , . . . , m } and server p o ols in J = { 1 , . . . , J } . F orming their own queue, customers of eac h class are serv ed according to a First-Come-First-Served (FCFS) service discipline. W e assume throughout the pap er that customers do not abandon. F or all i ∈ I , let J ( i ) denote the subset of server p o ols that can serve customer class i . On the other hand, for all j ∈ J , let I ( j ) b e the subset of customer classes that can b e served by server p o ol j . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 5 W e form a bipartite undirected graph G = ( I ∪ J , E ) with a set of edges defined by E = { ( i, j ) ∈ I × J : j ∈ J ( i ) } , and use the notation i ∼ j , if ( i, j ) ∈ E , and i  j , otherwise. W e assume that the graph G is a tree. W e define R G + : =  ξ = [ ξ ij ] ∈ R m × J + : ξ ij = 0 for i  j  , (2.1) and analogously define R G , Z G + , and Z G . In each server po ol j , w e let N n j b e the num ber of serv ers, and assume that the servers are statistically iden tical. F or eac h i ∈ I , class i customer arrives according to a Poisson pro cess with arriv al rate λ n i > 0. These customers are served at an exp onen tial rate µ n ij > 0 at server p o ol j if j ∈ J ( i ), and µ n ij = 0 otherwise. Finally , w e assume that the arriv al and service pro cesses of all classes are mutually indep endent. W e study these netw orks in the Halfin–Whitt regime, which in volv es the following assumption on the parameters. There exist p ositiv e constants λ i and ν j , nonnegativ e constants µ ij , with µ ij > 0 for i ∼ j and µ ij = 0 for i  j , and constants ˆ λ i , ˆ µ ij and ˆ ν j , suc h that the following limits exist as n → ∞ : λ n i − nλ i √ n → ˆ λ i , √ n ( µ n ij − µ ij ) → ˆ µ ij , and N n j − nν j √ n → ˆ ν j . (2.2) The parameters λ i and µ ij are the limiting arriv al and service rates, ν j is the limiting service capacit y in p o ol j in the fluid scale, while the parameters ˆ λ i , ˆ µ ij and ˆ ν j are the asso ciated limits in the diffusion scale. An additional standard assumption referred to as the c omplete r esour c e p o oling condition [ 11 , 33 ] concerns the fluid scale equilibrium, and is stated as follows. The linear program (LP) given by Minimize max j ∈J X i ∈I ( j ) ξ ij , sub ject to X j ∈J ( i ) µ ij ν j ξ ij = λ i ∀ i ∈ I , (2.3) has a unique solution ξ ∗ = [ ξ ∗ ij ] ∈ R G + satisfying X i ∈I ξ ∗ ij = 1 , ∀ j ∈ J , and ξ ∗ ij > 0 for all i ∼ j . (2.4) W e define x ∗ ∈ R m , and z ∗ ∈ R G + b y x ∗ i = X j ∈J ξ ∗ ij ν j , and z ∗ ij = ξ ∗ ij ν j . (2.5) The quantit y ξ ∗ ij represen ts the fraction of serv ers in p o ol j allo cated to class i in the fluid equilib- rium, x ∗ i represen ts the total num b er of class i customers in the system, and z ∗ ij is the n umber of class i customers in p o ol j . Note that the constraint in ( 2.3 ) is the rate balance equation for eac h class i with allo cations in each service p o ol j . Also, ρ j := P i ξ ij can b e interpreted as the traffic in tensity in p o ol j , hence the condition in ( 2.4 ) implies that eac h p o ol is critically loaded. F or each i ∈ I and j ∈ J , we let X n i = { X n i ( t ) : t ≥ 0 } denote the total n umber of class i customers in the system (b oth in service and in queue), Z n ij = { Z n ij ( t ) , t ≥ 0 } the n umber of class i customers currently b eing served in p o ol j , Q n i = { Q n i ( t ) , t ≥ 0 } the num b er of class i customers in the queue, and Y n j = { Y n j ( t ) , t ≥ 0 } the num b er of idle servers in serv er p o ol j . Let X n = ( X n i ) i ∈I , Y n = ( Y n j ) j ∈J , Q n = ( Q n i ) i ∈I , and Z n = ( Z n ij ) i ∈I , j ∈J . The pro cess Z n is the scheduling control. W e hav e clearly the following b alanc e e quations Q n i ( t ) : = X n i ( t ) − X j ∈J Z n ij ( t ) , i ∈ I , Y n j ( t ) : = N n j − X i ∈J Z n ij ( t ) , j ∈ J , (2.6) 6 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Dropping the explicit dep endence on n for simplicity , let ( x, z ) ∈ Z m + × Z G + denote a state-action pair. W e rewrite (2.6) as q i ( x, z ) : = x i − X j ∈J z ij , i ∈ I , y j ( z ) : = N n j − X i ∈J z ij , j ∈ J . (2.7) and define the action sp ac e Z n ( x ) b y Z n ( x ) : =  z ∈ Z G + : q i ( x, z ) ∧ y n j ( z ) = 0 , q i ( x, z ) ≥ 0 , y n j ( z ) ≥ 0 ∀ ( i, j ) ∈ E  . Note that this space consists of work-conserving actions only . It should b e noted here that there is an abuse of notation in (2.7) . The quantities q i ( x, z ) and y j ( z ) still represen t the num b er of class i customers in the queue and the num b er of idle serv ers in p o ol j resp ectively . Equation (2.7) is used to sho w the dep endence on x and z through the balance equations. 2.2. Diffusion scaling. With ξ ∗ ∈ R G + the solution of the (LP), we define the cen tering quantities of the diffusion-scaled pro cesses ¯ z n ∈ R G + and ¯ x n ∈ R m b y ¯ z n ij : = 1 n ξ ∗ ij N n j , ¯ x n i : = X j ∈J ¯ z n ij , (2.8) and ˆ X n i ( t ) : = 1 √ n  X n i ( t ) − n ¯ x n i  , ˆ Q n i ( t ) : = 1 √ n Q n i ( t ) , ˆ Z n ij ( t ) : = 1 √ n  Z n ij ( t ) − n ¯ z n ij  , ˆ Y n j ( t ) : = 1 √ n Y n j ( t ) . (2.9) Using (2.6) and (2.9) , these ob ey the c enter e d b alanc e e quations ˆ X n i ( t ) = ˆ Q n i ( t ) + X j ∈J ( i ) ˆ Z n ij ( t ) ∀ i ∈ I , ˆ Y n j ( t ) + X i ∈I ( j ) ˆ Z n ij ( t ) = 0 ∀ j ∈ J . (2.10) W e introduce suitable notation in the diffusion scale as follows (see [ 7 , Definition 2.3]). F or x ∈ Z m + and z ∈ Z n ( x ), w e define ˆ x n : = x − n ¯ x n √ n , ˆ z n : = z − n ¯ z n √ n , (2.11) and let S n denote the state space in the diffusion scale, that is, S n : =  ˆ x ∈ R m : √ n ˆ x + n ¯ x n ∈ Z m +  . (2.12) It is clear that the diffusion-scaled work-conserving action space ˆ Z n ( ˆ x ) takes the form ˆ Z n ( ˆ x ) : =  ˆ z : √ n ˆ z + n ¯ z n ∈ Z n ( √ n ˆ x + n ¯ x n )  , ˆ x ∈ S n . Recall that a sc heduling p olicy is called stationary Marko v if Z n ( t ) = z ( X n ( t )) for some function z : Z m + → Z G + , in whic h case w e iden tify the p olicy with the function z . Under a stationary Mark ov p olicy , X n is Mark ov with controlled generator L n z f ( x ) : = X i ∈I λ n i  f ( x + e i ) − f ( x )  + X j ∈J ( i ) µ n ij z ij  f ( x − e i ) − f ( x )  ! (2.13) UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 7 for f ∈ C ( R m ) and x ∈ Z m + . Let ` n = ( ` n 1 , . . . , ` n m ) T b e defined by ` n i : = 1 √ n  λ n i − X j ∈J ( i ) µ n ij ξ ∗ ij N n j  . (2.14) By (2.8) , the assumptions on the parameters in (2.2) and (2.3) , we hav e ` n i − − − → n →∞ ` i : = ˆ λ i − X j ∈J ( i ) ˆ µ ij z ∗ ij − X j ∈J ( i ) µ ij ξ ∗ ij ˆ ν j , with z ∗ as in (2.5) . Let ` : = ( ` 1 , . . . , ` m ) T . Note that ` n i and ` i can b e regarded as the deficit or surplus in the n umber of serv ers (of order O ( √ n ) allocated to class i in the diffusion scale. Note also that ` n i and ` i app ears as constants in the drift of the diffusion-scaled and diffusion limit pro cesses; see (2.25) and (2.30) W e drop the dep endence on n in the diffusion-scaled v ariables in order to simplify the notation. A work-conserving stationary Marko v p olicy z , that is a map z : Z m + → Z G + suc h that z ( x ) ∈ Z n ( x ) for all x ∈ Z m + , gives rise to a p olicy ˆ z : S n → R G , with ˆ z ( ˆ x ) ∈ ˆ Z n ( ˆ x ) for all ˆ x ∈ S n , via (2.11) (and vice-v ersa). Using (2.9) , (2.13) , and (2.14) and rearranging terms, the controlled generator of the corresp onding diffusion-scaled pro cess can b e written as b L n ˆ z f ( ˆ x ) = X i ∈I λ n i n d f  ˆ x ; 1 √ n e i  + d f  ˆ x ; − 1 √ n e i  n − 1 − X i ∈I b n i ( ˆ x, ˆ z ) d f  ˆ x ; − 1 √ n e i  n − 1 / 2 , ˆ x ∈ S n , ˆ z ∈ ˆ Z n ( ˆ x ) , (2.15) where d f is giv en by d f ( x ; y ) : = f ( x + y ) − f ( x ) , x, y ∈ R m , and the ‘drift’ b n = ( b n 1 , · · · , b n m ) T is giv en by b n i ( ˆ x, ˆ z ) : = ` n i − X j ∈J ( i ) µ n ij ˆ z ij , ˆ z ∈ ˆ Z n ( ˆ x ) , i ∈ I . (2.16) Abusing the notation for ˆ x ∈ S n and ˆ z ∈ ˆ Z n ( ˆ x ), we define (compare with (2.10) ) ˆ q n i ( ˆ x, ˆ z ) : = ˆ x i − X j ∈J ( i ) ˆ z ij , i ∈ I , ˆ y n j ( ˆ z ) : = − X i ∈I ( j ) ˆ z ij , j ∈ J , (2.17) and ˆ ϑ n ( ˆ x, ˆ z ) : = h e, ˆ q n ( ˆ x, ˆ z )  ∧ h e, ˆ y n ( ˆ z )  . (2.18) Recall (2.9) . The parameter ˆ ϑ n can therefore b e regarded as the scaled minim um of the total n umber of customers in the queues and the total n umber of idle servers. By (2.17) , we hav e  e, ˆ q n ( ˆ x, ˆ z )  = ˆ ϑ n ( ˆ x, ˆ z ) + h e, ˆ x i + , and  e, ˆ y n ( ˆ z )  = ˆ ϑ n ( ˆ x, ˆ z ) + h e, ˆ x i − (2.19) for all ˆ x ∈ S n and ˆ z ∈ ˆ Z n ( ˆ x ). Define the ( m − 1) and ( J − 1) simplexes ∆ c : = { u ∈ R m : u ≥ 0 , h e, u i = 1 } , and ∆ s : = { u ∈ R J : u ≥ 0 , h e, u i = 1 } , (2.20) and let ∆ : = ∆ c × ∆ s . By (2.19) , there exists u = ( u c , u s ) ∈ ∆ such that ˆ q n ( ˆ x, ˆ z ) =  ˆ ϑ n ( ˆ x, ˆ z ) + h e, ˆ x i +  u c , and ˆ y n ( ˆ z ) =  ˆ ϑ n ( ˆ x, ˆ z ) + h e, ˆ x i −  u s . (2.21) Let D : = n ( α, β ) ∈ R m × R J : P m i =1 α i = P J j =1 β j o . 8 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG As sho wn in [ 10 , Prop osition A.2], there exists a unique linear map Φ = [Φ ij ] : D → R G solving X j ∈J ( i ) Φ ij ( α, β ) = α i ∀ i ∈ I , and X i ∈I ( j ) Φ ij ( α, β ) = β j ∀ j ∈ J . (2.22) The solutions Φ ij corresp ond to the resource allo cations ˆ z ij in the diffusion scale, see ( 2.23 ). Since  ˆ x − ˆ q n ( ˆ x, ˆ z ) , − ˆ y n ( ˆ z )  ∈ D by (2.7) and (2.17) , using the linearity of the map Φ and (2.21) and (2.22) , it follo ws that ˆ z = Φ  ˆ x − ˆ q n ( ˆ x, ˆ z ) , − ˆ y n ( ˆ z )  = Φ  ˆ x − h e, ˆ x i + u c , −h e, ˆ x i − u s  − ˆ ϑ n ( ˆ x, ˆ z ) Φ( u c , u s ) . (2.23) W e describ e an imp ortan t prop ert y of the linear map Φ which we need later. Consider the matrices B n 1 ∈ R m × m and B n 2 ∈ R m × J defined b y X j ∈J ( i ) µ n ij Φ ij ( α, β ) =  B n 1 α + B n 2 β  i , ∀ i ∈ I , ∀ ( α, β ) ∈ D . (2.24) It is clear that for B n 1 to b e a nonsingular matrix the basis used in the representation of the linear map Φ should b e of the form D =  α, ( β ) − j  , j ∈ J , where ( β ) − j = { β ` , ` 6 = j } . Since Φ has a unique representation in terms of such a basis, and since B n i , i = 1 , 2, are determined uniquely from Φ b y (2.24) , abusing the terminology , we refer to suc h an D as a basis for B n i , i = 1 , 2. In [ 5 , Lemma 4.3], the follo wing prop erty is asserted: Given any ˆ ı ∈ I , there exists an ordering of { α i , i ∈ I } with α ˆ ı the last element, and ˆ  ∈ J , such that the matrix B n 1 is lo wer diagonal with p ositiv e diagonal elements with resp ect to this ordered basis  α, ( β ) − ˆ   . F or more details, we refer the reader to [ 5 , Section 4.1]. In view of (2.23) and (2.24) , for any ˆ z ∈ e Z n ( ˆ x ) with ˆ x ∈ S n , there exists u = u ( ˆ x, ˆ z ) ∈ ∆ such that the drift b n in (2.16) tak es the form b n ( ˆ x, ˆ z ) = ` n − B n 1  ˆ x − h e, ˆ x i + u c  + B n 2 u s h e, ˆ x i − + ˆ ϑ n ( ˆ x, ˆ z )  B n 1 u c + B n 2 u s  . (2.25) 2.3. Join t and system-wide work conserv ation. W e start with the following definition. Definition 2.1. W e sa y that an action ˆ z ∈ ˆ Z n ( ˆ x ) is jointly work c onserving (JWC), if ˆ ϑ n ( ˆ x, ˆ z ) = 0. Recall that a work conserving p olicy refers to an action in whic h a server is idle if and only if there is no customer waiting in the queue that this server can serve. A jointly w ork conserving action k eeps all servers busy unless all queues are empt y . Let ˆ ϑ n ∗ ( ˆ x ) : = min ˆ z ∈ ˆ Z n ( ˆ x ) ˆ ϑ n ( ˆ x, ˆ z ) , ˆ x ∈ S n , and e Z n ( ˆ x ) : =  ˆ z ∈ ˆ Z n ( ˆ x ) : ˆ ϑ n ( ˆ x, ˆ z ) = ˆ ϑ n ∗ ( ˆ x )  , ˆ x ∈ S n . W e refer to e Z n ( ˆ x ) as the system-wide work c onserving (SWC) action set at ˆ x . A stationary Mark ov sc heduling p olicy ˆ z is called SWC if ˆ z ( ˆ x ) ∈ e Z n ( ˆ x ) for all ˆ x ∈ S n . W e let e Z n denote the class of all such p olicies. Since z and ˆ z are related by (2.11) , abusing this terminology , we also refer to a Mark ov p olicy z : Z m + → Z G + as SW C, if it satisfies z ( x ) − n ¯ z n √ n ∈ e Z n  x − n ¯ x n √ n  , and w e write z ∈ e Z n . W e recall [ 11 , Lemma 3] which states that there exists M 0 > 0 such that the collection of sets ˘ X n defined b y ˘ X n : =  ˆ x ∈ S n : k ˆ x k 1 ≤ M 0 √ n  , (2.26) UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 9 has the following prop ert y . If ˆ x ∈ ˘ X n , then for any pair ( ˆ q , ˆ y ) suc h that √ n ˆ q ∈ Z m + , √ n ˆ y ∈ Z J + , and satisfying h e, ˆ q i ∧ h e, ˆ y i = 0 , h e, ˆ x − ˆ q i = h e, − ˆ y i , and ˆ y j ≤ N n j , j ∈ J , it holds that Φ( ˆ x − ˆ q , − ˆ y ) ∈ ˆ Z n ( ˆ x ). It follows from this lemma and Definition 2.1 that if ˆ x ∈ ˘ X n , then the actions in e Z n ( ˆ x ) are JWC. R emark 2.1 . Using (2.19) , we know that under the JW C condition  e, ˆ q n ( ˆ x, ˆ z )  = h e, ˆ x i + , and  e, ˆ y n ( ˆ z )  = h e, ˆ x i − . Rewriting (2.21) under the JWC condition, w e therefore hav e the following u c i =    ˆ q n i ( ˆ x, ˆ z ) h e, ˆ q n ( ˆ x, ˆ z )  , if h e, ˆ q n ( ˆ x, ˆ z )  > 0 , e 1 , otherwise and u s j =    ˆ y n j ( ˆ z ) h e, ˆ y n ( ˆ z )  , if h e, ˆ y n ( ˆ z )  > 0 , e 1 , otherwise . Hence, one can see that the con trol u c i represen ts the prop ortion of the total queue length in the net work at queue i , while u s j represen ts the prop ortion of the total num b er of idle servers in the net work at p o ol j . Recall from Definition 2.1 that a JWC action keeps all serv ers busy unless all queues are empt y . It is clear that this cannot be alw ays enforced in multiclass multi-po ol net works ov er the entire state space. A SW C control enforces the complementarit y b etw een o verall queue length and service. 2.4. The diffusion limit. The diffusion appr oximation or diffusion limit of the queueing mo del describ ed ab ov e is an m -dimensional sto chastic differen tial equation (SDE) of the form (see [ 10 ] and [ 11 , Section 2.5]) d X t = b ( X t , U t ) d t + σ ( X t ) d W t , X 0 = x ∈ R m . (2.27) Here, { W t } t ≥ 0 is a standard m -dimensional Brownian motion, and the control U t tak es v alues in the set ∆ = ∆ c × ∆ s defined in (2.20) . The drift b can b e derived as follo ws. Recall R G in ( 2.1 ). F or u = ( u c , u s ) ∈ ∆ , let b Φ[ u ] : R m → R G b e defined by b Φ[ u ]( x ) : = Φ  x − ( e · x ) + u c , − ( e · x ) − u s  , (2.28) with Φ as defined in (2.22) . Then the drift b takes the form b i ( x, u ) = ` i − X j ∈J ( i ) µ ij b Φ ij [ u ]( x ) . (2.29) By [ 5 , Lemma 4.3], we also know that (2.29) can b e expressed as b ( x, u ) = ` − B 1  x − h e, x i + u c  + B 2 u s h e, x i − , (2.30) where B 1 ∈ R m × m is a low er diagonal matrix with p ositive diagonal elements, and B 2 ∈ R m × J . Of course B i in (2.30) and B n i in (2.25) , i = 1 , 2, ha ve the same functional form with resp ect to { µ ij } and { µ n ij } , resp ectiv ely . The diffusion matrix σ ∈ R m × m is constan t, and a : = σσ T = diag(2 λ 1 , . . . , 2 λ m ) . In addition, for f ∈ C 2 ( R m ), w e define L u f ( x ) : = 1 2 trace  a ∇ 2 f ( x )  +  b ( x, u ) , ∇ f ( x )  , (2.31) 10 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG with ∇ 2 f denoting the Hessian of f . In the follo wing, we describ e the drift of the diffusion limit of the net works under consideration while sho wing some examples in Figure 1 . 2.4.1. Networks with a dominant server p o ol. This net work has one non-leaf server no de, which, without loss of generality , we lab el as j = 1. As in Subsection 2.1 , the customer no des are denoted b y I = { 1 , 2 , . . . , m } , and the server no des by J = { 1 , 2 , . . . , J } . Recall that J ( i ) is the collection of sever nodes connected to customer i . Owing to the tree structure of the net work, serv er 1 ∈ J ( i ) for all i ∈ { 1 , 2 , . . . , m } . Let J 1 ( i ) : = J ( i ) \ { 1 } for all i ∈ I . Recall the form of the drift in (2.29) . Using (2.22) , it is simple to sho w that the matrix b Φ ij [ u ] for this netw ork is giv en by b Φ ij [ u ]( x ) =        x i − h e, x i + u c i + P j ∈J 1 ( i ) h e, x i − u s j for j = 1 , −h e, x i − u s j for j ∈ J 1 ( i ) , 0 otherwise. (2.32) Using (2.32) , the drift takes the follo wing simple form: b i ( x, u ) = ` i − µ i 1  x i − u c i h e, x i +  + X j ∈J 1 ( i ) µ i 1  η ij − 1  u s j h e, x i − , i ∈ I , (2.33) with η ij : = µ ij µ i 1 for j ∈ J 1 ( i ) and i ∈ I . Note that B 1 = diag ( µ 11 , . . . , µ m 1 ), and so ` = − % m B 1 e , where % is given by (2.35) . W e define ¯ η : = max i ∈I max j ∈J 1 ( i ) η ij , and η : = min i ∈I min j ∈J 1 ( i ) η ij . (a) Generalized ‘N’ Netw ork (b) Generalized ‘M’ Netw ork sq u a r e - cu st o m e r cl a ss, ci r cl e - se r ve r p o o l , t h e so l i d ci r cl e i s t h e d o m i n a t i n g se r ve r p o o l (c) Netw ork with a dominant server p ool Figure 1. Examples of multiclass multi-po ol netw orks with a dominant p ool (square–customer classes, circle–server p o ols, solid circle–the dominant server p o ol) 2.4.2. Networks with class-dep endent servic e r ates. W e consider in this part arbitrary tree netw orks where the service rates are dictated by the customer t yp e; namely µ ij = µ i for all ( i, j ) ∈ E . Recall the definition in (2.28) . Using (2.22) and (2.29) , the drift of this net work takes the form b i ( x, u ) = ` i − X j ∈J ( i ) µ ij b Φ ij [ u ]( x ) = ` i − µ i  x i − u c i h e, x i +  , ∀ i ∈ I . (2.34) Note then that B 1 = diag( µ 1 , . . . , µ m ) and B 2 = 0. W e remark that b oth classes of netw orks hav e one common feature: the matrix B 1 is diagonal. UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 11 2.5. A necessary and sufficient condition for uniform stability. W e define the spare capacit y (or the safet y staffing) for the n th system (prelimit) and the diffusion limit by % n : = −  e ( B n 1 ) − 1 , ` n  , and % : = −  eB − 1 1 , `  , (2.35) resp ectiv ely . Note, of course, that % n → % as n → ∞ by (2.2) . Recall the expressions of ` n and ` in (2.14) . Recall that ` n i and ` i can b e regarded as the deficit or surplus in the num b er of servers (of order O ( √ n ) allocated to class i in the diffusion scale. Hence % n and % can b e regarded as the optimal reallo cation of the capacity fluctuations (p ositiv e or negative) of order √ n when each serv er p o ol emplo ys a square-ro ot staffing rule. W e summarize the main results of the pap er. Theorem 2.1. Consider a network with a dominant server p o ol, or with class-dep endent servic e r ates. Then the c onditions % > 0 and % n > 0 ar e ne c essary and sufficient for the uniform stability of the limiting diffusion and the diffusion-sc ale d queueing pr o c esses, r esp e ctively. Mor e pr e cisely: (i) if % < 0 , the pr o c ess { X t } t ≥ 0 in (2.27) is tr ansient under any stationary Markov c ontr ol. In addition, if % = 0 , then { X t } t ≥ 0 c annot b e p ositive r e curr ent. (ii) if % n < 0 , the pr o c ess { X n t } t ≥ 0 is tr ansient under any stationary Markov sche duling p olicy. In addition, if % n = 0 , then { X n t } t ≥ 0 c annot b e p ositive r e curr ent. (iii) if % > 0 , the pr o c esses { X t } t ≥ 0 ar e uniformly exp onential ly er go dic over stationary Markov c ontr ols. (iv) if % n > 0 , the pr o c esses { X n t } t ≥ 0 ar e uniformly exp onential ly er go dic over SWC sche duling p olicies, and the invariant distributions have exp onential tails. P arts (i) and (ii) of Theorem 2.1 follow from Prop ositions 3.1 and 3.2 , resp ectiv ely . P art (iii) follo ws from Theorem 4.2 , and part (iv) from Theorem 5.1 . 3. Two proper ties of the sp are cap acity and transience In the first part of this section we pro ve the results in Theorem 2.1 (i) and (ii). It is imp ortant to note that for the mo dels in Subsections 2.4.1 and 2.4.2 w e hav e 1 + h e, B − 1 1 B 2 u s i > 0 . (3.1) Note that for netw orks with a dominan t server p o ol, w e hav e 1 = u s 1 + P i ∈I P j ∈J 1 ( i ) u s j . Hence 1 + h e, B − 1 1 B 2 u s i = 1 + X i ∈I X j ∈J 1 ( i ) ( η ij − 1) u s j = 1 − X i ∈I X j ∈J 1 ( i ) u s j + X i ∈I X j ∈J 1 ( i ) η ij u s j = u s 1 + X i ∈I X j ∈J 1 ( i ) η ij u s j > 0 . Note also that for netw orks with a class dep enden t service rate we hav e established that B 2 = 0. Hence 1 + h e, B − 1 1 B 2 u s i = 1 > 0. Prop osition 3.1. Supp ose that % = −h e, B − 1 1 ` i < 0 . Then the pr o c ess { X t } t ≥ 0 in (2.27) is tr ansient under any stationary Markov c ontr ol. In addition, if % = 0 , then { X t } t ≥ 0 c annot b e p ositive r e curr ent. Pr o of. Recall that the first and second order deriv ativ es of the h yp erb olic tangent function are tanh 0 ( x ) = 1 cosh 2 ( x ) ; tanh 00 ( x ) = − 2 tanh( x ) cosh 2 ( x ) 12 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Let H ( x ) : = tanh  β h e, B − 1 1 x i  , with β > 0. Then trace  a ∇ 2 H ( x ))  = β 2 tanh 00  β h e, B − 1 1 x i    σ T B − 1 1 e   2 . W e hav e L u H ( x ) = 1 2 trace  a ∇ 2 H ( x )  +  b ( x, u ) , ∇ H ( x )  = − β 2 tanh  β h e, B − 1 1 x i  cosh 2  β h e, B − 1 1 x i  | σ T B − 1 1 e | 2 + β cosh 2  β h e, B − 1 1 x i    e, B − 1 1 `  + h e, x i −  1 +  e, B − 1 1 B 2 u s    . (3.2) Th us, for 0 < β < h e, B − 1 1 ` i | σ T B − 1 1 e | − 2 , w e obtain L u H ( x ) > 0 by (3.1) . Therefore, { H  X t  } t ≥ 0 is a b ounded submartingale, so it conv erges almost surely . Since X is irreducible, it can b e either recurren t or transient. If it is recurren t, then H should b e constant a.e. in R m , which is not the case. Thus X is transient. W e now turn to the case where % = 0. Supp ose that the pro cess { X ( t ) } t ≥ 0 (under some stationary Mark ov control) has an in v arian t probabilit y measure π (d x ). It is w ell kno wn that π m ust ha ve a p ositiv e density . Let h 1 ( x ) and h 2 ( x ) denote resp ectively the first and the second terms on the righ t hand side of (3.2) . Applying Itˆ o’s form ula to (3.2) with X 0 = x as in (2.27) , we obtain E π  H ( X t ∧ τ r )  − H ( x ) = X i =1 , 2 E π " Z t ∧ τ r 0 h i ( X s )d s # , (3.3) where τ r denotes the first exit time from B r , r > 0. Note that h 1 ( x ) is bounded and h 2 ( x ) is non-negativ e. Thus using dominated and monotone con vergence, w e can tak e limits in (3.3) as r → ∞ for the terms on the right side to obtain Z R m H ( x ) π (d x ) − H ( x ) = t X i =1 , 2 Z R m h i ( x ) π (d x ) , t ≥ 0 . Since H ( x ) is b ounded, w e can divide b oth sides b y t and β and take the limit as t → ∞ to get Z R m β − 1 h 1 ( x ) π (d x ) + Z R m β − 1 h 2 ( x ) π (d x ) = 0 . (3.4) Since β − 1 h 1 ( x ) tends to 0 uniformly in x as β & 0, the first term on the left hand side of (3.4) v anishes as β & 0. Ho wev er, since β − 1 h 2 ( x ) is b ounded aw a y from 0 on the op en set { x ∈ R m : h e, x i − > 1 } , this contradicts the fact that π (d x ) has full supp ort.  R emark 3.1 . In the pro of of Prop osition 3.1 , the function H ( x ) is a b ounded test function which w as c hosen so that L u H ( x ) ≥ 0. Assume that X t is recurrent. Since H ( X t ) conv erges a.s. (b eing a b ounded submartingale), and since X t ‘visits every op en neighborho o d’ in R m with probability 1, it follo ws that H must b e a constan t function which is a contradiction. The formalism b ehind the ab o ve argument is as follows: Let τ b e the first hitting time to the unit ball B centered at x = 0. If E x denotes the exp ectation op erator on the canonical space of the Marko v pro cess { X t } t ≥ 0 with initial condition X 0 = x , then by Dynkin’s formula we obtain E x [ H ( X τ )] ≥ H ( x ). If the pro cess is recurrent, then of course P x ( τ < ∞ ) = 1, which giv es E x [ H ( X τ )] ≤ sup y ∈ B H ( y ). Th us sup y ∈ B H ( y ) ≥ H ( x ) for all x ∈ R m whic h is not true. Moving the cen ter of the the ball B to an arbitrary p oin t z , and denoting it as B ( z ), w e similarly hav e sup y ∈ B ( z ) H ( y ) ≥ H ( x ) ∀ x, z ∈ R m . This implies that H ( x ) = constant whic h is a contradiction and hence X is transient. UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 13 Prop osition 3.2. Supp ose that % n < 0 . Then the state pr o c ess { X n t } t ≥ 0 of the n th system is tr ansient under any stationary Markov sche duling p olicy. In addition, if % n = 0 , then { X n t } t ≥ 0 c annot b e p ositive r e curr ent. Pr o of. The pro of mimics that of Prop osition 3.1 . W e apply the function H in that pro of to the op erator b L n ˆ z in (2.15) , and use the identit y H  x ± 1 √ n e i  − H ( x ) ∓ 1 √ n ∂ x i H ( x ) = 1 n Z 1 0 (1 − t ) ∂ x i x i H  x ± t √ n e i  d t , (3.5) to express the first and second order incremental quotien ts, together with (2.25) whic h implies that  b n ( ˆ x, ˆ z ) , ∇ H ( ˆ x )  = β cosh 2  β h e, ( B n 1 ) − 1 ˆ x i    e, ( B n 1 ) − 1 ` n  +  ˆ ϑ n ( ˆ x, ˆ z ) + h e, ˆ x i −   1 +  e, ( B n 1 ) − 1 B n 2 u s    . The rest follo ws exactly as in the pro of of Prop osition 3.1 .  3.1. Spare capacit y and av erage idleness. It is shown in [ 4 , 8 ] that if the diffusion limit of the ‘V’ netw ork with no abandonment has a inv ariant distribution π under some stationary Marko v con trol, then % represents the ‘av erage idleness’ of the system, that is, % = R R m h e, x i − π (d x ). In calculating this av erage for multi-po ol netw orks, idle servers are not w eighted equally across differen t p o ols and the term  e, B − 1 1 B 2 u s ( x )  app ears in the expression, see (3.6) . It is imp ortant to note that only the control on the idleness allo cations among server p o ols u s app ears in the identit y , and the con trol comp onent u c do es not. Theorem 3.1. Consider a network with a dominant server p o ol, or with class-dep endent servic e r ates, and supp ose that % > 0 . L et π u denote the invariant invariant pr ob ability me asur e c orr e- sp onding to a stationary Markov c ontr ol u ∈ U sm , whose existenc e fol lows fr om The or em 2.1 ( iii ) . Then % = Z R m  1 +  e, B − 1 1 B 2 u s ( x )   h e, x i − π u (d x ) . (3.6) Pr o of. The pro of is similar to that of [ 8 , Corollary 5.1], but more inv olv ed for the multiclass mul ti- p o ol netw orks. W e first recall some definitions and notations. Let χ r ( t ), ˘ χ r ( t ), r > 1 b e smo oth, conca ve and conv e x functions, resp ectively , defined by χ r ( t ) = ( t , t ≤ r − 1 , r − 1 2 , t ≥ r , and ˘ χ r ( t ) = ( t , t ≥ 1 − r , 1 2 − r , t ≤ − r . Let g r ( x ) = ˘ χ r  h e, B − 1 1 x i  , and f r ( x ) = χ r  g r ( x )  . A straigh tforward calculation shows that h b ( x, u ) , ∇ f r ( x ) i = h 1 ( x ) + h 2 ( x ) , 1 2 trace  a ( x ) ∇ 2 f r ( x )  = h 3 ( x ) + h 4 ( x ) , where h 1 ( x ) : = − %χ 0 r  f ( x )  ˘ χ 0 r  h e, B − 1 1 x i  , h 2 ( x ) : =  1 +  e, B − 1 1 B 2 u s  χ 0 r  f ( x )  ˘ χ 0 r  h e, B − 1 1 x i  h e, x i − , h 3 ( x ) : = 1 2 χ 00 r  f ( x )   ˘ χ 0 r  h e, B − 1 1 x i  2   σ T B − 1 1 e   2 , h 4 ( x ) : = 1 2 χ 0 r  f ( x )  ˘ χ 00 r  h e, B − 1 1 x i    σ T B − 1 1 e   2 . 14 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG W e note that g r ( x ) is p ositiv e and bounded b elo w aw a y from 0, and f r ( x ) is smo oth, b ounded, and has b ounded deriv ativ es. Also note that h i , i = 1 , 2 , 3, are b ounded, and h 2 is nonnega- tiv e. Therefore, if { X ( t ) } t ≥ 0 is p ositiv e recurrent with an inv arian t probability measure π u (d x ), a straigh tforward application of Itˆ o’s formula shows that π u  L u f r  = 0. Therefore, w e obtain π u ( − h 1 ) = π u ( h 2 ) + π u ( h 3 ) + π u ( h 4 ) . (3.7) By the definition of χ r and ˘ χ r , it is straightforw ard to v erify that lim r →∞ π u ( h 3 ) = lim r →∞ π u ( h 4 ) = 0 . (3.8) In addition, using dominated conv ergence theorem lim r →∞ π u ( h 1 ) = − % , lim r →∞ π u ( h 2 ) = Z R m  1 +  e, B − 1 1 B 2 u s  h e, x i − π u (d x ) . (3.9) Com bining (3.7) – (3.9) , w e obtain (3.6) .  R emark 3.2 . Using Theorem 3.1 and the drift in Subsections 2.4.1 and 2.4.2 , it is easy to v erify that in the case of a netw ork with a dominan t server p o ol w e hav e % = Z R m " 1 + X i ∈I X j ∈J 1 ( i )  µ ij µ i 1 − 1  u s j # h e, x i − π u (d x ) , whereas in the case of a netw ork with class-dep enden t service rates, (3.6) tak es the form % = Z R m h e, x i − π u (d x ) . 4. Uniform exponential ergodicity of the diffusion limit In this section w e sho w that if % > 0 then the diffusion limit of a netw ork with a dominant serv er p o ol, or with class-dep endent service rates, is uniformly exp onen tially ergo dic and the in v arian t distributions ha ve exp onential tails. W e start b y reviewing the notion of uniform exp onential er go dicity for a controlled diffusion. W e do this under fairly general hypotheses. Consider a controlled diffusion pro cess X = { X t , t ≥ 0 } whic h takes v alues in the m -dimensional Euclidean space R m , and is gov erned by the Itˆ o equation d X t = b  X t , v ( X t )  d t + σ ( X t ) d W t . (4.1) All random pro cesses in (4.1) live in a complete probabilit y space (Ω , F , P ). The pro cess W is a m -dimensional standard Wiener pro cess indep endent of the initial condition X 0 . The function v maps R m to a compact, metrizable set U and is Borel measurable. The collection of such functions comprising of the set of stationary Marko v con trols is denoted by U sm . The parameters of the equation ( 4.1 ) satisfy the following: (1) L o c al Lipschitz c ontinuity: The functions b : R m × U → R m and σ : R m → R m × m are con tinuous, and satisfy | b ( x, u ) − b ( y , u ) | + k σ ( x ) − σ ( y ) k ≤ C R | x − y | ∀ x, y ∈ B R , ∀ u ∈ U . for some constan t C R > 0 dep ending on R > 0. (2) Affine gr owth c ondition: F or some C 0 > 0, we hav e sup u ∈ U h b ( x, u ) , x i + + k σ ( x ) k 2 ≤ C 0  1 + | x | 2  ∀ x ∈ R m . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 15 (3) Nonde gener acy: F or each R > 0, it holds that m X i,j =1 a ij ( x ) ξ i ξ j ≥ C − 1 R | ξ | 2 ∀ x ∈ B R , and for all ξ = ( ξ 1 , . . . , ξ m ) T ∈ R m , where a = 1 2 σσ T . It is w ell kno wn that, under hypotheses (1)–(2), (4.1) has a unique strong solution whic h is also a strong Marko v pro cess for any v ∈ U sm [ 22 ]. W e let E v x denote the exp ectation op erator on the canonical space of the pro cess controlled by v , with initial condition X 0 = x . Let τ ( A ) denote the first exit time from the set A ∈ R m . W e say that the process { X t } t ≥ 0 is uniformly exp onential ly er go dic if for some ball B ◦ there exist δ ◦ > 0 and x ◦ ∈ ¯ B c ◦ suc h that sup v ∈ U sm E v x ◦ [e δ ◦ τ ( B c ◦ ) ] < ∞ . W e let b A denote the op erator b A φ ( x ) : = 1 2 trace  a ( x ) ∇ 2 φ ( x )  + max u ∈ U  b ( x, u ) , ∇ φ ( x )  , x ∈ R m , for φ ∈ C 2 ( R m ). F or a lo cally b ounded, Borel measurable function f : R m → R , which is b ounded from b elo w in R m , i.e., inf R m f > −∞ , we define the generalized principal eigenv alue of b A + f by Λ( f ) : = inf n λ ∈ R : ∃ ϕ ∈ W 2 ,m loc ( R m ) , ϕ > 0 , b A ϕ + ( f − λ ) ϕ ≤ 0 a.e. in R m o , where W 2 ,m loc is a lo cal Sob olev space. W e ha ve the follo wing e quiv alent c haracterizations of uniform exp onen tial ergo dicit y . This is a straightforw ard extension of [ 8 , Theorem 3.1] for con trolled diffu- sions, and is stated without pro of. Recall that a map f : R m → R is called c o er cive , or inf-c omp act , if inf B c r f → ∞ as r → ∞ . Theorem 4.1. The fol lowing ar e e quivalent. ( a ) F or some b al l B ◦ ther e exists δ ◦ > 0 and x ◦ ∈ ¯ B c ◦ such that sup v ∈ U sm E x ◦ [e δ ◦ τ ( B c ◦ ) ] < ∞ . ( b ) F or every b al l B ther e exists δ > 0 such that sup v ∈ U sm E v x [e δ τ ( B c ) ] < ∞ for al l x ∈ B c . ( c ) F or every b al l B , ther e exists a c o er cive function V ∈ W 2 ,p loc ( R m ) , p > d , with inf R m V ≥ 1 , and p ositive c onstants κ 0 and δ such that b A V ( x ) ≤ κ 0 1 B ( x ) − δ V ( x ) ∀ x ∈ R m . (4.2) ( d ) Equation (4.1) is r e curr ent, and Λ( 1 B c ) < 1 for every b al l B . R emark 4.1 . Recall (1.2) . Let P v t ( x, d y ) denote the transition probability of { X t } t ≥ 0 in (4.1) under a con trol v ∈ U sm . It is w ell known that (4.2) implies that there exist constants γ and C γ whic h do not dep end on the control v c hosen, such that   P v t ( x, · ) − π v ( · )   V ≤ C γ V ( x ) e − γ t , ∀ x ∈ R m , ∀ t ≥ 0 , (4.3) where π v denotes the in v arian t probability measure of { X t } t ≥ 0 under the con trol v . In addition, for an y control v ∈ U sm w e hav e Z R m V ( x ) π v (d x ) ≤ κ 0 δ . (4.4) In particular for V ( x ) b eing an exp onen tial function, the moment generating function of { X t } t ≥ 0 is finite. 16 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG 4.1. A class of intrinsic Lyapuno v functions for the queueing net w ork mo del. As seen in Theorem 4.1 , uniform exp onen tial ergo dicit y is equiv alen t to the F oster–Lyapuno v inequality in (4.2) . In establishing this prop erty for the diffusion limit of sto c hastic netw orks, a prop er choice of a Ly apunov function is of tantamoun t imp ortance. W e first describ e an in trinsic class of such functions. W e fix a conv ex function ψ ∈ C 2 ( R ) with the prop erty that ψ ( t ) is constant for t ≤ − 1, and ψ ( t ) = t for t ≥ 0. This is defined by ψ ( t ) : =        − 1 2 , t ≤ − 1 , ( t + 1) 3 − 1 2 ( t + 1) 4 − 1 2 t ∈ [ − 1 , 0] , t t ≥ 0 . F or ε > 0 w e define ψ ε ( t ) : = ψ ( εt ) , Th us ψ ε ( t ) = εt for t > 0. A simple calculation also shows that ψ 00 ε ( t ) ≤ 3 2 ε 2 . Supp ose that B 1 = diag( ˜ µ 1 , . . . , ˜ µ m ). Using the function ψ ε in tro duced ab o ve, we let Ψ( x ) : = X i ∈I ψ ε ( x i ) ˜ µ i , (4.5) with ε : = % 3 m X i ∈I λ i (3 ˜ µ i + 2) ˜ µ 2 i ! − 1 . (4.6) W e also define V 1 ( x ) : = exp  θ Ψ( − x )  , V 2 ( x ) : = exp  Ψ( x )  , and V ( x ) : = V 1 ( x ) + V 2 ( x ) , (4.7) with θ a p ositiv e constant. As a result of fixing the v alue of ε in (4.6) , Ψ dep ends only on the parameter θ . This simplifies the statemen ts of the results in the rest of the pap er. W e review some useful prop erties of the function ψ ε . Note that for ε > 0 we hav e ψ 0 ε ( t ) t =        0 εt ≤ − 1 , εt ( εt + 1) 2 ( − 2 εt + 1) εt ∈ [ − 1 , 0] , εt εt ≥ 0 . The minim um of ψ 0 ε ( t ) t when εt ∈ [ − 1 , 0] is − 3(39+55 √ 33) 4096 ≥ − 1 2 for εt = − 1+ √ 33 16 . Therefore X i ∈I ψ 0 ε ( x i ) x i = X i ∈I ψ 0 ε ( x i ) x i 1 { εx i ≥ 0 } + X i ∈I ψ 0 ε ( x i ) x i 1 { εx i ∈ [ − 1 , 0] } ≥ ε X i ∈I x i 1 { x i ≥ 0 } − 1 2 X i ∈I 1 { εx i ∈ [ − 1 , 0] } ≥ ε k x + k 1 − m 2 , (4.8) and similarly − P i ∈I ψ 0 ε ( − x i ) x i ≥ ε k x − k 1 − m 2 , where I : = { 1 , 2 , . . . , m } . Note also that − X i ∈I ψ 0 ε ( − x i ) x i ≤ ε h e, x i ≤ X i ∈I ψ 0 ε ( x i ) x i . The function V in (4.7) , scaled by the parameter θ which are suppressed in the notation, is our c hoice of a Lyapuno v function when B 1 is a diagonal matrix. It is constructed in an in tricate manner in order to capture b oth the total w orkload (using Ψ( x )) on the p ositive half-space and the idleness UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 17 (using Ψ( − x )) on the negativ e half-space. As one cannot simply use x + or x − , w e must construct mollified smo oth functions for them. In addit ion, one must recognize that the effects of the workload and idleness on the system are not the same (not symmetric), so we also introduce a parameter θ in the definition of V 1 ( · ). The reader will notice the similarities in (4.7) and [ 4 , Definition 2.2]. Ho wev er the function V used in this pap er is the sum of the tw o exp onen tial functions V 1 and V 2 , whereas their pro duct is used in [ 4 , Lemma 2.1]. As will b e seen later, in the c ase of multiclass m ulti-p o ol netw orks, the analysis is considerably more complex. In Subsection 4.3.1 for example, the proofs required V 1 ≥ V 2 2 on a subset of the state space while requiring the opposite inequalit y on its complemen t in order to establish the F oster-Lyapuno v inequalit y in (4.15) . See the discussions follo wing (4.16) and (4.22) . This is mainly the reason b ehind using the sum instead of the pro duct when defining the function V . In the follo wing subsections, we establish the uniform exp onential ergo dicity of the netw orks under consideration. T o help with the exp osition, w e study the ‘N’ net work in detail and then pro ceed to the more general net works with a dominant server p o ol, and netw orks with class- dep enden t service rates. In establishing the desired drift inequalities, w e often partition the space appropriately , and fo cus on the subsets of the partition. The following cones app ear quite often in the analysis. F or δ ∈ [0 , 1], w e define the cones K + δ : =  x ∈ R m : h e, x i ≥ δ k x k 1  , K − δ : =  x ∈ R m : h e, x i ≤ − δ k x k 1  . (4.9) It is clear that K + 0 ( K − 0 ) corresp onds to the nonnegativ e (nonp ositiv e) canonical half-space, and K + 1 ( K − 1 ) is the nonnegative (nonp ositive) closed orthant. The follo wing identities are very useful. h e, x + i = 1 ± δ 2 k x k 1 , h e, x − i = 1 ∓ δ 2 k x k 1 for x ∈ ∂ K ± δ , δ ∈ [0 , 1] . (4.10) In addition, it is straightforw ard to sho w that X i ∈I ψ ε ( x i ) ≤ X i ∈I ψ ε ( − x i ) if x ∈ K − 0 . (4.11) Also, the follo wing inequality is true in general for an y I 0 ⊂ I . X i ∈I 0 ψ 0 ε ( x i ) x i − ε X i ∈I 0 x i = X x i < 0 , i ∈I 0  ψ 0 ε ( x i ) − ε ) x i ≥ 0 . (4.12) R emark 4.2 . There is an imp ortan t scaling of the drift whic h we employ . Note that if we let ζ = % m e + B − 1 1 ` , with % as in (2.35) , then a mere translation of the origin of the form ˜ X t = X t + ζ results in a diffusion of with the same drift as (2.29) , except that the vector ` gets replaced b y ` = − % m B 1 e . Therefore, we may assume without any loss of generality that the drift in (2.30) takes the form b ( x, u ) = − % m B 1 e − B 1  x − h e, x i + u c  + B 2 u s h e, x i − . (4.13) 4.2. The F oster–Ly apunov inequality. Recall the definition of the op erator L u in (2.31) . W e start with the following simple assertion. Lemma 4.1. L et V b e the function in (4.7) with ε as in (4.6) , and θ ≥ θ 0 wher e θ 0 is a c onstant. Supp ose that for any δ ∈ (0 , 1) ther e exist p ositive c onstants c 0 and c 1 such that the drift b in (4.13) satisfies h b ( x, u ) , ∇ V ( x ) i ≤ c 0 − εc 1 k x k 1 V ( x ) ∀ ( x, u ) ∈  K + δ  c × ∆ , (4.14a) h b ( x, u ) , ∇ V 2 ( x ) i ≤ − %ε m V 2 ( x ) ∀ ( x, u ) ∈ K + δ × ∆ . (4.14b) 18 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Then, ther e exists a c onstant C 0 such that L u V  x  ≤ C 0 − %ε 3 m V ( x ) ∀ ( x, u ) ∈ R m × ∆ . (4.15) Pr o of. A straightforw ard calculation, using the fact that ψ 00 ε ( t ) ≤ 3 2 ε 2 , sho ws that 1 2 trace  a ∇ 2 V 2 ( x )  ≤ ε 2 X i ∈I λ i (3 µ i + 2) 2 µ 2 i V 2 ( x ) ∀ x ∈ R m . Therefore, the c hoice of ε in (4.6) implies that 1 2 trace  a ∇ 2 V 2 ( x )  ≤ %ε 4 m V 2 for all x ∈ R 2 , and th us L u V 2 ( x ) ≤ − 3 %ε 4 m V 2 ( x ) ∀ ( x, u ) ∈ K + δ × ∆ (4.16) b y (4.14b) . Since | x + | ≥ 1+ δ 1 − δ | x − | for all x ∈ K + δ , we may select δ sufficiently close to 1 suc h that V 2 ≥ V 2 1 on K + δ ∩ K c r for some r > 0. Since V 2 has exp onential growth in k x k 1 on K + δ and L u V 1 ( x ) ≤ C (1 + | x | 1 ) V 1 ( x ) on K + δ × ∆ , it then follo ws that (4.15) holds on K + δ × ∆ by (4.16) . It is also clear that (4.15) also holds on  K + δ  c × ∆ by (4.14a) . This completes the pro of.  W e now hav e the following result. Theorem 4.2. Assume that % > 0 . L et V b e the function in (4.7) with ε as in (4.6) , and θ ≥ θ 0 wher e θ 0 is a c onstant. Then the diffusion limit of any network with a d ominant server p o ol or with class-dep endent servic e r ates is uniformly exp onential ly er go dic and the invariant distributions have exp onential tails. In p articular, ther e exists C 0 such that L u V  x  ≤ C 0 − %ε 3 m V ( x ) ∀ ( x, u ) ∈ R m × ∆ . (4.17) In p articular,   P v t ( x, · ) − π v ( · )   V ≤ C γ V ( x ) e − γ t , ∀ x ∈ R m , ∀ t ≥ 0 , (4.18) wher e P v t ( x, d y ) denote the tr ansition pr ob ability of { X t } t ≥ 0 under v = ( u c , u s ) and π v denotes the invariant pr ob ability me asur e of { X t } t ≥ 0 under the c ontr ol v . Pr o of. In Lemmas 4.3 and 4.4 in the section which follows, we establish (4.14a) and (4.14b) for these net works. Thus the pro of of (4.17) follows directly from Lemma 4.1 .  4.3. Three tec hnical lemmas. In this section, w e pro ve (4.14a) and (4.14b) for the netw orks under consideration which implies Theorem 4.2 . Even though the ‘N’ netw ork is a sp ecial case of the netw orks with a dominan t serv er po ol, w e first establish the result for this net w ork in Lemma 4.2 as understanding the results in R 2 will definitely help the reader in understanding the equations in Lemmas 4.3 and 4.4 . 4.3.1. The c ase of t he ‘N’ network. Here, m = 2, and the matrices B i , i = 1 , 2, in (2.30) are giv en b y B 1 = µ 11 0 0 µ 21 ! , and B 2 = 0 µ 12 − µ 11 0 0 ! . Th us, using (4.13) , the drift b : R 2 → R 2 for the ‘N’ netw ork is giv en by b ( x, u ) = − % 2 µ 11 µ 21 ! − µ 11 0 0 µ 21 !  x − h e, x i + u c  + ( µ 12 − µ 11 ) u s 2 0 ! h e, x i − . (4.19) Note that for the ‘N’ net work, w e ha ve Ψ( x ) = ψ ε ( x 1 ) µ 11 + ψ ε ( x 2 ) µ 21 b y (4.5) . Recall the definition of the cub e K r in (1.1) . W e ha ve the following lemma that verifies the drift inequalities (4.14a) and (4.14b) for the ‘N’ netw ork. UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 19 Lemma 4.2. Consider an ‘N’ network satisfying % > 0 . L et δ ∈ (0 , 1) , θ ≥ θ 0 : = 2( η ∨ η − 1 ) , with η : = µ 12 µ 11 , and V ( x ) b e as in (4.7) . Then, (4.14a) and (4.14b) hold with m = 2 . Pr o of. T o simplify the notation w e define F i ( x, u ) : = 1 V i ( x )  b ( x, u ) , ∇ V i ( x )  , i = 1 , 2 . (4.20) W e use (4.19) , and apply (4.8) and the inequalities % 2 P i ∈I ψ 0 ε ≤ %ε , and ψ 0 ε ( − x 1 )( η − 1) u s 2 h e, x i ≤ − ε (1 − η ) + h e, x i ≤ ε (1 − η ) + k x − k 1 ∀ ( x, u ) ∈ K − 0 × ∆ . to obtain 1 θ F 1 ( x, u ) = % 2 X i ∈I ψ 0 ε ( − x i ) + X i ∈I ψ 0 ε ( − x i ) x i − ψ 0 ε ( − x 1 )( η − 1) u s 2 h e, x i − ≤ 1 + %ε − ε k x − k 1 + ε (1 − η ) + k x − k 1 ≤ (1 + %ε ) − ε ( η ∧ 1) k x − k 1 ≤ (1 + %ε ) − ε 2 ( η ∧ 1) k x k 1 ∀ ( x, u ) ∈ K − 0 × ∆ . (4.21) Similarly , we hav e F 2 ( x, u ) = − % 2 X i ∈I ψ 0 ε ( x i ) − X i ∈I ψ 0 ε ( x i ) x i + ψ 0 ε ( x 1 )( η − 1) u s 2 h e, x i − ≤ ε  1 + ( η − 1) +  k x k 1 ≤ ε ( η ∨ 1) k x k 1 ∀ ( x, u ) ∈ K − 0 × ∆ . (4.22) Note that, due to (4.11) and the choice of θ , w e ha ve V 1 ≥ V 2 2 on K − 0 . Thus, since V 1 has exp onential gro wth in k x k 1 on K − 0 , com bining (4.21) and (4.22) and c ho osing an appropriate cube K r , w e obtain h b ( x, u ) , ∇ V ( x ) i ≤  θ (1 + %ε ) − ε 4 ( η ∧ 1) k x k 1  V ( x ) ∀ ( x, u ) ∈ ( K − 0 \ K r ) × ∆ . (4.23) W e contin ue with estimates on K + 0 . A straigh tforward calculation shows that 1 θ F 1 ( x, u ) = % 2 X i ∈I ψ 0 ε ( − x i ) + X i ∈I ψ 0 ε ( − x i ) x i − X i ∈I ψ 0 ε ( − x i ) u c i h e, x i F 2 ( x, u ) = − % 2 X i ∈I ψ 0 ε ( x i ) − X i ∈I ψ 0 ε ( x i ) x i + X i ∈I ψ 0 ε ( x i ) u c i h e, x i ∀ ( x, u ) ∈ K + 0 × ∆ . Again using (4.8) , we hav e 1 θ F 1 ( x, u ) ≤ 1 + %ε − ε k x − k 1 ∀ ( x, u ) ∈ K + 0 × ∆ . (4.24) W e break the estimate of F 2 in t wo parts. First, for any δ ∈ (0 , 1), using (4.8) , we obtain F 2 ( x, u ) ≤ − %ε 2 + 1 − ε k x + k 1 + ε h e, x i ≤ − %ε 2 + 1 − ε k x − k 1 ∀ ( x, u ) ∈  K + 0 \ K + δ  × ∆ . (4.25) Com bining (4.24) and (4.25) , we get h b ( x, u ) , ∇ V ( x ) i ≤  θ (1 + %ε ) − ε (1 − δ ) 2 k x k 1  V ( x ) ∀ ( x, u ) ∈  K + 0 \ K + δ  × ∆ , (4.26) and θ ≥ 1, where we use the fact that k x − k 1 ≥ 1 − δ 2 k x k 1 on K + 0 \ K + δ b y (4.10) . Thus (4.14a) follo ws b y (4.23) and (4.26) . 20 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Next, using (4.8) , we hav e F 2 ( x, u ) ≤ − %ε 2 ∀ ( x, u ) ∈ K + δ × ∆ , and this completes the pro of.  4.3.2. The c ase of networks with a dominant server p o ol. Consider the class of netw orks describ ed in Subsection 2.4.1 . W e hav e the following lemma. Lemma 4.3. Consider a network with a dominant server p o ol, such that % > 0 . L et δ ∈ (0 , 1) , θ ≥ θ 0 : = 2 max i µ i 1 min i µ i 1 , and V ( x ) b e as in (4.7) . Then, (4.14a) and (4.14b) hold. Pr o of. The metho d we follo w is analogous to the pro of of Lemma 4.2 . Recall the definitions in (4.20) . A straigh tforward calculation using (2.33) sho ws that 1 θ F 1 ( x, u ) = % m X i ∈I ψ 0 ε ( − x i ) + X i ∈I ψ 0 ε ( − x i )  x i − u c i h e, x i +  − h e, x i − X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( − x i )( η ij − 1) u s j , F 2 ( x, u ) = − % m X i ∈I ψ 0 ε ( x i ) − X i ∈I ψ 0 ε ( x i )  x i − u c i h e, x i +  + h e, x i − X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( x i )( η ij − 1) u s j . Let η : = min ij η ij , and ¯ η : = max ij η ij . Noting that −h e, x i − X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( − x i )( η ij − 1) u s j ≤ ε (1 − η ) + h e, x i − , it is easy to verify using (4.9) and (4.10) that 1 θ F 1 ( x, u ) ≤ %ε + m 2 − ε k x − k 1 + ε (1 − η ) + k x − k 1 ≤ %ε + m 2 − ε (1 − δ ) 2 ( η ∧ 1) k x k 1 ∀ ( x, u ) ∈  K + δ  c × ∆ . (4.27) Note that the drift equations on K + 0 are similar to those of the ‘N’ model, with the only exception that % 2 is replaced b y % m , and the sum ranges from i = 1 , . . . , m instead of i = 1 , 2. Hence, we obtain F 2 ( x, u ) ≤ − %ε m + m 2 − ε k x + k 1 + ε h e, x i ≤ − %ε m + m 2 − ε k x − k 1 ≤ − %ε m + m 2 − ε (1 − δ ) 2 k x k 1 ∀ ( x, u ) ∈  K + 0 \ K + δ  × ∆ , (4.28) and F 2 ( x, u ) ≤ − %ε m − ε h e, x i + ε h e, x i ≤ − %ε m ∀ ( x, u ) ∈ K + δ × ∆ , (4.29) for an y δ ∈ (0 , 1). The c hoice of θ implies that V 1 ≥ V 2 2 on K − 0 . Thus (4.14a) holds by (4.27) and (4.28) , while (4.29) is equiv alent to (4.14b) .  4.3.3. The c ase of networks with class-dep endent servic e r ates. Consider the class of netw orks de- scrib ed in Subsection 2.4.2 . Suc h netw orks hav e a limiting diffusion with the same drift structure studied in [ 4 ], and that pap er sho ws that when % > 0, then the diffusion (and the prelimit) is uniformly exp onen tially ergo dic in the presence or absence of abandonment. How ever, the pro of of uniform exp onential ergo dicit y of the prelimit for mo dels with class-dep endent service rates do es not see m to carry through with the Lyapuno v function used in [ 4 ]. Th us for the sake of proving the result for the n th system in Section 5 , we adopt here the Lyapuno v function in (4.7) . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 21 Lemma 4.4. Consider a network satisfying µ ij = µ i for al l i ∈ I , and % > 0 . L et δ ∈ (0 , 1) , and θ ≥ θ 0 : = 2 µ max µ min . Then, (4.14a) and (4.14b) hold. Pr o of. A simple calculation using (2.34) sho ws that 1 θ F 1 ( x, u ) = % m X i ∈I ψ 0 ε ( − x i ) + X i ∈I ψ 0 ε ( − x i )  x i − u c i h e, x i +  , F 2 ( x, u ) = − % m X i ∈I ψ 0 ε ( x i ) − X i ∈I ψ 0 ε ( x i )  x i − u c i h e, x i +  . (4.30) Using (4.30) , we obtain 1 θ F 1 ( x, u ) ≤ %ε + m 2 − ε 2 k x k 1 ∀ ( x, u ) ∈ K − 0 × ∆ , Therefore, (4.14a) holds on K − 0 × ∆ by this inequality and the choice of θ . On K + 0 × ∆ , the equations in (4.30) are identical to the corresp onding ones for a netw ork with a dominan t server p o ol, for whic h the result has already b een established in Lemma 4.3 . This completes the pro of.  5. Uniform exponential ergodicity of the n th system In this section we show that if % n > 0 then the prelimit of a netw ork with a dominant server p o ol, or with class-dep endent service rates, is uniformly exp onen tially ergo dic and the in v arian t distributions ha ve exp onential tails. Recall that { ˜ µ i , i ∈ I } are the elemen ts of the diagonal matrix B n 1 in (2.25) . Throughout this section V denotes the function in (4.7) , with ε giv en by ε = ε n : = % n 3 m X i ∈I 1 n λ n i (3 ˜ µ n i + 2) ( ˜ µ n i ) 2 ! − 1 exp − 1 √ n X i ∈I 1 ˜ µ n i ! . (5.1) Recall the definition of the op erator b L n ˆ z in (2.15) , and the definitions of S n and e Z n ( ˆ x ) in (2.12) and Definition 2.1 . W e start with the following simple assertion. Lemma 5.1. L et V b e the function in (4.7) with ε as in (5.1) , and θ fixe d at some value. Supp ose that for any δ ∈ (0 , 1) ther e exist p ositive c onstants c 0 and c 1 such that the drift b n in (2.25) satisfies  b n ( ˆ x, ˆ z ) , ∇ V ( ˆ x )  ≤ c 0 − εc 1 k x k 1 V ( ˆ x ) ∀ ˆ x ∈ S n \ K + δ , ∀ ˆ z ∈ e Z n ( ˆ x ) ,  b n ( ˆ x, ˆ z ) , ∇ V 2 ( ˆ x )  ≤ − % n ε n 2 m V 2 ( ˆ x ) ∀ ˆ x ∈ S n ∩ K + δ , ∀ ˆ z ∈ e Z n ( ˆ x ) . (5.2) Then, ther e exists a c onstant b C 0 such that b L n ˆ z V  ˆ x  ≤ b C 0 − % n ε n 4 m V ( ˆ x ) ∀ ˆ x ∈ S n , ∀ ˆ z ∈ e Z n ( ˆ x ) . (5.3) Pr o of. A simple calculation shows that Z 1 0 (1 − t ) ∂ x i x i V 2  ˆ x ± t √ n e i  d t ≤ ε 2 n 2 X i ∈I (3 ˜ µ n i + 1) (2 ˜ µ n i ) 2 ! exp 1 √ n X i ∈I 1 ˜ µ n i ! V 2 ( ˆ x ) . Th us, using (3.5) to express the first and second order incremental quotients in (2.15) , w e obtain b L n ˆ z V 2  ˆ x  ≤ % n ε n 4 m V 2 ( ˆ x ) +  b n ( ˆ x, ˆ z ) , ∇ V 2 ( ˆ x )  . The rest follo ws as in the pro of of Lemma 4.1 by selecting δ sufficien tly close to 1.  22 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG R emark 5.1 . Recall (2.25) . In direct analogy to Remark 4.2 , if we let ζ n = % n m e + ( B n 1 ) − 1 ` n , with % n as in (2.35) , then a mere translation of the origin of the form ˜ X n = ˆ X n + ζ n results in a diffusion of with the same drift as (2.25) , except that the vector ` n gets replaced b y ` n = − % n m ( B 1 ) n e . Therefore, w e may assume without any loss of generalit y that the drift in (2.25) tak es the form b n ( ˆ x, ˆ z ) = − % n m B n 1 e − B n 1  ˆ x − h e, ˆ x i + u c  + B n 2 u s h e, ˆ x i − + ˆ ϑ n ( ˆ x, ˆ z )  B n 1 u c + B n 2 u s  . (5.4) Note that this centering has the effect of translating the ‘equilibrium’ allo cations z n ij giv en in (2.8) . Since this translation is of O ( n − 1 / 2 ), it has no effect on the results for large n . Ho w ever, in the interest of pro viding precise estimates we calculate the new v alues of z n ij . Note that h e, ζ n i = 0, and recall the map Φ in (2.22) . Let ˇ z n ij = Φ( ζ n , 0). Then, the cen tering of ˆ x that results in (5.4) is giv en by (compare with (2.8) ) ¯ z n ij = 1 n ξ ∗ ij N n j + ˇ z n ij √ n , ¯ x n i : = X j ∈J ¯ z n ij . (5.5) Throughout this section the family { z n ij ( i, j ) ∈ E } is as given in (5.5) . R emark 5.2 . Recall Definition 2.1 and (2.26) . Let ˆ z ∈ e Z n ( ˆ x ). Then, ˆ ϑ n ( ˆ x, ˆ z ) = ˆ ϑ n ∗ ( ˆ x ) = 0 for all ˆ x ∈ ˘ X n , and in view of (5.4) , for any ˆ x ∈ ˘ X n , there exists u = u ( ˆ x, ˆ z ) ∈ ∆ such that b n ( ˆ x, ˆ z ) = − % n m B n 1 e − B n 1  ˆ x − h e, ˆ x i + u c  + B n 2 u s h e, ˆ x i − . (5.6) In view of Lemma 5.1 , and using Lemmas 4.2 to 4.4 , it is clear that F oster–Lyapuno v equations for L u carry ov er to analogous equations for b L n ˆ z on ˘ X n uniformly ov er SW C p olicies. Ho wev er, even though ˘ X n fills the whole space as n → ∞ , b and b n differ in functional form when ˆ ϑ n ( ˆ x, ˆ z ) 6 = 0, and this makes the stabilit y analysis of m ulticlass m ulti-p ool net w orks muc h harder than the ‘V’ netw ork studied in [ 4 ]. Notation 5.1. Let ε n and ¯ z n ij as in (5.1) and (5.5) , resp ectively . F or a net work with a dominan t serv er p o ol as in Subsection 4.3.2 define n 0 : = max  n ∈ N : 1 √ n ≥ ε n min i ∈I ¯ z n i 1  , (5.7) while for net work with class-dep enden t service rates, we let n 0 : = max  n ∈ N : 1 √ n ≥ ε n 2 m min i ∼ j ¯ z n ij  . (5.8) Since { ε n } and { ¯ z n i 1 } are b ounded aw ay from 0 by the con vergence of the parameters in (2.2) , the n umber n 0 is finite. The next theorem is the main result for the uniform exp onen tial ergo dicity of the prelimit pro cesses. Recall (2.35) . Notice the similarities b et ween the results in Theorem 5.1 for the diffusion- scaled pro cesses and the results in Theorem 4.2 and Remark 4.1 for the diffusion limit of the net works under consideration. Theorem 5.1. Assume that % n > 0 , and let n 0 b e as in Notation 5.1 . Then the pr elimit dynamics of any network with a dominant server p o ol or with class-dep endent servic e r ates ar e uniformly exp onential ly er go dic and the invariant distributions have exp onential tails for al l n > n 0 . In p articular, due to the c onver genc e of the p ar ameters, ther e exists b C 0 indep endent of n such that b L n ˆ z V  ˆ x  ≤ b C 0 − % n ε n 4 m V ( ˆ x ) ∀ ˆ x ∈ S n , ∀ ˆ z ∈ e Z n ( ˆ x ) . (5.9) wher e V and ε ar e as in (4.7) and (5.1) . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 23 In addition, with P n, ˆ z t and π n ˆ z denoting r esp e ctively the tr ansition pr ob ability and the stationary distribution of ˆ X n ( t ) under a p olicy ˆ z ∈ e Z n , ther e exist p ositive c onstants γ and C γ not dep ending on n ≥ 0 or ˆ z , such that   P n, ˆ z t ( ˆ x, · ) − π n ˆ z ( · )   V ≤ C γ V ( ˆ x ) e − γ t , ∀ ˆ x ∈ X n , ∀ t ≥ 0 . (5.10) Pr o of. In Lemmas 5.4 and 5.5 in the section which follows, we establish (5.2) for these netw orks. Th us the pro of of (5.9) follows directly from Lemma 5.1 . Since the process ˆ X n is irreducible and ap erio dic under an y stationary Mark o v sc heduling ˆ z ∈ e Z n (see Definition 2.1 ), a con vergence prop erty completely analogous to (4.3) follows from (5.9) . The pro of of this fact is identical to [ 4 , Theorem 2.1(b)].  R emark 5.3 . Using (4.18) and (5.10) , it is clear that under an y sc heduling p olicy ˆ z ∈ e Z n with a corresp onding control v , the stationary distribution of the diffusion-scaled pro cess ˆ X n ( t ) conv erges to that of the limiting diffusion { X t } t ≥ 0 for the t wo classes of net works, that is, π n ˆ z ( · ) → π v ( · ) , as n → ∞ . (5.11) That is, the interc hange of limits prop ert y holds. 5.1. F our tec hnical lemmas. In this section, we establish the technical results used in the pro of of Theorem 5.1 . Let e X n : =  ˆ x ∈ S n : ˆ ϑ n ∗ ( ˆ x ) 6 = 0  , (5.12) with ˆ ϑ n ∗ as in Definition 2.1 . As seen in Subsection 2.3 , the set ˘ X n in (2.26) is contained in S n \ e X n . In establishing (5.2) on S n \ e X n , the results in Section 4 pa ve the wa y , since the drift of of the con trolled generator b L n ˆ z o ver the class of SWC stationary Marko v p olicies e Z n (see (5.6) ) has the same functional form as the drift of the diffusion in (2.31) . So it remains to establish (5.2) in e X n . W e start by establishing a b ound for ˆ ϑ n in (2.18) o ver all SWC p olicies. As done earlier in the interest of notational econom y , we suppress the dep endence on n in the diffusion scaled v ariables ˆ x n and ˆ z n in (2.11) . Lemma 5.2. Ther e exists a numb er κ n ◦ < 1 dep ending only on the p ar ameters of the network such that ˆ ϑ n ( ˆ x, ˆ z ) = ˆ ϑ n ∗ ( ˆ x ) ≤ κ n ◦  k ˆ x + k 1 ∧ k ˆ x − k 1  ∀ ˆ z ∈ e Z n ( ˆ x ) . (5.13) In addition, due to the c onver genc e of the p ar ameters in (2.2) , such a c onstant κ ◦ < 1 may b e sele cte d which do es not dep end on n . Before pro ceeding to the pro of of Lemma 5.2 , we provide an in terpretation of (5.13) . Recall from (2.10) that ˆ x i = 1 √ n ( x i − P j ∈J ξ ∗ ij N n j ), and observ e that ˆ x i is p ositive if the total n umber of class i customers exceeds the total num b er of servers assigned to class i in the fluid equilibrium  P j ∈J ξ ∗ ij N n j  and is negative otherwise. This means that the quantit y k ˆ x + k 1 ∧ k ˆ x − k 1 on the right hand side of (5.13) represen ts the minim um of the total num b er of customers in the queues and the total n umber of idle servers if the serv ers are assigned to customer classes according to the fluid equilibrium. Recall also that ˆ ϑ n ∗ represen ts the minimum of the total num b er of customers in the queues and the total num b er of idle serv ers under a SWC p olicy . The result in Lemma 5.2 is no w clear, that is, the minimum of the total n umber of customers in the queues and idle servers is smaller under a SWC p olicy compared to a p olicy that assigns servers according to the fluid equilibrium. 24 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Pr o of. Let ˆ x ∈ e X n , ˆ z ∈ e Z n ( ˆ x ), and define e J : = ( j ∈ J : X j ∈J ( i ) ˆ z ij < 0 ) , and e I : =  i ∈ I : ( i, j ) ∈ E for some j ∈ e J  , and e E : =  ( i, j ) ∈ E : ( i, j ) ∈ e I × e J  . W ork conserv ation implies that x n i = P j ∈J ( i ) z n ij for all i ∈ e I . Let ˆ ı ∈ I b e suc h that ˆ q ˆ ı > 0, and consider the unique path (since the graph of the netw ork is a tree) connecting ˆ ı and e I , that is, a path ˆ ı → j 1 → i 1 → j 2 . . . → j k → ˜ ı , with j ` ∈ J \ e J for ` = 1 , . . . , k , i ` ∈ I \ e I for ` = 1 , . . . , k − 1, and ˜ ı ∈ e I . W e claim that z n ˜ ı,j k = 0, or equiv alently , that ˆ z ˜ ı,j k = − ¯ z n ij / √ n , with ¯ z n ij as defined in (5.5) . If not, then w e can mov e a job of class ˜ ı from p o ol j k to some p o ol in e J , and pro ceeding along the path to place one additional job from class ˆ ı into service, thus con tradicting the hypothesis that ˆ z ∈ e Z n ( ˆ x ). Remo ving all such paths, w e are left with a strict subnetw ork (p ossibly disconnected) G ◦ =  I ◦ ∪ J ◦ , E ◦ ), with I ◦ ⊃ e I , J ◦ ⊃ e J , and E ◦ : =  ( i, j ) ∈ E : ( i, j ) ∈ I ◦ × J ◦  , suc h that x n i = X j ∈J ( i ) ∩J ◦ z n ij , ∀ i ∈ I ◦ . (5.14) Let E 0 ◦ : =  I ◦ × ( J \ J ◦ )  ∩ E . By (5.14) we ha ve X ( i,j ) ∈E 0 ◦ z n ij = 0 . Th us we hav e k ˆ x − k 1 ≥ − X i ∈I ◦ ˆ x i = √ n X ( i,j ) ∈E 0 ◦ ¯ z n ij − X ( i,j ) ∈E ◦ ˆ z ij ≥ √ n X ( i,j ) ∈E 0 ◦ ¯ z n ij + ˆ ϑ n ∗ ( ˆ x ) (5.15) b y the construction ab o ve. By (5.15) , we obtain ˆ ϑ n ∗ ( ˆ x ) ≤ P ( i,j ) ∈E ◦ ˆ z ij √ n P ( i,j ) ∈E 0 ◦ ¯ z n ij − P ( i,j ) ∈E ◦ ˆ z ij X i ∈I ◦ ˆ x i ≤ − P ( i,j ) ∈E ◦ ¯ z n ij P ( i,j ) ∈E ◦ ¯ z n ij + P ( i,j ) ∈E 0 ◦ ¯ z n ij X i ∈I ◦ ˆ x i . (5.16) Similarly , k ˆ x + k 1 ≥ X i ∈I \I ◦ ˆ x i ≥ √ n X ( i,j ) ∈E 0 ◦ ¯ z n ij + ˆ ϑ n ∗ ( ˆ x ) . (5.17) Using the b ound ˆ ϑ n ∗ ( ˆ x ) ≤ √ n P ( i,j ) ∈E ◦ ¯ z n ij w e obtain from (5.17) that ˆ ϑ n ∗ ( ˆ x ) ≤  P i ∈I \I ◦ ˆ x i − √ n P ( i,j ) ∈E 0 ◦ ¯ z n ij  ∧ √ n P ( i,j ) ∈E ◦ ¯ z n ij P i ∈I \I ◦ ˆ x i X i ∈I \I ◦ ˆ x i ≤ P ( i,j ) ∈E ◦ ¯ z n ij P ( i,j ) ∈E ◦ ¯ z n ij + √ n P ( i,j ) ∈E 0 ◦ ¯ z n ij X i ∈I \I ◦ ˆ x i . (5.18) It should b e now clear how to construct κ n ◦ . F or any given subset J 0 ( J , let I J 0 : = ∪ j ∈J 0 I ( j ) , E J 0 : =  ( i, j ) ∈ E : ( i, j ) ∈ I J 0 × J 0  , UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 25 and E 0 J 0 : =  I J 0 × ( J \ J 0 )  ∩ E , and define κ n ◦ : = max J 0 ( J P ( i,j ) ∈E J 0 ¯ z n ij P ( i,j ) ∈E J 0 ¯ z n ij + P ( i,j ) ∈E 0 J 0 ¯ z n ij . Then the result clearly follows from (5.16) and (5.18) since k ˆ x − k 1 ≥ − X i ∈I ◦ ˆ x i , and k ˆ x + k 1 ≥ X i ∈I \I ◦ ˆ x i . This completes the pro of.  R emark 5.4 . W e pro vide an example of a SW C p olicy for the ‘N’ netw ork with tw o classes of customers and tw o serv er p o ols where class 1 can b e served by b oth serv er p o ols while class 2 can only served by p o ol 2 to explain the analysis in the pro of of Lemma 5.2 . The SW C p olicy is given b y z 11 ( x ) = x 1 ∧ N n 1 z 12 ( x ) = ( ( x 1 − N n 1 ) + ∧ ξ ∗ 12 N n 2 if x 2 ≥ ξ ∗ 22 N n 2 ( x 1 − N n 1 ) + ∧ ( N n 2 − x 2 ) otherwise, z 22 ( x ) = ( x 2 ∧ ξ ∗ 22 N n 2 if x 1 ≥ N n 1 + ξ ∗ 12 N n 2 x 2 ∧  N n 2 − ( x 1 − N n 1 ) +  otherwise. This is a priorit y p olicy in whic h customer class 1 prefers p o ol 1 o ver p o ol 2. This means that customer class 1 can use servers in p o ol 2 only if there are no idle servers in p o ol 1. Recall from (2.9) and (2.10) that ˆ q i ≥ 0 for all i ∈ I . Recall also that w e are only considering w ork-conserving p olicies and that ˆ ϑ n represen ts the minimum of the total n umber of customers in the queues and the total num b er of idle servers. F or the ‘N’ net work with tw o classes of customers and t wo server p o ols, w e hav e the following cases: Case 1 : ˆ q 1 > 0 and ˆ q 2 ≥ 0, in which case we hav e no idle servers and hence ˆ ϑ n ∗ ( ˆ x ) = 0. Case 2 : ˆ q 1 = 0 and ˆ q 2 = 0. Again, this is a trivial case and means that we ha ve no customers in the queues whic h implies that ˆ ϑ n ∗ ( ˆ x ) = 0. Case 3 : ˆ q 1 = 0 and ˆ q 2 > 0. This is actually the case that requires some analysis. Here again we ha ve the following cases: a : ˆ y 1 = 0 which means that there are no idle servers in p o ol 1 and hence ˆ ϑ n ∗ ( ˆ x ) = 0. This is b ecause ˆ q 2 > 0 which implies that ˆ y 2 = 0 under w ork conserv ation. b : ˆ y 1 > 0 which means that there are some idle servers in p o ol 1. This is the case analyzed in the pro of of Lemma 5.2 . (Observ e the following notation in the pro of: e J = { 1 } , e I = { 1 } , ˆ ı = 2, the path is class 2 → p o ol 2 → class 1, and G ◦ = { (1 , 1) } .) Note that since we are using a SWC p olicy , this means that no class 1 customers are b eing served in p o ol 2 b ecause otherwise, one can mov e a class 1 customer from p o ol 2 to p ool 1 and get a smaller ˆ ϑ n . Hence, ˆ ϑ n ∗ ( ˆ x ) = 1 √ n  ( N n 1 − x 1 ) ∧ ( x 2 − N n 2 )  , where ( N n 1 − x 1 ) is the total n umber of idle server and ( x 2 − N n 2 ) is the total n umber of customers in the queues. This is b ecause the only p o ol having idle servers is p o ol 1 and the only class having customers in the queue is class 2 ( ˆ q 1 = 0). 26 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Recall that ˆ x i = 1 √ n ( x i − P j ∈J ξ ∗ ij N n j ). This means in this case that ˆ x + 1 = 0 ; ˆ x − 1 = 1 √ n ( N n 1 + ξ ∗ 12 N n 2 − x 1 ) ˆ x + 2 = 1 √ n ( x 2 − ξ ∗ 22 N n 2 ) ; ˆ x − 2 = 0 . Therefore, w e hav e the following equation k ˆ x + k 1 ∧ k ˆ x − k 1 = 1 √ n  ( N n 1 + ξ ∗ 12 N n 2 − x 1 ) ∧ ( x 2 − ξ ∗ 22 N n 2 )  . Note also that (2.4) ( ξ ∗ 12 + ξ ∗ 22 = 1) implies that N n 1 − x 1 ≤ x 2 − N n 2 ⇐ ⇒ N n 1 + ξ ∗ 12 N n 2 − x 1 ≤ x 2 − ξ ∗ 22 N n 2 . It is now clear that there exists a constant κ n ◦ < 1 suc h that ˆ ϑ n ∗ < κ n ◦  k ˆ x + k 1 ∧ k ˆ x − k 1  where κ n ◦ = N n 1 − x 1 N n 1 + ξ ∗ 12 N n 2 − x 1 ∨ x 2 − N n 2 x 2 − ξ ∗ 22 N n 2 . This completes the analysis of all the cases. R emark 5.5 . It is easy to see that the estimates of the b ounds on ˆ ϑ n can b e impro ved. It is clear from (5.15) and (5.16) , that ˆ ϑ n ∗ ( ˆ x ) ≤ − κ n 0 X i ∈I 0 ˆ x i ! ∧ X i ∈I \I ◦ ˆ x i ! ∀ ˆ x ∈ e X n . Also, since there can b e at most P j ∈J N n j idle serv ers, it follows that e κ n ◦ ∈ (0 , 1), such that − X i ∈I 0 ˆ x i ≥ e κ n ◦ k ˆ x − k 1 ∀ ˆ x ∈ e X n , where the constan t e κ n ◦ ∈ (0 , 1) can b e selected as e κ n ◦ : = X j ∈J N n j ! − 1 min ( i,j ) ∈ E ξ ∗ ij N n j . Due to the con vergence of the parameters in (2.2) , e κ n ◦ is b ounded aw ay from 0 uniformly in n ∈ N . Note that for the ‘N’ netw ork this translates to e κ n ◦ = N n 1 ∧ ξ ∗ 12 N n 2 ∧ ξ ∗ 22 N n 2 N n 1 + N n 2 . Ev en though the ‘N’ netw ork is a sp ecial case of netw orks with a dominant server p o ol w e first establish the result for this netw ork in Lemma 5.3 in order to exhibit with simpler c alculations ho w Lemma 5.2 is applied. Throughout the pro ofs of Lemmas 5.3 to 5.5 we use the functions (compare with (4.20) ) F n i ( ˆ x, ˆ z ) : = 1 V i ( x )  b n ( ˆ x, ˆ z ) , ∇ V 1 ( ˆ x )  , i = 1 , 2 , and let n 0 b e as in Notation 5.1 . Moreo ver, w e suppress the dep endence on n in the v ariables ˆ q n , ˆ y n , and ˆ ϑ n in (2.17) and (2.18) , and from ε n in (5.1) . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 27 5.1.1. The diffusion-sc ale of the ‘N’ network. W e recall here [ 27 ]. Recall also that w e label the non-leaf server no de as j = 1 without loss of generality and hence w e present Stolyar’s work in [ 27 ] using our notation. In this work, Stoly ar considers the ‘N’ netw ork with O ( √ n ) safety staffing in p o ol 1, under the priorit y discipline that class 2 has priority in p o ol 1 and class 1 prefers p o ol 2, and shows tightness of the inv ariant distributions. First note that for an y stationary Marko v sc heduling p olicy z , suc h that class 2 has priorit y in p o ol 1 we hav e z n 21 ( x ) = x n 2 ∧ N n 1 , and it is clear that such a p olicy is SWC. The same applies to Mark ov p olicies under which class 1 prefers p o ol 2 (here z n 12 ( x ) = x n 1 ∧ N n 2 ). As a result, SWC p olicies are more general than the particular p olicy considered in [ 27 ]. Recall that the matrices B n 1 and B n 2 in the drift (5.4) are given by B n 1 = µ n 11 0 0 µ n 21 ! , and B n 2 = 0 µ n 12 − µ n 11 0 0 ! (5.19) It is also worth noting here, that the spare capacity % n of the n th system is giv en by % n = − 1 √ n  λ n 1 µ n 11 + λ n 2 µ n 21 − µ n 12 N n 2 + µ n 11 ξ ∗ 11 N n 1 µ n 11 − µ n 21 ξ ∗ 21 N n 1 µ n 21  , with ξ ∗ 11 = λ 1 − µ 12 ν 2 µ 11 ν 1 , and ξ ∗ 21 = λ 2 µ 21 ν 1 . This is clear by (2.14) , (2.35) , and (5.19) , together with [ 6 , Equation (2)]. W e let η n : = µ n 12 µ n 11 . Lemma 5.3. Consider the ‘N’ network, and assume that % n > 0 . Then for any θ ≥ θ n 0 : = µ n 11 ∨ µ n 21 µ n 11 ∧ µ n 21 , and δ ∈ (0 , 1) , ther e exist p ositive c onstants c 0 and c 1 such that (5.2) holds for al l n ≥ n 0 . Pr o of. Recall here that j = 1 is the non-leaf server po ol. As discussed earlier, it suffices to establish (5.2) in e X n . It is clear that ˆ x − 1 = ˆ y 2 + √ n ¯ z 11 , and ˆ x 2 = ˆ q 2 + √ n ¯ z 11 for all ˆ z ∈ ˆ Z n ( ˆ x ) and ˆ x ∈ e X n , with e X n as defined in (5.12) . Hence u c 1 = 0, u c 2 = 1, u s 2 = 1, and u s 1 = 0. Note here that SW C p olicies are interpreted through u c and u s as follo ws: idle servers are only allow ed in p o ol 1 and customer queues are only allow ed in class 2 which means that class 1 m ust use all the servers in p o ol 2 b efore using servers in p o ol 1. Also, b y the definitions of ψ ε and n 0 w e hav e ψ 0 ε ( ˆ x 1 ) = ψ 0 ε ( − ˆ x 2 ) = 0 ∀ ˆ x ∈ e X n , ∀ n ≥ n 0 . (5.20) By (5.4) and (5.19) , we ha ve 1 θ F n 1 ( ˆ x, ˆ z ) = % n 2 X i ∈I ψ 0 ε ( − ˆ x i ) + X i ∈I ψ 0 ε ( − ˆ x i )  ˆ x i − u c i h e, ˆ x i +  − ψ 0 ε ( − ˆ x 1 )( η n − 1) h e, ˆ x i − − ˆ ϑ n  ψ 0 ε ( − ˆ x 2 ) + ψ 0 ε ( − ˆ x 1 )( η n − 1)  , (5.21) F n 2 ( ˆ x, ˆ z ) = − % n 2 X i ∈I ψ 0 ε ( ˆ x i ) − X i ∈I ψ 0 ε ( ˆ x i )  ˆ x i − u c i h e, ˆ x i +  + ψ 0 ε ( ˆ x 1 )( η n − 1) h e, ˆ x i − + ˆ ϑ n  ψ 0 ε ( ˆ x 2 ) + ψ 0 ε ( ˆ x 1 )( η n − 1)  . (5.22) 28 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG Using the fact that ˆ x − 1 − ˆ x 2 = ˆ y 2 − ˆ q 2 , and ˆ ϑ n = ˆ q 2 when h e, ˆ x i ≤ 0, we obtain from (5.20) and (5.21) that 1 θ F n 1 ( ˆ x, ˆ z ) = % n 2 ψ 0 ε ( − ˆ x 1 ) − ψ 0 ε ( − ˆ x 1 ) ˆ x − 1 − ψ 0 ε ( − ˆ x 1 )( η n − 1) h e, ˆ x i − − ψ 0 ε ( − ˆ x 1 )( η n − 1) ˆ ϑ n = ε  % n 2 − ˆ x − 1 − ( η n − 1)  ˆ x − 1 − ˆ x 2  − ( η n − 1) ˆ ϑ n  = ε  % n 2 − ˆ x − 1 − ( η n − 1) ˆ y 2  ≤ ε  % n 2 − ˆ x − 1 + (1 − η n ) + ˆ x 1  ≤ ε  % n 2 − ( η n ∧ 1) ˆ x − 1  ≤ ε % n 2 − ε 2 ( η n ∧ 1) k ˆ x k 1 ∀ ( ˆ x, ˆ z ) ∈  e X n ∩ K − 0  × ˆ Z n ( ˆ x ) , ∀ n ≥ n 0 . (5.23) Similarly , from (5.22) , w e obtain F n 2 ( ˆ x, ˆ z ) = − % n 2 X i ∈I ψ 0 ε ( ˆ x i ) − X i ∈I ψ 0 ε ( ˆ x i ) ˆ x i + ψ 0 ε ( ˆ x 1 )( η n − 1) h e, ˆ x i − + ˆ ϑ n  ψ 0 ε ( ˆ x 2 ) + ψ 0 ε ( ˆ x 1 )( η n − 1)  ≤ 1 − ε k ˆ x + k 1 + ε κ 0  k ˆ x − k 1 ∧ k ˆ x + k 1  ≤ 1 − ε (1 − κ 0 ) k ˆ x + k 1 ∀ ( ˆ x, ˆ z ) ∈  e X n ∩ K − 0  × ˆ Z n ( ˆ x ) , ∀ n ≥ n 0 , (5.24) where w e also use (4.8) and Lemma 5.2 . W e contin ue with the estimate on K + 0 . W e hav e 1 θ F n 1 ( ˆ x, ˆ z ) ≤ % n 2 ψ 0 ε ( − ˆ x 1 ) − ψ 0 ε ( − ˆ x 1 ) ˆ x − 1 − ψ 0 ε ( − ˆ x 1 )( η n − 1) ˆ ϑ n = ε  % n 2 − ˆ x − 1 − ( η n − 1) ˆ ϑ n  ≤ ε  % n 2 − ( η n ∧ 1) ˆ x − 1  ∀ ( ˆ x, ˆ z ) ∈  e X n ∩ K + 0  × ˆ Z n ( ˆ x ) , ∀ n ≥ n 0 , (5.25) where in the last inequality we also use Lemma 5.2 . W e break the estimate of F n 2 in tw o parts. First, using ( 4.8 ), ( 4.10 ), ( 5.20 ) and Lemma 5.2 , we obtain F n 2 ( ˆ x, ˆ z ) ≤ − % n ε 2 − ε ˆ x 2 + ε h e, ˆ x i + ε κ 0  ˆ x − 1 ∧ ˆ x 2  ≤ − % n ε 2 − ε (1 − κ 0 ) ˆ x − 1 ≤    − % n ε 2 − ε (1 − δ ) 2 (1 − κ 0 ) k ˆ x k 1 for ( ˆ x, ˆ z ) ∈  e X n ∩ ( K + 0 \ K + δ )  × ˆ Z n ( ˆ x ) − % n ε 2 for ( ˆ x, ˆ z ) ∈  e X n ∩ K + δ  × ˆ Z n ( ˆ x ) . (5.26) Th us, (5.2) follows by (4.10) and (5.23) – (5.26) . This completes the pro of.  5.1.2. The diffusion sc ale of networks with a dominant p o ol. W e describ e these net works exactly as in Subsection 4.3.2 where the dominant server p o ol is j = 1. W e first note that the spare capacity % n of the n th system is giv en by % n = − 1 √ n X i ∈I λ n i µ n i 1 − X i ∈I X j ∈J ( i ) µ n ij µ n i 1 ξ ∗ ij N n j ! , UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 29 where ξ ∗ ij satisfies X j ∈J ( i ) µ ij ξ ∗ ij ν j = λ i . This is again due to (2.4) , (2.14) , (2.34) , and (2.35) . Recall from (5.4) that the drift reduces to the following form: b n i ( ˆ x, ˆ z ) = − % n m µ n i 1 − µ n i 1  ˆ x i − u c i h e, ˆ x i +  + X j ∈J 1 ( i ) µ n i 1  η n ij − 1  u s j h e, ˆ x i − + ˆ ϑ n µ n i 1 u c i + X j ∈J 1 ( i ) µ n i 1  η n ij − 1  u s j ! , i ∈ I , (5.27) with η n ij : = µ n ij µ n i 1 for j ∈ J 1 ( i ) : = J ( i ) \ { 1 } and i ∈ I . In analogy to Subsection 4.3.2 , we define W e define ¯ η n : = max i ∈I max j ∈J 1 ( i ) η n ij , and η n : = min i ∈I min j ∈J 1 ( i ) η n ij . Lemma 5.4. Consider a network with a dominant server p o ol, and assume % n > 0 . Then for any θ ≥ θ n 0 : = 2 max i µ n i 1 min i µ n i 1 , and δ ∈ (0 , 1) , ther e exist p ositive c onstants c 0 and c 1 such that (5.2) holds for al l n ≥ n 0 . Pr o of. Supp ose ˆ x ∈ e X n . A simple calculation using (5.27) shows that 1 θ F n 1 ( ˆ x, ˆ z ) = % n m X i ∈I ψ 0 ε ( − ˆ x i ) + X i ∈I ψ 0 ε ( − ˆ x i )  ˆ x i − u c i h e, ˆ x i +  − X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( − ˆ x i )  η n ij − 1  u s j h e, ˆ x i − − ˆ ϑ n X i ∈I ψ 0 ε ( − ˆ x i ) u c i + X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( − ˆ x i )  η n ij − 1  u s j ! , (5.28) and F n 2 ( ˆ x, ˆ z ) = − % n m X i ∈I ψ 0 ε ( ˆ x i ) − X i ∈I ψ 0 ε ( ˆ x i )  ˆ x i − u c i h e, ˆ x i +  + X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( ˆ x i )  η n ij − 1  u s j h e, ˆ x i − + ˆ ϑ n X i ∈I ψ 0 ε ( ˆ x i ) u c i + X i ∈I X j ∈J 1 ( i ) ψ 0 ε ( ˆ x i )  η n ij − 1  u s j ! . (5.29) By (5.28) we obtain 1 θ F n 1 ( ˆ x, ˆ z ) ≤ % n ε + m 2 − ε X i ∈I ˆ x − i + ε  1 − η n  + h e, ˆ x i − + ε  1 − η n  + ˆ ϑ n ≤ % n ε + m 2 − ε X i ∈I ˆ x − i + ε  1 − η n  + X i ∈I  ˆ x − i − ˆ x + i  + k ˆ x − k 1 ∧ k ˆ x + k 1 ! ≤ % n ε + m 2 − ε  η n ∧ 1  k ˆ x − k 1 ≤ % n ε + m 2 − ε (1 − δ ) 2  η n ∧ 1  k ˆ x k 1 ∀ ( ˆ x, ˆ z ) ∈  e X n \ K + δ  × ˆ Z n ( ˆ x ) , ∀ n ∈ N , (5.30) 30 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG where w e used (4.8) in the first inequalit y , Lemma 5.2 in the second, and (4.10) in the fourth. Next, we estimate a b ound for F n 2 ( ˆ x, ˆ z ). Recall the definitions of I ◦ , J ◦ , E ◦ , and E 0 ◦ in the pro of of Lemma 5.2 . Since x ∈ e X n , w e ha ve u c i = 0 for all i ∈ I ◦ , and u s j = 0 for all j ∈ J c ◦ . SW C p olicies are interpreted here as follows: customer classes must use all the servers in the leaf p o ols av ailable to them b efore using servers in p ool j = 1. Additionally , ˆ x i ≤ − P ( i,j ) ∈E 0 ◦ ¯ z ij for i ∈ I ◦ , whic h implies that ψ 0 ε ( x i ) = 0 for all i ∈ I ◦ and n > n 0 , by Notation 5.1 . Hence, since P i ∈I c ◦ ˆ x i > 0, where I c ◦ ≡ I \ I ◦ , w e hav e X i ∈I ψ 0 ε ( ˆ x i )  ˆ x i − u c i h e, ˆ x i +  ≥ X i ∈I c ◦ ψ 0 ε ( ˆ x i ) ˆ x i − ε X i ∈I c ◦ ˆ x i − ε X i ∈I ◦ ˆ x i ≥ − ε X i ∈I ◦ ˆ x i (5.31) b y (4.12) . Using (5.29) together with Remark 5.5 and (5.31) , we obtain F n 2 ( ˆ x, ˆ z ) ≤ − % n ε m + ε ˆ ϑ n + ε X i ∈I ◦ ˆ x i + ˆ ϑ n X i ∈I ◦ X j ∈J 1 ( i ) ψ 0 ε ( ˆ x i )  η n ij − 1  u s j ! ≤ − % n ε m + ε (1 − κ n ◦ ) X i ∈I ◦ ˆ x i ≤    − % n ε m − ε (1 − δ ) 2 (1 − κ n ◦ ) e κ n ◦ k ˆ x k 1 for ( ˆ x, ˆ z ) ∈  e X n ∩ ( K + 0 \ K + δ )  × ˆ Z n ( ˆ x ) − % n ε 2 for ( ˆ x, ˆ z ) ∈  e X n ∩ K + δ  × ˆ Z n ( ˆ x ) , (5.32) for all n ≥ n 0 . Thus, the result follows b y (5.30) and (5.32) , noting also that the choice of θ implies that V 1 ≥ V 2 2 on K − 0 .  5.1.3. The diffusion-sc ale of networks with class-dep endent servic e r ates. Recall from Subsection 4.3.2 that the drift in (5.4) reduces to b n ( ˆ x, ˆ z ) = − % n m B n 1 e − B n 1  ˆ x − h e, ˆ x i + u c  + ˆ ϑ n ( ˆ x, ˆ z ) B n 1 u c . (5.33) where B n 1 = diag( µ n 1 , . . . , µ n m ). Thus, the spare capacit y % n is giv en by % n = − 1 √ n  X i ∈I λ n i µ n i − X i ∈I X j ∈J ( i ) ξ ∗ ij N n j  . Lemma 5.5. Supp ose that µ n ij = µ n i , for al l i ∈ I , and % n > 0 . Then, for any θ ≥ θ n 0 : = 2 µ n max µ n min , and δ ∈ (0 , 1) , the c onclusions of L emma 5.4 fol low. Pr o of. Supp ose ˆ x ∈ e X n . A simple calculation using (5.33) shows that 1 θ F n 1 ( ˆ x, ˆ z ) = % n m X i ∈I ψ 0 ε ( − ˆ x i ) + X i ∈I ψ 0 ε ( − ˆ x i )  ˆ x i − u c i h e, ˆ x i +  − ˆ ϑ n X i ∈I ψ 0 ε ( − ˆ x i ) u c i , (5.34) F n 2 ( ˆ x, ˆ z ) = − % n m X i ∈I ψ 0 ε ( ˆ x i ) − X i ∈I ψ 0 ε ( ˆ x i )  ˆ x i − u c i h e, ˆ x i +  + ˆ ϑ n X i ∈I ψ 0 ε ( ˆ x i ) u c i . (5.35) By (5.34) , we obtain 1 θ F n 1 ( ˆ x, ˆ z ) ≤ % n ε + m 2 − ε k ˆ x − k 1 ≤ % n ε + m 2 − ε (1 − δ ) 2 k ˆ x k 1 ∀ ( ˆ x, ˆ z ) ∈  e X n \ K − δ  × ˆ Z n ( ˆ x ) , ∀ n ≥ n 0 . In computing the analogous b ound to (5.32) , there is a difference here. It is not the case here that ψ 0 ε ( x i ) = 0 for all i ∈ I ◦ and n > n 0 . UNIFORM ST ABILITY OF LARGE-SCALE P ARALLEL SER VER NETWORKS 31 So instead, recalling that u c i = 0 for all i ∈ I ◦ , and since ˆ x ∈ K + 0 , w e write − X i ∈I ψ 0 ε ( ˆ x i )  ˆ x i − u c i h e, ˆ x i +  + ˆ ϑ n X i ∈I ψ 0 ε ( ˆ x i ) u c i ≤ − X i ∈I ψ 0 ε ( ˆ x i ) ˆ x i + ε h e, ˆ x i + ε ˆ ϑ n ≤ ε ˆ ϑ n − X ` ∈I ◦ ˆ x − ` ! − X i ∈I ◦ ψ 0 ε ( ˆ x i ) ˆ x i − X i ∈I c ◦ ψ 0 ε ( ˆ x i ) ˆ x i − ε X i ∈I c ◦ ˆ x i ! . (5.36) The third term on the right-hand side is nonp ositive by (4.12) . W e also ha ve ˆ ϑ n − X ` ∈I ◦ ˆ x − ` = − √ n X ( i,j ) ∈E 0 ◦ ¯ z n ij , (5.37) and − X i ∈I ◦ ψ 0 ε ( ˆ x i ) ˆ x i ≤ X ˆ x − i ≤ 1 2 m √ n min i ∼ j ¯ z n ij ˆ x − i ≤ 1 2 √ n min i ∼ j ¯ z n ij . (5.38) Therefore, b y (5.35) , (5.36) – (5.38) and Remark 5.5 , we obtain F n 2 ( ˆ x, ˆ z ) ≤ ε 2 ˆ ϑ n − X ` ∈I ◦ ˆ x − ` ! ≤ ε 2 (1 − κ n ◦ ) e κ n ◦ k ˆ x − k 1 ∀ ( ˆ x, ˆ z ) ∈  e X n ∩ K + 0 ) × ˆ Z n ( ˆ x ) . The rest follo ws as in Lemma 5.4 .  A cknowledgment This work is supp orted in part by the Army Research Office through grant W911NF-17-1-0019, and in part by NSF grants DMS-1715210, CMMI-1635410, and DMS/CMMI-1715875, and in part b y the Office of Na v al Research through grant N00014-16-1-2956 and was appro ved for public release under DCN #43-5454-19. References [1] Z. Aksin, M. Armony, and V. Mehrotra, The mo dern c al l c enter: a multi-disciplinary p ersp e ctive on op er ations management rese ar ch , Pro duction and Op erations Management 16 (2007), 665–688. [2] A. Arapostathis, A. Bisw as, and G. P ang, Er go dic c ontr ol of multi-class M / M / N + M queues in the Halfin-Whitt r e gime , Ann. Appl. Probab. 25 (2015), no. 6, 3511–3570. MR 3404643 [3] A. Arap ostathis, V. S. Bork ar, and M. K. Ghosh, Er go dic c ontr ol of diffusion pr o cesses , Encyclop edia of Math- ematics and its Applications, vol. 143, Cambridge Universit y Press, Cambridge, 2012. MR 2884272 [4] A. Arapostathis, H. Hmedi, and G. Pang, On uniform exp onential er go dicity of Markovian multiclass many-server queues in the Halfin–Whitt r e gime , Math. Op er. Res. (2020). to app ear. [5] A. Arap ostathis and G. Pang, Er go dic diffusion c ontr ol of multiclass multi-p o ol networks in the Halfin-Whitt r e gime , Ann. Appl. Probab. 26 (2016), no. 5, 3110–3153. MR 3563203 [6] A. Arap ostathis and G. Pang, Infinite-horizon average optimality of the N-network in the Halfin-Whitt r e gime , Math. Op er. Res. 43 (2018), no. 3, 838–866. MR 3846075 [7] A. Arap ostathis and G. Pang, Infinite horizon asymptotic aver age optimality for lar ge-sc ale p ar al lel server net- works , Sto c hastic Pro cess. Appl. 129 (2019), no. 1, 283–322. MR 3906999 [8] A. Arap ostathis, G. Pang, and N. Sandri´ c, Er go dicity of a L´ evy-driven SDE arising fr om multiclass many-server queues , Ann. Appl. Probab. 29 (2019), no. 2, 1070–1126. MR 3910024 [9] M. Armony, S. Israelit, A. Mandelbaum, Y. N. Marmor, Y. Tseytlin, and G. B. Y om-T ov, On p atient flow in hospitals: a data-b ase d queueing-scienc e persp e ctive , Sto c h. Syst. 5 (2015), no. 1, 146–194. MR 3442392 [10] R. Atar, A diffusion mo del of sche duling c ontr ol in queueing systems with many servers , Ann. Appl. Probab. 15 (2005), no. 1B, 820–852. MR 2114991 32 HASSAN HMEDI, ARI ARAPOST A THIS, AND GUODONG P ANG [11] R. Atar, Sche duling c ontr ol for queueing systems with many servers: asymptotic optimality in he avy tr affic , Ann. Appl. Probab. 15 (2005), no. 4, 2606–2650. MR 2187306 [12] S. Borst, A. Mandelbaum, and M. I. Reiman, Dimensioning lar ge c al l c enters , Oper. Res. 52 (2004), no. 1, 17–34. MR 2066238 [13] A. Brav erman, J. G. Dai, and M. Miyaza wa, He avy tr affic appr oximation for the stationary distribution of a gener alize d Jackson network: the BAR appr o ach , Sto ch. Syst. 7 (2017), no. 1, 143–196. MR 3663340 [14] L. Bro wn, N. Gans, A. Mandelbaum, A. Sak ov, H. Shen, S. Zeltyn, and L. Zhao, Statistic al analysis of a telephone c al l c enter: a queueing-scienc e p ersp e ctive , J. Amer. Statist. Asso c. 100 (2005), no. 469. MR 2166068 [15] A. Budhira ja and C. Lee, Stationary distribution c onver gence for gener alize d Jackson networks in he avy tr affic , Math. Op er. Res. 34 (2009), no. 1, 45–56. MR 2542988 [16] A. B. Dieker and X. Gao, Positive r e curr enc e of pie c ewise Ornstein-Uhlenb e ck pr o c esses and c ommon quadr atic Lyapunov functions , Ann. Appl. Probab. 23 (2013), no. 4, 1291–1317. MR 3098433 [17] D. Gamarnik and A. L. Stolyar, Multiclass multiserver queueing system in the Halfin-Whitt he avy tr affic r e gime: asymptotics of the stationary distribution , Queueing Syst. 71 (2012), no. 1-2, 25–51. MR 2925789 [18] D. Gamarnik and A. Zeevi, V alidity of he avy tr affic ste ady-state appr oximation in gener alize d Jackson networks , Ann. Appl. Probab. 16 (2006), no. 1, 56–90. MR 2209336 [19] N. Gans, G. Ko ole, and A. Mandelbaum, T elephone c al l centers: T utorial, review and r ese ar ch pr osp e cts , Manu- facturing and Service Op erations Management 5 (2003), 79–141. [20] O. Garnett, A. Mandelbaum, and M. I. Reiman, Designing a c al l c enter with imp atient customers , Man ufacturing and Service Op erations Management 4 (2002), no. 3, 208–227. [21] I. Gurvich, V alidity of he avy-tr affic ste ady-state appr oximations in multiclass queueing networks: the case of queue-r atio disciplines , Math. Op er. Res. 39 (2014), no. 1, 121–162. MR 3173006 [22] I. Gy¨ ongy and N. Krylov, Existenc e of str ong solutions for Itˆ o’s sto chastic e quations via appr oximations , Probab. Theory Related Fields 105 (1996), no. 2, 143–158. MR 1392450 [23] S. Halfin and W. Whitt, He avy-tr affic limits for queues with many exp onential servers , Op er. Res. 29 (1981), no. 3, 567–588. MR 629195 [24] H. Hmedi, A. Arap ostathis, and G. Pang, On system-wide safety staffing of lar ge-sc ale p ar al lel server networks , Op erations Research (2022). [25] P . Shi, M. Chou, J. G. Dai, D. Ding, and J. Sim, Mo dels and insights for hospital inp atient op er ations: Time- dep endent e d b o ar ding time , Management Science 62 (2016), no. 1, 1–28. [26] A. L. Stolyar, Diffusion-sc ale tightness of invariant distributions of a lar ge-sc ale flexible servic e system , Adv. in Appl. Probab. 47 (2015), no. 1, 251–269. MR 3327324 [27] A. L. Stolyar, Tightness of stationary distributions of a flexible-server system in the Halfin-Whitt asymptotic r e gime , Sto ch. Syst. 5 (2015), no. 2, 239–267. MR 3442428 [28] A. L. Stolyar and E. Y udovina, Tightness of invariant distributions of a lar ge-sc ale flexible servic e system under a priority discipline , Sto ch. Syst. 2 (2012), no. 2, 381–408. MR 3354771 [29] A. L. Stoly ar and E. Y udovina, Systems with lar ge flexible server p o ols: instability of “natur al” load b alancing , Ann. Appl. Probab. 23 (2013), no. 5, 2099–2138. MR 3134731 [30] M. v an der Bo or, S. C Borst, J. S. v an Leeuw aarden, and D. Mukherjee, Scalable lo ad b alancing in networked systems: A survey of r e c ent advanc es , ArXiv e-prints 1806.05444 (2018), av ailable at 1806.05444 . [31] J. S. V an Leeuw aarden, B. W. Mathijsen, and B. Zwart, Ec onomies-of-sc ale in r esour c e sharing systems: tutorial and p artial r eview of the QED he avy-tr affic r e gime , ArXiv e-prin ts 1706.05397 (2017), av ailable at https: //arxiv.org/abs/1706.05397 . [32] W. Whitt, Understanding the efficiency of multi-server servic e systems , Management Science 38 (1992), no. 5, 708–723. [33] R. J. Williams, On dynamic sche duling of a par al lel server system with c omplete r esour c e p o oling , Fields Inst. Comm un., vol. 28, Amer. Math. So c., Providence, RI, 2000. MR 1788708 [34] H.-Q. Y e and D. D. Y ao, Diffusion limit of fair r esour c e c ontr ol—stationarity and interchange of limits , Math. Op er. Res. 41 (2016), no. 4, 1161–1207. MR 3544792 [35] H.-Q. Y e and D. D. Y ao, Justifying diffusion appr oximations for multiclass queueing networks under a moment c ondition , Ann. Appl. Probab. 28 (2018), no. 6, 3652–3697. MR 3861823