A Reduced Order Model for Monitoring Aeroengines Condition in Real-Time

1. Introduction

The health condition of a mechanical system determines its performance. Thus, monitoring such condition is essential for various tasks, including control, maintenance, and safety. In particular, predictive maintenance (also called condition-based maintenance) requires monitoring the system condition for diagnosis [1,2,3,4], to anticipate when it should be repaired/restored before, perhaps dangerous, breakdown occurs. Thus, predictive maintenance (i) reduces cost [5] because maintenance tasks are carried out when and where they are truly needed, not just by schedule; and (ii) increases safety by avoiding dangerous events during operation. This type of maintenance has received an increasing attention in the literature during the last decades [6,7,8], due to its interest in a variety of industrial sectors, such as the aerospace [9,10], manufacturing [11], and railway [12] industries.

Although the methods developed in the paper apply to more general systems, they will be illustrated (and their performance tested) considering a turbofan aeroengine for commercial aviation. This system is representative, in terms of both difficulties and opportunities, of many other complex mechanical systems.

Concentrating on aeroengines, during the last decade, an increasing interest has arisen in the main industrial companies in the field, on developing health monitoring tools to improve their predictive maintenance capability. Many of these tools are either data-driven or model-based [13,14]; besides, there are hybrid methods, such as the physics informed neural networks [15].

Purely data-driven methods [16] rely on just data and are usually based on machine learning and statistical algorithms. Their online operation is fairly inexpensive, but require a very computationally expensive offline training, which is performed using former experience in the type of engine that is being analyzed. They lack flexibility, since they do not give good results for other engines (or modified engines) not accounted for in the training stage.

Model-based methods are more flexible and precise because they explicitly take the aeroengine physics into account, using a nonlinear engine model. Such model computes the sensors outputs in terms of several degradation parameters that depend on the engine condition (health status) and the operating regime. The operating regime is determined by the flight altitude and Mach number, and the turbine inlet temperature. The latter is very important since it determines the engine thrust (and affects the engine lifetime) in aircraft aeroengines, and the engine power in related systems, such as naval gas turbines and aeroderivatives. Concerning the aeroengine degradation parameters, Stamatis et al. [17] proposed flow capacities and adiabatic efficiencies of several engine components, which are generally acepted nowadays. The main turbofan manufactures have developed their own models for the engines they produce. Unfortunately, these engine models are confidential and not available in the literature. In this paper, we shall use an engine model based on PROOSIS [18], which is able to simulate different aeroengine types. It must be noted that this (and any other) engine model, uses discrete empirical data obtained by actual measurements in the main engine components, and is based on an iterative method. The latter exhibits convergence difficulties when trying to simulate cases that are outside the operating regime where the empirical data have been acquired. A review of model-based methods can be found in our recent work [19], where a methodology was developed for aeroengine diagnosis. The present paper is a step further in the sense that, while diagnosis was performed in [19] using a PROOSIS-based engine model, here we shall construct and use a convenient reduced order model (ROM), combined with gradient-like condition monitoring tools. This means that the methodology to be developed in this paper shares the advantages of the data-driven and model-based methods, namely fast online operation and precise, robust results. In fact, we shall use the experience gained in [19] to directly select the most convenient gradient-like method, which is an appropriately adapted global, constraint Newton-based method. The offline construction of the ROM requires: (i) running a full engine model a very large number of times, to obtain data that are collected in a (high-order) tensor; this may be quite computationally expensive. And (ii) applying a tensor decomposition, which instead is fairly inexpensive. The online operation of the ROM is quite fast, permitting real-time results. The word real-time is somewhat ambiguous since it has several possible meanings, depending of the context. To be precise, in the present context we define two relevant meanings that will be used along the paper: (a) monitoring data in the aircraft cockpit requires a very small computational time, not larger than, say, 0.5 CPU seconds, which is enough for human eye detection and will be called continuous real-time; on the other hand, (b) since the flight time is, at least, 30 minutes, a computational time of up to a few CPU minutes will be called in this paper in-flight real-time.

Specifically, the ROM will be constructed by applying a tensor extension of standard singular value decomposition (SVD), called higher order singular value decomposition (HOSVD), combined with one-dimensional interpolation, to a set of data organized in tensor form. This combination is a higher dimensional extension of the celebrated, very useful method [20,21], which combines proper orthogonal decomposition (or standard SVD) and interpolation. The required data will be obtained using an engine model. In other words, a surrogate of the engine model will be constructed that gives the sensors outcomes in terms of the engine operating regime and the degradation parameters. Similar ROMs based on HOSVD plus interpolation have been already used in several fields, including generation of aerodynamic databases [22], continuous real-time control of reciprocating engines [23], and aircraft conceptual design [24,25]. HOSVD-based methods also permit constructing tools for recovering lost data [26] and filtering large errors [27] in multidimensional databases, which are relevant in condition monitoring tasks, including:

• Recovering sensors data not provided by the engine model, due to convergence difficulties. Also, recovering data at places in the engine not accessible to sensor measurements.
• Correcting wrong sensors data due to ill-functioning of the sensors themselves or the engine instrumentation, which may produce large errors.

However, these applications are beyond the scope of this paper and will be considered elsewhere. Here, we shall only consider small errors in the sensors data that somewhat mimic experimental noise, see below.

Using the surrogate engine model, two condition monitoring tools will be developed that compute both the engine degradations and the turbine inlet temperature, $T_{4t}$, at which the sensors data have been acquired:

• The condition monitoring tool 1 is based on a linear approximation and thus requires small degradations around a baseline state, which are expected in each individual flight. Its online operation is extremely fast, namely it can operate in continuous real-time. It precisely computes the turbine inlet temperature $T_{4t}$ and gives a first approximation of the engine degradation. This tool has been obtained as an unexpected byproduct, when checking the importance of nonlinear effects.
• The condition monitoring tool 2 is the counterpart, using the ROM surrogate, of a very efficient tool developed in [19] for the full engine model. This tool operates in in-flight real-time and is designed for non-small degradations (up to $\sim 2$%), for which nonlinear effects must be accounted for. Let us mention here that degradations of 2% are fairly high and usually require performing maintenance tasks on the degraded engine component [28,29,30] to avoid future dangerous events. On the other hand, the very small degradations that are expected during each flight accumulate in subsequent flights. In this sense, if an accumulated degradation becomes larger than 2 %, then it must be monitored since it could reach dangerous values as time proceeds. Such monitoring can be performed using the method developed in [19], which accurately performs diagnosis for very large degradations. In this sense, the tool 2 developed in this paper is complementary to the methodology presented in [19].

The raw (dimensional) sensors data account for properties such as temperatures, pressures, rotational speeds, and fuel consumption, which take very disparate values. This could produce ill-functioning of condition monitoring tools. Thus, the sensors outcomes need to be appropriately nondimensionalized/scaled.

On the other hand, in practice, the outcomes of actual sensors mounted in the engine are noisy. In order to take this into account, small random noise with zero mean will be added to the engine model outcomes, concluding that the developed tools are robust in connection with noise.

Computations will be performed using standard (uncompiled) MATLAB in a standard PC, with a microprocessor Intel Core I7-6500U at 2.5 GHz, with 16 Gb RAM memory. The CPU times reported along the paper for reference could be significantly reduced using more advanced, customized software and hardware. Quantitative results will be given in the form of tables, which is convenient to facilitate the reader to reproduce results.

Against the above background, the remainder of the paper is organized as follows. The specific aeroengine that will be used to test the performance of the methodology is presented in Section 2, where the convenient degradation parameters and sensors are described and the output sensors data are appropriately scaled. The required data processing methods, standard SVD and HOSVD, are briefly described in Section 3, where the detailed construction of the surrogate engine model and the aforementioned condition monitoring tools 1 and 2 are also addressed. Specific condition monitoring results are given in Section 4 and some concluding remarks, in Section 5.

2. Test case to evaluate the performance of the methodology: the CFM56 aeroengine.

As in [19], we consider a CFM56 aeroengine, which is representative of those widely used in commercial aviation nowadays [31]. It is the two-spool turbofan engine sketched in Figure 1.

In this figure, the main engine components are indicated, as well as the engine locations where several relevant pressures P, temperatures T, rotational speeds, $N_L$ and $N_H$, and the fuel consumption $W_f$ can be measured. As can be seen, the engine contains both a low-pressure spool and a high-pressure spool, which rotate concentrically at different speeds. The high-pressure compressor (HPC), the high-pressure turbine (HPT), and the combustor constitute the high-pressure part of the engine, also known as the core of the engine. The fan, the low-pressure compressor (LPC), also known as the booster, and the low-pressure turbine (LPT) are mechanically connected and constitute the low-pressure part of the engine.

For this engine, condition monitoring will be performed in the following ranges of the flight altitude h, the Mach number M, and the turbine inlet temperature $T_{4t}$

$$31,000\mathrm{ft}\le h\le 39,000\mathrm{ft},0.8\le M\le 0.86,1,350\mathrm{K}\le T_{4t}\le 1,550\mathrm{K},$$

(1)

which cover typical cruise conditions in commercial aviation.

The engine is characterized by the ten degradation parameters, whose deviation from their baseline values measure the engine degradation. They are measured in %, which represents an automatic scaling, and are nonnegative because the engine is degraded, not upgraded. In the present paper, as anticipated, they are allowed to vary in the range

$$0\le x_j\le 2,\mathrm{for}j=1\dots ,10.$$

(2)

The ten degradation parameters give the efficiencies and adiabatic capacities of various engine stations, as displayed in Table 1.

Our condition monitoring tools will use sensors data to compute the ten degradations and the turbine inlet temperature, which amount to eleven unknowns. After the analysis in [19], we shall consider eleven sensors that are conveniently uncorrelated and described in Table 2, second column.

The average values of the sensors and the associated standard deviations, displayed in the third and fourth columns of the table, have been computed from the outcomes using the full engine model at the 45 combinations of the following three values of the altitude h, three values Mach number M, and five values of the turbine inlet temperature $T_{4t}$

$$\begin{array}{cc}& h=31,000,35,000,\mathrm{and}39,000\mathrm{ft},\hfill \end{array}$$

(3)

$$\begin{array}{cc}& M=0.8,0.83,\mathrm{and}0.86,\hfill \end{array}$$

(4)

$$\begin{array}{cc}& T_{4t}=1,350,1,400,1,450,1,500,\mathrm{and}1,550.\hfill \end{array}$$

(5)

Note that these 45 combinations cover well the range defined in (1) in which the engine it to be inspected, and are necessarily computed because they will be precisely the discrete values of these parameters to construct the ROM. As can be seen in Table 2, the dimensional values of the sensors outcomes differ among each other by several orders of magnitude and, also, the standard deviations are dissimilar, since there is a factor of $\sim 7$ between the smallest and the largest. Thus, it is convenient to scale them in such a way that the scaled values are comparable and vary in comparable ranges. As thoroughly explained in [19], a convenient scaling meeting these requirements is performed as follows. The $N_{sensors}$ (=11 in the present case) scaled sensor outcomes, denoted as $y^{scal.}\left(1\right),\dots ,y^{scal.}\left(N_{sensors}\right)$, are calculated as

