Adam Wolfe Gordon

Paul Lu
email: paullu@cs.ualberta.ca

Elastic Phoenix: Malleable MapReduce for Shared-Memory Systems

Keywords: MapReduce, parallel programming, malleable, shared-memory systems

We present the design, implementation, and an evaluation of Elastic Phoenix. Based on the original Phoenix from Stanford, Elastic Phoenix is also a MapReduce implementation for shared-memory systems. The key new feature of Elastic Phoenix is that it supports malleable jobs: the ability add and remove worker processes during the execution of a job. With the original Phoenix, the number of processors to be used is fixed at start-up time. With Elastic Phoenix, if more resources become available (as they might on an elastic cloud computing system), they can be dynamically added to an existing job. If resources are reclaimed, they can also be removed from an existing job. The concept of malleable jobs is well known in job scheduling research, but an implementation of a malleable programming system like Elastic Phoenix is less common. We show how dynamically increasing the resources available to an Elastic Phoenix workload as it runs can reduce response time by 29% compared to a statically resourced workload. We detail the changes to the Phoenix application programming interface (API) made to support the new capability, and discuss the implementation changes to the Phoenix code base. We show that any additional run-time overheads introduced by Elastic Phoenix can be offset by the benefits of dynamically adding processor resources.

Introduction

In many programming systems, the number of processes, threads, or workers is fixed at start-up time. After all, creating and initializing per-thread data structures is most easily handled at the beginning of the job. Furthermore, the classic batch schedulers (e.g., Load Sharing Facility (LSF), Portable Batch System (PBS), Sun Grid Engine (SGE), LoadLeveler) either do not support a dynamically changing allocation of processors in the middle of running a job, or support it awkwardly. However, with the emergence of elastic cloud computing systems (e.g., Amazon EC2 [START_REF]Amazon Elastic Compute Cloud (Amazon EC2)[END_REF], Eucalyptus [START_REF] Nurmi | The Eucalyptus Open-Source Cloud-Computing System[END_REF]) comes a new resource provisioning model where it is feasible to allocate more processors on demand. With additional processors, a MapReduce job might finish with a lower response time. But, what programming model or system can support malleable jobs [START_REF] Feitelson | Theory and practice in parallel job scheduling[END_REF] that have the ability to add (and remove) processor resources during job execution? Custom applications can, of course, create and destroy threads on demand. But, one of the most appealing aspects of the MapReduce programming model [START_REF] Dean | MapReduce: Simplified data processing on large clusters[END_REF] is the separation of the application code from the management of processes and data, for a large class of data-parallel problems. Extending a MapReduce system to support malleable execution seems like a natural thing to do, with minimal additional complexity for the application developer. And, although the benefits of running malleable jobs on elastic cloud systems depends greatly on the specific application, the option of dynamically adding processors to a job is desirable, as long as the overheads are reasonable and the benefits are significant.

We present the design, implementation, and an evaluation of Elastic Phoenix, a MapReduce implementation for shared-memory multiprocessors. With a few minor changes to the original Phoenix [14] application programming interface (API) (Table 1) and some changes to MapReduce applications (Table 2), Elastic Phoenix maintains the original data-parallel programming model and has reasonable overheads. In our empirical evaluation, we show the overheads to be small (6% to 8%) for many applications (Figure 2). And, in a synthetic workload experiment (Figure 4), we show that Elastic Phoenix can improve the response time of the workload by up to 29%, over an equivalent workload where processor allocations can change only at job start-up time.

Phoenix [14] is a MapReduce framework written in C for shared-memory computers. Unlike other MapReduce frameworks, such as Google's MapReduce [START_REF] Dean | MapReduce: Simplified data processing on large clusters[END_REF] and Hadoop [START_REF]Apache Hadoop[END_REF], Phoenix is not distributed: a Phoenix application runs entirely in one address space on a single host, using only the local filesystem and memory. But, we argue that as the core count of individual servers increases over time, and if they have sufficient I/O capacity, the desire to run a MapReduce application on a single, sharedmemory server will also grow.

One limitation of the original Phoenix system is that, because all worker threads run in the same process, an application cannot take advantage of additional resources that become available after a job is started. That is, if a Phoenix application is started with four threads, and four more processors become idle while the application runs, there is no way to use these additional processors to accelerate the Phoenix application. In this way, we say that Phoenix is not malleable, or in the terminology of cloud computing, Phoenix is not elastic.

Elastic Phoenix is a modification to the original Phoenix framework that allows worker threads to be added or removed while a job is running, which makes our new system elastic and malleable. 1

Background

Phoenix [14] is a MapReduce framework designed for shared-memory systems. A Phoenix application is written in C and linked with the Phoenix library into a single executable. A Phoenix job runs in a single process, within which there are multiple worker threads and a single master thread. As in other MapReduce frameworks, execution proceeds in four stages: split, map, reduce, and merge. Input data is usually stored in a file on the local filesystem. The split stage is performed by the master thread, generating map tasks. The map tasks are executed by the worker threads, generating intermediate key-value pairs in memory. This intermediate data is consumed by the reduce tasks, which are also executed by the worker threads. The final data emitted by the reduce tasks is merged in multiple rounds by the worker threads. The merged data can be written to a file or processed further by the master thread.

The main function of a Phoenix application is provided by the application, not by the framework. The application configures a MapReduce job by filling in a data structure which is passed to the framework. The application must provide a map function. Most applications provide a reduce function, although an identity reduce function is provided by the framework. The application can optionally provide a splitter function for dividing input data (the default splitter treats the input data as an uninterpreted byte array and splits it into even parts); a combiner that merges intermediate values emitted by a map thread; a partitioner that divides the intermediate data among reduce tasks; and a locator function that helps the framework assign each input split to an appropriate map thread based on the location of the data. The application must also provide a comparator for the application-defined intermediate data keys.

Phoenix creates a worker thread pool, which is used for map, reduce, and merge tasks during execution of a job. Tasks are placed in a set of queues, one for each worker thread, by the master thread. The worker threads in the pool get work by pulling tasks from their queues. If a worker's queue runs out of tasks, it will steal work from other threads to avoid sitting idle.

Map tasks emit intermediate key-value pairs into a set of arrays. Each map thread has an output array for each reduce task (the number of reduce tasks is defined by the application). A reduce task processes the data produced for it by each map thread, emitting final data into another array of key-value pairs. These final data arrays are merged pairwise until a single, sorted output array is produced.

To avoid confusion, we will henceforth refer to Phoenix as "original Phoenix," in contrast to our Elastic Phoenix framework. We use the name Phoenix on its own only for discussion applicable to both frameworks.

Related Work

Malleable Programming Environments. Malleable applications are supported by the Message Passing Interface (MPI) through its dynamic process management capability [START_REF] Gropp | Dynamic process management in an MPI setting[END_REF]. However, this technique is only applicable to applications developed directly on top of MPI, which offers a low level of abstraction.

At a higher level, support for malleability has been implemented in the Process Checkpointing and Migration (PCM) library, which is built using MPI [START_REF] El Maghraoui | Dynamic Malleability in Iterative MPI Applications[END_REF]. PCM is similar to Elastic Phoenix in that it adds malleability to an existing API. However, unlike Elastic Phoenix, the PCM implementation uses explicit splitting and merging of MPI processes; that is, the application developer is required to reason about malleability in order to support it.

An unrelated project, also called Phoenix [START_REF] Taura | Phoenix: a parallel programming model for accommodating dynamically joining/leaving resources[END_REF], provides a non-MPI-based messagepassing API that supports malleability by allowing nodes in a cluster to join and leave a computation while it runs. Unlike Elastic Phoenix, it is not intended for shared-memory systems.

MapReduce Improvements. Other research groups have also tried to reduce the response time of MapReduce jobs by modifying the framework or semantics of the programming model. The Hadoop Online Prototype (HOP) [START_REF] Condie | MapReduce online[END_REF] improves response time of jobs in the Hadoop framework by modifying the flow of data during execution. In HOP, intermediate key-value pairs are sent directly from map workers to reduce workers, before a map task is completed. The framework can also produce partial results by running the reducer on partial map output, though the reducer must be executed over all the map output to produce a final result. By better pipelining and overlapping the work between the map and reduce stages, response time can be improved.

Verma et al. [START_REF] Verma | Breaking the MapReduce Stage Barrier[END_REF] present a Hadoop-based MapReduce framework in which reduce tasks not only receive but consume map output as it is produced, eliminating the effective barrier between the map and reduce stages. This technique changes the semantics of the MapReduce model, since the reducer no longer has guarantees about the ordering of its input data. Application developers are thus required to store and sort the data if necessary. Evaluation shows that if careful and efficient techniques are used to handle intermediate data, this modification to MapReduce can improve response time for certain classes of problems.

Chohan et al. [START_REF] Chohan | See spot run: using spot instances for mapreduce workflows[END_REF] explore the advantages of using Amazon EC2 Spot Instances (SI), inexpensive but potentially transient virtual machines, to accelerate Hadoop jobs. This is similar to our work in that it provides a type of malleability to a MapReduce framework: a job can be accelerated by adding worker nodes in SIs. However, because of the nature of Hadoop, termination of a spot instance is expensive, since completed work often needs to be re-executed. Elastic Phoenix does not encounter this expense in reducing the number of workers, since the data produced by all completed tasks is in shared memory and will not be lost. We note, though, that our current implementation is not suitable for the semantics of SIs, in which the worker node is terminated without warning, since an Elastic Phoenix worker requires a signal and a little bit of extra time to clean up when it is removed.

