Urtzi Ayesta

Tejas Bodas
email: tejaspbodas@gmail.com.

Ina Maria Verloop

On a unifying product form framework for redundancy models

In this paper, we present a unifying analysis for redundancy systems with cancel-on-start (c.o.s.) and cancel-on-complete (c.o.c.) with exponentially distributed service requirements. With c.o.s. (c.o.c.) all redundant copies are removed as soon as one of the copies starts (completes) service. As a consequence, c.o.s. does not waste any computing resources, as opposed to c.o.c..

 where state-dependent departure rates are permitted. This allows us to show that the stationary distribution for both the c.o.c. and c.o.s. models have a product form. We give a detailed first-time analysis for c.o.s and derive a closed form expression for important metrics like mean number of jobs in the system, and probability of waiting. We also note that the c.o.s. model is equivalent to Join-Shortest-Work queue with redundancy (JSW(d)). In the latter, an incoming job is dispatched to the server with smallest workload among d randomly chosen ones. Thus, all our results apply mutatis-mutandis to JSW(d).

Comparing the performance of c.o.s. with that of c.o.c. with i.i.d copies gives the unexpected conclusion (since c.o.s. does not waste any resources) that c.o.s. is worse in terms of mean number of jobs. As part of ancillary results, we illustrate that this is primarily due to the assumption of i.i.d copies in case of c.o.c. (together with exponentially distributed requirements) and that such assumptions might lead to conclusions that are qualitatively different from that observed in practice.

Introduction

Using redundancy to minimize latency in parallel server systems has become very popular in recent years [START_REF] Ganesh Ananthanarayanan | Reining in the outliers in map-reduce clusters using mantri[END_REF][START_REF] Vulimiri | Low latency via redundancy[END_REF][START_REF] Joshi | Queues with redundancy: Latency-cost analysis[END_REF][START_REF] Nihar | When do redundant requests reduce latency?[END_REF][START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]. While there are several variants of a redundancy-based system, the general notion of redundancy is to create multiple copies of the same job that will be sent to a subset of servers. By allowing for redundant copies, the aim is to minimize the system latency by exploiting the variability in the queue lengths of the different queues. Several recent works, both empirically [START_REF] Ganesh Ananthanarayanan | Reining in the outliers in map-reduce clusters using mantri[END_REF][START_REF] Vulimiri | Low latency via redundancy[END_REF][START_REF] Ganesh Ananthanarayanan | Effective straggler mitigation: Attack of the clones[END_REF][START_REF] Dean | The tail at scale[END_REF] and theoretically [START_REF] Joshi | Queues with redundancy: Latency-cost analysis[END_REF][START_REF] Nihar | When do redundant requests reduce latency?[END_REF][START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Lee | The mds queue: Analysing the latency performance of erasure codes[END_REF][START_REF] Lee | On scheduling redundant requests with cancellation overheads[END_REF][START_REF] Gardner | A better model for job redundancy: Decoupling server slowdown and job size[END_REF], have provided indications that redundancy can help in reducing the response time of a system.

Broadly speaking, depending on when replicas are deleted, we can consider two classes of redundancy systems: cancel-on-start (c.o.s) and cancel-on-completion (c.o.c). In redundancy systems with c.o.c, once one of the copies has completed service, the other copies are deleted and the job is said to have received service. On the other hand, in redundancy systems with c.o.s, copies are deleted as soon as one copy starts being served. From a practical point of view, if servers have similar computing capacity the c.o.s system is preferable: both configurations require the same amount of signaling among servers, but the c.o.s system does not waste any computation resources. It should be noted that most of the recent literatures have focused on systems with c.o.c, and c.o.s. has remained largely elusive to exact analysis.

In a recent series of papers, Gardner et al. [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] have provided a thorough analysis of redundancy systems with c.o.c in a system with K servers each with their own queue. In the redundancy-d model of [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF], redundant copies of an arriving job are sent to d ≤ K homogeneous servers chosen uniformly at random. Under the additional assumptions that service times are exponentially distributed and that the redundant copies are independent and identically distributed (we call this the independence assumption), Gardner et al. have shown that the steady-state distribution has a product form. Further, they derive the mean response time in closed form and it is shown that the stability region, i.e., the arrival rate for which steady state exists, does not decrease in d. The latter seems counter-intuitive, because having multiple copies of the same job should imply more work for the servers and as a consequence a reduction of the stability region. The reason that this does not happen is due to the independence assumption made in [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]: the exponential assumption together with i.i.d. copies creates the situation where the effective departure rate of a job is d times faster than that in a single server. One of the objectives in this paper is to further assess the impact of this independence assumption on the performance.

In this paper, we provide, to the best of our knowledge, the first analysis on a redundancy-d system with c.o.s. for 1 ≤ d ≤ K. We assume exponentially distributed service times, but copies do not need to be i.i.d. This is because for c.o.s., whether the copies are i.i.d or not is irrelevant for the analysis and hence we do not need to make the independence assumption for c.o.s.. We adopt as a benchmark the same underlying multi-server topology as the one of [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]. This model is very convenient, everything is symmetric and homogeneous and this permits to isolate the impact of the redundancy scheme. We first show that c.o.s is equivalent to a system with a single central queue and multi-type jobs and multiple servers, as analyzed in Visschers et. al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. This allows us to conclude that the steady-state distribution of c.o.s is of product-form. In addition, we obtain an expression for the generating function for the number of waiting jobs in the system and for the mean number of jobs in the system.

We then show that redundancy-d with c.o.s is equivalent to Join-Shortest-Work queue with redundancy (JSW(d)). In the latter, an incoming job is dispatched to the server with smallest workload among d randomly chosen ones. The redundancy model hence represents a method of implementing JSW(d) without requiring to know the residual work in each of the queues. In addition, performance measures obtained for the c.o.s model carry over to JSW (d).

We then extend and adapt the model of Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] in such a way that, using a new state descriptor (different from [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]), the redundancy systems with c.o.c and i.i.d. copies can be analyzed.

We believe that our results open an important avenue for research, as it establishes the link between redundancy models and the central queue model as analyzed in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. A summary of our main contribution is given below:

1. We obtain a unifying approach to derive the stationary distribution for both c.o.s. and c.o.c. model, which is of product form (Proposition 1 and Proposition 7).

2. We provide the first exact analysis for the c.o.s model and derive the generating function for the number of waiting jobs and its mean (in Proposition 2). A key interpretation from the generating function is the fact that the number of waiting jobs is a mixture of sums of geometric random variables.

3. By allowing redundant copies that cancel-on start, we achieve the performance of JSW(d) without knowing the job sizes and the workload in the servers (Proposition 3). Thereby we also provide the first exact analysis for the performance under JSW(d) (Corollary 1). [START_REF] Joshi | Queues with redundancy: Latency-cost analysis[END_REF]. In Proposition 6, we generalize the modeling framework of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] and use this generalization to provide an alternative analysis of redundancy systems with c.o.c.

5. Numerically we observe that the impact of the modeling assumption of independent copies is non-negligible.

The rest of the paper is organized as follows. In the next section, we look at some of the related work and discuss the model and preliminaries in Section 3. We analyze the c.o.s. model in Section 4, its equivalence to JSW(d) in Section 5 and discuss some asymptotic regimes in Section 6. We give an alternate product form for the c.o.c. model in Section 7 followed by some numerical results in Section 8.

Related work

The main motivation to investigate redundancy models comes from empirical evidence suggesting that redundancy can help improve the performance in real-world applications. For example Vulimiri et al. [START_REF] Vulimiri | Low latency via redundancy[END_REF] illustrate the advantages of redundancy in a Domain Name System (DNS) query network where a host computer can query multiple DNS servers simultaneously to resolve a name. Deal et al. [START_REF] Dean | The tail at scale[END_REF] note that Google's big table services use redundancy in order to improve latency.

Most of the literature deals with the c.o.c model with exponentially distributed service times and the independence assumption. Gardner et al. [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] have provided a comprehensive analysis of redundancy c.o.c queueing models. In [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF] the authors consider a class-based model where redundant copies of an arriving job type are dispatched to a type-specific subset of servers, and show that the steady state distribution has a product form. In [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF], the previous result is applied to analyze a multi-server model with homogeneous servers where incoming jobs are dispatched to randomly selected d servers. Bonald et al. [START_REF] Bonald | Balanced fair resource sharing in computer clusters[END_REF][START_REF] Bonald | Performance of balanced fairness in resource pools: A recursive approach[END_REF] have shown that the c.o.c. model under the independence assumption and the balance fairness allocation yields the same steady-state distribution. An important result in the c.o.c. model with exponential service times and i.i.d. redundant copies is that the necessary condition for stability, i.e. that the arrival rate per server is smaller than the service rate, is also sufficient, see both [START_REF] Bonald | Balanced fair resource sharing in computer clusters[END_REF] and [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]. An important observation made in [START_REF] Bonald | Balanced fair resource sharing in computer clusters[END_REF] is that the c.o.c model is a special case of the Order Independent queue [START_REF] Ae Krzesinski | Order independent queues[END_REF], which provides a direct path to derive the steady-state distribution. However we note that the departure rates from a state under the c.o.s model fail to satisfy the order independence condition, which in turn implies the necessity of a different approach from that used in [START_REF] Bonald | Balanced fair resource sharing in computer clusters[END_REF] to analyze c.o.s. We give more details in Section 7. In a recent paper [START_REF] Adan | FCFS parallel service systems and matching models[END_REF] A more general redundancy model is the so-called (n, k, r) model proposed and studied by Joshi et al. [START_REF] Joshi | Efficient redundancy techniques for latency reduction in cloud systems[END_REF]. In this model there are n servers, r redundant copies and a job is said to be served if k copies are served. For this model, in [START_REF] Joshi | Efficient redundancy techniques for latency reduction in cloud systems[END_REF], the authors analyze the impact of the service time distribution on the performance. Shah et al. [START_REF] Nihar | When do redundant requests reduce latency?[END_REF] have showed that under i.i.d memoryless service the average response time reduces as r increases, and that r = n minimizes the average latency among all possible redundancy policies.

To the best of our knowledge, the only results available for the c.o.s model have appeared in [START_REF] Joshi | Efficient redundancy techniques for latency reduction in cloud systems[END_REF]Section 4] and [START_REF] Gardner | A better model for job redundancy: Decoupling server slowdown and job size[END_REF]Section 9]. In both cases, the analysis is restricted to the particular case in which redundancy copies are sent to all the servers.

