Decentralized Robust Direct MRAC via E-Modification for the Pose of a Quadrotor UAV

1. Introduction

When the system output is influenced by uncontrolled inputs, and not only by the control input, such inputs create disturbances which may affect the system in unpredictable ways. In general, external perturbations may influence the output persistently although these may be present for a limited time.

In this work, in order to assure stability and to improve the performance when considering presence of external perturbations and modeling errors in the adaptive control of the pose for a quadrotor UAV, a modification to the parametric estimation algorithm is examined, which can leads to what is known as robust learning algorithms. A robust learning algorithm is characterized by preserving, within some specifications of design, stability properties under external perturbations and modeling errors [1,2,3]. It is well known that conventional adaptive (update) laws show parameter drift in presence of external perturbations, unmodeled dynamics, measurement noise and time variations. The parameter drift is a phenomenon in which the parametric estimates drift from their ideal values, possibly to infinity . As a result of the on-line parametric estimator, when any value for the parameters does not exist from a matching condition, in the attempt to find for a function an exact matching, it brings parameter drift. Two alternatives can be addressed in order to evade parameter drift. When the training error is too small, the first one, the so-called dead-zone modification [4], consists into indirectly not performing parametric adaptation. In the last one, whether the well-known e-modification [5,6], $\sigma$-modification [7], also known as leakage factor [1,8], and projection algorithms [1,8,9,10,11], the on-line parametric estimator is modified such that the parametric estimates are constrained from the drift to infinity.

The design of stable adaptive controllers is the main issue when trying with asymptotic model following. Unbounded solutions could emerge when trying with unmodeled dynamics and time-varying parameters as perturbations. From systems under the effects of the presence of perturbations, exact model following is no longer possible [12]. To this last eventuality, certain performance conditions must be satisfied through the controller design task having as the main goal the boundedness for all the signals from the closed-loop system.

It is well known that in the presence of measurement noise the parametric estimates do not converge to their ideal values if the input signal does not satisfy the persistency of excitation property, in other words, if it is not sufficiently rich, which could originate that these may also drift. The use of sufficiently rich inputs can prevent the parameter drift, although when such input signal is not available then update law modifications, namely, dead-zone modification, $\sigma$-modification, and e-modification, could be an alternative option.

The e-modification has been developed [5] in order to overcome a limitation with the $\sigma$-modification in that it can achieve asymptotic tracking under certain conditions while improving robustness. The aim of the e-modification is to reduce the damping term proportionally to the tracking error norm. The damping term is reduced to zero as the tracking error tends to zero, giving back the asymptotic tracking property of the ideal case with MRAC. Nevertheless, the e-modification only achieves bounded tracking in accordance with the stability analysis, so, the ideal property of the MRAC is not preserved. Moreover, asymptotic tracking is not attained in general with increased robustness.

An adaptive controller is said to be robust if it guarantees signal boundedness in the presence of bounded perturbations, unmodeled dynamics, and time delay [1,2,4,7,13].

In this work, a decentralized direct MRAC for the pose of a quadrotor UAV with e-modification in the update law, in order to evade parameter drift, here is proposed. The rotational dynamics from the quadrotor UAV is parameterized in the form of a decentralized model with matched uncertainties via control channels. This parametrization is preferred due to the matched uncertainties supply an extra realism to the ideal system. Also, MRACs are designed from each subsystem that constitutes the whole rotational system. Simulation results validate the proposal.

The manuscript is organized as follows. In Section 2 is described the nonlinear dynamic model of the quadrotor. The parameterized model for the quadrotor is given in Section 3. In Section 4, the design of the decentralized robust MRAC is developed. In Section 5, simulation results about the validation of the proposal are discussed. The conclusion is drawn in Section 6.

2. Quadrotor’s Dynamics

In this work, a four-rotor helicopter having one pair of rotors rotating clockwise and the remaining pair of rotors rotating counterclockwise is considered. Each rotor comprises a DC motor, a rotor blade and a drive gear. By increasing the speed of the rear rotor while simultaneously reducing the forward rotor by the same amount the forward motion is fulfilled. Backward, leftward, and rightward motion are fulfilled in the same manner. By accelerating the two clockwise turning rotors while decelerating the counterclockwise ones is performed the yaw motion.

The dynamics from a quadrotor is described through the equations of motion derived from the rigid body theory with six degrees of freedom [14], which can be decoupled into kinematic and dynamic equations [15]. The dynamic equations

$$\begin{array}{cc}\hfill F_{\mathrm{ext}}& =m_q{\dot{V}}_b+\omega \times m_q{V}_b,\hfill \end{array}$$

(1)

$$\begin{array}{cc}\hfill T_{\mathrm{ext}}& =J\dot{\omega }+\omega \times J\omega ,\hfill \end{array}$$

