Lionel Lacassagne

Daniel Etiemble

Ali Hassan Zahraee

Alain Dominguez
email: alain.dominguez@intel.com

Pascal Vezolle
email: pascal.vezolle@fr.ibm.com

High Level Transforms for SIMD and low-level computer vision algorithms

Keywords: High Level Transforms, SIMD, Intel SSE & XeonPhi, IBM Altivec, ARM Neon, code optimization, 2D stencil, low-level computer vision and image processing algorithms

This paper presents a review of algorithmic transforms called High Level Transforms for IBM, Intel and ARM SIMD multi-core processors to accelerate the implementation of low level image processing algorithms. We show that these optimizations provide a significant acceleration. A first evaluation of 512-bit SIMD Xeon-Phi is also presented. We focus on the point that the combination of optimizations leading to the best execution time cannot be predicted, and thus, systematic benchmarking is mandatory. Once the best configuration is found for each architecture, a comparison of these performances is presented. The Harris points detection operator is selected as being representative of low level image processing and computer vision algorithms. Being composed of five convolutions, it is more complex than a simple filter and enables more opportunities to combine optimizations. The presented work can scale across a wide range of codes using 2D stencils and convolutions.

Introduction

Graphic Processing Units (GPU) are efficient for High Performance Computing (HPC) [START_REF] Volkov | qr and cholesky factorizations using vector capabilities of gpus[END_REF] where the operations involved lend themselves to massive parallelization [START_REF] Garland | Parallel computing experiences with cuda[END_REF]. Some papers claim "orders-ofmagnitude performance increase" versus General Purpose Processors (GPP). Recently, papers from Intel [START_REF] Lee | Debunking the 100x gpu vs. cpu myth: an evaluation of throughput computing on cpu and gpu[END_REF] and IBM [START_REF] Bordawekar | Can cpus match gpus on performance with productivity?[END_REF] claim that GPPs can match GPUs if optimizations are applied. Some papers propose a fair benchmarking, by optimizing as much as possible the implementations on these architectures and compare them rigorously to find out the applications that really achieve significant speedups (like n-body [START_REF] Arora | Direct n-body kernels for multicore platforms[END_REF] or stencil [START_REF] Datta | Stencil computation optimization and auto-tuning on state-of-the-art multicore architectures[END_REF]).

The aim of this paper is to present some high level transforms (HLT) for SIMD applied to low-level image processing and computer vision algorithms. We focus on the fact that combining SIMD with OpenMP is not enough to reach and sustain a high level of performance. Optimizing memory accesses is an issue and HLT should be cache-aware and also external-memory-aware when data do not fit in the caches. In our study, the Harris operator [START_REF] Harris | A combined corner and edge detector[END_REF] for point of interest detection is chosen. Widely used for image stabilization, velocity analysis or visual tracking, this operator is also a representative example of the regular low-level image processing algorithms class. As it is composed of 8 operators (Fig. 1), it enables more opportunities for optimizations and parallelization than a unique convolution kernel. HLT can be also efficiently applied to any code using 2D stencils.

The paper is organized as follows: the first section details the Harris operator, the software optimizations and HLT that can be applied to the Harris detector. The second section presents the targeted SIMD machines (with Altivec, SSE or Neon instruction set extensions) and a multistep benchmark that evaluates the impact of these optimizations. A first evaluation of XeonPhi is also presented. Then performance of GPPs and GPUs is compared. The Harris detector computation relies on a set of point-topoint operations such as products and additions, along with (3 × 3) convolutions (the Sobel gradients and Gaussian smoothing). Given an input image I, the first derivatives Ix and Iy of Sobel gradients are computed by GradX and GradY. The cross-products Ixx = Ix × Ix, Ixy = Ix × Iy and Iyy = Iy × Iy performed by the Mul operator. These cross-products are smoothed by a Gaussian filter (Gauss) and finally, "coarsity" K is a linear combination of the smoothed cross-products. A producer-consumer model is also provided under each figure with explicit input and output patterns: (3 × 3) → (1 × 1) for Grad and Gauss and (1 × 1) → (1 × 1) for Mul and coarsity. Moreover, considering the four stages of its classical computation, the so called Nopipe version (Fig. 1) is an interesting candidate for higher level optimizations such as fusion of operators.

Harris detector and HLT

Optimizing a code is a two-step process. First, the code is optimized to reach the highest level of performance in a stressless context, when data fit in the cache. Next, performance is sustained in a stressed context when data do not fit in the cache. We first describe the classical compilation optimizations like loop-unrolling, scalarization and reduction for High Performance Computation. These techniques are then adapted to the low level operators in image processing: some high level algorithmic transformations named Halfpipe and Fullpipe are introduced. These optimizations are then improved by considering the application domain: both the separability of the 2D filters and the convolutions overlapping are considered and combined to loop-unrolling and scalarization to add reduction in Halfpipe and Fullpipe transforms. Their complexity is finally evaluated in terms of arithmetic operations and memory access. The second part presents an advanced memory layout transform, circular buffers and modular addressing associated with another level of operator pipelining to optimize spatial and temporal data locality. Finally SIMD reduction is detailed.

Classical transformations, decomposition and reduction

B3×3 = 1 16 " 1 2 1 2 4 2 1 2 1 # = 1 16
" 1 2 1 # * [ 1 2 1 ] (1)
Algorithm 1: 1-pass implementation of the 3 × 3 binomial filter with the 2D-filter corresponding to equation ( 1)

for i = 1 to n -1 do 1 for j = 1 to n -1 step 3 do 2 a 0 ← X(i -1, j -1), b 0 ← X(i -1, j), c 0 ← 3 X(i -1, j + 1) a 1 ← X(i, j -1), b 1 ← X(i, j), c 1 ← X(i, j + 1) 4 a 2 ← X(i + 1, j -1), b 2 ← X(i + 1, j), c 2 ← 5 X(i + 1, j + 1) s ← 1a 0 +2b 0 +1c 0 +2a 1 +4b 1 +2c 1 +1a 2 +2b 2 +1c 2 6 Y (i, j) ← s/16 7
To facilitate the presentation of different algorithms, we assume that computations use local variables and that memory accesses correspond to transfers between image pixels and local variables. This assumption fits with RISC processors, for which the computations operate on register values and memory accesses correspond to LOAD and STORE instructions. In that case, local variables are held in processor registers. The situation is different for CISC processors using the IA-32 or Intel 64/AMD64 instruction sets for which memory operands are used by computing instructions. In that case, due to the limited numbers of registers, some local variables will be held in the memory hierarchy. However, this difference does not change the impact of the optimizations that are presented. In the rest of the paper, we present the transformations assuming a RISC instruction set.

The classical software optimizations for optimizing compilers [START_REF] Allen | Optimizing compilers for modern architectures: a dependence-based approach[END_REF] aim to improve the operation of the processor pipeline. Given a 3 × 3 binomial filter aka Gaussian filter (Eq. 1), one can perform scalarization (to put data into registers), and register rotation to avoid reloading data from an iteration of the filter to another. The 9 LOADs of the original algorithm (Alg. 1) are replaced by only 3 (Alg. 2). It is an efficient optimization, considering that many algorithms are memory-bound, but the algorithm complexity remains the same.

In this presentation of optimization techniques, we assume general filters with unknown coefficients: this is why we indicate multiplications by 1, 2 or 4, that are transformed into a set of additions (4x = t+t with t = x+x, i.e. the strength reduction optimization) or replaced by a shift (4x = x<<2) for integer computations. Note also that the apron processing problem is not considered in order to simplify the explanation.

Algorithm 2: 1-pass implementation of the 3 × 3 binomial filter with one 2D-filter Register Rotation

