Prescribed performance control of morphing aircraft based on switched nonlinear systems and reinforcement learning

Abstract

The prescribed performance control of a morphing aircraft with variable sweep wings is investigated based on switched nonlinear systems and reinforcement learning. Switched nonlinear systems in lower triangular form are first adopted to describe the longitudinal altitude motion, and an error transformation is applied to handle the prescribed performance bound. Then, the designed controller is divided into the basic part and supplementary part. For the basic part, the backstepping method with involvement of the modified dynamic surface control technique is utilized to avoid the “explosion of complexity” problem. Improved disturbance observers inspired from the idea of extended state observer are then designed to estimate the disturbances and combined with radial basis function neural networks to develop the common virtual control laws. Moreover, by using the error variables defined in the backstepping design, the reinforcement learning–based supplementary part controller is devised with the critic-action neural network structure, which can adjust its parameters online and further decrease the altitude tracking error. It is proved that all signals of the closed-loop system are uniformly ultimately bounded, and the prescribed performance bound for convergence of the altitude tracking error can be satisfied. Finally, comparative simulations demonstrate the effectiveness of the proposed control approach.

Introduction

The morphing aircraft can alter its aerodynamic configuration to expand mission capabilities and obtain optimal flight performance.^1–3 The change in aircraft geometrical shape resolves performance conflicts between high-speed and low-speed regimes. However, the large variations in mass distribution and applied aerodynamic forces during the morphing process may cause unstability of the aircraft, together with internal uncertainties and external disturbances. Hence, compared with traditional flight vehicles, the control of morphing aircraft is more challenging and has been explored by the worldwide researchers in recent years.^4,5

The extant research on the morphing aircraft control mainly focuses on the longitudinal dynamics, where several representative control methods have been investigated based on linear parameter varying (LPV) model or general nonlinear systems. An inner loop linear quadratic controller is combined with outer loop gain-scheduled controller in Yue et al. ⁶ to guarantee stability. The linear matrix inequalities and sliding mode–based approaches are investigated in He et al. ⁷ and Wen et al. ⁸ for the LPV model, respectively. A high-order chained differentiator is combined with neural networks (NNs) to determine the control input in Wu et al., ⁹ where a coordinate transformation is adopted to obtain normal form of the aircraft motion model. The adaptive dynamic surface control is proposed in Wu et al. ¹⁰ to deal with input–output constraints via a barrier Lyapunov function, where system uncertainties are approximated by NNs. The constraints on active region of the NNs are further relaxed in Wu et al. ¹¹ by a smooth switching between the neural controller and robust controller. He et al. ¹² adopt NNs to deal with the model uncertainties and design a disturbance observer to attenuate the effects of the disturbances of a flapping wing micro aerial vehicle.

It is noted that the aforementioned LPV model–based approaches may degrade the control performance when the nonlinearities are severe. Meanwhile, the approximation of the variations in mass distribution and aerodynamic forces is difficult by a single NN in Wu et al.^9–11 due to the time-varying characteristics. In contrast, switched systems provide a powerful tool for modeling of morphing aircraft dynamics, where the switching law describes the morphing process, and the subsystems correspond to the specific configurations along the transition trajectory. The longitudinal dynamics of a morphing aircraft with variable sweepback wings is approximated by uncertain switched linear systems in Wang et al., ¹³ where a robust state feedback controller is designed to guarantee finite-time boundedness of the closed-loop system. Switched LPV models in the continuous domain are established in Jiang et al., ¹⁴ and the smooth switching control problem is investigated via overlapped scheduled parameter subsets. The non-fragile $H_{\infty }$ controller is proposed in Cheng et al. ¹⁵ for the switched LPV models in the discrete domain and is expanded to finite-time $H_{\infty }$ control in Cheng et al. ¹⁶ Nevertheless, the methods in previous studies^14–16 require linearization of the nonlinear dynamics, and switched nonlinear systems can be further adopted to enlarge the modeling capability, which contributes to reducing the conservativeness of controller design.

Considerable progress has been made in stability analysis and controller design for switched nonlinear systems during the past decade,^17–21 where the studied control schemes mainly deal with stability of the closed-loop system and do not address the performance of tracking errors. Recently, the prescribed performance control proposed in Bechlioulis and Rovithakis^22–24 and developed for various classes of nonlinear systems in previous studies^25–27 has been extended to switched nonlinear systems. The main idea of prescribed performance control is coordinate transformation, and the convergence of tracking error can then be guaranteed with prescribed decay rate, maximum overshoot and steady value. The prescribed performance control of switched nonlinear systems in non-strict-feedback form is investigated in Li et al., ²⁸ where semi-globally uniform ultimate boundedness is achieved. A prescribed performance adaptive controller is developed in Li and Xiang ²⁹ to deal with the input saturation and unmodelled dynamics of a class of switched nonlinear systems. Mode-dependent adaptive laws are devised in Zhai et al. ³⁰ to reflect the switchings between subsystems, which is less conservative compared with the common adaptive laws. However, to our best knowledge, the related research on control of morphing aircraft via switched nonlinear systems has not been fully investigated yet. Besides, the extended state observer (ESO)-based control approach provides a universal way of dealing with disturbances,^31,32 where the disturbances are regarded as an extended state estimated by the ESO and compensated for by the controller. It is noteworthy that the research on ESO for switched nonlinear systems is also rare.

However, reinforcement learning has been studied to find optimal control solutions for various classes of systems over the last few decades.^33–35 The critic-action structure is commonly utilized in the reinforcement learning, where two agents called critic agent and action agent are adopted to estimate performance index and revise the control policy, respectively. The reinforcement learning control problems of nonlinear systems including partially and fully unknown dynamics are addressed in literature,^36–38 where NNs are adopted to implement the algorithm. The reinforcement learning control has also been applied to robot manipulators, ³⁹ hypersonic vehicles, ⁴⁰ autonomous underwater vehicles ⁴¹ and so on. In particular, it is worth pointing out that Mu et al. ⁴⁰ combine sliding mode control and reinforcement learning to solve the tracking problem of hypersonic vehicles, where a supplementary strategy is generated by the reinforcement learning approach to further reduce tracking errors and improve control performance. In view of the online learning and adaption merit, it is desirable to integrate the reinforcement learning into morphing aircraft controller design.

Motivated by the above considerations, this paper investigates prescribed performance control of morphing aircraft based on switched nonlinear systems and reinforcement learning, where the prescribed performance bound is handled by error transformation. The main contributions of our work are summarized as follows:

1. The longitudinal altitude motion of morphing aircraft is modeled as switched nonlinear systems in lower triangular form, and it provides an efficient way to cope with the variations in mass distribution and applied aerodynamic forces.

2. The designed controller is split into two parts: the basic part and supplementary part. For the basic part, the backstepping approach with integration of the dynamic surface control technique is exploited to eliminate the “explosion of complexity” problem. The disturbance observers inspired from the idea of ESO are designed and combined with the radial basis function (RBF) NNs to generate the common virtual control laws.

3. The critic-action NN structure of the reinforcement learning scheme is adopted to devise the supplementary part on the basis of tracking errors defined in the backstepping design. The weight updating algorithm of the established critic-action NN structure and related convergence analysis are provided. The improved control performance via generating a supplementary strategy is validated by comparative simulation studies.

The remaining part of this study is organized as follows. Section “Model description” states the modeling of morphing aircraft by switched nonlinear systems. Section “Controller design” presents the basic part of the controller based on the switched nonlinear systems and the supplementary part of the controller based on reinforcement learning. Comparative simulation studies are addressed in section “Numerical simulation,” and final conclusion is drawn in section “Conclusion.”

Notations

The notations in this paper are standard. $\mathbb{R}$ and ${\mathbb{R}}^n$ denote the set of real numbers and n-dimensional real vectors, respectively. ${\mathit{I}}_n$ refers to the n-dimensional identity matrix. The maximum and minimum eigenvalues of a matrix $\mathit{P}$ are denoted by ${\lambda }_{\max }(\mathit{P})$ and ${\lambda }_{\max }(\mathit{P})$ , respectively. $\parallel ·\parallel \mathrm{represents}$ the Euclidean vector norm.

Model description

The morphing aircraft with variable sweep wings considered in this study is similar with the research in Wu et al.^9–11 and the longitudinal motion model is formulated as ¹⁰

