Two-Stage Hyperelliptic Kalman Filter-Based Hybrid Fault Observer for Aeroengine Actuator under Multi-Source Uncertainty

Abstract

An aeroengine faces multi-source uncertainty consisting of aeroengine epistemic uncertainty and the control system stochastic uncertainty during operation. This paper investigates actuator fault estimation under multi-source uncertainty to enhance the fault diagnosis capability of aero-engine control systems in complex environments. With the polynomial chaos expansion-based discrete stochastic model quantification, the optimal filter under multi-source uncertainty, the Hyperelliptic Kalman Filter, is proposed. Meanwhile, by treating actuator fault as unknown input, the Two-stage Hyperelliptic Kalman Filter (TSHeKF) is also proposed to achieve optimal fault estimation under multi-source uncertainty. However, considering that the biases of the model are often fixed for the individual, the TSHeKF-based fault estimation is robust and leads to inevitable conservativeness. By adding the additional estimation of the unknown deviation in state function caused by probabilistic system parameters, the hybrid fault observer (HFO) is proposed based on the TSHeKF and realizes conservativeness-reduced estimation for actuator fault under multi-source uncertainty. Numerical simulations show the effectiveness and optimality of the proposed HFO in state estimation, output prediction, and fault estimation for both single and multi-fault modes, when considering multi-source uncertainty. Furthermore, Monte Carlo experiments have demonstrated that the HFO-based optimal fault estimation is less conservative and more accurate than the Two-stage Kalman Filter and TSHeKF, providing better safety and more reliable aeroengine operation assurance.

1. Introduction

With the continuous development of aerospace technology, pursuing better performance also requires higher robustness and reliability of aeroengines [1]. In addition, the safety of aeroengine operation can be significantly assured by the contribution of fault diagnosis systems in detecting, isolating, and estimating faults [2,3]. Meanwhile, considering that uncertainty is an essential feature of complex systems [4], it can be classified as stochastic and epistemic uncertainty [5]. According to [6], the control system stochastic uncertainty and aeroengine epistemic uncertainty form the multi-source uncertainty affecting aeroengine operation. However, the traditional design process of the fault diagnosis system is based on a deterministic model and ignores the widespread uncertainty effects [2]. With the development of UQ, fault diagnosis of aeroengine control systems under multi-source uncertainty is gradually becoming a research hotspot, and fruitful results have been achieved in sensor fault diagnosis [6] and estimation [7]. But more research on actuator fault estimation is needed. Therefore, to improve the performance and effectiveness of the diagnostic system to ensure the safety of aeroengine operation, it is of great importance to consider multi-source uncertainty while designing the actuator fault estimation system.

Due to the higher interpretability and reliability of model-based methods, fruitful research results of actuator fault estimation are achieved by designing appropriate filters or observers [8,9]. The LKF is the most classical state-estimation method for linear dynamic systems, and Simon et al. designed a Kalman filter bank for the fault diagnosis of the aeroengine control system in 2003 [10,11]. It is widely known that the health management of aeroengines is mainly achieved by treating the unknown health parameters input as state and estimating them with the corresponding ASKF [12,13]. The actuator fault can also be treated as unknown input and is estimated in the same way as the health parameters. However, as the number of states increases, the computational effort of the ASKF increases geometrically [14]. To reduce computation ford while guaranteeing state estimation accuracy, TSKF, composed of two sub-filters, is proposed by introducing U-V transformations [14]. By extending the sub-filter of the TSKF from LKF to the extended Kalman filter [15] or unscented Kalman filter [16], the enhanced TSKF can realize more accurate fault estimation in a more extensive range. However, the above methods need to consider uncertainties sufficiently, and more advanced fault estimation methods are required for more complex environments [17,18]. Under non-Gaussian noise conditions, Lu et al. combined the particle swarm algorithm with LKF to propose the Kalman particle filter for state estimation and used it for aeroengine health management [19]. With the development of UQ, Combastel first proposed the ZKF and achieved optimal estimation of the state motion boundary [20]. And in 2021, Zhang et al. used ZKF for fault interval estimation [21]. However, the ZKF-based estimation is an interval estimation method and cannot characterize the distribution of states. PCE, the most popular way for UQ today, can adequately characterize the statistical properties of the state [22] and is increasingly used by scholars for state estimation [23,24], fault diagnosis [6] and performance monitoring [25].

In summary, research on fault estimation of aeroengine actuators under multi-source uncertainty needs to be revised, especially considering both stochastic processes and systems. On the foundation of PCE-based UQ of the discrete stochastic model, HeKF is proposed in this paper and realizes optimal filtering under multi-source uncertainty. Further, the TSHeKF is proposed and achieves optimal actuator fault estimation under multi-source uncertainty. Moreover, by adding estimation of the unknown state deviation, the HFO is established on the basis of the TSHeKF; it provides fault estimation with better stability and reliability and is less conservative.

This paper is organized as follows. The mathematical description of aeroengine under multi-source uncertainty is presented in Section 2. PCE-based discrete stochastic model quantification and the further optimal filtering under multi-source uncertainty are described in Section 3. The TSHeKF is proposed in Section 4, and the HFO is designed and realized conservativeness-reduced fault estimation. Numerical and MC simulations are presented in Section 5. Some brief conclusions are shown in Section 6.

Notation: In this article, $ℝ$ represents the set of real numbers, and ${ℝ}^n$ and ${ℝ}^{n\times m}$ denote real vectors and real matrices of appropriate dimensions shown in superscript, respectively. $O$ and $I$ represent the null matrix and identity square matrix, respectively; dimensions are also shown in the corresponding superscript. The bounded basic function set is represented by $\mho$, and $L^2$ represents the $L2$ space. For a given real vector $Y\in {ℝ}^N$, its square is written as $Y^2$ and is defined as follows. Using $soe$ to represent the sum of each element, for a given matrix $A\in {ℝ}^{N\times M}$, $soe(A)$ has the following expression.

(1)$\begin{array}{l}soe(A)={\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^M{a}_{ij}}}\\ Y^2=Y{Y}^T\end{array}$

2. Mathematical Description of the Problem

From [6], the multi-source uncertainty faced by aeroengines consists of the stochastic uncertainty of control system and the epistemic uncertainty of aeroengine. Based on uncertainty modeling of one type or a batch of aeroengines [26], the linear dynamical system considers the actuator fault with independent probabilistic parameters, as shown in Equation (2):