for i = 1 to n -1 do 1 j ← 1 2 a 0 ← X(i -1, j -1), b 0 ← X(i -1, j) 3 a 1 ← X(i, j -1), b 1 ← X(i, j) 4 a 2 ← X(i + 1, j -1), b 2 ← X(i + 1, j) 5 for j = 1 to n -1 step 3 do 6 c 0 ← X(i -1, j + 1) 7 c 1 ← X(i, j + 1) 8 c 2 ← X(i + 1, j + 1) 9 s ← 1a 0 +2b 0 +1c 0 +2a 1 +4b 1 +2c 1 +1a 2 +2b 2 +1c 2 10 Y (i, j) ← s/16 11 a 0 ← b 0 , b 0 ← c 0 // Rot 12 a 1 ← b 1 , b 1 ← c 1 // Rot 13 a 2 ← b 2 , b 2 ← c 2 // Rot 14
Taking into account the application domain, the 2D-filter can be replaced by two 1D-filters (Eq. 1, right part). This algorithmic transformation reduces complexity and number of memory accesses, but requires two passes on the image, which can generates cache misses, when the image is too large to entirely fit in the cache.

So we need to introduce another optimization to combine the two 1D-filters with a single pass, to factor the computations and reduce the number of memory accesses simultaneously. First, the result of the first 1D-filter is stored in a register (Alg. 3, line 5, 6 and 11). This transformation is called a reduction. In our case, it is a column-wise reduction. Then the second 1D-filter is directly applied to the reduced values (Alg. 3, line 12). Algorithm 3: 1-pass implementation of the 3 × 3 binomial filter with two 1D-filters, with Register Rotation and reduction

for i = 1 to n -1 do 1 a 0 ← X(i -1, j -1), b 0 ← X(i -1, j) 2 a 1 ← X(i, j -1), b 1 ← X(i, j) 3 a 2 ← X(i + 1, j -1), b 2 ← X(i + 1, j) 4 ra ← 1a 0 + 2a 1 + 1a 2 // Red part #1 5 r b ← 1b 0 + 2b 1 + 1b 2 // Red part #1 6 for j = 1 to n -1 do 7 c 0 ← X(i -1, j + 1) 8 c 1 ← X(i, j + 1) 9 c 2 ← X(i + 1, j + 1) 10 rc ← 1c 0 + 2c 1 + 1c 2 // Red part #1 11 s ← 1ra + 2r b + 1rc // Red part #2 12 Y (i, j + 0) ← s/16 13 ra ← r b , r b ← rc // Rot of reduced values 14
For a general 1D k-tap filter, there are still (k -1) copies for register rotation. They can be entirely removed by loop unrolling. The order of the unrolling is usually chosen by the compiler, which may lead to suboptimal unroll. But in signal or image processing, the unrolling order has not to be chosen by heuristics: the optimal order leading to a perfect unroll is equal to the filter order. In our case, the smallest unrolling order is k = 3.

Operator fusion: the Halfpipe and Fullpipe transforms

The fusion of several operators is much more than a simple loop fusion (aiming to improve data locality). By pipelining operators and storing first results in registers instead of memory, the transform avoids any intermediate memory access. For that aim, each operator is described by the producer-consumer model with a consumption pattern and a production pattern. Such a model is derived from Synchronous Data Flow [START_REF] Lee | Multidimensional streams rooted in dataflow[END_REF]. The only condition needed to fuse two operators (in the sense of mathematical composition of functions f • g) is that the patterns must be either similar or adaptive: the output pattern of the first operator should be the same as the input pattern of the second one. In its Nopipe version (Fig. 1), the Harris detector is composed of four computation stages: computation of the gradients, product of the gradients, smoothing of the products and computation of the point coarsity. This computation chain requires the access to eight intermediate arrays. It is therefore possible to pipeline the Sobel and Mul operators on the one hand, and the Gaussian and coarsity operators on the other hand, since their consumption and production patterns are compatible. The transformation is called Halfpipe2 (Fig. 2, top) and reduces the number of memory accesses (Tab. 1).

It is possible to entirely pipeline the operators and suppress any intermediate memory access, resulting in the Fullpipe version (Fig. 2, bottom). In that case, it is necessary to adapt the consumption and production patterns: to produce (1 × 1) points in the output, it is necessary to consume (5 × 5) points and to perform the operations Sobel+Mul (3 × 3) times. A reduction along the columns is applied to optimize both the number of memory accesses and the number of computations.

Note that the Fullpipe version is much more complex than the other ones with up to 3 times more computations compared to the Nopipe version (Tab. 1). Indeed, without any intermediate image to store the common computations (which were performed 9 times by the Gaussian operators), these computations have to be reperformed.

The Halfpipe1 version is designed to balance Halfpipe2 and Fullpipe. Figure 2 shows that Halfpipe2 has the lowest complexity, and the highest number of memory access. But once optimized (with reduction), Halfpipe2 and Halfpipe1 versions have similar complexities (Tab. 1) and memory accesses, while Fullpipe has a far lower number of memory accesses. If we consider the arithmetic intensity AI (ratio between computations and memory accesses) as a metric of performance, Fullpipe is particularly interesting from a parallelization perspective since its arithmetic intensity is higher than Halfpipe and Nopipe versions. Fullpipe version is compute bound while the two other versions are memory bound with different levels of memory stress.

Multithreading & cache overflow