$$\{\begin{array}{c}\stackrel{·}{h}=V\sin \gamma \hfill \\ \multicolumn{1}{c}{\stackrel{·}{\gamma }=\frac{L+T\sin \alpha -F_{\mathrm{Ikz}}}{\mathrm{mV}}-\frac{g}{V}\cos \gamma }\\ \multicolumn{1}{c}{\stackrel{·}{\alpha }=-\frac{L+T\sin \alpha -F_{\mathrm{Iz}}}{\mathrm{mV}}+\frac{g}{V}\cos \gamma +q}\\ \multicolumn{1}{c}{\stackrel{·}{q}=\frac{-{\stackrel{·}{I}}_{y}q-S_{x}g\cos \theta +M_A+T{Z}_T+M_{\mathrm{Iy}}}{I_y}}\end{array}$$

(1)

where $h$ , $\gamma$ , $\alpha$ , and $q$ represent the altitude, flight path angle, angle of attack, and pitch rate, respectively. $m$ and $I_y$ denote the aircraft mass and inertia, $T$ is the engine thrust, $Z_T$ is the thrust moment arm, $g$ is the acceleration of gravity, $L$ and $M_A$ are the lift and pitch moments, and $F_{\mathrm{Ikz}}$ , $F_{\mathrm{Iz}}$ , and $M_{\mathrm{Iy}}$ are inertial forces and moment related to morphing process which can be modeled as

$$\begin{array}{c}{F}_{\mathrm{Iz}}=F_{\mathrm{Ikz}}=S_x(\stackrel{·}{q}\cos \alpha -q^2\sin \alpha )+2{\stackrel{·}{S}}_{x}q\cos \alpha +{\stackrel{··}{S}}_x\sin \alpha \hfill \\ M_{\mathrm{Iy}}=S_x(\stackrel{·}{V}\sin \alpha +V\stackrel{·}{\alpha }\cos \alpha -\mathrm{Vq}\cos \alpha )s\hfill \end{array}$$

where $S_x$ denotes the static moment distribution in x-axis of the body frame. The related definitions are given as follows

$$\begin{array}{c}\multicolumn{1}{c}{S_x(\xi )\approx 2m_1{r}_{1x}+m_3{r}_{3x}}\\ \multicolumn{1}{c}{L}& \approx L_{\alpha }(\xi )\alpha +L_0(\xi )=\overline{q}{S}_w(\xi )(C_L^{\alpha }(\xi )\alpha +C_L^0(\xi ))\hfill \\ \multicolumn{1}{c}{M_A}& \approx M_{A\alpha }(\xi )\alpha +M_{A0}(\xi )+M_{A{\delta }_e}(\xi ){\delta }_e+M_{\mathrm{Aq}}(\xi )\frac{q{c}_A(\xi )}{2V}\hfill \\ =\overline{q}{S}_w(\xi )c_A(\xi )(C_M^{\alpha }(\xi )\alpha +C_M^0(\xi )+C_M^{{\delta }_e}(\xi ){\delta }_e+\frac{C_M^{q}q{c}_A(\xi )}{2V})\hfill \end{array}$$

where $m_1$ and $m_3$ denote the mass of aircrafts’ wing and fuselage, $r_{1x}$ and $r_{3x}$ denote the position of corresponding components in the body frame, $\overline{q}=(1/2){\rho }_h{V}^2$ denotes the dynamic pressure, ${\rho }_h$ denotes the air density, $\xi$ denotes the sweep angle, $S_w(\xi )$ denotes the wing surface, and ${\delta }_e$ denotes the elevator angle. The other related definitions can be found in Wu et al.^9–11

Assumption 1

The flight path angle $\gamma$ is small, that is, $\sin \gamma \approx \gamma$ .^10,11 From Equation (1), when the sweep angle takes values within the available range, the related altitude dynamics can be written as

$$\{\begin{array}{c}\stackrel{·}{h}=V\gamma \hfill \\ \multicolumn{1}{c}{\stackrel{·}{\gamma }=f_{2,\sigma (t)}+g_{2,\sigma (t)}\theta +d_{\gamma }}\\ \multicolumn{1}{c}{\stackrel{·}{\theta }=q}\\ \multicolumn{1}{c}{\stackrel{·}{q}=f_{4,\sigma (t)}+g_{4,\sigma (t)}{\delta }_e+d_q}\end{array}$$

(2)

where $f_{2,\sigma (t)}=(L_{0,\sigma (t)}-L_{\alpha ,\sigma (t)}\gamma )/\mathrm{mV}-(g/V)\cos \gamma$ , $g_{2,\sigma (t)}=L_{\alpha ,\sigma (t)}/\mathrm{mV}$ , $d_{\gamma }=(T\sin \alpha -F_{\mathrm{Ikz}})/\mathrm{mV}+\Delta d_{\gamma }$ , $f_{4,\sigma (t)}=(T{Z}_T+M_{A\alpha ,\sigma (t)}\alpha +M_{A0,\sigma (t)}/I_{y,\sigma (t)})+(M_{\mathrm{Aq},\sigma (t)}q{c}_{A,\sigma (t)}/2I_{y,\sigma (t)}V)$ , $g_{4,\sigma (t)}=M_{A{\delta }_e,\sigma (t)}/I_{y,\sigma (t)}$ , $d_q=(M_{\mathrm{Iy}}-{\stackrel{·}{I}}_{y}q-S_{x}g\cos \theta /I_y)+\Delta d_q$ , and $\sigma (t):[0,+\infty )\to {\mathit{S}}_M$ represent the switching signal related to the sweep angle $\xi$ , and ${\mathit{S}}_M=\{1,2,\dots ,M\}$ denotes the finite index set of the aircraft configurations for certain sweep angles. $\Delta d_{\gamma }$ and $\Delta d_q$ denote the equivalent disturbances caused by model uncertainties and exogenous perturbations, respectively.

In this paper, the RBF NNs are applied to Equation (2) for obtaining common virtual control laws in the backstepping design procedure. For any continuous function $f(\mathit{Z}):{\mathbb{R}}^n\to \mathbb{R}$ defined in a compact set ${\Omega }_Z\subset {\mathbb{R}}^n$ with $n$ being the input dimension, there exist an RBF NN and a positive constant $\overline{\omega }$ such that

$$f(\mathit{Z})={\mathit{W}}^{\mathrm{T}}\mathit{S}(\mathit{Z})+\omega ,|\omega |⩽\overline{\omega }$$

(3)

where $\mathit{Z}\in {\Omega }_Z\subset {\mathbb{R}}^n$ denotes the NN’s input, $\omega$ is the approximation error, and $\mathit{W}=[w_1,w_2,\dots ,w_{\overline{l}}{]}^{\mathrm{T}}\in {\mathbb{R}}^{\overline{l}}$ denotes the ideal weight vector with $\overline{l}$ being the NN node number, which is defined as

$$\mathit{W}:=\arg \underset{\hat{\mathit{W}}\in {\mathbb{R}}^{\overline{l}}}{\min }\{\underset{\mathit{Z}\in {\Omega }_Z}{sup}|f(\mathit{Z})-{\hat{\mathit{W}}}^{\mathrm{T}}\mathit{S}(\mathit{Z})|\}$$

(4)

where $\mathit{S}(\mathit{Z})={[s_1(\mathit{Z}),s_2(\mathit{Z}),\dots ,s_{\overline{l}}(\mathit{Z})]}^{\mathrm{T}}$ represents the basis function vector in Gaussian function form, that is

$$s_i(\mathit{Z})=\exp [-\frac{{(\mathit{Z}-{\mathbf{\mu }}_i)}^{\mathrm{T}}(\mathit{Z}-{\mathbf{\mu }}_i)}{{\zeta }_i^2}],i=1,2,\dots ,\overline{l}$$

(5)

with $m_i=[{\mu }_{i1},{\mu }_{i2},\dots ,{\mu }_{\mathrm{in}}{]}^{\mathrm{T}}$ and ${\zeta }_i$ denoting the center of the receptive field and the width of the Gaussian function, respectively. Moreover, we review the following lemma for dealing with bounds on estimated values of the weight vector.

Lemma 1

The following inequality is satisfied for any $z\in \mathbb{R}$ and any positive constant $ϵ$ ⁴²

$$0|z|-z\tan h\left(\frac{z}{ϵ}\right)\chi ϵ,\chi =0.2785$$

(6)

Controller design

The designed controller of this paper comprised two parts, that is, ${\delta }_e={\delta }_{e1}+{\delta }_{e2}$ , where ${\delta }_{e1}$ and ${\delta }_{e2}$ denote the basic and supplementary parts, respectively. The structure of the controller is depicted in Figure 1. As will be stated in section “Basic part design,” the basic part is developed by the backstepping method, where improved disturbance observers and dynamic surfaces are devised to make the controller less conservative. Specifically, the coupling terms related to the tracking errors are introduced in the disturbance observers and dynamic surfaces to deal with cross-product terms of the time derivative of the Lyapunov function. Meanwhile, the supplementary part is generated by reinforcement learning algorithm with involvement of the critic-action NN structure to further expedite the control performance.

Figure 1.

Structure diagram of the proposed controller.

Basic part design

The following assumptions are made for the switched nonlinear systems (2).

Assumption 2

The reference altitude $h_r$ and its first two-order time derivatives are continuous and bounded.

Assumption 3

The velocity $V$ and its time derivative are continuous and bounded.

Assumption 4

The time derivatives of disturbances $d_{\gamma }$ and $d_q$ are bounded, that is, there exist positive constants $N_{\gamma }$ and $N_q$ such that $|{\stackrel{·}{d}}_{\gamma }|⩽N_{\gamma }$ and $|{\stackrel{·}{d}}_q|⩽N_q$ .

Assumption 5

For any $i\in {\mathit{S}}_M$ , $g_{2i}$ and $g_{4i}$ have fixed signs and satisfy $g_{2m}⩽|g_{2i}|⩽g_{2M}$ and $g_{4m}⩽|g_{4i}|⩽g_{4M}$ , where $g_{2m}$ , $g_{2M}$ , $g_{4m}$ , and $g_{4M}$ are positive constants.

Remark 1

Assumption 2 is rather general in the literature on aircraft longitudinal motion control.^43,44 Assumption 3 is reasonable since the thrust force, drag force, and inertial force caused by morphing are indeed continuous and bounded. Assumption 4 is commonly adopted in extant references on disturbance observer design^31,32,45 and facilitates the convergence analysis of the disturbance observer. Assumption 5 guarantees that the designed control laws are nonsingular ²¹ and can be satisfied by the model used in the numerical simulation.

The basic part of the controller is designed by the backstepping technique, which is described as follows.

• Step 1. Define $\tilde{h}=h-h_r$ and according to the idea in Bechlioulis and Rovithakis, ²⁴ the prescribed performance of altitude tracking error $\tilde{h}$ can be determined by the following inequality

$$-\underset{¯}{\varsigma }\rho (t)<\tilde{h}(t)<\overline{\varsigma }\rho (t),\forall t0$$

(7)

where $\rho (t)=({\rho }_0-{\rho }_{\infty }){\mathrm{e}}^{-\mathrm{at}}+{\rho }_{\infty }$ , $\underset{¯}{\varsigma }$ , $\overline{\varsigma }$ , ${\rho }_0$ , ${\rho }_{\infty }$ , and $a$ are positive constants and ${\rho }_0$ satisfies $-\underset{¯}{\varsigma }{\rho }_0<\tilde{h}(0)<\overline{\varsigma }{\rho }_0$ . It is not difficult to see that $\tilde{h}(t)$ is less than $\max \{\underset{¯}{\varsigma }{\rho }_0,\overline{\varsigma }{\rho }_0\}$ if Equation (7) holds. To facilitate the controller design, the constrained error condition (7) is converted into the equivalent boundedness requirement by introduction of a transformed error $\vartheta (t)$

$$\tilde{h}(t)=\rho (t)\mathcal{P}(\vartheta (t))$$

(8)

where $\mathcal{P}(\vartheta )=(\overline{\varsigma }{\mathrm{e}}^{\vartheta }-\underset{¯}{\varsigma }{\mathrm{e}}^{-\vartheta })/e^{\vartheta }+e^{-\vartheta }$ is a smooth and strictly increasing function. Differentiating $\vartheta (t)$ with respect to time yields

$$\begin{array}{c}\multicolumn{1}{c}{\stackrel{·}{\vartheta }(t)=\frac{1}{2}(\frac{1}{\underset{¯}{\varsigma }+\mathcal{P}(·)}+\frac{1}{\overline{\varsigma }-\mathcal{P}(·)})\stackrel{·}{\mathcal{P}}(·)}\\ =\Xi (t)(\stackrel{·}{\stackrel{~}{h}}(t)-\frac{\stackrel{·}{\rho }(t)\tilde{h}(t)}{\rho (t)})\hfill \end{array}$$

(9)

with $\Xi (t)=1/2\rho (t)((1/\underset{¯}{\varsigma }+\mathcal{P}(·))+(1/\overline{\varsigma }-\mathcal{P}(·)))$ . Notice that $\vartheta (t)=(1/2)\ln (\underset{¯}{\varsigma }/\overline{\varsigma })$ when $\tilde{h}(t)=0$ , hence we define $e_1=\vartheta -(1/2)\ln (\underset{¯}{\varsigma }/\overline{\varsigma })$ , and the time derivative of $e_1$ reads

$${\stackrel{·}{e}}_1=\Xi (\stackrel{·}{\stackrel{~}{h}}-\frac{\stackrel{·}{\rho }\tilde{h}}{\rho })=\Xi (V\gamma -{\stackrel{·}{h}}_r-\frac{\stackrel{·}{\rho }\tilde{h}}{\rho })$$

(10)

It is shown in Bechlioulis and Rovithakis ²⁴ that the prescribed performance of altitude tracking error $\tilde{h}$ is satisfied if $e_1$ is bounded. Using Equation (10), the virtual control law ${\gamma }_r$ is selected as

$${\gamma }_r=\frac{1}{V}(\frac{-k_1{e}_1}{\Xi }+{\stackrel{·}{h}}_r+\frac{\stackrel{·}{\rho }\tilde{h}}{\rho })$$

(11)

where $k_1$ is a positive constant. Introducing a new state variable ${\overline{\gamma }}_r$ , which is obtained by the following modified first-order filter with a time constant ${\tau }_1$

$${\tau }_1{\stackrel{·}{\overline{\gamma }}}_r+{\overline{\gamma }}_r={\gamma }_r-{\tau }_1\Xi V{e}_1$$

(12)

Define $e_2=\gamma -{\overline{\gamma }}_r$ and ${\tilde{\gamma }}_r={\overline{\gamma }}_r-{\gamma }_r$ , and it then follows that

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{e}}_1=\Xi (V{\gamma }_r+V{e}_2+V{\tilde{\gamma }}_r-{\stackrel{·}{h}}_r-\frac{\stackrel{·}{\rho }\tilde{h}}{\rho })}\\ =-k_1{e}_1+\Xi V{e}_2+\Xi V{\tilde{\gamma }}_r\hfill \end{array}$$