Regarding the independence assumption, Gardner et al. [START_REF] Gardner | A better model for job redundancy: Decoupling server slowdown and job size[END_REF] propose a more accurate model in which the service time of a redundant copy is decoupled into two components, one related to the inherent job size of the task, and the other related to the server's slowdown. As written in the abstract of [START_REF] Gardner | A better model for job redundancy: Decoupling server slowdown and job size[END_REF] "The independent runtime model is unrealistic, and has led to theoretical results which can be at odds with computer systems implementation results." One of our contributions is to provide some numerical examples of this statement, and to provide an exact analysis for a redundancy model (c.o.s.) that does not need the independence assumption.

Redundancy models have also been studied in other application domains. For example an analysis of redundant models in the context of storage systems is provided in [START_REF] Lee | The mds queue: Analysing the latency performance of erasure codes[END_REF], and Ananthanarayanan et al. [START_REF] Ganesh Ananthanarayanan | Reining in the outliers in map-reduce clusters using mantri[END_REF][START_REF] Ganesh Ananthanarayanan | Effective straggler mitigation: Attack of the clones[END_REF] have studied the benefits of redundancy in straggler mitigation in data centers and cloud computing.

Model, Preliminaries and Notation

In this section, we will first describe the redundancy-d model with c.o.s. and c.o.c., give some preliminaries, as well as discuss the impact of the independence assumption by looking at the case d = K. This is followed by a brief description of the multi-type job and multi-server model of Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] which forms the basis for our analysis of the two redundancy models.

3. Under cancel-on-start (c.o.s.), once any of the copies is taken for service, the remaining copies are canceled immediately. Further, on arrival of a job, if more than one of its d feasible servers are idle, then the job is served immediately at one of the idle feasible servers based on a predefined assignment rule and the remaining copies are canceled. For this c.o.s. model, a natural assignment rule is the uniform assignment rule, that is, choose an idle feasible server uniformly at random.

Under cancel-on-complete (c.o.c.), when one of the d copies of a job completes its service (hence the name cancelon-complete), the other copies are removed. In addition, when referring to c.o.c., we assume that the service requirements of the copies of the same job are independent. We have referred to this assumption as the independence assumption.

Before giving the preliminaries, we first illustrate the impact of the independence assumption for the particular case of d = K.

Example 1. Case d = K: impact of independence assumption: When redundant copies are independent, the c.o.c. model is equivalent to a single-server queue with arrival rate λ and service rate µK. 1 In particular, this means that the mean number of jobs in the system is ρ/(1-ρ) where ρ := λ Kµ . However, if instead we assume that all copies have exactly the same service requirement, for d = K, the c.o.c. model reduces to K completely synchronized single-server queues, each having an arrival rate of λ and a service rate of µ. 2 The stability condition is now λ < µ, i.e., ρ < 1 K , and thus considerably smaller than the stability condition under the independence assumption (ρ < 1). In addition, without the independence assumption, the mean number of jobs increases significantly to ρK/(1 -Kρ).

Turning now to the c.o.s. model, when d = K, with and without the independence assumption, this model is equivalent to an M/M/K queue (see Section 9 for details). We can thus conclude that with independent copies, the c.o.c. model significantly outperforms the c.o.s. model. However, when this assumption is removed, the performance of c.o.c. degrades severely, as the stability condition is more stringent than for c.o.s. This illustrates well the unfair advantage due to the independence assumption.

Recall that redundancy models were created to improve the performance by sending copies to d servers simultaneously. However, for c.o.c. without the independence assumption, implementing full redundancy (d = K) gives a much worse performance than when there is no redundancy (d = 1). This is because d = 1 corresponds to Bernoulli routing in which case we have K independent M/M/1 queues with rates λ/K and µ. This has a better performance as compared to c.o.c. with d = K which results in K synchronized single-servers with rates λ and µ. This inefficiency comes from the fact that under c.o.c. without the independence assumption, all servers work on copies that have the exact same size. Under c.o.s., such inefficiency does not arise. This therefore shows the importance to study c.o.s. when one assumes identical copies.

Preliminaries

We define the total load as ρ = λ Kµ and assume ρ < 1. For both c.o.s. and for c.o.c., ρ < 1 is also the stability condition, that is, the condition under which a stationary distribution exists. This was proved in [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF] for c.o.c. For c.o.s. it follows as well: by symmetry, each server serves a fraction 1/K of the jobs. Since the service rate of a server is µ, each server is stable if λ/K < µ, i.e., ρ < 1.

Mean number of busy servers: Let p(i) denote the probability that there are i servers busy. Using the rate balance principle, we have (both for c.o.s. and for c.o.c.)

λ = µ K i=0 ip(i). (1)
Hence, the mean number of busy servers in the redundancy-d system is

K i=0 ip(i) = λ µ = ρK.

Probability of an idle server:

We now identify the probability of an arbitrary server (say server 1) being idle. Let P 1 denote the probability that server 1 is idle. Given i busy servers, let q i denote the probability that server 1 is idle. We have,

P 1 = K i=0 p(i)q i = K i=0 p(i) K-1 i K i
= 1 -ρ, from Eq. (1).

Central queue architecture

An important point to note is that both these redundancyd models have a central queue architecture with multiple servers. Refer to [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF]

3 4 1 3 4 1 2 (1,2) (1,3) (1,3) (2,3) (2,3) (2,4) (2,4) 1 2 3 4 (1,2) (1,2) 2 (1,2) (1,3) (2,3) (2,4) c.o.c c.o.s (1,3) (2,3) (2,4) (1,4) (1,4) (1,4)
(1,4) It is this central queue architecture for the two redundancy models, that enables us to offer a unifying framework for their analysis. This unifying analysis is based on the multitype job and server model of Visschers et al [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] which also has a central queue architecture. We use the model of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] to analyze redundancy-d with c.o.s. in detail in Section 4. In Section 7, we give an alternative analysis for redundancy-d with c.o.c. model that is based on an extension, together with an adaptation of the state descriptor, of Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. Thus, our extension of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] provides a comprehensive analysis for redundancy systems with c.o.s. and c.o.c..

Multi-type job and server model of [1]

In this section we present the multi-type job and multi-type server model and the results as presented in Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. This will form the basis in order to prove the steady-state probabilities for our redundancy models.

Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] consider a central queue system with K possibly heterogeneous servers and multi-type jobs. Let C denote the set of job types and M = {m 1 , . . . , m K } denote the set of servers respectively. Let µ i denote the service rate of server m i . For a ∈ C, S a denotes the set of feasible servers that can serve jobs of type a. Similarly, associated with a server M ∈ M is a set of job types it can serve, denoted henceforth by C(M ). The jobs of type a ∈ C arrive according to a Poisson process with rate λ a and the total arrival rate λ = a∈C λ a . All jobs have an exponential service requirement with unit mean and the service rate of server m i ∈ M is denoted by µ i for 1 ≤ i ≤ K. Each arriving job waits in a central queue and is served in a first-in first-out (FIFO) basis. Each job can only be served by one of its feasible server. If an arriving job finds that more than one of its feasible servers are idle, then an assignment rule decides which of the idle feasible server to be assigned.

State space representation

An appropriate Markovian descriptor used in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] to describe the system is of the type (n i , M i , n i-1 , M i-1 , . . . , n 1 , M 1 ) which denotes states with i busy servers (denoted by M 1 , . . . , M i ) and n j waiting jobs between servers M j and M j+1 for 1 ≤ j ≤ i -1. In this state space representation, jobs and active servers are arranged in a FIFO basis from right to left. Therefore all n j jobs have arrived before n k jobs where 1 ≤ j < k ≤ i. Note that a job is waiting in the central queue only if all of its feasible servers are busy (serving jobs that came before it). Therefore n 1 denotes the number of those jobs (who have arrived before n j jobs where j > 1) that have to wait since they have server M 1 as their only feasible server which happens to be busy. Let us denote the set of such job types by U ({M 1 }) where clearly U ({M 1 }) = {a ∈ C : S a = M 1 }. Similarly, n j represent jobs (that arrived after n j-1 but before n j+1 jobs) that have to wait because their feasible servers are busy. The feasible servers for these n j jobs must clearly be a subset of the active servers {M 1 , M 2 , . . . , M j } ahead of it (otherwise the job would not have been waiting if any of its feasible server was idle). Therefore for 1 ≤ j ≤ i, we have U ({M 1 , M 2 , . . . , M j }) := {a ∈ C : S a ⊆ {M 1 , . . . M j }} as the set of possible job types for n j waiting jobs in the representation. (See [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] for more details.)

Denote the state space of the Markov chain by S and let any generic state s ∈ S be of the type s = (n i , M i , . . . , n 1 , M 1 ). Define λ M ({M 1 , . . . , M i }) as the activation rate of server M, when the set of active servers is {M 1 , . . . , M i }. In state s = (n i , M i , . . . , n 1 , M 1 ), this is the transition rate from state s to state (0, M, s). This activation rate λ M ({M 1 , . . . , M i }) depends on the assignment rule defined for the model, which determines to which idle feasible machine (if any) a job is routed. By considering assignment rules that satisfy the following assignment condition, Visschers et al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] obtain a product form stationary distribution.

Assignment Condition. For i = 1 . . . K, and for every subset {M i , . . . , M 1 } of M of size i, the following holds:

i j=1 λ Mj ({M 1 , . . . M j-1 }) = i j=1 λ Mj ({ M1 , . . . Mj-1 }) (2)
for every permutation M1 , . . . Mi of M 1 , . . . M i .

We now recall the following theorem from [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] that provides a product form stationary distribution.

Theorem 1. For the multi-type job and server model with an assignment rule satisfying the assignment condition, the steady state probability for all states s = (n i , M i , . . . , n 1 , M 1 ) ∈ S, is given by

π(s) = α i ni . . . α 1 n1 Π λ ({M 1 , . . . , M i }) Π µ (M i , . . . , M 1 ) π(0) (3)
where

Π λ ({M 1 , . . . , M i }) = i j=1 λ Mj ({M 1 , . . . M j-1 }), Π µ (M i , . . . , M 1 ) = i j=1 µ {M1,...Mj } , α j = λ U ({M1,...,Mj }) µ {M1,...Mj } , µ {M1...Mi} = i j=1 µ Mj and λ U ({M1,M2,...,Mj }) = a∈U ({M1,M2,...,Mj }) λ a .
4 Redundancy-d with c.o.s.: An exact analysis

In this section, we provide the first exact analysis for the redundancy-d model with c.o.s.. For this model, we will say that all jobs that choose the same d servers are said to be of the same type. In all, there are K d job types and the arrival rate of any job type is

λ type = λ ( K d )
.

