Fault Detection and Isolation for Linear Discrete-Time Systems Using ...
Joint 48th IEEE Conference on Decision and Control and 28th Chinese Control Conference Shanghai, P.R. China, December 16-18, 2009
ThBIn1.10
Fault detection and isolation for linear discrete-time systems using input/output measurement analysis Ze Zhang, Imad M. Jaimoukha Abstract— In this work, a new approach to fault detection and isolation is proposed for linear discrete-time systems subject to faults and bounded additive disturbances using a system model, together with input/output measurements over a finite estimation horizon. We use the proposed method to compute upper and lower bounds on the faults. These bounds are the tightest possible given the system model, input/output measurements and the bounds on the initial state and disturbances. Linear programming optimization techniques are used to obtain the bounds. Moreover, a subsequent-state-estimation technique, together with an estimation horizon update procedure are given to allow the fault detection and isolation process to be repeated in an iterative procedure. Finally, the approach is verified using a numerical example.
I. I NTRODUCTION Model–based fault detection and isolation (FDI) for complex dynamic systems has received extensive attention both in literature and practice due to reliability, security and fault tolerance considerations. One of the most common approaches, the observer-based FDI approach, designs a dedicated observer to generate residual signals to provide fault signatures. See [13], [10], [7], [9], [11] for more details. The observer effectively cancels the (nominal) process dynamics and is sensitive only to disturbances, plant/ model mismatch (often recast as disturbances) and faults. The filter design objective is then to reduce the sensitivity to disturbances and/or plant/model mismatch for a given sensitivity to faults. There are two main approaches for achieving this objective, namely, perfect and approximate disturbance decoupling. The former aims to decouple the residual signal from disturbances exactly, while in the latter, the transfer matrix from disturbances to the residual signal is required to be small in either the H2 or H∞ norm sense. Moreover, in either case, the transfer matrix function from faults to residual is required to be diagonal for fault isolation. Both approaches have been given much attention and implemented using several methods. Patton and Chen [12] proposed left and right eigenvalue assignment method. Zhang and Ding [15] used a periodic parity space for linear discrete-time system. Chen, Patton and Zhang [1] solved the robust FDI problem by using an unknow input observer with disturbances decoupled in state estimation error. See also [4], [2] for other approaches. By using these methods, the decoupling problems can be The authors are with the Control and Power Group, Department of Electrical and Electronic Engineering, Imperial College, London. Email: [email protected] and [email protected]. Corresponding author: Ze Zhang, Control and Powe Group, Department of Electrical and Electronic Engineering, Imperial College, London SW7 2BT, UK. Fax: +44 207 594 6282. Telephone: +44 207 594 6279. Email: [email protected]
978-1-4244-3872-3/09/$25.00 ©2009 IEEE
transformed to a sensitivity optimization problem, which seeks to increase the sensitivity of the residual to faults and simultaneously reduces the sensitivity to disturbances and plant/model mismatch. Recently, linear matrix inequality (LMI) approaches have received increasing attention since they provide numerically attractive techniques for formulating robust FDI decoupling problems. Hou and Patton give a realization of fault detection observer design based on the bounded real lemma [6]. Furthermore, a new performance index which introduces a reference residual model is formulated by Zhong, Ding, Lam and Wang using LMI techniques in [16]. As illustrated above, most of existing FDI schemes propose generating a residual that attempts to reproduce the fault signal. In this work, we propose a new approach to FDI for linear discrete-time systems subject to bounded additive disturbances and faults using a dynamic system model as well as input/output measurements over a finite estimation horizon. By using the proposed method, upper and lower bounds on the faults are computed. Linear Programming (LP) optimization techniques are employed to obtain the bounds. Furthermore, we also propose a subsequent-stateestimation technique together with an estimation horizon update procedure which allow the on-line fault detection process to be operated in an iterative procedure. This work is organized as follows. After defining our notations, Section II defines the problem settings. Section III gives the problem formulation and derives upper and lower bounds on the fault signal in the form of solutions to linear programming problem. In Section IV, the method for subsequent-state-estimation as well as an estimation horizon update procedure is described, which complete the description of our fault detection algorithm. Examples are given to demonstrate the effectiveness of the proposed scheme in Section V. Finally, Section VI summarizes this work. The notation we use is fairly standard. The set of real n×m matrices is denoted by R n×m . For A ∈ Rn×m we use the notation AT to denote the transpose. The i–th eigenvalue of A ∈ Rn is denoted by λi (A). For a symmetric matrix A ∈ Rn×n , A 0 (A 0) denotes that A is positive semidefinite (negative semidefinite), that is, λ i (A) ≥ 0, ∀ i (λi (A) ≤ 0, ∀ i). For vectors, x, y ∈ Rn , x < y and x ≤ y denote element-wise inequalities. The n × n identity matrix is denoted as In and the n × m null matrix is denoted as 0n,m with the subscripts dropped if they can be inferred from context.
4908
ThBIn1.10 Bf
II. FDI PROBLEM SETTING
Consider a linear time-invariant (LTI) discrete-time system subject to disturbances and process, actuator and sensor faults of the form xk+1
=
Axk + Bd dk + Bf fk + Bu uk ,
yk
=
Cxk + Dd dk + Df fk + Du uk ,
(1)
for k ∈ N , where N := {0, 1, . . . , N − 1} is the estimation horizon, xk ∈ Rn , uk ∈ Rnu and yk ∈ Rny are the state, input and output vectors, respectively, and d k ∈ Rnd and fk ∈ Rnf are the disturbance and fault vectors, respectively. Here, Bf ∈ Rn×nf and Df ∈ Rny ×nf are the component and instrument fault distribution matrices, respectively, while Bd ∈ Rn×nd and Dd ∈ Rny ×nd are the corresponding disturbance distribution matrices [3], [5]. For simplicity, model uncertainties are assumed to be recast as disturbances [12]. Assume that upper and lower bounds x0 , x0 and dk , dk on the initial state and disturbances, respectively, are available and have of the form x0 ≤ x0 ≤ x0 ,
uN −1
d0 d1 d= . .. dN −1 x0 x1 x= . ..
yN −1
, d =
d0 d1 .. .
dN −1
xN −1
dN −1
f0 f1 , f = .. . fN −1
.
,
x=
AN −1
+
x0 + 0 Bd .. .
AN −2 Bd
0 Bu .. . AN −2 Bu
0 0 .. . AN −3 Bu
Bd
0 0 .. . AN −3 Bd
Dyu
+
+
C CA .. .
··· ··· .. .
0 Du .. .
Du CBu .. .
0 0 .. .
x0 + u N −2 N −3 CA Bu CA Bu · · · Du
CAN −1
Dyd
Dd CBd .. .
0 Dd .. .
CAN −2 Bd
CAN −3 Bd
··· ··· .. .
0 0 .. .
· · · Dd
Dyf
Df CBf .. .
0 Df .. .
CAN −2 Bf
CAN −3 Bf
··· ··· .. .
d
0 0 .. .
· · · Df
f
(3)
Thus, the output over the estimation horizon is given as y
= Cω ω + Dyf f + Dyu u,
(4)
where ω ∈ W := {ω : ω ≤ ω ≤ ω}.
,
III. FAULT SIGNAL BOUNDS
Bu
AN −3 Bf
To simplify the notation, define
x0 x0 x0 ,ω = . Cω = Cy0 Dyd , ω = ,ω = d d d
(2)
A
I A .. .
y =
and yk are
Using an iterative computation, the system algebraic formulation is described as
AN −2 Bf
··· 0 ··· 0 .. f .. . . ··· 0
Cy0
dk ≤ dk ≤ dk .
Assume also that the input and output signals u k available for all k ∈ N . For the estimation horizon, let u0 y0 d0 u1 y1 d1 u = . , y = . , d = . .. .. ..
+
0 0 .. .
0 Bf .. .
··· 0 ··· 0 . u .. . .. ··· 0
As mentioned in Section I, fault detection in this work is achieved by obtaining upper and lower bounds on the fault signal over the time horizon. This section gives the formulation of the problem of obtaining upper and lower bounds. Then, an approach for computing these bounds is proposed by solving an LP optimization. The cases for upper and lower bounds are dealt with separately. Firstly, we describe the problem as follows. Problem 3.1: Let the system dynamics be as defined in (1). Assume that the upper and lower bounds on x 0 and d are known. Find upper and lower bounds on f , f and f defined by
··· 0 ··· 0 . d .. . .. ··· 0
4909
min fi f i ≥ eTi f ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given max fi f i ≤ eTi f ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given
(5)
(6)
ThBIn1.10 where, f i and f i are the upper and lower bounds on e Ti f , respectively, ei is a vector with only the i th element 1 and others zero, and where i = 1, 2, . . . , N × n f . The following theorem shows that evaluating the upper bound on each f i can be obtained by solving an LP problem. Note that we use boldface to denote variables in the optimization. Theorem 3.1: Let all data be as defined in Section II. Let {v ∈ Rm : v ≥ 0}. Then eTi f ≤ f i for i = Rm + := ×nd 1, 2, . . . , N × nf , if there exist µi , µi ∈ Rn+N and + N ×ny such that µi ∈ R µi − µi + CωT µi = 0 T Dyf µi = e i (7) T f i ≥ µi ω + µi T (y − Dyu u) − µi T ω Proof. A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , eTi f − f i
=
−µi T (ω − ω) − µi T (ω − ω) T
−µi (y − Cω ω − Dyf f − Dyu u)
− ωT f T
ω 1 Li (f i , µi , µi , µi ) f 1
Proof. A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , f i − eTi f
0 0 1 T 2 µi Dyf
⋆ ⋆ T T T µi ω − µi ω − µi (y − Dyu u) + f i
(8)
where ⋆ denotes terms readily inferred from symmetry. It follows that a sufficient condition for e Ti f ≤ f i , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is Li (f i , µi , µi , µi ) 0
µi T (ω − ω) + µi T (ω − ω)
+ ωT f T
(9)
That the condition is necessary follows from Farkas Lemma [8]. Here, it can be seen that evaluating f i subject to (9) is equivalent to solve the following LP problem: µi
T T T µ (10) fi = min ω −ω (y −Dyu u) i T µi = ei Dyf µi µi −µi +CωT µi = 0 µi ≥ 0, µi ≥ 0 Therefore, the upper bound on the fault can be computed by solving (10). 2 Similarly, the theorem below is proposed to obtain the lower bound. Theorem 3.2: Let all data be as defined in Section II. Let Rm {v ∈ Rm : v ≤ 0}. Then f i ≤ eTi f for − :=
ω 1 Li (f i , µi , µi , µi ) f 1
where Li (f i , µi , µi , µi ) = 0 0 T 1 1 1 T T µ µ − 2 i 2 i + 2 µi Cω
0 0 1 T 2 µi Dyf
− 21 eTi
⋆ ⋆ T T T µi ω − µi ω − µi (y − Dyu u) + f i
(12)
where ⋆ denotes terms readily inferred from symmetry. It indicates that a sufficient condition for f i ≤ eTi f , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is µi , µi ≤ 0,
− 21 eTi
=
+µi T (y − Cω ω − Dyf f − Dyu u)
where Li (f i , µi , µi , µi ) = 0 0 T 1 1 1 T T µ µ − 2 i 2 i + 2 µi Cω
µi , µi ≥ 0,
i = 1, 2, . . . , N × nf , if and only if there exist µ i , µi ∈ ×nd Rn+N and µi ∈ RN ×ny such that − µi − µi + CωT µi = 0 T µi = e i Dyf (11) T f i ≤ µi ω + µi T (y − Dyu u) − µi T ω
Li (f i , µi , µi , µi ) 0
(13)
And the necessary condition follows from Farkas Lemma [8]. Here, it is noted that evaluating f i subject to (13) is equivalent to solve the following LP problem: µi
T T T µ −ω max (y −D u) (14) fi = ω yu i T µi = e i Dyf µi µi −µi +CωT µi = 0 µi ≥ 0, µi ≥ 0 Hence, the lower bound on the fault can be computed by solving (14). 2 Remark 3.1: From the proofs of Theorems 3.1 and 3.2, it can be seen that obtaining the tightest bounds on f is equivalent to solving 2(N × n f ) LP problems. Moreover, since the LP problems are relatively easy to solve, the propose method can be implemented online. IV. S UBSEQUENT STATE ESTIMATION As shown in Section II, the initial state x 0 is assumed to be in a given interval. After obtaining the bounds of each f i within the estimation horizon from 0 to N − 1, it is required to predict the interval of the next state x 1 to carry on the proposed method for the next horizon which is from 1 to N . In this section, a modified version of the method in the last section is given to obtain upper bound and lower bounds on the subsequent state.
4910
ThBIn1.10 Let x1i represent the i th element of x1 and x1i , x1i stand for tightest upper and lower bounds on x 1i , where i = 1, 2, . . . , n, given the model, the measurements over the estimation horizon and the apriori bounds on x 0 and d. Define
Bu ABu . . . AN −2 Bu AN −1 Bu , Buu =
Bd ABd . . . AN −2 Bd AN −1 Bd , Bud =
Bf ABf . . . AN −2 Bf AN −1 Bf . Buf = Inu 0 . . . 0 Ind 0 . . . 0 0 0 ... 0 0 0 ... 0 Iu = . .. . . .. , Id = .. .. . . .. , .. . . . . . . . 0 0 ... 0 0 0 ... 0 Inf 0 . . . 0 0 0 ... 0 If = . .. . . .. . .. . . . 0
0 ... 0
Next, we give the problem formulation for obtaining upper and lower bounds on x 1 . Problem 4.1: Let the system dynamics be as defined in Section II. Assume that the upper and lower bounds on x 0 and d are known. Find upper and lower bounds x 1i and x1i , for i = 1, 2, . . . , n, defined by x1i , (15) min x1i ≥ eTi x1 ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given max (16) x1i . x1i ≤ eTi x1 ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given Next, the upper and lower bounds on x 1i are obtained separately. Firstly, the upper bound on each x 1i can be evaluated by solving the LP problem in the following theorem. Theorem 4.1: Let all data be as defined in Section II. Then eTi x1 ≤ x1i for i = 1, 2, . . . , n, if and only if there exist ×nd and µi ∈ RN ×ny such that µi , µi ∈ Rn+N +
T µi − µi + CωT µi − A Bud Id ei = 0 T T (17) Dyf µi = If Buf ei x ≥ µ Tω+µ T(y −D u)− µ Tω+eT B I u 1i i yu uu u i i i .
2
where Lxi (x1i , µi , µi , µi ) = 0 0
T 1 T T T +µ C (µ −µ i ω−ei A Bud Id ) i 2 i
⋆ ⋆ µi Tω−µi Tω−µi T(y −Dyu u) −eTiBuuIu u+x1i
µi , µi ≥ 0,
x1i = uT Iu Buu ei + min T T Dyf µi = If Buf ei
T T µi−µi+Cω µi= A Bud Id ei µi ≥ 0, µi ≥ 0
.
2 The proof is omitted since it can be obtained just by following the lines in the proof of Theorem 4.1. Similarly as given in (21), the lower bound on x 1i can be computed by solving the following LP problem:
x1i = uT Iu Buu ei + T T µi = If Buf ei Dyf
T T µi −µi+Cω µi= A Bud Id ei µi ≥ 0, µi ≥ 0
(18)
−µi T (y − Cω ω − Dyf f − Dyu u)
− ωT f T
ω
1 Lxi (x1i , µi , µi , µi ) f 1
µi
ωT −ωT (y −Dyuu)T µi µi
Therefore, the upper bound on the fault can be computed by solving (21). 2 Similarly, it is proposed the following theorem which solves the lower bound on each x 1i . Theorem 4.2: Let all data be as defined in Section II. Then x1i ≤ eTi x1 for i = 1, 2, . . . , n, if and only if there exist ×nd µi , µi ∈ Rn+N and µi ∈ RN ×ny such that −
T µi − µi + CωT µi − A Bud Id ei = 0 T T (22) Dyf µi = If Buf ei x ≤ µ Tω+µ T(y −D u)− µ Tω+eT B I u i yu i 1i i uu u i
max
= −µi T (ω − ω) − µi T (ω − ω)
(21)
A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , eTi x1 − x1i
Lxi (x1i , µi , µi , µi ) ≥ 0 (20)
That the condition is necessary follows from Farkas Lemma [8]. Note that evaluating x 1i subject to (20) is equivalent to solve
= eTi Ax0 + eTi Buu Iu u + eTi Bud Id d +eTi Buf If f − x1i
(19)
where ⋆ denotes terms readily inferred from symmetry. It follows that a sufficient condition for e Ti x1 ≤ x1i , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is
Proof. From (1), it can be seen that eTi x − x1i
0 0 1 T T (µ D −e i yf i Buf If ) 2
ωT −ω T
µi
(y −Dyuu)T µi µi
(23) Remark 4.1: Note that in evaluating x 1 and x1 , we use all the information available over the whole estimation horizon. Remark 4.2: From the demonstration of the approach above, it can be seen that there is no filter or observer design involved, fault detection is achieved only by computing the
4911
ThBIn1.10 upper and lower bounds on the fault signal. This means that observer design is possibly no longer necessary and fault detection and isolation can be implemented on-line based on the input and output data. 2 Remark 4.3: Traditional observer-based FDI schemes only use the current state and data to generate the residual which reflects the fault information. However, the current state or data can be influenced by system parameter changes or other perturbations, which possibly causes a false alarm. In the proposed approach, the fault information is analyzed within the estimation horizon N , which provides more evidence to validate the existence of the fault. And the horizon N can be chosen as needed according to control process speed and available computing resources. 2 Remark 4.4: Although the bounds on f can be directly obtained from solving the primal Problem 3.1, we chose to present the bounds in terms of the solution of the dual problem in Theorems 3.1 and 3.2, since this format is more convenient for subsequent analysis of model uncertainty, which involves quadratic costs and constraints, which is under active research by the authors. 2 As demonstrated above, the upper and lower bounds on faults and subsequent state can simply be obtained by solving LP problems. Next, a procedure is given to illustrate the steps by which the proposed scheme is applied to real systems.
0.2666 0.0365 0.2666 0.0365 1.7629 −3.2664 , Bf = 1.7629 −3.2664 , Bu = −2.3152 1.7209 −2.3152 1.7209 −0.6083 0.4660 −0.6083 0.4660 1 0 0 0 −0.0069 0.0026 0 1 0 0 −0.0688 0.3896 C = 0 0 1 0 , Bd = 0.4433 0.0358 , 0 1 1 1 0.1144 0.0063 0 0.2 0 0.1 , Dd = Du = 04×2 , Df = 04×2 . 0.3 0 0 0 The same faults and disturbances are used as given in T T [14]: f (k)= 1.2 0.8 , for k ≥ 20 and f (k)= 0 0
T elsewhere. Also, d(k) = 0.5 cos(k)e−0.05k 0.5 sin(k) for all k. Generally, FDI schemes give better results as the number of outputs increases. However, to demonstrate the effectiveness of the proposed method, the system model is modified by removing the first output. Taking N = 10 and using the procedure given above, we obtained the following figures of upper and lower bounds on the faults from the sampling instant 12 to 21 and 15 to 24. From the graphs 1.5 Lower bound Fault1 Upper bound
1.2
1
Magnitude
1) Set i = 0. Obtain the input and output data u i and yi . Compute the upper and lower bounds on the potential fault, f0 , . . . , fi . 2) Take i = i + 1, compute the upper and lower bounds on fi until i = N − 1. 3) Set j = 1. 4) Compute upper lower bounds on x j . 5) Obtain uj+N −1 , and yj+N −1 , compute upper and lower bounds on f j , fj+1 , . . . , fj+N −1 , and xj+1 . 6) Set j = j + 1 and go to Step 4.
0.5
0
−0.5 12
13
14
15
16
17
18
19
20
21
24
25
20
21
Time
Fig. 1.
Fault1 from the instant 12 to 21
1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
Steps 1 and 2 correspond to an increasing horizon until N is reached. Remark 4.5: Although we have presented our procedure for linear time-invariant systems, they can be easily modified for linear time varying systems.
0.5
0
−0.5 16
17
18
19
20
21
22
23
Time
Fig. 2.
Fault1 from the instant 16 to 25
1.5 Lower bound Fault2 Upper bound 1
V. N UMERICAL EXAMPLES To illustrate the effectiveness of the proposed scheme, we consider a discrete-time model of a vertical takeoff and landing aircraft in the vertical plane used in [14], where the states are the horizontal velocity (knot), vertical velocity (knot), pitch rate (degree/second) and pitch angle (degree), respectively. The actuator inputs are collective pitch control and longitudinal cyclic pitch control, respectively. The parameters of the linearized model are given as 0.9828 0.0083 −0.0454 −0.2461 0.0117 0.5813 −0.3898 −1.6662 , A = 0.0458 0.1274 0.8230 0.4803 0.0117 0.0358 0.4433 1.1361
Magnitude
0.8
0.5
0
−0.5 12
13
14
Fig. 3.
15
16
17
18
19
Time
Fault2 from the instant 12 to 21
shown in Figure 1, 2, 3 and 4 and the values of the bounds, it can be seen at the 21st sampling instant that both the upper and lower bounds were greater than zero at the 20th sampling instant, which indicates fault occurrence in the system. The band between the upper and lower bounds on the fault is narrow and the fault is well-located in it. Note that since, in this example Df = 0, the bounds on the faults at the last sampling instant of each estimation horizon are not reliable.
4912
ThBIn1.10 1.5
1.5 Lower bound Fault2 Upper bound
1
1
0.8
0.8
Magnitude
Magnitude
Lower bound Fault2 Upper bound
0.5
0
−0.5 16
0.5
0
17
18
19
20
21
22
23
24
−0.5 16
25
17
18
Time
Fig. 4.
Fig. 8.
Also, using Theorem 4.1 and 4.2, it can be shown that state x16 is within the lower and upper bounds:
T x16 = −0.9649 0.7665 −0.1595 −0.3022
T x16 = 0.2508 0.8626 −0.1283 −0.2579
T x16 = 1.4179 0.8627 −0.0081 −0.2064
Furthermore, the proposed method can give even better result when it is applied using all four outputs as in [14]. Figures 5, 6, 7 and 8 illustrate the FDI results with the same faults and disturbances. It shows that both the upper and lower bounds are tight so that the faults are accurately detected at the 21st sampling instant. 1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
0.5
0
−0.5 12
13
14
15
16
17
18
19
20
21
24
25
20
21
Time
Fig. 5.
Fault1 from the instant 12 to 21
1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
0.5
0
−0.5 16
17
18
19
20
21
22
23
Time
Fig. 6.
Fault1 from the instant 16 to 25
1.5 Lower bound Fault2 Upper bound
1
Magnitude
0.8
0.5
0
−0.5 12
13
14
15
16
17
18
19
Time
Fig. 7.
19
20
21
22
23
24
25
Time
Fault2 from the instant 16 to 25
Fault2 from the instant 12 to 21
VI. S UMMARY In this work, we propose a new approach to fault detection and isolation for linear discrete-time systems subject to faults and bounded additive disturbances using a system model, together with input/output measurements over a finite estimation horizon. First, we use the proposed method to compute upper and lower bounds on the faults. These bounds
Fault2 from the instant 16 to 25
are the tightest possible based on the system model, output measurements and the bounds on the initial state and disturbances. Linear programming optimization techniques are used to obtain the bounds. Furthermore, a subsequent-stateestimation technique, together with an estimation horizon update procedure are given to allow the fault detection and isolation process to be repeated in an iterative procedure. Finally, a numerical example demonstrates the effectiveness of the proposed scheme. R EFERENCES [1] J. Chen, P.R. Patton, and H.Y. Zhang. Design of unknown input observers and robust fault detection filters. Int. J. Control, 63(1):85– 105, 1996. [2] R.K. Douglas and J.L. Speyer. An H∞ bounded fault detection filter. In Proc. Amer. Control Conf., pages 86–90, Seattle, Washington, June 1995. IEEE Press, New York. [3] P.M. Frank. Fault diagnosis in dynamic systems using analytical and knowledge–based redundancy: A survey and some new results. Automatica, 26(3):459–474, 1990. [4] P.M. Frank and X. Ding. Frequency domain approach to optimally robust residual generation and evaluation for model–based fault diagnosis. Automatica, 30(5):789–804, 1994. [5] P.M. Frank and X. Ding. Survey of robust residual generation and evaluation methods in observer–based fault detection systems. J. Proc. Cont., 7(6):403–424, 1997. [6] M. Hou and R.J. Patton. An LMI approach to H− /H∞ fault detection observers. In UKACC International Conference on Control. ’96, pages 305–310, England, 1996. [7] R. Isermann. Model-based fault-detection and diagnosis - status and applications. Annual Reviews in Control, 29:71–85, 2005. [8] R. Jagannathan and S. Schaible. Duality in generalized fractional programming via farka’s lemma. Journal of optimization theory and applications, 41(3), 1983. [9] I. M. Jaimoukha, Z. Li, and E. Mazars. Fualt isolation filter with linear matrix inequality solution to optimal decoupling. In Proc. Amer. Control Conf., pages 2339–2344. IEEE Press, New York, June 2006. Minneapolis, Minnesota. [10] I.M. Jaimoukha, Z. Li, and V. Papakos. A matrix factorization solution to the H− /H∞ fault detection problem. Automatica, 42(11):1907– 1913, 2006. [11] Z. Li and I. M. Jaimoukha. Observer-based fault detection and isolation filter design for linear time-invariant systems. Int. J. Control, 82:171–182, 2009. [12] R.J. Patton and J. Chen. On eigenstructure assignment for robust fault diagnosis. Int. J. Robust & Nonlinear Control, 10(14):1193–1208, 2000. [13] R.J. Patton, P.M. Frank, and R.N. Clark (eds). Issues of Fault Diagnosis for Dynamic Systems. Springer, 2000. [14] H. Wang and G. H. Yang. A finite frequency domain approach to fault detection for linear discrete-time systems. Int. J. Control, 81:1162– 1171, 2008. [15] P. Zhang and S. X. Ding. Distribution decoupling in fault detection of linear periodic systems. Automatica, 43:1410–1417, 2007. [16] M. Zhong, S.X. Ding, J. Lam, and H. Wang. An LMI approach to design robust fault detection filter for uncertain LTI systems. Automatica, 39(3):543–550, 2003.
4913
ThBIn1.10
Fault detection and isolation for linear discrete-time systems using input/output measurement analysis Ze Zhang, Imad M. Jaimoukha Abstract— In this work, a new approach to fault detection and isolation is proposed for linear discrete-time systems subject to faults and bounded additive disturbances using a system model, together with input/output measurements over a finite estimation horizon. We use the proposed method to compute upper and lower bounds on the faults. These bounds are the tightest possible given the system model, input/output measurements and the bounds on the initial state and disturbances. Linear programming optimization techniques are used to obtain the bounds. Moreover, a subsequent-state-estimation technique, together with an estimation horizon update procedure are given to allow the fault detection and isolation process to be repeated in an iterative procedure. Finally, the approach is verified using a numerical example.
I. I NTRODUCTION Model–based fault detection and isolation (FDI) for complex dynamic systems has received extensive attention both in literature and practice due to reliability, security and fault tolerance considerations. One of the most common approaches, the observer-based FDI approach, designs a dedicated observer to generate residual signals to provide fault signatures. See [13], [10], [7], [9], [11] for more details. The observer effectively cancels the (nominal) process dynamics and is sensitive only to disturbances, plant/ model mismatch (often recast as disturbances) and faults. The filter design objective is then to reduce the sensitivity to disturbances and/or plant/model mismatch for a given sensitivity to faults. There are two main approaches for achieving this objective, namely, perfect and approximate disturbance decoupling. The former aims to decouple the residual signal from disturbances exactly, while in the latter, the transfer matrix from disturbances to the residual signal is required to be small in either the H2 or H∞ norm sense. Moreover, in either case, the transfer matrix function from faults to residual is required to be diagonal for fault isolation. Both approaches have been given much attention and implemented using several methods. Patton and Chen [12] proposed left and right eigenvalue assignment method. Zhang and Ding [15] used a periodic parity space for linear discrete-time system. Chen, Patton and Zhang [1] solved the robust FDI problem by using an unknow input observer with disturbances decoupled in state estimation error. See also [4], [2] for other approaches. By using these methods, the decoupling problems can be The authors are with the Control and Power Group, Department of Electrical and Electronic Engineering, Imperial College, London. Email: [email protected] and [email protected]. Corresponding author: Ze Zhang, Control and Powe Group, Department of Electrical and Electronic Engineering, Imperial College, London SW7 2BT, UK. Fax: +44 207 594 6282. Telephone: +44 207 594 6279. Email: [email protected]
978-1-4244-3872-3/09/$25.00 ©2009 IEEE
transformed to a sensitivity optimization problem, which seeks to increase the sensitivity of the residual to faults and simultaneously reduces the sensitivity to disturbances and plant/model mismatch. Recently, linear matrix inequality (LMI) approaches have received increasing attention since they provide numerically attractive techniques for formulating robust FDI decoupling problems. Hou and Patton give a realization of fault detection observer design based on the bounded real lemma [6]. Furthermore, a new performance index which introduces a reference residual model is formulated by Zhong, Ding, Lam and Wang using LMI techniques in [16]. As illustrated above, most of existing FDI schemes propose generating a residual that attempts to reproduce the fault signal. In this work, we propose a new approach to FDI for linear discrete-time systems subject to bounded additive disturbances and faults using a dynamic system model as well as input/output measurements over a finite estimation horizon. By using the proposed method, upper and lower bounds on the faults are computed. Linear Programming (LP) optimization techniques are employed to obtain the bounds. Furthermore, we also propose a subsequent-stateestimation technique together with an estimation horizon update procedure which allow the on-line fault detection process to be operated in an iterative procedure. This work is organized as follows. After defining our notations, Section II defines the problem settings. Section III gives the problem formulation and derives upper and lower bounds on the fault signal in the form of solutions to linear programming problem. In Section IV, the method for subsequent-state-estimation as well as an estimation horizon update procedure is described, which complete the description of our fault detection algorithm. Examples are given to demonstrate the effectiveness of the proposed scheme in Section V. Finally, Section VI summarizes this work. The notation we use is fairly standard. The set of real n×m matrices is denoted by R n×m . For A ∈ Rn×m we use the notation AT to denote the transpose. The i–th eigenvalue of A ∈ Rn is denoted by λi (A). For a symmetric matrix A ∈ Rn×n , A 0 (A 0) denotes that A is positive semidefinite (negative semidefinite), that is, λ i (A) ≥ 0, ∀ i (λi (A) ≤ 0, ∀ i). For vectors, x, y ∈ Rn , x < y and x ≤ y denote element-wise inequalities. The n × n identity matrix is denoted as In and the n × m null matrix is denoted as 0n,m with the subscripts dropped if they can be inferred from context.
4908
ThBIn1.10 Bf
II. FDI PROBLEM SETTING
Consider a linear time-invariant (LTI) discrete-time system subject to disturbances and process, actuator and sensor faults of the form xk+1
=
Axk + Bd dk + Bf fk + Bu uk ,
yk
=
Cxk + Dd dk + Df fk + Du uk ,
(1)
for k ∈ N , where N := {0, 1, . . . , N − 1} is the estimation horizon, xk ∈ Rn , uk ∈ Rnu and yk ∈ Rny are the state, input and output vectors, respectively, and d k ∈ Rnd and fk ∈ Rnf are the disturbance and fault vectors, respectively. Here, Bf ∈ Rn×nf and Df ∈ Rny ×nf are the component and instrument fault distribution matrices, respectively, while Bd ∈ Rn×nd and Dd ∈ Rny ×nd are the corresponding disturbance distribution matrices [3], [5]. For simplicity, model uncertainties are assumed to be recast as disturbances [12]. Assume that upper and lower bounds x0 , x0 and dk , dk on the initial state and disturbances, respectively, are available and have of the form x0 ≤ x0 ≤ x0 ,
uN −1
d0 d1 d= . .. dN −1 x0 x1 x= . ..
yN −1
, d =
d0 d1 .. .
dN −1
xN −1
dN −1
f0 f1 , f = .. . fN −1
.
,
x=
AN −1
+
x0 + 0 Bd .. .
AN −2 Bd
0 Bu .. . AN −2 Bu
0 0 .. . AN −3 Bu
Bd
0 0 .. . AN −3 Bd
Dyu
+
+
C CA .. .
··· ··· .. .
0 Du .. .
Du CBu .. .
0 0 .. .
x0 + u N −2 N −3 CA Bu CA Bu · · · Du
CAN −1
Dyd
Dd CBd .. .
0 Dd .. .
CAN −2 Bd
CAN −3 Bd
··· ··· .. .
0 0 .. .
· · · Dd
Dyf
Df CBf .. .
0 Df .. .
CAN −2 Bf
CAN −3 Bf
··· ··· .. .
d
0 0 .. .
· · · Df
f
(3)
Thus, the output over the estimation horizon is given as y
= Cω ω + Dyf f + Dyu u,
(4)
where ω ∈ W := {ω : ω ≤ ω ≤ ω}.
,
III. FAULT SIGNAL BOUNDS
Bu
AN −3 Bf
To simplify the notation, define
x0 x0 x0 ,ω = . Cω = Cy0 Dyd , ω = ,ω = d d d
(2)
A
I A .. .
y =
and yk are
Using an iterative computation, the system algebraic formulation is described as
AN −2 Bf
··· 0 ··· 0 .. f .. . . ··· 0
Cy0
dk ≤ dk ≤ dk .
Assume also that the input and output signals u k available for all k ∈ N . For the estimation horizon, let u0 y0 d0 u1 y1 d1 u = . , y = . , d = . .. .. ..
+
0 0 .. .
0 Bf .. .
··· 0 ··· 0 . u .. . .. ··· 0
As mentioned in Section I, fault detection in this work is achieved by obtaining upper and lower bounds on the fault signal over the time horizon. This section gives the formulation of the problem of obtaining upper and lower bounds. Then, an approach for computing these bounds is proposed by solving an LP optimization. The cases for upper and lower bounds are dealt with separately. Firstly, we describe the problem as follows. Problem 3.1: Let the system dynamics be as defined in (1). Assume that the upper and lower bounds on x 0 and d are known. Find upper and lower bounds on f , f and f defined by
··· 0 ··· 0 . d .. . .. ··· 0
4909
min fi f i ≥ eTi f ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given max fi f i ≤ eTi f ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given
(5)
(6)
ThBIn1.10 where, f i and f i are the upper and lower bounds on e Ti f , respectively, ei is a vector with only the i th element 1 and others zero, and where i = 1, 2, . . . , N × n f . The following theorem shows that evaluating the upper bound on each f i can be obtained by solving an LP problem. Note that we use boldface to denote variables in the optimization. Theorem 3.1: Let all data be as defined in Section II. Let {v ∈ Rm : v ≥ 0}. Then eTi f ≤ f i for i = Rm + := ×nd 1, 2, . . . , N × nf , if there exist µi , µi ∈ Rn+N and + N ×ny such that µi ∈ R µi − µi + CωT µi = 0 T Dyf µi = e i (7) T f i ≥ µi ω + µi T (y − Dyu u) − µi T ω Proof. A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , eTi f − f i
=
−µi T (ω − ω) − µi T (ω − ω) T
−µi (y − Cω ω − Dyf f − Dyu u)
− ωT f T
ω 1 Li (f i , µi , µi , µi ) f 1
Proof. A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , f i − eTi f
0 0 1 T 2 µi Dyf
⋆ ⋆ T T T µi ω − µi ω − µi (y − Dyu u) + f i
(8)
where ⋆ denotes terms readily inferred from symmetry. It follows that a sufficient condition for e Ti f ≤ f i , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is Li (f i , µi , µi , µi ) 0
µi T (ω − ω) + µi T (ω − ω)
+ ωT f T
(9)
That the condition is necessary follows from Farkas Lemma [8]. Here, it can be seen that evaluating f i subject to (9) is equivalent to solve the following LP problem: µi
T T T µ (10) fi = min ω −ω (y −Dyu u) i T µi = ei Dyf µi µi −µi +CωT µi = 0 µi ≥ 0, µi ≥ 0 Therefore, the upper bound on the fault can be computed by solving (10). 2 Similarly, the theorem below is proposed to obtain the lower bound. Theorem 3.2: Let all data be as defined in Section II. Let Rm {v ∈ Rm : v ≤ 0}. Then f i ≤ eTi f for − :=
ω 1 Li (f i , µi , µi , µi ) f 1
where Li (f i , µi , µi , µi ) = 0 0 T 1 1 1 T T µ µ − 2 i 2 i + 2 µi Cω
0 0 1 T 2 µi Dyf
− 21 eTi
⋆ ⋆ T T T µi ω − µi ω − µi (y − Dyu u) + f i
(12)
where ⋆ denotes terms readily inferred from symmetry. It indicates that a sufficient condition for f i ≤ eTi f , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is µi , µi ≤ 0,
− 21 eTi
=
+µi T (y − Cω ω − Dyf f − Dyu u)
where Li (f i , µi , µi , µi ) = 0 0 T 1 1 1 T T µ µ − 2 i 2 i + 2 µi Cω
µi , µi ≥ 0,
i = 1, 2, . . . , N × nf , if and only if there exist µ i , µi ∈ ×nd Rn+N and µi ∈ RN ×ny such that − µi − µi + CωT µi = 0 T µi = e i Dyf (11) T f i ≤ µi ω + µi T (y − Dyu u) − µi T ω
Li (f i , µi , µi , µi ) 0
(13)
And the necessary condition follows from Farkas Lemma [8]. Here, it is noted that evaluating f i subject to (13) is equivalent to solve the following LP problem: µi
T T T µ −ω max (y −D u) (14) fi = ω yu i T µi = e i Dyf µi µi −µi +CωT µi = 0 µi ≥ 0, µi ≥ 0 Hence, the lower bound on the fault can be computed by solving (14). 2 Remark 3.1: From the proofs of Theorems 3.1 and 3.2, it can be seen that obtaining the tightest bounds on f is equivalent to solving 2(N × n f ) LP problems. Moreover, since the LP problems are relatively easy to solve, the propose method can be implemented online. IV. S UBSEQUENT STATE ESTIMATION As shown in Section II, the initial state x 0 is assumed to be in a given interval. After obtaining the bounds of each f i within the estimation horizon from 0 to N − 1, it is required to predict the interval of the next state x 1 to carry on the proposed method for the next horizon which is from 1 to N . In this section, a modified version of the method in the last section is given to obtain upper bound and lower bounds on the subsequent state.
4910
ThBIn1.10 Let x1i represent the i th element of x1 and x1i , x1i stand for tightest upper and lower bounds on x 1i , where i = 1, 2, . . . , n, given the model, the measurements over the estimation horizon and the apriori bounds on x 0 and d. Define
Bu ABu . . . AN −2 Bu AN −1 Bu , Buu =
Bd ABd . . . AN −2 Bd AN −1 Bd , Bud =
Bf ABf . . . AN −2 Bf AN −1 Bf . Buf = Inu 0 . . . 0 Ind 0 . . . 0 0 0 ... 0 0 0 ... 0 Iu = . .. . . .. , Id = .. .. . . .. , .. . . . . . . . 0 0 ... 0 0 0 ... 0 Inf 0 . . . 0 0 0 ... 0 If = . .. . . .. . .. . . . 0
0 ... 0
Next, we give the problem formulation for obtaining upper and lower bounds on x 1 . Problem 4.1: Let the system dynamics be as defined in Section II. Assume that the upper and lower bounds on x 0 and d are known. Find upper and lower bounds x 1i and x1i , for i = 1, 2, . . . , n, defined by x1i , (15) min x1i ≥ eTi x1 ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given max (16) x1i . x1i ≤ eTi x1 ω≤ω≤ω y = Cω ω + Dyf f + Dyu u y, u given Next, the upper and lower bounds on x 1i are obtained separately. Firstly, the upper bound on each x 1i can be evaluated by solving the LP problem in the following theorem. Theorem 4.1: Let all data be as defined in Section II. Then eTi x1 ≤ x1i for i = 1, 2, . . . , n, if and only if there exist ×nd and µi ∈ RN ×ny such that µi , µi ∈ Rn+N +
T µi − µi + CωT µi − A Bud Id ei = 0 T T (17) Dyf µi = If Buf ei x ≥ µ Tω+µ T(y −D u)− µ Tω+eT B I u 1i i yu uu u i i i .
2
where Lxi (x1i , µi , µi , µi ) = 0 0
T 1 T T T +µ C (µ −µ i ω−ei A Bud Id ) i 2 i
⋆ ⋆ µi Tω−µi Tω−µi T(y −Dyu u) −eTiBuuIu u+x1i
µi , µi ≥ 0,
x1i = uT Iu Buu ei + min T T Dyf µi = If Buf ei
T T µi−µi+Cω µi= A Bud Id ei µi ≥ 0, µi ≥ 0
.
2 The proof is omitted since it can be obtained just by following the lines in the proof of Theorem 4.1. Similarly as given in (21), the lower bound on x 1i can be computed by solving the following LP problem:
x1i = uT Iu Buu ei + T T µi = If Buf ei Dyf
T T µi −µi+Cω µi= A Bud Id ei µi ≥ 0, µi ≥ 0
(18)
−µi T (y − Cω ω − Dyf f − Dyu u)
− ωT f T
ω
1 Lxi (x1i , µi , µi , µi ) f 1
µi
ωT −ωT (y −Dyuu)T µi µi
Therefore, the upper bound on the fault can be computed by solving (21). 2 Similarly, it is proposed the following theorem which solves the lower bound on each x 1i . Theorem 4.2: Let all data be as defined in Section II. Then x1i ≤ eTi x1 for i = 1, 2, . . . , n, if and only if there exist ×nd µi , µi ∈ Rn+N and µi ∈ RN ×ny such that −
T µi − µi + CωT µi − A Bud Id ei = 0 T T (22) Dyf µi = If Buf ei x ≤ µ Tω+µ T(y −D u)− µ Tω+eT B I u i yu i 1i i uu u i
max
= −µi T (ω − ω) − µi T (ω − ω)
(21)
A manipulation verifies the following identity which is valid for all µi , µi ∈ Rn+N ×nd , µi ∈ RN ×ny , eTi x1 − x1i
Lxi (x1i , µi , µi , µi ) ≥ 0 (20)
That the condition is necessary follows from Farkas Lemma [8]. Note that evaluating x 1i subject to (20) is equivalent to solve
= eTi Ax0 + eTi Buu Iu u + eTi Bud Id d +eTi Buf If f − x1i
(19)
where ⋆ denotes terms readily inferred from symmetry. It follows that a sufficient condition for e Ti x1 ≤ x1i , ω ≤ ω ≤ ω and y = Cω ω + Dyf f + Dyu u is
Proof. From (1), it can be seen that eTi x − x1i
0 0 1 T T (µ D −e i yf i Buf If ) 2
ωT −ω T
µi
(y −Dyuu)T µi µi
(23) Remark 4.1: Note that in evaluating x 1 and x1 , we use all the information available over the whole estimation horizon. Remark 4.2: From the demonstration of the approach above, it can be seen that there is no filter or observer design involved, fault detection is achieved only by computing the
4911
ThBIn1.10 upper and lower bounds on the fault signal. This means that observer design is possibly no longer necessary and fault detection and isolation can be implemented on-line based on the input and output data. 2 Remark 4.3: Traditional observer-based FDI schemes only use the current state and data to generate the residual which reflects the fault information. However, the current state or data can be influenced by system parameter changes or other perturbations, which possibly causes a false alarm. In the proposed approach, the fault information is analyzed within the estimation horizon N , which provides more evidence to validate the existence of the fault. And the horizon N can be chosen as needed according to control process speed and available computing resources. 2 Remark 4.4: Although the bounds on f can be directly obtained from solving the primal Problem 3.1, we chose to present the bounds in terms of the solution of the dual problem in Theorems 3.1 and 3.2, since this format is more convenient for subsequent analysis of model uncertainty, which involves quadratic costs and constraints, which is under active research by the authors. 2 As demonstrated above, the upper and lower bounds on faults and subsequent state can simply be obtained by solving LP problems. Next, a procedure is given to illustrate the steps by which the proposed scheme is applied to real systems.
0.2666 0.0365 0.2666 0.0365 1.7629 −3.2664 , Bf = 1.7629 −3.2664 , Bu = −2.3152 1.7209 −2.3152 1.7209 −0.6083 0.4660 −0.6083 0.4660 1 0 0 0 −0.0069 0.0026 0 1 0 0 −0.0688 0.3896 C = 0 0 1 0 , Bd = 0.4433 0.0358 , 0 1 1 1 0.1144 0.0063 0 0.2 0 0.1 , Dd = Du = 04×2 , Df = 04×2 . 0.3 0 0 0 The same faults and disturbances are used as given in T T [14]: f (k)= 1.2 0.8 , for k ≥ 20 and f (k)= 0 0
T elsewhere. Also, d(k) = 0.5 cos(k)e−0.05k 0.5 sin(k) for all k. Generally, FDI schemes give better results as the number of outputs increases. However, to demonstrate the effectiveness of the proposed method, the system model is modified by removing the first output. Taking N = 10 and using the procedure given above, we obtained the following figures of upper and lower bounds on the faults from the sampling instant 12 to 21 and 15 to 24. From the graphs 1.5 Lower bound Fault1 Upper bound
1.2
1
Magnitude
1) Set i = 0. Obtain the input and output data u i and yi . Compute the upper and lower bounds on the potential fault, f0 , . . . , fi . 2) Take i = i + 1, compute the upper and lower bounds on fi until i = N − 1. 3) Set j = 1. 4) Compute upper lower bounds on x j . 5) Obtain uj+N −1 , and yj+N −1 , compute upper and lower bounds on f j , fj+1 , . . . , fj+N −1 , and xj+1 . 6) Set j = j + 1 and go to Step 4.
0.5
0
−0.5 12
13
14
15
16
17
18
19
20
21
24
25
20
21
Time
Fig. 1.
Fault1 from the instant 12 to 21
1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
Steps 1 and 2 correspond to an increasing horizon until N is reached. Remark 4.5: Although we have presented our procedure for linear time-invariant systems, they can be easily modified for linear time varying systems.
0.5
0
−0.5 16
17
18
19
20
21
22
23
Time
Fig. 2.
Fault1 from the instant 16 to 25
1.5 Lower bound Fault2 Upper bound 1
V. N UMERICAL EXAMPLES To illustrate the effectiveness of the proposed scheme, we consider a discrete-time model of a vertical takeoff and landing aircraft in the vertical plane used in [14], where the states are the horizontal velocity (knot), vertical velocity (knot), pitch rate (degree/second) and pitch angle (degree), respectively. The actuator inputs are collective pitch control and longitudinal cyclic pitch control, respectively. The parameters of the linearized model are given as 0.9828 0.0083 −0.0454 −0.2461 0.0117 0.5813 −0.3898 −1.6662 , A = 0.0458 0.1274 0.8230 0.4803 0.0117 0.0358 0.4433 1.1361
Magnitude
0.8
0.5
0
−0.5 12
13
14
Fig. 3.
15
16
17
18
19
Time
Fault2 from the instant 12 to 21
shown in Figure 1, 2, 3 and 4 and the values of the bounds, it can be seen at the 21st sampling instant that both the upper and lower bounds were greater than zero at the 20th sampling instant, which indicates fault occurrence in the system. The band between the upper and lower bounds on the fault is narrow and the fault is well-located in it. Note that since, in this example Df = 0, the bounds on the faults at the last sampling instant of each estimation horizon are not reliable.
4912
ThBIn1.10 1.5
1.5 Lower bound Fault2 Upper bound
1
1
0.8
0.8
Magnitude
Magnitude
Lower bound Fault2 Upper bound
0.5
0
−0.5 16
0.5
0
17
18
19
20
21
22
23
24
−0.5 16
25
17
18
Time
Fig. 4.
Fig. 8.
Also, using Theorem 4.1 and 4.2, it can be shown that state x16 is within the lower and upper bounds:
T x16 = −0.9649 0.7665 −0.1595 −0.3022
T x16 = 0.2508 0.8626 −0.1283 −0.2579
T x16 = 1.4179 0.8627 −0.0081 −0.2064
Furthermore, the proposed method can give even better result when it is applied using all four outputs as in [14]. Figures 5, 6, 7 and 8 illustrate the FDI results with the same faults and disturbances. It shows that both the upper and lower bounds are tight so that the faults are accurately detected at the 21st sampling instant. 1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
0.5
0
−0.5 12
13
14
15
16
17
18
19
20
21
24
25
20
21
Time
Fig. 5.
Fault1 from the instant 12 to 21
1.5
1.2
Lower bound Fault1 Upper bound
Magnitude
1
0.5
0
−0.5 16
17
18
19
20
21
22
23
Time
Fig. 6.
Fault1 from the instant 16 to 25
1.5 Lower bound Fault2 Upper bound
1
Magnitude
0.8
0.5
0
−0.5 12
13
14
15
16
17
18
19
Time
Fig. 7.
19
20
21
22
23
24
25
Time
Fault2 from the instant 16 to 25
Fault2 from the instant 12 to 21
VI. S UMMARY In this work, we propose a new approach to fault detection and isolation for linear discrete-time systems subject to faults and bounded additive disturbances using a system model, together with input/output measurements over a finite estimation horizon. First, we use the proposed method to compute upper and lower bounds on the faults. These bounds
Fault2 from the instant 16 to 25
are the tightest possible based on the system model, output measurements and the bounds on the initial state and disturbances. Linear programming optimization techniques are used to obtain the bounds. Furthermore, a subsequent-stateestimation technique, together with an estimation horizon update procedure are given to allow the fault detection and isolation process to be repeated in an iterative procedure. Finally, a numerical example demonstrates the effectiveness of the proposed scheme. R EFERENCES [1] J. Chen, P.R. Patton, and H.Y. Zhang. Design of unknown input observers and robust fault detection filters. Int. J. Control, 63(1):85– 105, 1996. [2] R.K. Douglas and J.L. Speyer. An H∞ bounded fault detection filter. In Proc. Amer. Control Conf., pages 86–90, Seattle, Washington, June 1995. IEEE Press, New York. [3] P.M. Frank. Fault diagnosis in dynamic systems using analytical and knowledge–based redundancy: A survey and some new results. Automatica, 26(3):459–474, 1990. [4] P.M. Frank and X. Ding. Frequency domain approach to optimally robust residual generation and evaluation for model–based fault diagnosis. Automatica, 30(5):789–804, 1994. [5] P.M. Frank and X. Ding. Survey of robust residual generation and evaluation methods in observer–based fault detection systems. J. Proc. Cont., 7(6):403–424, 1997. [6] M. Hou and R.J. Patton. An LMI approach to H− /H∞ fault detection observers. In UKACC International Conference on Control. ’96, pages 305–310, England, 1996. [7] R. Isermann. Model-based fault-detection and diagnosis - status and applications. Annual Reviews in Control, 29:71–85, 2005. [8] R. Jagannathan and S. Schaible. Duality in generalized fractional programming via farka’s lemma. Journal of optimization theory and applications, 41(3), 1983. [9] I. M. Jaimoukha, Z. Li, and E. Mazars. Fualt isolation filter with linear matrix inequality solution to optimal decoupling. In Proc. Amer. Control Conf., pages 2339–2344. IEEE Press, New York, June 2006. Minneapolis, Minnesota. [10] I.M. Jaimoukha, Z. Li, and V. Papakos. A matrix factorization solution to the H− /H∞ fault detection problem. Automatica, 42(11):1907– 1913, 2006. [11] Z. Li and I. M. Jaimoukha. Observer-based fault detection and isolation filter design for linear time-invariant systems. Int. J. Control, 82:171–182, 2009. [12] R.J. Patton and J. Chen. On eigenstructure assignment for robust fault diagnosis. Int. J. Robust & Nonlinear Control, 10(14):1193–1208, 2000. [13] R.J. Patton, P.M. Frank, and R.N. Clark (eds). Issues of Fault Diagnosis for Dynamic Systems. Springer, 2000. [14] H. Wang and G. H. Yang. A finite frequency domain approach to fault detection for linear discrete-time systems. Int. J. Control, 81:1162– 1171, 2008. [15] P. Zhang and S. X. Ding. Distribution decoupling in fault detection of linear periodic systems. Automatica, 43:1410–1417, 2007. [16] M. Zhong, S.X. Ding, J. Lam, and H. Wang. An LMI approach to design robust fault detection filter for uncertain LTI systems. Automatica, 39(3):543–550, 2003.
4913
