Journal of Control Science and Engineering

Journal of Control Science and Engineering / 2008 / Article
Special Issue

Robustness Issues in Fault Diagnosis and Fault Tolerant Control

View this Special Issue

Research Article | Open Access

Volume 2008 |Article ID 274313 |

M. Benini, M. Bonfè, P. Castaldi, W. Geri, S. Simani, "Design and Analysis of Robust Fault Diagnosis Schemes for a Simulated Aircraft Model", Journal of Control Science and Engineering, vol. 2008, Article ID 274313, 18 pages, 2008.

Design and Analysis of Robust Fault Diagnosis Schemes for a Simulated Aircraft Model

Academic Editor: Jakob Stoustrup
Received30 Mar 2007
Revised27 Sep 2007
Accepted20 Dec 2007
Published30 Mar 2008


Several procedures for sensor fault detection and isolation (FDI) applied to a simulated model of a commercial aircraft are presented. The main contributions of the paper are related to the design and the optimisation of two FDI schemes based on a linear polynomial method (PM) and the nonlinear geometric approach (NLGA). The FDI strategies are applied to the aircraft model, characterised by tight-coupled longitudinal and lateral dynamics. The robustness and the reliability properties of the residual generators related to the considered FDI techniques are investigated and verified by simulating a general aircraft reference trajectory. Extensive simulations exploiting the Monte Carlo analysis tool are also used for assessing the overall performance capabilities of the developed FDI schemes, in the presence of turbulence, measurement, and model errors. Comparisons with other disturbance-decoupling methods for FDI based on neural networks (NNs) and unknown input kalman filter (UIKF) are finally reported.

1. Introduction

Increasing demands on reliability for safety critical systems such as aircraft or spacecraft require robust control and fault diagnosis capabilities as these systems are potentially subjected to unexpected anomalies and faults in actuators, input-output sensors, components, or subsystems. Consequently, fault diagnosis capabilities and requirements for aerospace applications have recently been receiving a great deal of attention in the research community [1, 2]. A fault diagnosis system needs to detect and isolate the presence and location of the faults, on the basis also of the control system architectures. Development of appropriate techniques and solutions for these tasks are known as the fault detection and isolation (FDI) problem. There are, broadly speaking, two main approaches for addressing the FDI problem, namely, hardware-based and model-based techniques [3, 4]. A common and important approach in model-based techniques is known as the residual-based method. A number of researchers have developed residual-based methods for dynamic systems such as the parity space [5], state estimation [6], unknown input observer (UIO), Kalman filters (KFs) [3], and parameter identification [6]. Intelligent techniques [7] can be also exploited. Furthermore, the Massoumnia's geometric method [8] was successfully extended to nonlinear systems [9, 10]. A crucial issue with any FDI scheme is its robustness properties and a viable procedure for practical application of FDI techniques is really necessary. Moreover, robust FDI for the case of aircraft systems and applications is still an open problem for further research.

The first part of this work deals with the residual generator design for the FDI of input-output sensors of a general aviation aircraft subject to turbulence, wind gust disturbances, and measurement noises. The developed PM scheme belongs to the parity space approach [5] and it is based on an input-output polynomial description of the system under diagnosis. In particular, the use of input-output forms allows to easily obtain the analytical description for the disturbance-decoupled residual generators. These dynamic filters, organised into bank structures, are able to achieve fault isolation properties. An appropriate choice of their parameters allows to maximise robustness with respect to both measurement noise and modelling errors, while optimising fault sensitivity characteristics. The development of NLGA methodology is based on the works by De Persis and Isidori [10]. It was shown that the problem of the FDI for nonlinear systems is solvable if and only if there is an unobservability distribution that leads to a quotient subsystem which is unaffected by all faults but one. If such a distribution exists, an appropriate coordinate transformations in the state space can be exploited for designing a residual generator only for the observable subsystem. This technique was applied for the first time to a vertical takeoff and landing (VTOL) aircraft with reference to a reduced-order model [11]. The NLGA residual generators have been designed in order to be analytically decoupled from the vertical and lateral components of the wind (gusts and turbulence). Moreover, a new full analytical developed mixed / optimisation is proposed in order to design the NLGA residual generators so that a good tradeoff between the fault sensitivity and the robustness with respect to measurements and model errors is achieved. The designed residual generators have been tested on a PIPER PA-30 aircraft flight simulator that was implemented in Matlab-Simulink environments. With respect to the related works by the same authors [12, 13], the main contribution of this paper regards the enhancement in the designs of the proposed FDI schemes. Moreover, the final performances have been evaluated by adopting a typical aircraft reference trajectory embedding several steady-state flight conditions, such as straight flight phases and coordinated turns. Comparisons with different disturbance-decoupling methods for FDI based on neural networks (NNs) and unknown input Kalman filter (UIKF) have been also provided. Finally, extensive experiments exploiting Monte Carlo analysis are used for assessing the overall capabilities of the developed FDI methods, in the presence of uncertainty, measurement, and modelling errors.

2. Aircraft Model Overview

This section recalls briefly the description of the monitored aircraft whose main parameters and variables are reported in Table 1.

𝛼 Angle of attack
𝛽 Angle of sideslip
𝑝 𝜔 Roll rate
𝑞 𝜔 Pitch rate
𝑟 𝜔 Yaw rate
𝜙 Bank angle
𝜃 Elevation angle
𝜓 Heading angle
𝑛 𝑒 Engine shaft angular rate
[ 𝐼 𝑥 0 𝐼 𝑥 𝑧 0 𝐼 𝑦 0 𝐼 𝑥 𝑧 0 𝐼 𝑧 ] Inertia moment matrix
𝑉 True airspeed (TAS)
𝛿 𝑒 Elevator deflection angle
𝛿 𝑎 Aileron deflection angle
𝛿 𝑟 rudder deflection angle
𝛿 t h Throttle aperture percentage
𝐻 Altitude
𝛾 Flight path angle
𝑚 Airplane mass
𝜔 𝑢 , 𝜔 𝑣 , 𝜔 𝑤 Wind gust components

The considered aircraft simulation model consists of a PIPER PA-30, based on the classical nonlinear 6 degrees of freedom (DoF) rigid body formulation [14] whose motion occurs as a consequence of applied forces and moments (aerodynamic, propulsive, and gravitational). A set of local approximations for these forces has been computed and scheduled depending on the values assumed by true airspeed (TAS), curvature radius, flight path angle, altitude, and flap deflection. In this way, it is possible to obtain a mathematical model for each flight condition. This model is suitable for a state-space representation, as it can be made explicit. The parameters in the analytic representation of the aerodynamic actions have been obtained from wind tunnel experimental data, and the aerodynamic actions are expressed along the axes of the wind reference system. It should be observed that aerodynamic forces and moments are not implemented by the classical linearised expressions (stability derivatives) but by means of cubic splines approximating the nonlinear experimental curves. The nonlinear 6 DoF model has been completed by means of the PIPER PA-30 propulsion system consisting of two 4-pistons aspirated engines, with the throttle valve aperture 𝛿th as input and the overall thrust intensity as output. The overall simulation model, used to perform all the following tests, consists of the aircraft 6 DoF flight dynamics and the engine model completed with the model of input-output sensors, the servo actuators, the atmosphere turbulence Dryden description, the wind gust disturbances, and a classical autopilot. Moreover, the sensor models embed all the possible sources of disturbance (calibration and alignment errors, scale factor, white and coloured noises, limited bandwidth, g-sensitivity, gyro drift, etc.).

The linear model used by the proposed PM FDI approach described in Section 3 embeds the linearisation both of the 6 DoF model and of the propulsion system as follows: ̇𝑥(𝑡)=𝐴𝑥(𝑡)+𝐵𝑐(𝑡)+𝐸𝑑(𝑡)(1) with 𝑥(𝑡)=[Δ𝑉(𝑡)Δ𝛼(𝑡)Δ𝛽(𝑡)Δ𝑝𝜔(𝑡)Δ𝑞𝜔(𝑡)Δ𝑟𝜔(𝑡)Δ𝜙(𝑡)Δ𝜃(𝑡)Δ𝜓(𝑡)Δ𝑛𝑒(𝑡)]𝑇,𝑐(𝑡)=[Δ𝛿𝑒(𝑡)Δ𝛿𝑎(𝑡)Δ𝛿𝑟(𝑡)Δ𝛿th(𝑡)]𝑇,𝑑(𝑡)=[𝑤𝑢(𝑡)𝑤𝑣(𝑡)𝑤𝑤(𝑡)]𝑇,(2) where Δ denotes the variations of the considered variables while 𝑐(𝑡) and 𝑑(𝑡) are the control inputs and the disturbances, respectively. The disturbance contribution of the wind gusts as air velocity components, 𝑤𝑢, 𝑤𝑣, and 𝑤𝑤, along body axes was also considered. The output equation associated with the model (1) is of the type 𝑦(𝑡)=𝐶𝑥(𝑡), where the rows of 𝐶 correspond to rows of the identity matrix, depending on the measured variables.