(13)

Choosing the following Lyapunov function candidate

$$L_{V_1}=\frac{1}{2}{e}_1^2+\frac{1}{2}{\tilde{\gamma }}_r^2$$

(14)

The time derivative of $L_{V_1}$ satisfies

$$\begin{array}{c}{\stackrel{·}{L}}_{V_1}=e_1(-k_1{e}_1+\Xi V{e}_2+\Xi V{\tilde{\gamma }}_r)\hfill \\ +{\tilde{\gamma }}_r(-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r-\Xi V{e}_1-{\stackrel{·}{\gamma }}_r)\hfill \\ =-k_1{e}_1^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2+\Xi V{e}_1{e}_2\hfill \\ -{\tilde{\gamma }}_r(\frac{\partial {\gamma }_r}{\partial e_1}{\stackrel{·}{e}}_1+\frac{\partial {\gamma }_r}{\partial {\stackrel{·}{h}}_r}{\stackrel{··}{h}}_r+\frac{\partial {\gamma }_r}{\partial V}\stackrel{·}{V}+\frac{\partial {\gamma }_r}{\partial \rho }\stackrel{·}{\rho }+\frac{\partial {\gamma }_r}{\partial V}\stackrel{·}{V}+\frac{\partial {\gamma }_r}{\partial \stackrel{·}{\rho }}\stackrel{··}{\rho })\hfill \\ =-k_1{e}_1^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2+\Xi V{e}_1{e}_2\hfill \\ -{\tilde{\gamma }}_r{B}_1(e_1,e_2,{\tilde{\gamma }}_r,h_r,{\stackrel{·}{h}}_r,{\stackrel{··}{h}}_r,V,\stackrel{·}{V},\rho ,\stackrel{·}{\rho },\stackrel{··}{\rho })\hfill \end{array}$$

(15)

with $B_1(·)$ being a continuous function corresponding to the time derivative of ${\gamma }_r$ .

• Step 2. Taking the time derivative of $e_2$ comes to

$${\stackrel{·}{e}}_2=\stackrel{·}{\gamma }-{\stackrel{·}{\overline{\gamma }}}_r=f_{2,i}+g_{2,i}\theta +d_{\gamma }+\frac{1}{{\tau }_1}{\tilde{\gamma }}_r+\Xi V{e}_1$$

(16)

where $i\in {\mathbf{S}}_M$ . The improved disturbance observer of the $\gamma$ loop is designed as follows

$$\{\begin{array}{c}\stackrel{·}{\hat{\gamma }}=f_{2,i}+g_{2,i}\theta +{\hat{d}}_{\gamma }+\frac{l_{1\gamma }}{{\epsilon }_{\gamma }}(\gamma -\hat{\gamma })+{\epsilon }_{\gamma }{l}_{3\gamma }{e}_2\hfill \\ \multicolumn{1}{c}{{\stackrel{·}{\hat{d}}}_{\gamma }=\frac{l_{2\gamma }}{{\epsilon }_{\gamma }^2}(\gamma -\hat{\gamma })+l_{4\gamma }{e}_2}\end{array}$$

(17)

where $\hat{\gamma }$ and ${\hat{d}}_{\gamma }$ are the states of the disturbance observer; $l_{1\gamma }$ , $l_{2\gamma }$ , and ${\epsilon }_{\gamma }$ are positive constants related to the observer gain; and $l_{3\gamma }$ and $l_{4\gamma }$ are constants to be specified later.

With the disturbance estimation ${\hat{d}}_{\gamma }$ , for any given constant ${\overline{\delta }}_{1,i}>0$ , an RBF NN is used to approximate the nonlinear term $f_{2,i}+{\hat{d}}_{\gamma }+(1/{\tau }_1){\tilde{\gamma }}_{\mathrm{r}}+2\Xi V{e}_1$ such that

$$\begin{array}{c}{f}_{2,i}+{\hat{d}}_{\gamma }+\frac{1}{{\tau }_1}{\tilde{\gamma }}_{\mathrm{r}}+2\Xi V{e}_1={\mathit{W}}_{1,i}^{\mathrm{T}}{\mathit{S}}_1({\mathit{X}}_1)+{\delta }_{1,i}({\mathit{X}}_1)\hfill \\ |{\delta }_{1,i}({\mathit{X}}_1)|{\overline{\delta }}_{1,i}\hfill \end{array}$$

(18)

where ${\mathit{W}}_{1,i}$ denotes the ideal weight vector, ${\mathit{X}}_1=[V,\gamma ,e_1,{\tilde{\gamma }}_r,\rho {]}^{\mathrm{T}}$ , and ${\delta }_{1,i}$ denotes the approximation error. Substituting Equation (18) into Equation (16) gives

$${\stackrel{·}{e}}_2={\mathit{W}}_{1,i}^T{\mathit{S}}_1(X_1)+{\delta }_{1,i}+g_{2,i}\theta -\Xi V{e}_1+d_{\gamma }-{\hat{d}}_{\gamma }$$

(19)

By Equation (19), it can be verified that

$$\begin{array}{cc}\multicolumn{1}{c}{e_2({\mathit{W}}_{1,i}^{\mathrm{T}}{\mathit{S}}_1({\mathit{X}}_1)+{\delta }_{1,i})}& ⩽|e_2|{\overline{\Theta }}_1(\parallel {\mathit{S}}_1({\mathit{X}}_1)\parallel +1)\hfill \\ \multicolumn{1}{c}{}& ⩽|e_2|{\overline{\Theta }}_1{Y}_1({\mathit{X}}_1)\hfill \end{array}$$

(20)

with ${\overline{\Theta }}_1=\max \{\parallel {\mathit{W}}_{1,i}\parallel ,{\overline{\delta }}_{1,i},i\in {\mathit{S}}_M\}$ , $Y_1({\mathit{X}}_1)=\parallel {\mathit{S}}_1({\mathit{X}}_1)\parallel +1$ . Applying Lemma 1 to Equation (20) yields

$$\begin{array}{c}{e}_2({\mathit{W}}_{1,i}^{\mathrm{T}}{\mathit{S}}_1({\mathit{X}}_1)+{\delta }_{1,i})⩽{\overline{\Theta }}_1{e}_2{Y}_1\tan h\left(\frac{e_2{Y}_1}{{\epsilon }_1}\right)\hfill \\ +{\overline{\Theta }}_1\chi {ϵ}_1\hfill \end{array}$$

(21)

where ${ϵ}_1$ is a positive design parameter. Let ${\theta }_r$ denote the virtual control law of $e_2$ and pass through a modified first-order filter with a time constant ${\tau }_2$ to obtain a new state variable ${\overline{\theta }}_r$

$${\tau }_2{\stackrel{·}{\overline{\theta }}}_r+{\overline{\theta }}_r={\theta }_r-{\tau }_2{g}_{2,i}{e}_2$$

(22)

Define $e_3=\theta -{\overline{\theta }}_r$ , ${\tilde{\theta }}_r={\overline{\theta }}_r-{\theta }_r$ , and ${\Theta }_1=g_{2m}^{-1}{\overline{\Theta }}_1$ and choose the Lyapunov function candidate as

$$L_{V_2}=L_{V_1}+\frac{1}{2}{e}_2^2+\frac{1}{2}{\tilde{\theta }}_r^2+\frac{g_{2m}}{2r_1}{\tilde{\Theta }}_1^2$$

(23)

where ${\tilde{\Theta }}_1={\hat{\Theta }}_1-{\Theta }_1$ and ${\hat{\Theta }}_1$ denotes the estimation of ${\Theta }_1$ and $r_1$ is a positive constant. The time derivative of $L_{V_2}$ is given by

$$\begin{array}{c}{\stackrel{·}{L}}_{V_2}⩽{\stackrel{·}{L}}_{V_1}+e_2{g}_{2,i}({\theta }_r+e_3+{\tilde{\theta }}_r)\hfill \\ +e_2({\overline{\Theta }}_1{Y}_1\times \tan h\left(\frac{e_2{Y}_1}{{\epsilon }_1}\right)-\Xi V{e}_1+d_{\gamma }-{\hat{d}}_{\gamma })\hfill \\ +{\tilde{\theta }}_r(-\frac{1}{{\tau }_2}{\tilde{\theta }}_r-g_{2,i}{e}_2-{\stackrel{·}{\theta }}_r)+\frac{g_{2m}}{r_1}{\tilde{\Theta }}_1{\stackrel{·}{\hat{\Theta }}}_1+{\overline{\Theta }}_1\chi {ϵ}_1\hfill \\ ⩽-k_1{e}_1^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2-\frac{1}{{\tau }_2}{\tilde{\theta }}_r^2-{\tilde{\gamma }}_r{B}_1(·)-{\tilde{\theta }}_r{\stackrel{·}{\theta }}_r+e_2{g}_{2,i}{e}_3\hfill \\ +e_2(g_{2m}{\Theta }_1{Y}_1\tan h\left(\frac{e_2{Y}_1}{1}\right)+g_{2,i}{\theta }_r)\hfill \\ +e_2(d_{\gamma }-{\hat{d}}_{\gamma })+\frac{g_{2m}}{r_1}{\tilde{\Theta }}_1{\stackrel{·}{\hat{\Theta }}}_1+{\overline{\Theta }}_1{\chi }_1\hfill \end{array}$$

(24)

From Equation (24), the virtual control law ${\theta }_r$ is designed as

$${\theta }_{\mathrm{r}}=-\mathrm{sign}(g_{2i})(\frac{k_2}{g_{2m}}{e}_2+{\hat{\Theta }}_1{Y}_1\tan h\left(\frac{e_2{Y}_1}{{ϵ}_1}\right))$$

(25)

where $k_2$ is a positive constant and $\mathrm{sign}(·)$ denotes the signum function. Furthermore, choosing the update law of ${\hat{\Theta }}_1$ as follows

$${\stackrel{·}{\hat{\Theta }}}_1=r_1{e}_2{Y}_1\tan h\left(\frac{e_2{Y}_1}{1}\right)-{\kappa }_1{\hat{\Theta }}_1$$

(26)

with ${\kappa }_1$ being a positive constant and ${\hat{\Theta }}_1(0)$ being set to satisfy ${\hat{\Theta }}_1⩽0$ . Applying Equations (25) and (26) to Equation (24) implies that

$$\begin{array}{c}{\stackrel{·}{L}}_{V_2}⩽-k_1{e}_1^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2-\frac{1}{{\tau }_2}{\tilde{\theta }}_r^2-{\tilde{\gamma }}_r{B}_1(·)-{\tilde{\theta }}_r{\stackrel{·}{\theta }}_r+e_2{g}_{2,i}{e}_3\hfill \\ +e_2(g_{2m}{\Theta }_1{Y}_1\tan h\left(\frac{e_2{Y}_1}{{ϵ}_1}\right)-\frac{|g_{2,i}|}{g_{2m}}{k}_2{e}_2\hfill \\ -|g_{2,i}|{\hat{\Theta }}_1{Y}_1\tan h\left(\frac{e_2{Y}_1}{{ϵ}_1}\right))+e_2(d_{\gamma }-{\hat{d}}_{\gamma })\hfill \\ +\frac{g_{2m}}{r_1}{\tilde{\Theta }}_1(r_1{e}_2{Y}_1\tan h\left(\frac{e_2{Y}_1}{{ϵ}_1}\right)-{\kappa }_1{\hat{\Theta }}_1)+{\overline{\Theta }}_1\chi {ϵ}_1\hfill \\ ⩽-k_1{e}_1^2-k_2{e}_2^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2-\frac{1}{{\tau }_2}{\tilde{\theta }}_r^2-{\tilde{\gamma }}_r{B}_1(·)\hfill \\ -{\tilde{\theta }}_r{B}_2(e_1,e_2,e_3,{\tilde{\gamma }}_r,{\tilde{\theta }}_r,\gamma -\hat{\gamma },d_{\gamma }-{\hat{d}}_{\gamma },{\hat{\Theta }}_1,h_r,{\stackrel{·}{h}}_r,{\stackrel{··}{h}}_r,\hfill \\ V,\stackrel{·}{V},\rho ,\stackrel{·}{\rho },\stackrel{··}{\rho })-\frac{g_{2m}{\kappa }_1}{r_1}{\tilde{\Theta }}_1{\hat{\Theta }}_1+e_2{g}_{2,i}{e}_3+e_2(d_{\gamma }-{\hat{d}}_{\gamma })\hfill \\ +{\overline{\Theta }}_1\chi {ϵ}_1\hfill \end{array}$$

(27)

where $B_2(·)$ denotes a continuous function related to the time derivative of ${\theta }_r$ .

• Step 3. Differentiating $e_3$ with respect to time yields

$${\stackrel{·}{e}}_3=\stackrel{·}{\theta }-{\stackrel{·}{\overline{\theta }}}_r=q+\frac{1}{{\tau }_2}{\tilde{\theta }}_r+g_{2,i}{e}_2$$

(28)

Let $q_r$ denote the virtual control law of $e_3$ and pass through a modified first-order filter with a time constant ${\tau }_3$ to obtain a new state variable ${\overline{q}}_r$

$${\tau }_3{\stackrel{·}{\overline{q}}}_r+{\overline{q}}_r=q_r-{\tau }_3{e}_3$$

(29)

Defining $e_4=q-{\overline{q}}_r$ and ${\tilde{q}}_r={\overline{q}}_r-q_r$ and using $L_{V_3}=L_{V_2}+(1/2)e_3^2+(1/2){\tilde{q}}_r^2$ as a Lyapunov function candidate, we obtain

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{L}}_{V_3}={\stackrel{·}{L}}_{V_2}+e_3(q_{\mathrm{r}}+e_4+{\tilde{q}}_r+\frac{1}{{\tau }_2}{\tilde{\theta }}_r+g_{2,i}{e}_2)}\\ +{\tilde{q}}_r(-\frac{1}{{\tau }_3}{\tilde{q}}_r-e_3-{\stackrel{·}{q}}_r)\hfill \\ ={\stackrel{·}{L}}_{V_2}+e_3(q_r+e_4+\frac{1}{{\tau }_2}{\tilde{\theta }}_r+g_{2,i}{e}_2)\hfill \\ +{\tilde{q}}_r(-\frac{1}{{\tau }_3}{\tilde{q}}_r-{\stackrel{·}{q}}_r)\hfill \end{array}$$