(2)$\left\{\begin{array}{c}\dot{x}(t)=\tilde{A}x(t)+\tilde{B}u(t)+Ef(t)+M{w}^m(t)\hfill \\ y(t)=\tilde{C}x(t)+\tilde{D}u(t)+Ff(t)+Nv(t)\hfill \\ u(t)=m^u(t)+f[m^u(t)]w^u(t)\hfill \end{array}\right.$

for $t\ge 0$ with state functions $x:ℝ\to {ℝ}^{N_{state}}$; similarly, input functions $u:ℝ\to {ℝ}^{N_{in}}$ and output function $y:ℝ\to {ℝ}^{N_{out}}$. $w^m:ℝ\to {ℝ}^{N_m}$, $v:ℝ\to {ℝ}^{N_v}$ and $w^u:ℝ\to {ℝ}^{N_u}$ are continuous independent disturbance functions, which represent system noise, measurement noise and disturbance of actuator, respectively. $m^u(t)$ is the expectation of $u(t)$ and $f[m^u(t)]$ is the deterministic function which characterize the intensity of actuator stochastic uncertainty. $\tilde{A}$, $\tilde{B}$, $\tilde{C}$, and $\tilde{D}$ denote the control matrices of appropriate dimensions with probabilistic system parameters which satisfy the controllability and the observability. $f:ℝ\to {ℝ}^{N_f}$ is the additive actuator fault vector; $E$ and $F$ are fault-drive matrices with constant parameters. $M$ and $N$ are constant drive matrices; especially when the disturbance of each state or output is independent, $M$ and $N$ become diagonal matrices. Use $w(t)$ to represent $w^u(t)$, $w^m(t)$, and $v(t)$, which is a continuous zero-mean random process with independent increments. According to [27], $w(t)$ can only be expressed as $dB(t)/dt$, where $B(t)$ is the Brownian motion. Using $\Delta t$ as the sampling time, the linear discrete stochastic dynamical system with actuator fault is shown in Equation (3):

(3)$\left\{\begin{array}{c}{x}_{k+1}={\tilde{A}}^d{x}_k+{\tilde{B}}^d{u}_k+{\tilde{E}}^d{f}_k+M{w}_k^m\hfill \\ y_k={\tilde{C}}^d{x}_k+{\tilde{D}}^d{u}_k+F^d{f}_k+N{v}_k/\sqrt{\Delta t}\hfill \\ u_k=m_k^u+M_k^u{w}_k^u/\sqrt{\Delta t}\hfill \end{array}\right.$

where $w_k^u=B^u(k+1)-B^u(k)~N(0_{N_u},\Delta t{I}_{N_u})$; similarly, $w_k^m~N(0_{N_m},\Delta t{I}_{N_m})$ and $v_k~N(0_{N_v},\Delta t{I}_{N_v})$ obey the mutually independent multi-source Gaussian distribution. $M_k^u$ is used to substitute $f(m_k^u)$; the superscript $d$ represents the discretization of the system matrix, which is obtained by Equation (4), where $p$ is the expansion order.

(4)$\begin{array}{ccc}{\tilde{A}}^d=e^{\tilde{A}\Delta t}\approx {\displaystyle \sum _{i=0}^p\frac{1}{i!}{\tilde{A}}^i{(\Delta t)}^i}& {\tilde{C}}^d=\tilde{C}& {\tilde{E}}^d=\widetilde{IA}\times E\\ {\tilde{B}}^d={\displaystyle {\int }_0^{\Delta t}{e}^{\tilde{A}t}dt}\tilde{B}=\widetilde{IA}\times \tilde{B}& {\tilde{D}}^d=\tilde{D}& F^d=F\end{array}$

Temporarily ignoring the parameter uncertainty of the system matrix in Equation (3), the actuator failure can be optimally estimated based on either the ASKF [13] or TSKF [15] by treating $f_k$ as an augmented state. Moreover, the LKF, which is the basis for the ASKF and TSKF, needs to be replaced by optimal filtering considering multi-source uncertainty. However, since the above model is for a type or a batch of aeroengines rather than for a specific engine, the fault estimation based on Equation (3) is also optimal for a range of aeroengines, and robust for an individual one. In order to improve the accuracy of fault estimation for an individual aeroengine, the adaptiveness of the fault estimation needs to be enhanced. The linear discrete stochastic system for one specific aeroengine is shown in Equation (5):

(5)$\left\{\begin{array}{c}{x}_{k+1}=({\overline{A}}^d+\Delta A^d)x_k+({\overline{B}}^d+\Delta B^d)u_k+({\overline{E}}^d+\Delta E^d)f_k+M{w}_k^m\hfill \\ y_k=({\overline{C}}^d+\Delta C^d)x_k+({\overline{D}}^d+\Delta D^d)u_k+F^d{f}_k+N{v}_k/\sqrt{\Delta t}\hfill \\ u_k=m_k^u+M_k^u{w}_k^u/\sqrt{\Delta t}\hfill \end{array}\right.$

where the overscore indicates the expectation of the stochastic matrix below, and $\Delta$ represents the differences of the following control matrix between the aimed aeroengine and the expectation of the corresponding type or batch of aeroengines, which can also be seen as the uncertainty of the discrete model. It is obvious that both the state function and output function are affected. In order to improve the adaptiveness of the state estimation, the effect of model uncertainty on the state and output needs to be quantified. $d_k^x=\Delta A^d{x}_k+\Delta B^d{u}_k+\Delta B_f^d{f}_k$ and $d_k^y=\Delta C^d{x}_k+\Delta D^d{u}_k$ are used to represent the bias of state and output caused by model uncertainty, respectively. By treating $d_k^x$ as unknown input, the effect of model uncertainty on the state can be estimated by the TSKF optimally. Moreover, since $d_k^y$ is additive to the output and does not affect the state, it is impossible to estimate $d_k^y$ by considering it as the augmented state. Considering that the optimal output estimation of an individual aeroengine is unachievable, the optimal output estimation for a range of aeroengines is the best available option. Therefore, the conservativeness-reduced actuator fault estimation system for individual aeroengine is designed optimally based on Equation (6).

(6)$\left\{\begin{array}{c}{x}_{k+1}={\overline{A}}^d{x}_k+{\overline{B}}^d{u}_k+{\overline{E}}^d{f}_k+d_k^x+M{w}_k^m\hfill \\ y_k={\tilde{C}}^d{x}_k+{\tilde{D}}^d{u}_k+F^d{f}_k+N{v}_k/\sqrt{\Delta t}\hfill \\ u_k=m_k^u+M_k^u{w}_k^u/\sqrt{\Delta t}\hfill \end{array}\right.$

To sum up, when considering multi-source uncertainty, TSKF, as the classical fault-estimation method when treating unknown inputs, is no longer feasible. TSKF, therefore, needs to be adapted from the optimal filtering under multi-source uncertainty to obtain optimal actuator fault estimation. However, directly using the optimal state estimation and output prediction for a batch of aeroengines would result in an unnecessary loss of accuracy by completely ignoring the differences between individuals. Because of the existence of $\Delta C$ and $\Delta D$, the optimal output prediction for individual aeroengine does not exist. Therefore, by combining the adaptive state estimation for individual aeroengine with the optimal output estimation for the corresponding range of aeroengines, the accuracy of actuator fault estimation is maximized under multi-source uncertainty. With the optimal estimation of $d_k^x$, the conservativeness-reduced estimation of $f_k$ is achieved based on the proposed HFO, which is the main contribution of this paper.

3. Optimal Filtering under Multi-Source Uncertainty

In this section, based on PCE, the statistical property of control matrices with independent probabilistic system parameters is quantified. Furthermore, the optimal state estimation of linear discrete stochastic dynamical system is provided. In addition, the actuator fault is temporarily not discussed.

3.1. PCE-Based UQ of Discrete Stochastic Model

3.1.1. Preliminaries of PCE

Considering a stochastic model $Y=f(\xi )$, with the probabilistic input $\xi$, which belongs to the probability space $(\mathsf{\Omega },\mathrm{A},\mathsf{\Gamma })$ [27]. For stochastic analysis, the expansion of the model $f(\xi )$ and the statistical property of $Y$ are as follows:

(7)$\begin{array}{l}Y(\xi )\approx {\displaystyle \sum _{k=0}^d{c}_k{P}^{(k)}(\xi )}\\ E(Y)=c_{0}V(Y)={\displaystyle \sum _{k=1}^d{c}_k^2}\end{array}$

where $d$ is the expansion order, $c_k$ is the expansion coefficient, which can be obtained by Galerkin projection [28,29], numerical integration [30], or collocation [31,32]. $P^{(k)}(\xi )$ is the polynomial of orthogonal basis $\left\{P^{(0)},\dots P^{(d)}\right\}$ with the respect to the measure $\mathsf{\Gamma }$. $E(Y)$ and $V(Y)$ are the expectation and variance of $Y$, which can be directly estimated by $c_k$.

Most realistic power systems combine multi random input and multi output, especially regarding the aeroengine. For $Y={[Y_1,\dots Y_m]}^T=f({\xi }_1,\dots {\xi }_N)$ with $N$ dimensional inputs and $m$ dimensional outputs, as an extension of one-dimensional PCE, the statistical property of $Y$ can be quantized as follows:

(8)$\begin{array}{l}Y({\xi }_1,\dots ,{\xi }_N)\approx {\displaystyle \sum _{k=0}^M{c}_k{\mathsf{\Phi }}_k({\xi }_1,\dots ,{\xi }_N)}\\ E(Y)=c_0\in {ℝ}^{m\times 1}\mathrm{cov}(Y_p,Y_q)={\displaystyle \sum _{k=1}^M{(c_k)}_p{(c_k)}_q}\\ \mathsf{\Psi }(Y)={\displaystyle \sum _{k=1}^M{c}_k{c}_k^T}\in {ℝ}^{m\times m}E(Y^2)=E{(Y)}^2+\mathsf{\Psi }(Y)\end{array}$

where ${\mathsf{\Phi }}_k({\xi }_1,\dots {\xi }_N)$ is the multi-variate polynomial orthonormal basis with $M$ as dimension. From [33], $M$ is determined by the dimension of the random inputs and the expansion order. In order to harmonize with one-dimensional PCE, ${\mathsf{\Phi }}_0({\xi }_1,{\xi }_2,\dots {\xi }_N)$ is defined as multi-dimensional basis with zero-order expansion [33]. $\mathrm{cov}(Y_p,Y_q)$ is the correlation between $pth$ and $qth$ output, and $\mathsf{\Psi }(Y)$ is the correlation matrix of the random vector.

To quantify the statistical property of the discrete model with probabilistic system parameters, the output of $f({\xi }_1,\dots {\xi }_n)$ is further extended to $Y\in {ℝ}^{m\times n}$. It is obvious that the expansion coefficients have the same size as the output. Moreover, considering the matrix computation demand in this paper, the correlation of the random matrix is indicated by the correlation matrix between row vectors or column vectors. The statistical property of $Y$ is expressed by ${\left\{c_k\right\}}_{k=0-M}$ as follows:

(9)$\begin{array}{l}E(Y)=c_0\in {ℝ}^{m\times n}\\ {\mathsf{\Psi }}_{ij}(Y)=\mathrm{cov}(Y_i,Y_j)={\displaystyle \sum _{k=1}^M{(c_k)}_i{(c_k)}_j^T}\in {ℝ}^{n\times n}\\ E(Y_i{Y}_j^T)=E{(Y)}_{i}E{(Y)}_j^T+{\mathsf{\Psi }}_{ij}(Y)\\ {\mathsf{\Psi }}_{ij}^{\prime }(Y)=\mathrm{cov}(Y_{i^{\prime }},Y_{j^{\prime }})={\displaystyle \sum _{k=1}^M{(c_k)}_{i^{\prime }}^T{(c_k)}_{j^{\prime }}}\in {ℝ}^{m\times m}\\ E(Y_{i^{\prime }}^T{Y}_{j^{\prime }})=E{(Y)}_{i^{\prime }}^{T}E{(Y)}_{j^{\prime }}+{\mathsf{\Psi }}_{ij}^{\prime }(Y)\end{array}$

where $E(Y)$ is the expectation, ${\mathsf{\Psi }}_{ij}(Y)$ is the correlation matrix between $ith$ and $jth$ columns, and ${\mathsf{\Psi }}_{ij}^{\prime }(Y)$ is based on row vectors. The subscript $i$ and $i^{\prime }$ represent the order of the column and row, respectively.

The statistical property of the product between the vector or matrix with probabilistic system parameters is discussed. With the assumption of mutual independence, the expectation of the product is simply equal to product of the individual expectation. Moreover, the correlation of the product is expressed in Equation (10) by extending the formula for variance of multiplied random variables [27]. Here, $\ast$ is Hadamard product, $a,b\in {ℝ}^{n\times 1}$ are random vectors with statistical property shown in Equation (8), and $A\in {ℝ}^{m\times n}$, $B\in {ℝ}^{n\times r}$ are random matrixes with statistical property obey Equation (9).

(10)$\begin{array}{cc}\hfill V(a^{T}b)=& soe[E(a^2)\ast E(b^2)-E{(a)}^2\ast E{(b)}^2]\hfill \\ \hfill \mathsf{\Psi }{(Ab)}_{lm}=& soe[E(A_{l^{\prime }}^T{A}_{m^{\prime }})\ast E(b^2)-E{(A_{l^{\prime }})}^{T}E(A_{m^{\prime }})\ast E{(b)}^2]\hfill \\ \hfill =& E(A_{l^{\prime }})\mathsf{\Psi }(b)E(A_{m^{\prime }}^T)+soe[{{\mathsf{\Psi }}^{\prime }}_{lm}(A)\ast E(b^2)]\hfill \\ \hfill {\mathsf{\Psi }}_{ij}{(AB)}_{lm}=& soe[E(A_{l^{\prime }}^T{A}_{m^{\prime }})\ast E(B_i{B}_j^T)-E{(A_{l^{\prime }})}^{T}E(A_{m^{\prime }})\ast E(B_i)E{(B_j)}^T]\hfill \\ \hfill {\mathsf{\Psi }}_{ij}^{\prime }{(AB)}_{lm}=& soe[E(A_{i^{\prime }}^T{A}_{j^{\prime }})\ast E(B_l{B}_m^T)-E{(A_{i^{\prime }})}^{T}E(A_{j^{\prime }})\ast E(B_l)E{(B_m)}^T]\hfill \end{array}$

3.1.2. Statistical Property of Discrete Stochastic Model

The statistical property of the stochastic model in Equation (3) is quantified based on multi-variate PCE. Considering the mutual independence between different outputs, the probabilistic parameters in $C$ and $D$ can be treated by row. Additionally, the multi-variate polynomial orthogonal basis of the dynamic model has the following form:

(11)$\begin{array}{cc}\tilde{A}\to \left\{{\mathsf{\Phi }}_{0-M_A}^{\tilde{A}}(a_{11}\dots a_{NN})\right\}={\left\{{\mathsf{\Phi }}_j^{\tilde{A}}\right\}}_{j=0-M_A}& \tilde{B}\to \left\{{\mathsf{\Phi }}_{0-M_A}^{\tilde{B}}(b_{11}\dots b_{N{N}_{in}})\right\}={\left\{{\mathsf{\Phi }}_j^{\tilde{B}}\right\}}_{j=0-M_B}\\ {\tilde{C}}_{i^{\prime }}\to \left\{{\mathsf{\Phi }}_{0-M_C}^{{\tilde{C}}_{i^{\prime }}}(c_{i1}\dots c_{iN})\right\}={\left\{{\mathsf{\Phi }}_j^{{\tilde{C}}_{i^{\prime }}}\right\}}_{j=0-M_C}& {\tilde{D}}_{i^{\prime }}\to \left\{{\mathsf{\Phi }}_{0-M_D}^{{\tilde{D}}_{i^{\prime }}}(d_{i1}\dots d_{i{N}_{in}})\right\}={\left\{{\mathsf{\Phi }}_j^{{\tilde{D}}_{i^{\prime }}}\right\}}_{j=0-M_D}\end{array}$

where $i^{\prime }$ represent the $ith$ row of the matrix, and $j$ is the sort order with the max number of $M_A$, $M_B$, $M_C$, and $M_D$ for each multi-variate polynomial orthonormal basis. Based on Equation (4), the discrete stochastic model is expressed by the random linear model and is determined by the relevant expansion coefficients. From [30], based on numerical integration, the expansion coefficient is calculated as follows:

(12)$\begin{array}{c}{c}_k^{{\tilde{A}}^d}={\displaystyle \int \dots }{\displaystyle {\int }_S{e}^{\tilde{A}\Delta t}{\mathsf{\Phi }}_k^{\tilde{A}}d\varphi (a_{11}\dots a_{NN})}\in {ℝ}^{N\times N}\hfill \\ c_k^{\widetilde{IA}}={\displaystyle \int \dots }{\displaystyle {\int }_S{\displaystyle {\int }_0^{\Delta t}{e}^{\tilde{A}t}dt}{\mathsf{\Phi }}_k^{\tilde{A}}d\varphi (a_{11}\dots a_{NN})}\in {ℝ}^{N\times N}\hfill \\ c_k^{\tilde{B}}={\displaystyle \int \dots }{\displaystyle {\int }_S\tilde{B}{\mathsf{\Phi }}_k^{\tilde{B}}d\varphi (b_{11}\dots b_{N{N}_{in}})}\in {ℝ}^{N\times N_{in}}\hfill \\ c_k^{{\tilde{C}}_{i^{\prime }}^d}={\displaystyle \int \dots }{\displaystyle {\int }_S{\tilde{C}}_{i^{\prime }}^d{\mathsf{\Phi }}_k^{{\tilde{C}}_{i^{\prime }}}d\varphi (c_{i1}\dots c_{iN})\in {ℝ}^{1\times N}}\hfill \\ c_k^{{\tilde{D}}_{i^{\prime }}^d}={\displaystyle \int \dots }{\displaystyle {\int }_S{\tilde{D}}_{i^{\prime }}^d{\mathsf{\Phi }}_k^{{\tilde{D}}_{i^{\prime }}}d\varphi (d_{i1}\dots d_{i{N}_{in}})}\in {ℝ}^{1\times N_{in}}\hfill \end{array}$

where $k$ is the sort order and $S$ is support of the arbitrary measure $\varphi$ [30,33]. Based on Equation (9), the statistical property of ${\tilde{A}}^d$ is quantified by ${\left\{c_k^{{\tilde{A}}^d}\right\}}_{k=1-M_A}$. By respectively combining $c_k^{{\tilde{C}}_{i^{\prime }}^d}$ and $c_k^{{\tilde{D}}_{i^{\prime }}^d}$, the expansion coefficient matrix ${\left\{c_k^{{\tilde{C}}^d}\right\}}_{k=1-M_C}$ and ${\left\{c_k^{{\tilde{D}}^d}\right\}}_{k=1-M_D}$, which determines the statistical property of ${\tilde{C}}^d$ and ${\tilde{D}}^d$, is shown in Equation (13).

(13)$\begin{array}{cc}{c}_k^{{\tilde{C}}^d}=[c_k^{{{\tilde{C}}^d}_1};c_k^{{{\tilde{C}}^d}_2};\dots c_k^{{{\tilde{C}}^d}_{N_{out}}}]\in {ℝ}^{N_{out}\times N}& c_k^{{\tilde{D}}^d}=[c_k^{{{\tilde{D}}^d}_1};c_k^{{{\tilde{D}}^d}_2};\dots c_k^{{{\tilde{D}}^d}_{N_{out}}}]\in {ℝ}^{N_{out}\times N_{in}}\end{array}$

From Equation (4), considering $E$ is a constant matrix, the randomness of ${\tilde{E}}^d$ is entirely caused by $\widetilde{IA}$, and the expansion coefficient matrix ${\left\{c_k^{{\tilde{E}}^d}\right\}}_{k=1-M_A}$, which expresses the statistical property of ${\tilde{E}}^d$, has the following form:

(14)$c_k^{{\tilde{E}}^d}=c_k^{\widetilde{IA}}{\tilde{E}}^d\in {ℝ}^{N\times N_f}$

Finally, since ${\tilde{B}}^d=\widetilde{IA}\times \tilde{B}$, the statistical property of ${\tilde{B}}^d$ can only be quantified based on Equation (10) with the statistical property of $\widetilde{IA}$ and $\tilde{B}$, which can be obtained by ${\left\{c_k^{\widetilde{IA}}\right\}}_{k=1-M_A}$ and ${\left\{c_k^{\tilde{B}}\right\}}_{k=1-M_A}$ and has the following form.

(15)$\begin{array}{c}\hspace{1em}E({\tilde{B}}^d)=E(\widetilde{IA})E(\tilde{B})=c_0^{\widetilde{IA}}{c}_0^{\tilde{B}}=c_0^{{\tilde{B}}^d}\in {ℝ}^{N\times N_{in}}\hfill \\ {\mathsf{\Psi }}_{ij}{({\tilde{B}}^d)}_{lm}=soe[{\displaystyle \sum _{k=0}^{M_A}{(c_k^{\widetilde{IA}})}_{l^{\prime }}^T{(c_k^{\widetilde{IA}})}_{m^{\prime }}}\ast {\displaystyle \sum _{k=0}^{M_B}{(c_k^{\tilde{B}})}_i{(c_k^{\tilde{B}})}_j^T}]-soe[{(c_0^{\widetilde{IA}})}_{l^{\prime }}^T{(c_0^{\widetilde{IA}})}_{m^{\prime }}\ast {(c_0^{\tilde{B}})}_i{(c_0^{\tilde{B}})}_j^T]\hfill \\ {\mathsf{\Psi }}_{ij}^{\prime }{({\tilde{B}}^d)}_{lm}=soe[{\displaystyle \sum _{k=0}^{M_A}{(c_k^{\widetilde{IA}})}_{i^{\prime }}^T{(c_k^{\widetilde{IA}})}_{j^{\prime }}}\ast {\displaystyle \sum _{k=0}^{M_B}{(c_k^{\tilde{B}})}_l{(c_k^{\tilde{B}})}_m^T}]-soe[{(c_0^{\widetilde{IA}})}_{i^{\prime }}^T{(c_0^{\widetilde{IA}})}_{j^{\prime }}\ast {(c_0^{\tilde{B}})}_l{(c_0^{\tilde{B}})}_m^T]\hfill \end{array}$

3.2. Hyperelliptic Kalman Filter

For discrete stochastic system under multi-source uncertainty in Equation (3), the optimal state estimation is obtained by HeKF [6] and has the following form.

(16)$\left\{\begin{array}{c}{\widehat{x}}_{k|k-1}={\tilde{A}}^d{\widehat{x}}_{k-1|k-1}+{\tilde{B}}^d{u}_{k-1}+M{w}_{k-1}^m\hfill \\ {\widehat{x}}_{k|k}={\widehat{x}}_{k|k-1}+K_k(y_k-{\tilde{C}}^d{\widehat{x}}_{k|k-1}-{\tilde{D}}^d{u}_k-N{v}_k/\sqrt{\Delta t})\hfill \end{array}\right.$

where ${\widehat{x}}_{k|k-1},{\widehat{x}}_{k|k}\in {ℝ}^N$ are the one-step prediction and the optimal estimation of the state, which are both random processes due to the multi-source uncertainty. Based on Equation (8), the expectation and the correlation matrix together express the statistical property of the random vector, and $E({\widehat{x}}_{k-1|k-1})$ and $\mathsf{\Psi }({\widehat{x}}_{k-1|k-1})$ are assumed to be known. $K_k$ is the constant estimated gain matrix, which is optimized later. The statistical property of ${\widehat{x}}_{k|k-1}$ can be quantified by Equation (17).

(17)$\begin{array}{l}E({\widehat{x}}_{k|k-1})={\overline{A}}^d{\widehat{x}}_{k-1|k-1}+{\overline{B}}^d{m}_{k-1}^u\\ \mathsf{\Psi }({\widehat{x}}_{k|k-1})=\mathsf{\Psi }({\tilde{A}}^d{\widehat{x}}_{k-1|k-1})+\mathsf{\Psi }({\tilde{B}}^d{u}_{k-1})+\mathsf{\Psi }(M{w}_{k-1}^m)\end{array}$

where the $E$, $\mathsf{\Psi }$, and ${\mathsf{\Psi }}^{\prime }$ for both ${\tilde{A}}^d$ and ${\tilde{B}}^d$ are quantified in Section 3.1; the expectation and correlation matrix for $u_k$ is presented in Equation (18).

(18)$E(u_k)=m_k^u\mathsf{\Psi }(u_k)=M_k^u{(M_k^u)}^T$

Next, the statistical property of ${\widehat{x}}_{k|k}$ is quantified as follows:

(19)$\begin{array}{l}E({\widehat{x}}_{k|k})=E({\widehat{x}}_{k|k-1})+K_k[y_k-{\overline{C}}^{d}E({\widehat{x}}_{k|k-1})-{\overline{D}}^d{m}_k^u]\\ \mathsf{\Psi }({\widehat{x}}_{k|k})=\mathsf{\Psi }[(I_N-K_k{\tilde{C}}^d){\widehat{x}}_{k|k-1}]+K_k{S}_k{K}_k^T\end{array}$

where $S_k=\mathsf{\Psi }({\tilde{D}}^d{u}_k)+\mathsf{\Psi }(N{v}_k/\sqrt{\Delta t})$, consistent with the case in Equation (17), the statistical properties of all the random parameters in Equation (19) can be quantified. $E({\widehat{x}}_{k|k})$ is the optimal estimation of the state, and $\mathsf{\Psi }({\widehat{x}}_{k|k})$ represent the disturbance intensity. Therefore, in order to obtain the optimal state estimation, the optimal gain matrix $K_k$ is required to minimize $\mathsf{\Psi }({\widehat{x}}_{k|k})$, since $\mathsf{\Psi }({\widehat{x}}_{k|k})$ is a real symmetric positive definite matrix, and the diagonalization, $\overline{\mathsf{\Psi }}({\widehat{x}}_{k|k})=U_k\mathsf{\Psi }({\widehat{x}}_{k|k})U_k^T$, is introduced based on the transfer matrix $U_k$ which is determined by $\mathsf{\Psi }({\widehat{x}}_{k|k})$. Naturally, the trace of $\overline{\mathsf{\Psi }}({\widehat{x}}_{k|k})$ can be used as the evaluation of the disturbance intensity. Using $J_k$ as the optimization objective, it can be expressed as $J_k=tr[\overline{\mathsf{\Psi }}({\widehat{x}}_{k|k})]$. For calculation convenience, $\overline{\mathsf{\Psi }}({\widehat{x}}_{k|k})$ can be simplified as follows based on Equation (10).

(20)$\begin{array}{l}\mathsf{\Psi }({\widehat{x}}_{k|k})=(I_N-K_k{\overline{C}}^d)\mathsf{\Psi }({\widehat{x}}_{k|k-1}){(I_N-K_k{\overline{C}}^d)}^T+K_k(S_k+R_k)K_k^T\\ {(R_k)}_{lm}=soe[{\mathsf{\Psi }}_{lm}^{\prime }({\tilde{C}}^d)\ast E({\widehat{x}}_{k|k-1}^2)]\end{array}$

It is obvious that $K_k$ is the only unknown parameter in $J_k$. Since $\mathsf{\Psi }({\widehat{x}}_{k|k})$ is a real symmetric positive definite matrix, based on the calculation rules, $J_k$ is minimum only when $\partial J_k/\partial K_k=0$, and the corresponding $K_k$ becomes the optimal gain matrix. Considering $U_k$ is not affect by $K_k$, and thus, $K_k$ satisfying $\partial tr[\mathsf{\Psi }({\widehat{x}}_{k|k})]/\partial K_k=0$, also satisfies $\partial J_k/\partial K_k=0$, and $K_k$ has the following form:

(21)$\begin{array}{l}\partial tr[\mathsf{\Psi }({\widehat{x}}_{k|k})]/\partial K_k=2[L_k{K}_k^T-{\overline{C}}^d\mathsf{\Psi }({\widehat{x}}_{k|k-1})]\\ L_k={\overline{C}}^d\mathsf{\Psi }({\widehat{x}}_{k|k-1}){({\overline{C}}^d)}^T+S_k+R_k\\ K_k=\mathsf{\Psi }({\widehat{x}}_{k|k-1}){({\overline{C}}^d)}^T{L}_k^{-1}\end{array}$

Moreover, combining Equations (19) and (21), $\mathsf{\Psi }({\widehat{x}}_{k|k})$ can be simplified as follows:

(22)$\mathsf{\Psi }({\widehat{x}}_{k|k})=(I_N-K_k{\overline{C}}^d)\mathsf{\Psi }({\widehat{x}}_{k|k-1})$

Equations (17), (19), (21) and (22) together achieved the optimal state estimation considered multi-source uncertainty; using ${\widehat{y}}_k\in {ℝ}^{N_{out}}$ as the optimal prediction of output, the statistical property of ${\widehat{y}}_k$ can be expressed as Equation (23), where the statistical properties of all random variables have been quantified, and the calculation rules are presented in Equation (10).

(23)$\begin{array}{l}E({\widehat{y}}_k)={\overline{C}}^{d}E({\widehat{x}}_{k|k})+{\overline{D}}^d{m}_k^u\\ \mathsf{\Psi }({\widehat{y}}_k)=\mathsf{\Psi }({\tilde{C}}^d{\widehat{x}}_{k|k})+\mathsf{\Psi }({\tilde{D}}^d{u}_k)+\mathsf{\Psi }(N{v}_k/\sqrt{\Delta t})\end{array}$

4. Conservativeness-Reduced Fault Estimation under Multi-Source Uncertainty

In this section, the traditional TSKF-based fault estimation ideas are briefly described. Meanwhile, the TSKF is extended to the TSHeKF to provide optimal fault estimation of fault considered multi-source uncertainty. Moreover, by adding quantification of $d_k^x$, the HFO is established on the basis of the TSHeKF, and the conservativeness-reduced fault estimation is achieved.

4.1. TSKF-Based Fault Estimation

By treating fault as unknown input, the TSKF is widely used in fault estimation [15,16], and the data flow of TSKF-based fault optimal estimation is shown in Figure 1. Where $x_k^a$ is the augmented state, which is composed of state and fault, ${\widehat{x}}_{k|k}^a$ is the ASKF-based optimal estimation. ${\widehat{x}}_{k|k}^{\overline{a}}$ is the augmentation of ${\widehat{x}}_{k|k}$ and ${\widehat{f}}_{k|k}$, which are estimated by TSKF, and ${\widehat{f}}_{k|k}$ is the estimation of fault. By introducing the potential U-V transformation, ${\widehat{x}}_{k|k}^{\overline{a}}$ can reach the same estimation accuracy as $x_k^a$ with less computational cost. The design process of the TSKF can be divided into four parts, namely, ASKF, U-V transformation, decoupling, and TSKF, as shown in Figure 1, where ASKF presents the optimal estimation, which is also the aim of TSKF. U-V transformation is hypothetically present and is quantified, while the decoupling part establishes the conversion relationship between ${\widehat{x}}_{k|k}^a$ and ${\widehat{x}}_{k|k}^{\overline{a}}$. By repeating U-V transformation, ${\widehat{x}}_{k|k}^a$ is eliminated; moreover, the parts bounded by the dotted line in Figure 1 only contribute to the design process of the TSKF and are not involved its calculation. With the refresh of the U-V transformation matrix, the optimal state estimation ${\widehat{x}}_{k|k}^o$ can be formed by the linear combination of ${\widehat{x}}_{k|k}$ and ${\widehat{f}}_{k|k}$, and ${\widehat{f}}_{k|k}$ becomes the optimal fault estimation. More detailed TSKF design processes can be found in [14,34]. However, due to the lack of consideration of multi-source uncertainty in the design process, the TSKF-based fault estimation needs to be improved. By combining with the optimal filtering in Section 3.2, the TSHeKF is proposed in the subsequent study and realizes optimal fault estimation considering multi-source uncertainty.

4.2. TSHeKF-Based Optimal Estimation

Considering the multi-source uncertain environment, the basic filter of the TSKF needs to be replaced by HeKF proposed in Section 3.2. Meanwhile, in this section, the optimal fault estimation for one type or batch of aeroengines is proposed based on Equation (3). Same as the TSKF, the TSHeKF is designed in the same four-step process which shown in Figure 1.

• Augmented State estimator

Considering the dynamics of $f_k$ as $f_{k+1}=f_k+w_k^f$, where $w_k^f\sim N(O^{N_f\times 1},Q_k^f)$ and is independent from $w_k^m$, $v_k$, $w_k^u$, $E(w_k^m{w}_k^{f T})=Q_k^{xf}$. Treating $f_k$ as the augmented state, the optimal augmented state estimator under multi-source uncertainty is designed as follows:

(24)$\left\{\begin{array}{l}{\widehat{x}}_{k|k-1}^a={\tilde{A}}^a{\widehat{x}}_{k-1|k-1}^a+{\tilde{B}}^a{u}_{k-1}+M^a{w}_{k-1}^a\\ {\widehat{x}}_{k|k}^a={\widehat{x}}_{k|k-1}^a+K_k^a(y_k-{\tilde{C}}^a{\widehat{x}}_{k|k-1}^a-{\tilde{D}}^d{u}_k-N{v}_k/\sqrt{\Delta t})\end{array}\right.$

where

(25)$\begin{array}{l}{x}_{(.)}^a=\left[\begin{array}{c}{x}_{(.)}\\ f_{(.)}\end{array}\right],w_k^a=\left[\begin{array}{c}{w}_k^m\\ w_k^f\end{array}\right],{\tilde{A}}^a=\left[\begin{array}{cc}{\tilde{A}}^d& {\tilde{E}}^d\\ O^{N\times N_f}& I^{N_f}\end{array}\right],{\tilde{B}}^a=\left[\begin{array}{c}{\tilde{B}}^d\\ O^{N_f\times N_{in}}\end{array}\right],M^a=\left[\begin{array}{cc}M& O^{N\times N_f}\\ O^{N_f\times N}& I^{N_f}\end{array}\right],\\ {\tilde{C}}^a=\left[\begin{array}{cc}{\tilde{C}}^d& F^d\end{array}\right]\end{array}$

Obviously, Equation (24) is a typical filter as in Equation (16), and presents the optimal state estimation ${\widehat{x}}_{k|k}^a={[{\widehat{x}}_{k|k}^T,{\widehat{f}}_{k|k}^T]}^T$, which is only used as the control group during the design process of the TSHeKF.

• U-V transformations

The filtering process of the sub-filter is indicated by an adding an upper underline in this paper, and the one-step prediction and the optimal state estimation have been augmented as ${\overline{x}}_{k|k-1}^a={[{\overline{x}}_{k|k-1}^T,{\overline{f}}_{k|k-1}^T]}^T$ and ${\overline{x}}_{k|k}^a={[{\overline{x}}_{k|k}^T,{\overline{f}}_{k|k}^T]}^T$, the same as the optimal gain matrix ${\overline{K}}_k^a={[{\overline{K}}_k^{x T},{\overline{K}}_k^{f T}]}^T$. By introducing $U_k$ and $V_k$ as the gain of fault to state, ${\overline{x}}_{k|k}^a$ meets the estimation accuracy as ${\widehat{x}}_{k|k}^a$, and the U-V transformation for the TSHeKF is as follows:

(26)$\left\{\begin{array}{c}{\widehat{x}}_{k|k-1}^a=T(U_k){\overline{x}}_{k|k-1}^a\\ {\widehat{x}}_{k|k}^a=T(V_k){\overline{x}}_{k|k}^a\end{array}\right.$

where

(27)$\begin{array}{cc}T(U_k)=\left[\begin{array}{cc}{I}^N& U_k\\ O^{N_f\times N}& I^{N_f}\end{array}\right],& T(V_k)=\left[\begin{array}{cc}{I}^N& V_k\\ O^{N_f\times N}& I^{N_f}\end{array}\right]\end{array}$

Obviously, $T(U_k)+T(V_k)=T(U_k+V_k)$, and the inverse of $T(U_k)$ is $T(-U_k)$. By repeating the U-V transformation twice, the estimation of augmented state can be removed, and Equation (26) becomes:

(28)$\left\{\begin{array}{c}{\overline{x}}_{k|k-1}^a=T(-U_k)[{\tilde{A}}^{a}T(V_{k-1}){\overline{x}}_{k-1|k-1}^a+{\tilde{B}}^a{u}_{k-1}+M^a{w}_{k-1}^a]\hfill \\ {\overline{x}}_{k|k}^a=T(U_k-V_k){\overline{x}}_{k|k-1}^a+{\overline{K}}_k^a[y_k-{\tilde{C}}^{a}T(U_k){\overline{x}}_{k|k-1}^a-{\tilde{D}}^d{u}_k-N{v}_k/\sqrt{\Delta t}]\hfill \end{array}\right.$

Additionally, the statistical property of ${\overline{x}}_{k|k-1}^a$ and ${\overline{x}}_{k|k}^a$ is quantified below:

(29)$\begin{array}{c}E({\overline{x}}_{k|k-1}^a)=T(-U_k)[{\overline{A}}^{a}T(V_{k-1})E({\overline{x}}_{k-1|k-1}^a)+{\overline{B}}^a{m}_{k-1}^u]\hfill \\ \mathsf{\Psi }({\overline{x}}_{k|k-1}^a)=T(-U_k)\{\mathsf{\Psi }[{\tilde{A}}^{a}T(V_{k-1}){\overline{x}}_{k-1|k-1}^a]+{\overline{Q}}_{k-1}^a\}T{(-U_k)}^T\hfill \\ E({\overline{x}}_{k|k}^a)={\overline{K}}_k^a[y_k-{\overline{C}}^{a}T(U_k)E({\overline{x}}_{k|k-1}^a)-{\overline{D}}^d{m}_{k-1}^u]+T(U_k-V_k)E({\overline{x}}_{k|k-1}^a)\hfill \\ {\overline{K}}_k^a=T(U_k-V_k)\mathsf{\Psi }({\overline{x}}_{k|k-1}^a){[{\overline{C}}^{a}T(U_k)]}^T{({\overline{L}}_k^a)}^{-1}\hfill \\ \mathsf{\Psi }({\overline{x}}_{k|k}^a)=[T(U_k-V_k)-{\overline{K}}_k^a{\overline{C}}^{a}T(U_k)]\mathsf{\Psi }({\overline{x}}_{k|k-1}^a)T{(U_k-V_k)}^T\hfill \end{array}$

where $\ast$ represents the transpose of the corresponding part of partitioned matrix.

(30)$\begin{array}{c}\mathsf{\Psi }[{\tilde{A}}^{a}T(V_{k-1}){\overline{x}}_{k-1|k-1}^a]=\left[\begin{array}{cc}{\mathsf{\Psi }}_{AVx}^x& {\mathsf{\Psi }}_{AVx}^{xf}\\ \ast & {\mathsf{\Psi }}_{AVx}^f\end{array}\right]=\left[\begin{array}{cc}\begin{array}{l}\mathsf{\Psi }({\tilde{A}}^d{\overline{x}}_{k-1|k-1})\\ +\mathsf{\Psi }({\tilde{E}}^d{\overline{f}}_{k-1|k-1})\\ +\mathsf{\Psi }({\tilde{A}}^d{V}_{k-1}{\overline{f}}_{k-1|k-1})\end{array}& ({\overline{A}}^d{V}_{k-1}+{\overline{E}}^d)\mathsf{\Psi }({\overline{f}}_{k-1|k-1})\\ \ast & \mathsf{\Psi }({\overline{f}}_{k-1|k-1})\end{array}\right]\hfill \\ {\overline{Q}}_{k-1}^a=\mathsf{\Psi }({\tilde{B}}^a{u}_{k-1}+M^a{w}_{k-1}^a)=\left[\begin{array}{cc}{\overline{Q}}_{k-1}^x& {\overline{Q}}_{k-1}^{xf}\\ \ast & {\overline{Q}}_{k-1}^f\end{array}\right]=\left[\begin{array}{cc}\Delta tM{Q}_{k-1}^m{M}^T+\mathsf{\Psi }({\tilde{B}}^d{u}_{k-1})& Q_{k-1}^{xf}\\ \ast & Q_{k-1}^f\end{array}\right]\hfill \\ {\overline{L}}_k^a={\overline{C}}^{a}T(U_k)\mathsf{\Psi }({\widehat{x}}_{k|k-1}^a)T{(U_k)}^T{({\overline{C}}^a)}^T+S_k+{\overline{R}}_k^a\hfill \\ {\overline{R}}_k^a=\mathsf{\Psi }({\tilde{C}}^d{\overline{x}}_{k|k-1})+\mathsf{\Psi }({\tilde{C}}^d{U}_k{\overline{f}}_{k|k-1})+F\mathsf{\Psi }({\overline{f}}_{k|k-1})F^T\hfill \end{array}$

• Decoupling

By setting the correlation of state and fault in $\mathsf{\Psi }({\overline{x}}_{k|k-1}^a)$ and $\mathsf{\Psi }({\overline{x}}_{k|k}^a)$ to zero, the decoupling process is achieved and obtained:

(31)$\begin{array}{l}0={\overline{Q}}_{k-1}^{xf}-U_k{\overline{Q}}_{k-1}^f+{\overline{U}}_k\mathsf{\Psi }({\overline{f}}_{k-1|k-1})-U_k\mathsf{\Psi }({\overline{f}}_{k-1|k-1})\\ 0=U_k-V_k-{\overline{K}}_k^x{H}_k\end{array}$

where ${\overline{U}}_k={\overline{A}}^d{V}_{k-1}+{\overline{E}}^d$, $H_k={\overline{C}}^d{U}_k+F$. Moreover, the optimal sub-filter of fault has the following form:

(32)$\left\{\begin{array}{l}{\overline{f}}_{k|k-1}={\overline{f}}_{k-1|k-1}+w_{k-1}^f\\ {\overline{f}}_{k|k}={\overline{f}}_{k|k-1}+{\overline{K}}_k^f(y_k-{\tilde{C}}^d{\overline{x}}_{k|k-1}-{\tilde{D}}^d{u}_k-H_k{\overline{f}}_{k|k-1}-N{v}_k/\sqrt{\Delta t})\end{array}\right.$

Based on Equations (17), (19), (21), and (22), the statistical property of the optimal fault estimation is quantified as follows:

(33)$\begin{array}{l}E({\overline{f}}_{k|k-1})=E({\overline{f}}_{k-1|k-1})\\ \mathsf{\Psi }({\overline{f}}_{k|k-1})=\mathsf{\Psi }({\overline{f}}_{k-1|k-1})+Q_{k-1}^f\\ E({\overline{f}}_{k|k})=E({\overline{f}}_{k|k-1})+{\overline{K}}_k^f[y_k-E({\tilde{C}}^d{\overline{x}}_{k|k-1}+{\tilde{D}}^d{u}_k+H_k{\overline{f}}_{k|k-1})]\\ {\overline{K}}_k^f=\mathsf{\Psi }({\overline{f}}_{k|k-1})H_k^T{({\overline{L}}_k^f)}^{-1}\\ \mathsf{\Psi }({\overline{f}}_{k|k})=(I^{N_f}-{\overline{K}}_k^f{H}_k)\mathsf{\Psi }({\overline{f}}_{k|k-1})\end{array}$

where ${\overline{L}}_k^f=\mathsf{\Psi }({\tilde{C}}^d{\overline{x}}_{k|k-1})+H_k\mathsf{\Psi }({\overline{f}}_{k|k-1})H_k^T+S_k$. By adding the compensation of the unknown fault on the state $\Delta x$, which is a random vector, the optimal sub-filter of state is expressed in Equation (34).

(34)$\left\{\begin{array}{l}{\overline{x}}_{k|k-1}={\tilde{A}}^d{\overline{x}}_{k-1|k-1}+{\tilde{B}}^d{u}_{k-1}+M{x}_{k-1}^m+\Delta x_{k-1}\\ {\overline{x}}_{k|k}={\overline{x}}_{k|k-1}+{\overline{K}}_k^x(y_k-{\tilde{C}}^d{\overline{x}}_{k|k-1}-{\tilde{D}}^d{u}_k-N{v}_k/\sqrt{\Delta t})\end{array}\right.$

With ${\overline{L}}_k^x=E(C)\mathsf{\Psi }({\overline{x}}_{k|k-1})E(C^T)+S_k+{\overline{R}}_k^x$ and ${({\overline{R}}_k^x)}_{lm}=soe[{\mathsf{\Psi }}^{lm}(C)\ast E({\overline{\widehat{x}}}_{k|k-1}^2)]$, the statistical property of the sub-filter-based state estimation is quantified as follows:

(35)$\begin{array}{l}E({\overline{x}}_{k|k-1})={\overline{A}}^{d}E({\overline{x}}_{k|k-1})+{\overline{B}}^d{m}_{k-1}^u+E(\Delta x_{k-1})\\ \mathsf{\Psi }({\overline{x}}_{k|k-1})=\mathsf{\Psi }({\tilde{A}}^d{\overline{x}}_{k|k-1})+\mathsf{\Psi }({\tilde{B}}^d{u}_{k-1})+\mathsf{\Psi }(M{x}_{k-1}^m)+\mathsf{\Psi }(\Delta x_{k-1})\\ E({\overline{x}}_{k|k})=E({\overline{x}}_{k|k-1})+{\overline{K}}_k^x[y_k-{\overline{C}}^{d}E({\overline{x}}_{k|k-1})-{\overline{D}}^d{m}_k^u]\\ {\overline{K}}_k^x=\mathsf{\Psi }({\overline{x}}_{k|k-1}){({\overline{C}}^d)}^T{({\overline{L}}_k^x)}^{-1}\\ \mathsf{\Psi }({\overline{x}}_{k|k})=[I^N-{\overline{K}}_k^x{\overline{C}}^d]\mathsf{\Psi }({\overline{x}}_{k|k-1})\end{array}$

The sub-filters are all established, and based on Equations (29), (31), (33), and (35), the transformation matrices are quantified as:

(36)$\begin{array}{l}{U}_k={\overline{U}}_k+(Q_{k-1}^{xf}-{\overline{U}}_k{Q}_{k-1}^f)\mathsf{\Psi }{({\overline{f}}_{k|k-1})}^{-1}\\ V_k=U_k-{\overline{K}}_k^x{H}_k\end{array}$

Meanwhile, the statistical property of the random vector $\Delta x$ can be expressed based on Equations (31), (33), and (35) as follows:

(37)$\begin{array}{l}E({\Delta }_{k-1})=({\overline{U}}_k-U_k)E({\overline{f}}_{k-1|k-1})\\ \mathsf{\Psi }({\Delta }_{k-1})=({\overline{U}}_k-U_k)\mathsf{\Psi }({\overline{f}}_{k-1|k-1})U_k^T\end{array}$

• Optimal estimation under multi-source uncertainty

Finally, the design of the TSHeKF is completed based on the above steps and is optimal under multi-source uncertainty compared to standard TSKF. The optimal estimation of state and actuator fault is as follows:

(38)$\left\{\begin{array}{c}{\widehat{x}}_{k|k-1}={\overline{x}}_{k|k-1}+U_k{\overline{f}}_{k|k-1}\hfill \\ {\widehat{x}}_{k|k}={\overline{x}}_{k|k}+V_k{\overline{f}}_{k|k}\hfill \\ {\overline{x}}_{k|k-1}={\tilde{A}}^d{\overline{x}}_{k-1|k-1}+{\tilde{B}}^d{u}_{k-1}+M{x}_{k-1}^m+\Delta x_{k-1}\hfill \\ {\overline{x}}_{k|k}={\overline{x}}_{k|k-1}+{\overline{K}}_k^x(y_k-{\tilde{C}}^d{\overline{x}}_{k|k-1}-{\tilde{D}}^d{u}_k-N{v}_k/\sqrt{\Delta t})\hfill \\ {\overline{f}}_{k|k-1}={\overline{f}}_{k-1|k-1}+w_{k-1}^f\hfill \\ {\overline{f}}_{k|k}={\overline{f}}_{k|k-1}+{\overline{K}}_k^f(y_k-{\tilde{C}}^d{\overline{x}}_{k|k-1}-{\tilde{D}}^d{u}_k-H_k{\overline{f}}_{k|k-1}-N{v}_k/\sqrt{\Delta t})\hfill \end{array}\right.$

where ${\widehat{x}}_{k|k}$ is the optimal estimation of state when considering multi-source uncertainty in Equation (3), and the statistical property of ${\widehat{x}}_{k|k}$ is quantified based on Equations (33), (35), (36), and (38) as follows. Meanwhile, the optimal output prediction can also be presented based on ${\widehat{x}}_{k|k}$. Especially since the different sensors are mutually independent, $\mathsf{\Psi }({\widehat{y}}_k)$ is a diagonal matrix with variances of the outputs as the diagonal elements.

(39)$\begin{array}{l}E({\widehat{x}}_{k|k})=E({\overline{x}}_{k|k})+V_{k}E({\overline{f}}_{k|k})\\ \mathsf{\Psi }({\widehat{x}}_{k|k})=\mathsf{\Psi }({\overline{x}}_{k|k})+V_k\mathsf{\Psi }({\overline{f}}_{k|k})V_k^T\\ E({\widehat{y}}_k)={\overline{C}}^{d}E({\widehat{x}}_{k|k})+{\overline{D}}^d{m}_k^u+FE({\overline{f}}_{k|k})\\ \mathsf{\Psi }({\widehat{y}}_k)=\mathsf{\Psi }({\tilde{C}}^d{\widehat{x}}_{k|k})+\mathsf{\Psi }({\tilde{D}}^d{u}_k)+F\mathsf{\Psi }({\overline{f}}_{k|k})F^T+N{N}^T\end{array}$

4.3. Conservativeness-Reduced Fault Estimation Based on HFO

The TSHeKF-based fault estimation in the previous section takes into account the multi-source uncertainty shown in Equation (3), which means it is optimal for the corresponding type or batch of aeroengines, but still retains conservativeness for the individual. Therefore, in order to improve the accuracy of estimation and reduce the conservativeness of the fault estimation system, multiple TSHeKFs are used for estimation of $d_k^x$ and $f_k$, respectively, to achieve optimal fault estimation for individual aeroengines. The structure of the HFO is shown in Figure 2. The design process of the conservativeness-reduced fault estimation for individual aeroengines is based on Equation (6).

Where $y_k$ is the measurement of sensors, $m_k^u$ is the expectation of the control system output and can also be regarded as the expected input of the aeroengine. ${\widehat{f}}_{k|k}$ and ${\widehat{d}}_{k|k}^x$ are the TSHeKF-based unknown input estimations of fault and bias of state function, respectively. Moreover, ${\widehat{x}}_{k|k}^f$ and ${\widehat{x}}_{k|k}^d$ are both optimal state estimations of the corresponding TSHeKF, and together contribute to optimal state estimation of state under multi-source uncertainty. ${\widehat{y}}_k$ is the output prediction, ${\widehat{f}}_k$ is the estimation of actuator fault, and ${\widehat{d}}_k^x$ is the optimal quantification of the state function bias $d_k^x$, which are all random vectors with the statistical property in Equation (8). Compared with the fault estimation in Section 4.2, the HFO-based fault estimation can accommodate the effect caused by model uncertainty in state function by adding the estimation of $d_k^x$, which improves the adaptiveness and reduces the conservativeness. Since the output function bias $d_k^y$ cannot be quantified, the state optimal estimation process and output optimal prediction process are faced with multi-source uncertainty. Meanwhile, based on Equation (6), with the estimation of $d_k^x$, there is no model uncertainty with the one-step prediction. Therefore, during the design process of the HFO, the one-step prediction of each filter follows the LKF, and the estimation of state and output should be handled by the proposed TSHeKF. The data flow of the HFO is shown in Figure 3, where the gain matrices $V_0^f$ and $V_0^d$ and the statistical properties of ${\overline{f}}_{0|0}$, ${\overline{d}}_{0|0}$, ${\overline{x}}_{0|0}^f$, ${\overline{x}}_{0|0}^d$ are assumed to be known.

5. Simulation

Performance of the proposed HFO is verified by numerical simulation under the MATLAB platform in this section. Meanwhile, the TSKF- and TSHeKF-based actuator fault estimation is used as a control group.

5.1. Linear Model under Multi-Source Uncertainty

The structure and the section definition of the biaxial separated exhaust turbofan engine, which is taken as an example, is shown in Figure 4. The nonlinear and uncertainty modeling processes are shown in [3,35], and the linearization processes are clearly described in [36]. The rotational speed of low-pressure rotor $N_L$ and high-pressure rotor $N_H$ together form the state of the linear system. Meanwhile, the control variables of the system are composed of fuel flow $Wf$, variable stator vane ($VSV$), and variable bleed valve ($VBV$). In contrast, the actuator fault modes discussed in this paper correspond to the control variables. Concurrently, the output involves $N_L$ and $N_H$, the total pressure on the outlet of the low-pressure compressor $P_{25}$ and high-pressure compressor $P_3$, and the total temperature on the outlet of the high-pressure compressor $P_3$, high-pressure turbine $T_{45}$, and low-pressure turbine $T_5$. Considering the ground state (H = 0 km, Ma = 0), with the working condition of $PCNL=80$ (Percentage of $N_L$), the parameter of the above control and output variables are shown in Table 1 and Table 2. To facilitate the evaluation of the performance of actuator fault estimation, the fault drive matrices are normalized in this paper, as shown in Table 3. As for uncertainty modeling of state-space system, this paper follows [6], which assumes that the probabilistic system parameters obey Gaussian distribution, and the statistical property of control matrices of the state-space system for each element is shown in Table 4.

5.2. Optimal Filtering under Multi-Source Uncertainty

In this section, the performances of the TSKF, TSHeKF, and HFO are compared under fault-free conditions, and numerical simulation is carried out under multi-source uncertainty, as shown in Table 1, Table 2, Table 3 and Table 4. The simulation process is given as follows. The simulation duration is $60 \mathrm{s}$ with a sampling time of $20 \mathrm{ms}$, and a $10\%$ step of fuel flow is given from $13 \mathrm{s}$ to $45 \mathrm{s}$. The comparative results of the TSKF, TSHeKF, and HFO-based state estimation, with and without accounting for probabilistic system parameters, are shown in Figure 5. The HFO-based optimal output prediction and estimation of state function bias are shown in Figure 6 and Figure 7, respectively.

As can be seen from Figure 5, when the model uncertainty is not considered, all of the three filters can realize precise state estimation. However, when considering systems with probabilistic parameters, neither the TSKF nor the TSHeKF can achieve accurate tracking of states. Due to the additional estimation of state function bias in Figure 7, only the HFO enables optimal state estimation under multi-source uncertainty. As it is well known that the probability of a Gaussian random variable takes value in the range of $(u-\sigma ,u+\sigma )$ and is approximately equivalent to $68.27\%$, it is referred to in this paper as σ-percentage for simplicity of description. Where both the states and outputs follow the Gaussian distribution, and also considering the mutually independence between different outputs, the $\sigma$ for outputs at time $k$ are the corresponding diagonal elements of $V[\widehat{y}(k)]$. For the output prediction in Figure 6, the σ-percentage for the shown outputs are $66.44\%$, $67.41\%$, $70.74\%$, and $67.74\%$. It is clear that the HFO-based predictions of σ-percentage of the outputs are basically equal to the theoretical value, which further illustrates that the HFO-based state estimation is optimal. As for Figure 7, since increment of state and input is around zero while the aeroengine is at the steady-state point, state function bias is only valued when system deviates from the steady-state. It is clear that after the step of fuel flow, $d_k^x$ gradually stabilizes to another steady-state, and the HFO realized accurate tracking of unknown input $d_k^x$ lays the foundation for optimal filtering under multi-source uncertainty.

To further validate the performance of different filters under multi-source uncertainty, the ${10}^5$ times MC simulation is designed, and the result is shown in Table 5. Firstly, from the aspect of mean relative error, the HFO-based output prediction is undoubtedly the best performer, due to the increased estimation of $d_k^x$ compared to the other two methods. Additionally, the TSHeKF slightly outperforms the standard TSKF, which is attributable to the proposed optimal filtering under multi-source uncertainty in Section 3.2. By adding the PCE-based UQ of the discrete stochastic model, compared with the TSKF, the TSHeKF achieves optimality under more uncertainties and proposes a more robust filter for the individual. Aiming at the minimum relative error, which means the best prediction of filtering, the prediction accuracies of the three methods are very close. This, again, demonstrates that in Figure 5c,d, all three methods achieve high confidence filtering when ignoring the model uncertainty. Finally, by comparing the maximum estimation error longitudinally, it is obvious that the standard TSKF is the least robust due to the fewest considerations of uncertainty in the design process. Meanwhile, compared to the TSHeKF, the proposed HFO greatly reduces the conservativeness by adding the estimation of state function bias which increases the adaptability to individuals. Therefore, the HFO-based state estimation balances both multi-source uncertainty for one type or batch of aeroengines and the model uncertainty of an individual aeroengine, and provides a theoretical basis for fault optimal estimation under multi-source uncertainty.

5.3. Actuator Optimal Fault Estimation

In this section, the accuracies of fault estimation based on the TSKF, TSHeKF, and HFO are evaluated firstly by keeping the control process of numerical simulation in Section 5.2 unchanged and adding fault inputs. Numerical simulations of step fault and ramp fault of fuel flow are used for verification, and the estimation results of different method are shown in Figure 8 and Figure 9, respectively. The step fault is of $10\%$ and remains between $10 \mathrm{s}$ and $50 \mathrm{s}$; meanwhile, the ramp fault starts from $20 \mathrm{s}$ and maintains at $20 \mathrm{s}$, with a slope of $0.01/\mathrm{s}$.

As can be clearly seen from Figure 8 and Figure 9, the HFO provides optimal fault estimation when considering multi-source uncertainty. Meanwhile, errors in the TSKF or TSHeKF-based fault estimation are inevitably large due to the inability to achieve accurate tracking of states in the presence of probabilistic system parameters. To further illustrate the effectiveness of the proposed HFO under multi-source uncertainty and to verify its robustness to input noise, ${10}^4$ MC simulations are designed for each noise intensity of fuel flow, and the relative errors in the estimation of 10% step faults for $Wf$ are shown in Table 6.

The following conclusions can be drawn from the data in Table 6. Firstly, as the noise intensity of the actuator output increases, the accuracy of fault estimation based on different methods decreases accordingly. Additionally, the three methods are ranked in the descending order of mean relative errors as the HFO, TSHeKF, and TSKF, which is also the order of the methods’ fault estimation performances under multi-source uncertainty. Same as the comparison of the output prediction accuracy, since the multi-source uncertainty is considered in the design process of the TSHeKF, the TSHeKF-based estimation outperformed that of the TSKF. Moreover, on the foundation of the TSHeKF, the HFO enhances the adaptiveness for individual aeroengines, thus enabling greater reduction in conservativeness of fault estimation and realizing adaptive optimal estimation under multi-source uncertainty. Secondly, by comparing the minimum mean estimation errors of the different methods, it can be seen that the precisions of the HFO and TSHeKF are comparable and slightly less than the precision of the TSKF. This is because the situation with minimal mean estimation error can be seen as a system without probabilistic parameters, and aiming at this situation, the TSKF is certainly optimal. Due to the adaptation of the HFO and TSHeKF to multi-source uncertainty, especially in the presence of model uncertainty, the estimation error increases in the above situation. Lastly, from the aspect of maximum estimation error and standard deviation, the TSHeKF provides more robust fault estimation than the TSKF. Moreover, compared to the TSHeKF, the HFO is able to adaptively follow the state function bias for each individual and can provide more accurate fault estimation. Therefore, the TSKF-based fault estimation is optimal for individual aeroengines when the model is completely known. Since uncertainty is widespread and unavoidable, the TSHeKF-based estimation is optimal for the corresponding type or batch but robust for the individual. Finally, the HFO further reduces the conservativeness by adding an estimation process of the state function bias. Considering the output function bias is not estimable, when considering multi-source uncertainty, the HFO-based fault estimation is optimal while satisfying achievability.

Next, the effectiveness of HFO-based multi-fault estimation is demonstrated in Figure 10, Figure 11, Figure 12, Figure 13 and Figure 14. Each figure corresponds to a specific multi-actuator fault mode. The multi-fault modes are set by considering the components, type, sequence of occurrence, and duration of the fault, and following five fault mode were established for simulation verification.

• $F{M}_1$: $10\%$ step fault of $Wf$ start from $20 \mathrm{s}$, with $20\%$ step fault of $VSV$ start from $30 \mathrm{s}$; both faults last for $20 \mathrm{s}$.

• $F{M}_2$: $15\%$ step fault of $VSV$ and $20\%$ step fault of $VBV$; both last from $20 \mathrm{s}$ to $50 \mathrm{s}$.

• $F{M}_3$: a ramp fault of $Wf$ with slope of $0.005/\mathrm{s}$ from $10 \mathrm{s}$ to $50 \mathrm{s}$, and a $18\%$ step fault of $VSV$ start at $20 \mathrm{s}$ and last for $20 \mathrm{s}$.

• $F{M}_4$: ramp faults of $VSV$ and $VBV$, with slope of $0.008/\mathrm{s}$ from $10 \mathrm{s}$ to $30 \mathrm{s}$ for $VSV$ and $0.006/\mathrm{s}$ as slope starts from $20 \mathrm{s}$ and stays at $30 \mathrm{s}$ for $VBV$.

• $F{M}_5$: ramp faults of $Wf$ and $VBV$ within $20 \mathrm{s}$ to $50 \mathrm{s}$ with slopes of $0.004/\mathrm{s}$ and $0.005/\mathrm{s}$, respectively.

It is clear that the proposed HFO realizes optimal fault estimation with high accuracy whether the two faults occur simultaneously or are staggered, or whether they are step faults or ramp faults or both.

To sum up, compared with the commonly used fault estimation method, the HFO-based optimal fault estimation system has higher estimation accuracy, less conservativeness, and higher reliability under multi-source uncertainty.

6. Conclusions

This paper presents an innovative HFO for handling aeroengine actuator fault estimation under multi-source uncertainty consisting of system noise, measurement noise, input noise, and probabilistic system parameters. On the foundation of the PCE-based UQ of the discrete stochastic model, the optimal filtering under multi-source uncertainty is realized. Additionally, the TSHeKF-based fault estimation is also proposed by replacing the basic filter of the TSKF with HeKF, realizing optimal filtering for the linear stochastic system under multi-source uncertainty, and leading to the robust estimation for individuals. Since the bias in output function is relatively small and not estimable, by adding an estimate of state function bias, the HFO provides the best possible fault estimation for individuals and significantly reduces the conservativeness. Numerical simulations demonstrate the high accuracy of the HFO-based state estimation, output prediction, and unknown deviation estimation. Meanwhile, the effectiveness of the HFO in handling single and multiple faults has also been proven. Furthermore, MC simulations demonstrate that, compared with estimations based on the TSKF and TSHeKF, the proposed HFO has higher accuracy and better adaptiveness with guaranteed robustness when considering multi-source uncertainty. This enormously improves the fault diagnosis performance of the engine control system in complex environments and better ensures the safety and reliability of the aeroengine.

Author Contributions

Conceptualization, Y.W. and R.-Q.S.; methodology, R.-Q.S.; software, Y.W.; validation, Y.W., R.-Q.S. and L.-F.G.; formal analysis, Y.W.; investigation, R.-Q.S.; resources, L.-F.G.; data curation, Y.W.; writing—original draft preparation, Y.W.; writing—review and editing, R.-Q.S.; visualization, R.-Q.S.; supervision, R.-Q.S.; project administration, L.-F.G.; funding acquisition, L.-F.G. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data presented in this study are available on request from the corresponding author. The data are not publicly available due to privacy and ethical concerns.

Conflicts of Interest

The authors declare no conflicts of interest.

Abbreviation

ASKF	Augmented State Kalman Filter
HFO	Hybrid Fault Observer
HeKF	Hyperelliptic Kalman Filter
LKF	Linear Kalman filter
MC	Monte Carlo
PCE	Polynomial Chaos Expansion
TSKF	Two-stage Kalman Filter
TSHeKF	Two-stage Hyperelliptic Kalman Filter
UQ	Uncertainty Quantification
ZKF	Zonotope Kalman Filter

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. TSKF-based optimal fault estimation.

Figure 2. Structure of HFO.

Figure 3. Data flow of HFO.

Figure 4. Biaxial separated exhaust turbofan engine.

Figure 5. State estimation with (without) probabilistic system parameters.

Figure 6. Output prediction of HFO (partial).

Figure 7. HFO-based optimal estimation of [Forumla omitted. See PDF.].

Figure 8. Estimation of fuel flow step fault.

Figure 9. Estimation of fuel flow ramp fault.

Figure 10. HFO-based actuator multi-fault estimation of [Forumla omitted. See PDF.].

Figure 11. HFO-based actuator multi-fault estimation of [Forumla omitted. See PDF.].

Figure 12. HFO-based actuator multi-fault estimation of [Forumla omitted. See PDF.].

Figure 13. HFO-based actuator multi-fault estimation of [Forumla omitted. See PDF.].

Figure 14. HFO-based actuator multi-fault estimation of [Forumla omitted. See PDF.].

Table 1

Turbofan engine control variables.

Control Variable	Relative Error (%)	Nominal Value
$Wf$	2.0	0.75 $\mathrm{kg}/\mathrm{s}$
$VSV$	5.0	−13.62°
$VBV$	5.0	0.52°

Table 2

Turbofan engine output variables.

Output Variable	Relative Error (%)	Nominal Value
$N_L$	0.5	3190 $\mathrm{RPM}$
$N_H$	0.5	$\mathrm{10,699} \mathrm{RPM}$
$P_{25}$	0.5	$181.17 \mathrm{kPa}$
$P_3$	0.5	$1697.41 \mathrm{kPa}$
$T_3$	2.0	$696.64 \mathrm{K}$
$T_{45}$	2.0	$875.02 \mathrm{K}$
$T_5$	2.0	$666.65 \mathrm{K}$

Table 3

State-space model of turbofan engine normalized fault drive matrix.

Drive Matrix	Value
$E$	$\left|\begin{array}{ccc}830.45& -183.72& -61.49\\ 2315.91& 262.64& 28.66\end{array}\right|$
$F$	$\left|\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 1.35& 9.40& -11.07\\ 297.12& -73.92& -48.20\\ 24.37& -5.44& -2.23\\ 250.96& 10.74& 7.16\\ 168.41& 11.64& 7.56\end{array}\right|$

Table 4

State-space model of turbofan engine.

Control Matrix	Expectation	Standard Deviation
$A$	$\left|\begin{array}{cc}-2.22& 0.68\\ -0.35& -1.44\end{array}\right|$	$\left|\begin{array}{cc}11.09& 3.41\\ 1.75& 7.21\end{array}\right|$(${10}^{-2}$)
$B$	$\left|\begin{array}{ccc}1100.82& 13.48& -118.44\\ 3069.93& -19.27& 55.20\end{array}\right|$	$\left|\begin{array}{ccc}33.02& 0.67& 5.92\\ 92.10& 0.96& 2.76\end{array}\right|$
$C$	$\left|\begin{array}{cc}1& 0\\ 0& 1\\ 0.096& -0.033\\ 0.368& 0.270\\ 0.022& 0.049\\ -0.040& -0.015\\ -0.067& -0.025\end{array}\right|$	$\left|\begin{array}{cc}0& 0\\ 0& 0\\ 0.48& 0.16\\ 1.84& 1.35\\ 0.11& 0.24\\ 0.20& 0.07\\ 0.34& 0.12\end{array}\right|$(${10}^{-2}$)
$D$	$\left|\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 1.79& -0.69& -21.33\\ 393.86& 5.42& -92.83\\ 32.30& 0.40& -4.30\\ 332.67& -0.799& 13.79\\ 223.24& -0.85& 14.56\end{array}\right|$	$\left|\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0.04& 0.10& 1.71\\ 7.88& 0.81& 7.43\\ 0.65& 0.06& 0.34\\ 6.65& 0.12& 1.10\\ 4.46& 0.13& 1.16\end{array}\right|$

Table 5

MC simulations for output prediction.

Method	Index	Relative Error of Output Prediction (%)
		${\mathit{N}}_{\mathit{L}}$	${\mathit{N}}_{\mathit{H}}$	${\mathit{P}}_{25}$	${\mathit{P}}_3$	${\mathit{T}}_3$	${\mathit{T}}_{45}$	${\mathit{T}}_5$
TSKF	Min	0.382	0.385	0.389	0.389	1.513	1.530	1.519
	Mean	0.507	0.408	0.412	0.412	1.585	1.619	1.627
	Max	2.679	0.945	0.562	0.580	1.656	3.150	3.012
TSHeKF	Min	0.376	0.377	0.381	0.376	1.512	1.531	1.520
	Mean	0.493	0.399	0.404	0.398	1.584	1.609	1.617
	Max	2.511	0.892	0.527	0.511	1.655	3.084	2.989
HFO	Min	0.367	0.375	0.359	0.370	1.512	1.530	1.516
	Mean	0.386	0.395	0.383	0.391	1.584	1.591	1.598
	Max	0.404	0.419	0.407	0.412	1.656	1.668	1.689

Table 6

MC simulations of fault estimation under different $Wf$ noise.

Method	Relative Error ofFault Estimation		Input Noise Intensity
			2%	4%	6%	8%	10%
TSKF	Mean	Min	0.056	0.066	0.074	0.081	0.089
		Mean	0.183	0.224	0.278	0.329	0.373
		Max	1.518	1.927	2.485	3.025	3.471
	Standard deviation		0.160	0.196	0.245	0.290	0.326
TSHeKF	Mean	Min	0.055	0.070	0.081	0.094	0.108
		Mean	0.178	0.191	0.211	0.231	0.248
		Max	1.624	1.703	1.860	2.010	2.127
	Standard deviation		0.157	0.159	0.168	0.177	0.184
HFO	Mean	Min	0.054	0.060	0.068	0.077	0.092
		Mean	0.086	0.102	0.121	0.142	0.166
		Max	0.218	0.2288	0.2498	0.287	0.350
	Standard deviation		0.023	0.022	0.023	0.027	0.032

References

1. Lv, C.; Chang, J.; Bao, W.; Yu, D. Recent research progress on airbreathing aero-engine control algorithm. Propuls. Power Res.; 2022; 11, pp. 1-57. [DOI: https://dx.doi.org/10.1016/j.jppr.2022.02.003]

2. Rath, N.; Mishra, R.K.; Kushari, A. Aero engine health monitoring, diagnostics and prognostics for condition-based mainte-nance: An overview. Int. J. Turbo Jet Eng.; 2024; 40, pp. s279-s292. [DOI: https://dx.doi.org/10.1515/tjj-2022-0020]

3. Gou, L.; Sun, R.; Han, X. FDIA System for Sensors of the Aero-Engine Control System Based on the Immune Fusion Kalman Filter. Math. Probl. Eng.; 2021; 2021, 6662425. [DOI: https://dx.doi.org/10.1155/2021/6662425]

4. Roy, C.J.; Oberkampf, W.L. A comprehensive framework for verification, validation, and uncertainty quantification in sci-entific computing. Comput. Methods Appl. Mech. Eng.; 2011; 200, pp. 2131-2144. [DOI: https://dx.doi.org/10.1016/j.cma.2011.03.016]

5. Oberkampf, W.L.; Helton, J.C.; Joslyn, C.A.; Wojtkiewicz, S.F.; Ferson, S. Challenge problems: Uncertainty in system response given uncertain parameters. Reliab. Eng. Syst. Saf.; 2004; 85, pp. 11-19. [DOI: https://dx.doi.org/10.1016/j.ress.2004.03.002]

6. Sun, R.Q.; Han, X.B.; Chen, Y.X.; Gou, L.F. Hyperelliptic Kalman filter-based aeroengine sensor fault FDIA system under multi-source uncertainty. Aerosp. Sci. Technol.; 2023; 132, 108058. [DOI: https://dx.doi.org/10.1016/j.ast.2022.108058]

7. Liu, Z.; Wang, Z.; Wang, Y.; Ji, Z. Sensor fault estimation based on the constrained zonotopic Kalman filter. Int. J. Robust Nonlinear Control; 2021; 31, pp. 5984-6006. [DOI: https://dx.doi.org/10.1002/rnc.5629]

8. Gou, L.; Shen, Y.; Zheng, H.; Zeng, X. Multi-Fault Diagnosis of an Aero-Engine Control System Using Joint Sliding Mode Observers. IEEE Access; 2020; 8, pp. 10186-10197. [DOI: https://dx.doi.org/10.1109/ACCESS.2020.2964572]

9. Hamayun, M.T.; Edwards, C.; Alwi, H. Design and Analysis of an Integral Sliding Mode Fault-Tolerant Control Scheme. IEEE Trans. Autom. Control; 2012; 57, pp. 1783-1789. [DOI: https://dx.doi.org/10.1109/TAC.2011.2180090]

10. Pourbabaee, B.; Meskin, N.; Khorasani, K. Sensor Fault Detection, Isolation, and Identification Using Multi-ple-Model-Based Hybrid Kalman Filter for Gas Turbine Engines. IEEE Trans. Control Syst. Technol.; 2016; 24, pp. 1184-1200. [DOI: https://dx.doi.org/10.1109/TCST.2015.2480003]

11. Simon, D. A comparison of filtering approaches for aircraft engine health estimation. Aerosp. Sci. Technol.; 2008; 12, pp. 276-284. [DOI: https://dx.doi.org/10.1016/j.ast.2007.06.002]

12. Jin, P.; Lu, F.; Huang, J.; Kong, X.; Fan, M. Life cycle gas path performance monitoring with control loop parameters uncertainty for aeroengine. Aerosp. Sci. Technol.; 2021; 115, 106775. [DOI: https://dx.doi.org/10.1016/j.ast.2021.106775]

13. Lu, F.; Jiang, C.; Huang, J.; Qiu, X. A multi-rate sensor fusion approach using information filters for estimating aero-engine performance degradation. Chin. J. Aeronaut.; 2019; 32, pp. 1603-1617. [DOI: https://dx.doi.org/10.1016/j.cja.2019.04.024]

14. Hsieh, C.S.; Chen, F.C. Optimal solution of the two-stage Kalman estimator. IEEE Trans. Autom. Control; 1999; 44, pp. 194-199. [DOI: https://dx.doi.org/10.1109/9.739135]

15. Hajiyev, C.; Cilden-Guler, D.; Hacizade, U. Two-Stage Kalman Filter for Fault Tolerant Estimation of Wind Speed and UAV Flight Parameters. Meas. Sci. Rev.; 2020; 20, pp. 35-42. [DOI: https://dx.doi.org/10.2478/msr-2020-0005]

16. Chen, X.; Sun, R.; Wang, F.; Song, D.; Jiang, W. Two-stage unscented Kalman filter algorithm for fault estimation in space-craft attitude control system. IET Control Theory A; 2018; 12, pp. 1781-1791. [DOI: https://dx.doi.org/10.1049/iet-cta.2017.1369]

17. Patton, R.J.; Chen, J. On eigenstructure assignment for robust fault diagnosis. Int. J. Robust Nonlinear Control; 2000; 10, pp. 1193-1208. [DOI: https://dx.doi.org/10.1002/1099-1239(20001215)10:14<1193::AID-RNC523>3.0.CO;2-R]

18. Zhou, D.; Ma, S.; Chen, Y.; Wei, T.; Zhang, H.; Wei, F. A gas path fault diagnostic model for gas turbine based on deep belief network with prior information. Proceedings of the 2018 IEEE International Conference on Prognostics and Health Management (ICPHM); Seattle, WA, USA, 11–13 June 2018; pp. 1-6.

19. Lu, F.; Huang, J.; Lv, Y. Gas Path Health Monitoring for a Turbofan Engine Based on a Nonlinear Filtering Approach. Energies; 2013; 6, pp. 492-513. [DOI: https://dx.doi.org/10.3390/en6010492]

20. Combastel, C. Zonotopes and Kalman observers: Gain optimality under distinct uncertainty paradigms and robust convergence. Automatica; 2015; 55, pp. 265-273. [DOI: https://dx.doi.org/10.1016/j.automatica.2015.03.008]

21. Zhang, W.; Wang, Z.; Guo, S.; Shen, Y. Interval estimation of sensor fault based on zonotopic Kalman filter. Int. J. Control; 2021; 94, pp. 1641-1650. [DOI: https://dx.doi.org/10.1080/00207179.2019.1662489]

22. Xiu, D. Numerical Methods for Stochastic Computations: A Spectral Method Approach; Princeton University Press: Princeton, NJ, USA, 2010; 125.

23. Dutta, P.; Bhattacharya, R. Nonlinear estimation with polynomial chaos and higher order moment updates. Proceedings of the 2010 American Control Conference; Baltimore, MD, USA, 30 June–2 July 2010; pp. 3142-3147.

24. Dutta, P.; Bhattacharya, R. Nonlinear Estimation of Hypersonic Flight Using Polynomial Chaos. Proceedings of the AIAA Guidance, Navigation, and Control Conference; Toronto, ON, Canada, 2–5 August 2010.

25. Sun, R.-Q.; Gou, L.-F.; Liu, Z.-Y.; Han, X.-B. Three-stage hyperelliptic Kalman filter for health and performance monitoring of aeroengine under multi-source uncertainty. Int. J. Engine Res.; 2023; 25, pp. 557-572.

26. Ding, S.X. Model-Based Fault Diagnosis Techniques: Design Schemes, Algorithms and Tools; Springer: New York, NY, USA, 2013.

27. Øksendal, B. Stochastic Differential Equations: An Introduction with Applications; Spring: New York, NY, USA, 2003.

28. Xia, B.; Wu, H.; Qiu, Y.; Lou, B.; Song, Y. A Galerkin Method-Based Polynomial Approximation for Parametric Problems in Power System Transient Analysis. IEEE Trans. Power Syst.; 2019; 34, pp. 1620-1629. [DOI: https://dx.doi.org/10.1109/TPWRS.2018.2879367]

29. Manfredi, P.; De Zutter, D.; Ginste, D.V. On the relationship between the stochastic Galerkin method and the pseudo-spectral collocation method for linear differential algebraic equations. J. Eng. Math.; 2018; 108, pp. 73-90. [DOI: https://dx.doi.org/10.1007/s10665-017-9909-7]

30. Oladyshkin, S.; Nowak, W. Data-driven uncertainty quantification using the arbitrary polynomial chaos expansion. Reliab. Eng. Syst. Saf.; 2012; 106, pp. 179-190. [DOI: https://dx.doi.org/10.1016/j.ress.2012.05.002]

31. Shin, Y.; Xiu, D. Nonadaptive quasi-optimal points selection for least squares linear regression. SIAM. J. Sci. Comput.; 2016; 38, pp. A385-A411. [DOI: https://dx.doi.org/10.1137/15M1015868]

32. Constantine, P.G.; Gleich, D.F.; Iaccarino, G. Spectral methods for parameterized matrix equations. SIAM J. Matrix Anal. Appl.; 2010; 31, pp. 2681-2699. [DOI: https://dx.doi.org/10.1137/090755965]

33. Oladyshkin, S.; Nowak, W. Incomplete statistical information limits the utility of high-order polynomial chaos expansions. Reliab. Eng. Syst. Saf.; 2018; 169, pp. 137-148. [DOI: https://dx.doi.org/10.1016/j.ress.2017.08.010]

34. Hsieh, C.S. Robust two-stage Kalman filters for systems with unknown inputs. IEEE Trans. Autom. Control; 2000; 45, pp. 2374-2378. [DOI: https://dx.doi.org/10.1109/9.895577]

35. Montomoli, F.; Carnevale, M.; D’Ammaro, A.; Massini, M.; Salvadori, S. Uncertainty Quantification in Computational Fluid Dynamics and Aircraft Engines; Springer: Cham, Switzerland, 2019.

36. Gou, L.; Zeng, X.; Wang, Z.; Han, G.; Lin, C.; Cheng, X. A Linearization Model of Turbofan Engine for Intelligent Analysis Towards Industrial Internet of Things. IEEE Access; 2019; 7, pp. 145313-145323. [DOI: https://dx.doi.org/10.1109/ACCESS.2019.2945337]