In
(n i , M i , n i-1 , M i-1 , . . . , n 1 , M 1 )
to represent any state s ∈ S (following [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], see Section 3.2). Recall that active servers and jobs are arranged in a FIFO basis from right to left. Further, M j in the representation denotes that server M j is currently serving a job. For a given state s, n j (j ≥ d) represents jobs that have their d feasible servers to be a subset of {M j , . . . , M 1 }. In fact, from the discussion in Section 3.2, the type of such jobs belongs to the set U ({M 1 , M 2 , . . . , M j }). An important point to note is that since each job type has d feasible servers, it can never happen that there are jobs in the central queue waiting to be served and that there are less than d busy servers in the system. Therefore, n 1 = . . . = n d-1 = 0.

Recall that to apply Theorem 1, we need to consider an assignment rule that satisfies the assignment condition. We consider the uniform assignment rule: if an arriving job finds α of its d (α ≤ d) servers busy, then the probability of choosing each idle server for service is 1 d-α . For the uniform assignment rule, we prove that it satisfies the assignment condition. The proof is in the Appendix.

Lemma 1. The uniform assignment rule satisfies the assignment condition given by Eq. (2).

Product form stationary distribution

In this section we provide the steady state distribution for the Markovian descriptor s = (n i , M i , . . . , n 1 , M 1 ), which turns out to have a product form. We obtain this by applying Theorem 1 to the c.o.s. model. (See the Appendix for the proof.) Proposition 1. The steady state distribution for any state s = (n i , M i , . . . , n 1 , M 1 ) ∈ S is given by

π(s) =      G i (K, d) π(0) i!µ i for i < d and n 1 = . . . = n d-1 = 0 ri n i . . . r d n d G i (K, d) π(0) i!µ i for i ≥ d and n 1 = . . . = n d-1 = 0 0 elsewhere (4)
where π(0) is the probability of an empty system and where

G i (K, d) = i j=1 G j (K, d) G j (K, d) = λ K d min(j-1,d-1) α=max(0,j-1+d-K) j-1 α K-j d-α-1 d -α and
r i = λ i d iµ K d .
The stationary distribution for any state s can be written in an alternative form as follows. Define

Ḡj (K, d) = G j (K, d) • K λ = d K-1 d-1 min(j-1,d-1) α=max(0,j-1+d-K) j-1 α K-j d-α-1 d -α , Ḡi (K, d) = i j=1 Ḡj (K, d) and ri = ri ρ = ( i-1 d-1 ) ( K-1 d-1 )
. The stationary distribution for any states s with n 1 = . . . = n d-1 = 0 can now be expressed as

π(s) = Ḡi (K, d) ρ i π(0) i! for i < d rni i . . . rn d d Ḡi (K, d) π(0) i! ρ (i+ i j=d nj ) for i ≥ d. (5)
Note that i + i j=d n j represents the total number of jobs in the system. From (5), one directly concludes that the stationary distribution is a function only of ρ, K and d.

It is interesting to point out that the stationary distribution does not depend on the identity of the servers that are active since the servers are assumed to be homogeneous. Hence the stationary probabilities for states with the same number of active servers (i) and same number of waiting (n i ) jobs between servers are the same.

Normalization constant and probability of busy servers

We denote by p(i) the stationary probability that the redundancy-d model with c.o.s. has i busy servers, for 1 ≤ i ≤ K. This metric would be useful in obtaining the normalizing constant π(0), which is also the probability that the system is empty, i.e., there are no busy servers (p(0) = π(0)).

Noting that

p(i) = s∈Si π(s) (6)
where S i = {s ∈ S : exactly i servers are busy}, we have the following lemma. The proof can be found in the Appendix.

Lemma 2. The probability that i servers are busy is given by

p(i) = p(i)π(0) ( 7
)
where

p(i) = K i Ḡi (K, d)ρ i , for i < d 1 1-ri . . . 1 1-r d K i Ḡi (K, d)ρ i , for i ≥ d. (8)
Further, π(0) is given by

π(0) = 1 + K i=1 p(i) -1 . (9)
A related performance metric of interest is P w (j), the probability that an arriving job sees j of its d feasible servers as busy. When a job arrives, conditioned on i busy servers, the probability that j of its servers are already busy is

( K-j i-j ) ( K i )
, i ≥ j. Therefore we have

P w (j) = K i=j p(i) K-j i-j K i . ( 10
)
In particular, the probability a job has to wait, denoted by P w , is given by P w (d).

Distribution of number of jobs

An important performance metric for the c.o.s. system is the number of waiting jobs. In this section, we obtain the expression for the probability-generating function for the number of waiting jobs in the system. An inversion of this transform gives us the distribution of the number of waiting jobs. Define p(i, m) as the probability that the system has i busy servers and m waiting jobs in the system. When i < d, there are no jobs waiting in the central queue, hence, we have p(i, m) = p(i) for m = 0 0 elsewhere.

For i ≥ d, it follows from Proposition 1 that

p(i, m) = π(0) K i Ḡi (K, d)ρ i l i (m)
where

l i (m) = {ni...n d : i l=d n l =m} r ni i . . . r n d d .
Let Q denote the random variable corresponding to the number of waiting jobs in the c.o.s. system. The probability that there are m waiting jobs (Q = m) is given by p(m) = K i=1 p(i, m). Using the above expressions for p(i, m), we can derive the distribution for the number of waiting jobs. Its proof can be found in the Appendix.

Proposition 2. The P.G.F. for the number of waiting jobs is given by

E(z Q ) = d-1 i=0 p(i) + K i=d p(i)   i j=d Geom rj (z)   ( 11
)
where Geom rj (z) = 1 -r j 1 -r j z and p(i) is given by Eq. ( 7) and [START_REF] Dean | The tail at scale[END_REF]. The expected number of waiting jobs in the system is given by

E(Q) = K i=d p(i)   i j=d r j 1 -r j   ( 12
)
and the expected number of jobs in the system is given by ). When d = K, redundancy-K c.o.s. model is equivalent to an M/M/K multiserver queue (see Section 9) whereas the JSW (K) policy is simply the traditional JSW (join shortest work) policy. But the JSW policy is known to be equivalent to the M/M/K queue [START_REF] Harchol-Balter | Performance Modeling and Design of Computer Systems: Queueing Theory in Action[END_REF]. While the equivalence is now apparent for d = 1 and d = K, we give a formal proof of this fact for any d in this section. The reasoning is based on a sample-path argument, hence, the equivalence holds for generally distributed service requirements and heterogeneous servers. Proposition 3. Assume generally distributed service requirements and heterogeneous servers. For any given sample-path realization, a given job will be served under both JSW(d) and redundancy-d with c.o.s. in the same server. As a result, the following performance measures are the same under both models:

E(N ) = E(Q) + ρK. Remark 1.
• Joint probability of servers being busy or idle.

• Delay distribution of a job.

• Total number of (waiting) jobs in the system.

Proof. For a given realization, we couple the arrivals of jobs, the d servers sampled upon arrival, and their service requirements. For redundancy-d with c.o.s., a job cannot be overtaken by copies that arrive after him in any of its d sampled servers (since it has copies in all d FIFO-queues). It therefore follows directly that upon arrival of the job, one can determine the server in which the job will be served. In order to describe this server, we define the effective workload, V red i (t) to denote the work present in server i that will be served in this server. This ignores the workload due to those copies in server i that will not be served in the server. As such, for redundancy-d with c.o.s. an arriving job will be served at the server with smallest workload V red i (t), among the d servers.

Under JSW(d), an arriving job is sent to the server with smallest workload W JSQ i (t), among the d servers. Hence, in order to prove that a given job will be served in the same server under both models we are left to show that W JSQ i (t) = V red i (t). This can be seen as follows. First note that the (effective) workload decreases at rate 1 whenever the process is positive. Now assume Now from the equivalence between JSW(d) and redundancy-d with c.o.s., we have the following corollary that to the best of our knowledge provides performance metrics for JSW(d) that have not been obtained before.

W JSQ i (t) = V red i (t)
Corollary 1. For a JSW(d) system with Poisson arrivals and exponential service requirements, the P.G.F. for the number of waiting jobs is given by Eq. [START_REF] Lee | On scheduling redundant requests with cancellation overheads[END_REF] ad the expected number of waiting jobs in JSW(d) is given by Eq. [START_REF] Gardner | A better model for job redundancy: Decoupling server slowdown and job size[END_REF]. The expected number of jobs in JSW(d) is given by E(N ) = E(Q)+ρK.

6 Asymptotic regimes for c.o.s

In this section, we consider the redundancy-d system with cancel-on-start under limiting regimes. We first consider the heavy-traffic regime where the number of servers is fixed and the traffic intensity approaches 1. In the second regime, we scale the number of servers and the arrival rate while keeping the traffic intensity unchanged, that is, the mean-field regime.

Heavy-traffic regime

For the heavy-traffic analysis, we keep the number of servers, K, fixed and let λ µ ↑ K, so that ρ ↑ 1. When d = 1, we have a system of K independent M/M/1 queues with parameters λ/K and µ. Hence, after scaling by 1 -ρ, the total number of jobs in the system converges, as ρ ↑ 1, to the sum of K exponentially distributed random variables with mean 1. In the result below, we derive the distribution of the total scaled number of waiting jobs for d ≥ 2. for any s > 0. In words, this implies that (1 -ρ)Q converges to an exponential random variable with unity parameter, as ρ ↑ 1.

Proof. From Eq. (8) we have p(i) = O(ρ i ) for i < d, while p(K) = O(ρ K /(1 -ρ)). Since ρ ↑ 1, this implies that lim ρ↑1 p(K) = 1. Now, substituting z = e -s(1-ρ) in Eq. ( 11), together with lim ρ↑1 p(K) = 1 and taking the limit ρ → 1 gives

lim ρ→1 E[e -s(1-ρ)Q ] = lim ρ→1   K j=d Geom rj (z)   .
For i < K, we have

lim ρ→1 r i = ri < 1, since d ≥ 2.
Further, lim ρ→1 r K = rK = 1. This implies that lim ρ→1 Geom rj (e -(1-ρ)s ) = 1 for j = K. For the case j = K, using L'Hopital's rule it can be seen that lim ρ→1 Geom rK (e -(1-ρ)s ) = 1 1+s . This completes the proof.

The above result shows that there is a large drop in the performance when passing from d = 1 to d = 2. Further, for d ≥ 2, the performance is independent on the parameter d.

Mean-field regime