(2)

are derived from the center of mass, with $m_q$ denoting the quadrotor mass, $V_b$ is for the velocity in the body frame, $\omega$ is the angular rate of the quadrotor, $J=\mathrm{diag}\{J_x,J_y,J_z\}$ is the inertia matrix, the external force $F_{\mathrm{ext}}\in {\mathbb{R}}^3$ takes into account the total thrust, the quadrotor weight, and the aerodynamic force, whereas the external torque $T_{\mathrm{ext}}\in {\mathbb{R}}^3$ considers the aerodynamic moment vector and the difference of torque and thrust wielded by the two pairs of rotors. From a Z-Y-X sequence, consider the rotation matrix

$$\mathcal{R}\left(\begin{array}{ccc}\psi ,\hfill & \theta ,\hfill & \varphi \hfill \end{array}\right)=\left(\begin{array}{ccc}\mathrm{c}\theta \mathrm{c}\psi & \mathrm{c}\psi \mathrm{s}\theta \mathrm{s}\varphi -\mathrm{c}\varphi \mathrm{s}\psi & \mathrm{c}\varphi \mathrm{c}\psi \mathrm{s}\theta +\mathrm{s}\varphi \mathrm{s}\psi \\ \mathrm{c}\theta \mathrm{s}\psi & \mathrm{s}\theta \mathrm{s}\varphi \mathrm{s}\psi +\mathrm{c}\varphi \mathrm{c}\psi & \mathrm{c}\varphi \mathrm{s}\theta \mathrm{s}\psi -\mathrm{c}\psi \mathrm{s}\varphi \\ -\mathrm{s}\theta & \mathrm{c}\theta \mathrm{s}\varphi & \mathrm{c}\theta \mathrm{c}\varphi \end{array}\right)$$

(3)

with $\psi$, $\theta$, and $\varphi$ the Euler angles denoting yaw, pitch, and roll, respectively, where the $sin(·)$ and $cos(·)$ functions are abbreviated by $\mathrm{s}(·)$ and $\mathrm{c}(·)$, respectively. The angular velocity vector and the vector of Euler angles are related by

$$\left(\begin{array}{c}\dot{\psi }\\ \dot{\theta }\\ \dot{\varphi }\end{array}\right)=\mathcal{M}\begin{array}{ccc}(\psi ,& \theta ,& \varphi )\omega \end{array}$$

(4)

where $\psi \in [-\pi ,\pi ),\theta \in (-\pi /2,\pi /2),\varphi \in (-\pi /2,\pi /2)$, and