On the other hand, regarding the NLGA FDI scheme described in Section 4, it requires a nonlinear input affine system [10], but the adopted simulation model of the aircraft does not fulfil this requirement. For this reason, the following simplified aircraft model is used: ̇𝑉=(𝐶𝐷0+𝐶𝐷𝛼𝛼+𝐶𝐷𝛼2𝛼2)𝑚𝑉2++𝑔sin𝛼cos𝜃cos𝜙cos𝛼sin𝜃cos𝛼𝑚𝑡𝑝𝑉𝑡0+𝑡1𝑛𝑒𝛿th+𝑤𝑣𝐶sin𝛼,̇𝛼=𝐿0+𝐶𝐿𝛼𝛼𝑚𝑔𝑉+𝑉cos𝛼cos𝜃cos𝜙+sin𝛼sin𝜃+𝑞𝜔+sin𝛼𝑚𝑡𝑝𝑉2𝑡0+𝑡1𝑛𝑒𝛿th+cos𝛼𝑉𝑤𝑣,̇𝛽=(𝐶𝐷0+𝐶𝐷𝛼𝛼+𝐶𝐷𝛼2𝛼2)sin𝛽+𝐶𝑌𝛽𝛽cos𝛽𝑚𝑉+𝑔cos𝜃sin𝜙𝑉+𝑝𝜔sin𝛼𝑟𝜔+cos𝛼cos𝛼sin𝛽𝑚𝑡𝑝𝑉2𝑡0+𝑡1𝑛𝑒𝛿th+1𝑉𝑤,̇𝑝𝜔=𝐶𝑙𝛽𝛽+𝐶𝑙𝑝𝑝𝜔𝐼𝑥𝑉2+𝐼𝑦𝐼𝑧𝐼𝑥𝑞𝜔𝑟𝜔+𝐶𝛿𝑎𝐼𝑥𝑉2𝛿𝑎,̇𝑞𝜔=𝐶𝑚0+𝐶𝑚𝛼𝛼+𝐶𝑚𝑞𝑞𝜔𝐼𝑦𝑉2+𝐼𝑧𝐼𝑥𝐼𝑦𝑝𝜔𝑟𝜔+𝐶𝛿𝑒𝐼𝑦𝑉2𝛿𝑒+𝑡𝑑𝐼𝑦𝑡𝑝𝑉𝑡0+𝑡1𝑛𝑒𝛿th,̇𝑟𝜔=𝐶𝑛𝛽𝛽+𝐶𝑛𝑟𝑟𝜔𝐼𝑧𝑉2+𝐼𝑥𝐼𝑦𝐼𝑧𝑝𝜔𝑞𝜔+𝐶𝛿𝑟𝐼𝑧𝑉2𝛿𝑟,̇𝜙=𝑝𝜔+𝑞𝜔sin𝜙+𝑟𝜔̇cos𝜙tan𝜃,𝜃=𝑞𝜔cos𝜙𝑟𝜔𝑞sin𝜙,̇𝜓=𝜔sin𝜙+𝑟𝜔cos𝜙cos𝜃̇𝑛𝑒=𝑡𝑛𝑛3𝑒+𝑡𝑓𝑛𝑒𝑡0+𝑡1𝑛𝑒𝛿th,(3) where 𝐶() are the aerodynamic coefficients; 𝑡() are the engine parameters; and 𝑤𝑣, 𝑤𝑙 are the vertical and lateral wind disturbance components. In particular, the model of (3) has been obtained on the basis of some assumptions. In particular, the expressions of aerodynamic forces and moments have been represented by means of series expansions in the neighbourhood of the steady-state flight condition, then only the main terms are considered. The engine model has been simplified by linearising the power with respect to the angular rate behaviour in the neighbourhood of the trim point. The second-order coupling between the longitudinal and lateral-directional dynamics have been neglected. The 𝑥-body axis component of the wind has been neglected. In fact, the aircraft behaviour is much more sensitive to the 𝑦-body and 𝑧-body axis wind components. Finally, the rudder effect in the equation describing the 𝛽 dynamics has been neglected. It is worth noting that Section 5 has shown that the designs and the simulations of the NLGA residual generators are robust with respect to the last approximation. In fact, the model of the 𝛽 dynamics will never be used.

3. PM Residual Generators

Let us consider the input-output representation of a continuous-time, time-invariant linear dynamic system affected by faults and disturbances in the form 𝑃(𝑠)𝑦(𝑡)=𝑄𝑐(𝑠)𝑐(𝑡)+𝑄𝑑(𝑠)𝑑(𝑡)+𝑄𝑓(𝑠)𝑓(𝑡),(4) where 𝑦(𝑡)𝑚 is the output vector, 𝑐(𝑡)𝑙𝑐 is the input vector, 𝑑(𝑡)𝑙𝑑 is the disturbance vector, and 𝑓(𝑡)𝑙𝑓 is the fault vector; 𝑃(𝑠), 𝑄𝑐(𝑠), 𝑄𝑑(𝑠), and 𝑄𝑓(𝑠) are known polynomial matrices of proper dimensions.

Models of type (4) can be frequently found in practice by applying well-known physical laws to describe the input-output dynamical links of various systems. Algorithms to transform multivariable state-space models to equivalent multiple-input-multiple-output (MIMO) polynomial representations and vice versa are available [15]. Suitable software routines for multivariable system transformations have been implemented by the authors in the Matlab environment. In fact, the Matlab software for state-space and transfer function conversions is not able to manage directly MIMO models, since they are considered as concatenations of single-input-single-output (SISO) systems.

An important aspect of the residual generator design concerns the decoupling properties of the disturbance 𝑑(𝑡). The decoupling can be obtained premultiplying all the terms of (4) by the matrix 𝐿(𝑠)𝒩𝑙(𝑄𝑑(𝑠)), that is, the left null-space of the matrix 𝑄𝑑(𝑠): 𝐿(𝑠)𝑃(𝑠)𝑦(𝑡)𝐿(𝑠)𝑄𝑐(𝑠)𝑐(𝑡)=𝐿(𝑠)𝑄𝑓(𝑠)𝑓(𝑡).(5) Hence, the residual generator for the system of (4) is represented by 𝑅(𝑠)𝑟(𝑡)=𝐿(𝑠)𝑃(𝑠)𝑦(𝑡)𝐿(𝑠)𝑄𝑐(𝑠)𝑐(𝑡)=𝐿(𝑠)𝑄𝑓(𝑠)𝑓(𝑡),(6) where it is assumed that 𝑟(𝑡) and 𝐿(𝑠) is a polynomial row vector. The polynomial 𝑅(𝑠) can be arbitrarily selected among the polynomials with degree greater than or equal to 𝑛𝑟, where 𝑛𝑟 is the maximum row-degree of the pair {𝐿(𝑠)𝑃(𝑠),𝐿(𝑠)𝑄𝑐(𝑠)}. Moreover, if all the roots of 𝑅(𝑠) lie in the open left-half 𝑠-plane, it assures the stability of the filter of (6). Without loss of generality, it is assumed that 𝑅(0)=1. Remark 1. If the matrix 𝑄𝑑(𝑠) is of full-column rank (i.e., rank𝑄𝑑(𝑠)=𝑙𝑑), 𝒩𝑙(𝑄𝑑(𝑠)) has dimension 𝑚𝑙𝑑. Therefore, a polynomial matrix 𝐵(𝑠), whose rows represents a minimal polynomial basis of 𝒩𝑙(𝑄𝑑(𝑠)), has 𝑚𝑙𝑑 rows and 𝑚 columns.

This work is focused on the problem of detecting and isolating additive faults acting on the input and output sensors of the monitored system. If the input-output measurements are modelled by the relations of (7): 𝑐(𝑡)=𝑐(𝑡)+𝑓𝑐𝑦(𝑡),(𝑡)=𝑦(𝑡)+𝑓𝑜(𝑡),(7) the system of (4) becomes 𝑦𝑃(𝑠)(𝑡)𝑓𝑜(𝑡)=𝑄𝑐𝑐(𝑠)(𝑡)𝑓𝑐(𝑡)+𝑄𝑑(𝑠)𝑑(𝑡),(8) under the assumptions that 𝑄𝑓(𝑠)𝑓(𝑡)=[𝑄𝑐(𝑠),𝑃(𝑠)][𝑓𝑇𝑐(𝑡),𝑓𝑇𝑜(𝑡)]𝑇. Thus, the residual generator of (6) is written as 𝑅(𝑠)𝑟(𝑡)=𝐿(𝑠)𝑃(𝑠)𝑦(𝑡)𝐿(𝑠)𝑄𝑐(𝑠)𝑐(𝑡)=𝐿(𝑠)𝑃(𝑠)𝑓𝑜(𝑡)𝐿(𝑠)𝑄𝑐(𝑠)𝑓𝑐(𝑡).(9)Remark 2. The residual generator described by (7) and (9) can be seen as an errors-in-variables (EIV) model [16] with respect the input and output variables, as the measurements that feed the residual function are affected by additive faults. This description highlights the importance of the residual generator in the form of (9).Remark 3. The diagnostic capabilities of the residual generator of (6) strongly depend on the choice of the terms 𝐿(𝑠) and 𝑅(𝑠). This paper proposes a method for the design of these polynomials, under the assumption that 𝑓(𝑡) is a scalar and, consequently, 𝑄𝑓(𝑠) is a vector. The rationale of this assumption is commented in Section 3.2 where the fault isolation method is proposed.