For a multi-core SIMD processor, maximum performance is indeed reached when all cores are running. That means that evaluation of HLT for SIMD should be done when the code is multi-threaded. This corresponds to the configuration that maximizes the stress on the external memory bus. OpenMP is used for this purpose. The necessary code modification is very light. It mainly consists of pragmas that are simply put before the code sections to be parallelized. Considering a SPMD parallelization with a sub-band calculation, the modifications consist in parallelizing the outer loops, and privatizing variables to avoid serialization of accesses. The parallelization of a for loop only adds the directive (#pragma omp parallel for) to the sequential code. Such parallelization may have a super-linear speedup. Let us introduce the measurement unit used for the analysis and the benchmarking: cpp (cycle per pixel). It is the number of clock cycles normalized by the number of pixels processed.

Figure 3 presents the cpp execution times of our benchmarks according to the size of the working set for different configurations of processors and caches. Let us first consider (Fig. 3a). On the horizontal axis, a point corresponds to the size of the working set, i.e. to the size of the data to be processed. As cpp is the number of clocks per pixel, a horizontal cpp curve means that the execution time, including computation and memory accesses, is exactly proportional to the data size. As a matter of fact, the actual cpp curve shows two plateaux and a transition area from one plateau to the other. This transition corresponds to the switch between a state for which the working set can be held in the cache to another one for which the working set is larger than the cache size. We call this phenomenon a cache overflow. This phenomenon is widely detailed in [START_REF] Drepper | what every programmer should know about memory[END_REF].

Using more processors reduces execution time, even when the Amdhal law forbids a perfect speedup (Fig. 3c). Adding caches (or using a larger cache) would shift right the cpp curve and results are improved (Fig. 3b). For multi-core processors, the impact resulting from the combination of these two modifications is clearly visible (Fig. 3d). Multithreading with multi-core processors increases the available cache size and can be seen as a mean to postpone cache overflow. The performance decrease is more important for larger working sets. Note that the memory optimizations can lead to a super-linear speedup. For a given working set, it appears when the basic version (a) does not fit in the cache while the optimized one (d) does.

Data interleaving memory optimization

The data interleaving consists in replacing the accesses to several different arrays by accesses to one single array containing the same data [START_REF] Truong | Improving cache behavior of dynamically allocated data structures[END_REF]. This is the SoA-AoS transformation (Structure of Arrays vs Array of Structures). For the Harris detector, the arrays that are produced at the same computation stages are interleaved, namely Ix et Iy for the Sobel operator, Ixx, Ixy and Iyy for Mul and Sxx, Sxy et Syy for the Gaussian operators. This optimization is particularly efficient when the number of arrays to manipulate is lower than the associativity of the memory caches, therefore avoiding systematic cache misses. For scalar computations, the interleaving is 1:1 and for SIMD the interleaving is 4:4 (that is 1 vector of 4 points for 1 vector of 4 points) to avoid adding SIMD instructions to deinterleave the data. For example we have a full SIMD register of Ix then a full SIMD register of Iy.

Circular buffers with modular addressing optimization

If array interleaving helps avoid systematic cache eviction, cache overflow for large data sets can not be avoided. The circular buffers with modular addressing optimization are designed to address this problem.

Circular buffers are widely used in embedded applications (typically signal processing on DSPs) since they reduce the memory required to store intermediate results and improve the performance Circular buffers (CB) can be transparently introduced into user code thanks to Iliffe pointer (offset addressing) [START_REF] Iliffe | The use of the genie system in numerical calculation[END_REF] popularized by Numerical Recipes in C matrix [START_REF] Press | Numerical Recipes in C book set: Numerical Recipes in C: The Art of Scientific Computing[END_REF] (Fig. 4, left). So, for a given point of coordinates (i, j), T[i] is pointing to the row (i mod k) instead of (i) and j is the offset within the row (Fig. 4, center). Such a spatial locality optimization can also be used in a multithreaded context using OpenMP, but with some restrictions. First the number of threads should be known before the memory allocation, because the number of sets of circular buffers should be equal to p, the number of created threads. Hence, the number of threads cannot be dynamically changed (except if a buffer reallocation is performed). Second, the loop space should be split into p blocks of contiguous lines: OpenMP parallel macro should be tuned by specifying the thread policy with schedule static. Finally, p sets of k circular buffers should be allocated for an execution on p processors with a k × k convolution. Because of the overlapping of addressing due to the convolution kernel, the sets of circular buffers do not define a bijective application. On the right part of figure 4, the lines addressing the internal aprons (in dark gray) are duplicated while the external aprons (in light gray) are pointing to the first or last line of the block.

In the Halfpipe case (Halfpipe1 and Halfpipe2), rather than applying the Sobel and Mul operators to the entire image, this algorithm is first executed on the first two rows of the image (prolog). Once arrived at the third row, the first two operators have produced enough data so that the two subsequent operators can start their own work. Therefore, the two groups of two operators are chained together. If the Halfpipe transformation can be viewed as the pipeline of two operators, the use of circular buffers with modular addressing can be viewed as the pipeline of two tiles of data with cacheblocking. In the following, we call this optimization mod.

HLT adaptation to SIMD computations

Optimizing a convolution with SIMD instructions relies on the optimization of memory accesses and the management of unaligned vectors. For example, for the 1D horizontal Gaussian filter, X(i -1) + 2X(i) + X(i + 1) implies the access to the left unaligned and right unaligned data. The naive way is to perform an unaligned load using mm loadu ps instruction when available (only on SSE, not on Altivec nor Neon). It is not efficient as 3 LOADs are required and cannot be combined with register rotation. The efficient way is to build unaligned vectors from aligned vectors (Alg. 4, lines 10, and vec right encapsulate the instructions: vec sld for Altivec, vextq f32 for Neon and mm shuffle ps for SSE and SSE2 extension or mm alignr epi32 for SSSE3+ extension.

Algorithm 4: 1-pass implementation of the 3 × 3 binomial filter with one 2D-filter with Register Rotation

for i = 1 to n -1 do 1 j ← 1 2 a 0 ← X(i -1, j -1), b 0 ← X(i -1, j) 3 a 1 ← X(i, j -1), b 1 ← X(i, j) 4 a 2 ← X(i + 1, j -1), b 2 ← X(i + 1, j) 5 for j = 1 to n -1 do 6 c 0 ← X(i -1, j + 1) 7 c 1 ← X(i, j + 1) 8 c 2 ← X(i + 1, j + 1) 9 a 0 ← vec left(a 0 , b 0 ), c 0 ← vec right(b 0 , c 0 ) 10 a 1 ← vec left(a 1 , b 1 ), c 1 ← vec right(b 1 , c 1 ) 11 a 2 ← vec left(a 2 , b 2 ), c 2 ← vec right(b 2 , c 2 ) 12 s ← 1a 0 +2b 0 +1c 0 +2a 1 +4b 1 +2c 1 +1a 2 +2b 2 +1c 2 13 Y (i, j) ← s/16 14 a 0 ← b 0 , b 0 ← c 0 // Rot of 1st line 15 a 1 ← b 1 , b 1 ← c 1 // Rot of 2nd line 16 a 2 ← b 2 , b 2 ← c 2 // Rot of 3rd line 17
Such a technique is adaptable to the scalar column-wise reduction optimization (Alg. 3). Instead of using the load-permute-add sequence, we can swap the permute (4, lines 10, 11 and 12) and add steps: reductions are computed on aligned data for the vertical 1D filter, then the unaligned vectors are built (Alg. 5, lines 12 and 13). As the vec left and vec right instructions are applied to reduced registers, fewer are needed: from 6 in the version without reduction down to 2 in the version with reduction. That computation scheme makes the SIMD reduction possible and very efficient as it can be combined with register rotation or loop-unrolling.

Algorithm 5: 1-pass implementation of the 3 × 3 binomial filter with two 1D-filter, with Register Rotation and reduction Three SIMD extensions are evaluated : SSE to SSE4 for Intel, Altivec (VMX) for IBM and Neon for ARM (Tab. 2). For each architecture the company's compiler is used : icc, xlc and armcc. To simplify the benchmark, and also allow a fair comparison with Neon, only the 128-bit SIMD extension is evaluated. Even with such a limitation, the hardware of 256-bit AVX Sandy/IvyBride processors can pair two identical 128-bit instructions into a 256-bit instruction. The second restriction is that FMA (Fused Multiply-Add) is not used, as this instruction is not present on all architectures. But the main point limiting the performance is the bandwidth, not the lack of FMA. In order to get a multi-architecture code, both kernel computations and load/store instructions are rewritten with macros. A header file holds the translation for each SIMD instruction set. This feature ensures that the same algorithm is executed on all architectures. We plan to rewrite the code with Boost.SIMD soon [START_REF] Estérie | Boost. simd: generic programming for portable simdization[END_REF]. Let us first compare the impact of HLT without and with reduction and focus on Penryn and Nehalem processors. Before cache overflow (Tab. 3), Penryn and Nehalem behavior and cpp are very similar: red provides a speedup in the range of [×1.3 : ×2.3] and the HLT total gains (Nopipe versus Halpipe1 are also close: ×2.5 and ×2.9. After cache overflow, red is inefficient for Halfpipe versions. As Fullpipe version is memory bound, the cpp values remain unchanged and red oprimization still provides a speedup of ×1.4. Because of the difference in memory bus (FSB for Penryn, QPI for Nehalem) the impact of cache overflow on Nopipe version is very different: the Penryn slowdown is 16.6 whereas Nehalem one is only 3.1. The total HLT gain increases a little for Nehalem : ×3.6 and a lot for Penryn ×14.9 ! One important point is that the best transform changes with the data set size. Before cache overflow the best transform is the Halfpipe1 and after, it is the Fullpipe version. This result is different for IBM and ARM processors where Halfpipe is always the best optimization (see next paragraph).

for i = 1 to n -1 do 1 a 0 ← X(i -1, j -1), b 0 ← X(i -1, j) 2 a 1 ← X(i, j -1), b 1 ← X(i, j) 3 a 2 ← X(i + 1, j -1), b 2 ← X(i + 1, j) 4 ra ← 1a 0 + 2a 1 + 1a 2 5 r b ← 1b 0 + 2b 1 + 1b 2 6 for j = 1 to n -1 do 7 c 0 ← X(i -1, j + 1) 8 c 1 ← X(i, j + 1) 9 c 2 ← X(i + 1, j + 1)

Impact of vectorization

Concerning the vectorization, C99 source code is used with the restrict qualifier to enforce the fact that there is no pointer aliasing. For 16, 32 or 64 byte alignments, mm malloc is used for Intel, posix memalign for ARM. For IBM, the alignment of classical malloc is the size of Altivec registers (when activated). Under these conditions, the relative vectorization speedup goes from 0% (cannot vectorize due to data dependency) for red versions up to 90% (quite perfect vectorization) for the version without reduction. In fact, column-wise reduction is too restrictive for the compiler. Regarding the vectorization capability of compilers, it appears that only HTL without red optimization is vectorized with a performance very close to the hand-coded SIMD version. If data do not fit in the cache, there is no need to apply red optimization nor to SIMDize the code. Otherwise, code SIMDization combined with HLT with reduction is efficient. So the "ninja gap" [START_REF] Satish | Can traditional programming bridge the ninja performance gap for parallel computing applications[END_REF] between the basic version (vectorized and parallelized by the compiler -we assume efficient compiler parallelization) and the best handcrafted version is equal to the total HLT gain here. The impact of data interleaving is not detailed in the paper. The speedup is about 20% when there is no reduction and falls to a few percents when reduction is active. The impact of manual loop-unrolling compared to unrolling by the compiler is not evaluated either. Compiler unrolling option is set for all compilers and architectures.

Table 5. Xeon Phi performance for HLT+red optimizations

The codes have been easily ported on the Xeon Phi: Compilation tools and computing model are same as classic processors. As our codes use macro to handle 3 SIMD dialects, the only modification was to replace mm intrinsics prefix by mm512 , and set the SIMD cardinal to 16 instead of 4. FMA is not used in order to be consistent with other processors. Because its arithmetic intensity is very high, Halfpipe1+red is faster than Halfpipe2+red. For the same reason, Fullpipe+red is faster to Halfpipe1+red and is actually the fastest version. The mod transform is currently inefficient. More investigations are required with Vtune to understand the problem and fix it. Anyway compared to one IvyBridge Xeon, the Phi has a better cpp and -due to CPU frequency differencean execution time that is 25% longer. When the bandwidth (Tab. 5) exceeds stream triad performance, it is not an error but the evidence that some data still remain in the caches.

Conclusion for SIMD processors

As previously said, the first step to accelerate a code is to optimize when data fit in the cache. This is the goal of HLT with red optimization. While the second step is to sustain the performance after cache overflow: that is the goal of the additional mod optimization. Table 4 provides cpp for HLT+red optmization before overflow and HLT+red+mod after overflow. In fact mod only provides an additional speedup of ×2, but it prevents performance deceleration. Without it, the computing performance delivered for big images (2048 × 2048) by the evaluated parallel machines is only 16.3 % down to 1.1 % of the performance delivered for small images (512 × 512). HLT are mandatory.

It appears that Intel IvyBridge and IBM Power7 have very close cpp performance, despite the architecture differences: 2 × 12 cores versus 4 × 8, 256-bit AVX versus 128-bit Altivec. Power7 is faster with red optimization and IvyBridge is faster for mod optimization. From an execution point of view they have the same execution time 0.23 ms. We can also notice the impact of HLT transforms: the higher AI the high impact. Even if latest processors have a better memory architecture -that will decrease the impact of HLT compared to previous processor -it is still very important: ×7.7 for Power7 and ×35.3 for IvyBridge. For Cortex, impact should increase with latest quad-core processor. In order to make a fair comparison, we also applied these HLT to Nvidia GPU with CUDA 4. There are three kinds of memory: Global memory that is shared by all threads of all SM processors, texture memory that provides hardware bilinear interpolation and shared memory that is private to a SM processor. Texture usage reduces both the complexity and the number of memory accesses of the convolution. For example X

Multi-core SIMD GPP versus GPU

[i-1]+2*X[i]+X[i+1] is equal to 0.5*(X[i-0.5]+X[i+0.5]
). 2D-convolutions only require 4 LOADs instead of 9, and 3 ADDs instead of 8 ADDs and 5 MULs. For shared memory, the Nvidia example "convolution2D" is specialized to a 3 × 3 convolution wit loop unwinding. As GPU are very sensitive to tile size, an exhaustive search of the best tile has been done for all benchmarked GPU (GTX 285, GTX 480, GTX 580, Quadro 4000). Usually 16 × 8 is a nice size for many algorithms, but the best one provide a boost up to 50 %.

Profiling helps to explain the Shared performance: the occupancy is 25 % and there are 41 registers per thread. The reason is that the GPU uses a lot of Shared memory. Two other metrics are interesting: the ipc (instruction per cycle) equals 1.85 (for a max of 2) and the L1 hit rate is 99.6 %. That clearly means that we are very close to peak performance, with a quite optimal pipeline feed. Note that the other tested configurations have a better occupancy, but are slower. So, we are confident that this configuration is globally optimal, and that "Fermi has better performance at lower occupancy" [START_REF] Volkov | Better performance at lower occupancy[END_REF] (if there is enough computations). We can see (Tab. 6) that HLT transforms are also efficient for GPU: ×5.6. Moreover, the Fullpipe versions (with texture memory and interpolation, or with shared memory) are faster than Halfpipe (version that requires two kernels instead of one for the Fullpipe). That focuses on the importance of avoiding communication and synchronization on GPU (Halfpipe requires two syncthreads instead of one for Fullpipe).

To perform a State-of-the Art comparison (Tab. 6), GTX Titan and K40 performance have been estimated according to the clock frequency and cores numbers. We do not take into account PCI transfer duration for Xeon Phi nor for GPU. If we compare the performance without any optimization, (Nopipe column), the GPUs are faster than the GPPs. But if we compare the most optimized version, GPPs match GPUs performance thanks to mod optimization. Only Xeon Phi can comes close to GPUs without the mod optimization.

Main conclusion and future works

In this paper, we have presented algorithmic High Level Transforms for SIMD General Purpose Processors applied to low-level computer vision algorithms. HLT have a major impact on performance, and make the difference. They can be applied to any code using 2D stencils or convolution and so, can scale across a wide range of codes

The combination of operator fusion/pipelining with algorithmic reduction and circular buffers with modular addressing provide huge speedups: from ×6.1 for dual Cortex-A15 up to ×89.3 for Penryn (that has a small bandwidth). For State-of-the-Art processors with a more important bandwidth like IvyBridge and Power7, the respective speedups are ×35.3 and ×7.7. Xeon Phi is easy to program thanks to SIMD intrinsics and OpenMP and match other processors. These benchmark optimizations have clarified several important points. First of all, SIMDization is really effective on the SIMD multi-core machines and is the only way to match GPU performance. Secondly, considering the fusion of operators, the Halfpipe and Fullpipe transforms provide consequent additional gains. Since these transformations remain out of reach even for the best current compilers [START_REF] Pouchet | Iterative optimization in the polyhedral model: part ii, multidimensional time[END_REF], the manual coding is fully justified. Third, the algorithmic reduction the use of modular addressing is beyond the scope of compilers as such transforms modify the code semantic and should be also handed-coded.

Depending on the user skills and the required level of performance, there are finally two choices. On one hand, the cost-effective implementation is to combine the compiler auto-vectorization with openMP. It may be sufficient for some applications with soft realtime constraints -especially if data fit in the caches. On another hand, the most-effective version is to combine SIMDization with reduction and modular addressing. The circular buffers enforce spatial and temporal locality that are mandatory to sustain performance for processors that have a cache overflow. It is a significant modification, but it is a worthwhile modification as it ensures a high and quite constant level of performance, without cache overflow for any realistic image size.

In future works, we will finalize the on-going benchmarks on Intel Xeon Phi and Nvidia Kepler GPUs. The general-purpose code will be specialized to take into account the specificities of processors (larger SIMD, FMA, for GPP, CUDA 5 and communica-tions for GPU) in order to provide more accurate results. More important, we also plan to implement the algorithmic reduction into a DSL (Domain Specific Language) above Boost.SIMD for signal/image/stencil operators. Concerning circular buffers and modular addressing, the memory layout modification should be also embedded into a tool that capture such abstraction to combine it with the semantic of image and multidimensional stencil operators.

Figure 1 .

 1 Figure 1. Harris detector, Nopipe version

Figure 2 .

 2 Figure 2. Halfpipe and Fullpipe transforms of Harris operator, from the lowest to the highest complexity

Figure 3 .

 3 Figure 3. Expected transformations impact: evolution of the operating point according to the number & size of caches, cache overflow (in gray) and number of processors

Figure 4 .

 4 Figure 4. Iliffe matrix with 0, 1 or 2 sets of three circular buffers for an image made of 6 rows and a border of 1, for a 3 × 3 convolution.

10 rc ← 1c 0 + 2c 1 + 1c 2 11 r 12 s 13 Y

 10111213 a ← vec left(ra, r b ), r c ← vec right(r b , rc) ← 1r a + 2r b + 1r c (i, j) ← 14 ra ← r b , r b ← rc // Rot of reduced values 15 3. Benchmarking 3.1 Targeted SIMD processors

Table 1 .

 1 Algorithmic

	version	MUL + ADD	LOAD + STORE	AI
		without reduction
	Nopipe	5 + 44 = 49	48 + 9 = 57	0.9
	Halfpipe 2	5 + 44 = 49	36 + 4 = 40	1.2
	Halfpipe 1	29 + 44 = 73	27 + 3 = 30	2.4
	Fullpipe	29 + 124 = 153	25 + 1 = 26	5.9
		with reduction
	Nopipe+red	5 + 27 = 32	21 + 9 = 30	1.1
	Halfpipe 2 +red	5 + 27 = 32	12 + 4 = 16	2.0
	Halfpipe 1 +red	11 + 27 = 38	9 + 3 = 12	3.2
	Fullpipe+red	29 + 82 = 111	5 + 1 = 6	18.5

complexity, number of memory accesses and arithmetic intensity for Harris with HLT

Table 2 .

 2 Main characteristics of the evaluated machines.

	processor	nb cores	freq (GHz)	perf. GFlops	BW GB/s	AI ratio
	Cortex A9 OMAP4	1 × 2	1.2	4.8	1.2	4.0
	Cortex A15 Exynos5	1 × 2	1.7	13.6	5.8	2.3
	PowerPC 970MP	2 × 2	2.5	40.0	5.4	7.4
	Power6	2 × 2	4.0	64.0	15.1	4.2
	Power7	4 × 8	3.8	486	265	1.8
	Penryn X3370	2 × 4	3.0	96.0	15	6.4
	Nehalem X5550	2 × 4	2.67	85.1	22	3.9
	Westmere X5680	2 × 6	3.33	159.8	25	6.4
	IvyBridge E5-2697v2	2 × 12	2.7	518.4	92	5.6
	Xeon Phi	1 × 61	1.33	1298	170	7.6
	3.2 Halfpipe and Fullpipe impact with reduction
	3.2.1 Impact of reduction
			Penryn		Nehalem
	HLT	without red	with red	gain	without red	without red	gain
		before cache overflow (in 512 × 512)
	Nopipe	3.76	2.62	×1.4	3.40	3.36	×1.0
	Halfpipe 2	2.96	2.04	×1.5	1.86	1.56	×1.2
	Halfpipe 1	3.42	1.50	×2.3	1.76	1.17	×1.5
	Fullpipe	5.75	4.10	×1.4	3.92	2.92	×1.3
	total gain		×2.5			×2.9
		after cache overflow (in 2048 × 2048)
	Nopipe	62.50	62.50	×1.0	10.50	10.50	×1.0
	Halfpipe 2	27.10	26.90	×1.0	4.47	4.43	×1.0
	Halfpipe 1	19.10	19.10	×1.0	3.42	3.30	×1.4
	Fullpipe	5.80	4.20	×1.4	3.99	2.95	×3.6
	total gain		×14.9			×3.6

Table 3 .

 3 Impact of reduction for Penryn and Nehalem

Table 6 .

 6 Global comparison of State-Of-The Art GPP and GPU for 2048 × 2048 images (H and F stand for Halpipe and Fullpipe transforms, R and M for red and mod optimizations).

	execution time (ms) of SIMD multi-core and many-core GPP
		Nopipe	H+R	F+R	H+R+M	gain
	Cortex-A15	84.08	34.2	60.3	13.69	×6.1
	IvyBridge	8.23	1.01	1.26	0.23	×35.3
	Xeon Phi	3.12	1.10	0.57	-	×5.5
	Power 7	1.79	0.44	1.56	0.23	×7.7
		execution time (ms) of GPU
		Nopipe	H	F	F
		Global	Tex	Tex	Shared	gain
	GTX 580	6.52	2.24	1.40	1.16	×5.6
	GTX Titan (est.)	2.29	0.79	0.49	0.41	×5.6
	K40 (est.)	2.41	0.83	0.52	0.43	×5.6

Impact of Halfpipe and Fullpipe transforms

The three processor families, have the same kind behavior. First, all Halfpipe versions are memory-bound: a cache overflow occurs when the bandwidth is too much stressed. Secondly, Halfpipe1+red is always faster than Halfpipe2+red: even for ARM and Power processors that have a smaller AI, it is more efficient to reduce the number of memory accesses than the algorithmic complexity. The four generations of Intel processors have very close behavior (Fig. 5). The only differences are the scale (related to core number and bandwidth of the processor) of the graph and the data set size when a cache overflow occurs. After overflow, Fullpipe+red is faster than Halfpipe1+red. For IBM processors (Fig. 6), depending on AI, Fullpipe+red is faster (for PowerPC) than Halfpipe2+red, but never outperforms Halfpipe1+red. ARM Cortex-A9 and A15 processors (Fig. 7) have also a similar behavior. But the scales are very different: the Cortex-A15 has around a twice smaller cpp than Cortex-A9. It has a globally improved architecture compared to A9. Moreover A15 can execute one SIMD Neon instruction every cycle instead of one instruction every two cycles for A9. To address cache overflows and to prolongate HLT performance after overflow, mod optimization has been implemented. If a cache overflow should occur, mod will only postpones it to bigger data set, as only 3p × n × sizeof(float) bytes should stay in the cache instead of n × n × sizeof(float) for each image (3 being the size of the convolution, p the number of threads and n the image size). That makes the Halfpipe1+red+mod alway the best version as it requires only two arrays (Ix and Iy) to stay in the cache, instead of three (Ixx, Ixy and Iyy) for Halfpipe2+red+mod (Fig. 8).

Impact of circular buffers and modular addressing

When there is no cache overflow or before it happens, the mod optimization provides an extra speedup on red version close to ×2. After cache overflow, the mod speedup is in the range [×1.8 : ×6.3] (Tab. 4). The higher AI the higher speedup: the max is reach for Westmere that has the same bandwidth than Nehalem, but 50% more cores. As no compiler can vectorize scalar red + mod versions, it is mandatory to write them manually, but it's worth it. Finally, compared to the Nopipe SIMDized+OpenMPized version, the total gain of HLT is in the range [×6.


A Philippot

P Marangé

F Gellot

J F Pétin

B Riera

Fault Tolerant Control for manufacturing discrete systems by filter and diagnoser interactions

The paper deals with an online safety mechanism to define interactions between a diagnoser and a control filter for fault tolerant control of manufacturing discrete systems. The diagnoser observes the plant behavior whereas the control filter ensures the safety from the controller. This online interaction is based by events communication where the control law is never reconfigured. The proposed approach is applied to CISPI platform from the CRAN laboratory (Research Center for Automatic Control of Nancy).

INTRODUCTION

Engineering systems become more and more complex and consequently, faults are more and more present and cause undesired behaviors. Diagnosis information can lead the user in its decision for maintenance or reconfiguration (Nke and Lunze, 2011), but can also allow fault tolerant control. The aim of diagnosis approaches is to detected and isolated with certainty a fault. After this step, it is necessary to reconfigure the controller in order to guarantee the dependability and safety but also to propose a Fault Tolerant Control (FTC) in a degraded mode [START_REF] Blanke | Diagnosis and Fault-Tolerant Control[END_REF], (Paoli et al., 2011[START_REF] Brown | A Prognostic Health Management Based Framework for Fault-Tolerant Control[END_REF].

Ensuring safety of manufacturing system control is currently based on two complementary approaches: control design activities with the objective to avoid unexpected behaviors and safe design activities by the development of online barriers.

First one, we focus on the control design activities with the objective to avoid unexpected behavior. Two main approaches are suggested in this way [START_REF] Faure | Methods for safe control systems design and implementations[END_REF]: (i) control validation and verification (V&V) (Roussel and Faure, 2002), (ii) Supervisory Control Theory (SCT) based on synthesis controller (Ramadge and Wonham, 1989), that enables automatic generation of the controller from the specification, and the uncontrolled behavior of the plant. Most of the time, those designing approaches make two strong assumptions: the behavior of plant devices is not faulty and the designed control is exactly the same as the program that is implemented on the control devices (i.e. code generation deviations or code modifications by maintenance agents are not considered). These assumptions being not realistic in practice, a second approach complements the safe design activities by the development of online barriers like diagnosis or filtering control. Diagnosis of manufacturing systems aims at detecting unsafe behavior of the plant and localizing the components that are involved in the behavioral deviation (Sampath, 1995). Control filtering aims at avoiding that a PLC program provokes plant damages, whatever the PLC program (Marangé, 2008, Riera et al., 2012). The filter is placed between the controller and the plant and inhibits potential dangerous evolutions by checking a set of safety constraints. Nevertheless, the diagnosis and the filter are formally built from models of process behavior. Consequently, hypothesis that the information from the process is correct is made. At least, if the plant situation is unknown, automatic procedures implemented by control filtering and diagnosis may be not efficient. This case generally requires the intervention of human expert to analyze the unknown situation of the plant, and to take emergency decision to drive back the plant in acceptable states.

The aim of this paper is to propose an approach of FTC where diagnosis provides information about the plant to the filter; and vice-versa. Control laws are never reconfigured but the system must always be in safety situation thanks to the filter even in case of plant fault. Models of the plant devices behavior as well as the control rules can be described as Discrete Event Systems (DES), i.e., dynamical systems with discrete state spaces and event-driven transitions [START_REF] Cassandras | Introduction to discrete event systems[END_REF]. The proposed approach provides similar results in term of detection to classical approaches (Sampath, 1995[START_REF] Debouk | Coordinated decentralized protocols for failure diagnosis of discrete events systems[END_REF], Wang et al., 2007 …) but it continues to improve the safety even in presence of faults thanks to the control filter.

The paper is organized as follows. In section 2, the fault tolerant control architecture proposed is presented with a diagnosis and a filtering control sub-sections. A benchmark is studied with results in section 3 before to conclude and propose some future works.

FTC ARCHITECTURE

From the previous discussion, diagnosis approaches make hypothesis that controller information is safe whereas filtering controller approaches are supposed free of faults. The figure 1 presents the FTC architecture. Control law, diagnoser and filter are present in a Remote Terminal Unit (RTU) as a Programmable Logic Controller (PLC) for example. The diagnoser does not use directly the orders sent by the controller but the orders validated by the filter, which set to allows to guarantee the orders correctness. Also, the filter confirms orders according to the plant information (value of sensors/actuators) and the plant state defined by the diagnoser. User can send requests but also have situation awareness thanks to filter and diagnoser.

Diagnoser

In industrial processes, a manufacturing system is a functional chain composed of a controller that emits signals to a plant and receives sensor values. This exchange between controller and plant represents the only observable information available online. Since a diagnoser is defined as an observer of the system, it is necessary to use this information to rebuild behaviors through models.

From literature (Sampath, 1995, Qiu, 2005), centralized approaches appear as unthinkable for large and complex systems. As manufacturing system is composed of mechanical components (actuators/sensors), a methodology to obtain a decentralized diagnosis approach, as [START_REF] Debouk | Coordinated decentralized protocols for failure diagnosis of discrete events systems[END_REF], Wang et al., 2007, Kan et al., 2010) A local diagnoser is a special case of an observer that carries fault information by means of labels attached to states. These labels indicate the types of faults that have been occurred. A local diagnoser is considered as an extended automaton: D i = (X i ∪ X DFi , Z io , δ i , x i0 , T i , l i ) where:

• X i is the set of normal states of NBM i ,
• X DFi is the set of faulty states,

• Z io is the set of observable events by the PE i ,

• δ i : X i × Z i * → X i ∪ X DFi
is the transition function with the expected (δ ei ) and unexpected (δ ui ) functions from a state, • x i0 is the initial state,

• T i is the set of interval time where transition functions are expected between [t min , t max ],

• l i is the set of decision functions of the local diagnoser D i with l i (x) the decision function of the state x which can be one or more fault labels {Fj}. The sets of failure events corresponding to partitions, noted Π f .

Indeed, the methodology is dependent of the control specification (step 2) and if the controller is not safe or if it changes, then diagnosers can return a bad decision in the first case or must be reconstructed in the second case. To have diagnosis independent from the control, diagnoser is obtained from the behavior of PE and the addition of the possible faulty events.

From decentralized diagnosers, a transition function δ i corresponds to a logical expression composed by all the events. It is possible to define all transition functions by the 2 n possibility (with n: number of events and intervals). However, the mechanical structure of components and the use of filters make it impossible some combinations. For example, only one interval time can be activate simultaneously, or thanks to the control filter, opposite orders cannot be sent. Consequently, the complexity depends on the granularity of the local models but also on the performance of the control filter. These diagnosers are independent of the controller specification in its structure thanks to the control filter but not in the definition of the set of interval T i .

The choice of an automaton to represent a local diagnoser permits to compose a library of commonly components. However, this model can be translated as Markov chain or Causal Temporal Signature under some hypothesis.

Control Filter

The control filtering consists in interlacing a filter between the plant and the control law to inhibit the evolutions that can lead the system to a dangerous situation for operators and production resources. This aim is to ensure that the controller outputs (Σ c ), are legal according to plant safety. It means that, for each new evolution of actuators output vector (at t), the filter verifies that these outputs are compatible with the plant state perceived by means of uncontrollable variables Σ uc (inputs sensors (at t, t-1, t-2…), previous outputs (at t-1, t-2…), observers (at t, t-1, t-2…)).

The filter is built according to a set of logical constraints that must be satisfied to let the outputs getting out of the control filter. It is based on the use of safety constraints, which act as logical guards placed at the end of the PLC program, and forbids sending unsafe controllable events to the plant (Marangé et al. 2008), (Riera et al., 2014). Constraints (or guards) are always modeled with the point of view of the control part (PLC), and it is assumed that the PLC scan time is sufficient to detect any changes of the input vector (synchronous operation, possible simultaneous changes of state of PLC inputs).

Safety constraints are expressed in the form of a logical monomial function (product of logical variables, as ∏ ) which must always be equal to 0 (FALSE) at each PLC scan time in order to guarantee the safety. It is considered in this work that the initial safe state for all the actuators (o k ) is defined to 0.

Initially, the constraints are defined in order to ensure a permissive control, and it is assumed that, with the filter, the system remains controllable. In other words, it is possible to design a controller which matches the specifications. For example, considering the previous hypothesis about the safe initial state, a filter which resets all outputs is safe but does not ensure the controllability. Some guards involve a single output at time t (simple safety constraints CSs), other constraints involve several outputs at time t (combined safety constraints CSc). Constraints require the knowledge of Σ c and Σ uc at the current time t and possibly previous times (presence of edge (t-1) for instance noted *). Hence, the filter requires a memory function.

The set of constraints CS is considered as necessary and sufficient to guarantee the safety. In this approach, it is assumed that safety constraints can always be represented as a monomial and depend on the uncontrollable and controllable variables (at t, t-1, t-2…). Filter stops has to stop the process in a safe situation if a safety constraint is not respected.

CSs and CSc can be represented respectively by equation ( 1) and equation ( 2) which are Boolean monomial functions and have always to be False at each PLC scan time. N CSs and N CSc are respectively the number of simple safety constraints and the number of combined safety constraints. N o is the number of outputs.

∀ ∈ 1, N , ∃! ∈ 1, N / = ∏( , Σ ) = 0
(1)

∀ ∈ 1, N , ∃! ( , !, … ) ∈ 1, N O $%&ℎ ≠ ! ≠ ⋯ / * = ∏( , + , … , Σ ) = 0
(2) There are 2 forms of Simple Safety Constraints CSs because they are expressed as a monomial function, and they only involve a single output at time t (equation ( 3) or ( 4)):

∀ ∈ 1, N , ∃! ∈ 1, N / = . ℎ -(Σ ./ )
(3) xor = 000. ℎ 1 (Σ ./ )

(4) These simple safety constraints (CSs) express the fact that if ℎ -(Σ ./ ) which is a monomial (product) function of only uncontrollable variables at t, is TRUE, o k must be necessarily FALSE (equation ( 3)) in order to keep the constraints equal to 0. If ℎ 1 (Σ ./ ) is TRUE, o k must be necessarily TRUE (equation ( 4)).

For each output, it is possible to write equation ( 5) corresponding to a logical OR of all simple safety constraints. (5) 5 = ( , Σ ./ ) is a logical ∑ ∏ function independent of the other outputs at t because only CSs are considered. 5 = ( , Σ ./ ) can be developed in equation ( 6) where 5 =- and 5 =1 are polynomial functions (sum of products, ∑ ∏ ) of uncontrollable variables. Equation ( 6) has always to be FALSE because all simple safety constraints must be FALSE at each PLC scan time. 5 = ( , Σ ./ ) = . 5 =-(Σ ./ ) + 000. 5 =1 (Σ ./ ) = 0 (6)

Taking into account all CSs; it is possible to write equation (7). The definition of constraints set is not formal and the filter robustness must be verified. In (Marangé, 2008) and (Riera et al., 2012), authors proposed to enrich this expert-based approach by a formal identification of the constraints set to ensure its completeness.

The use of this filter allows detecting errors resulting from the controller by making a hypothesis on the accuracy of the information resulting from the plant. Indeed, a fault on the plant can lead:

• Too much restriction: sensor information is going to be blocked in the most critical state and the constraint is not verified while the plant is not in a critical situation. • Too much tolerant: sensor information is going to be in the state which verifies all the time the constraint and thus the filter is going to allow to pass dangerous orders for the plant. This case is to be avoided. The consideration of diagnosis information allows to use the filter in degraded mode. For that purpose, the information resulting from the plant is added by taking into account diagnoser. When a failure arises on a sensor or an actuator, the filter constraints that contain the logical variables associated to the faulty devices becomes unreliable. Authorized signals may be forbidden, and, worse forbidden signals may be authorized. Consequently, the filter constraints must consider the occurrence of a fault or not.

For every fault partition, a flag is set to true when the diagnoser reaches a faulty decision state. This flag determines if the considered variable can be used into the filter constraint (flag=0), or if an equivalent reconstructed information must be used (flag=1). Only the sensor information can be reconstituted by using:

• the expert knowledge (timed or temporal model),

• redundant information or reconstruction logics.

The property defining the dangerous situation has been verified using a model-checker meaning that the filter delivers correct inhibition and authorization even in presence of device faults (with the assumption that the diagnoser is able to detect and localize the fault).

Moreover, as the control filter only concerns safety part and not the functional part, if the component is exchanged or replaced, only the set of constraints corresponding to this component must evolve. For industrial systems, establishment of a constraints library is feasible. In fact, constraints sets are defined for a sub-system of component interaction.

CASE STUDY

The approach is applied to the CISPI platform from the CRAN laboratory (figure 2). This platform implements hydraulic processes involving valves, pumps and tanks and various transmitters (flow, pressure…). Local controllers implement basic control loops and are involved in a global mode management control that enables concurrent access to devices for start, shutdown and normal operation procedures. To avoid damages and failures of the system, as well as the human operator's errors, this experimental platform promotes new forms of control organization that exploits the capacity ambient technologies (sensor network, PDA, mobile control...) to favor safe human/system interactions in any place, at any instant and for any plant operation.

Within the framework of this project, the control filter and diagnoser are implanted to bring a help during the supervisory control of the CISPI system. To illustrate the approach presented in this paper, an automatic valve is considerate. This valve can Independently of the control laws, the sub-system valve must always be in a safety mode. For this, an assumption is made that when a fault is on an actuator, all outputs must be deactivated by the filter. If a fault is on a sensor, the subsystem can be tolerant to this fault.

Diagnoser

From the illustrative example, the valve with sensors fsc and fso constitute one PE and it is possible to identify each faulty event by a label:

• Sensor fsc stuck to 0 (F1) or to 1 (F2) • Sensor fso stuck to 0 (F3) or to 1 (F4) • Valve stuck to fsc (F5) or fso (F6) position • Unexpected fsc (F7) or fso (F9) from 0 to 1 • Unexpected fsc (F8) or fso (F10) from 1 to 0 • Unexpected movement from fsc to fso (F11) or from fso to fsc (F12) • Valve blocked between fsc and fso (F13) Three fault partitions are defined belong to:

• Sensor fsc: Π fsc = {F1, F2, F7, F8} • Sensor fso: Π fso = {F3, F4, F9, F10} • Valve: Π Va = {F5, F6, F11, F12, F13}
With the consideration of the controller information, and thanks to the filter, the valve diagnoser is composed of 9 normal states and 16 abnormal states (Fig. 3) where:

• double circle is the initial state,

• 9 white states are the normal states, • 3 grey states noted F2, F7, F8 represent the abnormal states with detection and isolation of an abnormal behavior with certainty from Π fsc ,

• 3 grey states noted F4, F9, F10 represent the abnormal states with detection and isolation of an abnormal behavior with certainty from Π fso ,

• 4 grey states noted F5, F6, F11, F13 represent the abnormal states with detection and isolation of an abnormal behavior with certainty from Π Va ,

• 6 black states describe the detection of a fault but not the isolation (4 intermediate before isolation).

Figure 3. Valve Diagnoser

The reliability of sensors ensures to be into a safety mode (white states). However, after the detection and isolation of a fault (grey and black states), this diagnoser cannot be anew used. Indeed, it is not possible to rely on misinformation. That is why, it is necessary to preserve the state of the system until the fault is been corrected and reset.

Control Filter

Constraints take into account information of the diagnosers. Information used in the filter is noted X filter and diagnosis information is noted defX. The following flags are done:

• deffso for the partition of valve sensor fso,

• deffsc for the partition of valve sensors fsc,

• defV for the partition of valve actuator V,

To be tolerant on sensors' faults, an expert knowledge is used to estimate the plant information by temporal information. This knowledge can be optimally obtained by FMEA (Failure Mode and Effects Analysis) and so provide a reactivity of detection. For example, figure 4 shows equivalent information of fso and fsc sensors information from a learning chronogram where the estimated value of fso is given by a flag TON1 when an On Delay Timer is activated, and respectively a flag TON2 for the estimated value of fsc.

• for 5 = 1 by 5 B = CDE1

• for 5 = 1 by 5 B = CDE2 Consequently, orders must be deactivated by the filter even in case of faulty event by: G3+HIJ = KL5O 0000000 .

(9)

D G3+HIJ = KL5O 0000000 . D (10)
The set of constraints is defined as following. It is forbidden to maintain an order when the position valve is done (equations ( 11) & ( 12)). It is forbidden to deactivate an order until the ending position valve is not done (equations ( 13) & ( 14)). It is forbidden to activate an order until the starting position valve is not done (equations ( 15) & ( 16)).

For the combined safety constraint of equation ( 17), it is forbidden to activate orders C and O together: . D GT+HIJ 00000000 . 5 GT+HIJ 00000000000 = 0 (14)

1 = G3+HIJ . 5 G3+HIJ = 0 (11) Q = D G3+HIJ . 5 G3+HIJ = (12)

Results and Key Performance Indicators

A first analysis shows that the system is detectable in a bounded delay with certainty for the defined fault partitions. Indeed, all labels are represented in an abnormal state. However, the system is non-diagnosable with certainty. 10 labels on 13 possible are isolated with certainty (one unique label), 3 labels are with an ambiguity. For example, it is not possible to isolate with certainty states with labels {F1, F13} and {F3, F13}. Diagnostic Coverage (DC) is the ratio of the probability of detected dangerous failures (dd) to the probability all the dangerous failures (d). This meaning of the term DC is common to (ISO13849-1) and (IEC/EN 62061). For the valve, the DC is to 76.9%. The standard ISO13849-1 divides DC into four basic ranges: i) <60% = none, ii) 60% to <90% = low, iii) 90% to <99% = medium and iv) 99%+ = high. Consequently, another rule must be present to improve it and to guarantee complete diagnosability notion as defined in (Lin, 1994).