In this subsection we consider the mean-field regime and assume that the arrival rate into the system is λ = λK where λ < µ, and let K → ∞. We note that in a recent work Hellemans and van Houdt [START_REF] Hellemans | On the powerof-d-choices with least loaded server selection[END_REF] have analyzed JSW (d) in the meanfield regime. In particular, they obtain the explicit expression for the limiting workload and response time distribution for exponential job sizes. In the following, we analyze the meanfield limit for the metric P w (j) due to its meaningfulness for c.o.s.

Define ρ = λ µ . For JSW(d), asymptotic independence in the mean-field limit has been proved in [START_REF] Bramson | Asymptotic independence of queues under randomized load balancing[END_REF]. This, together with the equivalence between JSW(d) and red-d c.o.s. and the fact that a given server under red-d c.o.s. is busy with probability ρ (see preliminaries), gives the following result for the probability P w (j) that an arriving job sees j busy feasible servers.

Proposition 5. Let λ = K λ. We have

lim K→∞ P w (j) = ρj (1 -ρ) d-j d j .
Proof. In [20, Theorem 2.1], JSW(d) is considered in the mean-field regime. It is proved there that given a finite set of servers, these servers behave independent in the mean-field limit. Note that the joint distribution of servers being idle is the same for both JSW(d) and red-d c.o.s. (Proposition 3).

In addition, we showed that the probability a server is busy equals ρ. Hence, the probability that j servers are busy out of d, follows a binomial distribution with parameters d and ρ.

7 Redundancy with c.o.c.: an alternate product from

In this section, we give an alternate analysis for the c.o.c. redundancy models of [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] by generalizing the multi-type job and multi-server model of Visschers et. al. [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF].

In the c.o.c. model, multiple feasible servers can serve a job simultaneously and hence a direct application of Theorem 1 is not possible. Additionally, the service rate of a job in the central queue of the c.o.c. model depends on the service rate received by jobs ahead of it (in the central queue), since the job will be served simultaneously in all those servers that cannot serve jobs ahead of him. Therefore, the departure rate of a job is in fact state dependent.

In Section 7.1 we first present our generalization of [1] that allows for state dependent departure rates. This generalization shall allow us to provide a Markovian state space descriptor for the c.o.c. redundancy systems, whose stationary distribution is of product form, see Subsection 7.2.

A generalized multi-type job and server model

In the existing model of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], in any state s ∈ S, an active server M j ∈ M has a fixed service rate µ Mj . Associated with the model is what we call a standard departure set function µ {•} ({•} is used to denote a set) which is defined as

µ Y = t∈Y µ t for Y ⊆ M. (13)
Therefore the total departure rate in state s = (n i , M i , . . . , n 1 , M 1 ), is given by µ {M1...Mi} = i j=1 µ Mj . In this case, one can see that the departure rate from a server is independent of the machines that are active. We will show that it is in fact possible to generalize this and allow state dependent departure rates without affecting the nature of the results in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. We say that a departure set function µ {•} is a generalized departure set function if it satisfies the following condition.

Condition 1. For 1 ≤ j ≤ K and a subset {M 1 . . . M j } of M, µ {M1...Mj } -µ {M1...Mj-1} ≥ 0. ( 14
)
For any state s = (n i , M i , . . . , n 1 , M 1 ) and a generalized departure set function µ {•} define μMj (s) := µ {M1...Mj } -µ {M1...Mj-1} . μMj (s) now denotes the state dependent departure rate from server M j where 1 ≤ j ≤ i. From the definition of a generalized departure set function, it should be clear that we are restricting to state dependent departure rates of a specific type where the departure rate from an active server depends only on the active servers ahead of it. Of course, the total departure rate in state s = (n i , M i , . . . , n 1 , M 1 ) is µ {M1...Mi} where µ {M1...Mi} = i j=0 μMj (s) (noting that for j = 1, we have µ {M1...Mj-1} = µ {φ} = 0.) It should be noted that the set function defined by Eq. ( 13) is just a special case as it satisfies the Condition 1. Also note that the generalized departure set functions could be such that in state s, μMj (s) = 0 for some 1 ≤ j ≤ i. This implies that although server M j has picked a job it can serve, the configuration of servers in front of it is such that its current service rate is zero.

We shall henceforth refer to the multi-type job and server model with a generalized departure set function satisfying Condition 1 by the name generalized multi-type job and server model. This model is thus a generalization of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] that allows a more general departure rate functions. We now have the following proposition that upholds the validity of Theorem 1 for the generalized multi-type job and server model. Proposition 6. Consider a generalized multi-type job and server model with an assignment rule satisfying the assignment condition (Eq. ( 2)). Assuming the steady state exists, the steady state probability for any state s = (n i , M i , . . . , n 1 , M 1 ) ∈ S is given by

π(s) = α i ni . . . α 1 n1 Π λ ({M 1 , . . . , M i }) Π µ (M i , . . . , M 1 ) π(0) (15)
where

Π λ ({M 1 , . . . , M i }) = i j=1 λ Mj ({M 1 , . . . M j-1 }), Π µ (M i , . . . , M 1 ) = i j=1 µ {M1,...Mj } , α j = λ U ({M1,...,Mj }) µ {M1,...Mj } .
Proof. The proof of this proposition is similar to that of Theorem 2 of [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] (Theorem 1 in this paper). An outline of the proof is given in the appendix to point out to the changes due to the generalized departure rate functions.

(s) := µ (|F j (s)| -|F j-1 (s)|) ,
where F j (s) := j l=1 S c l . The total departure rate in state s is the sum of the departure rates at each O cj for 0 ≤ j ≤ i. This corresponds to the departure rate set function for the model and is given by

µ {Oc 1 ,...Oc j } = t∈ j l=1 Sc l µ = µ|F j (s)|, ( 16
)
which satisfies Condition 1 (Eq. ( 14)).

Redundancy

π(s) = π(0) i! i j=1 jλ type µ|F j (s)| nj +1 . (17)
The redundancy-d c.o.c. system has already been extensively studied by Gardner et al [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] and hence we do not analyze this any further. The main purpose of Proposition 7 is to illustrate the unifying framework for redundancy models via Visschers et al [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] and an adaptation of its generalization. To recall the main results from [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF], the normalizing constant π(0) in Proposition 7 is given by

π(0) = K i=d 1 - i-1 d-1 K-1 d-1 ρ ,
and the mean number of customers in the system E(N ) is

E[N ] = K i=d ρ ( K-1 d-1 ) ( i-1 d-1 ) -ρ . ( 18
)
Since our state space representation for analyzing c.o.c. is different from that used in [START_REF] Gardner | Queueing with redundant requests: exact analysis[END_REF][START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF], a direct comparison between Proposition 7 and [7, Theorem 2] is cumbersome. However one can compare certain states of the system to show that the stationary probabilities obtained from Proposition 7 and [7, Theorem 2] are indeed the same. For example, consider a state with m customers, all belonging to the same type. Such a state in our representation is denoted by (m -1, O c ),

hence i = 1 and F 1 (s) = d in Proposition 7, so π((m -1, O c )) = π(0) λ type dµ m .
This is the same as the steady state probability given in [7, Theorem 2]. Similarly, consider a state with m customers all from distinct types. In our representation, such a state is represented by (O cm , . . . , O c1 ) and from Proposition 7,

π((O cm , . . . , O c1 )) = π(0) λ type dµ m .
Again, this coincides with [7, Theorem 2].

Numerical Results

In this section we present numerical results to further analyze the redundancy-d model considered in this paper. In Section 8.1 we numerically compare the performance between redundancy-d systems with c.o.s. and c.o.c. and study the impact of independence assumption made in the analysis of c.o.c.. In Section 8.2, we provide some more numerical results for redundancy-d c.o.s. This is followed by numerical results for the mean-field regime in Section 8.3.

8.1 Comparing c.o.c. and c.o.s., impact of independence assumption

In this section we compare the performance of both systems and in particular we assess the impact of the independence assumption (as made in c.o.c.) on the performance. Recall that due to cancellation on start, the analysis for co.s. with or without the independence assumption remains the same. As we saw in Example 1, the c.o.c. model with independence assumption significantly outperforms the c.o.s. model for d = K. However, when copies are no longer i.i.d., the performance of c.o.c. degrades severely. In the rest of this section, we numerically analyze the impact of the independence assumption for arbitrary values of d.

We first compare the mean number of customers (E(N )) in the system under c.o.s. and c.o.c.. The mean number of customers E(N ) in the c.o.c. system is given by Eq. ( 18), as derived in [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]. Consider c.o.c. and c.o.s. systems with parameters K = 10, ρ = 0.35 and 0.7. For varying values of d, Fig. 2 compares E(N ) for the two systems. The case d = 1 for both models is equivalent to Bernoulli routing to K servers in which case we have E(N ) = Kρ 1-ρ . The case, d = K for the c.o.c. model coincides with an M/M/1 server with arrival rate of λ and service rate of Kµ and hence E(N ) = ρ 1-ρ . On the other hand, d = K for c.o.s. corresponds to an M/M/K system with arrival rate λ and K servers each with a service rate 1. Naturally as expected, E(N ) for the M/M/K system is larger than that of the M/M/1 system associated with To elaborate on this, in Table 1, we compare E(N ) for c.o.s. and c.o.c. when increasing the redundancy parameter from d = 1 to d = 2. We fix K = 10 and µ = 1 and take varying values of ρ. Note from the discussion after Proposition 1 that the stationary distribution for the c.o.s. system is merely a function of ρ, K and d. This is also the case for the c.o.c. system (see Proposition 7) and this justifies the choice of ρ as a parameter for this table. Table 1 indicates that E(N ) for c.o.c. system is lower than that of the c.o.s. system for various values of ρ. We also quantify the percentage reduction in E(N ) when increasing the redundancy parameter from d = 1 to d = 2.

A key observation from the table is the fact that for low values of ρ, the reduction with c.o.s. is smaller than that obtained with c.o.c.. This is due to the independence assumption: at low loads, an arriving job will likely find the system empty. Under c.o.c., the job will then get an instantaneous service rate of 2µ, which explains the reduction of the order of 50% in the table. However, with c.o.s., even if the system is empty, an arriving job will only get served in one server. However, for higher values of ρ, the gain obtained with c.o.s. becomes comparable to that obtained with c.o.c.. In this case, the common feature that helps improve the performance with both systems is that by sending redundant copies, they can more efficiently use the capacity of the system.

