Computation of a Reference Model for Robust Fault Detection and Isolation Residual Generation

This paper considers matrix inequality procedures to address the robust fault detection and isolation (FDI) problem for linear time-invariant systems subject to disturbances, faults, and polytopic or norm-bounded uncertainties. We propose a design procedure for an FDI filter that aims to minimize a weighted combination of the sensitivity of the residual signal to disturbances and modeling errors, and the deviation of the faults to residual dynamics from a fault to residual reference model, using the H∞-norm as a measure. A key step in our procedure is the design of an optimal fault reference model. We show that the optimal design requires the solution of a quadratic matrix inequality (QMI) optimization problem. Since the solution of the optimal problem is intractable, we propose a linearization technique to derive a numerically tractable suboptimal design procedure that requires the solution of a linear matrix inequality (LMI) optimization. A jet engine example is employed to demonstrate the effectiveness of the proposed approach.


INTRODUCTION
In the past decade, great attention has been devoted to the design of model-based fault detection systems and their robustness [1,2].With the rapid development of robust control theory and H ∞ optimization control techniques, more and more methods have been presented to solve the robust FDI problem.The H ∞ -filter is designed such that the H ∞norm of the estimation error is minimized (see [3][4][5] and the references therein).Some of the approaches used for this problem include frequency domain approaches [6], left and right eigenvector assignment [7], structure parity equation [8], and an unknown input observer with disturbances decoupled in the state estimation error [9].Recently developed LMI approaches offer numerically attractive techniques [10][11][12].
A reference residual model is an ideal solution for robust FDI under the assumption that there are no disturbances or model uncertainty.The idea is to design a filter for the uncertain system that approximates the solution given by the reference model [13].In [14], a new performance index is proposed using such a reference residual model.Frisk and Nielsen [15] give an algorithm to design a reference model and a robust FDI filter that fits into the framework of standard robust H ∞ -filtering relying on established and efficient methods.However, their framework consists in solving two optimization problems successively, which results in a suboptimal solution.
In this paper, we propose a performance index that captures the requirements of fault isolation and disturbance rejection as well as the design of the optimal reference model.The fault isolation performance is measured by the size of the deviation of the fault to residual dynamics from the reference dynamics model, while the disturbance rejection performance is measured by the size of the input to residual and disturbance to residual dynamics.In all cases, the H ∞ norm is used as a measure.The design of the optimal reference model is incorporated in the robust FDI framework.We consider systems subject to norm-bounded or polytopic uncertainties.For systems described by polytopic and unstructured norm-bounded uncertainties, we derive an optimal FDI filter obtained as the solution of a QMI optimization.For systems described by structured uncertainties, we derive a suboptimal QMI-based solution.Since the solution of a QMI optimization is, in general, intractable, we propose a linearization technique to derive a suboptimal design procedure that requires the solution of a numerically tractable LMI optimization.This note extends our work in [16] by proposing algorithms for the design of suboptimal reference dynamics.
The structure of the work is as follows.After defining the notation, we review filter-based FDI techniques for residual signal generation and give the problem formulation in Section 2. Section 3 presents a matrix inequality formulation for the FDI problem, and gives the solution and the design of the optimal reference model for both norm-bounded and polytopic uncertainties in a form of QMI's.Section 4 gives a suboptimal solution in both cases through the use of an algorithm that necessitates solving LMI's.Finally, a numerical example is presented in Section 5, and Section 6 summarizes our results.
The notation we use is fairly standard.The set of real n × m matrices is denoted by R n×m .For A ∈ R n×m , we use the notation A T to denote transpose.For A = A T ∈ R n×n , A 0 (A ≺ 0) denotes that A is positive (negative) definite, that is, all the eigenvalues of A are greater (less) than zero.The n × n identity matrix is denoted as I n and the n × m null matrix is denoted as 0 n,m with the subscripts occasionally dropped if they can be inferred from context.
Let L 2 be the set of square integrable functions.The and dependence on the variable s will be suppressed.For a stable transfer matrix G, we define In section 3, we use the following result.
Lemma 1 (see [17]).Let φ(s if and only if there exists symmetric P that satisfies the following linear matrix inequality:

PROBLEM FORMULATION
Consider a linear time-invariant dynamic system subject to disturbances, modeling errors and process, sensor and actu-ator faults modeled as where x(t) ∈ R n , u(t) ∈ R nu , and y(t) ∈ R ny are the process state and input and output vectors, respectively, and where d(t) ∈ R nd and f (t) ∈ R n f are the disturbance and fault vectors, respectively.Here, B f ∈ R n×n f and D f ∈ R ny ×n f are the component and instrument fault distribution matrices, respectively, while B d ∈ R n×nd and D d ∈ R ny ×nd are the corresponding disturbance distribution matrices [18].We consider two types of uncertainties: norm-bounded and polytopic uncertainties.In the case of norm-bounded uncertainties, where M o represents the nominal model, Δ H = Δ(I −HΔ) −1 , where and where , and H are known and constant matrices with appropriate dimensions.This linear fractional representation of uncertainty, which is assumed to be well-posed over Δ (i.e., det (I − HΔ) / = 0 for all Δ ∈ Δ), has great generality and is widely used in control theories.
In the case of polytopic uncertainties, where M i , i = 1, . . ., p, are known constant matrices with appropriate dimensions.A residual signal in an FDI system should represent the inconsistency between the system variables and the mathematical model.The objective is to design an FDI filter of the form Filter where x(t) ∈ R nk is the filter state and r(t) ∈ R n f is the residual signal.Figure 1 illustrates this filter in the robust residual generation scheme.By defining an augmented state z(t) = [ x(t) T x(t) T ] T the residual dynamics are given by ż(t) where T M r f , T M rd , and T M ru denote the dynamics from faults, disturbances, and inputs to the residual, respectively.Note that dependence on the uncertain data is indicated by a superscript M.
Ideally, the residual signal is required to be sensitive only to faults.This corresponds to T M rd = 0, T M ru = 0, and T M r f / = 0.For fault isolation, it is further required that the fault signature can be deduced from the residual.This corresponds to T M r f ∈ S, where S is a set of transfer matrices with a special structure that depends on the nature of the faults, for example, if the faults can occur simultaneously, S is the set of stable diagonal transfer matrices with nonzero diagonal en-tries.Unfortunately, characterizing a general structured set S is intractable, and we will assume that we can define a subset such that subsequent optimizations over the structured statespace data sets S A , S B , S C , and S D are tractable.For example, if S denotes the set of all stable n f × n f diagonal transfer matrices with nonzero diagonal entries, we may define S A as the set of all n f × n f diagonal matrices with negative entries, and S B , S C , and S D as the sets of all n f × n f diagonal matrices.An example of this simplification procedure is given in Section 5 below.We also replace the requirement of nonzero diagonal entries for S by a condition of the form Due to the presence of disturbances and modelling uncertainties, exact FDI is not possible.For robust FDI, we propose the following, more realistic, problem formulation.
Problem 1. Assume that the system dynamics in (4) are quadratically stable.For any γ > 0, find a stable fault reference dynamics to ensure the requirement of nonzero diagonal entries for S), and find a stable filter of the form given in (8), if it exists, such that the residual dynamics in (10) are quadratically stable and where M = M Δ for norm-bounded uncertainties and M = M ξ for polytopic uncertainties.
Recall that quadratic stability for the dynamics in ( 4) is equivalent to the existence of P = P T 0 such that A T P + PA ≺ 0 for all A (see [19] for more details).
A modified version of Problem 1 uses a weighted cost function, say, of the form where W f , W d , and W u are constant or frequency-dependent weighting functions that can vary the emphasis between fault detection (small T M rd ∞ and T M ru ∞ ) and fault isolation (small In the sequel, we assume that any such weighting functions are absorbed in the system data. Remark 1.The objective is to find the smallest γ for which (12) is satisfied.Indeed, by minimizing γ, we will ensure T M rd ∞ (which measures the disturbance rejection level), T M ru ∞ (which measures the effect of uncertainty), and which measures the deviation of the fault to residual dynamics from the reference dynamics, and hence is a measure of fault isolation capability) to be small.
A poorly chosen reference model can result in a residual generator with poor robustness.Here, we incorporate its design into our scheme so as to improve the robustness of the FDI filter.
In some approaches, a common assumption is that D f = 0 and/or CB f has full rank [20,21]; in others, the assumption D f has full rank and is widely used [22,23].Here, we do not impose any of these assumptions.Note, however, that if D f does not have full column rank, for example, if D f = 0, then this will have an adverse effect on the minimum val- . This would therefore limit the overall performance of the filter, which is measured by the value of γ.

