Nonlinear control study on aircraft engine using a generalized minimum variance approach

Abstract

With the requirements of the advanced propulsion systems for aircraft control performance increasing, it has attracted the attention to design more efficient controller that is based on the nonlinear model over the flight envelope directly. A nonlinear aircraft engine control method based on generalized minimum variance approach is studied. A nonlinear autoregressive with extra input aircraft engine model for a certain turbofan engine is built and is chosen as the control plant in this study. Then, a controller is designed for the identification model based on the nonlinear generalized minimum variance theory. And simulation results manifest the controller to satisfy the requirements of the tracking with the transient process and the regulation with the steady process.

Introduction

With the relevant technology fields developed in, for example, machinery, electronics, materials, and so on, the performance of aircraft engine is continuously improved, and the number of control variables is also increasing. ¹ In traditional hydraulic system, engine control laws were designed primarily using classical control design methods, namely, frequency response and/or time-domain methods. It was always single loop, which means one input variable and one output variable. ² However, today, the number of advanced propulsion system control variables has also been reached 15–25.^3,4 The F119 has eight control variables in 1990s, ⁵ and the PW1128 engine equipped in the F16-XL also has up to 24 control variables. ⁶ The linear control theories about multivariable are abundant, such as linear quadratic Gaussian (LQG) or loop transfer recovery (LTR)^7–9 and H_∞ robust control.^10–12 The nonlinear control such as sliding mode control,^13,14 fuzzy control,^15–17 adaptive control,^18–20 and neural networks^21,22 has also been investigated for engine control. It is reported that American engineers had researched the linear quadratic method to apply it for the engine control and had adopted the LQG or LTR method in the Design Method for Integrated Control Systems (DMICS) research programs. ²³ Apart from this, there is no other information reported, and the traditional approach gain scheduling is still the primary method to solve the problem in practical project. ²⁴

This phenomenon can be attributed to the following reasons. First, parts of the above advanced control algorithms are based on linear system theory that is extremely inefficient for such a strong nonlinear control of aircraft engine. Furthermore, due to the limited speed of engine control computation unit, most of the above methods could not achieve their owned performance. But, however, as the deepening study of nonlinear theory, breakthroughs come into reality in optical fiber technology, ²⁵ composites, ²⁶ and electronic equipment of heat resistant. ²⁷ And it makes adopting the nonlinear control algorithm as a trend in engine design.

In this article, with the analysis of aircraft engine characteristics, a structure of nonlinear autoregressive with extra input (NARX) model has been identified using the data of a turbofan engine, and nonlinear generalized minimum variance (NGMV) control method^28,29 is introduced to design the nonlinear control method for the turbofan engine.

This article is structured as follows: section “Engine model description” identifies the nonlinear model with the structure of NARX. Section “NGMV controller design” gives the complete control structure and formulates the solution by the method of NGMV. In section “Compare the performance of the controllers,” simulations are given to support the statements mentioned above, including the situation with disturbance and nondisturbance. Section “Conclusion” concludes this article.

Engine model description

It is difficult to obtain one analytical mathematical model which can describe the aircraft engine over the flight envelop. This article considers identifying model by the data from idle to take-off. And it is conducted under the condition of precision guaranteed.

NARX model is understood as an extension of the linear autoregressive with extra input (ARX) structure. It is a kind of identification method that combines ARX model and neural networks with multi-step prediction. ³⁰ NARX includes feedback connections enclosing several layers of the network. The defining equation of the NARX model is as follows

$$y_p(t)=f[y_1(t-1),\dots ,y_1(t-n_{a_{p1}}),\dots ,y_i(t-1),\dots ,y_i(t-n_{a_{pi}}),u_1(t-n_{k_{p1}}),\dots ,u_1(t-n_{k_{p1}}-n_{b_{p1}}+1),\dots ,u_j(t-n_{k_{pj}}),\dots ,u_j(t-n_{k_{pj}}-n_{b_{pj}}+1)]$$

(1)