(30)

Taking

$$q_r=-k_3{e}_3-\frac{1}{{\tau }_2}{\tilde{\theta }}_r$$

(31)

delivers

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{L}}_{V_3}={\stackrel{·}{L}}_{V_2}-k_3{e}_3^2-\frac{1}{{\tau }_3}{\tilde{q}}_r^2+e_3{e}_4+g_{2,i}{e}_2{e}_3}\\ -{\tilde{q}}_r{B}_3(e_1,e_2,e_3,e_4,{\tilde{\gamma }}_r,{\tilde{\theta }}_r,{\tilde{q}}_r,\gamma -\hat{\gamma },d_{\gamma }-{\hat{d}}_{\gamma },{\hat{\Theta }}_1,\hfill \\ h_r,{\stackrel{·}{h}}_r,{\stackrel{··}{h}}_r,V,\stackrel{·}{V},\rho ,\stackrel{·}{\rho },\stackrel{··}{\rho })\hfill \end{array}$$

(32)

with $B_3(·)$ being a continuous function related to the time derivative of $q_r$ .

• Step 4. The time derivative of $e_4$ is given by

$${\stackrel{·}{e}}_4=\stackrel{·}{q}-{\stackrel{·}{\overline{q}}}_r=f_{4,i}+g_{4,i}{\delta }_{e1}+d_q+\frac{1}{{\tau }_3}{\tilde{q}}_r+e_3$$

(33)

The improved disturbance observer for the $q$ loop is designed as