$$y_m^{scal.}=\frac{y_m^{dimens.}-y_m^{aver.}}{sd\left(m\right)·y_m^{aver.}},$$

(6)

where, for $m=1,\dots ,N_{sensors}$, $y_m^{dimens.}$ are the original dimensional sensor outcomes, while $y_m^{aver.}$ and $sd\left(m\right)$ are the associated averages and standard deviations displayed in Table 2. Sensors will always be scaled using (6) along the paper. This scaling can be computed for any aeroengine, with any value of $N_{sensors}$ and, mutatis mutandis, to the sensor outcomes in other mechanical systems.

3. Methods and tools

After briefly describing the well known SVD and HOSVD data processing tools, the construction of the surrogate model and the development of the already mentioned health monitoring tools 1 and 2 will be addressed.

3.1. Standard SVD and HOSVD

Standard SVD was introduced [33] by Beltrami and Jordan (in 1873-74) for square matrices, and further developed by Eckart and Young [34] for general, rectangular matrices. There are many versions of SVD [35]. In the present version, for a real, generally rectangular, $I\times J$ matrix $\mathbf{A}$, its SVD is given by

$$\mathbf{A}=\mathbf{U}\mathbf{S}{\mathbf{V}}^{\top }\mathrm{or}{A}_{ij}=\sum _{q=1}^Q{s}_q{U}_{iq}{V}_{jq},$$

(7)

where the superscript ${}^{\top }$ denotes the transpose and $\mathbf{U}$ and $\mathbf{V}$ are $I\times Q$ and $Q\times J$ matrices, respectively, with $Q=min\{I,J\}$; the columns of these matrices are mutually orthonormal with the usual Euclidean inner product. The matrix $\mathbf{S}$ is diagonal and its elements, known as the singular values and denoted as $s_n$, are real, non-negative, and sorted in decreasing order, namely $s_1\ge s_2\ge \dots \ge s_Q\ge 0$. The decomposition (7) can be calculated using the MATLAB command ‘svd’, option ‘economy’. SVD is also useful to define the condition number of a matrix $\mathbf{A}$, which is directly computed with the MATLAB command ‘cond’. It indicates the accuracy when solving linear systems whose coefficient matrix is $\mathbf{A}$; such accuracy worsens as the condition number increases [36].

The extension of SVD to higher than two order tensors (i.e., multi-dimensional databases) is highly non-trivial. There are many such extensions [37]. The most obvious one, known as the canonical decomposition, consists in re-organizing the tensor as a linear combination of rank-one tensors. The minimum number of rank-one tensors such that this decomposition is possible defines the rank of the tensor, whose computation for general tensors is an open problem still nowadays [38]. Because of this difficulty, other extensions have been developed. A very robust extension was introduced by Tucker [39] and more recently popularized by de Lathauwer et al. [40,41]. This extension, called higher order singular value decomposition (HOSVD), is also known as the Tucker decomposition. For a N-th order tensor $\mathbf{T}$ of size $I_1\times I_2\times \dots \times I_N$, its HOSVD rewrites the tensor components as

$$T_{i_1{i}_2\dots i_N}=\sum _{i_1^{\prime }=1}^{I_1}\sum _{i_2^{\prime }=1}^{I_2}\cdots \sum _{i_N^{\prime }=1}^{I_N}{S}_{i_1^{\prime }{i}_2^{\prime }\dots i_N^{\prime }}{U}_{i_1{i}_1^{\prime }}^1{U}_{i_2{i}_2^{\prime }}^2\cdots U_{i_N{i}_N^{\prime }}^N,$$

(8)

for $i_1=1,2,\dots ,I_1$, $i_2=1,2,\dots ,I_2$, …, $i_n=1,2,\dots ,I_N$. $S_{i_1^{\prime }{i}_2^{\prime }\dots i_N^{\prime }}$ are the components of a N-th order tensor $\mathbf{S}$, of size $I_1\times I_2\times \dots \times I_N$, and the square matrices ${\mathbf{U}}^1,{\mathbf{U}}^2,\dots ,{\mathbf{U}}^N$, of sizes $I_1,I_2,\dots I_N$, whose elements are $U_{i_1{n}_1}^1$, $U_{i_2{n}_2}^2\dots$, $U_{i_N{n}_N}^N$, are known as the mode matrices along the various dimensions of the tensor $\mathbf{T}$. Now, we note that the expansion appearing in (8) is a multilinear function of the components of the tensor $\mathbf{S}$ and the mode matrices, which can be considered as a tensor product of these. Thus, Equation (8), is rewritten as

$$\mathbf{T}=\mathbf{tprod}(\mathbf{S},{\mathbf{U}}^1,{\mathbf{U}}^2,\dots ,{\mathbf{U}}^N).$$

(9)

A MATLAB function to compute this tensor product can be found in [42].

Let us now explain briefly how HOSVD is constructed. To begin with, for each $n=1\dots ,N$, the tensor $\mathbf{T}$ is unfolded (which can be done using a MATLAB function given in [42]) into the matrix ${\mathbf{A}}^n$, whose columns are the fibers of $\mathbf{T}$ along the n-th tensor dimension. Then, SVD is applied to the matrix ${\mathbf{A}}^n$, which gives

$${\mathbf{A}}^n={\mathbf{U}}^n{\mathbf{S}}^n{\left({\mathbf{V}}^n\right)}^{\top }.$$

(10)

Here, the matrix ${\mathbf{U}}^n$ is precisely the mode matrix along the n-th dimension appearing in (9). Then, the tensor $\mathbf{S}$ is computed as