$$\mathcal{M}\begin{array}{ccc}(\psi ,& \theta ,& \varphi )\end{array}=\left(\begin{array}{ccc}\mathcal{0}& \mathrm{s}\varphi sec\theta & \mathrm{c}\varphi sec\theta \\ \mathcal{0}& \mathrm{c}\varphi & -\mathrm{s}\varphi \\ \mathcal{1}& \mathrm{s}\varphi tan\theta & \mathrm{c}\varphi tan\theta \end{array}\right.$$

(5)

is nonsingular inside the region where the Euler angles are allowed to take values.

The propeller dynamics is reduced to linear relationships between the force and moment of the propellers and the squared rotor angular rate given as

$$\begin{array}{cc}\hfill F_{zi}& =-b{\varpi }_i^2,\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill {\tau }_{zi}& =\kappa {\varpi }_i^2,\hfill \end{array}$$

(7)

with ${\varpi }_i$ the rotor angular rate, $\kappa$ and b denoting drag and thrust coefficients, respectively, where i is the motor index. All the motors are aligned with the vertical axis of the body frame, so, force and moments arise in the z-direction of the body frame, where $F_z$ is negative being consistent with the downward z-axis convention. The motor torque inputs are linked to the rigid body forces and moments through the mapping

$$\left(\begin{array}{c}F\\ {\tau }_x\\ {\tau }_y\\ {\tau }_z\end{array}\right)=\left(\begin{array}{cccc}-b& -b& -b& -b\\ 0& bl& 0& -bl\\ bl& 0& -bl& 0\\ \kappa & -\kappa & \kappa & -\kappa \end{array}\right)\left(\begin{array}{c}{\varpi }_1^2\\ {\varpi }_2^2\\ {\varpi }_3^2\\ {\varpi }_4^2\end{array}\right)$$

(8)

where l denotes the distance between the center of mass (center of gravity) and the rotors from the quadrotor. The resultant gyroscopic moment $M_g=\omega \times J_r{\varpi }_r$ is determined by the rotational velocity of the quadrotor $\omega$ and the sum of angular velocity of the rotors ${\varpi }_r=-{\varpi }_1+{\varpi }_2-{\varpi }_3+{\varpi }_4$. $J_r$ denotes the inertia of the rotating rotors, having the same parametric value for each one of the four rotors.

Asumming low speeds, the equations of motion for a quadrotor in terms of the Euler angles are then given by [16]

$$\begin{array}{cc}\hfill \ddot{z}& =-g+\left(\mathrm{c}\theta \mathrm{c}\varphi \right){\displaystyle \frac{F}{m_q}},\hfill \end{array}$$

(9)

$$\begin{array}{cc}\hfill \ddot{\theta }& =\left({\displaystyle \frac{J_z-J_x}{J_y}}\right)\dot{\varphi }\dot{\psi }+{\displaystyle \frac{J_r}{J_y}}{\varpi }_r\dot{\varphi }+{\displaystyle \frac{l}{J_y}}{\tau }_y,\hfill \end{array}$$

(10)

$$\begin{array}{cc}\hfill \ddot{\varphi }& =\left({\displaystyle \frac{J_y-J_z}{J_x}}\right)\dot{\theta }\dot{\psi }-{\displaystyle \frac{J_r}{J_x}}{\varpi }_r\dot{\theta }+{\displaystyle \frac{l}{J_x}}{\tau }_x,\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill \ddot{\psi }& =\left({\displaystyle \frac{J_x-J_y}{J_z}}\right)\dot{\theta }\dot{\varphi }+{\displaystyle \frac{1}{J_z}}{\tau }_z,\hfill \end{array}$$

(12)

where z is for the translational vertical (altitude) dynamics. $J_x$, $J_y$, and $J_z$ are moments of inertia in the direction of the three-dimensional Cartesian coordinates.

It can be seen that the rotational dynamics is intrinsically highly coupled by the angles of motion.

3. Decentralized Robust MRAC Design

3.1. Quadrotor’s Parametric Model

In this work, a decentralized robust MRAC for the pose of a quadrotor UAV is proposed.

Consider the set of state variables as given by

$$\begin{array}{cccc}{x}_1=\varphi ,\hfill & x_3=\theta ,\hfill & x_5=\psi ,\hfill & x_7=z,\hfill \\ x_2={\dot{x}}_1=\dot{\varphi },\hfill & x_4={\dot{x}}_3=\dot{\theta },\hfill & x_6={\dot{x}}_5=\dot{\psi },\hfill & x_8={\dot{x}}_7=\dot{z}.\hfill \end{array}$$

(13)

The description of the model (9)-(12) from the quadrotor in state–variable equations then takes the form

$$\dot{x}=\left(\begin{array}{c}{x}_2\\ \left({\displaystyle \frac{J_y-J_z}{J_x}}\right)x_4{x}_6-{\displaystyle \frac{J_r}{J_x}}{\varpi }_r{x}_4+{\displaystyle \frac{l}{J_x}}{\tau }_x\\ x_4\\ \left({\displaystyle \frac{J_z-J_x}{J_y}}\right)x_2{x}_6+{\displaystyle \frac{J_r}{J_y}}{\varpi }_r{x}_2+{\displaystyle \frac{l}{J_y}}{\tau }_y\\ x_6\\ \left({\displaystyle \frac{J_x-J_y}{J_z}}\right)x_2{x}_4+{\displaystyle \frac{1}{J_z}}{\tau }_z\\ x_8\\ -g+\left(\mathrm{c}{x}_1\mathrm{c}{x}_3\right){\displaystyle \frac{F}{m_q}}\end{array}\right).$$

(14)

Under the assumption that Euler angles change slightly during flying, the gyroscopic torques are neglected and the rotational subsystem can be partitioned in the form

$${\dot{x}}^1=\left(\begin{array}{c}{x}_2\\ \left({\displaystyle \frac{J_y-J_z}{J_x}}\right)x_4{x}_6+{\displaystyle \frac{l}{J_x}}{\tau }_x\end{array}\right),$$

(15)

$${\dot{x}}^2=\left(\begin{array}{c}{x}_4\\ \left({\displaystyle \frac{J_z-J_x}{J_y}}\right)x_2{x}_6+{\displaystyle \frac{l}{J_y}}{\tau }_y\end{array}\right),$$

(16)

$${\dot{x}}^3=\left(\begin{array}{c}{x}_6\\ \left({\displaystyle \frac{J_x-J_y}{J_z}}\right)x_2{x}_4+{\displaystyle \frac{1}{J_z}}{\tau }_z\end{array}\right),$$

(17)

with ${\dot{x}}^1={\left({\dot{x}}_1{\dot{x}}_2\right)}^T,{\dot{x}}^2={\left({\dot{x}}_3{\dot{x}}_4\right)}^T,$ and ${\dot{x}}^3={\left({\dot{x}}_5{\dot{x}}_6\right)}^T$.

In order to achieve the goal of pose control task for the quadrotor UAV via a decentralized MRAC design approach, (15) is parameterized in the form

$${\dot{x}}^1=A_1{x}^1+B_1{\mathsf{\Lambda }}_1\left(u^1+{\mathsf{\Theta }}_1{\mathsf{\Phi }}_1\left(x\right)\right)+{\zeta }_1$$

(18)

with $x^1={\left(x_1{x}_2\right)}^T$, $a_1=(J_y-J_z)/J_x$, $b_1=l/J_x$, ${\mathsf{\Lambda }}_1=b_1$, ${\mathsf{\Theta }}_1=a_1/b_1$, ${\mathsf{\Phi }}_1\left(x\right)=x_4{x}_6$, $u^1={\tau }_x$, and

$$\begin{array}{cc}\hfill A_1=\left[\begin{array}{cc}0& 1\\ 0& 0\end{array}\right],& B_1=\left[\begin{array}{c}0\\ 1\end{array}\right].\hfill \end{array}$$

Similarly, (16) takes the form

$${\dot{x}}^2=A_2{x}^2+B_2{\mathsf{\Lambda }}_2\left(u^2+{\mathsf{\Theta }}_2{\mathsf{\Phi }}_2\left(x\right)\right)+{\zeta }_2$$

(19)

with $x^2={\left(x_3{x}_4\right)}^T$, $a_2=(J_z-J_x)/J_y$, $b_2=l/J_y$, ${\mathsf{\Lambda }}_2=b_2$, ${\mathsf{\Theta }}_2=a_2/b_2,$${\mathsf{\Phi }}_2\left(x\right)=x_2{x}_6$, $u^2={\tau }_y$, $A_2=A_1$, and $B_2=B_1$. From (17),

$${\dot{x}}^3=A_3{x}^3+B_3{\mathsf{\Lambda }}_3\left(u^3+{\mathsf{\Theta }}_3{\mathsf{\Phi }}_3\left(x\right)\right)+{\zeta }_3$$

(20)

where $x^3={\left(x_5{x}_6\right)}^T$, $a_3=(J_x-J_y)/J_z$, $b_3=1/J_z$, ${\mathsf{\Lambda }}_3=b_3$, ${\mathsf{\Theta }}_3=a_3/b_3,$${\mathsf{\Phi }}_3\left(x\right)=x_2{x}_4$, $u^3={\tau }_z$, $A_3=A_1$, and $B_3=B_1$.

From the decentralized parametric model for the rotational dynamics of the quadrotor, it must be noticed that since the rotational dynamics is intrinsically highly coupled (interconnected) by the angles of motion, an interconnection matrix here is not defined in contrast with the decentralized adaptive control methodology from [7].

3.2. Robust MRAC Design from Unmatched External Perturbations

Consider a nonlinear MIMO system of the form [2]

$$\dot{x}=Ax+B\mathsf{\Lambda }(u+f\left(x\right)),$$

(21)

where $x\in {\mathbb{R}}^n$ is the system state, $B\in {\mathbb{R}}^{n\times m}$ is the known control matrix, $u\in {\mathbb{R}}^m$ is the control input, $A\in {\mathbb{R}}^{n\times n}$ and $\mathsf{\Lambda }\in {\mathbb{R}}^{m\times m}$ are unknown constant matrices. The pair ($A,B\mathsf{\Lambda }$) is controllable with $\mathsf{\Lambda }$ a diagonal matrix with elements ${\lambda }_i$ strictly positive. The uncertainty $\mathsf{\Lambda }$ is introduced for modeling errors or to model control failures. The unknown vector–function $f\left(x\right):{\mathbb{R}}^n\to {\mathbb{R}}^m$ represents the system matched uncertainty. $f\left(x\right)$ is a linear combination of N known locally Lipschitz–continuous basis functions ${\phi }_i\left(x\right)$ with unknown constant coefficients. So,

$$f\left(x\right)={\mathsf{\Theta }}^T\mathsf{\Phi }\left(x\right)$$

(22)

with $\mathsf{\Theta }\in {\mathbb{R}}^{N\times m}$ a constant matrix of the unknown coefficients and $\mathsf{\Phi }\left(x\right)={({\phi }_1\left(x\right){\phi }_2\left(x\right)\cdots {\phi }_N\left(x\right))}^T\in {\mathbb{R}}^N$ is the measurable signals (regressor) vector [3]. From (21), it is said that the uncertainties are matched in the sense that these enter to the system dynamics through control channels. So, the system controllability property is not affected as long as $\mathsf{\Lambda }$ is invertible. More even, the matched uncertainty assumption implies existence of at least one control solution capable of steering the system state along the desired trajectories [2]. The quantity (22) is called a parametric matched uncertainty since it appears in the range space of the control input matrix B. When a parametric uncertainty is matched, the control input can cancel out the uncertainty completely when the adaptation is perfect [13].

This study is centered in the design of a MIMO state feedback adaptive controller such that the system state x globally uniformly asymptotically follows the state $x_{\mathit{m}}\in {\mathbb{R}}^n$ of the reference model

$${\dot{x}}_{\mathit{m}}=A_{\mathit{m}}{x}_{\mathit{m}}+B_{\mathit{m}}r$$

(23)

where $A_{\mathit{m}}\in {\mathbb{R}}^{n\times n}$ is Hurwitz, $B_{\mathit{m}}\in {\mathbb{R}}^{n\times m}$, and $r\in {\mathbb{R}}^m$ is the bounded command vector. It is also required that during tracking, all signals in the closed–loop system remain uniformly bounded.

In general, there is no guarantee that ideal gains $K_x$ and $K_r$ satisfy the matching conditions

$$\begin{array}{cc}\hfill A+B\mathsf{\Lambda }{K}_x^T& =A_{\mathit{m}},\hfill \end{array}$$

(24)

$$\begin{array}{cc}\hfill B\mathsf{\Lambda }{K}_r^T& =B_{\mathit{m}},\hfill \end{array}$$

(25)

given $(A,B,\mathsf{\Lambda },A_{\mathit{m}},B_{\mathit{m}})$. Anyway, the pair ($A_{\mathit{m}},B_{\mathit{m}})$ can be selected such that (24)-() have one solution pair $(K_x,K_r)$ from the premise that the structure of A may be known.

Thus, given any bounded command, the control input needs to be chosen such that the tracking error

$$e=x-x_{\mathit{m}}$$

(26)

globally uniformly asymptotically tends to zero, i.e.,

$${\lim }_{t\to \infty }\parallel x-x_{\mathit{m}}\parallel =0.$$

(27)

From every subsystem (18),(19),(20) it is desirable that the closed-loop dynamics takes the form

$$\begin{array}{c}\hfill \dot{x}=A_{m}x+B\mathsf{\Lambda }(u+{\mathsf{\Theta }}^T\mathsf{\Phi }\left(x\right))+B_{m}r+\zeta \left(t\right)\end{array}$$

(28)

with $\zeta \left(t\right)$ a perturbation term. Thus, selecting the control input

$$u=-{\widehat{\mathsf{\Theta }}}^T\mathsf{\Phi }\left(x\right)$$

(29)

and substituting (29) in (28) results

$$\dot{x}=A_{m}x-B\mathsf{\Lambda }{\tilde{\mathsf{\Theta }}}^T\mathsf{\Phi }\left(x\right)+B_{m}r+\zeta \left(t\right)$$

(30)

with

$${\tilde{\mathsf{\Theta }}}^T={\widehat{\mathsf{\Theta }}}^T-{\mathsf{\Theta }}^T$$

(31)

the parametric estimation error, and

$$e=x-x_m$$

(32)

the tracking error equation.

Subtracting the reference model dynamics (23) from (28) yields the tracking error dynamics

$$\dot{e}=A_{m}e-B\mathsf{\Lambda }{\tilde{\mathsf{\Theta }}}^T\mathsf{\Phi }\left(x\right)+\zeta \left(t\right).$$

(33)

To guarantee uniform ultimate boundedness (UUB) of all signals in the closed-loop system consider the Lyapunov function candidate

$$V(e,\tilde{\mathsf{\Theta }})=e^{T}Pe+\mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T{\mathsf{\Gamma }}_{\mathsf{\Theta }}^{-1}\tilde{\mathsf{\Theta }}\mathsf{\Lambda }\right),$$