$$\{\begin{array}{c}\stackrel{·}{\hat{q}}=f_{4,i}+g_{4,i}{\delta }_{e1}+{\hat{d}}_q+\frac{l_{1q}}{{\epsilon }_q}(q-\hat{q})+{\epsilon }_q{l}_{3q}{e}_4\hfill \\ \multicolumn{1}{c}{{\stackrel{·}{\hat{d}}}_q=\frac{l_{2q}}{{\epsilon }_q^2}(q-\hat{q})+l_{4q}{e}_4}\end{array}$$

(34)

where $\hat{q}$ and ${\hat{d}}_q$ are the states of the disturbance observer; $l_{1q}$ , $l_{2q}$ , and ${\epsilon }_q$ are the positive constants related to the observer gain; and $l_{3q}$ and $l_{4q}$ are the constants to be specified later.

With the estimated value ${\hat{d}}_q$ , the nonlinear term $f_{4,i}+{\hat{d}}_q+(1/{\tau }_3){\tilde{q}}_r+2e_3$ is approximated by an RBF NN that satisfies

$$\begin{array}{c}{f}_{4,i}+{\hat{d}}_q+\frac{1}{{\tau }_3}{\tilde{q}}_{\mathrm{r}}+2e_3={\mathit{W}}_{2,i}^{\mathrm{T}}{\mathit{S}}_2({\mathit{X}}_2)+{\delta }_{2,i}\hfill \\ |{\delta }_{2,i}|⩽{\overline{\delta }}_{2,i}\hfill \end{array}$$

(35)

for any given constant ${\overline{\delta }}_{2,i}>0$ , where ${\mathit{W}}_{2,i}$ denotes the ideal weight vector, ${\mathit{X}}_2=[V,\gamma ,\theta ,q,{\hat{d}}_q,e_3,{\tilde{q}}_r{]}^{\mathrm{T}}$ ; ${\delta }_{2,i}$ denotes the approximation error, and it can be verified that

$${\stackrel{·}{e}}_4={\mathit{W}}_{2,i}^{\mathrm{T}}{\mathit{S}}_2({\mathit{X}}_2)+{\delta }_{2,i}+g_{4,i}{\delta }_{e1}-e_3+d_q-{\hat{d}}_q$$

(36)

Since

$$\begin{array}{c}{e}_4({\mathit{W}}_{2,i}^{\mathrm{T}}{\mathit{S}}_2({\mathit{X}}_2)+{\delta }_{2,i})|e_4|⩽{\overline{\Theta }}_2(\parallel {\mathit{S}}_2({\mathit{X}}_2)\parallel +1)\hfill \\ |e_4|{\overline{\Theta }}_2{Y}_2({\mathit{X}}_2)\hfill \end{array}$$

(37)

where ${\overline{\Theta }}_2=\max \{\parallel {\mathit{W}}_{2,i}\parallel ,{\overline{\delta }}_{2,i},i\in {\mathit{S}}_M\}$ and $Y_2({\mathit{X}}_2)=\parallel S_2({\mathit{X}}_2)\parallel +1$ . Applying Lemma 1 to Equation (37) results in

$$\begin{array}{c}{e}_4({\mathit{W}}_{2,i}^{\mathrm{T}}{\mathit{S}}_2({\mathit{X}}_2)+{\delta }_{2,i}){\overline{\Theta }}_2{e}_4{Y}_2\tan h\left(\frac{e_4{Y}_2}{{ϵ}_2}\right)\hfill \\ +{\overline{\Theta }}_2\chi {ϵ}_2\hfill \end{array}$$

(38)

where ${ϵ}_2$ denotes a positive constant. Define ${\Theta }_2=g_{4m}^{-1}{\overline{\Theta }}_2$ and choose the Lyapunov function candidate as

$$L_{V_4}=L_{V_3}+\frac{1}{2}{e}_4^2+\frac{g_{4m}}{2r_2}{\tilde{\Theta }}_2^2$$

(39)

where ${\tilde{\Theta }}_2={\hat{\Theta }}_2-{\Theta }_2$ , ${\hat{\Theta }}_2$ denotes the estimation of ${\Theta }_2$ , and $r_2$ is a positive constant. It follows that

$$\begin{array}{c}{\stackrel{·}{L}}_{V_4}⩽{\stackrel{·}{L}}_{V_3}+e_4{g}_{4,i}{\delta }_{e1}+e_4(g_{4m}{\Theta }_2{Y}_2\hfill \\ \times \tan h\left(\frac{e_4{Y}_2}{{\epsilon }_2}\right)-e_3+d_q-{\hat{d}}_q)+\frac{g_{4m}}{r_2}{\tilde{\Theta }}_2{\stackrel{·}{\hat{\Theta }}}_2\hfill \\ +{\overline{\Theta }}_2\chi {\epsilon }_2\hfill \end{array}$$

(40)

Consequently, the final basic control input ${\delta }_{e1}$ is

$${\delta }_{e1}=-\mathrm{sign}(g_{4i})(\frac{k_4}{g_{4m}}{e}_4+{\hat{\Theta }}_2{Y}_2\tan h\left(\frac{e_4{Y}_2}{{\epsilon }_2}\right))$$

(41)

with $k_4$ being a positive constant. The update law for ${\hat{\Theta }}_2$ is defined as

$${\stackrel{·}{\hat{\Theta }}}_2=r_2{e}_4{Y}_2\tan h\left(\frac{e_4{Y}_2}{{\epsilon }_2}\right)-{\kappa }_2{\hat{\Theta }}_2$$

(42)

where ${\kappa }_2$ is a positive constant and ${\hat{\Theta }}_2(0)$ is selected to satisfy ${\hat{\Theta }}_2⩽0$ . Invoking Equations (41) and (42) entails that

$$\begin{array}{c}{\stackrel{·}{L}}_{V_4}⩽{\stackrel{·}{L}}_{V_3}-e_3{e}_4+e_4(g_{4m}{\Theta }_2{Y}_2\tan h\left(\frac{e_4{Y}_2}{{ϵ}_2}\right)\hfill \\ -\frac{|g_{4,i}|}{g_{4m}}{k}_4{e}_4-|g_{4,i}|{\hat{\Theta }}_2{Y}_2\tan h\left(\frac{e_4{Y}_2}{{ϵ}_2}\right))\hfill \\ +\frac{g_{4m}}{r_2}{\tilde{\Theta }}_2(r_2{e}_4{Y}_2\tan h\left(\frac{e_4{Y}_2}{{ϵ}_2}\right)-{\kappa }_2{\hat{\Theta }}_2)\hfill \\ +e_4(d_q-{\hat{d}}_q)+{\overline{\Theta }}_2\chi {ϵ}_2\hfill \\ ⩽{\stackrel{·}{L}}_{V_3}-k_4{e}_4^2-\frac{g_{4m}{\kappa }_2}{r_2}{\tilde{\Theta }}_2{\hat{\Theta }}_2-e_3{e}_4+e_4(d_q-{\hat{d}}_q)\hfill \\ +{\overline{\Theta }}_2\chi {ϵ}_2\hfill \end{array}$$

(43)

Remark 2

It should be pointed out that the nonlinearities including the disturbances of the systems cannot be canceled out directly when designing the common virtual control laws since the systems are of switched type. By the adoption of NNs, the focus becomes that the upper bound of related weight vector which is estimated and adopted to generate the common virtual control laws. Specifically, the updating law of the estimation of the upper bound is to ensure the stability of closed-loop system via the property of inequalities. Such technique is generally adopted for the research on backstepping control of switched nonlinear systems. ²¹

Remark 3

We note that there are also several excellent studies on neural learning design. For instance, tracking control of a flexible hypersonic vehicle via a combination of neural approximation and disturbance estimation is addressed in Xu et al., ⁴³ where both the tracking error and prediction error are exploited to generate the weight updating law of NNs. The robust design is integrated into the updating laws for both the weight vector and approximation error of the adopted NNs in Xu ⁴⁶ to develop a direct neural controller. The research in Xu et al.^43,46 mainly adopts the general nonlinear systems in lower triangular form and there exists no system switching. As a result, the nonlinearity of the system can be canceled out directly via the employment of NNs as well as disturbance observers. However, the dynamics of morphing aircraft is described by switched nonlinear systems where the nonlinearities cannot be eliminated directly since these nonlinearities switch as the sweep angle of the morphing aircraft varies, which makes our work different from the research in Xu et al.^43,46

For the basic part controller, we have the following theorem:

Theorem 1

Consider the closed-loop switched nonlinear systems formed of the plant (2); the improved disturbance observers (17) and (34); the modified dynamic surfaces (12), (22), and (29); the parameter update laws (26) and (42); and the basic control law (41). Suppose that ${\delta }_{e2}=0$ , Assumptions 1–5 are satisfied, and initial values of the closed-loop system are bounded. Then, there exist the design parameters $k_j(j=1,2,3,4)$ , ${\tau }_j(j=1,2,3)$ , $l_{j\iota }$ , ${\epsilon }_{\iota }(j=1,2,3,4;\iota =\gamma ,q)$ , $r_j$ , ${ϵ}_j$ , and ${\kappa }_j(j=1,2)$ such that all signals of the closed-loop system remain bounded, and the altitude tracking error $\tilde{h}$ converges to a neighborhood of the origin with prescribed performance bound.

Proof

Define scaled estimation errors of the designed disturbance observers (17) and (34) as follows

$${\mathbf{\eta }}_{\iota }=[{\eta }_{\iota 1},{\eta }_{\iota 2}{]}^{\mathrm{T}},{\eta }_{\iota 1}=\frac{\iota -\hat{\iota }}{{\epsilon }_{\iota }},{\eta }_{\iota 2}=d_{\iota }-{\hat{d}}_{\iota },\iota =\gamma ,q$$

and the dynamics of the scaled estimation errors $h_{\iota }$ then become

$$\{\begin{array}{c}{\epsilon }_{\iota }{\stackrel{·}{\eta }}_{\iota 1}={\eta }_{\iota 2}-l_{1\iota }{\eta }_{\iota 1}-{\epsilon }_{\iota }{l}_{3\iota }{e}_2\hfill \\ \multicolumn{1}{c}{{\epsilon }_{\iota }{\stackrel{·}{\eta }}_{\iota 2}=-l_{2\iota }{\eta }_{\iota 1}-{\epsilon }_{\iota }{l}_{4\iota }{e}_2+{\epsilon }_{\iota }{\stackrel{·}{d}}_{\iota },\iota =\gamma ,q}\end{array}$$

(44)

Let ${\mathit{A}}_{\iota }=\left[\begin{array}{cc}-l_{1\iota }\hfill & 1\hfill \\ \multicolumn{1}{c}{-l_{2\iota }}& 0\hfill \end{array}\right](\iota =\gamma ,q)$ , since $l_{1\iota }$ and $l_{2\iota }$ are positive constants, ${\mathit{A}}_{\iota }$ is a Hurwitz matrix, and then for any ${\overline{k}}_{\iota }>0$ , there exists a positive definite matrix ${\mathit{P}}_{\iota }$ that satisfies ${\mathit{P}}_{\iota }{\mathit{A}}_{\iota }+{\mathit{A}}_{\iota }^{\mathrm{T}}{\mathit{P}}_{\iota }=-{\overline{k}}_{\iota }{\mathit{I}}_2$ . Using $L_{V_{\iota }}={\mathbf{\eta }}_{\iota }^{\mathrm{T}}{\mathit{P}}_{\iota }{\mathbf{\eta }}_{\iota }$ as a Lyapunov function candidate for Equation (44), we obtain

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{L}}_{V_{\iota }}=-\frac{{\overline{k}}_{\iota }}{{\epsilon }_{\iota }}\parallel {\mathbf{\eta }}_{\iota }{\parallel }^2-2e_2[l_{3\iota }{l}_{4\iota }]{\mathit{P}}_{\iota }{\mathbf{\eta }}_{\iota }+2{\mathbf{\eta }}_{\iota }^{\mathrm{T}}{\mathit{P}}_{\iota }{[0{\stackrel{·}{d}}_{\iota }]}^{\mathrm{T}}}\\ =-\frac{{\overline{k}}_{\iota }}{{\epsilon }_{\iota }}\parallel {\mathbf{\eta }}_{\iota }{\parallel }^2-2(p_{\iota 11}{l}_{3\iota }+p_{\iota 12}{l}_{4\iota }){\eta }_{\iota 1}{e}_2\hfill \\ -2(p_{\iota 12}{l}_{3\iota }+p_{\iota 22}{l}_{4\iota }){\eta }_{\iota 2}{e}_2+2{\mathbf{\eta }}_{\iota }^{\mathrm{T}}{\mathit{P}}_{\iota }{[0{\stackrel{·}{d}}_{\iota }]}^{\mathrm{T}}\hfill \end{array}$$

(45)

with $p_{\iota j\ell }$ being the $j$ , $\ell$ entry of ${\mathit{P}}_{\iota }$ $(j=1,2\mathrm{and}\ell =1,2)$ .

The constants $l_{3\iota }$ and $l_{4\iota }$ are chosen to satisfy