$$\mathbf{S}=\mathbf{tprod}(\mathbf{T},{\left({\mathbf{U}}^1\right)}^{\top },{\left({\mathbf{U}}^2\right)}^{\top },\dots ,{\left({\mathbf{U}}^N\right)}^{\top }).$$

(11)

With this, all ingredients appearing in the right hand side of (8) have been calculated.

3.2. Construction of the HOSVD-based surrogate engine model

To begin with, we construct the fourteenth order tensor

$$T_{i{j}_1{j}_2{j}_3{k}_1\dots k_{10}},$$

(12)

containing the (scaled as explained in Section 2) sensors data, obtained using the engine model, for the 45 combinations of the parameters defining the flight regime displayed in (3)-(5). In this tensor:

• For $i=1,\dots ,11$, the first index indicates the eleven sensors.
• For $j_1=1,2,3$, the second index corresponds to the three flight altitudes displayed in (3).
• For $j_2=1,2,3$, the third index is associated with the three Mach numbers appearing in (4).
• For $j_3=1,\dots ,5$, the fourth index indicates the five values of the turbine inlet temperature displayed in (5).
• The remaining ten indexes, $k_1,\dots ,k_{10}$, correspond to the ten degradations denoted as x and are allowed to take the following three values

  $$x=0,1,\mathrm{and}2\mathrm{quantified}\mathrm{in}\%.$$

  (13)

Applying HOSVD to the tensor (12), we obtain

$$\begin{array}{cc}& T_{i{j}_1{j}_2{j}_3{k}_1\dots k_{10}}=\sum _{j_1^{\prime },j_2^{\prime }=1}^3\sum _{j_3^{\prime }=1}^5\sum _{k_1^{\prime },\dots ,k_{10}^{\prime }=1}^3{S}_{i{j}_1^{\prime }{j}_2^{\prime }{j}_3^{\prime }{k}_1^{\prime }\dots k_{10}^{\prime }}{V}_{j_1{j}_1^{\prime }}^1{V}_{j_2{j}_2^{\prime }}^2{V}_{j_3{j}_3^{\prime }}^3{W}_{k_1{k}_1^{\prime }}^1\cdots W_{k_{10}{k}_{10}^{\prime }}^{10}.\hfill \end{array}$$

(14)

Now, we construct the surrogate engine model, which gives as outputs the 11 scaled sensor data, denoted as $y_1,\dots ,y_{11}$, in terms of the flight altitude, h, the Mach number, M, the turbine inlet temperature, $T_{4t}$, and the ten degradations, $x_1,\dots ,x_{10}$, for continuous values of these 13 input parameters. Thus, we consider the functions

$$y_i=F_i(h,M,T_{4t},x_1,\dots ,x_{10})\mathrm{for}i=1,\dots ,11.$$

(15)

These functions will apply for the input parameters in the ranges defined in (1)-(2), and are constructed combining the various ingredients in the right hand side of the tensor decomposition (14) and standard one-dimensional spline interpolation (performed with the MATLAB function ‘spline’), as

$$F_i(h,M,T_{4t},x_1,\dots ,x_{10})=\sum _{j_1^{\prime },j_2^{\prime }=1}^3\sum _{j_3^{\prime }=1}^5\sum _{k_1^{\prime },\dots ,k_{10}^{\prime }=1}^3{S}_{i{j}_1^{\prime }{j}_2^{\prime }{j}_3^{\prime }{k}_1\dots k_{10}^{\prime }}{v}_{j_1^{\prime }}^1\left(h\right)v_{j_2^{\prime }}^2\left(M\right)v_{j_3^{\prime }}^3\left(T_{4t}\right)w_{k_1^{\prime }}^1\left(x_1\right)\cdots w_{k_{10}^{\prime }}^{10}\left(x_{10}\right),$$

(16)

where the tensor components $S_{i{j}_1^{\prime }{j}_2^{\prime }{j}_3^{\prime }{k}_1\dots k_{10}^{\prime }}$ are as computed in the tensor decomposition (14) and:

• The functions $v_{j_1^{\prime }}^1$, $v_{j_2^{\prime }}^2$, and $v_{j_3^{\prime }}^3$ are constructed by interpolation in the three intervals defined in (1), based on the intermediate points defined in (3)-(5), where the functions take the values $V_{j_1{j}_1^{\prime }}^1$, $V_{j_2{j}_2^{\prime }}^2$, and $V_{j_3{j}_3^{\prime }}^3$, respectively, computed in the tensor decomposition (14).
• Likewise, the functions $w_{k_1^{\prime }}^1,w_{k_2^{\prime }}^2,\dots ,w_{k_{10}^{\prime }}^{10}$ are constructed by interpolation in the (common) interval (2), based on the intermediate points defined in (13), where the functions take the values $W_{k_1{k}_1^{\prime }}^1,\dots ,W_{k_{10}{k}_{10}^{\prime }}^{10}$ computed in the tensor decomposition (14).

Summarizing, the methodology requires first constructing offline (at the outset) the various ingredients appearing in tensor decomposition (14). For the specific case considered here, using the computer described at the end of the Introduction, the construction of the surrogate engine model takes 30 CPU seconds to perform the HOSVD (14), plus $\sim 3$ CPU days to obtain the required sensors data using the full PROOSIS-based engine model. The latter is a very large computational cost. However, such offline construction must be performed just once for the brand-new engine and, perhaps, once more after a very important maintenance task if the engine is strongly degraded and a new baseline state needs to be considered. Once the surrogate engine model has been constructed, its online operation only requires performing some algebraic operations and one-dimensional interpolations, as explained above. In particular, computing the functions (16) for each specific case in the aforementioned computer takes $\sim 0.35$ CPU seconds.