In Figure 3 we compare the mean sojourn time. We observe that almost up to ρ = 0.7, the mean sojourn time with c.o.c. is below 1. Since the mean service time here is 1 In Figure 4 we plot the metric P w , the probability that an arriving job has to wait with K = 10, ρ = 0.7 and for various values of d. We use Equation ( 5) from Gardner et. al. [START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] to obtain the numerical values for p(i) in c.o.c. The probability an arriving job has to wait under c.o.c. can then be calculated using Equation [START_REF] Lee | The mds queue: Analysing the latency performance of erasure codes[END_REF] with j = d. With c.o.s. the probability that an arriving customer has to wait is lower than that for the c.o.c. system. The probability to wait under In Figure 5 we plot the metric p(i), that is, the probability that i servers are busy, for parameters K = 10, ρ = 0.7 and for different values of d. We observe from Figure 5 that p(K) increases in d. This is because a larger d implies a better utilization of all the servers and increases the probability of all servers being busy. This is however not the case in general for p(i), for 2 ≤ i < K. In fact, p(i) is increasing with d for lower values of i, (typically i ≤ 5) and p(i) is decreasing with d for higher values of i ( 5 < i < 10.) The effect reverses again for i = K and we have p(K) increasing in d. To explain this, let us consider the case d = 10. When i < K, p(i) corresponds to having exactly i jobs present in the system. However, p [START_REF] Lee | The mds queue: Analysing the latency performance of erasure codes[END_REF] is the probability that the number of jobs is equal or larger than 10, that is, it corresponds to the tail of the distribution.

ρ d = 1 %
A similar reason holds for other values of d.

In Figure 6, we examine the metrics P w , 1 -p(0) and p(K) for the c.o.s. system with parameters K = 10 and ρ = 0.7 and for various values of d. While P w denotes the probability that an arriving customer has to wait, 1 -p(0) denotes the probability that the system has waiting customers. Note that p(0) = A justification for this is the fact that all servers being busy implies that an arriving customer has to wait. When d = K, an arriving customer having to wait implies that all servers are busy (M/M/K system) and hence the equality.

d-1 i=0 p(i) + K l=d p(l, 0) = π(0)(1 + K i=1 f i (K, d)).

Asymptotic regimes

In this section we focus on the mean-field regime. That is, we set λ = λK, for a fixed λ, and let K → ∞. In other words, the load per server is kept constant, while taking the number of servers to infinity. In Fig. 7 we plot P w as a function of K and for different values of d. As expected, due to the asymptotic independence of the servers, P w converges to ρd where ρ = λ µ (see Proposition 5). Recall that for d = K, the models JSW , redundancyd c.o.s., and M/M/K are all the same. It is known for the M/M/K (and hence for JSW ), see for example [START_REF] Halfin | Heavy-traffic limits for queues with many exponential servers[END_REF], that in this regime the probability of waiting and the mean waiting time vanish in the limit. The same holds for the well-known JSQ dispatching policy. In recent work, van der Boor et al. [START_REF] Van Der Boor | Scalable load balancing in networked systems: Universality properties and stochastic coupling methods[END_REF] have shown that JSQ(d(K)) with d = o(K) can be asymptotically optimal in the sense of giving probability of waiting and mean waiting time equal to 0 in the limit. The latter implies that JSQ(d(K)) yields the same asymptotic performance as JSQ, but with a considerably smaller amount of overhead in terms of signaling. Getting back to our model, an interesting question is then how big d needs to be, in order for c.o.s. and JSW (d), to be asymptotically optimal with respect to the waiting probability. In Figure 8 we consider the Finally, we look at the asymptotic value of E(N )/K for redundancy-d system with c.o.s and c.o.c. This is illustrated in Figure 9. We can derive a simple approximation for the mean sojourn time under c.o.s. and hence for the mean number of jobs in the system. If an incoming job does not need to wait, its sojourn time will be 1/µ. If it has to wait, then its sojourn time will be its waiting time plus its service time 1/µ. To calculate the waiting time in the latter case, we assume that there are no other jobs (with feasible servers in common) waiting in front of the job. Hence, the waiting time will be equal to the minimum of d exponential service times, i.e. 1/(µd). Now, recalling that the probability of waiting in the limit is ρd (Proposition 5), and using Little's law, the mean number of jobs for K sufficiently large can be approximated by:

E(N ) ≈ λK (1 -P w ) 1 µ + P w ( 1 µd + 1 µ ) = ρK 1 + ρd d .
For the values of Figure 9, for d = 2 the approximation yields 0.496 whereas the value from the plots are 0.504, and for d = 4 these values are 0.4646 and 0.4553, respectively. For sufficiently large K, we expect the approximation to get more accurate as the load decreases and as parameter d increases, since then the of waiting decreases. Indeed, as K → ∞ in the mean-field regime, the above approximation coincides with the first order approximation of the corresponding mean-field result for JSW (d) given by Theorem 5.4 in [START_REF] Hellemans | On the powerof-d-choices with least loaded server selection[END_REF].

Concluding remarks

We have obtained a unifying framework to investigate the steady state distribution of c.o.s. and c.o.c. redundancy models. In the literature there are already several studies of c.o.c., but to the best of our knowledge, we are the first to give an exact analysis of c.o.s. for any 1 ≤ d ≤ K. In addition, due to the equivalence between c.o.s. and JSW (d), our approach provides a first exact analysis of workload-based load balancing schemes.

There are several interesting research problems that stem from our work. An important extension is to investigate c.o.s. with heterogeneous servers and class-based models. Our approach based on the multi-type multi-server queue will still be valid, however, the routing discipline that will satisfy the assignment condition, see Equation ( 2), needs to be determined. In this paper we have mostly focused on the queue length and related performance metrics. Obviously, the waiting time and the sojourn time distribution are very relevant metrics as well. We plan to extend the central queue approach to study the waiting time distribution of both c.o.s. and c.o.c. models. Another research path pertains with a thorough analysis of c.o.s. model in the mean-field regime. The proof for asymptotic optimality of JSWd(K) is part of future work. Last but not least, the analysis of c.o.c. without the independence assumption is a very relevant problem, and we plan to study whether this model can be analyzed through a similar approach as the one developed in this paper.

Here j-1 α is the number of ways of choosing α servers out of {M 1 , . . . M j-1 }. For a particular choice of α servers say {M 1 , . . . M α }, and given M j , the number of types that have these particular α + 1 servers as feasible servers is K-j d-α-1 . This explains the value of |C α Mj ({M 1 , . . . M j-1 })|. Note additionally that this cardinality depends only on the number of busy servers j -1 and not on the label of the servers that are actually busy. This implies that for any permutation M1 , . . . Mi of M 1 , . . . M i and for 0 ≤ j ≤ i, we have

|C α Mj ({M 1 , . . . M j-1 })| = |C α Mj ({ M1 , . . . Mj-1 })|. (19)
Here { M1 , . . . Mj-1 } is the first j servers from the permutation M1 , . . . Mi of M 1 , . . . M i for 0 ≤ j ≤ i. Now note that

λ Mj ({M 1 , . . . M j-1 }) = min(j-1,d-1) α=max(0,j-1+d-K) c∈C α M j ({M1,...Mj-1}) λ type d -α = min(j-1,d-1) α=max(0,j-1+d-K) λ type |C α Mj ({M 1 , . . . M j-1 })| d -α = min(j-1,d-1) α=max(0,j-1+d-K) λ type |C α Mj ({ M1 , . . . Mj-1 })| d -α = λ Mj ({ M1 , . . . Mj-1 })
where { M1 , . . . Mj-1 } is the first j servers from the permutation M1 , . . . Mi of M 1 , . . . M i for 0 ≤ j ≤ i. Here the range of α in the index of the summation is so chosen that the quantity j-1 α K-j d-α-1 is property defined. The first equality above follows from the fact the each type c ∈ C α Mj ({M 1 , . . . M j-1 }) activates its idle feasible server M j with rate λtype d-α and the total activation rate for M j is obtained by summing over all such types and then over all α. The third equality follows from Eq. ( 19) and the last equality in fact follows from the first. From this, it is easy to see that for

i = 1 . . . K, i j=1 λ Mj ({M 1 , . . . M j-1 }) = i j=1 λ Mj ({ M1 , . . . Mj-1 }) (20)
for every permutation M1 , . . . Mi of M 1 , . . . M i . This completes the proof.

B: Proof of Proposition 1

For a redundancy-d system with c.o.s, first note that if less than d servers are busy, then the number of waiting copies (in the traditional parallel server representation) is zero. This is because a waiting copy implies that all the other d -1 copies are also waiting for service; this is possible when at least d servers are busy in the system. Therefore for any feasible state s = (n i , M i , . . . , n 1 , M 1 ) we have n 1 = n 2 = . . . -n d-1 = 0. Since the redundancy-d system with c.o.s can be seen as a multi-type job and server model and since the uniform assignment rule satisfies the assignment condition (Lemma 1), we can apply Proposition 1. The steady state distribution for any state s = (n i , M i , . . . , n 1 , M 1 ) with i < d and n 1 = . . . = n d-1 = 0 is now given by For the redundancy-d system with c.o.s, the departure rate from any state with active servers {M 1 , . . . M j } is given by µ{M 1 , . . . M j } = jµ and hence Π µ (M i , . . . , M 1 ) = µ i i!.

Recall from the preliminaries that λ U ({M1,...,Mj }) denotes the arrival rate of types that can be served only at servers {M 1 , . . . , M j }. The number of such types is j d and the arrival rate of each such type is λ type = λ ( K d )

. Therefore C: Proof of Lemma 2

When i < d it is easy to see that

p(i) = K i i!π(s)
where s is any state of the form s = (0, M i , 0, M i-1 , . . . , 0, M 1 ) (since the servers are homogeneous, π(s) for all such states is the same). The i busy servers can be chosen in K i ways and there are i! ways to arrange the servers leading to the above expression.