In the following, the freedom design in the selection of the rows of the polynomial matrix 𝐿(𝑠) is investigated when 𝑞=𝑚𝑙𝑑2. These degrees of freedom are used to optimise the sensitivity properties of 𝑟(𝑡) with respect to the fault 𝑓(𝑡), for example, by maximising the steady-state gain of the transfer function 𝐺𝑓(𝑠)=𝐿(𝑠)𝑄𝑓(𝑠)/𝑅(𝑠).

If 𝑏𝑖(𝑠) (𝑖=1,,𝑞) are the row vectors of the basis 𝐵(𝑠), 𝐿(𝑠) can be expressed as linear combination of these vectors: 𝐿(𝑠)=𝑞𝑖=1𝑘𝑖𝑏𝑖(𝑠),(10) where 𝑘𝑖 are real constants maximising: lim𝑠01[𝑅(𝑠)𝑞𝑖=1𝑘𝑖𝑏𝑖(𝑠)]𝑄𝑓(𝑠)=[𝑞𝑖=1𝑘𝑖𝑏𝑖(0)]𝑄𝑓(0)(11) with the constraint 𝑞𝑖=1𝑘2𝑖=1.(12) Under these assumptions, when the fault 𝑓(𝑡) is a step-function of magnitude 𝐹, the steady-state residual value is lim𝑡𝑟(𝑡)=lim𝑠0𝑠𝐿(𝑠)𝑄𝑓(𝑠)𝐹𝑅(𝑠)𝑠=[𝑞𝑖=1𝑘𝑖𝑏𝑖(0)]𝑄𝑓(0)𝐹.(13) If the following real vectors are defined as 𝑘𝑘=1𝑘2𝑘𝑞,𝑎=𝐵(0)𝑄𝑓𝑎(0)=1𝑎2𝑎𝑞,(14) the problem of the maximisation of the residual fault sensitivity can be recast as follows.

Proposition 1. Given the vector 𝑎, the vector 𝑘 that maximises the steady-state fault sensitivity, that is, the function 𝑊 given by the expression 𝑊=𝑎𝑇𝑘=𝑞𝑖=1𝑎𝑖𝑘𝑖,(15) under the constraint of (12), can be found by solving ̃𝑘=argmax𝑊(𝑘).(16)

The solution to the problem described by Proposition 1 can be derived as follows. The constraint of (12) describes a hypersphere, whilst the expression of the function of (15) is a hyperplane. The unknown coefficients 𝑘 must belong to both the hyperplane and the hypersphere. Therefore, the points of tangency between the hypersphere and the hyperplane represents the solutions that maximise or minimise 𝑊. As shown below, the solution of the problem described by Proposition 1 exists and is unique.

Proof. From (12), 𝑘1 is expressed as a function of 𝑘2, 𝑘3,,𝑘𝑞, and it is substituted into (15): 𝑊=𝑎11𝑘22𝑘23𝑘2𝑞+𝑎2𝑘2++𝑎𝑞𝑘𝑞.(17) By computing 𝑊=0, that is, 𝜕𝑊𝜕𝑘2=12𝑎12𝑘21𝑘22𝑘23𝑘2𝑞+𝑎2=0,𝜕𝑊𝜕𝑘3=12𝑎12𝑘31𝑘22𝑘23𝑘2𝑞+𝑎3=0,𝜕𝑊𝜕𝑘𝑞=12𝑎12𝑘𝑞1𝑘22𝑘23𝑘2𝑞+𝑎𝑞=0,(18) and squaring the expression, after algebraic manipulation: 𝑎22=𝑎22+𝑎21𝑘22+𝑎22𝑘23++𝑎22𝑘2𝑞,𝑎23=𝑎23𝑘22+𝑎23+𝑎21𝑘23++𝑎23𝑘2𝑞,𝑎2𝑞=𝑎2𝑞𝑘22+𝑎2𝑞𝑘23𝑎++2𝑞+𝑎21𝑘2𝑞,(19) an expression in the form of 𝐴𝑥=𝑏 is obtained, where 𝑎𝐴=22+𝑎21𝑎22𝑎22𝑎23𝑎23+𝑎21𝑎23𝑎2𝑞𝑎2𝑞𝑎2𝑞+𝑎21,𝑘𝑥=22𝑘23𝑘2𝑞𝑎,𝑏=22𝑎23𝑎2𝑞.(20) The unknown vector ̃𝑥, under the constraint of (12), can be expressed as follows: ̃𝑥=1𝑞1𝑖=1𝐴1𝑏𝑖𝐴1𝑏,(21) where (𝐴1𝑏)𝑖 is the 𝑖th element of the vector 𝐴1𝑏. The vector ̃𝑥 represents the squares of the solution of the problem of Proposition 1.
Let us indicate Ω the set of the vectors 𝑘 whose elements are the square roots of the elements of ̃𝑥. As every element can be taken both with signs “+” and “”, such vectors are 2𝑞. Therefore, the solution ̃𝑘 of Proposition 1 can be reformulated as ̃𝑘=argmax𝑘Ω𝑊(𝑘).(22)

Remark 4. The matrix 𝐴 can be expressed as 𝐴=𝐸+𝑎21𝐼𝑞1, where 𝐸 is a matrix with equal columns. If 𝑎10, this assumption guarantees the existence of 𝐴1, and consequently the existence and the uniqueness of the solution 𝐴1𝑏. Obviously, if 𝑎1=0 and 𝑎𝑗0, it is sufficient to express 𝑘𝑗 as function of the remaining variables and to reapply the same procedure.

Remark 5. The same solution can be found by maximising the function |𝑊|. Due to the symmetry properties, the maximisation of |𝑊| admits two solutions corresponding to the maximum and the minimum of the function 𝑊. Moreover, the choice of the quadratic constraint of (12) guarantees the unicity of the solution to the problem of Proposition 1.

Remark 6. The problem described by Proposition 1 could have been solved also in a numerical way, that is, by searching 𝑘 that maximises 𝑊 on the surface of the 𝑞-dimensional hypersphere. However, the computational cost of this numerical solution can be a drawback when 𝑞 is big.

3.1. PM Residual Design

Section 3 has shown how to maximise the steady-state gain of the continuous-time transfer function 𝐺𝑓(𝑠)=𝐿(𝑠)𝑄𝑓(𝑠)/𝑅(𝑠) trough a suitable choice of the real vector 𝑘 (i.e., ̃𝑘=𝑘). The design of the filter of (6) has been completed here by introducing a method for assigning both the zeros and the poles of the continuous-time transfer function 𝐺𝑓(𝑠). The zeros and poles location influences the transient characteristics (maximum overshoot, delay time, rise time, settling time, etc.) of the filter of (6). In many applications, these characteristics must be kept within tolerable or prescribed limits in order to guarantee good performances of the filter in terms, for example, of fault detection times and false-alarm probabilities.

Remark 7. When ̃𝑘𝑘=, the polynomial 𝐿(𝑠)𝑄𝑓(𝑠)=𝑘𝑇𝐵(𝑠)𝑄𝑓(𝑠) is fixed and no freedom degree is left to arbitrarily assign the zeros. In order to solve this problem, a polynomial vector 𝑘(𝑠) can be considered. Under this assumption, 𝐿(𝑠) still belongs to the subspace 𝒩𝑙(𝑄𝑑(𝑠)), where the terms 𝑘𝑖 are polynomial coefficients.

The previous consideration leads to introduce the polynomial 𝐸(𝑠)=𝑘𝑇(𝑠)𝐵(𝑠)𝑄𝑓(𝑠), where 𝑘(𝑠) is a 𝑞-dimensional polynomial vector whose 𝑖th element has the form 𝑘𝑖(𝑠)=𝑛𝑘𝑗=0𝑘𝑗𝑖𝑠𝑗.(23) The degree 𝑛𝑘 and the 𝑞×𝑛𝑘 coefficients 𝑘𝑗𝑖 are freedom design (𝑗0) that are exploited for obtaining the desired roots of the polynomial 𝐸(𝑠). However, in order to maximise the steady-state gain, as shown in Section 3, the following condition must hold: ̃𝑘𝑘(0)=𝑘=1𝑘2𝑘𝑞𝑘0𝑖=𝑘𝑖,𝑖=1,,𝑞.(24)

Definition 1. 𝐻(𝑠) is the reference polynomial whose roots are the zeros to be assigned: 𝐻(𝑠)=𝑛𝑗=0𝑗𝑠𝑗.(25) Since the constraint of (24) must hold, 𝑘𝐻(0)=𝑇𝐵(0)𝑄𝑓(0). Obviously, this assumption does not provide any restriction on the roots assignable. Under the previous considerations, the zero assignment and pole placement problem is formulated as follows.

Proposition 2. The degree 𝑛𝑘 and the coefficients 𝑘𝑗𝑖 have to be determined under the constraint of (24) in order to obtain 𝐸(𝑠)=𝐻(𝑠).