Now, in order to check how well the surrogate engine model approximates the original full engine model, we simulate 200 cases with the engine model, for randomly chosen flight altitudes, Mach numbers, turbine inlet temperatures, and degradations, in the ranges (1)-(2). The resulting sensors outcomes are compared with their counterparts obtained with the surrogate engine model in terms of the relative root mean square (RMS) and the maximum differences between both, defined as

$$\frac{\left|\right|{\mathbf{y}}^{\mathrm{full}\mathrm{model}}-{\mathbf{y}}^{\mathrm{surrogate}}{\left|\right|}_{\mathrm{RMS}}}{\left|\right|{\mathbf{y}}^{\mathrm{full}\mathrm{model}}{\left|\right|}_{\mathrm{RMS}}}\mathrm{and}\frac{\left|\right|{\mathbf{y}}^{\mathrm{full}\mathrm{model}}-{\mathbf{y}}^{\mathrm{surrogate}}{\left|\right|}_{\infty }}{\left|\right|{\mathbf{y}}^{\mathrm{full}\mathrm{model}}{\left|\right|}_{\mathrm{RMS}}},$$

(17)

respectively. Here, the sensors vector $\mathbf{y}$ collects the eleven sensor outcomes, the RMS norm $\left|\right|·{\left|\right|}_{\mathrm{RMS}}$ is the usual Euclidean norm divided by $\sqrt{N_{Sensors}}$, and $\left|\right|·{\left|\right|}_{\infty }$ is the maximum of the absolute values of the vector elements. The relative RMS and maximum differences are $\sim 0.03$ and $\sim 0.15$, respectively, which is acceptable in the context of this paper.

Let us now consider some preliminary ingredients that will be needed to develop the two condition monitoring tools. For convenience, the outcomes of the surrogate engine model are recalled, rewriting Equation (15) in vector form as

$$\mathbf{y}=\mathbf{F}(h,M,T_{4t},\mathbf{x}),$$

(18)

where the sensors vector $\mathbf{y}$ is as defined above and the degradations vector $\mathbf{x}$ collects the values of the ten degradations. Our condition monitoring methods will consider the vector equation

$$\mathbf{F}(h,M,T_{4t},\mathbf{x})={\mathbf{y}}^{\mathrm{measur}.},$$

(19)

where ${\mathbf{y}}^{\mathrm{measur}.}$ collects the sensors outcomes. The aim is to solve this vector equation and obtain $T_{4t}$ and $\mathbf{x}$ for given values of h, M, and ${\mathbf{y}}^{\mathrm{measur}.}$.

In the present test case, we have $N_{\mathrm{degrads}}=10$ degradations and $N_{\mathrm{sensors}}=11$ sensors, but the methodology developed below is readily extended to general values of $N_{\mathrm{degrads}}$ and $N_{\mathrm{sensors}}$, provided that $N_{\mathrm{sensors}}\ge N_{\mathrm{degrads}}+1$. In fact, if $N_{\mathrm{sensors}}>N_{\mathrm{degrads}}+1$, then some linear systems appearing below are over-determined and must be solved using the pseudo-inverse, as done in [19].

Since the temperature $T_{4t}$ needs to be computed along with the ten degradations, it is convenient to consider $T_{4t}$ as a new 0-th unknown, and scale it in such a way that the scaled value becomes of order one, which is done transforming $T_{4t}$ as

$$x_0=(T_{4t}-1450)/1450,$$

(20)

where the value 1450 K is taken as a typical value in cruise conditions, which usually range between 1400 K and 1500 K in commercial aviation.

Equation (19) will be approximately solved in the next two sections using gradient-like methods. These methods require the Jacobian of the left hand side of (19) with respect to the unknowns, which are the temperature $T_{4t}$ and the degradations. The corresponding Jacobian is defined as

$${\mathbf{J}}_{T_{4t},\mathbf{x}}=[{\mathbf{j}}_0,{\mathbf{j}}_1,\dots ,{\mathbf{j}}_{10}],$$

(21)

and approximated using forward finite differences. The first column, associated with $T_{4t}$, is computed as

$${\mathbf{j}}_0\simeq [\mathbf{F}(h,M,1450(1+x_0+\delta ),\mathbf{x})-\mathbf{F}(h,M,1450(1+x_0),\mathbf{x})]/\left(1450\delta \right).$$

(22)

Note that, although this is written in terms of the variable $x_0$ defined defined in (20), it gives an approximation of the derivative of the sensors outcomes with respect to $T_{4t}$ (instead of the derivative with respect to $x_0$), which will permit computing $T_{4t}$ with more accuracy than the degradations; see below. The columns for the degradations are approximated by

$${\mathbf{j}}_m\simeq [\mathbf{F}(h,M,1450(1+x_0),\mathbf{x}+\delta {\mathbf{I}}_m)-\mathbf{F}(h,M,1450(1+x_0),\mathbf{x})]/\delta ,$$

(23)

for $m=1,\dots ,10$, where ${\mathbf{I}}_m$ is the m-th column of the $10\times 10$ unit matrix.

In (22)-(23), the finite differences increment is selected as $\delta ={10}^{-8}$, although the methodology is insensitive to the chosen value of $\delta$, provided that it is conveniently small. This is contrast to what happened in the method developed in [19], where $\delta$ cannot be chosen too small due to round-off artifacts in the outcomes of the PROOSIS solver.