Table 1 presents a comparison between solutions with or without filter and/or diagnosers by simulation of the 13 faulty events under ProcesSim (http://processim.hecfh.be/). Thirteen scenarii have been exploited to obtain these results. With no filter, the valve system is under blocked behavior in 8 cases, into a degraded mode in 1 case and induces a defect situation for 4 cases. We can see that the tolerant situation disappear with the use of the filter only because its purpose is to ensure a safety behavior. When the FTC solution is used, the degraded mode is tolerant to 4 faulty events and above all, it decreases 2 cases of defect situations.

The proposed FTC approach has not been extended on all CISPI platforms yet. But a study has been done on a subsystem composed of 2 automatic valves, one pump and 2 tanks. Another point of view can be also to evaluate the steady state transition probabilities as a KPI. Indeed, a repetitive sequence of normal events can provide an indicator of the system behavior. For the moment, this remark is not treated in these works.

CONCLUSION

A Fault Tolerant Control approach is presented around an interaction between diagnosers and filtering control.

Diagnosis design is refined using enriched information from the real implemented control rules (control + violated constraints of the filter) while control filter benefits from using diagnose information to adapt its set of constraints according to reliable raw or constructed information.

In future works, when diagnosers detect a fault on a component or when the filter detects a mistake on the controller, a significant explanation must be given to a human operator to choose the best policy. A graduated explanation with potential consequences is to return. As last remark, the control filter has been implemented and extended to control design pattern on a real complex system called CellFlex at the University of Reims (www.univreims.fr/meserp/).

Symposium on Information

First

  Author et al. This is an open-access article distributed under the terms of the Creative Commons Attribution 3.0 United States License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

  Figure 1. FTC Architecture

  be closed or open by respectively C and O boolean signals, and two sensors for the open position (fso) and for the closed position (fsc) are present.

Figure

  Figure 2. CISPI Platform

Figure 4 .

 4 Figure 4. Reconstruction of sensors' information For a sensor fault, the plant information is replaced by temporal information: 5 G3+HIJ = KL55 0000000000 . 5 + KL55 . 5 B (7)

  and \ * . \ 0 represent respectively a rising and a falling edge of an order X.

  (2), pp233-263.BIOGRAPHIESPhilippot Alexandre was born in France in 1979. After a Master degree in Systems Optimization and Safety at the University of Reims Champagne-Ardenne in France, he achieved a PhD in Diagnosis of Discrete-events Systems at the University of Reims Champagne-Ardenne in France in 2006. He passed a PostDoc at the Ecole Normale Supérieure de Cachan (ENS-Cachan) in France (2007). Research interests: Discrete Event Systems, Fault diagnosis, Modeling, Supervisory Control Theory, Optimal Control, Manufacturing systems. Currently, he is associated professor at the University of Reims Champagne-Ardenne (URCA) and realises its research at CReSTIC Laboratory (research center in Sciences and Technologies on Information and Communication). Pascale Marangé was born in 1982 in France. After a Master degree at the University of Reims Champagne-Ardenne (URCA) in France, she achieved a PhD in Control synthesis of Discrete-events Systems at URCA in 2008. Currently, she is associated professor at Lorraine University and realizes its research at Nancy Research Center for Automatic Control (CRAN). Her research interests include Verification and Validation of PLC program, Dependability of Discrete-Event Systems, Control synthesis. François Gellot is associated professor at the University of Reims Champagne-Ardenne (URCA) France. He was born in France in 1965. He achieved a PhD in analyze and simulation of Petri nets (1994). Its research interests include supervisory control, discrete events systems modeling and Verification and Validation of PLC program. Jean-François Pétin is a Professor at Lorraine University and realizes its research at Nancy Research Center for Automatic Control (CRAN). Its research interest concerns Dependability of Discrete-Event Systems, Control synthesis, Verification & Validation. Bernard Riera received the Ph.D. degree in Automatic Control from the University of Valenciennes (UVHC), in 1993. He is a Professor of Control Engineering at the University of Reims Champagne-Ardenne (URCA) France, a Researcher at the CReSTIC (research center in Sciences and Technologies on Information and Communication) and associate director of the CReSTIC. Its research interests include supervisory control, supervisory support systems, discrete events and hybrid systems modeling.

  High Level Diagnoser from global specifications is done for uncertainty cases. Diagnosers are implemented as online observers in the PLC. User's decision is given thanks to the set of local labels.

	1. From the plant components, decomposition is made to
	obtain local models called Plant Elements (PEs). A PE
	describes all possible mechanical evolution of the
	component independently of the controller.
	2. From each PE, local desired behavior is extracted.
	Temporal information, obtained by excited events
	simulation, is added to enrich the model. The result is
	an automaton called Normal Behavior Model (NBM).
	3. The third step identifies, from each normal state of
	NBMs, faults which can occur and composes the
	abnormal model by adding of labeled states to obtain
	local diagnosers (D i ). Faults are grouped according to
	the failing component (sensor/actuator) into partitions.
	4. A
	, for
	manufacturing systems with discrete sensors and actuators
	has been developed in previous works (Philippot and Carré-
	Ménétrier, 2011). It is composed of 4 offline steps describe:

Table 1 : Comparison with and without FTC solution

 1

	Diag		No Diag	FTC (Diag
	No Filter	Filter	and Filter)
	Blocked	8	9	7
	Tolerant	1	0	4
	Defect	4	4	2