Dynamic surface tracking controller design for a constrained hypersonic vehicle based on disturbance observer

Abstract

The tracking control problem of a flexible air-breathing hypersonic vehicle subjects to aerodynamic parameter uncertainty and input constraint is investigated by combining nonlinear disturbance observer and dynamic surface control. To design controller simply, a control-oriented model is firstly derived and divided into two subsystems, velocity subsystem and altitude subsystem based on the engineering backgrounds of flexible air-breathing hypersonic vehicle. In every subsystem, compounded disturbances are included to consider aerodynamic uncertainty and the effect of the flexible modes. Then, disturbance observer is not only used to handle the compounded disturbance but also to handle the input constraint, where the estimation error converges to a random small region through appropriately choosing the observer parameters. To sequel, the disturbance observer–based robust control scheme and the disturbance observer-based dynamic surface control scheme are developed for the velocity subsystem and altitude subsystem, respectively. Besides, novel filters are designed to alleviate the problem of “explosion of terms” induced by backstepping method. On the basis of Lyapunov stability theory, the presented control scheme can assure that tracking error converges to an arbitrarily small neighborhood around zero by rigorous theoretical analysis. At last, simulation result shows the effectiveness of the presented control method.

Introduction

Air-breathing hypersonic vehicle (AHV) features the slender geometry and relative light weight, which induces that the vibration modes significantly affect the aerodynamic forces. And the vibration modes of hypersonic vehicle were analyzed in the study by Williams et al. ¹ It is a challenge to explicit the characteristics of flight dynamics of flight AHV (FAHV) for the aerothermodynamics effects of the hypersonic speed, strong interaction among the elastic airframe, the structure dynamics, and the propulsion system. It is also difficult to measure the atmospheric properties and aerodynamic characteristics in the flight envelope of FAHV. So, the problem about the robustness and stability of flight control system becomes the central issue of FAHV. ² The flight control system design for hypersonic vehicles is a challenging task, and it has received more and more attention. ^3,4 The longitudinal analytical model of FAHV proposed in the study by Bolender and Doman ⁵ has been used for the controller design. Moreover, backstepping control technique, ^6,7 sliding mode control, ⁸ and intelligent technique ^9,10 have been used for controller design of FAHV.

It is known that input constraint can be detrimental to the stability of the control system from theoretical aspects. It has been received more and more attention in both theory and engineering areas. ^{11 –14} Due to the safe operation of flight, it is important to design the flight control system with the consideration of input constraints. If it is not considered, then the structure of vehicle may change seriously when input saturation occurs, which may cause the aircraft to disintegrate. Much research has focused on the controller design of hypersonic vehicle with input constraint. ^12,15,16 Although backstepping control was applied to design control system for AHVs, ^6,7,17 its drawback is the problem of “explosion of terms.” So to address this problem, dynamic surface control (DSC) was presented ^17,18 and applied in flight control system design field. ^12,14,19 In the study by Pan and Yu and Pan et al., ^17,18 composite learning–based adaptive DSC is designed for strict-feedback nonlinear systems with unknown parameters, where the control performance can be improved from practical asymptotic stability to practical exponential stability. In the study by Xu, ¹⁴ the DSC scheme is developed for nonlinear transport aircraft model during the process of continuous heavy cargo airdrop in case of disturbance and actuator saturation. The disturbance is handled by the disturbance observer (DOB) in each step, and the auxiliary signal is used for tackling actuator saturation. Then, in the work by Waseem et al., ¹⁹ an integral term was applied into DSC design procedure to avoid a large initial control signal, where the magnitude and rate constraints of the actuator commands were taken into account to assure feasibility. A robust adaptive DSC scheme ¹² was proposed for AHV with input saturation and uncertainty. It should be pointed out that the characteristics of aerodynamics are difficult to be measured and estimated because of the integration of aircraft with the propulsion system and interactions between external flow fields and the internal. Besides, plant parameter variations and uncertainties are the main issues to be handled at the level of control design for the lack of a broad flight dynamics database. ²⁰ For a nonlinear system with time-varying external disturbance, to achieve satisfactory closed-loop control performance, the disturbance needs to be efficiently handled. And in order to improve the compensation ability for disturbance, the DOB can be employed to facilitate the direct adaptive control for the uncertain nonlinear system. The DOB method has been proved to be effective in compensating the effects of unknown external disturbances and model uncertainties in control systems and received a great deal of attention in control society. ²¹ Besides, the DOB control (DOBC) approach provides powerful abilities in improving robustness and disturbance rejection performance. ²² It provides a flexible way to handle unknown disturbances or uncertainties of nonlinear systems. ²³ The DOBC received more and more attention. ^{24 –27}

Based on the aforementioned discussion, the contribution of this article is to design a robust dynamic surface control scheme for FAHV stabilization that considers input constraint explicitly, guarantees accurate response, and obtains compensation for the effect of parametric uncertainty. Compared with our previous work, ^28,29 where the design procedure of controller is also according to DSC, the difference in this article is concluded as follows. Firstly, nonlinear DOB technique is not only employed to approximate compounded uncertainty but also to cope with input constraint. And the estimation error can converge to random small region around zero. From theoretical and simulation aspects, the parameters of DOB are determined more simply than that of radial basis function neural network (RBFNN) used in the study by Zong et al. ²⁹ The DOB structure is simple and the designed controller is more succinctly. To sequel, in the previous study, ^28,29 the additional system was designed to tackle input constraint, but in this article, there is no need to construct additional system at the level of control design. Therefore, the controller design and stability analysis are simpler than that of other studies. ^28,29 Finally, the “explosion of terms” problem inherent in backstepping control is avoided by developing the novel first-order filter which is superior to the first-order filter used in the work by Zong et al. ²⁹

Hypersonic vehicle model

The curve-fitted model (CFM) of the hypersonic vehicle is firstly shown in this section, and then the control-oriented model (COM) is given and used for the design of control system.

Curve-fitted model

In this article, an adopted longitudinal model is originated from a model proposed in the study by Sigthorsson and Serrani, ³⁰ and it is constructed for longitudinal dynamics of an FAHV. And the model is denoted as CFM, which is described as follows ³⁰

$\dot{V}=\frac{T \cos \alpha -D}{m}-g\text{sin}\gamma$

1

$\dot{h}=V\text{sin}\gamma$

2

$\dot{\gamma }=\frac{L+T\text{sin}\alpha }{mV}-\frac{g \cos \gamma }{V}$

3

$\dot{\alpha }=q-\dot{\gamma }$

4

$\dot{q}=\frac{M_{yy}}{I_{yy}}$

5

${\ddot{\mathit{\eta }}}_i=-2{\xi }_i{\omega }_i{\dot{\mathit{\eta }}}_i-{\omega }_i^2{\dot{\mathit{\eta }}}_i+N_i,i=1,2,3$

6

Here, velocity V, altitude h, flight path angle(FPA) γ, angle of attack (AOA) α, pitch rate q are rigid body states and ${\mathit{\eta }}_1,{\mathit{\eta }}_2,{\mathit{\eta }}_3,{\dot{\mathit{\eta }}}_1,{\dot{\mathit{\eta }}}_2,{\dot{\mathit{\eta }}}_3$ are flexible states. Besides, it is obvious from equations (1) to (6) that the thrust T, drag D, lift L, pitching moment M_yy, and three generalized forces N₁,N₂,N₃ are functions about the system states and inputs, and they must be simplified to make the model analytically tractable. Flexible states ${\mathit{\eta }}_1,{\mathit{\eta }}_2,{\mathit{\eta }}_3,{\dot{\mathit{\eta }}}_1,{\dot{\mathit{\eta }}}_2,{\dot{\mathit{\eta }}}_3$ are corresponding to the first three vibrational modes and η₁,η₂,η₃ are related to the deflections of the fore-body turn angle τ₁ and aft-body vertex angle τ₂, denoted by $\Delta {\mathit{\tau }}_1$ and $\Delta {\mathit{\tau }}_2$ , respectively. ξ_i and ω_i are, respectively, the damping ratio and natural frequency of the mass-normalized generalized coordinates of the flexible structure. The vector $\mathbf{\eta }={[{\mathit{\eta }}_1,{\mathit{\eta }}_2,{\mathit{\eta }}_3]}^T$ and the matrix $E_j\in R^{1\times 3}$ depend on the fuel level. They describe the linear relationship $\Delta {\mathit{\tau }}_j=E_j\mathbf{\eta },j=1,2$ . Here, the three elements of $E_j\in R^{1\times 3}$ are the ith mass-normalized mode shape derivative response to the node locations in the free beam. Moreover, if the mass of the aircraft decreases with fuel consumption, node locations and modal frequencies of the flexible modes change, ¹ whereas the damping ratio and the mode shapes themselves remain relatively constant. Therefore, different fuel levels determine different values of $E_j\in R^{1\times 3}$ and $\Delta {\mathit{\tau }}_j$ .