MATRIX INEQUALITY FORMULATION
In this section, we derive a matrix inequality formulation of Problem 1.The main idea is to express (12) in terms of QMIs using the bounded real lemma and change of variables techniques, and then to derive necessary and sufficient conditions for solvability.
The dynamics in ( 12) can be written as It follows from the bounded real lemma that there exists a stable filter of the form given in (8) such that ( 12) is satisfied if and only if there exists P c = P T c such that P c 0 and for all M (see [24,Theorem 3]).We deal separately with norm-bounded and polytopic uncertainties.

Solution with norm-bounded uncertainties
For norm-bounded uncertainties, we separate the terms involving modeling uncertainties from the other terms as where Using this separation, the inequality ( 15) can be rewritten as A calculation shows that T Δ c = FΔ H E + E T Δ H F T , where The next result uses the fact that Δ ∈ Δ to remove explicit dependence on Δ.
Lemma 2 (see [25]).Let Δ be as described in (6) and define the subspaces Let R = R T , F, E, and H be matrices with appropriate dimensions.We have det (I − HΔ) / = 0 and R + FΔ(I − HΔ) and G ∈ Γ such that S 0 and for some scalar λ ≥ 0. In this case, condition (23) is both necessary and sufficient.
When the uncertainty set is unstructured, then Lemma 2 implies that for some λ ≥ 0. Using a Schur complement argument shows that ( 19) is equivalent to where denotes terms readily inferred from symmetry.Next, we introduce a change of variables to linearize the above matrix inequality [26].Assume that n k = n ref + n, that is, the filter order is equal to the order of the system plus the order of the reference model.Let where X, Y , X, Y ∈ R nk×nk are symmetric matrices, and Then T c ≺ 0 if and only if T T T c T ≺ 0. From P c P −1 c = I, we have the following calculations, where boldface letters are used to indicate optimization variables: where we have defined If M and N are invertible, the variables A k , C k , B ky , B ku , D ky , D ku can be replaced by the new variables A, B y , B u , C, D y , D u without loss of generality.We can now rewrite (19) as which is nonlinear in the variables.Note that the nonlinearities involve the terms A ref and B ref only.The constraint P c 0 can be expressed as an LMI as follows: The constraint T ref − ≥ 1 can be expressed as a quadratic matrix inequality using the next lemma.
It is easy to show that φ(s) can be written as follows: The result then follows from Lemma 1.
The next lemma summarizes our result so far by giving necessary and sufficient conditions for the solution of Problem 1, in the case that the uncertainty set is unstructured, in the form of a QMI feasibility problem.Lemma 4. Assume that Δ is unstructured.Then there exist a stable filter of the form given in (8) and a stable fault reference dynamics model where S is defined in (11), such that T ref − ≥ 1, residual dynamics in (10) are quadratically stable, and (12)  Approximate solutions to these QMIs can be obtained by using algorithms with guaranteed global convergence [27,28], as well as local numerical search algorithms that converge (without a guarantee) much faster [29,30].A related discussion of the solution algorithms for QMIs can also be found in [31].In Section 4 below, we develop an alternative procedure for the approximate solution of these QMIs.
Remark 2. In the case that T ref is preassigned to a known value, (37) becomes linear and (39) becomes irrelevant, therefore, the optimal solution is given in a form of an LMI optimization.This case has been considered in [16].
When Δ is structured, we proceed as follows.By using (22) from Lemma 2, we get where Using a Schur complement argument and the expression of E and F in (20), we get As we did for unstructured uncertainties, we use the same matrix T = diag (Π 1 , I) to allow to change variables, it follows that T SG ≺ 0 if and only if T SG := TT SG T T ≺ 0. We multiply T SG by K = diag (I, S) and K T from left and right, respectively, to get T SG ≺ 0 if and only if Therefore, when Δ is structured, we have the following sufficient condition for solvability.
Lemma 5. Suppose that Δ has the structure defined in (6).Then there exist a stable filter of the form given in (8)