Each computation of the Jacobian ${\mathbf{J}}_{T_{4t},\mathbf{x}}$ requires applying the surrogate model $10+2=12$ times and takes $\sim 4.2$ CPU seconds, which scales with the computational cost of the surrogate model ($\sim 0.35$ CPU seconds). Concerning the condition number of ${\mathbf{J}}_{T_{4t},\mathbf{x}}$, it ranges between 65 and 70 (depending on the degradations), which is acceptable for the purpose of this work.

In the next two subsections, we develop the two aforementioned tools, one that is appropriate for small degradations, and another useful for larger degradations. In both cases, the aim is to obtain the temperature $T_{4t}$ and the degradation vector $\mathbf{x}$ for given values of the altitude h and the Mach number M, using the set of sensors data ${\mathbf{y}}^{\mathrm{measur}.}$. The turbine inlet temperature will be assumed to be in a range of $\pm 50$ K around a known reference value provided by the engine instrumentation. This range more than covers the accuracy of the engine instrumentation to estimate $T_{4t}$.

3.3. Condition monitoring tool 1: a purely linear method

In principle, this linear method should provide good results for small degradations and $T_{4t}$ close to a reference value, denoted as $T_{4t}^{\mathrm{ref}.}$. Such value will be taken as that provided by the engine instrumentation. The unknown vector, denoted as $\widehat{\mathbf{x}}$, will be of size 11. Its first component, ${\widehat{x}}_0$, is the scaled value of $T_{4t}$, defined in (20), and the remaining components are those of the degradations. Likewise, the reference value of this vector, which corresponds to the undegraded engine, is denoted as ${\widehat{\mathbf{x}}}^{\mathrm{ref}.}=({\widehat{x}}^{\mathrm{ref}.},{\widehat{x}}_1^{\mathrm{ref}.},\dots ,{\widehat{x}}_{11}^{\mathrm{ref}.})$, with the first component of this vector, ${\widehat{x}}_0^{\mathrm{ref}.}$, obtained from $T_{4t}^{\mathrm{ref}.}$ using (20) and the remaining ten components all equal to zero.

In the linear approximation, the unknown vector $\widehat{\mathbf{x}}$ is written as

$$\widehat{\mathbf{x}}={\widehat{\mathbf{x}}}^{\mathrm{ref}.}+\Delta \widehat{\mathbf{x}},$$

(24)

where the perturbation $\Delta \widehat{\mathbf{x}}$ is computed by solving the following linear problem

$${\mathbf{J}}^{\mathrm{ref}.}\Delta \widehat{\mathbf{x}}={\mathbf{y}}^{\mathrm{measur}.}-{\mathbf{y}}^{\mathrm{ref}.},$$

(25)

which results by linearizing (19) around the reference ${\widehat{\mathbf{x}}}^{\mathrm{ref}.}$. Thus, the Jacobian ${\mathbf{J}}^{\mathrm{ref}.}$ is as defined in (21), computed for the reference value of $T_{4t}$, $T_{4t}^{\mathrm{ref}.}$, and zero degradations, and ${\mathbf{y}}^{\mathrm{ref}}$ is the sensors vector for the reference vector $\widehat{\mathbf{x}}={\widehat{\mathbf{x}}}^{\mathrm{ref}}$.

Since, as anticipated, the $11\times 11$ square matrix ${\mathbf{J}}^{\mathrm{ref}.}$ is reasonably well conditioned, the unknown vector $\Delta \widehat{\mathbf{x}}$ is readily obtained by solving the linear system (25).

For the present test case, the offline preparation of the method requires computing the reference sensors vector ${\mathbf{y}}^{\mathrm{ref}.}$ (which is calculated running the surrogate engine model once) and, also, to compute the Jacobian matrix ${\mathbf{J}}^{\mathrm{ref}.}$. These take $\sim 4$ CPU seconds. The online operation of the tool is extremely fast, since it only requires solving the linear system (25), which consumes ${10}^{-3}$ CPU seconds.

It turns out (after a very large number of computations, omitted here for the sake of brevity) that the accuracy of this linear method is much higher for the computation of the turbine inlet temperature than for the degradations. This is because, as anticipated, in the computed Jacobian (21), the first row, given by (22), approximates the derivatives with respect to $T_{4t}$, which is much larger than the degradations. Moreover, the method computes $T_{4t}$ with great accuracy, while it does not compute the degradations well. Thus, nonlinearity plays an essential role in condition monitoring calculations of the engine degradation.

When analyzing how small the degradations should be for the good computation of $T_{4t}$, we have concluded that the method is able to accurately compute this temperature even for $O\left(1\right)$ degradations, in the whole range, between 0 and 2%, accounted for in the surrogate engine model. For these non-small values of the degradations, it may happen that some of the computed values of either the temperature $T_{4t}$ or the degradations are outside the ranges, defined in (1)-(2), in which the surrogate engine model applies. In this case, the ‘problematic’ values are set as equal to the nearest value of the boundary of the ranges.

This tool can be continuously used to compute $T_{4t}$ along the engine operation, but not to compute degradations accurately. These, instead, can be computed using the tool 2 developed in the next subsection. In this sense, the present tool 1 complements well the tool 2.

In any event, the present tool is able to compute, extremely fast (in $\sim {10}^{-3}$ CPU seconds, as anticipated), the temperature $T_{4t}$ very accurately. Thus, this tool can be used to compute this temperature in continuous real-time, showing its value in a monitor installed in the aircraft cockpit.

3.4. Condition monitoring tool 2: a global, constraint Newton-based method