where $y(t)$ is a 1-by-i vector and represents the output of the identification, $u(t)$ is 1-by-j vector and represents the input of the identification model, $f(·)$ is nonlinear function with limited number of inputs and outputs, $n_a$ is an i-by-i matrix and represents the input order of the model, $n_b$ is an i-by-j matrix and represents the input order of the model, $n_k$ is an i-by-j matrix and represents the input order of the model $y_p$ is the p-th element of the y(t).

The structure of the NARX model is shown in Figure 1, which is combined by regressors and nonlinear estimator in series connection; at the same time, the output is fed back to the regressors.

Figure 1.

NARX structure.

The engine inputs used in this study are fuel flow $W_f$ and nozzle area $A_8$ . The rotation speeds of the compressor and fan, $n_H$ and $n_L$ , respectively, are the output variables. The target prototype is a two-spool turbofan engine in this study. In order to make the research work easy, a well-tuned nonlinear component-level engine model for the turbofan engine mentioned above is employed, which can reflect the real engine work state at a very high accuracy above idle state. After normalizing the data, the model is identified as follows: Different types of NARX models are studied by choosing various $n_a$ , $n_b$ , and $n_k$ , as shown in Table 1. Multiple datasets from the nonlinear component-level model are utilized to identify the NARX model for engine, and 12 models are obtained, as shown in Table 1. Figure 2 shows the results of the identified model, where $m_1$ means the best model we acquire and zv is the data to validate the model. It is also the data collected during the similar simulation and are applied to compare with the identification model. Fit represents the model matching degree, so $m_1$ is selected as the control plant.

Table 1.

Identification results.

Model	$n_a$	$n_b$	$n_k$	$n_{H\mathrm{Fit}\%}$	$n_{L\mathrm{Fit}\%}$
m 1	[1 0; 0 1]	[1 2; 1 1]	[1 1; 0 1]	95.7	90.56
m 2	[1 0; 0 1]	[1 1; 1 1]	[1 2; 0 2]	95.6	87.08
m 3	[1 0; 0 1]	[1 2; 1 1]	[1 2; 0 2]	95.32	87.08
m 4	[1 0; 0 1]	[1 2; 1 1]	[1 2; 1 1]	95.32	90.59
m 5	[1 0; 0 1]	[1 3; 1 1]	[1 2; 1 2]	87.96	88.65
m 6	[1 1; 1 1]	[1 2; 1 1]	[1 2; 1 1]	92.71	86.20
m 7	[2 0; 0 2]	[1 0; 0 1]	[1 3; 3 1]	93.85	78.49
m 8	[3 0; 0 3]	[1 0; 0 1]	[1 4; 4 1]	90.91	83.59
m 9	[4 0; 0 4]	[1 2; 2 1]	[1 0; 0 1]	87.24	80.62
m 10	[5 0; 0 5]	[1 0; 0 1]	[1 0; 0 1]	85.46	81.71
m 11	[1 0; 0 1]	[1 1; 1 1]	[1 2; 1 5]	95.63	86.65
m 12	[2 1; 1 2]	[2 2; 2 2]	[1 2; 2 1]	78.32	79.69

Figure 2.

Results of the best identified NARX model: (a) $n_H$ and (b) $n_L$ .

NGMV controller design

It is essential to introduce the framework of the complete NGMV control system first, as shown in Figure 3.

Figure 3.

Framework of the NGMV control system.

A complete control system generally consists of the following types of components: a plant to be controlled, a number of sensors, and a controller (or calculator). ³¹ Considering the integrity of the complete system, here the reference signal model and the sensor model are given as equation (2)