Proof. In Section 3, the polynomial vector 𝑎(𝑠)=𝐵(𝑠)𝑄𝑓(𝑠) was defined. Its 𝑖th element is a known polynomial of a certain degree, 𝑛𝑎𝑖. If 𝑛𝑎 is defined as follows: 𝑛𝑎=max𝑖=1,,𝑞𝑛𝑎𝑖,(26) the 𝑖th element of 𝑎(𝑠) can be always written as a polynomial of degree 𝑛𝑎: 𝑎𝑖(𝑠)=𝑛𝑎𝑗=0𝑎𝑗𝑖𝑠𝑗(27) by imposing that 𝑎𝑗𝑖=0 when 𝑗>𝑛𝑎𝑖.
As 𝐸(𝑠)=𝑘𝑇(𝑠)𝑎(𝑠), by multiplying (23) and (27), it results 𝐸(𝑠)=𝑞𝑛𝑖=1𝑘+𝑛𝑎𝑗=0(𝛼+𝛽=𝑗𝑘𝛼𝑖𝑎𝛽𝑖)𝑠𝑗=𝑛𝑘+𝑛𝑎𝑗=0𝑒𝑗𝑠𝑗,(28) where 𝑒𝑗=𝑞𝑖=1𝛼+𝛽=𝑗𝑘𝛼𝑖𝑎𝛽𝑖.(29) Equations (28) and (29) assume that 𝑘𝛼𝑖=0 when 𝛼>𝑛𝑘 and 𝑎𝛽𝑖=0 when 𝛽>𝑛𝑎. Note that the coefficients 𝑒1,,𝑒𝑛𝑘+𝑛𝑎 depend on the freedom design 𝑘1𝑖,,𝑘𝑛𝑘𝑖. On the other hand, 𝑒0 is fixed as the coefficients 𝑘0𝑖 are assigned by (24).
Let us suppose that 𝑛𝑛𝑘+𝑛𝑎. By imposing 𝐸(𝑠)=𝐻(𝑠), from (29) and (25), the following expressions are computed: 𝑞𝑖=1𝛼+𝛽=𝑗(𝛼0)𝑘𝛼𝑖𝑎𝛽𝑖=𝑗𝑞𝑖=1𝑘0𝑖𝑎𝑗𝑖,𝑗=1,,𝑛𝑘+𝑛𝑎.(30) Equations (24) and (30) represent a linear system with 𝑛𝑘+𝑛𝑎 equations and 𝑞×𝑛𝑘 unknowns, which can be expressed in the classical form 𝐴𝑥=𝑏, where 𝑎𝐴=01𝑎0𝑞0000𝑎01𝑎0𝑞𝑎𝑛𝑎1𝑎𝑛𝑎𝑞00𝑎𝑛𝑎1𝑎𝑛𝑎𝑞𝑎00000001𝑎0𝑞0000𝑎𝑛𝑎1𝑎𝑛𝑎𝑞,𝑘𝑥=11𝑘1𝑞𝑘21𝑘2𝑞𝑘𝑛𝑘1𝑘𝑛𝑘𝑞,𝑏=1𝑞𝑖=1𝑘0𝑖𝑎1𝑖𝑛𝑎𝑞𝑖=1𝑘0𝑖𝑎𝑛𝑎𝑖𝑛𝑎+1𝑛𝑎+𝑛𝑘.(31) The degree 𝑛𝑘 of the polynomials 𝑘𝑖(𝑠) has to be chosen in order to obtain a solvable system (i.e., rank𝐴=rank[𝐴𝑏]).

In order to understand the proposed solution, the following points should be considered.

(i)The choice of 𝑛𝑘 must guarantee that the relations 𝑛𝑛𝑘+𝑛𝑎 are satisfied.(ii) When 𝑞2, the difference between the number of unknown terms and the number of equations, that is, (𝑞1)×𝑛𝑘𝑛𝑎, is greater than zero if 𝑛𝑘 is selected sufficiently high.(iii)Even if the system admits solutions, the inverse of the matrix 𝐴 may not exist; in such case there are infinite solutions and the one associated to the pseudoinverse of 𝐴, that is, 𝐴+𝑏 can be computed.

Remark 8. The use of a polynomial vector 𝑘(𝑠) instead of a real vector 𝑘 has the drawback of increasing the complexity of the residual generator. Many FDI applications require that 𝐺𝑓(𝑠)𝐺𝑓=(0)𝐹(0)𝐹(𝑠),(32) where 𝐹(𝑠) is an arbitrary polynomial. These cases do not require a 𝑘(𝑠) such that 𝐸(𝑠)=𝐺𝑓(0), but it is enough considering ̃𝑘𝑘= and imposing 𝐸𝑅(𝑠)=0(𝑠)𝐹(𝑠)𝐺𝑓(0)𝐹(0),(33) where 𝐸0̃(𝑠)=𝑘𝐵(𝑠)𝑄𝑓(𝑠). However, there is a restriction on the choice of 𝐹(𝑠). In fact, due to the realisability condition, deg{𝐹(𝑠)}>𝑛𝑟deg{𝐸0(𝑠)}. Moreover, the method cannot be applied if 𝐸0(𝑠) admits one or more roots with positive real part, as the residual generator would become unstable. These cases require an approximate solution.

Remark 9. This section is focused on the design of residual generators on the basis of a given reference function with disturbance-decoupling and fault sensitivity maximisation properties. The pole location influences the transient dynamics of the designed residual filters, while the steady-state properties depend on the PM residual design, as it maximises the residual steady-state values with respect to step faults affecting input and output sensors. The poles of the residual functions could be optimised with respect to both fault and disturbance terms, as shown, for example, in a work by the same authors [17].

3.2. PM Fault Isolation

This section addresses the design problem of residual generator banks for the isolation of faults affecting the input and the output sensors. This design is performed by using the disturbance-decoupling method suggested in Section 3.

To univocally isolate a fault concerning one of the output sensors, under the hypotheses that the input sensors and the remaining output sensors are fault-free, a bank of residual generator filters is used. The number of these generators is equal to the number 𝑚 of the system outputs, and the 𝑖th device (𝑖=1,,𝑚) is driven by all but the 𝑖th output and all the inputs of the system. In this case, a fault on the 𝑖th output sensor affects all but the 𝑖th residual generator.

In presence of a fault on the 𝑗th output sensor, the measured output 𝑦(𝑡) can be expressed as follows: 𝑦(𝑡)=𝑦(𝑡)+𝑓𝑜𝑗(𝑡),(34) with 𝑓𝑜𝑗(𝑡)=00𝑜𝑗(𝑡)00𝑇,(35) where 𝑜𝑗(𝑡) represents the 𝑗th output fault function.

In these conditions, the system of (4) becomes 𝑃(𝑠)𝑦(𝑡)𝑝𝑗(𝑠)𝑜𝑗(𝑡)=𝑄𝑐(𝑠)𝑐(𝑡)+𝑄𝑑(𝑠)𝑑(𝑡),(36) where 𝑝𝑗(𝑠) is the 𝑗th column of the matrix 𝑃(𝑠).

Let us indicate 𝐿𝑜𝑖(𝑠) a polynomial row vector belonging to the basis of the left null space of the matrix [𝑄𝑑(𝑠)𝑝𝑖(𝑠)]. The expression of the 𝑖th filter when a fault is acting on the 𝑗th output sensor is obtained by multiplying (36) by 𝐿𝑜𝑖(𝑠): 𝑅𝑜𝑖(𝑠)𝑟𝑜𝑖(𝑡)=𝐿𝑜𝑖(𝑠)𝑃𝑖(𝑠)𝑦𝑖(𝑡)𝐿𝑜𝑖(𝑠)𝑄𝑐=𝐿(𝑠)𝑐(𝑡)𝑜𝑖(𝑠)𝑝𝑗(𝑠)𝑜𝑗(𝑡)for𝑗𝑖,0for𝑗=𝑖,(37) where 𝑃𝑖(𝑠) is the matrix obtained by deleting from 𝑃(𝑠) the 𝑖th column, and 𝑦𝑖(𝑡) represents the (𝑚1)-dimensional vector obtained by deleting from 𝑦(𝑡) its 𝑖th component.

From the comparison between (37) and (6) with 𝑓(𝑡) if 𝑞=𝑚𝑙𝑑12, the methods shown in Sections 3 and 3.1 can be exploited for the design of the 𝑖th filter. In particular, the parameters of this filter can be properly chosen in order to optimise its performances when a fault is acting on the 𝑗th output sensor.

In more detail, as shown in Section 3, 𝐿𝑜𝑖(𝑠) is chosen to maximise the steady-state gain in the presence of the fault 𝑓𝑜𝑗(𝑡). Moreover, as shown in Section 3.1, 𝑅𝑜𝑖(𝑠) is chosen to obtain a fixed behaviour of the transfer function due to the fault 𝑓𝑜𝑗(𝑡).

It is worth noting that the similar design technique can be used for input sensor fault isolation.

The problem requirements determine the selection of the specific fault with respect to which the design depends. Most often in practice, it is important to obtain good performance with respect to all possible faults rather than optimal behaviour with respect to one specific fault. In this situation, a different design of the filter behaviour for each fault situation is needed.

4. NLGA Residual Generators

The considered NLGA to the FDI problem is suggested in [18] and formally developed in [10]. It consists in finding, by means of a coordinate change in the state space and in the output space, an observable subsystem which, if possible, is affected by the fault and not affected by disturbance. In this way, necessary and sufficient conditions for the FDI problem to be solvable are given. Finally, a residual generator can be designed on the basis of the model of the observable subsystem.