As with the linear method, the aim is to obtain the temperature $T_{4t}$ and the degradation vector $\mathbf{x}$ for given values of h and M, using the set of sensors data ${\mathbf{y}}^{\mathrm{measur}.}$. Also, the temperature is rescaled according to (20), the various Jacobians that will be needed below are computed as indicated in (21)-(23), and the unknown vector $\widehat{\mathbf{x}}$ collects both the temperature $T_{4t}$ and the degradations.

The method solves the nonlinear system (19) iteratively, beginning with a prescribed initial condition. At the k-th iteration, the new value of $\widehat{\mathbf{x}}$, ${\widehat{\mathbf{x}}}_{k+1}$, is computed in terms of the former value as

$${\widehat{\mathbf{x}}}_{k+1}={\widehat{\mathbf{x}}}_k+\Delta {\widehat{\mathbf{x}}}_k,$$

(26)

where $\widehat{\Delta }{\mathbf{x}}_k$ is given by the following linear system

$${\mathbf{J}}_k\Delta {\widehat{\mathbf{x}}}_k={\mathbf{y}}^{\mathrm{measur}.}-{\mathbf{y}}_k^{\mathrm{ref}.}.$$

(27)

Here, ${\mathbf{J}}_k$ is the Jacobian computed at ${\widehat{\mathbf{x}}}_k$. As in the purely linear method, it may happen that at some iteration, some of the computed values of either the temperature $T_{4t}$ or the degradations are outside the ranges defined in (1)-(2), in which the surrogate engine model applies. In this case, the ‘problematic’ values are set as equal to the nearest value of the boundary of the ranges.

The method is considered as converged if the following conditions hold

$$\parallel {\widehat{\mathbf{x}}}_{k+1}-{\widehat{\mathbf{x}}}_k{\parallel }_{\mathrm{RMS}}<{\epsilon }_1\mathrm{and}{∥{\mathbf{y}}^{\mathrm{measur}.}-{\mathbf{y}}_k∥}_{\mathrm{RMS}}<{\epsilon }_2,$$

(28)

for some small thresholds ${\epsilon }_1>0$ and ${\epsilon }_2>0$. Also, a maximum number of iterations is allowed. If, after these, conditions (28) are not met, the obtained solution is not accepted and the process is repeated. The latter situation almost never occurs if the maximum number of iterations is reasonable because, as it will be seen in Section 4, the process generally converges in less than 5 iterations.

Comparison of (25) and (27) shows that, conceptually, the present Newton-based method could be seen as iterating the linear method used in the condition monitoring tool 1. However, their computational costs are quite different. This is because, while the Jacobian ${\mathbf{J}}^{\mathrm{ref}.}$ is computed offline just once in the linear method, whose online operation only involves solving the system (25), which is extremely inexpensive ($\sim {10}^{-3}$ CPU seconds), each iteration of the present method requires computing online the updated Jacobian ${\mathbf{J}}_k$. In other words, the online computational cost of each iteration of the present Newton-based method equals the offline ($\sim 4$ CPU seconds) plus the online computational costs of the linear model.

It must be noted that, once the iterations are close to the right solution, since this is a second order method [43], convergence is very fast. Thus, performing an additional iteration, the error decreases dramatically.

4. Representative condition monitoring results

The condition monitoring tools developed in Section 3.3 and Section 3.4 are now applied to several representative cases, in which the exact values of the turbine inlet temperature and the degradations are selected randomly in the range (1)-(2). In principle, the results obtained using the tools 1 and 2 will be obtained for the following representative values of the altitude, the Mach number, and the reference value of $T_{4t}$

$$h=34,000\mathrm{feet},M=0.82,\mathrm{and}{T}_{4t}^{\mathrm{ref}.}=1,400\mathrm{K}.$$

(29)

However, the robustness in the results will always be checked simulating a large number of cases, selecting h, M, and $T_{4t}^{\mathrm{ref}.}$ randomly in the ranges (1).

4.1. Results using the condition monitoring tool 1

To begin with, we apply this tool to two representative cases, one with small degradations and another with $O\left(1\right)$ degradations. Results for these two cases are displayed in Table 3.

Note in this table that, in both cases, the error for ${\widehat{x}}_0$ (giving the temperature $T_{4t}$) is about three orders of magnitude smaller than for the degradations. This was anticipated in Section 3.3 and means that the computed temperature differs from its exact counterpart by 1450 times $\mathrm{Error}\left({\widehat{x}}_0\right)\sim$ 0.1 K in both cases. This error is extremely small. The errors in the degradations, instead, are only moderately small, meaning that the method does not compute degradations well, although the obtained values are useful, see below. It is interesting that the method computes four of the very small degradations as equal to zero, which is due to the fact that the purely linear approximations computed by the method were negative and, thus, the method has set them equal to zero, as explained in Section 3.3. Likewise, for $O\left(1\right)$ degradations, the method sets one of the degradations as equal to zero and another as equal to two, because the linearly computed degradations were negative and larger than two, respectively.

Concerning the computational cost, the offline computation of the Jacobian requires $\sim 4$ CPU seconds, while the online operation of the tool just takes $\sim {10}^{-4}$ CPU seconds.

The computations above have been repeated many times, obtaining results are consistent with those in Table 3.

Let us now check the robustness of the linear tool to noisy sensors data. To this end, we add random noise with zero mean of size 0.005 (i.e., 0.5%, which is a typical accuracy of the engine instrumentation) to the sensors data. In order to clean the effect of errors in the computation of $T_{4t}$ and the degradations, we repeat computations 20 times (adding independent random noise each time) and take the means of the obtained values of $T_{4t}$ and the degradations. Note that this can be implemented in practice by using moving means on the turbine inlet temperature and the sensor outcomes during the engine operation. The results obtained using this averaged linear method to the two cases considered in Table 3 are displayed in Table 4.