(34)

with ${\mathsf{\Gamma }}_{\mathsf{\Theta }}={\mathsf{\Gamma }}_{\mathsf{\Theta }}^T>0$ the adaptation gain matrix and $P=P^T>0$ the symmetric positive-definite solution from the Lyapunov equation. Its time derivative along the tracking error dynamics (33) results

$$\dot{V}(e,\tilde{\mathsf{\Theta }})=-e^{T}Qe-2e^{T}PB\mathsf{\Lambda }{\tilde{\mathsf{\Theta }}}^T\mathsf{\Phi }\left(x\right)+2e^{T}P\zeta \left(x\right)+2\mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T{\mathsf{\Gamma }}_{\mathsf{\Theta }}^{-1}\dot{\widehat{\mathsf{\Theta }}}\mathsf{\Lambda }\right),$$

(35)

which can be rewritten as

$$\dot{V}(e,\tilde{\mathsf{\Theta }})=-e^{T}Qe+2\mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\left\{{\mathsf{\Gamma }}_{\mathsf{\Theta }}^{-1}\dot{\widehat{\mathsf{\Theta }}}-\mathsf{\Phi }{e}^{T}PB\right\}\mathsf{\Lambda }\right)+2e^{T}P\zeta \left(t\right).$$

(36)

So, selecting the adaptive law