More precisely, the approach considers a nonlinear system model in the form ̇𝑥=𝑛(𝑥)+𝑔(𝑥)𝑐+(𝑥)𝑓+𝑝(𝑥)𝑑,𝑦=(𝑥)(38) in which the state vector 𝑥𝒳 (an open subset of 𝑛), 𝑐(𝑡)𝑐 is the control input vector, 𝑓(𝑡) is the fault, 𝑑(𝑡)𝑑 the disturbance vector (embedding also the faults which have to be decoupled), and 𝑦𝑚 the output vector; whilst 𝑛(𝑥), (𝑥), the columns of 𝑔(𝑥) and 𝑝(𝑥) are smooth vector fields; and (𝑥) is a smooth map.

Therefore, if 𝑃 represents the distribution spanned by the column of 𝑝(𝑥), the NLGA method can be devised as it follows: first, determine the largest observability codistribution contained in 𝑃, denoted with Ω [10].

If (𝑥)Ω, the design procedure can continue, otherwise, the fault is not detectable; whenever the previous condition is satisfied, it can be found a surjection Ψ1 and a function Φ1 fulfilling Ωspan{𝑑}=span{𝑑(Ψ1)} and Ω=span{𝑑(Φ1)}, respectively. The functions Ψ(𝑦) and Φ(𝑥) defined as Ψ(𝑦)=𝑦1𝑦2=Ψ1𝐻(𝑦)2𝑦,Φ(𝑥)=𝑥1𝑥2𝑥3=Φ1𝐻(𝑥)2Φ(𝑥)3(𝑥)(39) are (local) diffeomorphisms, where 𝐻2 is a selection matrix (i.e., a matrix in which any row has all 0 entries but one, which is equal to 1), Φ1(𝑥) represents the measured part of the state which is affected by 𝑓 and not affected by 𝑑, and Φ3(𝑥) represents the not measured part of the state which is affected by 𝑓 and by 𝑑.

In the new (local) coordinate defined previously, the system of (38) is described by the relations in the form ̇𝑥1=𝑛1𝑥1,𝑥2+𝑔1𝑥1,𝑥2)𝑐+1𝑥1,𝑥2,𝑥3̇)𝑓,𝑥2=𝑛2𝑥1,𝑥2,𝑥3+𝑔2𝑥1,𝑥2,𝑥3𝑐+2𝑥1,𝑥2,𝑥3𝑓+𝑝2𝑥1,𝑥2,𝑥3̇𝑑,𝑥3=𝑛3𝑥1,𝑥2,𝑥3+𝑔3𝑥1,𝑥2,𝑥3𝑐+3𝑥1,𝑥2,𝑥3𝑓+𝑝3𝑥1,𝑥2,𝑥3𝑑,𝑦1=𝑥1,𝑦2=𝑥2(40) with 1(𝑥1,𝑥2,𝑥3) not identically zero. Denoting 𝑥2 with 𝑦2 and considering it as an independent input, it can be singled out the 𝑥1-subsystem: ̇𝑥1=𝑛1(𝑥1,𝑦2)+𝑔1(𝑥1,𝑦2)𝑐+1(𝑥1,𝑦2,𝑥3)𝑓,𝑦1=(𝑥1),(41) which is affected by the single fault 𝑓 and decoupled from the disturbance vector. This subsystem has been exploited for the design of the residual generator for the FDI of the fault 𝑓.

As already described in Section 2, the proposed NLGA FDI scheme is designed on the basis of a model structure of the input affine type as described in [10]. For this reason, the aircraft simulation model has been simplified and the nonlinear model of (3) has been considered for the NLGA design.

Under these assumptions, by means of computations detailed in [19], the residual generators for detecting the faults affecting the aircraft input sensors are obtained. In particular, the residual generator for the elevator 𝑟𝛿𝑒(𝑡), with 𝑘𝛿𝑒>0, is described by the relation ̇𝜉1=𝑉2𝑚𝐶𝐷0+𝐶𝐷𝛼𝛼+𝐶𝐷𝛼2𝛼2+𝑉cos𝛼2𝑚𝐶𝐿0+𝐶𝐿𝛼𝛼sin𝛼𝑔sin𝜃+𝑉𝑞𝜔𝐶sin𝛼𝑚0+𝐶𝑚𝛼𝛼+𝐶𝑚𝑞𝑞𝜔𝑚𝑡𝑑𝑉2𝐼𝑧𝐼𝑥𝑚𝑡𝑑𝑝𝜔𝑟𝜔𝐶𝛿𝑒𝑚𝑡𝑑𝑉2𝛿𝑒+𝑘𝛿𝑒𝐼𝑉cos𝛼𝑦𝑚𝑡𝑑𝑞𝜔𝜉1,𝑟𝛿𝑒=𝐼𝑉cos𝛼𝑦𝑚𝑡𝑑𝑞𝜔𝜉1.(42) The aileron residual generator 𝑟𝛿𝑎(𝑡), with 𝑘𝛿𝑎>0, has the form ̇𝜉2=𝐶𝑙𝛽𝛽+𝐶𝑙𝑝𝑝𝜔𝐼𝑥𝑉2+𝐼𝑦𝐼𝑧𝐼𝑥𝑞𝜔𝑟𝜔+𝐶𝛿𝑎𝐼𝑥𝑉2𝛿𝑎+𝑘𝛿𝑎𝑝𝜔𝜉2,𝑟𝛿𝑎=𝑝𝜔𝜉2.(43) The rudder residual generator 𝑟𝛿𝑟(𝑡), with 𝑘𝛿𝑟>0, is written in the form ̇𝜉3=𝐶𝑛𝛽𝛽+𝐶𝑛𝑟𝑟𝜔𝐼𝑧𝑉2+𝐼𝑥𝐼𝑦𝐼𝑧𝑝𝜔𝑞𝜔+𝐶𝛿𝑟𝐼𝑧𝑉2𝛿𝑟+𝑘𝛿𝑟𝑟𝜔𝜉3,𝑟𝛿𝑟=𝑟𝜔𝜉3.(44) The throttle residual generator 𝑟𝛿th(𝑡), with 𝑘𝛿th>0, has the form ̇𝜉4=𝑡𝑛𝑛3𝑒+𝑡𝑓𝑛𝑒𝑡0+𝑡1𝑛𝑒𝛿th+𝑘𝛿th𝑛𝑒𝜉4,𝑟𝛿th=𝑛𝑒𝜉4.(45)Remark 10. It is worth observing that each residual generator is affected by a single input sensor fault and is decoupled from the wind components and the faults affecting the remaining input sensors. This feature can be obtained by defining a different 𝑝(𝑥) for each residual generator design [19]. In this way, the tuning of the residual generator gains 𝑘𝛿𝑒, 𝑘𝛿𝑎, 𝑘𝛿𝑟, and 𝑘𝛿𝑡 can be carried out independently. Finally, by a straightforward analysis, the positive sign of each gain is a necessary and sufficient condition for the asymptotic stability of the dynamics (42)–(45). A procedure optimising the tradeoff between the fault sensitivity and the robustness to the modelling errors and disturbances of the generic residual generator is proposed in the next section.

4.1. NLGA Robustness Improvement

The proposed NLGA-based scheme consists of two design steps: (1)the structural decoupling of critical disturbances (wind gust and turbulence) and critical modelling errors can be obtained as described in Section 4;(2)the nonlinear residual generators robustness is improved by minimising the effects of both noncritical disturbances and modelling errors, whilst maximising the fault effects on the residual signals. In order to apply the robustness improvement procedure presented in this section, the considered framework is restricted to suitable scalar components of the 𝑥1-subsystem (41). In particular, the vectors 𝑥1 and 𝑦1 are decomposed as follows: 𝑥1=𝑥11𝑥1c,𝑦1=𝑦11𝑦1c,(46) where 𝑥11, 𝑦11, and, correspondingly, it follows that 𝑛1𝑛()=11𝑛()1c(),𝑔1𝑔()=11𝑔()1c(),1()=11()1c()(47) Let us consider the following conditions: 𝑦11=11(𝑥11)𝑦1c=1c(𝑥1c)11()0,(48) where 11() is a smooth map and 1c() is an invertible smooth map. It is important to highlight that if the constraints (48) are satisfied, the decomposition (46)-(47) can always be applied to obtain the following 𝑥11-subsystem: ̇𝑥11=𝑛11𝑥11,𝑦1c,𝑦2+𝑔11𝑥11,𝑦1c,𝑦2𝑐+11𝑥11,𝑦1c,𝑦2,𝑥3𝑓,𝑦11=11𝑥11.(49) As can be seen in [19], the conditions (48) are satisfied for the considered aircraft application, hence, from now on, the scalar 𝑥11-subsystem (49) is referred to in place of the 𝑥1-subsystem (41).

It can be noted that the tuning of the residual generator gains, in the framework of the 𝑥11-subsystem (49), cannot be properly carried out. In fact, the critical disturbances are structurally decoupled but the noncritical ones are not considered. For this reason, to achieve robustness of the residual generators, the tuning of the gains is performed by embedding the description of the noncritical disturbances in the 𝑥11-subsystem as follows: ̇𝑥11=𝑛11𝑥11,𝑦1c,𝑦2+𝑔11𝑥11,𝑦1c,𝑦2𝑐+11𝑥11,𝑦1c,𝑦2,𝑥3𝑓+𝑒𝑥11,𝑦1c,𝑦2,𝑥3𝜁,𝑦11=𝑥11+𝜈,(50) where, to simplify the treatment without loss of generality (accordingly to the considered aircraft application), the state variable 𝑥11 is supposed to be directly measured. Moreover, the variable 𝜈 is the measurement noise on 𝑥11. Finally, the variable 𝜁 and the related scalar field 𝑒() represent the noncritical effects which have not been considered in the simplified aircraft model (3) used for the NLGA scheme.