$$\{\begin{array}{c}{x}_r(t+1)=A_r{x}_r\left(t\right)+B_r(u_r\left(t\right)+\omega \left(t\right))\hfill \\ \multicolumn{1}{c}{r\left(t\right)=C_r{x}_r\left(t\right)+D_r{u}_r\left(t\right)}\end{array}$$

(2)

where

$$A_r=\left[\begin{array}{cc}-3.33333\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill -3.33333\hfill \end{array}\right]$$

$$B_r=\left[\begin{array}{cc}1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 1\hfill \end{array}\right]$$

$$C_r=\left[\begin{array}{cc}3.3333\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 3.3333\hfill \end{array}\right]$$

$$D_r=\left[\begin{array}{cc}0\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 0\hfill \end{array}\right]$$

and $\omega (t)$ is the bounded disturbance.

As the control plant is a dual-input and dual-output model, the reference signal is ordered to match with it. It is impractical to obtain the ideal step signal in practice; therefore, the signals similar to the first-order inertia link are obtained instead of it.

With the characteristics of sensors, parameters are selected as equation (3)

$$\{\begin{array}{c}x(t+1)=\mathrm{Ax}\left(t\right)+B{u}_0\left(t\right)\hfill \\ \multicolumn{1}{c}{y\left(t\right)=\mathrm{Cx}\left(t\right)+D{u}_0\left(t\right)}\end{array}$$

(3)

where

$$A=\left[\begin{array}{cc}-1000\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 1000\hfill \end{array}\right]$$

$$B=\left[\begin{array}{cc}1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 1\hfill \end{array}\right]$$

$$C=\left[\begin{array}{cc}1000\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 1000\hfill \end{array}\right]$$

$$D=\left[\begin{array}{cc}0\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 0\hfill \end{array}\right]$$

And it is represented as $W_{\mathrm{ok}}$ .

In Figure 3, $W_{\mathrm{lk}}$ represents the identified nonlinear model hereinabove; here, we use $u_0=W_{\mathrm{lk}}(u)$ to describe it. The module of time delay is selected to be the unit delay $z^{-1}$ . w and v are the disturbances that refer to the white noise signal with zero mean and unit variance. This section designs the controller for the generalized controlled plant of the dashed box, as shown in Figure 3.

Based on the optimal control theory, the performance index is chosen as J, which can be presented as equation (4)

$$J=E\{{\varphi }_0^T\left(t\right){\varphi }_0\left(t\right)\}$$

(4)

As shown in equation (5), the ${\varphi }_0(t)$ in equation (4) is weighted by the feedback error $e(t)$ and control input $u(t)$ , where $P_c$ and $F_{\mathrm{ck}}$ are the weighted variables

$${\varphi }_0\left(t\right)=P_{c}e\left(t\right)+F_{\mathrm{ck}}u\left(t\right)$$

(5)

Then, according to the theory of the optimal estimation and considering that ${\hat{\varphi }}_0(t+k|t)$ and ${\tilde{\varphi }}_0(t+k|t)$ are orthogonal, ³² the following equation (6) will be held

$$\begin{array}{c}\multicolumn{1}{c}{J=E\left\{{\hat{\varphi }}_0^T(t+k|t){\hat{\varphi }}_0(t+k|t)\right\}}\\ \multicolumn{1}{c}{+E\left\{{\tilde{\varphi }}_0^T(t+k|t){\tilde{\varphi }}_0(t+k|t)\right\}}\end{array}$$

(6)

where ${\hat{\varphi }}_0(t)$ is estimated value and ${\tilde{\varphi }}_0(t)$ is estimation error. And the estimation error is independent of the input signal, so the minimum value J satisfies the following condition, ${\hat{\varphi }}_0(t+k|t)=0$ . And k is the delay step. Let $k=1$ , and considering the description in Figure 3, equation (5) can be presented as equation (7)

$$\begin{array}{c}\multicolumn{1}{c}{{\hat{\varphi }}_0(t+1|t)=P_c[\hat{r}(t+1|t)-C\hat{x}(t+1|t)-D{u}_o\left(t\right)]}\\ \multicolumn{1}{c}{+F_{\mathrm{ck}}u\left(t\right)}\end{array}$$

(7)

where $P_c$ and $F_{\mathrm{ck}}$ are the weighted variables and $\hat{r}(t+1|t)$ is the best estimated value of the reference signal, which can be presented as equation (8)

$$\hat{r}(t+1|t)=C_r{x}_r(t|t)$$

(8)

$\hat{x}(t+1|t)$ is the best estimated value of the sensor, and it satisfies, as shown in equation (9)

$$\hat{x}(t+1|t)=A\hat{x}(t|t)+z^{-1}B{u}_0\left(t\right)$$

(9)

In equation (9), $\hat{x}(t|t)$ could be obtained by the Kalman filter.

According to the theory of Kalman ³³ filter, as shown in equations (10) and (11)

$$\begin{array}{c}\multicolumn{1}{c}{\mathrm{Time}\mathrm{Update}:}\\ \multicolumn{1}{c}{\{\begin{array}{c}\hat{x}(t+1|t)=A\hat{x}(t|t)+B{u}_0(t-1)\hfill \\ \multicolumn{1}{c}{P(t|t-1)=\mathrm{AP}(t-1|t-1)A^T+Q}\end{array}}\end{array}$$

(10)

$$\begin{array}{c}\multicolumn{1}{c}{\mathrm{Measurement}\mathrm{Update}:}\\ \multicolumn{1}{c}{\{\begin{array}{c}K\left(t\right)=P(t|t-1)C^T\left(t\right){[R+C\left(t\right)P(t|t-1)C^T\left(t\right)]}^{-1}\hfill \\ \multicolumn{1}{c}{\hat{x}(t|t)=\hat{x}(t|t-1)+K\left(t\right)[y\left(t\right)-C\left(t\right)\hat{x}(t|t-1)]}\\ \multicolumn{1}{c}{P(t|t)=[I-K\left(t\right)C\left(t\right)]P(t|t-1)}\end{array}}\end{array}$$

(11)

where Q and R are the variances of the noise of the time update process w and the noise of the measurement update process v, respectively; $P(t)$ is the minimum mean square error matrix; and $K(t)$ is the Kalman gain.

From the condition in equation (8), substitute $K(t)$ into the following equation (12)

$$\hat{x}(t+1|t+1)=\hat{x}(t+1|t)+K\left(t\right)[e(t+1)-\hat{e}(t+1|t)]$$

(12)

where

$$\hat{e}(t+1|t)=C_{r}x(t|t)-C\hat{x}(t+1|t)-D{u}_o(t-k+1)$$

(13)

After substituting equations (8) and (13) into equation (12), equation (14) could be obtained

$$\hat{x}(t|t)=G_1+G_2{u}_0\left(t\right)$$

(14)

where

$$G_1={[I-(I+K\left(t\right)C)A{z}^{-1}]}^{-1}K\left(t\right)[e\left(t\right)-z^{-1}{C}_r{x}_r(t|t)]$$

(15)

$$\begin{array}{c}\multicolumn{1}{c}{G_2={[I-(I+K\left(t\right)C)A{z}^{-1}]}^{-1}}\\ \multicolumn{1}{c}{[K\left(t\right)D+(I+K\left(t\right)C)B{z}^{-1}]z^{-1}}\end{array}$$

(16)

Hence, if $J=J_{min}$ , then optimal controller $u_{\mathrm{opt}}(t)$ is deduced. Substituting equations (8), (9), and (14) into equation (7) which contains the unknown solved above, equation (17) holds

$$\begin{array}{c}\multicolumn{1}{c}{P_c\{C_r{x}_r(t|t)-\mathrm{CA}\hat{x}(t|t)-[C{z}^{-1}B+D]\left(W_{\mathrm{lk}}{u}_{\mathrm{opt}}\right)\left(t\right)\}}\\ \multicolumn{1}{c}{+\left(F_{\mathrm{ck}}{u}_{\mathrm{opt}}\right)\left(t\right)=0}\end{array}$$

(17)

To avoid the process of inversing the nonlinear module, the optimal controller is given by

$$\begin{array}{c}\multicolumn{1}{c}{u_{\mathrm{opt}}=-F_{\mathrm{ck}}^{-1}{P}_c\{C_r{x}_r(t|t)-\mathrm{CA}\hat{x}(t|t)}\\ \multicolumn{1}{c}{-[C{z}^{-1}B+D]\left(W_{\mathrm{lk}}{u}_{\mathrm{opt}}\right)\left(t\right)}\end{array}\}$$

(18)

where $F_{\mathrm{ck}}$ and $P_c$ are the weighted variables and taken as the regulating parameters of the controller freely by equation (18).

Compare the performance of the controllers

The simulations are carried out in this section based on the above theoretical analysis. To verify the effect of the controller, the NGMV controller and the traditional proportional–integral–derivative (PID) controller are compared.^34,35 To consider the more realistic case, bounded disturbance is injected to examine the robustness of the system.

Application without disturbance

To consider the stability of the control system, a proportional–integral (PI) controller parameter is chosen as $F_{\mathrm{ck}}$ to satisfy the stability. After tuning

$$F_{\mathrm{ck}}=\left[\begin{array}{cc}30+\frac{10}{\mathrm{s}}\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 10+\frac{15}{s}\hfill \end{array}\right]$$

A Butterworth filter is chosen to be $P_c$ . It is a filter with 10th order, and the standard cutoff frequency is 0.9. Furthermore, PID controller is chosen as

$$\mathrm{PID}=\left[\begin{array}{cc}60+\frac{30}{s}\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 20+\frac{50}{s}\hfill \end{array}\right]$$

to compare the effect of different control methods. The results are shown in Figures 4 and 5.

Figure 4.

Comparison of $n_L$ performance for different control methods: (a) $n_L$ performance from 0 to 30 s and (b) magnified view of the acceleration to take-off.

Figure 5.

Comparison of $n_H$ performance for different control methods: (a) $n_H$ performance from 0 to 30 s and (b) magnified view of the acceleration to take-off.

Figure 4(a) shows the tracking control of the rotation speed of the low pressure between idle (76.2%) and take-off (100%), and Figure 4(b) is the magnified view of the acceleration process from 85.6% to 100%. The two control methods reach the steady state in 2 s, and no large overshoot exists during the tracking process.

Figure 5 shows the performance of the $n_H$ . As shown in the graph, NGMV is slower than PID around the acceleration from 0 to the idle (85.6%). It is due to the reason of the identified model that the scope of the input is between idle and take-off. Besides this, both controllers track the reference signal well. Figure 5(b) is the magnified view of the acceleration process from 95% to 100%. The acceleration time is about 4 s. For all this, the superiority of NGMV is shown as the following comparison.

Figure 6 manifests the two control inputs by two control methods. On one hand, both fuel flow $W_f$ and nozzle area $A_8$ appear as high-frequency spikes by the PID control. It is hard (or even impossible) to implement by the real actuators. NGMV shows a very smooth control output than PID. On the other hand, from the comparison of the $W_f$ , its consumption is much less than the PID control. It is a very important economic significance for the aircraft engine.

Figure 6.

Comparison of the control inputs: (a) $W_f$ and (b) $A_8$ .

Application with bounded disturbance

Any realistic system consists of disturbances in one form or another, and accordingly, the corresponding control design must consider certain robustness. A certain amplitude and frequency of the disturbance signal are injected in the system to verify the system. In view of the normalization processing (i.e. the range of the variables is between 0 and 1) of the controlled plant, here the amplitude of the disturbance signal with 0.1 is considered as 10% disturbance of the system. Similarly, 0.2 represents 20%. Figures 7 and 8 reveal the effects on the $n_L$ and $n_H$ by different amplitude disturbances.

Figure 7.

Comparison of the $n_L$ with different amplitude disturbances: (a) $n_L$ performance from 0 to 30 s and (b) magnified view of the acceleration to take-off.

Figure 8.

Comparison of the $n_H$ with different amplitude disturbances: (a) $n_H$ performance from 0 to 30 s and (b) magnified view of the acceleration to take-off.

It is obvious that both $n_L$ and $n_H$ fulfill the requirements of tracking control when the disturbance is <20%. The steady error is less than about 1%, as shown in Figures 7 and 8. And $n_L$ is more sensitive than $n_H$ . Even though there exist some amplitude fluctuations, it is acceptable. It is predictable that the control input of the NGMV control is better than the PID control.

As a consequence, the system possesses some property of disturbance suppression. When the amplitude of disturbance is <20%, the outputs are acceptable. As the amplitude increases, the suppressed effect becomes weak.

Conclusion

In this article, NGMV is proposed to design a nonlinear control system for turbofan engine represented by an NARX model structure. This type of model suits perfectly the situation where it is time-consuming to model the system dynamics. The proposed method demonstrates that its fuel consumption is more efficient with a smooth modulation than the traditional PID-type controllers. The method also possesses robustness to disturbances which may upset conventional design. The proposed approach is designed for multivariable control systems, variable cycle engine even adaptive engine control system can be further studied by adopting this framework for the next step of investigation.