$$\dot{\widehat{\mathsf{\Theta }}}={\mathsf{\Gamma }}_{\mathsf{\Theta }}\mathsf{\Phi }\left(x\right)e^{T}PB$$

(37)

the time derivative (36) is then simplified to

$$\dot{V}(e,\tilde{\mathsf{\Theta }})=-e^{T}Qe+2e^{T}P\zeta \left(t\right)\le -{\lambda }_{\min }\left\{Q\right\}{\parallel e\parallel }^2+2\parallel e\parallel {\lambda }_{\max }\left\{P\right\}{\zeta }_{\max }$$

(38)

and, consequently, $\dot{V}<0$ outside of the set

$$E_0=\left\{(e,\tilde{\mathsf{\Theta }}):\parallel e\parallel \ge 2{\displaystyle \frac{{\lambda }_{\max }\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }=e_0\right\}.$$

(39)

The trajectories $e\left(t\right)$ of the error dynamics (33) enter a compact set $({\Omega }_0\supset E_0)\subset R^n$ in finite time and remain there for all future times [17]. Nevertheless, ${\Omega }_0$ is not compact in the $(e,\tilde{\mathsf{\Theta }})$ space. Certainly, because the parametric estimation error $\tilde{\mathsf{\Theta }}$ is not constrained at all, ${\Omega }_0$ is unbounded. Consequently, $\dot{V}$ can become positive inside ${\Omega }_0$ and the parametric error $\tilde{\mathsf{\Theta }}$ may diverge, even though the tracking error norm stays finite at all times [2]. This phenomenon is known as parameter drift induced by the perturbation $\zeta \left(t\right)$. Hence, this argument claims that the update law (37) is not robust to bounded perturbations.

3.3. Robust MRAC Design via e-Modification

To gain robustness, from (37) the update law

$$\dot{\widehat{\mathsf{\Theta }}}={\mathsf{\Gamma }}_{\mathsf{\Theta }}(\mathsf{\Phi }\left(x\right)e^{T}PB-\sigma \parallel e^{T}PB\parallel \widehat{\mathsf{\Theta }})$$

(40)

is selected, which includes an extra term adding a tracking error-dependent damping, which will tends to zero as the tracking error signal decreases, to the adaptive dynamics with $\sigma$, a design parameter, strictly positive constant.

Substituting (40) in (36), it yields to

$$\dot{V}(e,\tilde{\mathsf{\Theta }})=-e^{T}Qe-2\sigma \parallel e^{T}PB\parallel \mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\widehat{\mathsf{\Theta }}\mathsf{\Lambda }\right)+2e^{T}P\zeta \left(t\right).$$