Solution with polytopic uncertainties
In this section, we derive necessary and sufficient conditions for solvability of the robust FDI problem for a system subject to polytopic uncertainties, in the form of LMIs.Now, where ) are as defined in ( 17), but with superscript (o) replaced by (i).We assume that the polytopic system is quadratically stable.Recall that ( 12) is satisfied if and only if (15) is satisfied for all M. Now, Emmanuel Mazars et al.

(15) ⇐⇒
Noting that the change of variable defined in (36) is independent of M i , we can therefore use it in this scheme as well.Letting T = diag (Π 1 , I), it follows that where the (L i jk ) T s are as defined in ( 27)-( 33) and ( 34)-( 35), but with the nominal model data M o replaced by M i .Therefore, for polytopic uncertainties, we can derive necessary and sufficient conditions for solvability of Problem 1 in the form of a QMI feasibility problem as follows.

ROBUST FDI FILTER DESIGN USING LMIS
In this section, we give a suboptimal solution to Problem 1.An optimal solution would necessitate the solution of a quadratic matrix inequality and is in general intractable.Here, we propose a linearization procedure to derive an upper bound on the optimal solution that requires solving an LMI optimization problem.
The following general result demonstrates that if we have one feasible solution to a QMI optimization, then we can construct an LMI optimization problem whose solution gives an upper bound on the QMI problem. Then Then by using a Schur complement argument, we get and it follows that That is, Using (53) and noting that ξ(Q, R) 0 and A similar proof can be used to derive (51) and (52).
In order to simplify the subsequent analysis, we adopt the convention that variables appended with a subscript "x" denote feasible values of the variables for the QMIs (37) and (39).
In (37), the only nonlinear terms are We denote the matrix in (37) by T r f w , and set S r f w to be the matrix T r f w , with the nonlinear terms removed. Let and B refx ∈ R nref×n f be given.We use Lemma 7 to derive an LMI formulation.We can write T r f w as where where To linearize the matrix inequality in (39), we need Lemma 7 and the following lemma, whose proof is similar to that of Lemma 7 and is therefore dropped. Then , and D refx ∈ R n f ×n f be given.Using (51) from Lemmas 7 and 8, it is easy to get where where The next lemma summarizes the results of this section by giving a linearized formulation of the optimization problem defined in Lemma 4 using (60) and (64).
and D refx ∈ S D be given.Then there exist a stable filter of the form given in (8) and a stable fault reference dynamics model where S is defined in (11), such that T ref − ≥ 1, residual dynamics in (10) are quadratically stable, and (12) Remark 3.This scheme can also be applied to Lemmas 5 and 6 to obtain a suboptimal solution involving linear matrix inequalities.
Next, we need to choose the initial parameters (Z x , Y x , P refx , A refx , B refx , C refx , and D refx ) to reduce γ.The idea is to derive an algorithm where at each iteration, we solve the optimization problem given in Lemma 9, using the solution of this problem at the previous iteration, for initial parameters.The algorithm will use initial values Z init x , Y init x , P init refx , A init refx , B init refx , C init refx , and D init refx , which must guarantee that the LMI constraints (66) and (67) will be feasible.