Comparison with Table 3 shows that accuracies are comparable in both cases, showing that the simple averaging method filters errors well. In particular, the method computes $T_{4t}$ (which is the goal of this tool) with a very small error, $\sim 0.1$ K, while the degradations are not computed so well. Concerning the computational cost, this is twenty times larger than in the clean case, but still extremely small, namely $\sim {10}^{-3}$ CPU seconds. However, the above mentioned moving means are applied to the outcomes of the engine instrumentation, whose sampling frequency is $\sim 100$ Hz. Thus, in practice, results obtained by this method are given in $20\times 0.01=0.2$ seconds, which means that these results can be displayed in continuous real-time.

As in the unperturbed case, computations have been repeated many times, obtaining results that are systematically consistent with those in Table 4.

Summarizing the above:

• The linear tool is robust in connection with random noise added to the clean sensors data.
• Nonlinear effects cannot be ignored in the considered range of degradations.
• In spite of that, the tool is able to compute the temperature $T_{4t}$ very fast and very accurately. Thus, this tool can be used to obtain this temperature in continuous real-time, continuously showing its value in a monitor installed in the aircraft cockpit, as anticipated.

4.2. Results using the condition monitoring tool 2

A first test using tool 2 is presented in Table 5, where exactly the same two cases considered in Table 3 are addressed.

In both cases, the method requires just four iterations, which take $\sim 17$ CPU seconds. As can be seen, the errors in all degradations are extremely small, in spite of the fact that the initial condition is generally far from the exact solution. Also, errors in ${\widehat{x}}_0$ are about three orders of magnitude smaller than those for the degradations, meaning once more that the temperature is calculated with extreme precision.

As in the linear method, the robustness of the tool 2 to noise is now tested by adding random noise, of size $0.005$ with zero mean, to the sensors data and repeat computations twenty times. In each computation, we apply the Newton-based method, taking as initial conditions the outputs of the linear method displayed in Table 4. The obtained results for $T_{4t}$ and the degradations are averaged. Such averages can be implemented in practice using moving means in both $T_{4t}$ and the degradations. This ‘averaged Newton-based method’ requires twenty times the CPU time of the basic Newton-based method ($\sim 17$ CPU seconds), which is $\sim 5$ CPU minutes. Applying this averaged Newton-based method to the two cases considered in Table 5 gives the results displayed in Table 6.

Note that errors are much smaller than their counterparts using the averaged linear method, given in Table 4. In fact, both the temperature $T_{4t}$ and the degradations are computed with more than enough accuracy in practical terms, which permits discerning safely which degradations remain small from those that require attention.

In order to further check the performance of the averaged Newton-based method, a large number of cases have been run, always selecting randomly the altitude, Mach number, and exact values of the unknowns. In all cases, the method has converged to the exact values of the unknowns in a number of iterations not larger than five. Thus, the method is both robust and accurate enough in practical terms.

Summarizing the results obtained with the condition monitoring tool 2:

• The global, constraint Newton-based method is robust in connection with random noise added to the clean data.
• This tool takes nonlinear effects into account and, thus, it gives precise results for both the turbine inlet temperature and the degradation parameters.
• The tool gives results in in-flight real-time, meaning that it can be applied several times in each flight. This permits recalculating the engine condition when some degradations are not small enough.

5. Concluding remarks

A methodology has been developed for monitoring aeroengines condition. First, a surrogate engine model was constructed, by applying a tensor decomposition plus interpolation to a multidimensional dataset containing the aeroengine sensor outcomes computed by an engine model. Such surrogate model can advantageously substitute the full engine model for two main reasons. First, (i) the surrogate model does not give erroneous/spurious sensors outcomes, as the full engine model does due to convergence difficulties. Secondly, (ii) once the surrogate engine model has been constructed, it gives results in a computationally inexpensive way, so that it can be installed as digital twin, meaning that the full engine model is not needed anymore in the on board aircraft software. The full engine model needs only be used to update the surrogate engine model after a major maintenance, which seldom occurs along the engine lifetime.

Using the surrogate engine model, two gradient-like condition monitoring tools were constructed. The first tool was based on a linear approximation around the baseline state, and able to very precisely compute the turbine inlet temperature, $T_{4t}$, in continuous real-time. Thus, it is able to continuously show, very precisely, the instantaneous value of $T_{4t}$ in a monitor installed in the aircraft cockpit. This represents a considerable added value of the present methodology, since $T_{4t}$ is a very important property for various reasons, including the aeroengine control. However, $T_{4t}$ cannot be precisely computed by the current aeroengine instrumentation in civil aviation. Our tool, instead, computes this temperature very precisely and could well be incorporated to improve the engine instrumentation. On the other hand, even though this tool does not compute degradations accurately, it gives a first approximation of them that permits discerning, whether all degradations remain small or some of them have abnormally increased. In the latter case, the second engine monitoring tool must be used to precisely compute all degradations.

The second condition monitoring tool, instead, takes nonlinearity into account and precisely computes both the turbine inlet temperature and the degradation parameters. This tool operates in in-flight real-time, meaning that it can be used several times in each flight. In other words, this tool can be used to either safely ensure that degradations remain conveniently small, or to identify the in-flight evolution of some non-small degradations that must be taken care of before they degrade too much. In the later case, the tools developed in [19], which give precise results for very large degradations, can be used. Thus, both dangerous events and costly maintenance are avoided.

It is important that both condition monitoring tools are robust in connection with random noise of realistic size added to the sensors data. Besides the considered aeroengine, the developed methodology also applies to other aeroengines and, moreover, to many complex mechanical systems.