(41)

From (31), (41) can be rewritten as

$$\dot{V}(e,\tilde{\mathsf{\Theta }})=-e^{T}Qe-2\sigma \parallel e^{T}PB\parallel \mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\tilde{\mathsf{\Theta }}\mathsf{\Lambda }\right)-2\sigma \parallel e^{T}PB\parallel \mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\mathsf{\Theta }\mathsf{\Lambda }\right)+2e^{T}P\zeta \left(t\right).$$

(42)

By definition,

$$\mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\tilde{\mathsf{\Theta }}\mathsf{\Lambda }\right)=\sum _{i=1}^N\sum _{j=1}^m{\tilde{\mathsf{\Theta }}}_{ij}^2{\mathsf{\Lambda }}_{ii}\ge {\parallel \tilde{\mathsf{\Theta }}\parallel }_F^2{\mathsf{\Lambda }}_{\min }$$

(43)

with ${\mathsf{\Lambda }}_{\min }$ and $\parallel \tilde{\mathsf{\Theta }}{\parallel }_F^2={\sum }_{i=1}^N{\sum }_{j=1}^m{\tilde{\mathsf{\Theta }}}_{ij}^2$ the minimum diagonal element of $\mathsf{\Lambda }$ and the Frobenius norm, respectively. From the Schwarz inequality,