Approximations of the forces and moments can be expressed as follows ³⁰

$\begin{array}{l}T\approx \overline{q}\left[\varphi C_{T,\varphi }\left(\alpha ,\Delta {\mathit{\tau }}_1,M_{\infty }\right)+C_T\left(\alpha ,\Delta {\mathit{\tau }}_1,M_{\infty },A_d\right)\right],\\ D\approx \overline{q}S{C}_D\left(\alpha ,{\delta }_e,{\delta }_c,\Delta {\mathit{\tau }}_1,\Delta {\mathit{\tau }}_2\right),\\ L\approx \overline{q}S{C}_L\left(\alpha ,{\delta }_e,{\delta }_c,\Delta {\mathit{\tau }}_1,\Delta {\mathit{\tau }}_2\right),\\ M_{yy}\approx z_{T}T+\overline{q}S\overline{c}{C}_M\left(\alpha ,{\delta }_e,{\delta }_c,\Delta {\mathit{\tau }}_1,\Delta {\mathit{\tau }}_2\right),\\ N_i\approx \overline{q}{C}_{N_i}\left(\alpha ,{\delta }_e,{\delta }_c,\Delta {\mathit{\tau }}_1,\Delta {\mathit{\tau }}_2\right),i=1,2,3.\end{array}$

7

where $M_{\infty }=V/M_0$ is free-stream Mach number and M₀ is the speed of sound at a given altitude and temperature. $\overline{q}=0.5\rho (h)V^2$ is the dynamic pressure and ρ(h) is the altitude-dependent air density. The coefficients obtained from fitting the curves are given below