Algorithm 1. (1) Set initial values
( Algorithm 1 is convergent, possibly to a local minimum, in the sense that the quantity γ is nonincreasing after each iteration.This can be easily shown using (50) and (52) from Lemma 7, and (63) from Lemma 8. Remark 4. Lemmas 7, 8, and Algorithm 1 can also be applied to other problems involving QMIs and bilinear matrix inequalities (BMIs).The procedure potentially gives an improvement and seems to work well in practice.
Remark 5.In the case that we choose a diagonal structure for T ref , we may use as initial values.This will ensure that T ref is stable and We can solve the LMI optimization problem given in Lemma 4. This will give Z init x , Y init x , P init refx , and so forth.These initial values are not unique and can be chosen using other criteria.Remark 7. The requirement for T ref to be diagonal is to ensure fault isolability in the case that all faults may occur simultaneously.If we ignore disturbances and uncertainty, and assuming perfect fault isolation, then r = T ref f so that fault f i only affects residual r i .If none of the faults can occur simultaneously, then we need only to impose the constraint that T ref is upper (or lower) triangular.While it is not difficult to modify our analysis under these conditions, we only consider the case when T ref is diagonal so that all faults may occur simultaneously since our contribution is focussed on reducing the effect of disturbances and uncertainties under the most stringent fault scenarios.

NUMERICAL EXAMPLE
In order to illustrate the efficiency of Algorithm 1 and the importance of the choice of T ref , we consider a jet engine statespace model [32] supplied by NASA Glenn Research Center given as (71 The system is subject to three disturbances and three potential faults.Here, the setup is given by A square and diagonal structure for T ref is necessary to get fault isolation in the case that the faults occur at the same time. The Then, by following the first two steps of Algorithm 1, we get the optimal γ as and the values given below for Z, Y , and P refx , which will be used to initialize the main loop of the algorithm as follows.(77) Taking into account Remark 8, we implemented Algorithm 1 in Matlab to minimize γ.Table 1 shows the evolution of the optimization following Algorithm 1.
Algorithm 1 can clearly improve the result in a few iterations; γ is reduced by 33% compared to the one obtained with a fixed T ref .This shows that the choice of T ref is essential in our FDI scheme.We get Remark 9.In our numerical experimentation, other choices for T init ref have been used; however, all converged to the same solution but with different numbers of iterations.
In order to show that our filter is robust against disturbances and model uncertainties, we introduce a randomly generated Δ given by   as well as three disturbances.Simulated through MATLAB and SIMULINK, these disturbances are three-band-limited white noises with mean 0 and standard deviation 2. Faults f 1 and f 2 are both simulated by a unit positive jump connected from the 14th second.Fault f 3 , simulated by a soft bias (slope = 0.2), is connected from the 20th second.Figure 2 gives the residual responses before the algorithm, where each residual is dedicated to a particular fault, while Figure 3 gives the optimized residual response using our algorithm.
The lines opt 1 , opt 2 , and opt 3 represent the optimal trajectories that each residual must follow In both figures, each fault can be distinguished from the others and the disturbances; however, in Figure 3, the faults can be more easily distinguished and each residual follows its optimal trajectory (green line) with more accuracy.Furthermore, the disturbances are more attenuated compared to Figure 2, and the jumps that indicate faults are clearer in Figure 3 since their amplitudes are bigger and therefore allow a better fault detection using thresholds [33,34].The isolation performance is clearly effective as each fault produces a deviation of its residual only, without modifying the trajectory of the others.This example illustrates that the designed filter satisfies the performance requirement of FDI which is sufficiently robust against disturbances and modeling errors.

CONCLUSION
This paper has addressed the problem of fault detection and isolation for linear time-invariant systems subject to faults, disturbances, and model uncertainties.We proposed a performance index that captures the FDI requirements.Through QMI formulations, we gave the design of an optimal filter for polytopic or unstructured norm-bounded uncertainties, and a suboptimal filter for structured normbounded uncertainties.Suboptimality in the latter case is inherited from the bounded real lemma, which gives only sufficient LMI conditions for structured uncertainties (see Lemma 2).By allowing the reference model T ref to be variable in our formulation, we get its optimal design, which can be used in other schemes dedicated to fault isolation.The optimal design of this reference model is initially given in a form of a QMI optimization, then a suboptimal solution is obtained by using a linearization procedure which derives an upper bound on the optimal solution of the QMI formulation that requires the solution of an LMI optimization.Note, however, that we have no indication concerning the deviation of our design from the optimal filter.We have also illustrated the effectiveness of our algorithm using a jet-engine example.

Lemma 3 .
Let T ref be as defined above.Then T ref − ≥ 1 if and only if there exists P ref = P ref T ∈ R nref×nref such that and symmetric matrices P ref , Y, and Z such that (37), (38), and (39) are satisfied.If such variables exist, the filter dynamics are obtained by solving (36) where M and N are chosen such that NM T = I − YZ −1 .

Lemma 6 .
Let M = M ξ .Then there exist a stable filter of the form given in(8) and a stable fault reference dynamics modelT ref s = (A ref , B ref , C ref , D ref ) ∈ S, where S is defined in(11), such that T ref − ≥ 1, the residual dynamics in(10) are quadratically stable, and (12) is satisfied for M = M Δ if and only there existA ref ∈ S A , B ref ∈ S B , C ref ∈ S C , D ref ∈ S D , A, B y , B u , C, D y , D u , S ∈ Σ, G ∈ Γ,and symmetric matrices P ref , Y, and Z such that (38), (39), and T i pol ≺ 0 are satisfied for i = 1, . . ., p.
the value of Q c when Z and Y are replaced by Z x and Y x , respectively, let R cx denote the value of R c when A ref and B ref are replaced by A refx and B refx , respectively, and define E x = Q cx − R cx .Using (49) from Lemma 7, we get ) Find the solutions Z, Y, and P ref of the optimization derived from Lemmas 3 and 4, which is linear since the matrix inequalities (37), (38), and (39) become linear when T ref := T init ref is fixed.(The matrices Z, Y, and P ref always exist since the cost function in (12) is bounded because Δ is bounded).(3) Set Z init x := Z, Y init x := Y, and P init refx := P ref .(4) Start loop.
and D init refx = D ref , and go to Step 5. (7) End loop.

Remark 6 .
A more systematic technique for generating valid initial values is as follows: first, generate any B init ref , C init ref , D init ref , and a stable A init ref with the structure chosen, then, compute α= T init ref − .If α > 0, redefine the matrices as A init ref = A init ref /α, B init ref = B init ref /α, C init ref = C init ref /α,and D init ref = D init ref /α, which fulfill that the conditions T init ref are stable and T init ref − = 1.

3 Figure 2 :
Figure 2: Time response of the residual before the optimization.

Figure 3 :
Figure 3: Time response of the residual after the optimization.

Figure 3
Figure 3 also justifies the efficiency of Algorithm 1 to improve the design of the reference model and therefore the overall performance of our filter.
and symmetric matrices P ref , Y, and Z such that terms in B ref can be incorporated in C ref .It follows that we can set B ref = I.Therefore, the nonlinearity in (37) comes only from bilinear terms ZE 1 A ref E T 1 and YE 1 A ref E T 1 .