$$|\mathrm{tr}\left({\tilde{\mathsf{\Theta }}}^T\mathsf{\Theta }\mathsf{\Lambda }\right)|\le \parallel {\tilde{\mathsf{\Theta }}}^T{\mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F\le \parallel \tilde{\mathsf{\Theta }}{\parallel }_F{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F.$$

(44)

Replacing (43) and (44) into (42) yields

$$\begin{array}{cc}\hfill \dot{V}(e,\tilde{\mathsf{\Theta }})\le & -{\lambda }_{\min }{\left\{Q\right\}\parallel e\parallel }^2+2\parallel e\parallel {\lambda }_{\max }\left\{P\right\}{\zeta }_{\max }-2\sigma \parallel e\parallel \parallel PB\parallel \parallel \tilde{\mathsf{\Theta }}{\parallel }_F^2{\mathsf{\Lambda }}_{\min }\hfill \\ \hfill & +2\sigma \parallel e\parallel \parallel PB\parallel \parallel \tilde{\mathsf{\Theta }}{\parallel }_F{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F.\hfill \end{array}$$

(45)

By completion of squares, the first terms can be written as

$$\begin{array}{c}-{\lambda }_{\min }{\left\{Q\right\}\parallel e\parallel }^2+2\parallel e\parallel {\lambda }_{\max }\left\{P\right\}{\zeta }_{\max }\hfill \\ \hfill =-{\lambda }_{\min }\left\{Q\right\}{\left(\parallel e\parallel -{\displaystyle \frac{{\lambda }_{\max }\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }\right)}^2+{\displaystyle \frac{{\lambda }_{\max }^2\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }^2.\end{array}$$

(46)

Also, the remaining terms, by completion of squares, can be written as

$$\begin{array}{c}-2\sigma \parallel e\parallel \parallel PB\parallel \parallel \tilde{\mathsf{\Theta }}{\parallel }_F^2{\mathsf{\Lambda }}_{\min }+2\sigma \parallel e\parallel \parallel PB\parallel \parallel \tilde{\mathsf{\Theta }}{\parallel }_F{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F\hfill \\ \hfill =-2\sigma \parallel e\parallel \parallel PB\parallel {\mathsf{\Lambda }}_{\min }{\left(\parallel \tilde{\mathsf{\Theta }}{\parallel }_F-{\displaystyle \frac{1}{2}}{\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F}{{\mathsf{\Lambda }}_{\min }}}\right)}^2+{\displaystyle \frac{1}{2}}\sigma \parallel e\parallel \parallel PB\parallel {\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F^2{\parallel \mathsf{\Lambda }\parallel }_F^2}{{\mathsf{\Lambda }}_{\min }}}.\end{array}$$

(47)

Then, (45) can be rewritten as

$$\begin{array}{cc}\hfill \dot{V}(e,\tilde{\mathsf{\Theta }})\le & -{\lambda }_{\min }\left\{Q\right\}{\left(\parallel e\parallel -{\displaystyle \frac{{\lambda }_{\max }\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }\right)}^2+{\displaystyle \frac{{\lambda }_{\max }^2\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }^2\hfill \\ \hfill & -2\sigma \parallel e\parallel \parallel PB\parallel {\mathsf{\Lambda }}_{\min }{\left(\parallel \tilde{\mathsf{\Theta }}{\parallel }_F-{\displaystyle \frac{1}{2}}{\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F}{{\mathsf{\Lambda }}_{\min }}}\right)}^2+{\displaystyle \frac{1}{2}}\sigma \parallel e\parallel \parallel PB\parallel {\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F^2{\parallel \mathsf{\Lambda }\parallel }_F^2}{{\mathsf{\Lambda }}_{\min }}}.\hfill \end{array}$$

(48)

Consequently, $\dot{V}(e,\tilde{\mathsf{\Theta }})<0$ if

$$\parallel e\parallel >{\displaystyle \frac{1}{2}}\sigma \parallel PB\parallel {\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F^2{\parallel \mathsf{\Lambda }\parallel }_F^2}{{\lambda }_{\min }\left\{Q\right\}{\mathsf{\Lambda }}_{\min }}}+2{\displaystyle \frac{{\lambda }_{\max }\left\{P\right\}}{{\lambda }_{\min }\left\{Q\right\}}}{\zeta }_{\max }=c_1$$

(49)

or

$$\parallel \tilde{\mathsf{\Theta }}{\parallel }_F>\sqrt{{\displaystyle \frac{{\lambda }_{\max }\left\{P\right\}}{\sigma \parallel PB\parallel {\mathsf{\Lambda }}_{\min }}}{\zeta }_{\max }+{\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F^2{\parallel \mathsf{\Lambda }\parallel }_F^2}{4{\mathsf{\Lambda }}_{\min }^2}}}+{\displaystyle \frac{1}{2}}{\displaystyle \frac{{\parallel \mathsf{\Theta }\parallel }_F{\parallel \mathsf{\Lambda }\parallel }_F}{{\mathsf{\Lambda }}_{\min }}}=c_2$$

(50)

Hence, outside of the compact set

$$\Omega =\{(e,\tilde{\mathsf{\Theta }}):(\parallel e\parallel \le c_1)\wedge (\parallel \tilde{\mathsf{\Theta }}{\parallel }_F\le c_2\left)\right\}\subset R^n\times R^{N\times m},$$

(51)

closed and bounded, so, $\dot{V}(e,\tilde{\mathsf{\Theta }})<0$.

Thus, (51) proves UUB tracking of the external command r by the output $y\left(t\right)$. So, in the presence of bounded time-varying perturbations $\zeta \left(t\right)$ and parametric uncertainties $(\mathsf{\Lambda },\mathsf{\Theta })$, command tracking is accomplished. Moreover, the command tracking problem for the MIMO system dynamics (28) is solved, which can be formulated as is shown next.

Theorem 1. Given the MIMO dynamics (18) with control uncertainty ${\mathsf{\Lambda }}_1$ and matched unknown function $f_1\left(x\right)={\mathsf{\Theta }}_1^T{\mathsf{\Phi }}_1\left(x\right)$, the MRAC subsystem (29), (31), (32), (40) imposes uniform ultimate bounded tracking of the reference model (23), driven by any bounded time–varying command $r\left(t\right)$. Also, all signals from the respective closed–loop subsystem are UUB in time.

4. Simulation Results

Simulation results about the performance of the decentralized robust MRAC are argued. Both absence and presence of perturbations are considered when evaluating the proposal. A Gaussian noise source is considered as the source of perturbations. The command signal to be tracked is exacting and the simulation time is extended in order to search for the parameter drift. From [18,19], with $J_x=J_y=0.03\mathrm{kg}{\mathrm{m}}^2$, $J_z=0.04\mathrm{kg}{\mathrm{m}}^2$, $l=0.2\mathrm{m}$, $m=1.79\mathrm{kg}$, and selecting ${\omega }_n^2=32.653\mathrm{rad}/\mathrm{s}$ and $\xi =0.7$, the matching conditions (Eqs. (24), ()) are then satisfied to ideal gains $K_{x^1}={(-4.89795-1.2)}^T$ and $K_{r^1}=4.89795$ [20].

The symmetric positive-definite matrix

$$P=\left(\begin{array}{cc}20.5931& 0.0153\\ 0.0153& 0.6269\end{array}\right)$$

(52)

results from solving the Lyapunov equation

$$P{A}_{\mathit{m}}+A_{\mathit{m}}^{T}P=-Q$$

(53)

for

$$Q=\left(\begin{array}{cc}1& 0\\ 0& 10\end{array}\right).$$

(54)

The adaptation gain for the update law (40) is chosen as ${\mathsf{\Gamma }}_{{\mathsf{\Theta }}_1}=9500$. The same design procedure is followed for pitch and yaw angles subsystems, so, only the design procedure for roll angle stabilization is given. The same P matrix is obtained when considering the same values for ${\omega }_n$, $\xi$, and Q matrix. However, the remaining adaptation rates are chosen as ${\mathsf{\Gamma }}_{{\mathsf{\Theta }}_2}=9430$ and ${\mathsf{\Gamma }}_{{\mathsf{\Theta }}_3}=9500$.

It must be noticed that a PID controller is implemented for altitude control [21,22].

Simulation results from Figure 1, Figure 2, Figure 3, Figure 4 and Figure 5 show the performance of the decentralized robust MRAC in absence of perturbations when performing model reference tracking to the pose of the quadrotor UAV. In this work, it has been corroborated that in absence of perturbations the parametric estimates tend to their ideal values in contrast with that proposal in [23], where the parametric estimates from the adaptive law with $\sigma$-modification tend to the origin. Simulation results from Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13 show that the decentralized robust MRAC is capable to confronting the perturbations when performing the task of model reference following to the quadrotor UAV’s pose. Figure 14 and Figure 15 show the parametric estimates dynamics from the adaptive law with e-modification (40) but from different values for $\sigma$ when the quadrotor UAV is under the influence of external perturbations. It should be noticed that the results about the model reference following by the decentralized direct MRAC for these latter cases are not included since these are very similar to those shown from Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13.

The norms and bounds of the error signal from the decentralized system are exhibited from Table 1, Table 2 and Table 3 for different values of $\sigma$ in the update law (40). So, from Table 1, Table 2 and Table 3, it can be seen that the inequality (49) is satisfied. Then, from Theorem 1, the robust model reference tracking problem under external unmatched perturbations for the MIMO system is solved.

5. Conclusions

It can be seen that, from simulation results, the parametric estimates from the adaptive laws converge to their ideal values in absence of perturbations. By the contrary, in the presence of external perturbations the parametric estimates tend to the origin, i.e., the parametric estimates from the adaptive laws do not drift to infinity when the decentralized MRAC for the pose of the quadrotor UAV has to confronting external perturbations. From the above, it can be concluded that the decentralized MRAC is robust. Also, the parametric estimation errors (31) do not converge to the origin when the decentralized MRAC is under the influence of external perturbations. Moreover, it must be noticed that the quadrotor UAV arrives to the reference for the altitude in both scenarios.