Robust LPV model-based sensor fault diagnosis using relative fault sensitivity signature and residual directions approaches in a PEM fuel cell

Robust LPV Model-Based Sensor Fault Diagnosis Using Relative Fault Sensitivity Signature and Residual Directions Approaches in a PEM Fuel Cell S. de Lira and V. Puig and J. Quevedo Research group of Automatic Control Departament of Technical University of Catalonia (UPC) Pau Gargallo, 5, 08028 Barcelona Email: [email protected] A. Husar Institut de Robtica i Informtica Industrial (CSIC-IRI) Barcelona Technology Park Carrer Llorens i Artigas 4-6, Barcelona Email: [email protected] Abstract—In this paper, a model-based fault diagnosis methodology for PEM fuel cell systems is presented. The methodology is based on computing residuals using an LPV observer. Sensor fault detection faces the problem of robustness using adaptive thresholds generated with interval observer. Fault isolation is performed using the Euclidean distance between the observed relative residuals and theoretical relative sensitivities. To illustrate the results, a commercial fuel cell Ballard Nexa c is used in simulation where a set of typical fault scenarios have been considered. Finally, the diagnosis results corresponding to those fault scenarios are presented. It is remarkable that with this methodology it is possible to diagnose all the considered faults in contrast with other well known methodologies which use the classic binary signature matrix approach. Index Terms—Fault Detection; Fault Isolation, Fault Isolation, PEM Fuel Cell I. INTRODUCTION The energy generation systems based on fuel cells are complex since they involve thermal, fluidic and electrochemical phenomena. Moreover, they need a set of auxiliary elements (valves, compressor, sensors, regulators, etc.) to make the fuel cell working at the pre-established optimal operating point. For these reasons, they are vulnerable to faults that can cause a emergency shut down or a permanent damage of the fuel cell. To guarantee a safe operation of the fuel cell systems, it is necessary to use systematic techniques, like the recent methods of Fault Tolerant Control (FTC) in [1], which allow increasing the fault tolerance of this technology. The first task to achieve active tolerant control is based on the inclusion of a fault diagnosis system operating in real-time. The diagnosis system should not only allow the fault detection and isolation but also the fault magnitude estimation. In this paper, a model based fault diagnosis approach is proposed as a way to diagnose faults in fuel cell systems. The model-based fault diagnosis is based on comparing on-line the real behavior of the monitored system obtained by means of sensors with a predicted behavior obtained using a Linear Parameter Varying (LPV) dynamic model with an Luenberger observer scheme [2]. In case of a significant discrepancy (residual) is detected between the observer outputs and the measurements obtained by the sensors, the existence of a fault is assumed. Fault isolation is based on generating a set of residuals with the available sensors thanks to they present different directional sensitivity to the set of possible faults. The contributions of this paper are: first, the use of a LPV observer for fault detection and the second and the most important is the fact that dealing with fault relative sensitivity approach methodology for fault isolation is presented for the kind of faults that otherwise would not be separable using a classic fault isolation approach. II. FOUNDATIONS OF THE FAULT DIAGNOSIS METHODOLOGY The proposed methodology of fault diagnosis for fuel cell systems which is used in this paper is mainly based on classic FDI theory of model-based diagnosis described in [3], [4] and [5]. The task of fault diagnosis consists of determining the type of fault with as much details as possible (fault location, time and size). Thus, two subtasks can be considered: fault detection and fault isolation. The principle of model-based fault detection is based on checking the consistency of measured (yk) and predicted behaviors (ˆyk) by computing residuals rk. Those residuals are obtained from the discrepancy of measured input (uk) and outputs (yk) using the set of sensors installed and the analytical relations obtained by system modeling: rk=ψ(yk,uk)(1) where ψis the residual generator function that depends on the type of detection strategy used (parity equation [3] or observer [6]. At each time instance, k, the residual is compared with a threshold value that should be determined taking into noise and modelling uncertainty. Considering the whole set of residuals available, a set of fault indicators,φk=[φ1k,φ 2k...,φ nφ]are obtained as follows: φik =0if |ri|≤τi 1if |ri|>τ i(2) where τiis the threshold associated to the residual ri(k). Fault isolation consists of identifying the fault that affects the system. It is carried out by using the fault indicators, φ (generated by the detection module) and its relation with all the considered faults, fk={f1k,f 1k,...,f nk}. The method most often applied is based on the relation defined on the Cartesian product of the set of faults FSM ⊂φ×f, where FSM, known, the theoretical fault signature matrix [3]. An element FSMij of this matrix will be one, if a fault fjk is affected by the residual rik. In this case, the value of the fault indicator φi(k) must be equal to one when the fault appears in the monitored system. Otherwise, the element FSMij will be zero. Fault signature matrix can be obtained from the structural analysis of analytical relations coming from the model and the set of available sensors [1]. Sensitivity analysis is the key point in the evaluation of fault symptoms in order to obtain a final diagnosis, identification of significant and problematic components plays an important role. That sensitivity analysis require understanding the system structure and component relationships. According to [3], the sensitivity of the residual to a fault is given by Sf=∂r ∂f (3) which is a transfer function that describes the effect on the residual (r) of a given fault (f). Sensitivity provides quantitative information about the effect of the fault on the residual and qualitative information in their sense of variation (sign). The use of this information at the stage of diagnosis allows to separate faults that although show the same theoretical binary fault signature, they present different sensitivity values. III. FAULT DETECTION Many model-based fault detection techniques, are mostly based on linear models. However, fuel cells are inherently non-linear [7]. An attractive alternative to represent non-linear systems is to use techniques based on LPV models. The LPV approach is particularly appealing whenever non-linear plants can be modeled as a linear time-varying systems with on-line measurable state depending parameters as is the case of PEM fuel cells. A. Linear Parameter Varying Model A LPV system in discrete-time state space form with input and output sensor faults can be expressed as xk+1 =A(ϑk)xk+B(ϑk)uk+Fa(ϑk)fak +wk(4) yk=C(ϑk)xk+D(ϑk)uk+Fy(ϑk)fyk +vk where xk∈ nx,uk∈ nuand yk∈ nyare, respectively, the state, input, and output vectors, A,B,C,Dare the parameter model at ϑk,fy(k)∈ ny and fa(k)∈ nu represent faults in the system output sensors and actuators respectively being Fy(˜ ϑ)∈ ny×ny and Fa(˜ ϑ)∈ ny×nu their associated matrices, in faultless mode the faand fyare zeros. The process and measurement noises are wk∈ nx and vk∈ nv.ϑk∈ nθis the system vector of time-varying parameters that change with the operating point scheduled by some measured system variables pk(pk:=p(k)) that can be estimated using some known function: ϑk=f(pk)(5) The type of LPV systems considered in this paper assumes an affine dependence within the parameter vector Θspace: Θ=ϑk∈ nθ|ϑk≤ϑk≤ϑk There exists different ways to obtain a LPV model for a non-linear system. Some methods use directly the non-linear equations of the system to derive the LPV model (using for example a state transformation or the Jacobian linearization) [8]. Another kind of methods uses multi-model identification that consists basically in two different steps. First part, a set of LTI model is identified at different equilibrium points by classical methods (on-line or off-line). As second part of this methodology a multi-model is obtained by using an interpolation law that commutes the local LTI model according to the operating point[9]. B. Linear Parameter Varying Observer A LPV observer with Luenberger structure for the state estimation of the system described in Eq. (4) is given by ˆxk+1 =A(ϑk)ˆxk+B(ϑk)uk+(6) L(ϑk)(yk−ˆyk)+Fa(ϑk)fak ˆyk=C(ϑk)ˆxk+D(ϑk)uk+Fy(ϑk)fyk ˆxk+1 =A(ϑk)ˆxk+B(ϑk)uk+(7) L(ϑk)(yk−ˆyk)+Fa(ϑk)fak ˆyk=C(ϑk)ˆxk+D(ϑk)uk+Fy(ϑk)fyk where Lis the observer gain to be designed to guarantee stability for ϑk∈Θ. This gain is computed at operating point ϑkusing LMI formulation for Pole-Placement within a wide class of pole clustering regions that is founded in an extended Lyapunov Theorem (see [10]). The motivation for seeking pole clustering in specific regions inside the unitary circle is to obtain a fast observer dynamics for all considered operating points. C. LPV interval observer Definition 1.Consider the state estimator given by Eq. (4) in faultless mode, an initial compact set X0and a sequence of measured inputs (uj)k−1 0and outputs (yj)k 0. The exact uncertain estimated state set at time k is expressed by Xk={ˆxk:(ˆxj=A(ϑj−1)ˆxj−1+B(ϑj−1)uj−1 +wj−1+L(yj−1−ˆyj−1))k j=1, (ˆyj−1=Cˆxj−1+vj−1)k j=1|x0,ˆx0∈X0, (ϑj−1∈Θ,wj−1∈Wj−1,vj−1∈Vj−1)k j=1 (8) The uncertain state set described in Definition 1 at time kcan be computed as approximation of uncertain state set at time k−1. A exhaustive computing process could be presented if the exact set of estimated states is required. In order to reduce complexity a bounded set could be used such as a box (interval hull) or other geometric regions easy to compute. Before introducing such algorithm, an additional definition should be introduced. Definition 2.Consider the state estimator given by Eq. (4), the set of uncertain states at time k-1, Xk−1and the input/ouput values (uk−1,yk−1,yk). Then, the approximated set of estimated states at time k based on the measurements up to time k-1 is defined as Xe k={ˆxk:A(ϑk−1)ˆxk−1+B(ϑk−1)uk−1+wk−1 +L(yk−1−ˆyk−1),ˆyk−1=Cˆxk−1+vk−1 |ˆxk−1∈Xk−1,ϑ k−1∈Θ,¯wk−1∈Wk−1,vk−1∈Vk−1} (9) Analogously, considering a measurement equation in (4) the approximated set of estimated outputs Ye kcan be determined. Using previous definition, the set of estimated states (or outputs) introduced in Definition 1 will be approximated iteratively using zonotopes. From these zonotopes, an interval for each state variable can also be obtained by computing the interval hull of the zonotope. The sequence of interval hulls Xe kwith k∈[0,N]will be called the interval observer estimation of the system given by Eq. (4). Analogously, the sequence of interval hulls Ye kcan be obtained. Following previous idea, Algorithm 1 is proposed to determine an approximation of set of uncertain estimated states. Algorithm 1 Fault Detection using Interval observer 1: fault ←FALSE 2: k←0 3: Xe k⇐X0 4: while fault =FALSE do 5: Obtain input-output data {uk,yk} 6: Compute the approximated set of estimated states, Xe k 7: Compute the approximated set of estimated outputs, Ye k 8: Compute the interval hull of the approximated set of estimated states, Xe k=[xk,xk] 9: Compute the interval hull of the approximated set of estimated outputs, Ye k=yk,yk 10: if [yk]∩Ye k=∅then 11: fault ←TRUE 12: end if 13: k←k+1 14: end while D. Residual Generation The application of interval observers to fault detection involves in testing whether the measured output is consistent with the one given by the observer using a faultless model and parameter uncertainty. If an inconsistency is detected, the existence of a fault is proved. The consistency match is based on generating a residual by comparing the measurements of physical variables ykof the process with their estimation ˆyk provided by the observer: rk=yk−ˆyk(10) where rk∈ nyis the residual set. Then, the fault detection test consists in checking the satisfaction of: |rik|≤τik (11) where τik is the threshold associated to the residual rik. In case that this condition is not satisfied a fault can be indicated. IV. FAULT ISOLATION A. Fault Sensitivity Analysis The isolation approach presented in Section II uses a set of binary detection (Boolean) tests to compose the observed fault signature. However, the use of binary codification of the residual produces a lack of information that can lead to wrong diagnosis when applied to dynamic systems. This derives in some faults that are not isolable because they present the same theoretical binary fault signature [11]. To avoid this problem is possible to use additional information associated with the relationship between the residuals and faults, such as sign, sensitivity, and order or activation time, to improve the isolation results [11]. B. Proposed methodology In this work, a new method for fault diagnosis system design is proposed that exploits the information provided by the sensitivity in the case where the fault magnitude is unknown. 1) Sensitivity of the residual to an output sensor fault: Considering computational form of the residual generator using Eq. (11) and input-output format of plant model with q−shift expressed in terms of the effects caused by faults r(k)=ro(k)+(I−H(q−1,ϑ))(Gfa(q−1,˜ ϑ)fa(k) +Gfy(q−1,˜ ϑ)fy(k)) −Gfu(q−1,˜ ϑ)fu(k)(12) where ro(k)=−G(q−1,ϑ)uo(k)+(I−H(q−1,ϑ))yo(k)(13) Considering the residual internal form given by Eq. (12) and consideringthe fault residual sensitivity definition given by Eq. (3), it is possible to compute the sensitivity for the case of an output sensor fault, fy, which is given in q−transfer form by a transfer function matrix Sfy with size of ny ×ny, expressed as: Sfy(q−1,ϑ)=I+Cy(ϑ)(qI−Ay(ϑ))−1By(ϑ) Dy(ϑ) =I+Hy r(q−1,ϑ)−1Fy(ϑ)(14) where Ay(ϑ)=A(ϑ)−L(ϑ)C(ϑ),By=Fy(ϑ)fy, Cy(ϑ)=C(ϑ),Dy(ϑ)=D(ϑ)and Hy r(q−1,ϑ)= Cy(qI−Ay(ϑ)) By(ϑ)+Dy(ϑ). 2) Sensitivity of the residual to an input sensor fault: The residual sensitivity to an input sensor fault, fu, is given by the transfer function matrix Sfu which dimation is ny ×nu. The dynamic residual sensitivity is computed as: Sfu(q−1,ϑ)=I+Cu(ϑ)(qI−Au(ϑ))−1Bu(ϑ) =I+Hu r(q−1,ϑ)−1Fu(ϑ)(15) where Au(ϑ)=A(ϑ)−L(ϑ)C(ϑ),Bu=B(ϑ)− L(ϑ)D(ϑ),Cu(ϑ)=C(ϑ),Du(ϑ)=D(ϑ)and Hu r(q−1,ϑ)=Cu(qI−Au(ϑ)) Bu(ϑ)+Du(ϑ). In order to perform diagnosis, the algorithm uses a theoretical fault sensitivity matrix (FSMsensit), see Table I. Each value of this matrix, denoted as Srifj, contains the sensitivity of the residual rito the fault fj. Each value of this matrix, denoted as Srifj, contains the sensitivity of the residual rito the fault fj. ri/fjf1f2·· · fm r1S11 S12 ··· S1m r2S21 S22 ··· S2m . . .. . .. . .. . .. . . rnSn1Sn2··· Snm TABLE I THEORETICAL FAULT SENSITIVITYSIGNATURE MATRIX Although sensitivity depends on time in case of a dynamic system, here the steady-state value after a fault occurrence is considered as it is also suggested by [3]. In order to perform real time diagnosis, the observed sensitivity So rifjshould be computed using the current value of the residual ri(k)when a fault f(k)is detected. But, this requires the knowledge of the fault magnitude or an estimation of it. To solve this problem, this paper attempts to design the diagnosis using the concept of relative sensitivity rather than absolute sensitivity given by (3). The observed relative fault sensitivity is defined as Srel,o ri,rl,fj=Srifj Srlfj =ri(k)fj(k) rl(k)fj(k)=ri(k) rl(k)(16) where ri l=max i=1,...,m (|Srifm|)(17) The residual (rl) for each fjto be used as relative factor in Eq. (16), that guarantees the best isolation performance it is based on theoretical sensitivity value from one fault to another, see Eq. (17). Using the concept of relative sensitivity, the theoretical relative fault signature matrix FSMrel sensit presented in Table II is introduced. f1f2··· fm r2/rlSrel,t r2rl,f1Srel,t r2rl,f2 ··· Srel,t r2rl,fm r3/rlSrel,t r3rl,f1Srel,t r3rl,f2 ··· Srel,t r3rl,fm ··· ··· ··· ··· ··· rm/rlSrel,t rmrl,f1Srel,t rnrl,f2 ··· Srel,t rnrl,fm TABLE II THEORETICAL FAULT SIGNATURE MATRIX USING RELATIVE SENSITIVITY RESPECT TO rl The diagnostic Algorithm 2 computes in real-time the observed relative sensitivities (16) as a ratio of residuals providing a point vector in the relative sensitivities space. The vector generated will be compared with vectors of theoretical fault places stored into the relative sensitivity matrix FSMrel sensit. The theoretical fault signature vector with a minimum distance with respect to the fault observed vector is postulated as the possible fault as: min ds f1(k),...,d s fn(18) where the distance is calculated using the Euclidean distance between vectors ds fn(k)=sqrt ⎛ ⎜ ⎝Srel,o r2r1,f1(ϑk)−Srel,t r2r1,f1(ϑk)2 +... +Srel,o rmr1,fn(ϑk)−Srel,t rmr1,fn(ϑk)2⎞ ⎟ ⎠ (19) Algorithm 2 summarizes the fault isolation procedure. Algorithm 2 Fault Isolation 1: fault ←TRUE 2: k←k+1 3: while fault =TRUE do 4: Obtain input-output data {uk,yk} 5: Compute the approximated set of estimated outputs, using Algorithm 1 step: 6 to 9. 6: Compute the Theoretical Sensitivity Dynamic for process output ynat fj,St rlfj(ϑk) 7: Compute the Observed Sensitivity Dynamic for process output ynat fj,So rlfj(ϑk) 8: Compute the Euclidean distance between Theoretical and Observed Sensitivity, ds fn(k)= St rlfn(ϑk),So rlfn(ϑk) 9: Fault isolation →f 10: end while V. CASE STUDY A. Introduction In this paper, a PEM FC model developed by [12] is used as case study for fault diagnosis. This model has been modified in a calibration procedure for a comertial fuel cell, (Ballar Nexa c ) using lqs non linear data fitting from lab data. The overall FC system could be partitioned as: fuel cell stack and auxiliary components. The auxiliary components are compressor, supply and return manifold, cooling and humidifier. Stack voltage, anode and cathode flow, membrane hydration models belong to fuel cell stack subsystem. Figure 1 presents a conceptual diagram of the FC system. Fig. 1. PEM FC Process Block diagram B. PEM FC Dynamic Model The non-linear parameter model depends on state variables, which infers a high level of non-linear model, for that reason an observer approach is proposed. This model has the following: Characteristics: •The law of mass conservation is used for mass balance estimate. •Physics laws and some empirical equations are used. •The properties are based on electrochemical, thermodynamic and zero-dimensional fluid mechanics principles. Assumptions: •Cathode and Anode volume are taken as single volume. •Perfect temperature control at cooling system. •Stack temperature is constant. •The electrochemical reaction at membrane is performed instantaneously. The resulting dynamic model equation is described by ˙ωcp =1 Jcpωcp (Pcm −Pcp) ˙ Pom =RairTom Vom (Wca,o −Wom,o) ˙mom =Wca −Wom,o ˙ Pim =γRa Vim (WcpTcp −Wim,oTim) ˙mim =Wcp −Wsm,o ˙mH2=WH2,i −WH2,o −WH2,r ˙mw,an =Wvan ,i −Wvan,o −Wvmbr ˙mN2=WN2,i −WN2,o ˙mO2=WO2,i −WO2,o −WO2,r mstCst ˙ Tst =˙ Hreac −Pelec −˙ Qrad −˙ Qconv where the subindex i,o, represents inlet and outlet flow respectably and subindex an,rm,mbr,sm,cp means anode, return manifold, membrane, supply manifold and compressor, respectively. This non-linear model can be transformed into a LPV model in state space form by considering the following definition for •states: x=[ωcp Pom mom Pim mim mH2mw,an mN2 mO2Tst]T. •inputs: u=[Ist vcm]Twhere the scheduling variable is Ist. •outputs: y=[ωcp rO2vst Psm]and •perturbations: z=[Tamb]. The units of all these variables are in compatible magnitudes ( kRPM , gr , volt bar , amp ). C. LPV Model Analysis Since the model is a highly non-linear model, it is difficult to obtain a explicit dynamic model with independence of the parameter with the operating point, Ist. To face this problem the problem an LPV model is an alternative way to solve this problem. The structure of the model in LPV form has the following structure A= ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ a11 00a14 00 000 0 0a22 a23 a24 a25 00a28 a29 a210 0a32 0a34 a35 00a38 a39 a310 a41 00a44 a45 00a48 a49 a410 a51 00a54 00 0a58 a59 a510 000a64 0a66 a67 00a610 000a74 0a76 a77 00a710 0a82 0a84 a85 00a88 a89 a810 0a92 0a94 a95 00a98 a99 a910 00000a106 a107 0a109 0 ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ; C= ⎡ ⎢ ⎢ ⎣ c11 00000000 0 000c24 00000 0 000c34 c35 00c38 c39 c310 0000 0c46 c47 c48 c49 c410 ⎤ ⎥ ⎥ ⎦ ; B= ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ b11 b12 00 00 00 00 b71 0 b81 0 00 b91 0 b101 0 ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ;D= ⎡ ⎢ ⎢ ⎣ 00 00 d31 0 d41 0 ⎤ ⎥ ⎥ ⎦ subsectionFaulty PEM FC scenarios implementation In order to test the proposed methodology in the PEM FC system, a set of commons possible fault in the sensors was selected and implemented in simulation as a Fault Generator Block. Through this block it is possible to select and simulate the effects that cause each of the proposed faults over the Nexa c fuel cell. Each fault acts over a specific output sensor as additive fault. Note that interaction could appear in fuel cell dynamic when the fault acts over the system because since fault changes the process dynamics making the fault diagnosis process a non simple issue. The Table III describe the set of faults which were considered as additive and abrupt faults. D. Results Using Eq. (14) and evaluating it in steady state, then Sfy(∞) = lim q→1Sfy(q−1,ϑ)(20) =I+C(ϑ)(I−A(ϑ))−1L(ϑ)−1 Fy(ϑ) From Eq. (20) Table IV is computed ID Fault Description f1The speed sensor (ωcp) presents suddenly an offset value. step of 6units(kRPM). f2The supply manifold pressure (Psm) measurement suffers a suddenly offset. Step of 0.01(bar) f3Suddenly change in the Oxygen ratio measurement because of a sensor degradation. Ramp with slop of 0.5 f4The stack ltage (vst) measurement suffers a suddenly an offset . step of 4(Volts) TABLE III DESCRIPTIONOF THE ADDITIVE FAULT SCENARIOS IMPLEMENTED IN FGB. f1f2f3f4 r17.939483 5.76E−05 6.65E−05 0.002366 r2−2.83E−05 0.009996 −6.42E−06 −0.00019 r3−0.03957 0.009236 0.510927 0.373165 r40.265891 −0.07013 −0.08645 1.115084 TABLE IV THEORETICAL FAULT SENSITIVITYSIGNATURE MATRIX Based on the values computed in Table IV, its more clear if the range of study for fault sensitivity is [−1,1]. For that reason a scaling operation is introduced as Eq. (21), then Table V is obtained. rs i=ri n  i=1 r2 i (21) f1f2f3f4 rs 10.999 0.0008 0.000128 0.00201 rs 2 −3.6E−6 0.1399 −1.2E−5−0.00016 rs 3 −0.005 0.12928 0.9859 0.3173 rs 40.0334 −0.98169 −0.16682 0.9483 TABLE V THEORETICAL FAULT SENSITIVITYSIGNATURE MATRIX SCALED VALUES Using the values from Table V and selecting the ratio factor, rl, for each fjfault using the criteria introduced in Eq. (17). rj l=[r1,r 4,r 3,r 4] With the last information it is possible to compute the theoretical relative sensitivity matrix, which will be stored in the memory of the isolation algorithm. f1f2f3f4 r2/r1−3.6E−06 r1/r4−0.00082 r1/r30.00013 r1/r40.002122 r3/r1−0.00498 r2/r4−0.14253 r2/r3−1.3E−05 r2/r4−0.00017 r2/r10.03349 r3/r4−0.13169 r4/r3−0.16919 r3/r40.334651 TABLE VI THEORETICAL RELATIVE FAULT SENSITIVITY SIGNATURE MATRIX. In order to test the proposed methodology, an output sensor fault, f1is simulated. Using Eq. (14) and Table VI it is posible to know the Euclidian distance at each kfrom Eq. (19) for each fault in its space of residual ratios. 0 5 10 15 20 25 30 35 40 20 40 60 kRPM ωcp 0 5 10 15 20 25 30 35 40 5 10 15 Ist 0 5 10 15 20 25 30 35 40 1.02 1.04 1.06 bar Pim 0 5 10 15 20 25 30 35 40 5 10 15 Ist 0 5 10 15 20 25 30 35 40 0 5 10 ratio rO2 0 5 10 15 20 25 30 35 40 5 10 15 Ist 0 5 10 15 20 25 30 35 40 30 40 50 volts vst 0 5 10 15 20 25 30 35 40 5 10 15 Ist 0 5 10 15 20 25 30 35 40 −10 0 10 rωcp kRPM 0 5 10 15 20 25 30 35 40 −5 0 5x 10−3 rPsm bar 0 5 10 15 20 25 30 35 40 −5 0 5 rrO2 Time (Ts) ratio 0 5 10 15 20 25 30 35 40 −5 0 5 rvst Time (Ts) volts Schedule variable Threshold Plant output error Residual Threshold Fig. 2. Real and estimated output variable Figure 2 shows the scheduling variables, process outputs, estimated model output and the associated residual and thresholds. The plant suffers the fault f1at k= 80. Fault detection process in performed online. At time k=82the plant output, ωcp,cross one of the threshold, then the process isolation is performed. Figure 3 shows the Euclidean distance from observed relative residuals to relative sensitivity matrix already stored in the memory of the isolation algorithm. Note that computing of the dynamic theoretical relative sensitivity, improves the fault isolation process. Here the fault with the minimum distance is the fault presented in the ωcp sensor, f1. 0 5 10 15 20 25 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Time Distance from Observed relativity residuals to Theoretical relativity sensitivity matrix Fig. 3. Dynamic Euclidian Distance VI. CONCLUSION In this paper, a new LPV model-based fault diagnosis methodology based on the relative fault sensitivity has been presented and tested. An advantage of this new methodology is twofold: first, the variation of the dynamics with the operating point is considered by using an LPV observer when generating residuals. Second, a fault isolation algorithm based on the relative fault sensitivity concept is proposed. This method allows isolating faults that are not isolable considering only a binary (or a sign) fault signature matrix. To prove this methodology, a PEM fuel cell case study well-known in the literature has been used. The case study was modified to include a set of possible fault scenarios that try to reflect the most common faults. All the considered faults have been tested with the new diagnosis methodology, which has diagnosed correctly all of them. VII. ACKNOWLEDGMENTS This work has been partially funded by the grant CICYT HYFA DPI2008-01996 of Spanish Ministry of Education and by a grant from Consejo Nacional de Ciencia y Tecnologia (CONACyT), M´exico. The authors wish also to thank to the Instituto de Robotica i Informatica Industrial (IRI/UPC-CSIC) for its knowhow sharing, specially to J. Riera, M. Serra and D. Feroldi. REFERENCES [1] M. Blanke, M. Kinnaert, J. Lunze, and M. Staroswiecki, Diagnosis and Fault-Tolerant Control. Springer-Verlag Berlin Heidelberg, 2003. [2] D. Luenberger, “An introduction to observers.” IEEE Automatic Control, vol. 16, no. 6, pp. 596–602, 1971. [3] J. Gertler, Fault Detection and Diagnosis in Engineering Systems.New York: Marcel Dekker, 1998. [4] M. Staroswiecki and G. Comtet-Varga, “Analytical redundancy relations for fault detection and isolation algebraic dynamic systems,” Automatica, vol. 375, no. 5, pp. 687–699, 2001. [5] I. R., “Model-based fault-detection and diagnosis-status and applications,” Annual Reviews in Control, vol. 29, no. 1, pp. 71–85, 2005. [6] J. Chen and R. Patton, Robust Model-Based Fault Diagnosis for Dynamic Systems. Marcel Dekker, 1999. [7] C.Rayment and S. Sherwin, “Introduction to fuel cell technology.” in University of Notre Dame, 2003, department of Aerospace and Mechanical Engineering. [8] B. P., T. D., Y. Y., and P. R, “From differential equations to pdc controller design via numerical transformation,” Computers in Industry, vol. 51, pp. 281–297, 2003. [9] M.-S. R. and J. T. A., Multiple model approaches to modeling and control, 2nd ed. London, Taylor and Francis., 1997. [10] M. Chilali and P. Gahinet, “h∞design with pole placement constraints: an LMI approach,” IEEE Transactions on Automatic Control, vol. 41, no. 3, pp. 358–367, 1996. [11] V. Puig, J. Quevedo, T. Escobet, and J. Meseguer, “Towards a better integration of passive robust interval-based fdi algorithms,” 6th IFAC SAFEPROCESS, vol. 16, no. 5, 2006. [12] J. Pukrushpan, H. Peng, and A. Stefanopoulou, “Pem fc dynamic system,” Journal of Dynamics Systems, Measurement, and Control, vol. 126, pp. 14–25, 2002.