The Intermediate Storage System (ISS) [START_REF] Ko | Making cloud intermediate data fault-tolerant[END_REF] for Hadoop improves the response time of MapReduce jobs when faults occur. Hadoop normally stores map output data on local disk until it is required by a reduce worker, requiring the map task that produced the data to be re-executed if a fault occurs on the worker node. The ISS asynchronously replicates the intermediate data, saving the re-execution time if only one replica incurs a fault. However, this technique induces overhead in the absence of faults.

Design Goals

Our design goals for Elastic Phoenix were:

-Elasticity. The user should be able to add (or remove) worker threads to a running MapReduce job in order to accelerate it. The framework should not make assumptions about the stage of computation during which workers will be added.

-Design Flexibility. The implementation described in this paper uses the POSIX shared-memory mechanism, but our design assumes as little as possible about the mechanism being employed, such that it can be ported to other data-sharing mechanisms in the future. -API compatibility. It should take minimal effort to make an original Phoenix application work with Elastic Phoenix.

These design goals are sometimes conflicting. In particular, maintaining the API of the original Phoenix framework, which was not designed for elasticity, made it more difficult to achieve the other goals. The basic design change that enables Elastic Phoenix is separation of the master thread and worker threads into different processes, as illustrated in Figure 1. This presents a number of challenges, since original Phoenix was designed with a shared address space in mind. In addition, our implementation does not assume that all worker processes share a kernel or filesystem, allowing future versions to support distributed operation using another distributed data-sharing mechanism such as distributed shared memory or message passing.

Implementation