Recall that r i = ρr i . Since ρ < 1, and ri = ( i-1 d-1 ) ( K-

= K i 1 1 -r i . . . 1 1 -r d G i (K, d) π(0) µ i .
Therefore we have

p(i) = K i G i (K, d) π(0) µ i for i < d K i 1
1-ri . . .

r ni i l i-1 (m -n i ) = (h i ⋆ l i-1 )(m)
where h i (j) := r j i and ⋆ represents the convolution operator. We will now obtain the expression for P i (z), the P.G.F. for p(i, m). Let H i (z) and L i (z) be the P.G.F.'s for h i (•) and l i (•) respectively. When there are i < d servers busy, there are no jobs waiting. Hence, we directly have P i (z) = p(i), for i < d.

Let

f i (K, d) = K i Ḡi (K, d)ρ i .
For i ≥ d, we have p(i, m) = π(0)f i (K, d)l i (m). Further note when i ≥ d, we have

H i (z) = 1 1 -r i z L i+1 (z) = H i+1 (z)L i (z
where π(0) is given by Eq. ( 9). Now define Geom rj (z) := 1 -r j 1 -r j z .

Geom rj (z) denotes the P.G.F. of a geometric random variable X j with parameter r j . (X j has the distribution P (X j = n) = (1 -r j )r n j .). Together with Lemma 2, we then have for i ≥ d (31)

  , Adan et.al. show the equivalence of the c.o.c. model with two other parallel service models.

Figure 1 :

 1 Figure 1: A representation of redundancy-2 c.o.s. and c.o.c. system and their associated central queue. Jobs are denoted by their set of feasible servers and such jobs are arranged in the central queue in a FIFO manner. The choice of server in c.o.s. is by the uniform assignment rule.

5

  Equivalence of JSW(d) and redundancy-d with c.o.s. In this section, we will show the equivalence of the redundancy-d c.o.s. model with that of the join the shortest work among d servers (JSW(d)) policy. This is trivial to note for the case when d = 1 and d = K. When d = 1, the redundancy-1 c.o.s. model as well as the JSW(1) model are equivalent to probabilistic routing to K queues when the routing probabilities are Bernoulli ( 1 K

Proposition 4 .

 4 Assume d ≥ 2. Then lim ρ↑1 p(K) = 1 and lim ρ↑1 E[e -s(1-ρ)Q ] = 1 1 + s ,

Figure 2 :

 2 Figure 2: E(N ) for c.o.c. and c.o.s. for ρ = 0.35, 0.7 and K = 10.

Figure 3 :

 3 Figure 3: E(T ) for c.o.c. and c.o.s. with K = 10.

Figure 4 :

 4 Figure 4: P w for c.o.c. and c.o.s. with K = 10, ρ = 0.7.

Fig. 6 Figure 5 :

 65 Figure 5: p(i) for different d when K = 10, ρ = 0.7.

Figure 6 :

 6 Figure 6: P w , 1 -p(0) and p(K) for K = 10 and ρ = 0.7.

Figure 7 :Figure 8 :

 78 Figure 7: For λ = 0.5 and K → ∞, for a fixed d, we have lim K→∞ P w = ρd .

Figure 9 :

 9 Figure 9: Asymptotics for E(N )/K for c.o.s and c.o.c. system with λ = 0.45 and various values of d.

λ

  π(s) = Π λ ({M 1 , . . . , M i }) Π µ (M i , . . . , M 1 ) π(0)and for i ≥ d and n 1 = . . . = n d-1 = 0 byπ(s) = α i ni . . . α d n d Π λ ({M 1 , . . . , M i }) Π µ (M i , . . . , M 1 ) π(0).HereΠ λ ({M 1 , . . . , M i }) = i j=1 Mj ({M 1 , . . . M j-1 }) Π µ (M i , . . . , M 1 ) = i j=1 µ {M1,...Mj } α j = λ U ({M1,...,Mj }) µ {M1,...Mj } for 1 ≤ j ≤ i.

λ.λd

  U ({M1,...,Mj }) = λ Now using the fact thatλ Mj ({M 1 , . . . M j-1 }) λ ({M 1 , . . . , M i }) = i j=1 Mj ({M 1 , . . . M j-1 }) = i j=1 λ type min(j-1,d-1) α=max(0,j-1+d-K)Using these definitions, the statement of the theorem follows. This completes the proof.

1 1 - 1 .D: Proof of Proposition 2

 112 r d G i (K, d) d)ρ i and hence p(i) given in the statement now follows. Using the fact that K i=0 p(i) = 1 and p(i) = p(i)π(0), one can see that π(0) = 1 + K i=1 p(i) -We first analyze l i (m) and try to obtain a recursive form for it. First note that for i = d, we have l d (m) = r m d for m ≥ 0. Now for i > d, l i (m) = {ni...n d :

  ) andP i (z) = π(0)f i (K, d)L i (z), and since l d (m) = r m d , we have L d (z) = H d (z). We therefore have for i ≥ d, P i (z) = π(0)f i (K, d)

P

  i (z) = π(0)f i (K, d) redundancy-d c.o.c model, we have λ Oc j ({O c1 , . . . O cj-1 }) = λ cj = λ type . Further λ U ({Oc 1 ,...,Oc j }) = j l=1 λ c l = jλ type . Therefore we have, Π λ ({O c1 , . . . , O ci }) = i j=1 λ type , α j = jλ type µ|F j (s)| .After rearranging the terms and bit of algebra, the steady state distribution for state s = (n i , O ci , . . . , n 1 , O c1 ) is given by

  From the P.G.F above, one can conclude that the random variable Q is a mixture of sum of geometric random variables. This is also the case for redundancy-1 with c.o.c..The redundancy-K with c.o.s. is equivalent to an M/M/K system. This is easy to see from the central queue architecture; since each job can be served at any of the K servers, a server that becomes free ends up choosing the head of the queue customer for service. We verify in Appendix E that our expressions for c.o.s. for these special cases indeed coincide with the known results for these cases.

Remark 2. Special cases: We now discuss two special cases for the redundancy-d model with c.o.s. namely d = 1 and d = K. The case d = 1 or redundancy-1 c.o.s. corresponds to the system with Bernoulli routing of jobs to the K servers.

  is true at time t, and a new job arrives. (This is trivial to note for t = 0.) It is served in the server (among the d servers) with smallest value for W JSQ With the above proposition, redundancy-d with c.o.s. can be perceived as a method of implementing JSW(d) without requiring the knowledge of the residual work in each of the queue. But such an implementation of JSW(d) would only be possible if the underlying system allows for the possibility to use extra copies of the same job. Indeed, the power of allowing even a single extra copy per job may prove to be very beneficial in the absence of workload information at the queues.

	i respectively). Since W JSQ (t) (V red (t)), under JSW(d) (redundancy-d, i i (t) = V red i (t), it is the same server
	say server m, that has the smallest value in the two system
	and hence, the service requirement of this job is added both
	to W JSQ m	(t + ) and to V red m (t + ).
	Remark 3.

  7.[START_REF] Ganesh Ananthanarayanan | Reining in the outliers in map-reduce clusters using mantri[END_REF] An alternate product form for c.o.c. Proposition 6, we have generalized the results of[START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] by considering generalized departure set functions. This generalization can now be exploited to apply Proposition 6 to redundancy-d with c.o.c..For the redundancy-d c.o.c. model, one can view a type for a job as the choice of d servers that are chosen (see preliminaries). Since this choice is uniform over the set of all servers, each type is equally likely. Recall that λ type denotes the arrival rate of each type where λ type = λ where C denotes the set of all job types. For any type c ∈ C, let S c denote the set of d feasible servers for that type. For each type c ∈ C, we associate a label O c . This label will be used to identify the least recent (oldest) type c job present in the central queue.

	In ( K d )	. Also recall
	some earlier notation

A new state space representation for c.o.c.

We propose a new state space representation (n i , O ci , . . . , n 1 , O c1 ) to analyze the redundancy-d model of

[START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF]

. In such a state, the job at the head of line (of the central queue) is of type c 1 . Since the central queue has FIFO ordering from right to left, this is also the oldest type c 1 job in state s and is therefore indicated by the label O c1 . n 1 denotes the number of type-c 1 jobs that arrived after O c1 . These jobs were followed by a type-c 2 job, represented by O c2 . Now n 2 denotes the subsequent arrivals that are of either type c 1 or c 2 . In general, in state s and for 0 ≤ j ≤ i, O cj indicates the position in the central queue of the oldest type-c j job and n j denotes waiting jobs of type c ∈ j i=1 c i . Now in state s, the jobs that can receive service are the ones that are represented by the labels (O c ). This is because of the FIFO scheduling in the central queue and the fact that the waiting jobs (n j ) are not the oldest jobs in their type. In state s = (n i , O ci , . . . , n 1 , O c1 ), the departure rate of job O cj equals µ times the number of its feasible machines that are not used by the jobs O cj-1 , . . . , O c1 . Hence, it is given by μcj

  -d c.o.c. as generalized multi-type job and server model With a state space representation of s = (n i , O ci , . . . , n 1 , O c1 ), the c.o.c. model can now be viewed as a generalized multitype job and server model where O cj resembles a super-server of state dependent service rate µ (|F j(s)| -|F j-1 (s)|). Another artifact of this representation is that each super-server O cj can serve jobs only of type c j and therefore each type has a dedicated super-server that will only serve jobs of its type. From this interpretation, it should be clear that there is no need to specify any assignment rule because an arriving job is never faced with the prospect of being served by multiple super-servers (although he will be served by multiple servers that define the super-server). Due to this, we have λ Oc ({O c1 , . . . , O cj }) = λ c and hence the assignment condition is trivially satisfied. Proposition 6 can now be used to provide a new product form distribution for this redundancy model. This distribution is different from that of[START_REF] Gardner | Redundancy-d: The power of d choices for redundancy[END_REF] due to a different state space representation. At the end of this section we comment on the latter. We have the following proposition (proof in appendix).

Proposition 7. For the redundancy-d c.o.c. model, the steady state distribution for state s = (n i , O ci , . . . , n 1 , O c1 ) is given by

Table 1 :

 1 decrease in % decrease in E(N ) E(N ) for c.o.s. E(N ) for c.o.c. Reduction in the mean number of jobs with c.o.c. and c.o.s. for K = 10, when increasing the redundancy parameter from d = 1 to d = 2.

		with d = 2	with d = 2
	0.1 1.11	9.5	51.56
	0.2 2.5	18.17	53.3
	0.3 4.28	26.25	55.26
	0.4 6.66	33.94	57.5
	0.5 10	41.44	60.09
	0.6 15	48.95	63.17
	0.7 23.33 56.75	66.95
	0.8 40	65.28	71.79
	0.9 90	75.44	78.56

This can be seen as follows: starting from an empty queue, a first job will have copies in all K servers. A copy will hence finish at rate µK, after which a new job starts being served in all servers.

A job will be served simultaneously in all K servers. Since the service requirement is identical, all copies finish after an exponential time with mean µ.

Research partially supported by the French "Agence Nationale de la Recherche (ANR)" through the project ANR-15-CE25-0004 (ANR JCJC RACON)

Appendix A: Proof of Lemma 1

Assuming the uniform assignment rule, we will prove that for every subset {M i , . . . , M 1 } of M of size i and for every j (1 ≤ j ≤ i) we have λ Mj ({M 1 , . . . M j-1 }) = λ Mj ({ M1 , . . . Mj-1 }) for every permutation M1 , . . . Mi of M 1 , . . . M i . From Eq. ( 14), this will imply the assignment condition.

To characterize λ Mj ({M 1 , . . . , M j-1 }) (the activation rate of server M j , when the set of active servers is {M 1 , . . . , M j-1 }), we will require the number of types that can activate server M j when servers M 1 to M j-1 are busy. This can be obtained by first counting the number of types that can activate server M j while having exactly α of its d feasible servers already busy (and hence among M 1 to M j-1 ) and then summing over α. Now let C(m i ) denote the set of types which have m i as a feasible server. Likewise, C(M i ) would denote the same for a generic server M i . Let S c denote the set of d feasible servers for type c jobs. Now when the system has exactly j -1 busy servers denoted by {M 1 , . . . , M j-1 } and when M j is an idle server, define

as the set of types that have the idle server M j as feasible and have exactly α out of their d feasible servers in the set {M 1 , . . . , M j-1 }. Hence for any type c ∈ C α Mj ({M 1 , . . . , M j-1 }), there are exactly d -α idle servers (which include server M j ) and due to the uniform assignment rule, an arriving customer of this type will activate server M j with probability 1 d-α . Now it can be seen that

The P.G.F. for the number of waiting jobs is then given by

The required P.G.F. in the statement of the theorem now follows from Eq. [START_REF] Van Der Boor | Scalable load balancing in networked systems: Universality properties and stochastic coupling methods[END_REF].

The expected number of waiting jobs is given by

. Note that dPi(z)

(The second equality is because

The expected number of waiting jobs in the system (E(Q)) is hence given by

where p(i) is given by Eq. ( 7) and [START_REF] Dean | The tail at scale[END_REF]. Since the mean number of jobs at the servers is same as the mean number of busy servers we have E(N ) = E(Q) + Kρ.

E: Special cases; d = 1 and d = K

In this appendix, we will discuss the special cases for the c.o.s. model corresponding to d = 1 and d = K.

When d = 1, both c.o.s. and c.o.c. system are equivalent to a system with Bernoulli routing to the K servers. The queue length process of each server is independent and each server represents an M/M/1 queue with arrival rate λ K and service rate µ. It is easy to check from our analysis of Section 4 that when d = 1, we have

and r i = ρ. Using the notation that ρ = λ Kµ (the load per M/M/1 queue), the steady-state probability from Eq (4) for any state s = (n i , M i , . . . , n 1 , M 1 ) for i ≥ 0 is given by

p(i) is given by p(i)π(0) where from Eq. (8) we have p(i) = 9), the normalizing constant is obtained as

Since π(0) is also the probability that the system of K independent M/M/1 queues is empty, this justifies π(0) = (1 -ρ) K . Using this, p(i) can be simplified to p(i)

. Now using Proposition 2, it is easy to see that the P.G.F. for the number of waiting jobs in the redundancy-1 system is given by

Further, the mean number of waiting jobs and the total number of jobs in the system is given by

We now analyze the redundancy-K system with c.o.s.. The redundancy-K system is same as an M/M/K system with arrival rate of λ and K servers each with rate µ. For this case, substituting d = K in our analysis of Section 4 gives us

and r K = ρ. Again, by simple substitution, the probability of i busy servers p(i) is given by p(i) = π(0)p(i) where

The normalizing constant π(0) is given by

which coincides with the normalizing constant of an M/M/K system (see page 260, [START_REF] Harchol-Balter | Performance Modeling and Design of Computer Systems: Queueing Theory in Action[END_REF]). Note that p(K), the probability that all servers are busy corresponds to the blocking probability of the M/M/K system. In fact, we obtain

which as expected coincides with the celebrated Erlang-C formula. Based on substituting d = K in Proposition 2, we see that the P.G.F. for the number of waiting jobs in the redundancy-K system is given by

Further, the mean number of waiting jobs and total number of jobs in the system is given by E

These expressions verify with the corresponding quantities for the M/M/K system (see [START_REF] Harchol-Balter | Performance Modeling and Design of Computer Systems: Queueing Theory in Action[END_REF]).

F: Proof outline of Proposition 6

The proof outline is as follows. We first state the transition rates and the global balance equation for the generalized multi-type job and server model. Then, we identify the partial balance equations (as in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]) that need to be satisfied by our guess for the stationary distribution given in the statement of the theorem. Using the generalized departure set functions in the equations, we show that the stationary distribution from the theorem indeed satisfies the partial balance resulting in a product-form.

Transition rates and global balance

The transitions from an arbitrary state are of the following type.

Arrivals: When in state s = (n i , M i , . . . , n 1 , M 1 ), the total arrival rate is λ. An arriving job either finds no feasible server to be idle or is picked for service by a feasible idle server. In case of multiple feasible idle server, the assignment rule will determine which of the feasible and idle server will serve this arrival. It is asumed in the statement of the theorem that the assignment rule must satisfy the assignment condition. The arrival rate of jobs that find no feasible server to be idle is given by λ U ({M1,...,Mi}) . The arrival rate of jobs that activate a feasible server (according to the assignment rule) is λλ U ({M1,...,Mi}) respectively. As in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], we have