$$\{\begin{array}{c}{p}_{\iota 11}{l}_{3\iota }+p_{\iota 12}{l}_{4\iota }=0,\hfill \\ \multicolumn{1}{c}{p_{\iota 12}{l}_{3\iota }+p_{\iota 22}{l}_{4\iota }=\frac{1}{2}}\end{array}$$

(46)

then the derivative of $L_{V_{\iota }}$ along the trajectories of Equation (44) reads

$${\stackrel{·}{L}}_{V_{\iota }}=-\frac{{\overline{k}}_{\iota }}{{\epsilon }_{\iota }}\parallel {\mathbf{\eta }}_{\iota }{\parallel }^2-{\mathbf{\eta }}_{\iota 2}{e}_2+2{\mathbf{\eta }}_{\iota }^{\mathrm{T}}{\mathit{P}}_{\iota }[0{\stackrel{·}{d}}_{\iota }{]}^{\mathrm{T}}$$

(47)

For the closed-loop system, consider now the Lyapunov function candidate

$$L_V=L_{V_4}+L_{V_{\gamma }}+L_{V_q}$$

(48)

Invoking Equations (27), (32), and (43), it can be shown that

$$\begin{array}{c}{\stackrel{·}{L}}_{V_4}⩽-\sum _{j=1}^4{k}_j{e}_j^2-\frac{1}{{\tau }_1}{\tilde{\gamma }}_{\mathrm{r}}^2-\frac{1}{{\tau }_2}{\tilde{\theta }}_{\mathrm{r}}^2-\frac{1}{{\tau }_3}{\tilde{q}}_{\mathrm{r}}^2-{\tilde{\gamma }}_{\mathrm{r}}{B}_1(·)\hfill \\ -{\tilde{\theta }}_{\mathrm{r}}{B}_2(·)-{\tilde{q}}_{\mathrm{r}}{B}_3(·)+2g_{2,i}{e}_2{e}_3+e_2{\eta }_{\gamma 2}\hfill \\ \begin{array}{c}+e_4{\eta }_{q2}-\frac{g_{2m}{\kappa }_1}{r_1}{\tilde{\Theta }}_1{\hat{\Theta }}_1-\frac{g_{4m}{\kappa }_2}{r_2}{\tilde{\Theta }}_2{\hat{\Theta }}_2\hfill \\ \multicolumn{1}{c}{+\chi ({\overline{\Theta }}_1{ϵ}_1+{\overline{\Theta }}_2{ϵ}_2)}\end{array}\hfill \end{array}$$

(49)

By Young’s inequality, ⁴⁷ one gets

$$2g_{2,i}{e}_2{e}_3⩽g_{2M}(e_2^2+e_3^2)$$

(50)

Moreover, with the definition of ${\tilde{\Theta }}_j(j=1,2)$ , we can write

$$-{\tilde{\Theta }}_j{\hat{\Theta }}_j⩽-{\tilde{\Theta }}_j^2-{\tilde{\Theta }}_j{\Theta }_j⩽-\frac{1}{2}{\tilde{\Theta }}_j^2+\frac{1}{2}{\Theta }_j^2$$

(51)

Substituting Equations (50) and (51) into Equation (49) yields

$$\begin{array}{c}{\stackrel{·}{L}}_{V_4}⩽-k_1{e}_1^2-(k_2-g_{2M})e_2^2-(k_3-g_{2M})e_3^2-k_4{e}_4^2\hfill \\ -\frac{1}{{\tau }_1}{\tilde{\gamma }}_r^2-\frac{1}{{\tau }_2}{\tilde{\theta }}_r^2-\frac{1}{{\tau }_3}{\tilde{q}}_r^2-\frac{g_{2m}{\kappa }_1}{2r_1}{\tilde{\Theta }}_1^2-\frac{g_{4m}{\kappa }_2}{2r_2}{\tilde{\Theta }}_2^2\hfill \\ -{\tilde{\gamma }}_r{B}_1(·)-{\tilde{\theta }}_{\mathrm{r}}{B}_2(·)-{\tilde{q}}_r{B}_3(·)+e_2{\eta }_{\gamma 2}+e_4{\eta }_{q2}\hfill \\ +\frac{g_{2m}{\kappa }_1}{r_1}{\Theta }_1^2+\frac{g_{4m}{\kappa }_2}{2r_2}{\Theta }_2^2+\chi ({\overline{\Theta }}_1{ϵ}_1+{\overline{\Theta }}_2{ϵ}_2)\hfill \end{array}$$

(52)

By Assumptions 2 and 3, for given positive constants ${\nu }_j(j=1,2,\dots ,7)$ , the sets ${\Phi }_1=\{h_r^2+{\stackrel{·}{h}}_r^2+{\stackrel{··}{h}}_r^2⩽{\nu }_1\}$ , ${\Phi }_2=\{V^2+{\stackrel{·}{V}}^2⩽{\nu }_2\}$ , ${\Phi }_3=\{{\rho }^2+{\stackrel{·}{\rho }}^2+{\stackrel{··}{\rho }}^2⩽{\nu }_3\}$ , ${\Phi }_4=\{\sum _{j=0}^4{e}_j^2⩽{\nu }_4\}$ , ${\Phi }_5=\{{\tilde{\gamma }}_r^2+{\tilde{\theta }}_r^2+{\tilde{q}}_r^2⩽{\nu }_5\}$ , ${\Phi }_6=\{\parallel {\eta }_{\gamma }{\parallel }^2+\parallel {\eta }_q{\parallel }^2⩽{\nu }_6\}$ , and ${\Phi }_7=\{\sum _{j=1}^2{\tilde{\Theta }}_j^2⩽{\nu }_7\}$ are compact. It can then be verified that all arguments of the continuous functions $B_j(·)(j=1,2,3)$ are bounded on $\underset{\ell =1}{\overset{7}{\Pi }}{\Phi }_{\ell }$ , and there exist related positive constants $M_j$ such that $|B_j(·)|⩽M_j$ . Hence

$$\begin{array}{c}{\stackrel{·}{L}}_{V_4}⩽-k_1{e}_1^2-(k_2-g_{2M})e_2^2-(k_3-g_{2M})e_3^2-k_4{e}_4^2\hfill \\ -(\frac{1}{{\tau }_1}-\frac{1}{2}){\tilde{\gamma }}_r^2-(\frac{1}{{\tau }_2}-\frac{1}{2}){\tilde{\theta }}_r^2-(\frac{1}{{\tau }_3}-\frac{1}{2}){\tilde{q}}_r^2\hfill \\ -\frac{g_{2m}{\kappa }_1}{2r_1}{\tilde{\Theta }}_1^2-\frac{g_{4m}{\kappa }_2}{2r_2}{\tilde{\Theta }}_2^2+e_2{\eta }_{\gamma 2}+e_4{\eta }_{q2}\hfill \\ +\frac{1}{2}\sum _{j=1}^3{M}_j^2+\frac{g_{2m}{\kappa }_1}{r_1}{\Theta }_1^2+\frac{g_{4m}{\kappa }_2}{2r_2}{\Theta }_2^2\hfill \\ +\chi ({\overline{\Theta }}_1{ϵ}_1+{\overline{\Theta }}_2{ϵ}_2)\hfill \end{array}$$

(53)

Note that with ${\eta }_{\iota }\in {\Phi }_6$ , $2{\mathbf{\eta }}_{\mathbf{\iota }}^{\mathrm{T}}{\mathit{P}}_{\mathbf{\iota }}[0{\stackrel{·}{\mathbf{d}}}_{\mathbf{\iota }}{]}^{\mathrm{T}}⩽2{\mathbf{\lambda }}_{\max }({\mathit{P}}_{\mathbf{\iota }}){\mathbf{N}}_{\mathbf{\iota }}\sqrt{{\mathbf{\nu }}_6}$ $(\iota =\gamma ,q)$ . A combination of Equations (46) and (47) arrives at

$$\begin{array}{c}{\stackrel{·}{L}}_V⩽-k_1{e}_1^2-(k_2-g_{2M})e_2^2-(k_3-g_{2M})e_3^2-k_4{e}_4^2\hfill \\ -(\frac{1}{{\tau }_1}-\frac{1}{2}){\tilde{\gamma }}_r^2-(\frac{1}{{\tau }_2}-\frac{1}{2}){\tilde{\theta }}_r^2-(\frac{1}{{\tau }_3}-\frac{1}{2}){\tilde{q}}_r^2\hfill \\ -\frac{g_{2m}{\kappa }_1}{2r_1}{\tilde{\Theta }}_1^2-\frac{g_{4m}{\kappa }_2}{2r_2}{\tilde{\Theta }}_2^2-\frac{{\overline{k}}_{\gamma }}{{\epsilon }_{\gamma }}\parallel {\eta }_{\gamma }{\parallel }^2-\frac{{\overline{k}}_q}{{\epsilon }_q}\parallel {\eta }_q{\parallel }^2\hfill \\ +\frac{1}{2}\sum _{j=1}^3{M}_j^2+\frac{g_{2m}{\kappa }_1}{r_1}{\Theta }_1^2+\frac{g_{4m}{\kappa }_2}{2r_2}{\Theta }_2^2\hfill \\ +\chi ({\overline{\Theta }}_1{\epsilon }_1+{\overline{\Theta }}_2{ϵ}_2)+2({\lambda }_{\max }({\mathit{P}}_{\gamma })N_{\gamma }\hfill \\ +{\lambda }_{\max }({\mathit{P}}_q)N_q)\sqrt{{\nu }_6}\hfill \end{array}$$

(54)

Setting the parameters to satisfy

$$\begin{array}{c}{k}_1>0,k_2>g_{2M},k_3>g_{2M},k_4>0,0<{\tau }_1<2,\hfill \\ 0<{\tau }_2<2,0<{\tau }_3<2,{\kappa }_1>0,{\kappa }_2>0,{\epsilon }_{\gamma }>0,{\epsilon }_q>0\hfill \end{array}$$

(55)

Then

$${\stackrel{·}{L}}_V⩽-{\lambda }_V{L}_V+N_V$$

(56)

where ${\lambda }_V=min\{2k_1,2(k_2-g_{2M}),2(k_3-g_{2M}),2k_4,(2/{\tau }_1)-1,(2/{\tau }_2)-1,(2/{\tau }_3)-1,{\kappa }_1,{\kappa }_2,({\overline{k}}_{\gamma }/{\epsilon }_{\gamma }{\lambda }_{\max (P_{\gamma })}),({\overline{k}}_q/{\epsilon }_q{\lambda }_{\max (P_q)})\}$ , and

$$\begin{array}{c}\multicolumn{1}{c}{N_V=\frac{1}{2}\sum _{j=1}^3{M}_j^2+\frac{g_{2m}{\kappa }_1}{r_1}{\Theta }_1^2+\frac{g_{4m}{\kappa }_2}{2r_2}{\Theta }_2^2}\\ +\chi ({\overline{\Theta }}_1{ϵ}_1+{\overline{\Theta }}_2{ϵ}_2)+2({\lambda }_{\max }({\mathit{P}}_{\gamma })N_{\gamma }\hfill \\ +{\lambda }_{\max }({\mathit{P}}_q)N_q)\sqrt{{\nu }_6}\hfill \end{array}$$

Using the comparison principle in Equation (55), we have

$$L_V(t)(L_V(0)-\frac{N_V}{{\lambda }_V})e^{-{\lambda }_{V}t}+\frac{N_V}{{\lambda }_V}$$

(57)

It is straightforward to see that all the signals in the closed-loop system are bounded, and the transformed error $\vartheta$ is also bounded. Therefore, it can be concluded the altitude tracking error $\tilde{h}$ converges to a neighborhood of the origin with prescribed performance bound.

Remark 4

It is noteworthy that the improved disturbance observers (17) and (34) and the modified dynamic surfaces (12), (22), and (29) contain the coupling terms related to $e_j(j=1,2,3,4)$ , which constitutes the main difference between our method and those on ESO,^31,32 dynamic surface technique, ¹⁰ and the references therein. As can be seen in the proof of Theorem 1, the introduction of coupling terms related to $e_j(j=1,2,3,4)$ is mainly for the sake of canceling out the cross terms in the derivatives of Lyapunov function candidates and facilitates the stability analysis.

Supplementary part design

To expedite the control performance, the reinforcement learning approach is further applied to design the supplementary part controller, where both critic network and action network utilized are RBF NNs. The critic network is mainly for obtaining the estimation of a user-defined cost function, while the action network is mainly for generating the related control strategy which minimizes the estimated value of the cost function. The structure of the reinforcement learning scheme is shown in Figure 2. It is noted that the related error variables defined in the backstepping design are regarded as inputs of both the critic network and action network, and a supplementary control input is generated correspondingly.

Figure 2.

Structure diagram of the reinforcement learning scheme.

The input of the critic network is denoted by $x_c(t)=e_1(t)$ and the utility function is $\phi (x_c(t))={\overline{k}}_r{x}_c^2(t)$ , with ${\overline{k}}_r$ being a positive constant, and the related cost function is formulated as

$$J(x_c(t))=\underset{t}{\overset{\infty }{\int }}{e}^{-\beta (s-t)}\phi (x_c(s))\mathrm{d}s$$

(58)

where $\beta$ is a positive constant related to the discount factor. Since the exact value of $J(x_c(t))$ cannot be obtained directly, it is estimated by the critic network and then minimized by the action network. By definition, we can regard $J(x_c(t))$ as a continuous function of $x_c(t)$ . Let $J(x_c(t))={\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c(t))+{\sigma }_c(t)$ , where ${\sigma }_c(t)$ is the approximation error. The output of the critic network is the estimation of $J(x_c(t))$ , which is denoted by $\hat{J}(x_c(t))={\hat{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c(t))$ .

The error function for training the critic network is then constructed as

$$e_c=\phi (x_c)-\beta \hat{J}(x_c)+\nabla \hat{J}(x_c){\stackrel{·}{x}}_c$$

(59)

where $\nabla \hat{J}(x_c)$ denotes the gradient of $\hat{J}(x_c)$ with respect to $x_c$ . The related objection function for weight updating reads

$$E_c=\frac{1}{2}{e}_c^2$$

(60)

The weight updating of the critic network is to minimize $E_c$ in a gradient-descent way, that is

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{\hat{\mathit{W}}}}_c=-{\mu }_c\frac{\partial E_c}{\partial {\hat{\mathit{W}}}_c}}\\ =-{\mu }_c{e}_c\frac{\partial }{\partial {\hat{\mathit{W}}}_c}(\phi (x_c)-\beta \hat{J}(x_c)+\nabla \hat{J}(x_c){\stackrel{·}{x}}_c)\hfill \\ =-{\mu }_c{e}_c\frac{\partial }{\partial {\hat{\mathit{W}}}_c}(-\beta {\hat{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)+{\hat{\mathit{W}}}_c^{\mathrm{T}}\nabla {\mathit{S}}_c(x_c){\stackrel{·}{x}}_c)\hfill \\ =-{\mu }_c{e}_c(-\beta {\mathit{S}}_c(x_c)+\nabla {\mathit{S}}_c(x_c){\stackrel{·}{x}}_c)\hfill \\ =-{\mu }_c{e}_c\mathit{C}\hfill \end{array}$$

(61)

where ${\mu }_c>0$ denotes the learning rate, $\mathit{C}=-\beta {\mathit{S}}_c(x_c)+\nabla {\mathit{S}}_c(x_c){\stackrel{·}{x}}_c$ , and $\nabla {\mathit{S}}_c(x_c)$ denotes the gradient of ${\mathit{S}}_c(x_c)$ with respect to $x_c$ .

The action network takes ${\mathit{x}}_a=[e_1,e_2,e_3,e_4{]}^{\mathrm{T}}$ as input and generates the supplementary control input ${\delta }_{e2}$ , that is, ${\delta }_{e2}={\hat{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)$ . Let ${\mathit{W}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\sigma }_a$ denote the control input that minimizes the cost function $\hat{J}(x_c)$ , ${\tilde{\mathit{W}}}_a={\hat{\mathit{W}}}_a-{\mathit{W}}_a$ denote the weight error, and ${\sigma }_a$ be the approximation error. The error function for training the action network is designed as

$$e_a={\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a(\hat{J}(x_c)-U_c)$$

(62)

where ${\overline{k}}_a>0$ is a design parameter and $U_c$ is the desired ultimate cost objective. It can be seen from Equation (57) that we expect $J(x_c)$ to converge to zero and take $U_c=0$ . Similarly, the gradient-descent algorithm is applied to weight updating of the action network as follows

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{\hat{\mathit{W}}}}_a=-{\mu }_a\frac{\partial E_a}{\partial {\hat{\mathit{W}}}_a}}\\ =-{\mu }_a\frac{\partial E_a}{\partial e_a}\frac{\partial e_a}{\partial {\hat{\mathit{W}}}_a}\hfill \\ =-{\mu }_a({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a\hat{J}(x_c)){\mathit{S}}_a({\mathit{x}}_a)\hfill \end{array}$$

(63)

where $E_a=(1/2)e_a^2$ and ${\mu }_a>0$ denotes the learning rate. Since ${\tilde{W}}_a$ is unavailable, the weight updating law is changed to

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{\hat{\mathit{W}}}}_a=-{\mu }_a({\hat{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a\hat{J}(x_c)){\mathit{S}}_a({\mathit{x}}_a)}\\ =-{\mu }_a({\hat{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a{\hat{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)){\mathit{S}}_a(x_a)\hfill \end{array}$$

(64)

Remark 5

The reinforcement learning method adapts to the uncertain environment without model. Specifically, the reinforcement learning method is a data-driven framework for solving optimal control problems and the requirement of a model of the system dynamics or even the formulation of the reward function can be avoided. However, as stated in Buoniu et al., ³⁴ the model-free property of the reinforcement learning method does not mean that the knowledge on the system dynamics cannot be employed when the model is in fact available, no matter partially or fully. The model of the system dynamics can be exploited to develop reinforcement learning algorithms depending upon the specific conditions, and many related reinforcement learning approaches have been proposed for computing optimal control input of general affine nonlinear systems in theory.^35–37 To summarize, the reinforcement learning method is not restricted to the system model information, and the model can be appropriately adopted once available. In this paper, the reinforcement learning method without model information is devised to further decrease the altitude tracking error since the system model in the presence of the basic part controller is too complex and hence regarded as unknown.

Remark 6

For the supplement part, the RBF NNs are exploited to constitute the critic–action structure of the reinforcement learning algorithm. Both the critic network and action network can be of any appropriate type, such as the multi-layer perceptron structure with one hidden layer, ⁴⁰ and the adoption of RBF NNs is mainly for the convenience of weight updating and convergence analysis. The weight updating for the critic network and action network is mainly along the line of reinforcement learning approach in a gradient-descent way. It is noted that the specific choice and structure of the critic network and action network is not the focus of our work. We aim to show that the integration of reinforcement learning technique can further improve the control performance, which is demonstrated in the numerical simulation part.

For the supplementary part controller, we have the following theorem:

Theorem 2

Consider the critic–action NN structure established as Equations (61) and (64). Let ${\tilde{\mathit{W}}}_c={\hat{\mathit{W}}}_c-{\mathit{W}}_c$ , and then, the weight estimation errors ${\tilde{W}}_c$ and ${\tilde{W}}_a$ are uniformly ultimately bounded with the following condition

$${\overline{k}}_a⩽\min \{\frac{1}{4},\frac{{\mu }_c{({\tilde{\mathit{W}}}_c^{\mathrm{T}}C)}^2}{2{\mu }_a{({\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c))}^2}\}$$

(65)

Proof

By Equation (61), the dynamics of ${\tilde{\mathit{W}}}_c$ is formulated as

$${\stackrel{·}{\stackrel{~}{\mathit{W}}}}_c=-{\mu }_c\mathit{C}(\phi (x_c)-\beta \hat{J}(x_c)+\nabla \hat{J}(x_c){\stackrel{·}{x}}_c)$$

(66)

Notice that

$$0=\phi (x_c)-\beta J(x_c)+\nabla J(x_c){\stackrel{·}{x}}_c$$

(67)

Using $J(x_c)={\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)+{\sigma }_c$ in Equation (67) yields

$$\begin{array}{c}\multicolumn{1}{c}{\phi (x_c)=\beta ({\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)+{\sigma }_c)-\nabla ({\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)){\stackrel{·}{x}}_c}\\ ={\mathit{W}}_c^{\mathrm{T}}(\beta {\mathit{S}}_c(x_c)-\nabla {\mathit{S}}_c(x_c){\stackrel{·}{x}}_c)+\beta {\sigma }_c\hfill \\ =-{\mathit{W}}_c^{\mathrm{T}}\mathit{C}+\beta {\sigma }_c\hfill \end{array}$$

(68)

Substituting Equation (68) into Equation (66) delivers

$${\stackrel{·}{\stackrel{~}{\mathit{W}}}}_c=-{\mu }_c\mathit{C}{\mathit{C}}^{\mathrm{T}}{\tilde{\mathit{W}}}_c-\beta {\mu }_c{\sigma }_c\mathit{C}$$

(69)

Moreover, the dynamics of ${\tilde{\mathit{W}}}_a$ is given by

$${\stackrel{·}{\stackrel{~}{\mathit{W}}}}_a=-{\mu }_a({\hat{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a{\hat{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)){\mathit{S}}_a({\mathit{x}}_a)$$

(70)

Choosing the Lyapunov function candidate as

$$L_{\mathrm{ac}}=\frac{1}{2}{\tilde{\mathit{W}}}_c^{\mathrm{T}}{\tilde{\mathit{W}}}_c+\frac{1}{2}{\tilde{\mathit{W}}}_a^{\mathrm{T}}{\tilde{\mathit{W}}}_a$$

(71)

the derivative of $L_{\mathrm{ac}}$ along the trajectories of Equations (69) and (70) becomes

$$\begin{array}{c}\multicolumn{1}{c}{{\stackrel{·}{L}}_{\mathrm{ac}}={\tilde{\mathit{W}}}_c^{\mathrm{T}}{\stackrel{·}{\stackrel{~}{\mathit{W}}}}_c+{\tilde{\mathit{W}}}_a^{\mathrm{T}}{\stackrel{·}{\stackrel{~}{\mathit{W}}}}_a}\\ =-{\mu }_c{\tilde{\mathit{W}}}_c^{\mathrm{T}}(\mathit{C}{\mathit{C}}^{\mathrm{T}}{\tilde{\mathit{W}}}_c+\beta {\sigma }_c\mathit{C})\hfill \\ -{\mu }_a{\tilde{\mathit{W}}}_a^{\mathrm{T}}({\hat{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a)+{\overline{k}}_a{\hat{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)){\mathit{S}}_a(x_a)\hfill \\ =-{\mu }_c({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}{)}^2-\beta {\mu }_c{\sigma }_c{\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}-{\mu }_a({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a){)}^2\hfill \\ -{\mu }_a({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a))({\mathit{W}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a))-{\mu }_a{\overline{k}}_a({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a))\hfill \\ \times ({\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c))-{\mu }_a{\overline{k}}_a({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a))({\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c))\hfill \end{array}$$

(72)

By Young’s inequality, ⁴⁷ it follows that

$$\begin{array}{c}{\stackrel{·}{L}}_{\mathrm{ac}}⩽-\frac{{\mu }_c}{2}({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}+\beta {\sigma }_c{)}^2\hfill \\ -\frac{1}{2}({\mu }_c{({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C})}^2-{\mu }_a{\overline{k}}_a{({\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c))}^2)\hfill \\ -\frac{{\mu }_a(1-2{\overline{k}}_a)}{2}({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a){)}^2\hfill \\ +\frac{1}{2}({\mu }_c{\beta }^2{\overline{\sigma }}_c^2+{\mu }_a{({\mathit{W}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a))}^2\hfill \\ +{\mu }_a{\overline{k}}_a{({\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c))}^2)\hfill \end{array}$$

(73)

where ${\overline{\sigma }}_c$ denotes the upper bound of ${\sigma }_c$ . To show that ${\tilde{\mathit{W}}}_c$ and ${\tilde{\mathit{W}}}_a$ are bounded, we consider the following two cases:

• Case 1. ${\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)\ne 0$ . In this case, by setting ${\overline{k}}_a⩽\min \{1/4,({\mu }_c({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}{)}^2/2{\mu }_a({\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c){)}^2)\}$ , we can verify that

$$\begin{array}{c}{\stackrel{·}{L}}_{\mathrm{ac}}⩽-\frac{{\mu }_c}{4}({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}{)}^2-\frac{{\mu }_a}{4}({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a){)}^2+\Upsilon \hfill \\ ⩽-\frac{{\mu }_c}{4}({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}{)}^2-\frac{{\mu }_a}{4}({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a(x_a){)}^2+{\Upsilon }_{\max }\hfill \end{array}$$

(74)

where $\Upsilon =(1/2)({\mu }_c{\beta }^2{\overline{\sigma }}_c^2+{\mu }_a({\mathit{W}}_a^{\mathrm{T}}{\mathit{S}}_a({\mathit{x}}_a){)}^2+{\mu }_a{\overline{k}}_a({\mathit{W}}_c^{\mathrm{T}}{\mathit{S}}_c({\mathit{x}}_c){)}^2)$ and ${\Upsilon }_{\max }$ is a constant that satisfies $|\Upsilon |⩽{\Upsilon }_{\max }$ .

• Case 2. ${\tilde{\mathit{W}}}_c^{\mathrm{T}}{\mathit{S}}_c(x_c)=0$ . In this case, by setting ${\overline{k}}_a⩽(1/4)$ , we have

$${\stackrel{·}{L}}_{\mathrm{ac}}⩽-\frac{{\mu }_c}{2}({\tilde{\mathit{W}}}_c^{\mathrm{T}}\mathit{C}{)}^2-\frac{{\mu }_a}{4}({\tilde{\mathit{W}}}_a^{\mathrm{T}}{\mathit{S}}_a(x_a){)}^2+{\Upsilon }_{\max }$$

(75)

Summarizing the two cases, it can be verified that ${\stackrel{·}{L}}_{\mathrm{ac}}<0$ if $\parallel {\tilde{\mathit{W}}}_c\parallel >(4{\Upsilon }_{\max }/{\lambda }_{\min }(\mathit{C}{\mathit{C}}^{\mathrm{T}}))$ or $\parallel {\tilde{\mathit{W}}}_a\parallel >(4{\Upsilon }_{\max }/{\lambda }_{\min }({\mathit{S}}_a({\mathit{x}}_a){\mathit{S}}_a^{\mathrm{T}}({\mathit{x}}_a)))$ . A combination of the standard Lyapunov extension theorem ⁴⁸ and definition of $L_{\mathrm{ac}}$ entails that the weight estimation errors ${\tilde{W}}_c$ and ${\tilde{W}}_a$ are uniformly ultimately bounded, which completes the proof of Theorem 2.

Remark 7

The optimal cost function satisfies the Bellman equation according to optimal control theory, which has been addressed for general nonlinear systems in Buoniu et al. ³⁴ and the references therein. Since the exact solution of the Bellman equation is difficult to obtain, the NNs are generally adopted in the reinforcement learning technique to generate the approximate solution. It is noteworthy that the recursive relationship of the Bellman equation is critical for constructing the weight updating laws of the NNs, as can be seen from Equations (58) and (66).

Remark 8

We note that the stability analysis in the proof of Theorem 2 is mainly on the boundedness of weight estimation errors ${\tilde{W}}_c$ and ${\tilde{W}}_a$ . The convergence of the reinforcement learning method is generally presented by proving the boundedness of the weight estimation errors as done in Modares et al. ³⁵ and Mu et al. ⁴⁰ The stability proof of the closed-loop system with both ${\delta }_{e1}$ and ${\delta }_{e2}$ is very complicated since ${\delta }_{e1}$ and ${\delta }_{e2}$ are generated within different frameworks. We handle this issue in a similar way with the method in Mu et al., ⁴⁰ and this problem is worthy of further study.

Remark 9

We note that the prescribed performance technique adopted in this paper follows from Bechlioulis and Rovithakis, ²⁴ and many researchers have addressed prescribed performance control in recent years.^49–51 We do not deliberately emphasize on the specific type of prescribed performance control method. The focus of our work is a union of the prescribed performance control and the reinforcement learning approach, which is achieved by the basic part and supplementary part controllers, respectively. Nevertheless, we approve that the specific design of the prescribed performance control method is important to the final control performance, which will be the topic of our further research.

Numerical simulation

The effectiveness of the proposed control scheme is illustrated by a comparative simulation study with the method in Wu et al., ¹⁰ where the morphing aircraft model parameters are mainly taken from Wu et al. ⁹ and 20% uncertainties of aerodynamic coefficients are considered. Besides, we set $\Delta d_{\gamma }=2\sin (0.1t){\mathrm{e}}^{-0.1t}/(°/s)$ and $\Delta d_q=10\sin (0.2t){\mathrm{e}}^{-0.05t}/(°/{\mathrm{s}}^2)$ . The initial values of the states in Equation (2) are set as $[h_0,{\gamma }_0,{\alpha }_0,q_0]=[1000\mathrm{m},0°,0.4976°,0°/\mathrm{s}]$ . The reference altitude $h_r$ increases from $h_0$ to 1100 m through a transfer function $0.04/(s^2+0.4s+0.04)$ , and the sweep angle $\xi$ increases from $0°$ to $45°$ through a transfer function $0.01/(s^2+0.2s+0.01)$ . It is assumed that the system switching occurs every $5°$ within the $[0°,45°]$ range. The prescribed performance parameters are selected as $\underset{¯}{\varsigma }=1$ , $\overline{\varsigma }=1$ , ${\rho }_0=10$ , ${\rho }_{\infty }=2$ , and $a=0.1$ . For the basic part controller, $k_1=2$ , $k_2=4$ , $k_3=4$ , $k_4=6$ , ${\tau }_1=0.05$ , ${\tau }_2=0.05$ , and ${\tau }_3=0.05$ . The ESO parameters are taken as $l_{1\iota }=2$ , $l_{2\iota }=1$ , ${\epsilon }_{\iota }=0.05$ , ${\overline{k}}_{\iota }=1$ , and $\iota =\gamma ,q$ . The RBF NN ${\mathit{W}}_{1,i}^{\mathrm{T}}{\mathit{S}}_1({\mathit{X}}_1)$ has 100 nodes with the centers evenly spaced on $[25,30]\times [-0.3,0.3]\times [-1,1]\times [-0.3,0.3]\times [0,10]$ and the width is 1. The RBF NN ${\mathit{W}}_{2,i}^{\mathrm{T}}{\mathit{S}}_2({\mathit{X}}_2)$ has 150 nodes with the centers evenly spaced on $[25,30]\times [-0.3,0.3]\times [-0.3,0.3]\times [-0.3,0.3]\times [-5,5]\times [-0.3,0.3]\times [-0.3,0.3]$ , and the width is 1. The parameters for estimating ${\Theta }_1$ and ${\Theta }_2$ are set as $r_1=5$ , ${\kappa }_1=0.02$ , ${\epsilon }_1=0.5$ , $r_2=0.3$ , ${\kappa }_2=1$ , and ${\epsilon }_2=0.5$ , and the initial values are ${\hat{\Theta }}_1(0)=0$ and ${\hat{\Theta }}_2(0)=0$ .

For the reinforcement learning–based supplementary part, the critic network has 30 nodes with the centers evenly spaced on $[-0.3,0.3]$ and a width of 0.5, while the action network has 50 nodes with the centers evenly spaced on $[-0.3,0.3]$ and a width of 0.5. The other related parameters are taken as $\beta =0.1$ , ${\overline{k}}_r=10$ , ${\mu }_c=0.5$ , ${\overline{k}}_a=0.5$ , and ${\mu }_a=0.5$ . The initial values of both weights of the critic network and action network are set to zero, and the exploration noise is added in the first 3 s in the control input to enhance state exploration. The simulation results are shown in Figures 3–10, where the basic part control method is denoted by m1, the combination of basic part and supplementary part control method is denoted by m1 + RL, and the method in Wu et al. ¹⁰ is denoted by m2.

Figure 3.

Altitude tracking.

Figure 4.

Altitude tracking error.

Figure 5.

Altitude tracking error (m1 and m1 + RL).

Figure 6.

The control input ${\delta }_e$ .

Figure 7.

The response of ${\hat{\Theta }}_1$ .

Figure 8.

The response of ${\hat{\Theta }}_2$ .

Figure 9.

The responses of $d_{\gamma }$ and ${\hat{d}}_{\gamma }$ (m1).

Figure 10.

The response of $\parallel {\hat{\mathit{W}}}_c\parallel$ .

The altitude tracking responses of all three methods are depicted in Figures 3–5, and the boundedness of altitude tracking error can be guaranteed. However, the m2 method cannot ensure the prescribed performance bound, while the m1 and m1 + RL method can achieve the control goal. Besides, it is easy to see from Figure 5 that compared with the m1 method, the m1 + RL method has better transient performance and smaller steady error, which illustrates the effectiveness of the reinforcement learning.

The related control inputs are given in Figure 6, which are all less than $15°$ . It is shown that the related state variables are also bounded. Both the state variables and the control inputs of all three methods are of the same order. The responses of ${\hat{\Theta }}_1$ and ${\hat{\Theta }}_2$ are revealed in Figures 7 and 8, and they are bounded both for the m1 and m1 + RL method. The disturbance estimations in the $\gamma$ and $q$ loops can track the related disturbances with small errors even when the switching occurs for both the m1 method and m1 + RL method, which is depicted in Figure 9. In addition, the responses of $\parallel {\hat{\mathit{W}}}_c\parallel$ and $\parallel {\hat{\mathit{W}}}_a\parallel$ of the m1 + RL method are shown in Figures 10 and 11, and it can be seen that $\parallel {\hat{\mathit{W}}}_c\parallel$ and $\parallel {\hat{\mathit{W}}}_a\parallel$ stay bounded during the simulation.

Figure 11.

The response of $\parallel {\hat{\mathit{W}}}_c\parallel$ .

Remark 10

We have shown the advances of the work on reinforcement learning by comparative simulation studies. It is shown by the simulation results that the m1 method outperforms the m2 method in terms of prescribed performance of the altitude tracking, which indicates the advantage of the backstepping design via improved disturbance observer and dynamics surface technique. Furthermore, both the m1 method and m1 + RL method can ensure the prescribed performance bound, while the m1 + RL method has smaller tracking error when compared with the m1 method. Hence, the integration of reinforcement learning into the controller design can further improve control performance. It is noteworthy that the design philosophy of artificial intelligence community differs from that of control community although both of them contribute to the development of reinforcement learning greatly. ³⁴ The former cares about the convergence of the learning algorithm to the optimal solution, while the stability during the learning process is neglected. The latter focuses on the closed-loop stability, and the convergence to the optimal solution is not a priority. As a result, most of the extant reinforcement learning methods cannot be applied to systems with disturbances directly. For our work, it can be considered that the basic part controller is designed to provide a stable environment for the reinforcement learning scheme.

Conclusion

This paper addresses the prescribed performance control of the longitudinal altitude of a morphing aircraft based on switched nonlinear systems and reinforcement learning. Switched nonlinear systems in lower triangular form are first derived to describe the morphing aircraft dynamics, and the prescribed performance bound of the altitude tracking is considered by error transformation. Then, the controller is designed with integration of two parts: the basic part and supplementary part. The backstepping scheme is adopted to devise the basic part, and the modified dynamic surface control technique is involved to avoid the “explosion of complexity” problem. The common control laws are obtained by uniting improved disturbance observers and RBF NNs with consideration of the disturbances of the systems. Besides, a critic–action NN structure of the reinforcement learning is applied to develop the supplementary part controller via the error variables defined in the backstepping design. Comparative simulations clearly show that the altitude tracking error satisfies the prescribed performance bound by the proposed control method. In particular, the integration of reinforcement learning into the controller design can further improve control performance when compared with the controller which contains the basic part only.