3KRHQL[3URFHVV

Sharing Data

One challenge in designing Elastic Phoenix was how to get input data to the worker threads. The master thread splits the input data, but the worker threads then need to access this data for the map phase. A common idiom in original Phoenix applications is to mmap an input file, then have the splitter calculate offsets into the mapped file. Since Elastic Phoenix threads do not share an address space, this clearly will not work. In the interest of flexibility, we also assume that the worker processes cannot access the input file.

Our initial design, which made no changes to the original Phoenix API, required splitter functions to return pointers to the actual input data. The framework then copied this data into shared memory. This was a simple requirement, but one that did not exist in original Phoenix. As a result, some original applications returned pointers to structures that describe the input data, rather than pointers to the actual data itself. This design allows applications to continue using the mmap method described above, with most applications requiring no changes to the splitter code.

However, copying the input data into shared memory after it has been split is a performance bottleneck. Most application-provided splitters must examine at least some portion of the data in order to generate sane splits. Therefore, the copying strategy ends up reading in the data from disk during the split, then copying the data into the sharedmemory region, and possibly reading the data off disk again if the buffer cache has been recycled.

Our final design changes the API for splitter functions by adding an argument: a structure containing pointers to memory-allocation functions, which we require the splitter use to allocate space for the input data. The functions we pass to the splitter are custom memory-allocation functions that allocate space in the shared-memory region, keeping all allocation metadata in the region itself. The common design for splitters using this API is to allocate memory for a split, read data directly into this memory, then reduce the size of the split as necessary to reflect application-specific restrictions on the split data.

In addition to the input data, the master and the workers must share the intermediate and final data arrays. The intermediate data is straightforward to move into the shared-memory region, using the previously mentioned memory allocators to dynamically allocate the arrays. The only complication is that there is a set of intermediate data arrays for each worker thread, and in Elastic Phoenix the total number of worker threads is not known at initialization time. We deal with this by setting an upper bound on the number of map threads and pre-allocating space for the necessary pointers.

There is one additional issue with the final data array: it is returned to the application when the MapReduce job finishes. This causes a problem because many applications free this array when they are done processing the results it contains. This works in original Phoenix because the framework allocates the final data array with malloc. In Elastic Phoenix this array resides in the shared-memory region and is allocated with our custom allocator. We deal with this by adding a function to the API: map reduce cleanup. The cleanup function mirrors the framework's map reduce init function, and is called by the application when it is finished with the results of the MapReduce job.

Assigning Tasks

In original Phoenix, there is one task queue for each worker thread. The framework does allow threads to steal work from other queues when their own queue runs out, so pre-allocating a task queue for each of the maximum number of worker threads is one potential solution. However, given that Elastic Phoenix will rarely have the maximum number of worker threads, this would lead to the threads exhausting their own queues, then all sharing the other queues until they are exhausted. Given this, we opted for a simpler design with a single task queue. The task queue is pre-allocated in shared memory with a maximum number of tasks, and all worker threads take work from the queue. Since the task queue in original Phoenix was already thread-safe for dequeue operations, we re-used nearly all of its code.

The pre-allocation design causes an issue for the splitter: there is a maximum number of map tasks, so the splitter can only create so many splits. Elastic Phoenix uses the application-provided unit size estimate to request an appropriate number of data units for each split such that the task queue is filled but not exceeded. However, an application's estimated unit size will, in some cases, be inaccurate because the data is not uniform and not known beforehand. In this case, the framework may run out of space for map tasks before all the data has been split. We deal with this situation by requiring that splitters be "resettable". When the splitter is passed a special argument, it should thereafter operate as if it has not yet been called. When we run out of space for map tasks, we reset the splitter and re-split the data, requesting more units for each task. This is a rare case, as many applications provide good unit size estimates. We have observed that when re-splitting is required it often takes only one try, which tends to be fast due to filesystem caching.

Coordinating Workers

An Elastic Phoenix job is started by first running the master process, then running one or more worker processes. Coordinating the worker threads is somewhat different from original Phoenix, since the master does not create the worker threads.

We use two basic mechanisms for coordination in Elastic Phoenix: a worker counter and a barrier, both stored in shared memory. When the master starts, it executes the split phase, then waits for workers to join the computation. The first thing a worker process does when it starts is acquire a worker ID using the counter in shared memory. This is done using atomic operations to avoid the need for a lock.

Because we do not assume that all workers run under the same kernel, we cannot use a blocking barrier for synchronization between computation phases. Instead, we have implemented a barrier using spinlocks and polling. The barrier is used after each stage to synchronize the master and worker threads. The master thread reaches the barrier first, since it does no work in the map, reduce, and merge phases. To make the barrier code simpler, we give each worker process a constant number of threads; in our current implementation the default is four. In addition, we use a shared spinlock to prevent new workers from joining the computation at inconvenient times, such as after the first worker thread reaches the barrier (since this would change the number of threads expected, complicating the barrier code) or before the master has finished setting up the necessary shared environment.

Elastic Phoenix keeps track of the current phase of execution using a flag in shared memory. This allows workers to join during any stage of execution: split, map, reduce, or merge. A worker that joins during a later phase skips the previous phases entirely.

An additional coordination problem in Elastic Phoenix is that because we have attempted to preserve the original Phoenix API as much as possible, the application code is not intended to be executed in multiple processes. Many original Phoenix applications contain code that should be executed only once as part of initialization or cleanup. We add two optional application-provided functions to the API to deal with this: prep and cleanup. If an application may provides these functions, Elastic Phoenix guarantees that they will be executed only in the master process. The prep function is called before any other work is done, and can be used to, for example, open the input file and set up data for the splitter. The cleanup function is called by the map reduce cleanup function described earlier. These functions reduce the need for applications to explicitly manage concurrency, and allow our implementation to handle initialization and cleanup differently in future versions. One of our design goals, as described in Section 4, was to allow original Phoenix applications to be easily ported to Elastic Phoenix. In this section, we describe our experiences porting applications, and try to quantify the amount of change required. Table 1 lists the complete set of API changes.

Porting Phoenix Applications

The biggest change that is required in every application is to the splitter function. Since Elastic Phoenix requires that the splitter use a framework-provided memory allo-cator for allocating input splits, it has an additional argument. In some splitters, this is simply a matter of replacing calls to malloc with calls to the provided allocator, but most need to be restructured to read in a file incrementally instead of using mmap, as described in Section 5.1.

In Phoenix applications, most of the code in main serves one of four purposes: setting up the splitter; setting up the MapReduce job; processing the results; and cleaning up. The first can usually be moved verbatim into the prep function. The second is usually idempotent, operating only on local data structures, and can remain intact without harm. The third and fourth can usually be moved to the cleanup function to ensure they execute only in the master process. The original Phoenix framework includes seven sample applications: histogram, kmeans, linear regression, matrix multiply, PCA, string match, and word count. Of these, we have ported histogram, word count, and linear regression. PCA and k-means require multiple MapReduce jobs, a feature not currently supported in Elastic Phoenix (we discuss this limitation in Section 8). The matrix multiply included in original Phoenix is not written in the MapReduce style: the map tasks calculate the resultant matrix directly into a known output location, rather than emitting the calculated values as intermediate data. We wrote a new matrix multiply for original Phoenix that uses the MapReduce style, emitting a key-value pair for each entry in the resultant matrix, then ported it to Elastic Phoenix. String match also does not emit any intermediate or final data; in fact, it produces no results at all.

To quantify the amount of work required to convert an application to Elastic Phoenix, we analyzed the applications we converted. First, we calculated the lines of code for each application before and after conversion using SLOCCount [START_REF]SLOCCount[END_REF]. Second, we estimated how many lines were altered by taking a diff of the two source files and counting the number of altered lines in the Elastic Phoenix version. This does not entirely account for the fact that many of the changes simply involved moving code from one function to another. The results of this analysis are displayed in Table 2. Of course, lines of code are an incomplete and controversial metric for programming effort, but we include it as a data point for consideration.

Word count was the first sample application we ported, and we modified it incrementally during the development of Elastic Phoenix. Therefore, it underwent more change than would an application being ported directly to the final version of Elastic Phoenix. We ported histogram and linear regression to one intermediate version of Elastic Phoenix, then to the final version, so they display less change than word count. We developed our new version of matrix multiply for original Phoenix first, then we ported it directly to the final version of Elastic Phoenix. Therefore, matrix multiply provides the best indication of the effort needed to directly port an application. Anecdotally, a single developer wrote the new matrix multiply application for original Phoenix in several hours, and took only a few minutes to port it to Elastic Phoenix.

Evaluation

We evaluated Elastic Phoenix in two ways. First, we evaluated the overhead of using Elastic Phoenix by running our ported applications with both frameworks. Then, we explored the advantages of Elastic Phoenix by developing a workload of repeated MapReduce jobs on a system the availability of varies over time. All experiments were performed on an Intel Xeon X5550 (2.67 GHz) server with two sockets, eight cores, 16 hyperthreads, and 48 GB of RAM. We compiled original Phoenix with its default tuning parameters, and used the default four worker threads per process for Elastic Phoenix. The performance advantage of malleability is that MapReduce jobs can be accelerated as they run by adding threads to the computation. If the overhead of Elastic Phoenix is so large that, for example, using four threads in original Phoenix is faster than using eight threads in Elastic Phoenix, then this advantage disappears. In this section we uate the overhead using Elastic Phoenix, showing that for some applications it is possible to increase performance by adding threads.

To evaluate the overhead of our design, we ran four applications under both original and Elastic Phoenix. For each application, we used the largest input data we could given 4 GB of shared memory (we discuss these limits in Section 8). The applications we used and their input sizes are shown in Table 3. All times are averages over five runs, with cold filesystem caches.

We divide these applications into two classes, which we will discuss separately: I/Obound and CPU-bound. Histogram, word count, and linear regression are I/O-bound. This is due to the fact that the computation performed in the map and reduce stages is quite trivial, even with large input. The bottleneck, therefore, is reading and splitting the input data, which must be done serially in both original and Elastic Phoenix. Our matrix multiply application is CPU-bound because its input is generated in the splitter function, avoiding any file I/O.

Benchmark results for the I/O-bound applications are presented in Figure 2. Being I/O-bound, we do not expect any speedup as we add threads. We observe that Elastic Phoenix has little overhead compared to original Phoenix for the histogram and linear regression applications (8% and 6%, respectively, on average), showing some slowdown as we add more threads due to the additional coordination overhead. These applications have large input data, and generate intermediate data with a fixed number of keys.

In contrast, we observe significant overhead (1792% average) in the word count application. This overhead likely stems from the fact that the input size is small, and that the application produces a large amount of intermediate data with diverse keys, since many more words appear infrequently in the input than appear frequently. This causes a large amount of extra shared memory to be allocated but never used, as the framework allocates space for multiple values (10, by default) with the same key when the first value for a given key is emitted. In fact, word count, with only 10 MB of input data, can be run with at most 12 worker threads in Elastic Phoenix; with more threads, the small amount of additional per-thread allocation in shared memory causes it to exceed the 4 GB limit. Our shared-memory allocator is simple, and does not match the efficiency of the system's malloc, which is used in original Phoenix.

Benchmark results for our single CPU-bound application (i.e., matrix multiply) are shown in Figure 3. Here we observe significant speedup from four to eight threads, then moderate speedup as we increase to 12 and 16 threads. We theorize that the extreme speedup from four to eight threads is caused by the fact that our CPUs are quad-core, and with eight threads Phoenix uses both sockets. Original Phoenix pins threads to cores, spreading work out as much as possible, but we have configured it to use the same number of processors as threads, so the four threads are pinned onto one CPU. Using both sockets doubles the total usable cache size and the total memory bandwidth available. The decline in speedup beyond eight threads is likely due to the fact that we have only eight real CPU cores, and the hyperthreaded cores offer little advantage in an With four threads (i.e., a single worker process), Elastic Phoenix adds no overhead compared to original Phoenix. Elastic Phoenix adds a moderate amount of overhead when more than one process is used, due to the added overhead of coordinating threads in multiple processes, such as increased contention on spinlocks. Overall, the overhead averages 38%. Nonetheless, Elastic Phoenix with eight threads is significantly faster than original Phoenix with four threads, suggesting that the ability to add additional workers as a job runs can improve performance in some cases.

Elasticity Advantages

Elastic Phoenix can improve throughput and job latency by allowing idle processors to be used immediately by a running MapReduce job. A performance improvement is possible when the overhead of using Elastic Phoenix does not exceed the speedup provided by the additional threads. As we showed earlier, this is the case for some applications. In this section, we evaluate how much performance can be gained in such a situation. For this benchmark, job latency can be reduced by 29% by adding threads as processors become available.

Consider a scenario in which a sequence of MapReduce jobs is run repeatedly, starting again each time it finishes. For example, the sequence may consist of a job that creates a link graph from a set of hypertext pages, then a job that computes the most popular page. Assume the jobs run on a system where CPU resources become available unpredictably, as is typical in a multi-user, batch-scheduled environment. Original Phoenix always uses the number of CPUs that are available when it starts, while in Elastic Phoenix CPUs are added to the computation as they become available.

Our benchmark uses alternating linear regression and matrix multiply jobs with random-sized inputs of 250 MB to 500 MB for linear regression and 500 x 500 to 1000 x 1000 for matrix multiply. This sequence represents a typical workload of an I/Obound job that generates input for a compute-bound job, although we do not actually use the linear regression output as input for the matrix multiply. With one worker thread, the jobs take between 1 second and 5 minutes each. For this experiment, we compiled Elastic Phoenix with one thread per worker process, so that it can take advantage of a single processor becoming available.

We generated a set of random CPU availability changes, where every 5 seconds the number of CPUs available can change, with at least 2 and at most 8 CPUs always being available, since our host has 8 CPU cores (we do not employ the 8 hyperthreaded cores). The workload we generated contains 500 jobs in total: 250 linear regression jobs interleaved with 250 matrix multiply jobs. The results of the experiment are displayed in Figure 4. Original Phoenix completed the jobs in 5,451 seconds and Elastic Phoenix completed the jobs in 3,866 seconds, a 29% improvement. We see from these results that adding processors as they become available can drastically speed up a workload of this type.

Note that we have successfully tested the ability to dynamically remove workers from a running Elastic Phoenix job, for proper malleability, but we do not present a performance evaluation in this paper.

Limitations and Future Work

The design of Elastic Phoenix introduces a number of limitations not present in original Phoenix. Some of these are inherent to the design we have chosen; others are due to our current implementation and could be eliminated in the future without a significant design change. There are also some ways in which Elastic Phoenix could be improved, orthogonal to these limitations. Job Size. Currently, the size of the shared-memory region used with Elastic Phoenix limits the amount of input data that can be used. This is because the input data must all fit in shared memory, and will share the space with some of the intermediate data (input data is freed immediately after being used by a map task, but intermediate values are stored to shared memory as they are produced). Similarly, the intermediate data must share the space with some final data, so this data is also limited in size. The input size limit is not the same for all applications and cannot be calculated in advance. As mentioned in Section 7.1, the number of intermediate key-value as well as the diversity of the keys, affects the maximum input size.

Combiners. The original Phoenix framework allows applications to provide a combiner function, which performs a partial reduce on the map output of a single worker thread. This can speed up the reduce stage, and reduce the amount of space used for intermediate data. In particular, it helps applications that have a large number of repeated intermediate keys.

Elastic Phoenix currently does not support combiners. The application may specify a combiner, to keep the API consistent with original Phoenix, but it will not be used by the framework. This is a limitation of the current implementation and could be changed in the future. It is possible that this modification would help with the job size limitation by compacting intermediate data.

Multiple MapReduce Jobs. A common idiom in MapReduce applications, including

Phoenix applications, is to run multiple MapReduce jobs in sequence. For example, one might run a job to do a large computation, then another job to sort the results by value. In original Phoenix, multiple jobs can be performed without writing data to disk: the output data from one job becomes the input data to another job. That is, multiple calls are made to the Phoenix MapReduce scheduler within a single executable application.

The current Elastic Phoenix implementation does not support running multiple Map-Reduce jobs in one execution, due to the way in which Elastic Phoenix initializes data that will be shared among the master and worker processes. This limitation could be eliminated through some further design.

Number of Workers.

Elastic Phoenix pre-allocates a number of data structures that in original Phoenix are allocated based on the number of threads being used. In particular, there is one set of intermediate data queues and one final data queue for each worker thread, and the pointers for these are pre-allocated. They could be allocated dynamically as workers join the computation, but this would require using a linked list or similar structure, making lookups slower when inserting data.

In order to pre-allocate this data, we set a fixed limit on the number of worker tasks that can join a MapReduce job. The default limit in our current implementation is 32 threads, which allows for eight worker processes with four threads each. On any given system, the number of worker threads is naturally limited by the number of processors, so this is a minor limitation.

Other Future Work. Elastic Phoenix makes extensive use of spinlocks instead of blocking locks. This is because we designed Elastic Phoenix with flexibility in mind, and for a potential hybrid shared-memory and distributed-memory version of Elastic Phoenix where the worker processes may not share a kernel, making the use of standard blocking locks, such as pthreads mutexes and semaphores, impossible. One potential performance improvement would be to use blocking locks when they are available, or to exploit optimistic concurrency by making more extensive use of atomic operations. This would likely involve abstracting locking and atomic operations so that they could be implemented in an appropriate way for any given data sharing system.

The original Phoenix software implemented many optimizations to improve cache usage and memory locality for worker threads. In the process of building Elastic Phoenix, we have eschewed many of these optimizations, either because they do not apply in the case where worker threads are in different processes or because getting the initial version of Elastic Phoenix was easier without them. In the future, we may re-introduce some of these optimizations.

Concluding Remarks

We have presented the design, implementation, and an evaluation of Elastic Phoenix, a MapReduce implementation for shared-memory multiprocessors. Based on the Stanford Phoenix system, Elastic Phoenix adds the ability to support malleable jobs that can dynamically add processors to the computation, with the goal of reducing the turnaround time of the job or workload.

Using a synthetic workload (Figure 4), we have shown that Elastic Phoenix can improve the response time of the workload by up to 29%, over an equivalent workload where processor allocations can only change at job start-up time. Although the benefits of malleability or an elastic programming system will vary from application to application, Elastic Phoenix provides those benefits while maintaining the MapReduce programming model and much of the Phoenix API.

Fig. 1 .

 1 Fig. 1. High-level architecture of (a) original Phoenix and (b) Elastic Phoenix

Fig. 2 .

 2 Fig. 2. Overhead experiment results for I/O-bound applications

Fig. 3 .

 3 Fig. 3. Overhead experiment results for the CPU-bound matrix multiply application

Fig. 4 .

 4 Fig. 4. Job completion over time in the elasticity benchmark

Table 1 .

 1 API changes from original Phoenix to Elastic Phoenix

	Description		Original Phoenix	Elastic Phoenix
	Splitter signature	function	int splitter(void * , int, map_args_t * )	int splitter(void * , int, map_args_t * , splitter_mem_ops_t * )
	Splitter uses pro-vided allocator	free(out->data); out->data = malloc(out->length);	mem->free(out->data); out->data = mem->alloc(out->length);
				if(req_units < 0) {
	Splitter can be reset N/A	lseek(data->fd, 0, SEEK_SET); data->fpos = 0;
				}
	New API function cleanup	printf("Num. values: %d\n", mr_args.result->length); return(0);	printf("Number of values: %d\n", mr_args.result->length); map_reduce_cleanup(&mr_args); return(0);
			int main(int argc,	int prep(void * data,
			char ** argv) {	map_reduce_args_t * args) {
			...	struct stat s;
	New prep function application	struct stat s; stat("in", &s); data.fd = open("in");	stat("in", &s); ((app_data_t * )data)->fd = open("in");
			mr_args.data_size =	args->data_size =
			s.st_size;	s.st_size;
			...	}
			int main(int argc,	int cleanup(void * data) {
	New cleanup function application	char ** argv) { close(data.fd); ...	app_data_t * ad = (app_data_t * )data; close(ad->fd);
			...	}
	API init modifies argc and argv	int nargs = argc; char * fname = argv[1]; map_reduce_init();	map_reduce_init(&argc, &argv); int nargs = argc; char * fname = argv[1];

Table 2 .

 2 Porting MapReduce applications: original to Elastic Phoenix (lines of code)

	Application	Original Elastic Lines altered
	Histogram	245	285	110
	Word Count	234	202	114
	Linear Regression 230	243	97
	Matrix Multiply 168	185	48

Table 3 .

 3 Applications and inputs for evaluation

	Application	Description	Input Size
	Histogram	the frequency of each color value over pixels in a	1.4 GB
		bitmap image
	Word Count	the frequency of each word in a text file	10 MB
	Linear Regression a linear approximation of a set of points	500 MB
	Matrix Multiply Parallel matrix multiplication using row and column blocks 1000x1000

Elastic Phoenix is an open-source project and can be obtained from https://github. com/adamwg/elastic-phoenix.

Acknowledgements. This research was supported by a grant from the Natural Sciences and Engineering Research Council of Canada (NSERC).


Andrea Scherlozer

Mestaro Orsini

Sulvane Patole

Simulation and Numerical Analysis and Comparative Study of a PID Controller Based on Ziegler-Nichols and Auto Turning Method

Keywords: Control methods, Ziegler Nichols, Steady state error, Stability, PID controller, Simulation, Numerical study

whether they are published or not. The documents may come   L'archive ouverte pluridisciplinaire

Introduction:

Control systems has been used extensively during past decades. One of the most important and successful control systems is PID controller, specifically in industries where more than 90 percent of industrial controllers were PID family on 2002. The three term PID controller has three basic modes for controlling the characteristics of a second order system. The three basic modes are: the Proportional term, the Integral term as well as the Derivative term. The PID controller is basically used in controlling closed loop form of an open loop system. But it is broadly applicable since a PID controller relies only on the measured process variable, not on knowledge of the underlying process. The proportional term of the controller is in proportion to the error in the system as the name implies, the integral term is proportional to the integral of the past errors while the derivative term is proportional to the rate of change of the error [START_REF] Hunter | Automated process control systems: concepts and hardware[END_REF]Krishnaswamy, 2011). Where is the integral time constant and the derivative time constant. The proportional part acts on the present value of the error, the integral represents an average of past errors and the derivative can be interpreted as a prediction of future errors based on linear extrapolation.The controller can also be parameterized mathematically as:

( ) ( ( ) ( ) ( ) ) (1)
It is worthy of note that the control signal u is formed entirely from the error e. Overall in any system the Proportional term, the Integral term as well as the derivative term contributes in achieving a fast rise time, minimum overshoot, no oscillations and higher stability as well as no steady-state error (Zilouchian and Jamshidi, 2001;Owunna et al., 2016). The role of each of the terms in a PID controller are tabulated in Table 1. The proportional controller often times reduces the rise time and often reduces the steady state error of the system but never does it eliminate the error completely. The integral control on the other hand eliminates completely the steady state error of a system. A derivative control however increases the stability of the system, reduces the overshoot as well as improves the transient response of the system. It does not in any way alter the steady state error of the system. In this paper we compare our methods simulation result with the fuzzy system controllers implemented on aerial and surface vessels [START_REF] Abbasi | Controlling of quadrotor uav using a fuzzy system for tuning the pid gains in hovering mode[END_REF][START_REF] Yazdanpanah | Fuzzy LQR controller for heading control of an unmanned surface vessel[END_REF]. Our model has been developed based on the results of these papers and comparison with their model.

Type 1 Systems

Type 1 systems are systems that do not have any steady state error. Since these systems have no steady state error, the presence of integrators will be superfluous to the system. This is because the major function of integrator in controllers is to eliminate steady state errors which are absent in type 1 open loop systems [START_REF] Cooper | Car Cruise Control Loop Block Diagram[END_REF]. The significant challenge connected with the derivative controller noise issues and sensitivities that, a large frequency within a system associated with large changes in the system error may cause the derivative of the signal to amplify the signal significantly. Thus little levels of noise present in the system may cause the output of the system to increase greatly. In other words, the sensitivity of derivative controllers to noise may result in significant changes in the value of the output as a result of small level of noise in the system. In these circumstances, it is often sensible to use a PI controller or set the derivative action of a PID controller to zero. To eliminate/minimize this downside, an electronic signal filter may be enclosed within the loop. Electronic signal filters are unit electronic circuits that perform signal process functions, specifically supposed to get rid of unwanted signal components and/or enhance needed ones. Electronic filters can be: passive or active, analogue or digital, discrete-time (sampled) or continuous-time, linear or non-linear, etc. The most common types of electronic filters are linear filters, regardless of other aspects of their design.

Ziegler-Nichols Closed-loop Tuning Method

In the Ziegler-Nichols closed-loop tuning method, the ultimate gain and the ultimate period of oscillations is employed in calculating the needed which is the value of the proportional gain required for effective tuning of the system. The Ziegler-Nichols closed-loop tuning method is only applicable in closed-loop systems and cannot be applied in open loop systems. To determine the value of , the value of the proportional gain that will produce a steady oscillation in the system is first obtained. The gains for the integrator and the derivative controllers are initially set at zero for the procedure. When the systems oscillates steadily, the period of oscillation must therefore be obtained as it is required in calculating the integral and derivative times. The ultimate period is the time required to complete one full oscillation while the system is at steady state. To find the values of the PID parameters from the values of and period obtained, the following procedures must be adopted.

Closed Loop (Feedback Loop)

i.

The derivative controller gain and the integral gains must be set at zero. ii.

The proportional gain should then be varied till the system oscillates at constant amplitude. iii.

The values of and the period of the oscillation can then be recorded . To obtain the various values for the PID controllers Ziegler-Nichols equation must be used and the equation is presented in the Table 2. Closed Loop Systems-P Controller P controllers are often used in first order systems to stabilize unstable responses. P controller helps to majorly reduce the steady state error of the system. An increase in the proportional gain factor K of the P controller reduces the steady state error of the system (Ogata, 1997). It is however worthy of note that P controllers can reduce but not eliminate totally the steady state error of a system. As the proportional gain of the P controller increases, smaller amplitudes as well as smaller phase margin are introduced to the system. The dynamic of the system also becomes faster and the sensitivity of the system to noise reduces as the proportional gain increases. The system is applicable only in instances where the system can tolerate constant steady state error (Taeib and Chaari, 2015).

P-I Controller

P-I controller are used majorly in the elimination of steady state errors arising from P controllers. The PI controllers have a negative effect on the stability of a system as well as the response speed of the system. It is therefore important to note that P-I controllers are useful in systems where the response speed is insignificant. P-I controllers have no effect on the rise time and cannot eliminate oscillations in a system because they cannot accurately predict future errors within the system.

P-I-D Controller

With PID controllers, zero state errors are possible. The response of the system can be improved to achieve a fast response, oscillations in the system can be removed and the stability of the system can be improved. A derivative of the output response is often added to a PI controller to remove overshoot and oscillations in the system. PID controllers have the advantage of use in higher order systems.

P-D Controller

P-D controllers are used majorly to increase system stability as the controller is able to predict future error that can occur in the systems response. A derivative of the output response is often used instead of using the error in the signal to ensure there is no abrupt change in the value of the error of the signal (Padula and Visioli, 2011; Taeib and Chaari, 2015). The derivative controllers is often not used alone to prevent amplification of noise in the system.

Research Methodology

For optimum performance of control systems, the steady state errors must be eliminated to enable stability of the close loop systems which are basically controlled by PID. The closed loop transfer function was designed in Simulink and MATLAB M-file was generated to plot responses of the transfer function with different proportional gain. Ziegler-Nichols closed-loop tuning method and auto tuning system method to determine PID values and a MATLAB command was generated and simulated for both tuning methods. Result obtained from the MATLAB simulation was used to determine if the steady state error has been eliminated or not, and how effective each method is.

Steady State Error (P controllers)

For the given transfer function

( )

Recall that for a closed loop system,

* ( ) ( )+ (2)
For unit step input, ( )

Thus, * ( )

( )+ * ( ) + (4)
* ( ) + ( )

( ) ( ) (6)
For steady state error

( ) ( ) ( ) Therefore, * ( ) ( )+ * ( ) + [ ] ( ) ( )
For stability of the system, using Routh Hurwitz criterion, 1 5 0

For stability,

But

For , This implies that if the system will remain stable, there will remain within the system steady state error if a proportional controller is used. To establish the fact that there will remain within the system steady state error if a proportional controller is used, a MATLAB M-file was generated to plot responses of the transfer function with different proportional gain. The MATLAB code for the simulation is written as follows:

% The plot generated from the system shows plots of step response from the real transfer function and the closed loop systems with different proportional gains as well as the target response. Figure 1 shows the plots with different proportional gains, from the plots shown in Figure 1 it is obvious that with increasing value of proportional gain, the steady state gain remained.

PI Controllers

For PI controller with proportional gain and integral gain the closed loop transfer function becomes;

( )

[ ] ( ) (7) [ ] [ ] ( )
For steady state error

This implies that at any value of with proportional gain and integral gain , there will be no steady state error within the system. To establish the fact that there will no remains of any steady state error within the system if a proportional-Integrator controller is used, a MATLAB M-file was generated to plot responses of the transfer function with different integrator gains. The MATLAB code for the simulation is written as follows:

% To obtain the Unit-Step Response of the System num=[0 0 10]; % num_sys is the numerator of the system transfer function den= [1 5 6]; % den_sys is the denominator of the system transfer function % specify proportional and integrator controllers

Ziegler-Nichols Closed-loop Tuning Method

The derivative controller gain and the integral gains were set at zero. The proportional gain was varied, until a relatively stable system was obtained at Ku=8.5. The response obtained is shown in Figure 3. At this gain, the period of oscillation is obtained as 0.7. From Ziegler-Nichols method, the PID values are given as shown in Table 3. The PID Values if used in plotting the response of the system will give a situation as shown in Figure 4.

Figure 4: System response with Zeigler-Nichols Tuning Comparing the initial response and the final response, it is obvious that the steady state error is eliminated using the Ziegler-Nichols method of the analysis. A backdrop of this method is that, it is time consuming and may delay while entering into an unstable region for the system.

MATLAB PID Tuning

It is necessary to tune the system represented above using automatic Simulink tuning. The different steps used are as follows [START_REF] Nguyen | Simulink TUTORIAL[END_REF]; i. Design the closed loop transfer function in Simulink The design as generated from Simulink is shown in Figure 5. open the tune panel and increase the system response and increase the system robustness to highest

The result of the tuning is shown in Figure 6. The values obtained for the system is Kp=11.72, Ki=14.01, Kd=2.449 v.

At this value, the settling time reduced to 0.172s while the rise time reduced to 0.0909s.The system is very stable.

Comparing the tuned system to the initial transfer function response, Figure 7 was generated using MATLAB M-Codes

Conclusion

From the graph in Figure 8, it can be observed that both tuning methods removed the steady state errors within the system, thereby indicating stability of the close loop system. It is however easy to move the system in any direction with the auto tuning method than the Ziegler method, as the Ziegler method is a quick approximation of results.

Limitations of PID control

Although PID controllers can be used in many control situations with satisfactory performance, their performance in other applications may be relatively poor with no optimal performance. PID controllers, when used in cases that are non-linear may be unable to respond to the fluctuations in process behaviour and may ultimately lag in their response to large disturbances. To solve such discrepancies, a knowledge of the control system can help include a feedforward control allowing the PID controller to only deal with the steady state error [START_REF] Foley | A comparison of PID controller tuning methods[END_REF]. Another serious challenge with the use of PID controllers is that they are linear and symmetric. Their performance in non-linear systems is unpredictable. This means that overshoot cannot easily be corrected like in linear systems. In trying to reduce overshoot in non-linear systems, the performance of the system may be compromised. Another issue with PID systems is that the derivative term can amplify high frequencies in a system and for non-linear systems, any large frequencies may be amplified to cause large discrepancies in the output of the system.

  with gain the closed loop transfer function becomes;

  To obtain the Unit-Step Response of the System num=[0 0 10]; % num_sys is the numerator of the system transfer function den=[1 5 6]; % den_sys is the denominator of the system transfer function K=1,den); %Gs is the system transfer function %To plot the step response figure(1) %specifies figure number step(Gs, feedback (Gs*K,1),feedback (Gs*K2,1),feedback (Gs*K3,1)) hold on plot([-0.1,0,0,4],[0,0,1,1],'r'); axis([-0.1 4 0 1.8]); legend('Gs','K=1','K2=2','K3=5', 'target'); %specify title and grid grid on title ('Unit-Step Response of the system')

Figure 1 :

 1 Figure 1: Plot of Transfer Function with Proportional Controller

Figure 2 :

 2 Figure 2: Plot of Transfer Function with Proportional-Integrator Controller

Figure 3 :

 3 Figure 3: System response at Ku=8.5

Figure 5 :

 5 Figure 5: Simulink design for the closed loop transfer function ii. run the simulation iii.open the tune panel and increase the system response and increase the system robustness to highest

Figure 6 :

 6 Figure 6: Tuned Response with Auto-tuning iv.The values obtained for the system is Kp=11.72, Ki=14.01, Kd=2.449 v.At this value, the settling time reduced to 0.172s while the rise time reduced to 0.0909s.The system is very stable.

Figure 7 :Figure 8 :

 78 Figure 7: Tuned Response with Auto-tuning compared to

Table 1 : Roles of PID controller terms

 1

	Term	Response	Rise Time	Overshoot Settling Time S-S Error
	Proportional		Decrease	Increase	Small Change Decrease
	Integral	( )	Decrease	Increase	Increase	Eliminate
	Derivative		Small Change Decrease	Decrease	No Change

Table 2 : Ziegler-Nichols closed-loop tuning formula

 2

	Rule Name	Tuning Parameters
		Kp	Ki	Kd
	P	0.5 Ku
	PI	0.45Ku 1.2Kp/Pu
	PID	0.6 Ku 2Kp/Pu	KpPu/8

Table 3 : Ziegler-Nichols tuning parameters

 3

	Rule Name Tuning Parameters
		Kp	Ki	Kd
	P	4.25

Table 4 : Table showing the Results for both Tuning Methods

 4

	Ziegler-Nichols Auto tuning
	P Kp=5.1	Kp=11.72
	I Ki=24.28	Ki=14.01
	D Kd=0.74375	Kd=2.449

Kp=1; Ki=1; Ki2=2; Ki3=5; %To combine the P-I controllers K=pid(Kp,Ki); K2=pid(Kp,Ki2); K3=pid(Kp,Ki2); % To obtain the transfer function s=tf('s'); Gs=tf(num,den); %Gs is the system transfer function %To plot the step response figure [START_REF] Cooper | Car Cruise Control Loop Block Diagram[END_REF] %specifies figure number step(Gs, feedback (Gs*K,1),feedback (Gs*K2,1),feedback (Gs*K3,1))%to plot %all feedback functions hold on plot([-0.1,0,0,4],[0,0,1,1],'r'); axis([-0.1 4 0 1.8]); legend('Gs','Ki=1','Ki2=2','Ki3=5', 'target'); %specify title and grid grid on title ('Unit-Step Response of the system')

The plot generated from the system shows plots of step response from the real transfer function and the closed loop systems with different proportional-Integrator gains as well as the target response. From the plots shown in Figure 2, it is obvious that at whatever value of the integrator gain, the steady state error is always eliminated [START_REF] Messner | Control Tutorials for MATLAB and Simulink -Introduction: PID Controller Design[END_REF].

To compare the results from the Ziegler-Nichols closed-loop tuning method with the auto tuning system, a MATLAB m-file was written. The codes written for the comparism are given below: % To obtain the Unit-Step Response of the System num=[0 0 10]; % num_sys is the numerator of the system transfer function den= [1 5 6]; % den_sys is the denominator of the system transfer function % specify proportional and integrator controllers Kp=11.72; Ki=14.01; Kd=2.449; The result of the simulation is given in Figure 8


Klaus Szielasko

Meisam Amiri

Benjamin Ducharne

Bhaawan

Bhaawan Gupta

Magnetic Susceptibility Imaging as a New Approach towards Characterization and Testing of Para-and Diamagnetic Materials

published or not. The documents may come     Magnetic Susceptibility Imaging as a New Approach towards Characterization and Testing of Para-and Diamagnetic Materials

Introduction

Micromagnetic materials characterization is a well-known discipline of nondestructive evaluation and testing (NDE/NDT). It is based on the interaction of magnetic domain walls (Bloch walls) with microstructure, which is similar to the interaction of dislocations with microstructure [1 -3]. Every material interacts with a magnetic field to some extent, depending on its internal composition and structure. The interactions which define the material's behavior in a magnetic field take place on the micro-and nanometer scale and are partially quantum mechanical in nature. The magnetic susceptibility  describes how a material interacts with a magnetic field. It is assumed that a correlation exists between micro-or electromagnetic and mechanical properties in para-and diamagnetic materials as well. There is a lack of knowledge regarding the correlation between the mechanical properties and magnetic susceptibility of nonferromagnetic materials such as graphite, aluminum and plastics. This contribution presents an approach for magnetic susceptibility imaging and its applications in the field of NDE/NDT for non-ferromagnetic materials.

Experimental set-up

The novel magnetic force sensor includes a capacitive distance sensor in order to determine the deflection of a cantilever carrying a sensing magnet [START_REF] Gupta | [END_REF] (Fig. 1). This way, the magnetic force sensor can theoretically be used with common eddy current inspection equipment. Under the influence of the sample, the magnet is either repelled or attracted, which leads to a deflection of the cantilever. The output voltage change of the capacitive sensor electronics when approaching a material is proportional to the deflection of the cantilever. When investigating practical materials such as plastics or aluminum, the deflections are so small that the distance between sample and magnet changes to very small extent. In this case, the deflection (or voltage change) is almost proportional to the susceptibility of the material. Figure 2 depicts the diamagnetic behavior of a PVC sample. Many other non-ferromagnetic materials under test verified the efficiency of the sensor.

Conclusion

The principle is based on measuring the force acting upon a permanent magnet in the proximity of the material to be inspected. The sensitivity was sufficient for all kinds of plastics, metals, glass and adhesives investigated so far. Future work will focus on the improvement of lateral resolution and scanning speed, e.g. by increasing the magnetic field gradient and lowering the distance between magnet and sample. The improvement will also include the exclusion of cantilever to avoid possible vibrations. To minimize the cost, a strain gauge developed at LGEF INSA, with higher gauge factor will be used to test this principle. Possible applications are sensing and tracking of adhesives between opaque surfaces, material identification, inclusion detection and counterfeit identification.

Fig. 1 Fig. 2

 12 Fig. 1 Magnetic force sensor prototype (to the left) and schematic (to the right)


Anna Castañer

M Mercè Claramunt

Equilibrium distributions and discrete Schur-constant models

Keywords: Schur-constant property, discrete stationary-excess operator, discrete equilibrium distributions. Mathematics Subject Classification (2010): 60E05, 62H05

This paper introduces Schur-constant equilibrium distribution models of dimension n for arithmetic non-negative random variables. Such a model is defined through the (several orders) equilibrium distributions of a univariate survival function. First, the bivariate case is considered and analyzed in depth, stressing the main characteristics of the Poisson case. The analysis is then extended to the multivariate case. Several properties are derived, including the implicit correlation and the distribution of the sum.

Introduction

Schur-constant models for discrete survival data have been studied by several authors including [START_REF] Castañer | Discrete Schur-constant models[END_REF], [START_REF] Lefèvre | Markov property in discrete Schur-constant models[END_REF] and [START_REF] Ta | Some properties of bivariate Schur-constant distributions[END_REF]. In this case the discrete setting means that it is valued at N 0 = {0, 1, ...}. There is also a vaste literature regarding Schurconstant models for continuous survival data valued in R + ; let us mention, among others, [START_REF] Caramellino | Dependence and aging properties of lifetimes with Schurconstant survival function[END_REF], [START_REF] Chi | Decomposition of a Schur-constant model and its applications[END_REF] and [START_REF] Nelsen | Some properties of Schur-constant survival models and their copulas[END_REF]. [START_REF] Castañer | Discrete Schur-constant models[END_REF] discuss the properties of a Schur-constant vector (X 1 , . . . , X n ) generated by the survival function of an univariate random variable (r.v.). In that paper, it is shown that this survival function has to be n-monotone. The present paper aims to develop a class of Schur-constant models generated by an admissible univariate survival function. The admissibility condition is the existence of the (n -1) -th equilibrium distribution of the univariate r.v. Because of it, this new model is called Schur-constant (multivariate) equilibrium distribution model. Some properties of the bivariate continuous case have been previously studied by [START_REF] Nair | Modelling lifetimes with bivariate Schur-constant equilibrium distributions from renewal theory[END_REF]. For the sake of completeness, we include in this introduction some definitions and well-known results.

Let (X 1 , . . . , X n ) be a vector of n (≥ 2) arithmetic non-negative random variables, called lifetimes. It is said to have a Schur-constant joint survival function if for all (x 1 , . . . , x n ) ∈ N n 0 ,

P (X 1 ≥ x 1 , . . . , X n ≥ x n ) = S(x 1 + . . . + x n ), (1.1)
where S is an admissible function from N 0 to [0, 1]. In fact, S is the survival function of each of the marginal r.v. X i and it has to be n-monotone on N 0 . A function f (x):

N 0 → R is said to be n-monotone if it satisfies (-1) j ∆ j f (x) ≥ 0, for j = 0, . . . , n, (1.2)
where ∆ is the forward difference operator (i.e. ∆f (x) = f (x + 1) -f (x)) and ∆ j is its j -th iterated.

1

A general representation [START_REF] Castañer | Discrete Schur-constant models[END_REF] valid for any discrete Schur-constant model, putting a b = 0 when a < b, is:

S(x 1 + . . . + x n ) = E Z -(x 1 + . . . + x n ) + n -1 n -1 / Z + n -1 n -1 , (1.3)
where the variable Z is distributed as X 1 + . . . + X n , i.e., with a probability mass function (p.m.f.) given by

P (Z = z) = (-1) n ∆ n S(z) z + n -1 n -1 . (1.4)
As a particular case, for n = 2, P (Z = z) = ∆ 2 S(z)(z + 1).

(1.5)

The p.m.f. of any subvector in (X 1 , . . . , X n ) can be obtained from the j -th iterated forward difference of S, P (X 1 = x 1 , . . . , X j = x j ) = (-1) j ∆ j S(x 1 + . . . + x j ).

(1.6)

For bivariate Schur-constant vectors we can obtain an easy expression of the p.m.f. of the marginal. This is the next Lemma 1.1, which will be used in the next section of the paper.

Lemma 1.1. If the vector (X 1 , X 2 ) is Schur-constant, the probability mass function of X 1 fulfills

P (X 1 = x + 1) = P (X 1 = x) - P (Z = x) x + 1 , (1.7)
for x ≥ 0, with

P (X 1 = 0) = E 1 Z + 1 .
Proof. The survival function of X 1 , for n = 2, is obtained from (1.3) when x 2 = 0,

S(x) = E (Z -x + 1) + Z + 1 = ∞ z=x z -x + 1 z + 1 P (Z = z). (1.8)
From definition, P (X 1 = x) = -∆S(x). Then, considering (1.8),

P (X 1 = x) = ∞ z=x 1 z + 1 P (Z = z).
(1.9)

From (1.9), for x = 0, P (X 1 = 0) = E 1 Z+1 . Now, let us apply again the definition of the survival function so that ∆ 2 S(x) = P (X 1 = x) -P (X 1 = x + 1). Using this expression and (1.5), the desired result (1.7) is obtained.

On the basis of these definitions, we develop our analysis. The paper is organized as follows. In Section 2, we introduce the discrete Schur-constant bivariate equilibrium distribution model and provide its characteristics. We prove that in the Poisson bivariate Schur-constant equilibrium distribution model, Z is also Poisson distributed. In Section 3, we extend the analysis to the multivariate case. The paper ends with some concluding remarks.

Schur-constant bivariate equilibrium distribution model

In renewal theory, the asymptotic distribution of the age or residual life in a renewal process is known as the univariate equilibrium distribution [START_REF] Cox | Renewal Theory[END_REF]. Let X be a non negative discrete random variable representing the lifetime of a component with finite mean µ and survival function S(x) = P (X ≥ x), x ≥ 0. The (first order) equilibrium distribution is defined through the probability mass function, as

p * (x) = S(x + 1) µ , for x ∈ N 0 , (2.1)
or, via the survival function, as

S * (x) = ∞ h=x+1 S(h) µ , for x ∈ N 0 . (2.2) Expression (2.
2) is well defined as survival function as long as it is a decreasing function and S * (0) = 1. This latter condition can be easily checked taking into account that E(X) = ∞ h=1 S(h). Much attention has been paid to the equilibrium distributions (of first and higher orders) associated with a given distribution function, but most studies are for continuous univariate random variables (see for instance [START_REF] Deshpande | Aspects of positive ageing[END_REF], [START_REF] Sunoj | Characterizations of some continuous distributions using partial moments[END_REF] and the references therein). Regarding the discrete univariate case, two papers should be mentioned: [START_REF] Li | The equilibrium distribution of counting random variables[END_REF] and [START_REF] Willmot | Equilibrium compound distributions and stop-loss moments[END_REF]. In addition, [START_REF] Lefèvre | Stationary-excess operator and convex stochastic orders[END_REF] introduce a specific discrete version of the stationary-excess operator for discrete non-negative random variables. This discrete stationary-excess operator H maps any nonnegative discrete random variable X to an associated discrete non-negative random variable X H whose probability mass function is

P (X H = x) = P (X ≥ x + 1) E(X) , x ∈ N 0 .
Notice that the equilibrium distribution can be also obtained by applying this discrete stationaryexcess operator.

In this section, we discuss a (first order) bivariate equilibrium distribution of discrete random variables. It represents in fact a family of discrete Schur-constant distributions, that we call Schurconstant bivariate equilibrium distributions. Although it has been recently analized by [START_REF] Gupta | Bivariate equilibrium distribution and association measure[END_REF] and [START_REF] Nair | Modelling lifetimes with bivariate Schur-constant equilibrium distributions from renewal theory[END_REF] in a continuous setting, as far as we know, the bivariate discrete case has not previously been evaluated. The asymptotic joint distribution of age and residual life is defined by the survival function

G(x, y) = 1 µ ∞ h=x+y+1 S(h) = S * (x + y). (2.3)
Now we prove, in the next proposition, that this survival function corresponds to a Schur-constant survival function.

Proposition 2.1. The random vector (X * , Y * ) with survival function G given by (2.3) is Schurconstant and its marginal distribution is the equilibrium distribution of X given by (2.2).

Proof. It is enough to put y = 0 in expression (2.3) to check that X * is the equilibrium distribution of X. With respect to the Schur-constancy property, from (2.3), we obtain

∆S * (x) = - S(x + 1) µ ≤ 0, (2.4) and ∆ 2 S * (x) = - S(x + 2) µ + S(x + 1) µ = 1 µ [P (X ≥ x + 1) -P (X ≥ x + 2)] = P (X = x + 1) µ ≥ 0. (2.5)
Then, from (2.4) and (2.5), S * (x) is always 2-monotone although S(x) may not.

Now we derive several useful formulas. As (X * , Y * ) is Schur-constant,

P (X * = x, Y * = y) = ∆ 2 S * (x + y) = P (X = x + y + 1) µ , (2.6)
and the p.m.f of the sum Z = X * + Y * , by (1.5) and (2.5), is

P (Z = z) = P (X = z + 1)(z + 1) µ .
(2.7)

In fact, Z turns out to be the (first) length-biased type transform of X (as defined in [START_REF] Lefèvre | On multiply monotone distributions, continuous or discrete , with applications[END_REF]. The basic characteristics of (X * , Y * ) are available in terms of the expectations of simple functions of X and this is our next proposition.

Proposition 2.2. Let (X * , Y * ) be as in Proposition 2.1. Then,

µ * = E(X * ) = 1 2 E(X 2 ) µ -1 , (2.8) V (X * ) = 4µE(X 3 ) -3 E(X 2 ) 2 -µ 2 12µ 2 .
(2.9)

Proof. Considering (2.1), we have

µ * = E(X * ) = 1 µ ∞ x=0 xS(x + 1). (2.10) It is a simple exercise that ∆ x 2 -x 2 =
x, so that we can apply summation by parts to the last sum in (2.10). This leads to

∞ x=0 xS(x + 1) = 1 2 ∞ x=0 p(x + 1)(x + 1) 2 - ∞ x=0 p(x + 1)(x + 1) = 1 2 ∞ x=0 p(x)x 2 - ∞ x=0 p(x)x .
(2.11) Substituting this in (2.10), we obtain (2.8). Let us derive now a formula for the second ordinary moment of X * ,

E((X * ) 2 ) = ∞ x=0
x 2 S(x + 1) µ .

(2.12)

We know that ∆ x 3 3 -x 2 2 + x 6 = x 2 . Considering this, summation by parts yields

E((X * ) 2 ) = 1 µ E(X 3 ) 3 - E(X 2 ) 2 + µ 6 .
(2.13) By (2.8) and (2.13), we obtain formula (2.9). [START_REF] Castañer | Discrete Schur-constant models[END_REF] presented simple formulas for the Pearson correlation coefficient, ρ, between any two variables in a Schur-constant vector. One of them allows to compute ρ in terms of the expectation and the variance of one of the marginal random variables. But in these Schur-constant discrete bivariate equilibrium distributions, it is convenient to relate ρ to the main characteristics of the original variable X, the building block of the model.

Proposition 2.3. In terms of E(X * ) and V (X * ), ρ = V (X * ) -(E(X * )) 2 -E(X * ) 2V (X * ) , (2.14)
and, in terms of the ordinary moments of X,

ρ = 2µE(X 3 ) -3 E(X 2 ) 2 + µ 2 4µE(X 3 ) -3 (E(X 2 )) 2 -µ 2 .
(2.15)

Proof. Formula (2.14) is direct from Proposition 6.2. in [START_REF] Castañer | Discrete Schur-constant models[END_REF]. Let us now establish (2.15). We could evaluate E(X * Y * ) directly using (2.6). However, a simpler method consists in using (2.14) in which E(X * ) is substituted by (2.8) and V (X * ) is substituted by (2.9). Simple algebra yields the desired formula.

The corresponding formulas (2.8), (2.9) and (2.15) for continuous random variables can be found in [START_REF] Nair | Modelling lifetimes with bivariate Schur-constant equilibrium distributions from renewal theory[END_REF]. Now, we apply these results to two specific models generated by well-known distributions, the Poisson one and the geometric one.

Poisson Schur-constant bivariate equilibrium model

Let us consider that the Schur-constant equilibrium distribution model (X 1 , X 2 ) is generated by a r.v. X that follows a Poisson distribution with parameter λ; i.e., P (X = k) = e -λ λ k k! . This model is an interesting and particular one, because of its properties. Using (2.7), it is immediate to conclude that Z = X 1 + X 2 also follows a Poisson distribution of parameter λ,

P (Z = z) = λ z+1 e -λ (z + 1) λ(z + 1)! = λ z e -λ z! , z = 0, 1, . . . . (2.16)
Proposition 2.4. Let (X 1 , X 2 ) be a Schur-constant equilibrium distribution model generated by a r.v.

X and Z = X 1 + X 2 . Then, Z = d X if and only if X follows a Poisson distribution.

Proof. It is well known from [START_REF] Panjer | Recursive evaluation of a family of compound distributions[END_REF] that the Poisson distribution is the only one that fulfills a recursive equation on its probabilities of the form

P (X = k + 1) = P (X = k) b k + 1
, k = 0, 1, 2, . . . , (2.17) being b = λ. Then, substituting (2.17) in (2.7), we obtain

P (Z = z) = P (X = z), z = 0, 1, 2, . . . . (2.18)
For the other implication, if Z is Poisson distributed with parameter λ, from (2.7), P (X = z + 1) = P (Z=z)λ z+1 , which implies that X is Poisson distributed with the same parameter.

The distribution of X 1 can be found using the property that X 1 is the equilibrium distribution of

X. Since X is Poisson distributed, S(x) = ∞ h=x e -λ λ h h!
and then

P (X 1 = x) = 1 -Γ(x+1,λ) x! λ , x = 0, 1, . . . , (2.19)
being Γ(• , •) the incomplete Gamma function. The same result can be obtained using (1.7) and Proposition 2.4. The correlation coefficient between X 1 and X 2 in this model, computed using (2.15), is ρ = -λ 6+λ .

Geometric Schur-constant bivariate equilibrium model

Let us consider a geometric r.v. X with survival function S(x) = q x . From this survival function we build a Schur-constant equilibrium vector (X 1 , X 2 ). The r.v. X 1 is the equilibrium distribution of X.

It immediately follows that S * (x) = S(x) = q x . Then, X 1 and X 2 are independent r.v. [START_REF] Castañer | Discrete Schur-constant models[END_REF]. So, the Pearson correlation coefficient ρ computed from (2.14) is equal to 0.

Schur-constant multivariate equilibrium distribution model

In this section we generalize the results of the previous section. The first order equilibrium distribution of a discrete random variable X has been used in Section 2 to define a special class of bivariate Schurconstant models. In the present section, this analysis is extended to the general multivariate case through the use of the n -th order equilibrium distribution. Let X be a non negative discrete random variable with finite mean µ and survival function S(x). The second order equilibrium distribution can be defined by its survival function

S 2 * (x) = ∞ h=x+1 S * (h) µ 1:1 , for x ∈ N 0 , (3.1)
or, alternatively, by its p.m.f

p 2 * (x) = S * (x + 1) µ 1:1 , for x ∈ N 0 , (3.2)
being µ 1:1 the first order moment of X * (the first order equilibrium distribution of X). Likewise, we can define recursively the n -th order equilibrium distribution, (3.4) being µ i:1 < ∞ the first order moment of the i -th order equilibrium distribution of X, with µ 0:1 = µ, and S (0) * (x) = S(x).

S n * (x) = ∞ h=x+1 S (n-1) * (h) µ n-1:1 , for x ∈ N 0 , (3.3) and p n * (x) = S (n-1) * (x + 1) µ n-1:1 , for x ∈ N 0 ,
Following the same procedure than in Section 2, lets us define the random vector (X 1 , . . . , X n ) with survival function

P (X 1 ≥ x 1 , . . . , X n ≥ x n ) = S (n-1) * (x 1 + . . . + x n ).
(3.5)

Our next proposition is the parallel of Proposition 2.1 in the bivariate case. Now, we prove that the survival function (3.5) corresponds to a Schur-constant survival function.

Proposition 3.1. The vector (X 1 , . . . , X n ) with survival function (3.5) is Schur-constant and the marginal X 1 is the (n -1) -th order equilibrium distribution of X ((3.3) and (3.4)).

Proof. We know from Proposition 2.1 that it is true for n = 2, because S * (x) is 2-monotone. We assume now that it is true for n -1, i.e., S (n-2) * is (n -1)-monotone. From (3.3),

(-1) n ∆ n S (n-1) * (x) = (-1) n ∆ (n-1) ∆S (n-1) * (x) = (-1) (n-1) ∆ (n-1) S (n-2) * (x + 1) µ n-2:1 . (3.6)
Since S (n-2) * is (n -1)-monotone, expression (3.6) is equal or greater than zero, which in turn implies that S (n-1) * is n-monotone. Putting x 2 = . . . = x n = 0 in (3.5), P (X 1 ≥ x 1 ) = S (n-1) * (x 1 ), and thus, the marginal distribution X 1 is the (n -1) -th order equilibrium distribution of X.

A Schur-constant vector of n components can be defined using a survival function not necessarily n-monotone, using the previous proposition; we say then that such vector follows a Schur-constant multivariate equilibrium distribution model. By exchangeability, all of the X i 's have the same ordinary moments, and the same Pearson correlation coefficient. In this model we are able to obtain all of them from the ordinary moments of the original distribution X, as long as X i is the (n -1) -th order equilibrium distribution of X. Let us first establish the following lemma.

Lemma 3.2. The set of discrete functions F (x) such that ∆F (x) = x n-1 is given by P n (x) + C, being C a constant and P n (x) a polynomial of degree n, such that P n (x) = n r=0 a r (n)x r with a 0 (n) = 0 and a r (n), for r = 1, . . . , n, that can be computed recursively,

a n (n) = 1 n , a n-h+1 (n) = -1 n -h + 1 n r=n-h+2 r n -h a r (n), h = 2, . . . , n. (3.7)
Proof. Let us consider F (x) = P n (x) + C being P n (x) = n r=0 a r (n)x r whit a 0 (n) = 0, a polynomial of degree n on x. The difference of any of these functions is exactly the same and equals Proposition 3.3. The j -th ordinary moment of the i -th order equilibrium distribution, µ i:j = E X i * j , j = 1, 2, ..., i = 1, 2, ..., is given by the recurrence formula

∆P n (x) = n r=0 a r (n)∆x r = n r=0 a r (n) r s=0 r s x s -x r = n-1 r=0 a r (n) r s=0 r + 1 s x s = n-1 s=0 x s n-1 r=s r + 1 s a r+1 (n). ( 3
µ i:j = 1 µ i-1:1 j+1 r=1 µ i-1:r a r (j + 1), (3.11)
with a r (•) as in Lemma 3.2, and µ 0:r = E(X r ).

Proof. By definition,

µ i:j = 1 µ i-1:1 ∞ x=0
x j S (i-1) * (x + 1).

(3.12)

From Lemma 3.2, ∆ j+1 r=1 a r (j + 1)x r = x j , where the coefficients a r (•) are given by expression (3.7). Therefore, the sum in (3.12) can be computed by applying summation by parts,

∞ x=0 p (i-1) * (x) j+1 r=1 a r (j + 1)x r = j+1 r=1 a r (j + 1) ∞ x=0 p (i-1) * (x)x r .
(3.13)

In [START_REF] Li | The equilibrium distribution of counting random variables[END_REF] an equivalent expression relating the n -th factorial moments is found.

We focus now on the study of the Pearson correlation coefficient, ρ, in this Schur-constant multivariate equilibrium model. We see that ρ is a function of some ordinary moments of X.

3.4. Consider the vector (X 1 , . . . , X n ) survival function (3.5). In terms of expectation and variance of the (n -1)-th order equilibrium distribution of X, the Pearson correlation coefficient is given by

ρ = V (X (n-1) * ) -E X (n-1) * 2 -E(X (n-1) * ) 2V (X (n-1) * ) . (3.14)
Proof. It is direct from Proposition 6.2. in [START_REF] Castañer | Discrete Schur-constant models[END_REF].

Thanks to (3.11), an equivalent expression of ρ as a function of the ordinary moments of X can be deduced. For instance, let us make explicit the formula when n = 3 and n = 4. Let (X 1 , X 2 , X 3 ) be a Schur-constant multivariate equilibrium distribution model defined by X. The Pearson correlation coefficient can be computed as

ρ = 1 2 - µ -E(X 3 ) 2µ -3E(X 2 ) + E(X 3 ) 2µ 2 + 2(E(X 3 )) 2 + 3E(X 2 )(E(X 2 ) -E(X 4 )) + µ (-3E(X 2 ) -4E(X 3 ) + 3E(X 4 ))
.

If we consider now the Schur-constant multivariate equilibrium distribution model (X 1 , X 2 , X 3 , X 4 ),

ρ = 1 2 - 5 6µ -11E(X 2 ) + 6E(X 3 ) -E(X 4 ) 2µ -E(X 2 ) -2E(X 3 ) + E(X 4 ) 2(D -4µJ) , being D = 36µ 2 +65(E(X 2 )) 2 +20(E(X 3 )) 2 -70E(X 2 )E(X 4 )+5(E(X 4 )) 2 +24E(X 2 )E(X 5 )8E(X 3 )E(X 5 ), and J = 21E(X 2 ) + 8E(X 3 )-15E(X 4 ) + 4E(X 5 ).
As an example, Tables 2 and3 include the values of the Pearson correlation coefficient in a Schurconstant multivariate equilibrium distribution model built from a Poisson distributed random variable of mean λ, as a function of n, the number of elements of the random vector. One interesting property of these models is that the probability mass function of Z = X 1 +. . .+X n is easily calculated from the survival function of X. This is our next proposition. Proof. From Proposition 3.1, the survival function of X 1 is S (n-1) * (x), and from its definition (3.3),

∆ n S (n-1) * (x) = (-1) n P (X = x + n -1) µ • µ 1:1

Conclusions

Discrete Schur-constant models (of dimension n) can be generated by univariate survival functions that are n-monotone. This paper introduces a class of Schur-constant models that are generated by a wider class of univariate survival functions, i.e., those survival functions such that their (n -1) -th equilibrium distribution exists.

  we observe that ∆P n (x) is a new polynomial of degree n -1 and coefficients b s , that can be derived from the coefficients of the original polynomial, a r (n), and (3.9) and imposing the condition that ∆P n (x) = x n-1 , that is b n-1 = 1, and b s = 0 for s = 0, . . . , n -2, a linear system of n equations on a r (n) is defined. The value of a n = 1 n is obtained when s = n -1. When s = n -h, the corresponding equation is 0 a n-h+1 (n) from (3.10) the statement of the lemma is obtained.

Proposition 3. 5 .

 5 Consider the vector (X 1 , . . . , X n ) with survival function (3.5). The p.m.f. of Z = X 1 + . . . + X n isP (Z = z) = P (X = z + n -1) µ • µ 1:1 • • • µ (n-2

Table 1

 1 below gives the values of a r (n), for the first values of n.

Table 1 :

 1 Some coefficients a r (n)

	❍ ❍ ❍ ❍ ❍ ❍ r n	2	3	4	5	6	7	8	9	10
	1	-1/2	1/6	0 -1/30	0 1/42	0 -1/30	0
	2	1/2 -1/2	1/4	0 -1/12	0	1/12	0 -3/20
	3		1/3 -1/2	1/3	0 -1/6	0	2/9	0
	4			1/4	-1/2	5/12	0 -7/24	0	1/2
	5				1/5	-1/2	1/2	0 -7/15	0
	6					1/6 -1/2	7/12	0 -7/10
	7						1/7	-1/2	2/3	0
	8							1/8	-1/2	3/4
	9								1/9	-1/2
	10									1/10

Table 2 :

 2 ρ in a Poisson Schur-constant multivariate equilibrium model (X 1 , . . . , X n )

	n	2	3	4	5
	ρ	-λ 6 + λ	-λ 12 + 2λ	-λ 20 + 3λ	-λ 30 + 4λ

Table 3 :

 3 ρ in a Poisson Schur-constant multivariate equilibrium model (X 1 , . . . , X n ) for different values of λ

	❍ ❍ ❍ ❍ ❍ ❍