(23) Departures: On departure of a job in service, either the server serving this job become idle or the server serving this job picks some waiting job. These two cases are defined by the following two events that describe transition rates associated with a departure.

• { Rate P } : the rate of transition into state s due to a departure such that the server serving the departing job becomes idle.

• { Rate Q } : the rate of transition into state s due to a departure such that the server serving the departing job picks a waiting job.

The global balance equations for any state s = (n i , M i , . . . , n 1 , M 1 ) are as follows. When n i > 0, we have

When n i = 0, we have

We will now characterize Rate P and Q in detail. As in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], we denote the state (n i , M i , . . . , l, M, n k -l, M k , . . . , n

Let p M kl (s) (resp. q Mj kl (s)) denote the probability of transition from state insert M kl (s) (resp. swap Mj kl (s)) to state s due to departure of server M (resp. M j ). Therefore p M kl (s) is the probability that after departure of server M, the state insert M kl (s) is such that there are no waiting jobs in the system that are feasible to server M. Similarly, q Mj kl (s) is the probability that after departure of server M j , the state swap Mj kl (s) is such that after departure of server M j , the next feasible job for this server is the n j-1 + 1 st job after server M j-1 . Therefore we have,

The rates P and Q can be quantified as below:

where

Note the subtle difference between our definition of P M (s), Q Mj (s) given above and that in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF]. Unlike [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], we have subsumed the departure rates inside the definition of P M (s), Q Mj (s).

Partial Balance: Now as in the proof of Theorem 2, [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], the aim is to prove that the expression for π(s) in the statement of the theorem satisfies the following 4 partial balance equation. For n i > 0 π(s)µ {M1,...Mi} = λ U ({M1,...,Mi}) π(n i -1, M i , . . . , n 1 , M 1 ) (26) and for n i = 0,

Here we have made use of Eq (23) to replace λ in the global balance equations ( 24) and ( 25). Firstly, it is easy to see that Eq. ( 3) satisfies Eq. ( 26) and ( 27). Hence to complete the proof, it is sufficient to prove that Eq. ( 3) satisfies Eq. ( 28) and Eq. ( 29). Let us first prove that Eq. ( 3) satisfies Eq. ( 29). We proceed exactly as in the proof of Theorem 2 in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF] and divide each term in the summation of P M (s) by π(s)λ M ({M 1 , . . . , M i }). In the corresponding part of the proof of Theorem 2 in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], we have µ M instead of μM (insert M kl (s)). Proceeding in the same way and using the fact that μM (insert M kl (s)) = µ {M1,...M k ,M } -µ {M1,...M k } it can be shown that μM (insert M kl (s))p M kl (s)

A similar procedure for the boundary state gives us μM (s, 0, M )p M 1n1 (s)

Now as in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], this implies that

To prove, that Eq. (3) satisfies Eq. ( 28), as in the proof of Theorem 2 in [START_REF] Visschers | A product form solution to a system with multi-type jobs and multi-type servers[END_REF], divide each term in the summation of Summing over all j we have π(s)λ U ({M1,...,Mi}) = i j=1

which is same as Eq. ( 28). This completes the proof outline.

G: Proof of Proposition 7

From Section 3, recall that the steady state exists for any ρ <


Mamadou Mboup Projet

Cédric Join
email: cedric.join@cran.uhp-nancy.fr

Michel Fliess
email: michel.fliess@polytechnique.edu

A DELAY ESTIMATION APPROACH TO CHANGE-POINT DETECTION

Keywords: Jump parameter systems, Signal detection, Delay estimation, Symbol manipulation, Polynomial approximation

The change-point detection problem is cast into a delay estimation. Using a local piecewise polynomial representation and some elementary algebraic manipulations, we give an explicit characterization of a change-point as a solution of a given polynomial equation. A key feature of this polynomial equation is its coefficients being composed by short time window iterated integrals of the noisy signal. The so designed change-point detector shows good robustness to various type of noises.

INTRODUCTION

Given a piecewise regular signal, the purpose of this paper is to detect its discontinuities and estimate their locations. The problem is challenging especially for applications requiring on-line detection: the difficulties are stemming from corrupting noises which are blurring the discontinuities, and the combined need of fast calculations for real-time implementation and of reliable detection. A large amount of literature is devoted to these questions. We refer to the monograph of M. Basseville and V. Nikiforov [START_REF] Basseville | Detection of Abrupt Changes: Theory and Application[END_REF] for a nice tutorial presentation and many application examples. When a statistical model description of the observation data is available, then it is very common to approach the problem via hypothesis testing. And, the classical cumulative sum statistic [START_REF] Basseville | Detection of Abrupt Changes: Theory and Application[END_REF] (see [START_REF] Vaswani | Additive change detection in nonlinear systems with unknown change parameters[END_REF] for other statistics) is essential in most of the algorithms devised within this framework. Furthermore, since the wavelet transform coefficients of a signal inherently reveal the presence of irregularities, the combination wavelettest statistic has now became classical for change point detection [START_REF] Raimondo | A peaks over threshhold model for change-point detection by wavelets[END_REF]. Model selection is another common approach and therein, the Bayesian theory plays a important rôle. An example using penalized contrast may be found in [START_REF] Lavielle | Using penalized contrasts for change-point problem[END_REF] (see also [START_REF] Fan | Bayesian estimation of the number of change points in simple linear regression models[END_REF]). As opposed to the test statistic based methods which are sequential, the ones developed in this framework are rather global: the change points are simultaneously estimated. Clearly, this is not suited for on-line implementation. We may also mention other approaches such as the kernel based ones [START_REF] Gijbels | On the estimation of jump points in smooth curves[END_REF], [START_REF] Désobry | An online kernel change detection algorithm[END_REF]. These approaches are mainly local. The solution we are presenting here is both sequential and local. It is based on a direct estimation of the points of singularity of the signal's derivatives. Following this spirit, one immediate solution would be to use a signal derivative estimator, as those proposed in [START_REF] Mboup | A revised look at numerical differentiation with an application to nonlinear feedback control[END_REF], to locate the singularities. But, as we will shortly see, an explicit estimation of the derivative is not necessary. Using a local piecewise polynomial representation of the signal, we are able, in section 2, to cast the changepoint problem into a delay estimation. The proposed detection algorithm is next presented in section 3. It's principle is to express a change-point as a solution of a polynomial equation, the coefficients of which are composed by short time window iterated integrals of the noisy signal. These integrals lowpass filter the noise as witnessed (testified) by the extensive simulation results presented in section 4. Different types of singularities are considered therein1 : change-point in the mean, in the local slope, and a jump on the signal's second order derivative.