The following system, which is referred to as filter form, represents a generic scalar residual generator (based on the subsystem (50)) to which (42)–(45) belong as a particular case ̇𝜉𝑓=𝑛11𝑦11,𝑦1c,𝑦2+𝑔11𝑦11,𝑦1c,𝑦2𝑐+𝑘𝑓𝑦11𝜉𝑓,𝑟𝑓=𝑦11𝜉𝑓,(51) where the gain 𝑘𝑓 has to be tuned in order to minimise the effects of the disturbances 𝜁 and 𝜈 whilst maximise the effects of the fault 𝑓 on the residual 𝑟𝑓. The quantification both of the disturbances and of the fault effects on the residual can be obtained by defining the estimation error 𝑥𝑓=𝑥11𝜉𝑓,(52) which allows to write the following equivalent residual model: ̇̃𝑥𝑓=𝑛11𝑥11,𝑦1c,𝑦2𝑛11𝑦11,𝑦1c,𝑦2+𝑔11𝑥11,𝑦1c,𝑦2𝑐𝑔11𝑦11,𝑦1c,𝑦2𝑐+11𝑥11,𝑦1c,𝑦2,𝑥3𝑓+𝑒𝑥11,𝑦1c,𝑦2,𝑥3𝜁𝑘𝑓𝑥𝑓𝑘𝑓𝑟𝜈,𝑓=𝑥𝑓+𝜈.(53) In order to apply the effective mixed / approach [3, 20] to tune 𝑘𝑓, the system (53) has to be linearised in the neighbourhood of a stationary flight condition, as suggested in [2] with reference to the optimisation of nonlinear unknown input observers. It is worth observing that the considered aircraft application is characterised by small excursions of the state, input, and output variables with respect to their trim values 𝑥10, 𝑥30, 𝑐0, 𝑦10, and 𝑦20, hence the robustness of the nonlinear residual generator is achieved. The linearisation of (53) is the following: ̇̃𝑥𝑓=𝑘𝑓𝑥𝑓𝑘𝑓𝜈+mf+𝑞𝑟𝜁,𝑓=𝑥𝑓+𝜈,(54) where 𝑎=𝜕𝑛11()𝜕𝑥11|||(𝑥10,𝑦20),𝑏=𝑔11||()(𝑥10,𝑦20),𝑚=11||()(𝑥10,𝑦20,𝑥30)||,𝑞=𝑒()(𝑥10,𝑦20,𝑥30),𝑞𝜁=𝑞𝜁𝑎𝜈.(55) Now, it is important to note that in place of the residual generators in the filter form (51), the following observer form of the residual generators can be used: ̇𝜉𝑜=𝑛11𝜉𝑜,𝑦1c,𝑦2+𝑔11𝜉𝑜,𝑦1c,𝑦2𝑐+𝑘𝑜𝑦11𝜉𝑜,𝑟𝑜=𝑦11𝜉𝑜.(56) For the same reasons previously described, the estimation error 𝑥𝑜 is introduced: 𝑥𝑜=𝑥11𝜉𝑜,(57) hence ̇̃𝑥𝑜=𝑛11𝑥11,𝑦1c,𝑦2𝑛11(𝜉𝑜,𝑦1c,𝑦2+𝑔11(𝑥11,𝑦1c,𝑦2𝑐𝑔11(𝜉𝑜,𝑦1c,𝑦2𝑐+11𝑥11,𝑦1c,𝑦2,𝑥3𝑓+𝑒𝑥11,𝑦1c,𝑦2,𝑥3𝜁𝑘𝑜𝑥𝑜𝑘𝑜𝑟𝜈,𝑜=𝑥𝑜+𝜈.(58) The linearisation of (58) is ̇̃𝑥𝑜=𝑎𝑘𝑜𝑥𝑜𝑘𝑜𝑟𝜈+mf+𝑞𝜁,𝑜=𝑥𝑜+𝜈.(59) Both the linearised models (54) and (59) of the residual generators in the filter form and observer form, respectively, can be represented by the following general form: ̇𝐸̃𝑥=𝑎𝑘̃𝑥+1𝑘𝐸2𝜀+mf,𝑟=̃𝑥+𝐸2𝜀(60) with 𝐸1𝑒=[110] as well as the following positions: generalform̃𝑥𝜀𝑟𝑎𝑘𝑒11𝐸2lterform𝑥𝑓𝜁𝜈𝑟𝑓0𝑘𝑓𝑞01observerform𝑥𝑜𝜁𝜈𝑟𝑜𝑎𝑘𝑜𝑞.01(61) On the basis of (60) and (61), the mixed / [3, 20] procedure is developed for the robustness improvement of the residual generators both in the filter and observer form. Since the considered NLGA residual generators are scalar, the / procedure leads to a new analytical solution.

The following definition will be used throughout the section.

Definition 2. The norms and of a stable transfer function 𝐺 are defined as 𝐺=sup𝜔0𝜎𝐺,𝐺𝑗𝜔=𝜎𝐺𝑗0,(62) where 𝜎 represents the maximum singular value, whilst 𝜎 the minimum singular value. The problem of the tradeoff between disturbances robustness and fault sensitivity is stated as follows.

Problem 1 (Mixed / Residual Robustness Improvement). Given two scalars 𝛽>0 and 𝛾>0, find the set 𝒦 defined as: 𝐺𝒦={𝑘𝑎𝑘<0,𝑟𝜀𝐺<𝛾,𝑟𝑓>𝛽},(63) where 𝐺𝑟𝜀𝑠=𝑠𝑎+𝑘1𝐸1𝑘𝐸2+𝐸2𝐺,(64)𝑟𝑓𝑠=𝑠𝑎+𝑘1𝑚.(65) In order to obtain the analytical solution of Problem 1, the following propositions are given.

Proposition 3. For all 𝑘,(𝑎𝑘)<0, then 𝐺𝑟𝜀2𝑒=max{1,211+𝑎2𝑘𝑎2},(66)sup𝑘(𝑎𝑘)<0𝐺𝑟𝜀=+.(67)𝐺𝑟𝜀𝑠=𝑒11𝑠𝑎+𝑘𝑠𝑎𝑠𝑎+𝑘,(68)

Proof. From the definition (64) 𝜎𝐺𝑟𝜀𝑗𝜔2=𝑒211𝑘𝑎2+𝜔2+𝑎2+𝜔2𝑘𝑎2+𝜔2=𝑒211+𝑎2+𝜔2𝑘𝑎2+𝜔2(69) hence it is possible to write 𝐺𝑟𝜀2=sup𝜉0𝑒211+𝑎2+𝜉𝑘𝑎2+𝜉.(70) so that it follows 𝒦𝛾=𝐺𝑘𝑎𝑘<0,𝑟𝜀<𝛾,𝛾>1(71) From the last expression, it is straightforward to obtain (66) and (67).

Proposition 4. The set 𝑘>𝑘with𝑘=𝑎+𝑒211+𝑎2𝛾.(72) is given by 𝑒211+𝑎2𝑘𝑎2<𝛾2,(73)

Proof. By means of Proposition 3, it is possible to write 𝑘>𝑎+𝑒211+𝑎2𝛾.(74) which holds for 𝛾>1

Proposition 5. If {𝐺𝑟𝑓𝐺𝑟𝜀<𝛾}, then 𝐺0<𝑟𝑓<𝛽max𝛾where𝛽max𝛾=𝑚𝛾𝑒211+𝑎2.(75) is given by 𝐺𝑟𝑓(𝑠)=𝑚/(𝑠𝑎+𝑘)

Proof. From the definition (65), it results 𝑚>0 and assuming, without loss of generality, that 𝐺𝑟𝑓=𝑚/(𝑘𝑎), it follows 𝐺𝑟𝑓>𝛽. By imposing 𝛽>0 with 𝑘<𝑎+(𝑚/𝛽), the constraint 𝛽 has to hold. Then, by recalling the result of Proposition 4, the maximum feasible value of 𝐺𝑟𝜀<𝛾 fulfilling the constraint 𝑘𝑚=𝑎+𝛽max𝛾,(76) is given by 𝛽max𝛾=𝑚𝑘=𝑎𝑚𝛾𝑒211+𝑎2.(77) hence 𝛾>1

Theorem 1. Given 𝛽]0,𝛽max(𝛾)[ and 𝒦, the set 𝑘𝒦=𝑘𝑘,𝑘,𝑘𝑚=𝑎+𝛽max𝛾,𝑚𝑘=𝑎+𝛽.(78) fulfilling the constraints of Problem 1 is given by 𝐺𝐽=𝑟𝑓𝐺𝑟𝜀.(79)

Proof. The proof of the theorem is not reported, as it is straightforward from Propositions 3, 4, and 5.

Remark 11. Let us consider the following performance index to maximise 𝐺𝑟𝜀=1,𝑘>𝑎+𝑒211+𝑎2𝑒211+𝑎2𝑘𝑎,𝑎<𝑘(𝑎+𝑒211+𝑎2)(80) From (66), it follows 𝑚𝐽=𝑘𝑎,𝑘>(𝑎+𝑒211+𝑎2𝑚),𝑒211+𝑎2,𝑎<𝑘(𝑎+𝑒211+𝑎2).(81) hence 𝑚𝐽=<𝑚𝑘𝑎𝑒211+𝑎2,𝑘>(𝑎+𝑒211+𝑎2).(82) From (81), it can be observed that 𝐽 In this way, the maximum value of the performance index 𝐽max=𝑚𝑒211+𝑎2𝑘𝒦𝒥={𝑘𝑎<𝑘(𝑎+𝑒211+𝑎2)}.(83) is 𝐽 The method proposed in this paper guarantees the maximum value of the performance index 𝐺𝑟𝜀<𝛾 as well as the constraints 𝐺𝑟𝑓>𝛽 and 𝛽𝑚/𝑒211+𝑎2 if 𝛽𝑚/𝑒211+𝑎2.
In fact, from 𝐺𝑟𝑓=𝑚𝑚𝑘𝑎>𝛽𝑒211+𝑎2,(84) it follows 𝑘<(𝑎+𝑒211+𝑎2) hence 𝛽.
Finally, from (75) it is always possible to find a 𝑚𝑒211+𝑎2𝛽𝛽max(𝛾)𝛾>1.(85) such that 𝑘 On the basis of Theorem 1, 𝛾>1 can be designed by means of the following procedure.

Procedure 1. (1) Choose 𝛽max(𝛾) to obtain a desired level of disturbance attenuation.
(2) Compute 𝛽]0,𝛽max(𝛾)[ and choose 𝑘]𝑘,𝑘[ to obtain a desired level of fault sensitivity.
(3) Choose 𝑘=𝑎+𝑚/𝛽max(𝛾), with 𝑘=𝑎+𝑚/𝛽 and 𝑘.
(4) Apply the chosen gain 𝑘𝑓 to the 𝑘𝑜 of (51) or to the 4 of (56) if the NLGA residual generator is in the filter form or in the observer form, respectively.

5. FDI Performance Estimation

To show the diagnostic characteristics brought by the application of the proposed FDI schemes to general aviation aircrafts, some numerical results obtained in the Matlab and Simulink environments are reported. The final performances that are achieved with the developed FDI schemes are finally reported. These performances are evaluated by means of extensive simulations applied to the aircraft simulation model. This section presents also some comparisons of the developed PM and NLGA FDI strategies with NN and UIKF FDI schemes.

The designed PM residual generator filters are fed by the 𝑐(𝑡) component input vector 9 and the 𝑦(𝑡) component output vector 4 acquired from the simulation aircraft model previously described. In particular, a bank of 4 residual generator filters has been used to detect input sensor faults regarding the 4 input control variables. Moreover, in order to obtain the fault isolation properties, each residual generator function of the considered bank is fed by all but one of the 9 control input signals and by the 3 output variables. Obviously, the residual generator bank has been designed to be decoupled from the 𝑑(𝑡)=[𝑤𝑢(𝑡),𝑤𝑣(𝑡),𝑤𝑤(𝑡)]𝑇 component wind disturbance vector 50. As to the NLGA residual generator filters, the aircraft synthesis model of (3), adopted for the design, is simplified with respect to the simulation model. Analogously to the PM, the approximations of the NLGA synthesis nonlinear model are related to a particular steady-state flight condition. For this reason, the switching for the NLGA FDI scheme is also required when a generic reference trajectory is followed. Hence, it is important to evaluate the robustness characteristic of a single design of NLGA residual generators when a large set of flight conditions is dealt with.

It is worth noting that the aircraft reference trajectories are typically made up of a sequence of steady-state flight conditions, each one described by the associated input state output set point and the linearised model of (1). As a consequence, all the FDI linear techniques are usually implemented by switching among the residual generators related to the different steady-state flight conditions. The target of this work is to reduce the switching by adopting robust PM residual generators. In particular, the robustness is achieved by using the same residual generators for a large set of flight conditions. The chosen single steady-state flight condition for designing both of the PM and of the NLGA residual generators is a coordinated turn characterised as follows: (i)the true airspeed is 1000 m/s;(ii)the curvature radius is 0 m;(iii)the flight-path angle is 330;(iv)the altitude is 0 m;(v)the flap deflection is 𝑅(𝑠). This represents one of most general flight condition due to the coupling of the longitudinal and lateral dynamics. Moreover, it is used in simulation to highlight the performances of the proposed methods in the nominal flight condition.

Regarding the PM, the detection properties of the filters in terms of fault sensitivity and disturbance rejection can be achieved according to Section 3. The synthesis of the dynamic filters for FDI has been performed by choosing a suitable linear combination of residual generator functions. This choice has to maximise the steady-state gain of the transfer functions between input sensor fault signals. The roots of the 𝑐(𝑡) polynomial matrix have been optimised for maximising the fault detection promptness, as well as to minimise the occurrence of false alarms. In order to assess the PM diagnosis technique, different fault sizes have been simulated on each sensor. Single faults in the input sensors have been generated by producing abrupt (step) and ramp (slowly developing) faults in the input signals 𝑟𝜈𝜎𝑟𝑟(𝑡)𝑟+𝜈𝜎𝑟for𝑓(𝑡)=0,𝑟(𝑡)<𝑟𝜈𝜎𝑟or𝑟(𝑡)>𝑟+𝜈𝜎𝑟for𝑓(𝑡)0,(86). The residual signals indicate fault occurrence according to whether their values are lower or higher than the thresholds fixed in fault-free conditions. The residual processing methods can be based on simple residual geometrical analysis or comparison with fixed thresholds [3]. More complex residual evaluation can rely on statistical properties of the residuals and hypothesis testing [6], or based on adaptive thresholds, that is, the so-called threshold selector [21].

In this paper, the threshold test for FDI is performed with the logic described by (86): 𝑟(𝑡) that is, the comparison of 𝑟 with respect its statistical normal characteristics. 𝜎2𝑟 and 𝜈 are the normal values for the mean and variance of the fault-free residual, respectively. In order to separate normal from faulty behaviour, the tolerance parameter 𝜈3 (normally 𝜈) is selected and properly tuned. Hence, by a proper choice of the parameter 𝜈=4, a good tradeoff can be achieved between the maximisation of fault detection probability and the minimisation of false-alarm probability. In practice, the threshold values depend on the residual error amount due to measurement errors, linearised model approximations, and disturbance signals that are not completely decoupled.

Thus, in this case, a suitable value of 𝜈=4 for the computation of the positive and negative thresholds in (86) has been considered. To summarise the performances of the PM FDI scheme, the minimal detectable step faults on the various input sensors are collected in Table 2.

𝛿 𝑒 VariableFault sizeDelay time

Elevator 2 1 8 𝛿 𝑎  s
Aileron 3 6 𝛿 𝑟  s
Rudder 4 8 𝛿 t h  s
Throttle aperture % 2 1 5 % 𝜈 = 4  s

On the other hand, the minimal detectable ramp faults are reported in Table 3.

𝛿 𝑒 Variable Fault size Delay time

Elevator 0 . 1 1 / s 2 6 𝛿 𝑎  s
Aileron 0 . 5 0 / s 1 1 𝛿 𝑟  s
Rudder 0 . 4 9 / s 1 2 𝛿 t h  s
Throttle aperture % 0 . 1 3 / s 1 9 𝐾 𝛿 t h = 1  s

Concerning the NLGA, the synthesis of the residual generators has been performed by using filter gains that optimise the fault sensitivity and reduce as much as possible the occurrence of false alarms due to model uncertainties and to disturbances not completely decoupled. This robustness requirement has been fulfilled by designing the residual gains according to the Procedure 1. For example, with reference to the fourth residual generator, Procedure 1 has led to 𝛾=1.2 which satisfies the norm bounds 𝛽=400 and 𝑓20.05. This guarantees a good separation of the residual signal with 𝑑210 and 2, where 𝜈=8-norm is considered.

In order to assess the NLGA diagnosis technique, single step and ramp faults have been used. Moreover, also in this case the threshold values have been chosen in simulation according to (86). A suitable value of 𝜈=8 for the computation of the positive and negative thresholds in (86) has been considered. For what concern NLGA FDI scheme, the minimal detectable step faults on the various input sensors are summarised in Table 4.

𝛿 𝑒 Variable Fault size Delay time

Elevator 2 5 𝛿 𝑎  s
Aileron 2 3 𝛿 𝑟  s
Rudder 2 6 𝛿 t h  s
Throttle aperture % 6 3 % 𝜈 = 8  s

On the other hand, the minimal detectable input sensor ramp faults are reported in Table 5.

𝛿 𝑒 Variable Fault size Delay time

Elevator 0 . 2 1 / s 1 1 𝛿 𝑎  s
Aileron 0 . 4 5 / s 1 2 𝛿 𝑟  s
Rudder 0 . 3 2 / s 1 0 𝛿 t h  s
Throttle aperture % 0 . 1 5 / s 1 5 4  s

The minimal detectable step fault values in Tables 2, 3, 4 and 5 are expressed in the unit of measure of the sensor signals. The fault step sizes and ramp slopes are relative to the case in which the occurrence of a fault is detected and isolated as soon as possible. The detection delay times represent the worst case results, as they are evaluated on the basis of the time taken by the slowest residual function to cross the settled threshold. These experiments represent a further validation of the residual generator robustness with respect to the fault type, as the the residual function sensitivity was optimised only with respect to step faults.

As an example, Figure 1 shows the 𝑓𝑐1(𝑡) PM residual functions generated for the complete trajectory. On the basis of the fault-free and faulty conditions, this bank provides the correct isolation of the considered input sensor ramp fault.

The horizontal lines represent the levels of the fault-free thresholds that are settled according to test (86) with 𝑓𝑐(𝑡). The first residual function depicted in Figure 1, provides also the isolation of the fault 4 regarding the 1st input sensor.

The second example of Figure 2 shows the 𝜈=12 residual functions generated by the NLGA filter bank applied to the complete aircraft trajectory. The horizontal lines represent the thresholds with 𝛿𝑒. Note that, due to the NLGA design technique, only the 1st residual related to the 𝜈=12 signal of the filter bank is sensitive to a ramp fault affecting the 1st input sensor.

5.1. Reliability and Robustness Evaluation

In this section, the robustness characteristics of the proposed PM and NLGA FDI schemes have been evaluated and compared also with respect to the UIKF scheme [3] and the NN technique [7]. In particular, a bank of UIKF has been exploited for diagnosing faults of the monitored process. This technique seems to be robust with respect to the modelling uncertainties, the system parameter variations, and the measurement noise, which can obscure the performance of an FDI system by acting as a source of false faults. The procedure recalled here requires the design of a UIKF bank and the basic scheme is the standard one: a set of measured variables of the system is compared with the corresponding signals estimated by filters to generate residual functions. The diagnosis has been performed by detecting the changes of UIKF residuals caused by a fault. The FDI input sensor scheme exploits a number of KF equal to the number of input variables. Each filter is designed to be insensitive to a different input sensor of the process and its disturbances (the so-called unknown inputs). Moreover, the considered UIKF bank was obtained by following the design technique described in [3, Section 3.5, pages 99–105], whilst the noise covariance matrices were estimated as described in [22, Section 3.3, pages 70–74 and Section 4.6, pages 130-131]. Each of the 4 UIKF of the bank was decoupled from both one input sensor fault and the wind gust disturbance component, thus providing the optimal filtering of the input-output measurement noise sequences. On the other hand, a dynamic NN bank has been exploited in order to find the dynamic connection from a particular fault regarding the input sensors to a particular residual. In this case, the learning capability of NN is used for identifying the nonlinear dynamics of the monitored plant. The dynamic NN provides the prediction of the process output with an arbitrary degree of accuracy, depending on the NN structure, its parameters, and a sufficient number of neurons. Once the NN has been properly trained, the residuals have been computed as the difference between predicted and measured process outputs. The FDI is therefore achieved by monitoring residual changes. The NN learning is typically an offline procedure. Normal operation data are acquired from the monitored plant and are exploited for the NN training. Regarding the NN FDI method and according to a generalised observer scheme (GOS) [3], a bank of 15 time-delayed three-layers multilayer perceptron (MLP) NN with 25 neurons in the input layer, 1 neurons in the hidden layer, and 4 neuron in the output layer is implemented. Each NN was designed to be insensitive to each input sensor fault, and the NN were trained in order to provide the optimal output prediction on the basis of the training pattern and target sequences [7].

In the following of this section, the performances of the different FDI schemes have been evaluated by considering a more complex aircraft trajectory. This has been obtained by means of the guidance and control functions of a standard autopilot which stabilises the aircraft motion towards the reference trajectory as depicted in Figure 3. The reference trajectory is made up of 2 branches (2 straight flights and 2 turn flights) so that a closed path is obtained. It is worth observing that only 4 steady-state flight conditions are used to follow alternatively the =50 branches of the reference trajectory: (i)straight flight condition: true airspeed = [m/s]; radius of curvature =0; flight-path anglealtitude=330; =0 [m]; flap deflection=50;(ii)turn flight condition: true airspeed =1000 [m/s]; radius of curvature =0 [m]; flight-path angle =330; altitude =0 [m]; flap deflection 𝜈.

The reference turn flight condition is used to design the PM and the NLGA filters. The achieved results are reported in Tables 2 and 4, respectively. The performed tests represent a also a possible reliability evaluation of the considered FDI techniques. In fact, in this case the diagnosis requires that the residual generators are robust with respect to the flight conditions that do not match the nominal trajectory used for the design.

As an example, the fault-free and faulty residuals generated by the designed NN and UIKF banks are shown in Figures 4 and 5, respectively.

Table 6 summarises the results obtained by considering the observers and filters (corresponding to the PM, NLGA, UIKF, and NN) for the input sensor FDI whose parameters have been designed and optimised for the steady-state coordinated turn represented by the 2nd reference flight condition of the complete trajectory. Table 6 reports the performances of the considered FDI techniques in terms of the minimal detectable step faults on the various input sensors, as well as the corresponding parameters 𝜈 for the residual evaluation of (86). The mean detection delay is also reported in Table 6 in order to compare the effectiveness of the different FDI schemes.


4 1 2 9 5 𝛿 𝑒
4 3 4 3 𝛿 𝑎
5 3 5 4 𝛿 𝑟
5 3 4 4 𝛿 t h
7 1 0 % 1 1 % 1 2 % 2 6 %
Mean detection delay 2 5  s 3 1  s 2 7  s 𝜈  s

The choice of 𝜈 has been performed with reference to the particular flight conditions involved in the complete trajectory following. In particular, the selected value of 𝜈 for each FDI observer or filter represents a tradeoff between two objectives, that is, for increasing the residual fault sensitivity and promptness, as well as for minimising the occurrence of false alarms due to the switching among the reference flight conditions needed to stabilise the aircraft motion towards the reference trajectory. Table 6 shows how the proper design of the parameter 𝜈 allows to obtain good performances with all the considered FDI schemes, hence the robustness with respect to the proposed complete trajectory is always achieved.

It is worth noting that the NLGA has a theoretical advantage by taking into account the nonlinear dynamics of the aircraft. However, the behaviour of the related nonlinear residual generators is quite sensitive to the model uncertainties due to variation of the flight condition. In fact, the NLGA FDI scheme requires high values of 8 which have to be increased (from 12 to 𝜈 in this work) when the aircraft motion regarding the complete trajectory is considered in place of the nominal flight condition. In particular, even though the analysis was restricted just to the aircraft turn phase of the complete trajectory, a performance worsening would happen, since the steady-state condition (nominal flight condition) is quite far to be reached. However, the filter design based on the NLGA lead to a satisfactory fault detection, above all in terms of promptness. On the other hand, regarding the PM, it is rather simple to note the good FDI performances, even if optimisation stages can be required. The 𝑐(𝑡) values selected for the PM are lower, but the related residual fault sensitivities are even smaller. Similar comments can be made for the UIKF and NN techniques.

The simulation model applied to the complete trajectory is an effective way to test the performances of the proposed FDI methods with respect to modelling mismatch and measurement errors. The obtained results demonstrate the reliability of the PM-, NLGA-, UIKF-, and NN-based FDI schemes as long as proper design procedures are adopted.

5.2. Monte Carlo Analysis

In this section, further experiment results have been reported. They regard the performance evaluation of the developed FDI scheme with respect to uncertainty acting on the system. Hence, the simulation of different fault-free and faulty data sequences was performed by exploiting the aircraft Matlab-Simulink simulator and a Monte Carlo analysis implemented in the Matlab environment. The Monte Carlo tool is useful at this stage as the FDI performances depend on the residual error magnitude due to the system uncertainty, as well as the signal 𝑦(𝑡) and 𝑝 measurement errors. It is worth noting how the Monte Carlo simulations have been achieved by perturbing the parameters of the PM filter residuals by additive white Gaussian noises with standard deviation values equal to a fixed percentage 1000 of the element values. The same experiments have been performed by statistically varying the main parameters of the NLGA filters. In these conditions, the Monte Carlo analysis represents a further method for estimating the reliability and the robustness of the developed FDI schemes, when applied to the considered aircraft.

For robustness and reliability experimental analysis of the FDI schemes, some performance indices have been used. The performances of the FDI method are then evaluated on a number of Monte Carlo runs equal to 𝑟fa. This number of simulations is carried out to determine the indices listed below with a given degree of accuracy.

False-Alarm Probability (𝑟mf): the number of wrongly detected faults divided by total fault cases.

Missed-Fault Probability (𝑟td): for each fault, the total number of undetected faults, divided by the total number of times that the fault case occurs.

True Detection/Isolation Probability (𝑟ti, 𝜏md): for a particular fault case, the number of times it is correctly detected/isolated, divided by total number of times that the fault case occurs.

Mean Detection/Isolation Delay (𝜏mi, 𝑝=10): for a particular fault case, the average detection/isolation delay time.

These indices are hence computed for the number of Monte Carlo simulations and for each fault case. Table 7 summarises the results obtained by considering the PM dynamic filters for the input sensor FDI for a complete aircraft trajectory and with