$\begin{array}{l}{C}_{T,\varphi }=C_{T,\varphi }^{\alpha }\alpha +C_{T,\varphi }^{\alpha M_{\infty }^{-2}}\alpha M_{\infty }^{-2}+C_{T,\varphi }^{\alpha \Delta {\mathit{\tau }}_1}\alpha \Delta {\mathit{\tau }}_1+C_{T,\varphi }^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1{}^2}\Delta {\mathit{\tau }}_1{}^2+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_{T,\varphi }^0,\\ C_T=C_T^{A_d}{A}_d+C_T^{\alpha }\alpha +C_T^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+C_T^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_T^0,\\ C_D=C_D^{{(\alpha +\Delta {\mathit{\tau }}_1)}^2}{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2+C_D^{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}\left(\alpha +\Delta {\mathit{\tau }}_1\right)+C_D^{{\delta }_e^2}{\delta }_e^2+C_D^{{\delta }_e}{\delta }_e+C_D^{{\delta }_c^2}{\delta }_c^2+C_D^{{\delta }_c}{\delta }_c+C_D^{\alpha {\delta }_e}\alpha {\delta }_e+C_D^{\alpha {\delta }_c}\alpha {\delta }_c+C_D^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+C_D^0,\\ C_L=C_L^{\alpha }\alpha +C_L^{{\delta }_e}{\delta }_e+C_L^{{\delta }_c}{\delta }_c+C_L^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_L^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+C_L^0,\\ C_M=C_M^{\alpha }\alpha +C_M^{{\delta }_e}{\delta }_e+C_M^{{\delta }_c}{\delta }_c+C_M^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_M^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+C_M^0,\\ C_{N_i}=C_{N_i}^{\alpha }\alpha +C_{N_i}^{{\delta }_e}{\delta }_e+C_{N_i}^{{\delta }_c}{\delta }_c+C_{N_i}^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_{N_i}^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+C_{N_i}^0.\end{array}$

8

As shown in equations (7) and (8), the interaction between rigid body and flexible dynamics occurs through the aerodynamic forces. Besides, $\Delta {\mathit{\tau }}_1$ , $\Delta {\mathit{\tau }}_2$ , ${\delta }_e,{\delta }_c$ affect the aerodynamic forces/moments and the generalized elastic forces.

Control-oriented model

In the earlier subsection, the CFM is described and it is only used for simulation. From the study by Fiorentini et al., ³¹ a simplified model, called COM, is proposed for controller design. It is obtained from the CFM based on the neglect of the flexible dynamics for the measurements of the flexible states that are not assumed to be available for feedback. ³¹ Although the flexible dynamic effects are not taken into account directly during the controller design procedure, the flexible dynamics effects are regarded as the perturbations of COM and tested in the simulation. Then, the COM is constituted of five rigid body equations (1) to (5), and all the dominant features of the CFM are remained, which include flexible effects ( $\Delta {\mathit{\tau }}_1$ and $\Delta {\mathit{\tau }}_2$ ) and the couplings between propulsive and aerodynamic forces.

Here, all the uncertain aerodynamic parameters except $C_L^{{\delta }_e},C_L^{{\delta }_c}$ are considered. Besides, the uncertainty of interconnect gain k_ec between canard of deflection δ_c and elevator deflection δ_e is taken into account. That is to say, it is assumed that $k_{ec}=k_{ec0}\text{+}\Delta k_{ec}$ is uncertain, and $\Delta k_{ec}{\delta }_e$ is an uncertain term, which has a bounded magnitude and satisfies an integral quadratic constraint condition. ²⁸ The uncertainty of aerodynamic parameter $C_L^{\alpha }$ is denoted as $\Delta C_L^{\alpha }$ . Then, the expressions of the aerodynamic parameters can be obtained and they are written as follows

$C_L=C_L^{\alpha }\alpha +C_L^0+C_L^{{\delta }_c}\Delta k_{ec}{\delta }_e+C_L^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_L^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+\Delta C_L$

9

$C_D=C_D^{{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2}{\alpha }^2+C_D^{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}\alpha +C_D^0+C_D^{{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2}\left(2\alpha \Delta {\mathit{\tau }}_1+\Delta {\mathit{\tau }}_1{}^2\right)+C_D^{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}\Delta {\mathit{\tau }}_1+C_D^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+\Delta C_d$

10

$C_M=C_M^{\alpha }\alpha +\left(C_M^{{\delta }_e}+C_M^{{\delta }_c}{k}_{ec0}+C_M^{{\delta }_c}\Delta k_{ec}\right){\delta }_e+C_M^0+C_M^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_M^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2$

11

$C_{N_i}=C_{N_i}^{\alpha }\alpha +\left(C_{N_i}^{{\delta }_e}+C_{N_i}^{{\delta }_c}{k}_{ec0}+C_{N_i}^{{\delta }_c}\Delta k_{ec}\right){\delta }_e+C_{N_i}^0+C_{N_i}^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_{N_i}^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2$

12

with

$\Delta C_L=\Delta C_L^{\alpha }\alpha +\Delta C_L^0+\Delta C_L^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+\Delta C_L^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+C_L^{{\delta }_c}\Delta k_{ec}{\delta }_e$

13

$\Delta C_D=\Delta C_D^{{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2}{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2+\Delta C_D^{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}\left(\alpha +\Delta {\mathit{\tau }}_1\right)+\Delta C_D^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+\Delta C_D^0+\Delta C_d$

14

$\Delta C_d=\left[C_D^{{\delta }_e^2}+{\left(k_{ec0}+\Delta k_{ec}\right)}^2{C}_D^{{\delta }_c^2}\right]{\delta }_e^2+\left[C_D^{{\delta }_e}+C_D^{{\delta }_c}\left(k_{ec0}+\Delta k_{ec}\right)\right]{\delta }_e+\left[C_D^{\alpha {\delta }_e}+C_D^{\alpha {\delta }_c}\left(k_{ec0}+\Delta k_{ec}\right)\right]\alpha {\delta }_e$

where

$\begin{array}{l}\Delta C_d=\left[\Delta C_D^{{\delta }_e^2}+{\left(k_{ec0}+\Delta k_{ec}\right)}^2\Delta C_D^{{\delta }_c^2}\right]{\delta }_e^2+\left[\Delta C_D^{{\delta }_e}+\Delta C_D^{{\delta }_c}\left(k_{ec0}+\Delta k_{ec}\right)\right]{\delta }_e+\left[\Delta C_D^{\alpha {\delta }_e}+\Delta C_D^{\alpha {\delta }_c}\left(k_{ec0}+\Delta k_{ec}\right)\right]\alpha {\delta }_e,\\ \Delta C_M=\Delta C_M^{\alpha }\alpha +\left[\Delta C_M^{{\delta }_e}+\Delta C_M^{{\delta }_c}\left(k_{ec0}+\Delta k_{ec}\right)\right]{\delta }_e+\Delta C_M^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+\Delta C_M^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2+\Delta C_M^0,\\ \Delta C_{T,\varphi }=\Delta C_{T,\varphi }^{\alpha }\alpha +\Delta C_{T,\varphi }^{\alpha M_{\infty }^{-2}}\alpha M_{\infty }^{-2}+\Delta C_{T,\varphi }^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+\Delta C_{T,\varphi }^0+\Delta C_{T,\varphi }^{\alpha \Delta {\mathit{\tau }}_1}\alpha \Delta {\mathit{\tau }}_1+\Delta C_{T,\varphi }^{\Delta {\mathit{\tau }}_1{}^2}\Delta {\mathit{\tau }}_1{}^2+\Delta C_{T,\varphi }^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1,\\ \Delta C_T=\Delta C_T^{A_d}{A}_d+\Delta C_T^{\alpha }\alpha +\Delta C_T^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+\Delta C_T^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+\Delta C_T^0,\\ \Delta C_{N_i}=\Delta C_{N_i}^{\alpha }\alpha +\left(\Delta C_{N_i}^{{\delta }_e}+\Delta C_{N_i}^{{\delta }_c}{k}_{ec0}+\Delta C_{N_i}^{{\delta }_c}\Delta k_{ec}\right){\delta }_e+\Delta C_{N_i}^0+\Delta C_{N_i}^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+\Delta C_{N_i}^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2,i=1,2,3.\\ \Delta L\approx \overline{q}S\Delta C_L,\Delta D\approx \overline{q}S\Delta C_D,\Delta M_{yy}\approx z_T\Delta T+\overline{q}S\overline{c}\Delta C_M,\Delta T=\overline{q}\left(\varphi \Delta C_{T,\varphi }+\Delta C_T\right),\Delta N_i\approx \overline{q}\Delta C_{N_i}.\end{array}$

In equations (1) to (5), four control inputs are diffuser area ratio A_d (it is fitted as A_d = 1 in this study), canard deflection δ_c, fuel equivalence ratio (FER) ϕ, and elevator deflection δ_e. The outputs are chosen as velocity V and altitude h. And input constraints on FER and elevator deflection can be expressed as

$\text{sat}\left(\varphi \right)=\{\begin{array}{c}{\varphi }_{\text{max}},\begin{array}{c}{\varphi }_c\ge {\varphi }_{\text{max}}\end{array}\\ {\varphi }_c,\begin{array}{c}{\varphi }_{\text{min}}<{\varphi }_c<{\varphi }_{\text{max}}\end{array}\\ {\varphi }_{\text{min}},\begin{array}{c}{\varphi }_c\le {\varphi }_{\text{min}}\end{array}\end{array},\text{\hspace{0.17em}\hspace{0.17em}sat}\left({\delta }_e\right)=\{\begin{array}{c}{\delta }_e{}_{\text{max}},\begin{array}{c}{\delta }_{ec}\ge {\delta }_e{}_{\text{max}}\end{array}\\ {\delta }_{ec},\begin{array}{c}{\delta }_e{}_{\text{min}}<{\delta }_{ec}<{\delta }_e{}_{\text{max}}\end{array}\\ {\delta }_e{}_{\text{min}},\begin{array}{c}{\delta }_{ec}\le {\delta }_e{}_{\text{min}}\end{array}\end{array}$

15

where ϕ_c and δ_ec are two desired control inputs, which will be explained in the following section. ϕ_min and δ_{e min} are, respectively, the minimum values of FER and elevator deflection, ϕ_max and δ_{e max} are, respectively, the maximum values of FER and elevator deflection.

It is obvious from aerodynamic formulation (7) that thrust T affects velocity V and FER ϕ, and the velocity is mainly affected by FER. Meanwhile, the change of altitude h is dominantly influenced by the elevator deflection. Therefore, it is reasonable to decompose the equations of the COM into velocity subsystem and altitude subsystem. The former is composed of dynamic of velocity and the latter is comprised of the dynamics of altitude, FPA, AOA, and pitch rate. All the states of the rigid body system are assumed to be available, and the controller design only utilizes the feedback of these states. Because the measurements of the flexible states are not assumed to be available for feedback, ³⁰ the flexible states are treated as disturbances.

Based on the aforementioned discussion, equation (1) is expressed as

$\dot{V}=f_V+g_V{\varphi }_c+g_V\Delta \varphi +\Delta f_V$

16

where $\Delta \varphi =\text{sat}\left(\varphi \right)-{\varphi }_c$ . And the formulations of other functions in equation (16) are as follows

$f_V=\overline{q}{C}_{T,1}\cos \alpha /m-g\text{sin}\gamma -\overline{q}S\left(C_D^{{(\alpha +\Delta {\mathit{\tau }}_1)}^2}{\alpha }^2+C_D^{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}\alpha +C_D^0\right)/m$

$\begin{array}{l}\Delta f_V=\left[\Delta T+\overline{q}\varphi \left(C_{T,\varphi }^{\alpha \Delta {\mathit{\tau }}_1}\alpha \Delta {\mathit{\tau }}_1+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1{}^2}\Delta {\mathit{\tau }}_1{}^2+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1\right)+\overline{q}{C}_T^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1\right]\cos \alpha /m-\left[D-\overline{q}S\left(C_D^{{\left(\alpha +\Delta {\mathit{\tau }}_1\right)}^2}{\alpha }^2+C_D^{\alpha +\Delta {\mathit{\tau }}_1}\alpha +C_D^0\right)+\Delta D\right]/m,\\ g_V\text{=}\overline{q}{C}_{T,\varphi ,1}\cos \alpha /m,\end{array}$

with $C_{T,\varphi ,1}=C_{T,\varphi }^{\alpha }\alpha +C_{T,\varphi }^{\alpha M_{\infty }^{-2}}\alpha M_{\infty }^{-2}+C_{T,\varphi }^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+C_{T,\varphi }^0,C_{T,1}=C_T^{A_d}{A}_d+C_T^{\alpha }\alpha +C_T^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+C_T^0.$

By the same token, the dynamics of FPA, AOA, and pitch rate are rewritten as follows

$\dot{h}=V\text{sin}\gamma$

17

$\dot{\gamma }=f_{\gamma }+g_{\gamma }\alpha +\Delta f_{\gamma }$

18

$\dot{\alpha }=f_{\alpha }+g_{\alpha }q+\Delta f_{\alpha }$

19

$\dot{q}=f_q+g_q{\delta }_e{}_c+g_q\Delta {\delta }_e+\Delta f_q$

20

where $\Delta {\delta }_e=\text{sat}\left({\delta }_e\right)-{\delta }_e{}_c$ . And other functions in equations (18) to (20) are as follows.

$\begin{array}{l}{f}_{\gamma }\left(\gamma ,\alpha \right)=\overline{q}\varphi \text{sin}\alpha \left[\left(C_{T,\varphi }^{M_{\infty }^{-2}}{M}_{\infty }^{-2}+C_{T,\varphi }^0\right)\right]/\left(mV\right)-g\cos \gamma /\left(mV\right)+\overline{q}\text{sin}\alpha \left[\left(C_T^{A_d}{A}_d+C_T^{\alpha }\alpha +C_T^0+C_T^{M_{\infty }^{-2}}{M}_{\infty }^{-2}\right)+\overline{q}S{C}_L^0\right]/\left(mV\right),\\ \Delta f_{\gamma }\left(\gamma ,\alpha \right)=\left[\Delta T+\overline{q}\varphi \left(C_{T,\varphi }^{\alpha \Delta {\mathit{\tau }}_1}\alpha \Delta {\mathit{\tau }}_1+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1^2}\Delta {\mathit{\tau }}_1^2+C_{T,\varphi }^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1\right)+\overline{q}{C}_T^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1\right]\text{sin}\alpha /\left(mV\right)\text{+\hspace{0.05em}}\Delta L/\left(mV\right)+\overline{q}S\left(C_L^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_L^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2\right)/\left(mV\right),\\ g_{\gamma }\left(\alpha \right)=\left[\overline{q}S{C}_L^{\alpha }+\overline{q}\varphi \text{sin}\alpha \left(C_{T,\varphi }^{\alpha }+C_{T,\varphi }^{\alpha M_{\infty }^{-2}}{M}_{\infty }^{-2}\right)+\overline{q}\text{sin}\alpha C_T^{\alpha }\right]/\left(mV\right),\\ f_{\alpha }\left(\gamma ,\alpha \right)=-f_{\gamma }\left(\gamma ,\alpha \right),\Delta f_{\alpha }\left(\gamma ,\alpha \right)=-\Delta f_{\gamma }\left(\gamma ,\alpha \right),g_{\alpha }=1,T_1=\overline{q}\left(\varphi C_{T,\varphi ,1}+C_{T,1}\right),\\ f_q\left(\gamma ,\alpha ,q\right)=\left[z_T{T}_1+\overline{q}S\overline{c}\left(C_M^{\alpha }\alpha +C_M^0\right)\right]/I_{yy},g_q\left(\gamma ,\alpha ,q\right)=\overline{q}S\overline{c}\left(C_M^{{\delta }_e}+C_M^{{\delta }_c}{k}_{ec0}\right)/I_{yy},\\ \Delta f_q\left(\gamma ,\alpha ,q\right)=\left[z_T\Delta T+\overline{q}S\overline{c}\left(\Delta C_M+C_M^{\Delta {\mathit{\tau }}_1}\Delta {\mathit{\tau }}_1+C_M^{\Delta {\mathit{\tau }}_2}\Delta {\mathit{\tau }}_2\right)\right]/I_{yy}.\end{array}$

In equations (16) to (20), the terms $\Delta f_V$ , $\Delta f_{\gamma }$ , $\Delta f_{\alpha }$ , $\Delta f_q$ , $g_V\Delta \varphi ,g_q\Delta {\delta }_e$ are unknown, and they are handled through the DOB in the next section.

Based on earlier discussion, the control objective is to design FER and elevator deflection to achieve the stable tracking of velocity and altitude reference commands even with uncertain aerodynamic parameters and input constraints. The design procedure is shown in the following section.

DOB–DSC design

In this section, we will use the DSC approach and nonlinear DOB technique to design controllers for the velocity subsystem (16) and altitude subsystem (18) to (20). The detailed design procedure is described as follows.

Without loss of generality, the uncertainty in the COM is not too big from the background of hypersonic vehicle, so it is reasonable to make the following assumption 1 on the uncertain terms.

Assumption 1

For the uncertain terms $\Delta f_V+g_V\Delta \varphi ,\Delta f_{\gamma },\Delta f_q+g_q\Delta {\delta }_e$ , there are unknown positive constants ${\mu }_V,{\mu }_{\gamma },{\mu }_q$ such that $\left|\Delta f_V+g_V\Delta {\varphi }_V\right|\le {\mu }_V,\left|\Delta f_{\gamma }\right|\le {\mu }_{\gamma },\left|\Delta f_q+g_q\Delta {\delta }_e{}_q\right|\le {\mu }_q$ .

Velocity subsystem control design

In this subsection, the controller for velocity subsystem (16) is designed by combining DOB technique and robust control method. The DOB is firstly designed to handle lumped disturbance and then the controller is designed based on the designed DOB.

To systematically account for the input constraints during the controller design procedure, we will henceforth benefit from the method presented in the study by Zhou and Chen, ³² where the RBFNNs are used to approximate $\Delta \varphi =\text{sat}\left(\varphi \right)-{\varphi }_c$ . And in equation (16), the term $\Delta \varphi =\text{sat}\left(\varphi \right)-{\varphi }_c$ is considered as an uncertainty. So, in equation (16), the term $\Delta f_{V1}=\Delta f_V+g_V\Delta {\varphi }_V$ is the compounded uncertainty induced by aerodynamic uncertainty, the effects of flexible modes and input constraint. Furthermore, inspired by Sira-Ramírez et al., ³³ the following DOB is designed to estimate $\Delta f_{V1}$

$\begin{array}{l}\dot{\widehat{V}}=f_V+g_V{\varphi }_c+\Delta {\widehat{f}}_{V1}+{\mathit{\lambda }}_{13}\left(V-\widehat{V}\right)\\ \Delta {\dot{\widehat{f}}}_{V1}=z_{12}+{\mathit{\lambda }}_{12}\left(V-\widehat{V}\right)\\ {\dot{z}}_{12}=z_{13}+{\mathit{\lambda }}_{11}\left(V-\widehat{V}\right)\\ {\dot{z}}_{13}={\mathit{\lambda }}_{10}\left(V-\widehat{V}\right)\end{array}$

21

Define $e_V=\Delta f_{V1}-\Delta {\widehat{f}}_{V1}$ as the estimation error, and it converges to a small region around zero if the parameters ${\mathit{\lambda }}_{10},{\mathit{\lambda }}_{11},{\mathit{\lambda }}_{12},{\mathit{\lambda }}_{13}$ are chosen appropriately to make the coefficients of equation $e_{V1}^{\left(4\right)}\text{+}{\mathit{\lambda }}_{13}{e}_{V1}^{\left(3\right)}\text{+}{\mathit{\lambda }}_{12}{e}_{V1}^{\left(2\right)}\text{+}{\mathit{\lambda }}_{11}{e}_{V1}^{\left(1\right)}\text{+}{\mathit{\lambda }}_{10}\text{=}0$ satisfy the Hurwitz stability. ³³ And this small region is expressed as

$\left|e_V\right|\le {\epsilon }_V$

22

Here, ε_V > 0 is the random small parameter.

The velocity tracking error is defined as z_V = V − V_d, where V_d is the velocity reference command from equation (16), the velocity tracking error dynamics is

${\dot{z}}_V=f_V+g_V{\varphi }_c+\Delta f_{V1}-{\dot{V}}_d$

23

Based on equations (21) and (23), the controller of the velocity subsystem is designed as

${\varphi }_c=\left[-k_V{z}_V-f_V-\widehat{\Delta }{f}_{V1}+{\dot{V}}_d\right]/g_V$

24

Here, we adopt generalized proportional integral (GPI) observer which was proposed in the study by Sira-Ramírez et al. ³³ to estimate $\Delta f_{V1}$ If the disturbance $\Delta f_{V1}$ is estimated by the Extended state observer (ESO) in the study by Han, ³⁴ then the ESO is

$\{\begin{array}{l}{e}_{11}=\widehat{V}-V\\ \dot{\widehat{V}}=\Delta {\widehat{f}}_{V1}+f_V+g_V{\varphi }_c-{\beta }_{11}{e}_{11}\\ \Delta {\dot{\widehat{f}}}_{V1}=-{\beta }_{12}fal\left(e_{11},{\mathit{\lambda }}_1,{\epsilon }_1\right)\end{array}$

25

where $e_{11},\widehat{V},\Delta {\widehat{f}}_{V1}$ are the estimation errors of state V and outputs of observer, respectively, β₁₁,β₁₂ are the observer gains. The expression of function fal is

$fal\left(e_{11},{\mathit{\lambda }}_1,{\epsilon }_1\right)=\{\begin{array}{l}{\left|e_{11}\right|}^{{\mathit{\lambda }}_1}\text{sgn}\left(e_{11}\right),\left|e_{11}\right|>{\epsilon }_1\\ e_{11}/{\epsilon }_1^{1-{\mathit{\lambda }}_1},other\end{array}$

26

The main differences between ESO and GPI observer are as follows.

Remark 1

ESO is a state and disturbance observer, which has linear and nonlinear forms. For ESO, it only has a first-order extension, that is, one integrator to estimate the disturbance. ³⁴ But GPI observer is characterized by a higher dimensional extension, ³⁵ and thereby has a better estimation performance of disturbances than ESO. Besides, GPI observers are most naturally applicable to the control of perturbed differentially flat nonlinear systems with measurable flat outputs. ³⁶

If the disturbance is dealt with by the DOB presented in the study by Xu, ¹⁴ then the formulation of DOB is

$\{\begin{array}{l}\Delta {\widehat{f}}_{V1}={\widehat{z}}_1+k_{v1}{e}_v\\ {\dot{\widehat{z}}}_1=-k_{v1}\left(f_V+g_V{\varphi }_c+\Delta {\widehat{f}}_{V1}+k_{v1}{e}_v-{\dot{V}}_d\right)\end{array}$

27

with ϕ_c is the designed control input and e_v is the tacking error of velocity.

The comparison between the DOB of this article and that of the study by Xu ¹⁴ is shown as follows.

Remark 2

(1) Both two disturbance observers (21) and (27) can effectively estimate the unknown compounded disturbance where bounded approximation error is achieved. And the approximation error can be adjusted by tuning the observer gains. (2) The observer (21) is not only a state observer but also a disturbance observer, while (27) is only a disturbance observer. (3) The observer (27) is a first-order extension, but (21) is higher order extension, and thereby it has a better estimation performance to some extent. (27) has a simpler structure than that of (21). So it may have more application ability. Moreover, (27) is influenced by the designed input control ϕ_c more deeply than that of (21), since (27) is not only related to input control ϕ_c but also related to tracking error e_v.

To sequel, the stability of the velocity tracking error is analyzed.

Taking the tracking error of velocity into account, the Lyapunov function is chosen as

$Y_V=\frac{1}{2}{z}_V^2$

28

Based on equations (23) and (24), the derivative of Y_V with respect to time yields

$\begin{array}{l}{\dot{Y}}_V=z_V\left(f_V+g_V{\varphi }_c+\Delta f_{V1}-{\dot{V}}_d\right)\\ =-k_V{z}_V^2+z_V{e}_V\end{array}$

29

From equation (22), equality (29) becomes

${\dot{Y}}_V\le -\left(k_V-\frac{1}{2}\right)z_V^2+{\epsilon }_V^2=-2{\epsilon }_1{Y}_V+c_1$

30

where ${\epsilon }_1=\frac{2k_V-1}{4},c_1^{}={\epsilon }_V^2$ .

It can be obtained from (30) that the signals z_V,e_V are uniformly bounded. Standard arguments can be used to solve the inequality of inequality (30) as

$0\le Y_V\left(t\right)\le \frac{c_1}{2{\epsilon }_1}+\left(Y_V\left(0\right)-\frac{c_1}{2{\epsilon }_1}\right)\exp \left(-2{\epsilon }_{1}t\right),\forall t\ge 0$

31

From equations (28) and (31), we can obtain the following equation

$\left|z_V\right|\le \sqrt{\frac{c_1}{2{\epsilon }_1}+Y_V\left(0\right)}$

32

Altitude subsystem control design

The altitude subsystem is comprised of dynamics of h,γ,α,q. According to DSC technique, the design procedure starts from equation (18), and it consists of four steps. The brief design procedure is given as follows. In the first three steps, the virtual control inputs (here, states γ,α,q are regarded as virtual control inputs and denoted as u_h,u_γ,u_α, respectively) are developed. The control input elevator deflection δ_e is designed at last of this section. The detailed design procedure is shown as follows.

Step 1: Define the tracking error of altitude as s_h = h − h_d, where h_d is the reference command of altitude. From equation (2), the derivative of s_h with respect to time is

${\dot{s}}_h=V\text{sin\hspace{0.05em}\hspace{0.05em}\hspace{0.05em}\hspace{0.05em}\hspace{0.05em}\hspace{0.05em}}\gamma -{\dot{h}}_d$

33

Since in the cruise phase, FPA is very small, it is reasonable to consider that sinγ≈γ, and then the following virtual control input u_h is designed as follows

$u_h=\left(-k_h{s}_h\text{+}{\dot{h}}_d\right)/V$

34

where k_h > 0 is the parameter.

During the virtual control input u_γ design procedure, ${\dot{u}}_h$ is utilized. But ${\dot{u}}_h$ is difficult to be computed, since there is aerodynamic uncertainty in ${\dot{z}}_h$ . So, the following novel robust first-order filter is designed to estimate ${\dot{u}}_h$

${\dot{u}}_{h1}=-\left(u_h{}_1-u_h\right)/{\mathit{\tau }}_h-l_h\text{tanh}\left({\xi }_h\left(u_h{}_1-u_h\right)\right)$

35

Here, τ_h is the filter time parameter, and ${\xi }_h,l_h>0$ are constants. And the filter estimation error is defined as $e_{fh}=u_h{}_1-u_h$ .

The filter (35) is used to estimate ${\dot{u}}_h$ , and the estimation of ${\dot{u}}_h$ is applied for the controller design in the next step. It should be pointed out that if l_h is assumed to be zero, then the measurement noise can be reduced in the virtual control effort under the influence of this filter. It is known that the hyperbolic tangent function y = tanh(x) is smoother than the saturation function y = sat(x). Therefore, compared with the classical first-order filter in the DSC method, the fast transient response of each filter can be obtained.

Step 2: Define the first error surface as $s_{\gamma }=\gamma -u_{h1}$ , the derivative of s_γ with respect to time is

${\dot{s}}_{\gamma }=f_{\gamma }\left(\gamma ,\alpha \right)+g_{\gamma }\left(\alpha \right)\alpha +\Delta f_{\gamma }\left(\gamma ,\alpha \right)-{\dot{u}}_{h1}$

36

Then, the uncertain term Δf_γ is estimated by the following DOB

$\begin{array}{l}\dot{\widehat{\gamma }}=f_{\gamma }+g_{\gamma }\alpha +\Delta {\widehat{f}}_{\gamma }-{\dot{u}}_{h1}+{\mathit{\lambda }}_{23}\left(\gamma -\widehat{\gamma }\right)\\ \Delta {\dot{\widehat{f}}}_{\gamma }=z_{22}+{\mathit{\lambda }}_{22}\left(\gamma -\widehat{\gamma }\right)\\ {\dot{z}}_{22}=z_{23}+{\mathit{\lambda }}_{21}\left(\gamma -\widehat{\gamma }\right)\\ {\dot{z}}_{23}={\mathit{\lambda }}_{20}\left(\gamma -\widehat{\gamma }\right)\end{array}$

37

where $\Delta {\widehat{f}}_{\gamma }$ is the estimation of Δf_γ, and the estimation error is $e_{\gamma }=\Delta f_{\gamma }-\Delta {\widehat{f}}_{\gamma }$ . e_γ can converge to a small region around zero if the parameters ${\mathit{\lambda }}_{20},{\mathit{\lambda }}_{21},{\mathit{\lambda }}_{22},{\mathit{\lambda }}_{23}$ are chosen appropriately to make the coefficients of equation $e_{\gamma }^{\left(4\right)}\text{+}{\mathit{\lambda }}_{23}{e}_{\gamma }^{\left(3\right)}\text{+}{\mathit{\lambda }}_{22}{e}_{\gamma }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{21}{e}_{\gamma }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{20}\text{=}0$ satisfy the Hurwitz stability. And the small region is

$\left|e_{\gamma }\right|\le {\epsilon }_{\gamma }$

38

Here, ε_γ > 0 is the random small parameter.

The following virtual control input u_γ is designed as follows

$u_{\gamma }=\left(-k_{\gamma }{s}_{\gamma }-f_{\gamma }-\Delta {\widehat{f}}_{\gamma }+{\dot{u}}_{h1}\right)/g_{\gamma }$

39

where k_γ > 0 is the parameters to be designed.

A new variable $u_{\gamma 1}$ is obtained through the following first-order filter

${\dot{u}}_{\gamma 1}=-\left(u_{\gamma }{}_1-u_{\gamma }\right)/{\mathit{\tau }}_{\gamma }-l_{\gamma }\text{tanh}\left({\xi }_{\gamma }\left(u_{\gamma }{}_1-u_{\gamma }\right)\right)$

40

where τ_γ is the filter time constant, while l_γ,ξ_γ is the positive constant. And the filter estimation error is defined as $e_{f\gamma }=u_{\gamma 1}-u_{\gamma }$ .

• Step 3: The second error surface is defined as $s_{\alpha }=\alpha -u_{\gamma 1}$ , and the time derivative of s_α is

${\dot{s}}_{\alpha }=f_{\alpha }\left(\gamma ,\alpha ,q\right)+q+\Delta f_{\alpha }\left(\gamma ,\alpha ,q\right)-{\dot{u}}_{\gamma 1}$

41

In the above equation, the uncertain term Δf_α is estimated by the following DOB

$\begin{array}{l}\dot{\widehat{\alpha }}=f_{\alpha }+q+\Delta {\widehat{f}}_{\alpha }-{\dot{u}}_{\gamma 1}+{\mathit{\lambda }}_{33}\left(\alpha -\widehat{\alpha }\right)\\ \Delta {\dot{\widehat{f}}}_{\alpha }=z_{32}+{\mathit{\lambda }}_{32}\left(\alpha -\widehat{\alpha }\right)\\ {\dot{z}}_{32}=z_{33}+{\mathit{\lambda }}_{31}\left(\alpha -\widehat{\alpha }\right)\\ {\dot{z}}_{33}={\mathit{\lambda }}_{30}\left(\alpha -\widehat{\alpha }\right)\end{array}$

42

where $\Delta {\widehat{f}}_{\alpha }$ is the estimation of Δf_α, and the estimation error is defined as $e_{\alpha }=\Delta f_{\alpha }-\Delta {\widehat{f}}_{\alpha }$ . e_α can converge to a small region around zero if the parameters ${\mathit{\lambda }}_{30},{\mathit{\lambda }}_{31},{\mathit{\lambda }}_{32},{\mathit{\lambda }}_{33}$ are chosen appropriately to make the coefficients of equation $e_{\alpha }^{\left(4\right)}\text{+}{\mathit{\lambda }}_{33}{e}_{\alpha }^{\left(3\right)}\text{+}{\mathit{\lambda }}_{32}{e}_{\alpha }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{31}{e}_{\alpha }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{30}\text{=}0$ satisfy the Hurwitz stability. And the small region is

$\left|e_{\alpha }\right|\le {\epsilon }_{\alpha }$

43

Here, ε_α > 0 is the random small parameter.

Note that $f_{\gamma }=-f_{\alpha },\Delta f_{\alpha }=-\Delta f_{\gamma }$ , so we estimate Δf_α by DOB. Based on DOB (37) and the first-order filter (40), the virtual control input is designed as

$u_{\alpha }=-k_{\alpha }{s}_{\alpha }-f_{\alpha }+\Delta {\widehat{f}}_{\gamma }+{\dot{u}}_{\gamma 1}-g_{\gamma }{s}_{\gamma }$

44

A new variable u_α1 is obtained through the following first-order filter

${\dot{u}}_{\alpha 1}=-\left(u_{\alpha }{}_1-u_{\alpha }\right)/{\mathit{\tau }}_{\alpha }-l_{\alpha }\text{tanh}\left({\xi }_{\alpha }\left(u_{\alpha }{}_1-u_{\alpha }\right)\right)$

45

where τ_q is the filter time constant, while l_q,ξ_q is the positive constant. And the filter estimation error is defined as $e_{f\alpha }=u_{\alpha 1}-u_{\alpha }$ .

• Step 4: Define the third error surface as $s_q=q-u_{\alpha 1}$ , the time derivative of it along (20) is

${\dot{s}}_q=f_q+g_q{\delta }_{ec}+g_q\Delta {\delta }_e+\Delta f_q-{\dot{u}}_{\alpha 1}$

46

Then, the uncertain term Δf_α is estimated by the following DOB

$\begin{array}{l}\dot{\widehat{q}}=f_q+g_q{\delta }_{ec}+\Delta {\widehat{f}}_{q1}-{\dot{u}}_{\gamma 1}+{\mathit{\lambda }}_{43}\left(q-\widehat{q}\right)\\ \Delta {\dot{\widehat{f}}}_{q1}=z_{42}+{\mathit{\lambda }}_{42}\left(q-\widehat{q}\right)\\ {\dot{z}}_{42}=z_{43}+{\mathit{\lambda }}_{41}\left(q-\widehat{q}\right)\\ {\dot{z}}_{43}={\mathit{\lambda }}_{40}\left(q-\widehat{q}\right)\end{array}$

47

where $\Delta {\widehat{f}}_{q1}=g_q\Delta {\delta }_e+\Delta f_q$ and $\Delta {\widehat{f}}_{q1}$ is the estimation of $\Delta f_{q1}$ , the estimation error is defined as $e_q=\Delta f_{q1}-\Delta {\widehat{f}}_{q1}$ , e_q can converge to a small region around zero if the parameters ${\mathit{\lambda }}_{40},{\mathit{\lambda }}_{41},{\mathit{\lambda }}_{42},{\mathit{\lambda }}_{43}$ are chosen appropriately to make the coefficients of equation $e_{\alpha }^{\left(4\right)}\text{+}{\mathit{\lambda }}_{43}{e}_{\alpha }^{\left(3\right)}\text{+}{\mathit{\lambda }}_{42}{e}_{\alpha }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{41}{e}_{\alpha }^{\left(2\right)}\text{+}{\mathit{\lambda }}_{40}\text{=}0$ satisfy the Hurwitz stability. And this small region is

$\left|e_q\right|\le {\epsilon }_q$

48

here, ε_q > 0 is the random small parameter.

Based on equation (46) and the DOB (47), elevator deflection is designed as

${\delta }_{ec}=\left[-k_q{s}_q-f_q-\Delta {\widehat{f}}_q+{\dot{u}}_{\alpha 1}-s_{\alpha }\right]/g_q$

49

Remark 3

The differences of this article from the studies by Pan and Yu and Pan et al. ^17,18 are as follows. (1) In the studies by Pan and Yu and Pan et al., ^17,18 composite learning–based adaptive DSC is designed for strict-feedback nonlinear systems with unknown parameters. In this article, the DOB-based DSC is proposed for FAHV with aerodynamic uncertainty and input constraint. (2) In the study by Pan and Yu, ¹⁷ parameter convergence can be guaranteed by the interval excitation (IE) instead of persistent excitation (PE) condition and the control performance can be improved from practical asymptotic stability to practical exponential stability. In the study by Pan et al., ¹⁸ parametric uncertainties and external perturbations can be suppressed by the time-interval integral, whereas the DOB adopted in this article to estimate the lumped uncertainty includes parametric uncertainties and input constraint. Besides, the novel second-order filter is applied to estimate the time derivative of virtual control input to eliminate problem of explosion of terms. And the designed control scheme guarantees that the tracking error converges to arbitrarily small neighborhood around zero.

Compared with the existing results ^19,28,29 that also investigate the controller design of AHV based on DSC technique, this article contains the following different aspects. Firstly, the robustness of the scheme developed herein is shown through aerodynamic uncertainty, whereas the robustness of the designed scheme in the study by Waseem et al. ¹⁹ is evaluated through different fuel levels. In the study by Zong et al., ²⁸ the COM is decomposed into three subsystems, and the canard is introduced. Besides, in the previous studies, ^28,29 the additional system is designed to handle input constraint, whereas it does not need to construct additional system in this article. So, it simplifies the control system design procedure and stability analysis. Moreover, novel first-order filters that are different from the study by Zong et al. ²⁹ are developed to handle explosion of terms problem.

Remark 4

The difference of this article from the study by Xu et al. ¹² is as follows. Two subsystems don’t need to be transformed into linear parameterized form, and the uncertainty is tackled by the nonlinear DOB technique rather than adaptive law. Furthermore, input constraint is handled by taking it as a part of lumped uncertainty in this paper, whereas input constraint is handled by the application with compensation technique in the study by Xu et al. ¹² So, the control system design procedure and stability analysis are simplified to some extent.

Based on uniform approximation principle, neural network (NN) can approximate the randomly given smooth nonlinear function with random error accuracy. It has been widely employed in control design of uncertain nonlinear systems owing to its excellent approximation capabilities under certain conditions. However, for NN-based control, it is difficult to choose optimal weights for NNs, and the controlled performance relies on the online learning of NN. ³⁷ Moreover, the basis functions are generally not orthogonal or redundant, that is, the network representation is not unique and is probably not the most efficient one and the convergence of NNs may not be guaranteed. The training procedure may still be trapped in some local minima depending on the initial settings. ³⁸ It increases the complexity of the flight control system, and thereby it is inconvenient for engineering application. The main differences of DOBC approach from other robust control schemes can be concluded as following two aspects. ³⁹

Remark 5

Firstly, it can take the DOB compensation as a “patch” for existing controllers. The DOB-based compensation is added to improve the robustness and disturbance attenuation after the baseline controller is developed according to the existing procedures. So, there is no change to the baseline controller. That is to say, there is no need to employ a complete new and different control strategy that demands a new verification and certification process, and the verification of DOBC can be designed on the basis of the existing verification process to assure safety and reliability. Secondly, DOBC is not a worst-case-based design, while most of the existed robust controls approached are worst-case-based design which has been criticized as being “overconservative.” Promising robustness is achieved with the price of degraded nominal performance. For DOBC approach, the nominal performance of the baseline controller is recovered in the absence of disturbances or uncertainties.

Stability analysis

In what follows, stability analysis for altitude subsystem is presented on the basis of Lyapunov stability. The background of the Lyapunov analysis method is as follows. Stability theory plays a central role in systems theory and engineering. There are different kinds of stability problems that arise in the study of dynamical systems. Stability of equilibrium points is usually characterized in the sense of Lyapunov, a Russian mathematician and engineer, who laid the foundation of the theory, which now carries his name. ⁴⁰ An equilibrium point is stable if all solutions starting at nearby points stay nearby; otherwise, it is unstable. It is asymptotically stable if all solutions starting at nearby points not only stay nearby but also tend to the equilibrium point as time approaches infinity. In framework of Lyapunov analysis, appropriate Lyapunov function for the system is chosen firstly and then the time derivative of Lyapunov function is analyzed. In the studies by Liu et al., ^41,42 the discrete-time Lyapunov function is chosen for the stability analysis of discrete-time system. In the literature, ^{43 –45} for interconnected system, the Lyapunov function including adaptive estimation error and approximation error is chosen to analyze stability. And in the study by Chen et al., ⁴⁶ for nonlinear stochastic nonlinear pure feedback systems, using Lyapunov analysis, it is proven that all the signals of the closed-loop system are semiglobally uniformly ultimately bounded in probability and the system output tracks the reference signal to a bounded compact set. Then, in the study by Liu and Tong, ⁴⁷ for the nonlinear systems with full state constrained, barrier Lyapunov functions are adopted to guarantee that the states are not to violate their constraints.

The following theorem is given to show the stability of the closed-loop system.

Theorem considers the COM (18) to (20), when assumption holds, with application of virtual control inputs (34), (39), (44), and actual control input (49), first-order filters (35), (40), and (45), and nonlinear DOBs (37), (42), and (47). Then, if the design parameters are chosen appropriately, the algorithm can guarantee that the tracking error of altitude converges to an arbitrarily small neighborhood around zero.

Proof

Since DOB and filter are introduced to estimate the lumped uncertainty and time virtual control input, choose the following Lyapunov function

$Y=Y_h+Y_{\gamma }+Y_{\alpha }+Y_q+\frac{1}{2}{e}_{fh}^2+\frac{1}{2}{e}_{f\gamma }^2+\frac{1}{2}{e}_{f\alpha }^2$

50

with $Y_h=\frac{1}{2}{s}_h^2,Y_{\gamma }=\frac{1}{2}{s}_{\gamma }^2,Y_{\alpha }=\frac{1}{2}{s}_{\alpha }^2,Y_q=\frac{1}{2}{s}_q^2$ .

The time derivative of equation (50) is

$\dot{Y}={\dot{Y}}_h+{\dot{Y}}_{\gamma }+{\dot{Y}}_{\alpha }+{\dot{Y}}_q+e_{fh}{\dot{e}}_{fh}+e_{f\gamma }{\dot{e}}_{f\gamma }+e_{f\alpha }{\dot{e}}_{f\alpha }$

51

From equation (35) and filter error, if $l_h>{\left|{\dot{u}}_{\mathit{\lambda }}\right|}_{\text{max}}$ , from the study by Polycarpou and Ioannou, ⁴⁸ the term $e_{fh}{\dot{e}}_{fh}$ yields

$\begin{array}{l}{e}_{fh}{\dot{e}}_{fh}=-\frac{e_{fh}^2}{{\mathit{\tau }}_h}-l_h\text{tanh}\left({\xi }_h\left(e_{fh}\right)\right)e_{fh}-{\dot{u}}_h{e}_{fh}\\ \le -\frac{e_{fh}^2}{{\mathit{\tau }}_h}-\left(l_h\text{tanh}\left({\xi }_h\left(e_{fh}\right)\right)e_{fh}-{\dot{u}}_h\left|e_{fh}\right|\right)\\ \le -\frac{e_{fh}^2}{{\mathit{\tau }}_h}+\frac{k}{{\mathit{\lambda }}_h}\end{array}$

52

By the same token, if $l_i>{\left|{\dot{u}}_{\gamma }\right|}_{\text{max}}\left(i=\gamma ,\alpha \right)$ , $l_{\alpha }>{\left|{\dot{u}}_{\alpha }\right|}_{\text{max}}$ , then $e_{fi}{\dot{e}}_{fi}\le -e_{fi}^2/{\mathit{\tau }}_i+k/{\mathit{\lambda }}_i$ and equation (51) satisfies

$\dot{Y}\le {\dot{Y}}_h+{\dot{Y}}_{\gamma }+{\dot{Y}}_{\alpha }+{\dot{Y}}_q-\frac{e_{fh}^2}{{\mathit{\tau }}_h}-\frac{e_{f\gamma }^2}{{\mathit{\tau }}_{\gamma }}-\frac{e_{f\alpha }^2}{{\mathit{\tau }}_{\alpha }}+\frac{k}{{\mathit{\lambda }}_h}+\frac{k}{{\mathit{\lambda }}_{\gamma }}+\frac{k}{{\mathit{\lambda }}_{\alpha }}$

53

And in equation (53), ${\dot{Y}}_h=s_h{\dot{s}}_h$ , ${\dot{Y}}_{\gamma }=s_{\gamma }{\dot{s}}_{\gamma }$ , ${\dot{Y}}_{\alpha }=s_{\alpha }{\dot{s}}_{\alpha }$ , ${\dot{Y}}_q=s_q{\dot{s}}_q$ , then based on equations (33) and (34), ${\dot{Y}}_h$ yields

${\dot{Y}}_h\le -k_h{s}_h^2$

54

For the equation ${\dot{Y}}_{\gamma }=s_{\gamma }{\dot{s}}_{\gamma }$ , based on equation (36), the equation ${\dot{Y}}_{\gamma }=s_{\gamma }{\dot{s}}_{\gamma }$ becomes

${\dot{Y}}_{\gamma }=s_{\gamma }\left(f_{\gamma }+g_{\gamma }\left(s_{\alpha }\text{+}{e}_{f\gamma }\text{+}{u}_{\gamma }\right)+\Delta f_{\gamma }-{\dot{u}}_{h1}\right)$

55

where $s_{\alpha }=\alpha -u_{\gamma 1}$ is the second error surface.

Substituting equation (39) into (55), we have

${\dot{Y}}_{\gamma }\le -\left(k_{\gamma }-\frac{1}{2}\right)s_{\gamma }^2+\frac{1}{2}{e}_{\gamma }^2+g_{\gamma }{s}_{\gamma }{e}_{f\gamma }+g_{\gamma }{s}_{\gamma }{s}_{\alpha }$

56

For the equation ${\dot{Y}}_{\alpha }=s_{\alpha }{\dot{s}}_{\alpha }$ , based on equation (41), the equation ${\dot{Y}}_{\alpha }=s_{\alpha }{\dot{s}}_{\alpha }$ becomes

${\dot{Y}}_{\alpha }=s_{\alpha }\left(f_{\alpha }+s_q\text{+}{e}_{f\alpha }\text{+}{u}_{\alpha }+\Delta f_{\alpha }-{\dot{u}}_{\alpha 1}\right)$

57

where $s_q=q-u_{\alpha 1}$ is the third error surface.

Substituting equation (44) into it, we have

$\begin{array}{c}{\dot{Y}}_{\alpha }=-k_{\alpha }{s}_{\alpha }^2-s_{\alpha }{e}_{\alpha }+s_{\alpha }{e}_{f\alpha }-g_{\gamma }{s}_{\gamma }{s}_{\alpha }+s_{\alpha }{s}_q\\ \le -\left(k_{\alpha }-\frac{1}{2}\right)s_{\alpha }^2+\frac{1}{2}{e}_{\alpha }^2+\frac{1}{2}{e}_{f\alpha }^2-g_{\gamma }{s}_{\gamma }{s}_{\alpha }+s_{\alpha }{s}_q\end{array}$

58

For the equation ${\dot{Y}}_q=s_q{\dot{s}}_q$ , based on equation (46), it becomes

${\dot{Y}}_q = s_q\left(f_q+g_q{\delta }_{ec}+g_q\Delta {\delta }_e+\Delta f_q-{\dot{u}}_{\alpha 1}\right)$

59

From equations (47) and (49), the above equation yields

$\begin{array}{c}{\dot{Y}}_q = s_q\left(f_q+g_q{\delta }_{ec}+g_q\Delta {\delta }_e+\Delta f_q-{\dot{u}}_{\alpha 1}\right)\\ =-k_q{s}_q^2-s_q{e}_q-s_{\alpha }{s}_q\\ \le -\left(k_q-\frac{1}{2}\right)s_q^2+\frac{1}{2}{e}_q^2-s_{\alpha }{s}_q\end{array}$

60

Based on equations (38), (43), (47), (54), (56), (58), and (60), equation (53) satisfies

$\begin{array}{l}\dot{Y}\le -k_h{s}_h^2-\left(k_{\gamma }-1\right)s_{\gamma }^2-\left(k_{\alpha }-1\right)s_{\alpha }^2-\left(k_q-\frac{1}{2}\right)s_q^2\\ +\frac{1}{2}{e}_{\gamma }^2+\frac{1}{2}{e}_{\alpha }^2+\frac{1}{2}{e}_q^2+\frac{1}{2}{g}_{\gamma }^2{e}_{{}_{f\gamma }}^2+\frac{1}{2}{e}_{f\alpha }^2+\frac{1}{2}{e}_{fq}^2\\ -\frac{e_{fh}^2}{{\mathit{\tau }}_h}-\frac{e_{f\gamma }^2}{{\mathit{\tau }}_{\gamma }}-\frac{e_{f\alpha }^2}{{\mathit{\tau }}_{\alpha }}+k\left(\frac{1}{{\mathit{\lambda }}_h}+\frac{1}{{\mathit{\lambda }}_{\gamma }}+\frac{1}{{\mathit{\lambda }}_{\alpha }}\right)\\ \le -k_h{s}_h^2-\left(k_{\gamma }-1\right)s_{\gamma }^2-\left(k_{\alpha }-1\right)s_{\alpha }^2-\left(k_q-\frac{1}{2}\right)s_q^2\\ -\left(\frac{1}{{\mathit{\tau }}_h}-\frac{1}{2}{g}_{\gamma }^2\right)e_{fh}^2-\left(\frac{1}{{\mathit{\tau }}_{\gamma }}-\frac{1}{2}\right)e_{f\gamma }^2-\left(\frac{1}{{\mathit{\tau }}_{\alpha }}-\frac{1}{2}\right)e_{f\alpha }^2\\ +\frac{{\epsilon }_{\gamma }^2+{\epsilon }_{\alpha }^2+{\epsilon }_q^2}{2}+k\left(\frac{1}{{\mathit{\lambda }}_h}+\frac{1}{{\mathit{\lambda }}_{\gamma }}+\frac{1}{{\mathit{\lambda }}_{\alpha }}\right)\\ \le -2\epsilon Y+c\end{array}$

61

Here, $c=\frac{{\epsilon }_{\gamma }^2+{\epsilon }_{\alpha }^2+{\epsilon }_q^2}{2}+k\left(\frac{1}{{\mathit{\lambda }}_h}+\frac{1}{{\mathit{\lambda }}_{\gamma }}+\frac{1}{{\mathit{\lambda }}_{\alpha }}\right)$ and $\epsilon =\text{min}\left\{k_h,k_{\gamma }-1,k_{\alpha }-1,k_q-\frac{1}{2},\frac{1}{{\mathit{\tau }}_h}-\frac{1}{2},\frac{1}{{\mathit{\tau }}_h}-\frac{1}{2}{g}_{\gamma }^2,\frac{1}{{\mathit{\tau }}_{\alpha }}-\frac{1}{2}\right\}>0$ .

And the convergence domain of s_h is the following compact set

$R_h=\left\{s_h|\left|s_h\right|\le \sqrt{\frac{c}{2\epsilon }+Y\left(0\right)}\right\}$

62

From equation (62), as long as ε increases sufficiently large, the altitude tracking error s_h converges to arbitrarily small neighborhood around zero and the sufficiently large ε can be achieved via properly choosing the design parameters.

Simulation and analysis

Simulations are done to illustrate the control performance of the developed control scheme in this section. Stability of the flexible states can also be evaluated through simulation. The model parameters and the initial flight condition of the aircraft are referred to the study by Sigthorsson and Serrani. ³⁰ The fuel level is assumed at 50%. And 20% uncertainties on aerodynamic parameters are considered. The limits are set as ϕ_min = 0.1, ϕ_max = 1.2, δ_{e min} = −0.2618, δ_{e max} = −δ_{e min} to explore the capability of the designed control scheme in tackling the constraints.

The aggressive maneuver scenario simulation is considered to test the control performance of the designed control strategy. It is a climbing maneuver, where velocity reference command V_d and altitude reference command h_d are chosen independently. The velocity reference command is generated to make the vehicle accelerate from 7746 ft/s to 9746 ft/s, and the altitude reference command h_d is generated to make the aircraft climbs from 85,000 ft to105,000 ft. In order to obtain continuous and bounded derivatives of the desired flight attitudes, two-order reference model $G(s)=\frac{{\omega }_s^2}{s^2+2\xi {\omega }_s^{}s+{\omega }_s^2}$ is employed, where ω_s = 0.03 rad/s and ξ = 0.95 are, respectively, a natural frequency and a damping ratio. The adopted parameters of the control scheme are as follows: $k_V=15,k_h=\text{1, }{k}_{\gamma }=3,k_{\alpha }=8,k_q=35,$ ${\mathit{\lambda }}_{i0}=4,{\mathit{\lambda }}_{i1}=7,{\mathit{\lambda }}_{i2}={\mathit{\lambda }}_{i3}=8,\left(i=1,2,3,4\right),{\mathit{\tau }}_j=0.06,l_j=0.1,{\zeta }_j=0.4\left(j=h,\gamma ,\alpha \right)$ . The simulation result is shown in Figures 1 –4.

Figure 1.

Time response of velocity and altitude.

Figure 2.

Time response of FPA, AOA, pitch rate, and flexible states.

Figure 3.

Time response of control inputs.

Figure 4.

Time response of DOB estimation. DOB: disturbance observer.

As can be obtained from Figure 1 that though there are uncertainty and input constraint, under the designed robust DSC scheme, both velocity V and altitude h stably follow their reference commands V_d and h_d. Besides, the corresponding maximum absolute values of velocity tracking error and altitude tracking error are less than 0.4 ft/s and 4 ft, respectively. Besides, from Figure 2, during the flight phase, FPA γ, AOA α, pitch rate q change smooth, and three flexible states converge to their stable values despite excitations. As shown in Figure 3, all control inputs are in the saturation range under the proposed controller. And the saturations occur during the transient phase of the control inputs. It is clear from Figure 4 that the DOB estimates the uncertainty well.

Conclusions

An effective control strategy is designed for FAHV under aerodynamic uncertainty and input constraints based on DSC and DOB technique. The COM is firstly derived from CFM and is divided into two subsystems, which are velocity subsystem and altitude subsystem. A robust control scheme and the DSC scheme are, respectively, designed for the velocity subsystem and altitude subsystem, where input constraint and uncertainty are handled by the DOB. To handle “explosion of terms” problem, the novel first-order filters are employed. Both the tracking error and the disturbance observer error can converge to random small region around zero. At last, simulation results demonstrate the performance of the presented flight control strategy. In the further work, we will focus on fixed-time control for FAHV with input constraints and part states unmeasured.