PROBLEM STATEMENT

Signal description

Consider the following piecewise smooth signal model,

x(t) = K i=1 H(t -t i-1 ) f i (t -t i-1 ), (1)
where H(•) is the Heaviside function and where each f i (t) is a smooth segment. The only irregularities of x, commonly called change-points, occur at times t i , i = 1, . . . , K, where K is unknow. We set t 0 = 0. Based on the observation y(t) = x(t)+n(t), where n(t) is an additive noise corruption, we want to detect the change-points and estimate their locations t i . Let T be given and assume that there is at most one discontinuity point in each interval I T τ = (τ -T, τ ), τ T . In the sequel, we will set

x τ (t) = x(t + τ -T ), t ∈ [0, T ], τ T,
for the restriction of the signal in I T τ and we redefine the discontinuity point, say t τ , relatively to I T τ with: t τ = 0 if x τ (t) is smooth and 0 < t τ T otherwise. The problem now reduces to detect the nonsmoothness of x τ (t), for each τ T . For this task, we consider next a polynomial model for x τ (t).

Local polynomial model

By assumption, each sliding interval I T τ contains at most one change point. The corresponding signal x τ (t) then admits a representation of the form

x τ (t) = [1-H(t-t τ )]a(t)+H(t-t τ )b(t-t τ ) t ∈ [0, T ],
(2) where a(t) and b(t) are smooth and

x τ (t) = b(t) if t τ = 0.
Here the possible change-point t τ is viewed as a delay. Our detection method relies on the direct estimation of this delay. To proceed we need a model for a(t) and b(t) above. Now, a polynomial model is found adequate especially when the length of the interval I T τ , viz. T , is small. So henceforth, a(t) = p i=0 a i t i and b(t) = q j=0 b j t j are two polynomials of degree p and q respectively. Let us express Eq. ( 2) in the operational domain2 :

xτ (s) = p i=0 i!a i s i+1 + e -tτ s q j=0 j!b j s j+1 (3)
By using elementary differential algebraic operations on this expression and translating the result back in the time domain, one can show that t τ is identifiable with respect to x(t). These operations, which lead to the design of an estimator of t τ based on the noisy observation y(t), are exposed below.

DETECTOR SYNTHESIS AND IMPLEMENTATION

Change-point estimation

We start considering the simplest model for x τ (t), viz. a piecewise constant model. We thus have: a(t) ≡ a 0 and b(t) ≡ b 0 in (2). Eq. ( 3) then reduces to

xτ (s) = a 0 s + e -tτ s b 0 s . (4)
In the sequel, we will write xτ instead of xτ (s) to ease the notations. Now, we are going to show that t τ may be expressed as a function of x τ (t) only: we will say that t τ is identifiable with respect to xτ (see e.g. [START_REF] Fliess | An algebraic framework for linear identification[END_REF] for various types of identifiability). This will stem from the following manipulations. The first step is to eliminate the unknown coefficients a 0 and b 0 . For this, multiply both side of (4) by s: sx τa 0 = e -tτ s b 0 .

Noting that e -tτ s satisfies d ds + t τ e -tτ s = 0, we have

d ds + t τ (sx τ -a 0 ) = 0.
Differentiation with respect to s will eliminate a 0 . The delay t τ is therefore explicitly given by:

s d 2 ds 2 xτ + 2 d ds xτ + t τ (s d ds xτ + xτ ) = 0. (5)
Recall that, by the very classical rules of operational calculus, derivation with respect to s is the operational analogue of multiplication by -t: d ds xτ ¡-tx τ (t). And multiplication by s corresponds to time derivation: sx τ ¡ d dt x τ (t). While the first operation is easy to implement, numerical differentiation is known to be difficult and ill-conditioned. Let ν be a positive integer greater than the highest power of s in (5) and divide both members of (5) by a s ν . Then only negative powers of s will intervenes in the resulting equation. Replacing the unobserved signal x τ by its noisy observation counterpart y τ , we obtain the linear estimator t τ of t τ :

ŷ′′ τ s ν-1 + 2 ŷ′ τ s ν + t τ ŷ′ τ s ν-1 + ŷτ s ν = 0, (6)
where ŷ′ τ stands for d ds ŷτ . Recall that division by s k , k > 0 corresponds to k th -iterated time integration:

(k -1)!s -k x(s) ¡ t 0 (t -α) k-1 x(α)dα.
The estimate of t τ thus follows upon expressing (6) back in time domain, for ν 2:

T 0 (νt -T )(T -t) ν-2 y τ (t)dt t τ = T 0 {(ν + 1)t + 2T }(T -t) ν-2 ty τ (t)dt (7)
If the interval I T τ is devoid of a change point, then the right hand side and the bracketed term above must be zero up to the (output) noise level. In this situation, we will have t τ = 0. Note that noises are viewed highly as highly fluctuating phenomena (see [START_REF] Fliess | Analyse non standard du bruit[END_REF] for the mathematical details). They are attenuated by the iterated time integrals, which are simple examples of low-pass filters.

Synthesis of change point detector

The generalization of the preceding developments to higher order polynomial models is straightforward. To see this let us multiply both sides of (3) by s ℓ where ℓ △ = max{p, q} + 1. Then we may write: s ℓ xτ = A(s) + e -tτ s B(s) where A(s) and B(s) are two polynomials in s, with degree at most ℓ -1. Let us isolate B(s) as in B(s) = e tτ s (s ℓ xτ -A(s)) and annihilate it by d ℓ ds ℓ . We obtain:

ℓ i=0 ℓ i d ℓ-i ds ℓ-i (s ℓ xτ -A(s)) t i τ = 0.
Now, the unknown polynomial A(s) has to be eliminated. This is achieved by applying again d ℓ ds ℓ to the above expression. The result is next divided by s ν , ν > ℓ, to avoid time differentiation. Replacing xτ by its corresponding noisy observation ŷτ and returning to the time domain, we finally obtain, for each time instant τ , a polynomial expression D {p}{q} (τ, t τ ) of the form:

D {p}{q} (τ, t τ ) = T 0 y(τ -(T -t)) ℓ i=0 P i (t) t i τ dt. ( 8
)
Each P i (t) in ( 8) is a polynomial depending on the parameters p, q, ν and T and is defined as:

P i (t) = ℓ j=0 ℓ i 2ℓ -i j (-1) 2ℓ-i-j ℓ! (ℓ -j)!(j + ν -ℓ -1)! (T -t) j+ν-ℓ-1 (τ -(T -t)) 2ℓ-i-j
Clearly, this expression does not provide an estimator for t τ . However, it does provide us with a change-point detector, the implementation of which is described below.

Implementation

To see how to implement D {p}{q} (τ, t τ ) as a detector, let us set t τ to a fixed value, say t τ = rT for some r ∈ (0, 1). Consider the following scenarios:

• Case 1: I T τ does not contain a change-point. Then, the equality D {p}{q} (τ, rT ) = 0 holds except for the local mismodelling error and the noise contributions. This situation is depicted in figure 1-(a).

• Case 2: I T τ contains a change-point, located at t ⋆ τ = rt. In this case, D {p}{q} (τ, rT ) largely deviates from zero. The selected value for t τ does not correspond to a root of the polynomial in [START_REF] Mboup | A revised look at numerical differentiation with an application to nonlinear feedback control[END_REF], as shown in figure 1-(b).

• Case 3: I T τ contains a change-point at t ⋆ τ which coincides with rT . Then D {p}{q} (τ, rT ) must vanishes up to the noise and the local mismodelling errors. Figure 1-(c) illustrates this situation.

Remark 1 Choosing t τ = 0.5T seems to be the best choice especially for the detection of very close change points. According to these scenarios, a change-point is meant by the passage from case 2 to case 3, i.e. when starting beyong a threshold, |D {p}{q} (τ, t τ )| falls down and crosses zero. Since the characteristics of the noise are unknown, the selection of the threshold will subsequently be based on some heuristics. Also, the performance of the proposed detection method are investigated in the next section, through Monte Carlo simulation.

SIMULATION RESULTS

The change-point detector in ( 8) is now numerically implemented using the trapezoidal method, under various signal and noise settings. Very simple local polynomial models are used for all signal settings. Robustness to noise corruption is highlighted using several noise models with different powers: normal, uniform and Perlin noise 3 . Table 1, where the notation M {p}{q} refers to (3) with degrees p and q, summarizes the simulation results. For all simulations, the sampling period (without unit) is T e = 1 and the detector is inhibited during 10T e (10 samples) after each detection. For each setting, the distribution of the estimated change-points, computed over 100 Monte Carlo simulations, are shown below. We start with a piecewise constant signal. In view of the preceding experiments, the method presents a good robustness to noise. It also appears that the noise characterisctics do not have a significant effect on the detection. We consider now a piecewise non polynomial signal with 3 segments, through an additive white Gaussian noise (see fig-

Fig. 1 .

 1 Fig. 1. Detector

Fig. 2 .

 2 Fig. 2. Piecewise constant signal: changes in the mean

Fig. 3 .

 3 Fig. 3. Piecewise polynomial signal: order 1 and 2.

ure 4 -Fig. 4 .

 44 Fig. 4. Non polynomial signal

Table 1 .

 1 Simulation summarized

	The origi-

The first author would like to acknowledge M. Djafari (LSS, Gif-sur-Yvette, France) who kindly provided him with the test signals.

Algebraic manipulations of delays are easier via the formalism of operational calculus. See[START_REF] Fliess | Algebraic change point detection[END_REF] for more details.

Perlin's noises[START_REF] Perlin | An image synthesizer[END_REF] are fractal like and they are quite popular in computer graphics.