Generalized Proportional Integral Control for an Unmanned Quadrotor System

Abstract

In this article, a generalized proportional integral (GPI) control approach is presented for regulation and trajectory tracking problems in a nonlinear, multivariable quadrotor system model. In the feedback control law, no asymptotic observers or time discretizations are needed in the feedback loop. The GPI controller guarantees the asymptotically and exponentially stable behaviour of the controlled quadrotor position and orientation, as well as the possibilities of carrying out trajectory tracking tasks. The simulation results presented in the paper show that the proposed method exhibits very good stabilization and tracking performance in the presence of atmospheric disturbances and noise measurements.

1. Introduction

Unmanned aerial vehicles (UAVs) have attracted considerable interest for a wide variety of applications, including meteorological observation, fire monitoring and patrolling, to military purposes such as reconnaissance, monitoring and communication [1]. UAVs can be classified into four main categories based in their aerodynamic configuration [2]: (a) fixed-wing UAVs, which require a runway to land and take-off; (b) rotary-wing UAVs, which are able to land and take-off vertically and which have high manoeuvrability; (c) blimp UAVs, which appear like balloons or airships and which ensure lift by their helium-filled body; and finally (d) flapping-wing UAVs, which are inspired by flying insects and which can perform vertical take-off and landing. Nevertheless, rotary-wing platforms have captured a lot of attention in research projects because they present a number of advantages with regard to other UAV platforms. Their high manoeuvrability and ability to vertically take-off and land, as well as the capacity to fly in tough conditions to reach specified areas, make them ideal vehicles for these applications. Among rotary-wing UAVs, a new sub-classification can be considered [3]: (a) a single rotor is composed by a main rotor on the top and another rotor at the tail to achieve stability (a similar configuration to a helicopter); (b) a coaxial presents two propellers mounted in the same shaft rotating in opposite directions; (c) a quadrotor consists of four rotors fitted in a cross-like configuration; and (d) a multi-rotor consists of six or eight rotors. They are very agile and are able to fly even when a rotor does not work properly (there is redundancy due to the large number of rotors). We should note that, according to this classification, the quadrotor configuration has been the most widely used.

The quadrotor platform has four powerful rotors, each one of which has independent rotational speed, mounted in a square formation equidistant from the centre. The variation in the speed of the rotors generates the thrust and acceleration in the desired direction. Among the advantages of this platform are low cost, usability and ease of transportation, and it is also able to move at low speeds to ensure good quality images. However, despite the advantages of quadrotors with respect to other UAV platforms, the control of a quadrotor is a challenge due to its high manoeuvrability, its highly coupled multivariable dynamics of a nonlinear nature, and its underactuated condition, taking into account that it has six degrees of freedom (three for position and three for attitude) and only four rotors. For this reason, control techniques for these UAV platforms have witnessed rapidly expanding research to achieve not only autonomous hovering and orientation but also trajectory tracking [4]. Zuo [5] studied the command-filtered backstepping technique in order to stabilize a quadrotor's attitude without calculating analytically the pseudo-control signal derivative, as well as to decrease the dependent degrees of the analytic model. In Bouabdallah et al. [6], the application of two different approaches is presented, namely a PID approach assuming a simplified dynamics, and the LQ technique based on a more complete model. La Civita [7] proposed a robust control approach combined with linear rotorcraft models. Madani and Benallegue [8] presented a backstepping control strategy taking into consideration that the quadrotor can be divided into three subsystems: an under-actuated subsystem, a fully-actuated subsystem, and a propeller subsystem. Waslander et al. [9] presented two design approaches-integral sliding modes and reinforcement learning – for the altitude control loop. Both techniques resulted in stable controllers with similar response times, showing a significant improvement over linear controllers (which failed to stabilize the system adequately). Formentin and Lovera [10] developed a control law based on the differential flatness property of the position dynamics and the linearization of the system via feedforward and a passivity-based scheme. Furthermore, Gautam and Ha proposed in [11] a self-tuning fuzzy PID controller based on an EKF algorithm for the attitude and position control of a quadrotor. In a recent paper, Chen et al. described in [12] a reconfiguration control scheme for a quadrotor helicopter with actuator faults via adaptive control and combined multiple models. Sira-Ramrez studied in [13] an active disturbance rejection control scheme for efficient regulation and trajectory tracking tasks in a nonlinear, multivariable quadrotor system model. Escareño et al. [14] developed a nonlinear control strategy based on nested saturations that stabilizes the state of the quadrotor around the origin.

Recently, generalized proportional integral (GPI) controllers have demonstrated good performance in the control of nonlinear systems. GPI control has been found to present a better dynamic response than PID control in terms of the settling time while exhibiting a greater degree of robustness regarding disturbance rejection [15, 16]. GPI control sidesteps the need for traditional asymptotic state observers and directly proceeds to use, in a previously designed state feedback control law, structural state estimates in place of the actual state variables. These structural estimates are based on integral reconstructors and require only inputs, outputs and iterated integrals of such available signals. The effect of the neglected initial states is suitably compensated by means of a sufficiently large number of additional iterated integral output errors, integral input errors and control actions (see [17] for the relevant theoretical basis and [18–20] for the application of these ideas in diverse fields including laboratory experiments).

In this work, we extend the GPI control technique for both stabilization and trajectory tracking tasks of an unmanned quadrotor system. In the control law, neither asymptotic observers nor time discretizations are therefore needed in the feedback loop for the estimation of the states commonly required in traditional state-based feedback controllers for such systems.

The article is structured as follows. Section 2 presents the quadrotor model, and the problem to be solved is formulated. Section 3 establishes the flatness system of the quadrotor model. Next, the GPI controller to be used in the control of the unmanned quadrotor system is derived. In this section, it is proved that the GPI controller produces asymptotically, exponentially convergent tracking-error behaviour in relation to the origin of the coordinates in the error space. Section 4 depicts digital computer simulations showing the performance of the GPI controller and, finally, Section 5 presents the conclusions of the work.

2. Quadrotor Dynamics and Problem Formulation

A quadrotor is an underactuated aircraft with four rotors. The rotors are directed upwards and they are placed in a square formation at an equal distance from the centre of mass of the quadrotor, as shown in Figure 1.

Figure 1.

Quadrotor motion principle

In the quadrotor, there are four rotors with fixed angles that are basically the thrust generated by each propeller. The change of speed in the propellers modifies the lift forces. The pitch movement is obtained by increasing (reducing) the speed of the propeller (1) while reducing (increasing) the speed of the propeller (3). The roll movement is obtained similarly by increasing (reducing) the speed of the propeller (2) while reducing (increasing) the speed of the propeller (4). The yaw movement is achieved by increasing (decreasing) the speed between each pair of propellers.

The dynamic model of the quadrotor can be achieved through the Euler-Lagrange formalism in [22, 23],

$$m\ddot{x}=-u{S}_{\theta }$$

(1)

$$m\ddot{y}=u{C}_{\theta }{S}_{\varphi }$$

(2)

$$m\ddot{z}=u{C}_{\theta }{C}_{\varphi }-mg$$

(3)

$$\ddot{\psi }={\tau }_{\psi }$$

(4)

$$\ddot{\theta }={\tau }_{\theta }$$

(5)

$$\ddot{\varphi }={\tau }_{\varphi }$$

(6)

where S_θ ≡ sinθ, C_θ ≡ cosθ, S_φ ≡ sinφ, C_φ ≡ cosφ, m is the mass, g is the gravity acceleration, x and y are coordinates in the horizontal plane, z is the vertical position, the angles φ, θ and ψ express the independent orientation angles, u is defined as the total thrust, and τ_ψ, τ_θ and τ_φ denote the angular moments (yawing moment, pitching moment and rolling moment, respectively).

Assumption 1: The orientation angles θ and ϕ are upper- and lower-bounded in the intervals $-\frac{\pi }{2}<\varphi <\frac{\pi }{2}$ and $-\frac{\pi }{2}<\theta <\frac{\pi }{2}$ .

After defining the system dynamics, the problem formulation studied in this work is now stated:

Given a set of smooth reference trajectories $\left(x^{*},y^{*},z^{*},{\psi }^{*}\right)$ , devise a feedback control such that the horizontal coordinates x and y, the vertical position z, and the orientation variable ψ, asymptotically track the given references so that the tracking-error trajectories are ultimately confined to a small neighbourhood of the origin of the tracking-error phase space and the variables θ and φ are confined to move in the interval $(-\frac{\pi }{2},\frac{\pi }{2})$ .

3. Control Design

3.1 Flatness of the System

According to the theory of differential flatness [24], a dynamic system, $\dot{x}=f(x,u)$ , $y=h(x)$ , with $x\in R^n$ , $u\in R^m$ and $y\in R^p$ , is said to be differentially flat if there exist m differentially independent variables, called flat outputs (by differentially independent is meant that they are not related by any differential equation), which are functions of the state vector and, possibly, of a finite number of time derivatives of the state vector (i.e., derivatives of the inputs may be involved in their definition), such that all the system variables (states, inputs, outputs, and functions of these variables) can, in turn, be expressed as functions of the flat outputs and of a finite number of their time derivatives. We have the following proposition:

Proposition 1: The quadrotor model given in (1)–(6) is differentially flat, with a flat output vector given by (x, y, z, ψ), i.e., all the system variables in (1)–(6) can be differentially parameterized solely in terms of the flat output vector components x, y, z and ψ, and a finite number of their time derivatives. Their expressions are:

$$u=m\sqrt{{\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}$$

(7)

$$\varphi =\arctan \left(\frac{\ddot{y}}{\ddot{z}+g}\right)$$

(8)

$$\theta =-\arctan \left(\frac{\ddot{x}}{\sqrt{{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}}\right)$$

(9)

$$\dot{\varphi }=\frac{y^{(3)}(\ddot{z}+g)-\ddot{y}{z}^{(3)}}{{\ddot{y}}^2+{(\ddot{z}+g)}^2}$$

(10)

$$\begin{array}{l}\dot{\theta }=-\frac{1}{\left[{\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2\right]\sqrt{{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}}\\ \times \left[x^{(3)}\left[{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2\right]-\ddot{x}\left[\ddot{y}{y}^{(3)}+\left(\ddot{z}+g\right)z^{(3)}\right]\right]\end{array}$$

(11)

$${\tau }_{\psi }=\ddot{\psi }$$

(12)

$$\begin{array}{l}{\tau }_{\varphi }=\frac{y^{(4)}(\ddot{z}+g)-\ddot{y}{z}^{(4)}}{{\ddot{y}}^2+{(\ddot{z}+g)}^2}\\ -2\frac{\left(y^{(3)}(\ddot{z}+g)-\ddot{y}{z}^{(3)}\right)\left(\ddot{y}{y}^{(3)}+(\ddot{z}+g)z^{(3)}\right)}{{\left({\ddot{y}}^2+{(\ddot{z}+g)}^2\right)}^2}\end{array}$$

(13)

$$\begin{array}{l}{\tau }_{\theta }=\frac{\Pi \left(\ddot{x},\ddot{y},\ddot{z},x^{(3)},y^{(3)},z^{(3)},x^{(4)},y^{(4)},z^{(4)}\right)}{\left[{\ddot{x}}^2+{\ddot{y}}^2+{(\ddot{z}+g)}^2\right]\sqrt{{\ddot{y}}^2+{(\ddot{z}+g)}^2}}\\ +\frac{\left(x^{(3)}\left({\ddot{y}}^2+{(\ddot{z}+g)}^2\right)-\ddot{x}\left(\ddot{y}{y}^{(3)}+\left(\ddot{z}+g\right)z^{(3)}\right)\right)}{{\left({\ddot{x}}^2+{\ddot{y}}^2+{(\ddot{z}+g)}^2\right)}^2{\left({\ddot{y}}^2+{(\ddot{z}+g)}^2\right)}^{\frac{3}{2}}}\times \\ [2\left(\ddot{x}{x}^{(3)}+\ddot{y}{y}^{(3)}+(\ddot{z}+g)z^{(3)}\right)\left({\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2\right)\\ +\left[{\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2\right]\cdot \left[\ddot{y}{y}^{(3)}+\left({\ddot{z}}^2+g\right)z^{(3)}\right]]\end{array}$$

(14)

and

$$\begin{array}{l}\Pi \left(\ddot{x},\ddot{y},\ddot{z},x^{(3)},y^{(3)},z^{(3)},x^{(4)},y^{(4)},z^{(4)}\right)=-x^{(4)}\left[{\ddot{y}}^2+\left({\ddot{z}}^2+g\right)\right]\\ -x^{(3)}\left[\ddot{y}{y}^{(3)}+\left(\ddot{z}+g\right)z^{(3)}\right]+\ddot{x}\left[{\left(y^{(3)}\right)}^2+\ddot{y}{y}^{(4)}+{\left(z^{(3)}\right)}^2+\left(\ddot{z}+g\right)z^{(4)}\right]\end{array}$$

(15)

From expressions (7), (13) and (14), it can be seen that the relationship between the control input vector, (u,τ_φ,τ_θ,τ_ψ), and the flat output's highest derivatives, is not invertible. This reveals an obstacle in the input vector to achieve static feedback linearization, and points to the need for a second-order dynamic extension of the control input u in order to exactly linearize the system (see [25] for details on the use of dynamic feedback). This yields:

$$\begin{array}{l}\ddot{u}=m\frac{{\left(x^{(3)}\right)}^2+\ddot{x}{x}^{(4)}+{\left(y^{(3)}\right)}^2+\ddot{y}{y}^{(4)}+{\left(z^{(3)}\right)}^2\left(\ddot{z}+g\right)z^{(4)}}{\sqrt{{\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}}\\ -m\frac{{\left(\ddot{x}{x}^{(3)}+\ddot{y}{y}^{(3)}+\left(\ddot{z}+g\right)z^{(3)}\right)}^2}{{\left({\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2\right)}^{\frac{3}{2}}}\end{array}$$

(16)

Proof: By squaring expressions (1), (2) and (3), adding the resulting terms and rearranging, it follows that

$$u=m\sqrt{{\ddot{x}}^2+{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}$$

(17)

Now, from expressions (2) and (3), it is obtained that

$$\begin{array}{l}m\ddot{y}=u{C}_{\theta }{S}_{\varphi }\hfill \\ m(\ddot{z}+g)=u{C}_{\theta }{C}_{\varphi }\hfill \\ \hfill \end{array}\}\Rightarrow \varphi =\arctan \left(\frac{\ddot{y}}{\ddot{z}+g}\right)$$

(18)

and

$$S_{\varphi }=\frac{\ddot{y}}{\sqrt{{\ddot{y}}^2+{(\ddot{z}+g)}^2}};C_{\varphi }=\frac{\ddot{z}+g}{\sqrt{{\ddot{y}}^2+{(\ddot{z}+g)}^2}}$$

(19)

Then, by squaring solely expressions (2) and (3), adding the resulting terms and rearranging, yields the following result

$$u{C}_{\varphi }=m\sqrt{{\ddot{y}}^2+{(\ddot{z}+g)}^2}$$

(20)

Operating with expressions (1) and (19), we readily obtain

$$\theta =-\arctan \left(\frac{\ddot{x}}{\sqrt{{\ddot{y}}^2+{\left(\ddot{z}+g\right)}^2}}\right)$$

(21)

and

$$S_{\theta }=-\frac{\ddot{x}}{\sqrt{{\ddot{x}}^2+{\ddot{y}}^2+{(\ddot{z}+g)}^2}};C_{\theta }=\frac{\sqrt{{\ddot{y}}^2+{(\ddot{z}+g)}^2}}{\sqrt{{\ddot{x}}^2+{\ddot{y}}^2+{(\ddot{z}+g)}^2}}$$

(22)

Now, if the ψ angle is differentiated twice we arrive at

$${\tau }_{\psi }=\ddot{\psi }$$

(23)

On the other hand, if the expressions (1), (2) and (3) are differentiated with regard to the time and the terms are rearranged,

$$m{x}^{(3)}=-\dot{u}{S}_{\theta }-u{C}_{\theta }\dot{\theta }$$

(24)

$$m{y}^{(3)}=\dot{u}{C}_{\theta }{S}_{\varphi }-u{S}_{\theta }{S}_{\varphi }\dot{\theta }+u{C}_{\theta }{C}_{\varphi }\dot{\varphi }$$

(25)

$$m{z}^{(3)}=\dot{u}{C}_{\theta }{C}_{\varphi }-u{S}_{\theta }{C}_{\varphi }\dot{\theta }-u{C}_{\theta }{S}_{\varphi }\dot{\varphi }$$

(26)

and, expressed in matrix notation,

$$\left[\begin{array}{c}{x}^{(3)}\\ y^{(3)}\\ z^{(3)}\end{array}\right]=\underset{N^{-1}(u,\theta ,\varphi )}{\underbrace{\frac{1}{m}\left[\begin{array}{ccc}-S_{\theta }& -u{C}_{\theta }& 0\\ C_{\theta }{S}_{\varphi }& -u{S}_{\theta }{S}_{\varphi }& u{C}_{\theta }{C}_{\varphi }\\ C_{\theta }{C}_{\varphi }& -u{S}_{\theta }{C}_{\varphi }& -u{C}_{\theta }{S}_{\varphi }\end{array}\right]}}\left[\begin{array}{c}\dot{u}\\ \dot{\theta }\\ \dot{\varphi }\end{array}\right]$$

(27)

or

$$\left[\begin{array}{c}\dot{u}\\ \dot{\theta }\\ \dot{\varphi }\end{array}\right]=\underset{N(u,\theta ,\varphi )}{\underbrace{m\left[\begin{array}{ccc}-S_{\theta }& C_{\theta }{S}_{\varphi }& C_{\theta }{C}_{\varphi }\\ -\frac{C_{\theta }}{u}& -\frac{S_{\theta }{S}_{\varphi }}{u}& -\frac{S_{\theta }{C}_{\varphi }}{u}\\ 0& \frac{C_{\varphi }}{u{C}_{\theta }}& -\frac{S_{\varphi }}{u{C}_{\theta }}\end{array}\right]}}\left[\begin{array}{c}{x}^{(3)}\\ y^{(3)}\\ z^{(3)}\end{array}\right]$$

(28)

If the expression (26) is differentiated now with regard to the time, and the terms are rearranged taking into consideration that ${\tau }_{\theta }=\ddot{\theta }$ and ${\tau }_{\varphi }=\ddot{\varphi }$ , we obtain the following:

$$m{x}^{(4)}=-\ddot{u}{S}_{\theta }-2\dot{u}\dot{\theta }{C}_{\theta }+u{\dot{\theta }}^2{S}_{\theta }-u\ddot{\theta }{C}_{\theta }$$

(29)

$$\begin{array}{l}m{y}^{(4)}=\ddot{u}{C}_{\theta }{S}_{\varphi }-2\dot{u}\dot{\theta }{S}_{\theta }{S}_{\varphi }+2\dot{u}\dot{\varphi }{C}_{\theta }{C}_{\varphi }-2u\dot{\theta }\dot{\varphi }{S}_{\theta }{C}_{\varphi }\\ -u{\dot{\theta }}^2{C}_{\theta }{S}_{\varphi }-u{\dot{\varphi }}^2{C}_{\theta }{S}_{\varphi }-u\ddot{\theta }{S}_{\theta }{S}_{\varphi }+u\ddot{\varphi }{C}_{\theta }{C}_{\varphi }\end{array}$$

(30)

$$\begin{array}{l}m{z}^{(4)}=\ddot{u}{C}_{\theta }{C}_{\varphi }-2\dot{u}\dot{\theta }{S}_{\theta }{C}_{\varphi }-2\dot{u}\dot{\varphi }{C}_{\theta }{S}_{\varphi }+2u\dot{\theta }\dot{\varphi }{S}_{\theta }{S}_{\varphi }\\ -u{\dot{\theta }}^2{C}_{\theta }{C}_{\varphi }-u{\dot{\varphi }}^2{C}_{\theta }{C}_{\varphi }-u\ddot{\theta }{S}_{\theta }{C}_{\varphi }-u\ddot{\varphi }{C}_{\theta }{S}_{\varphi }\end{array}$$

(31)

and, expressed in matrix notation,

$$\begin{array}{l}\left[\begin{array}{c}{x}^{(4)}\\ y^{(4)}\\ z^{(4)}\end{array}\right]=N^{-1}(u,\theta ,\varphi )\left[\begin{array}{c}\ddot{u}\\ {\tau }_{\theta }\\ {\tau }_{\varphi }\end{array}\right]\\ +\underset{T(u,\varphi ,\theta ,\dot{u},\dot{\varphi },\dot{\theta })}{\underbrace{\frac{1}{m}\left[\begin{array}{c}-2\dot{u}\dot{\theta }{C}_{\theta }+u{\dot{\theta }}^2{S}_{\theta }\\ -2\dot{u}\dot{\theta }{S}_{\theta }{S}_{\varphi }+2\dot{u}\dot{\varphi }{C}_{\theta }{C}_{\varphi }-2u\dot{\theta }\dot{\varphi }{S}_{\theta }{C}_{\varphi }-u{C}_{\theta }{S}_{\varphi }\left({\dot{\theta }}^2+{\dot{\varphi }}^2\right)\\ -2\dot{u}\dot{\theta }{S}_{\theta }{C}_{\varphi }-2\dot{u}\dot{\varphi }{C}_{\theta }{S}_{\varphi }+2u\dot{\theta }\dot{\varphi }{S}_{\theta }{S}_{\varphi }-u{C}_{\theta }{C}_{\varphi }\left({\dot{\theta }}^2+{\dot{\varphi }}^2\right)\end{array}\right]}}\end{array}$$

(32)

or

$$\left[\begin{array}{c}\ddot{u}\\ {\tau }_{\theta }\\ {\tau }_{\varphi }\end{array}\right]=N(u,\theta ,\varphi )\left[\begin{array}{c}{x}^{(4)}\\ y^{(4)}\\ z^{(4)}\end{array}\right]+\underset{B(u,\varphi ,\theta ,\dot{u},\dot{\varphi },\dot{\theta })}{\underbrace{\left[\begin{array}{c}u\left({\dot{\theta }}^2+{\dot{\varphi }}^2{C}_{\theta }^2\right)\\ -\frac{2\dot{u}\dot{\theta }}{u}-{\dot{\varphi }}^2{S}_{\theta }{C}_{\theta }\\ -\frac{2\dot{u}\dot{\varphi }}{u}+2\dot{\theta }\dot{\varphi }\frac{S_{\theta }}{C_{\theta }}\end{array}\right]}}$$

(33)

Finally, the proof is completed after substituting (16), (18), (21) and (27) in (32). □

3.2 GPI Controller Design

As was described in Section 1, GPI control consists of the defective integral reconstruction of the state which, a priori, neglects the effects of unknown initial conditions and the effect of possible classical perturbation inputs (constant and low-order time polynomial errors, such as ramps and parabolic signals). They are based on the central observation that the states of observable linear systems may be integrally parameterized in terms of inputs and outputs alone (i.e., linear combinations of inputs, outputs and of a finite number of iterated integrals of signals). The errors of integral reconstruction are to be compensated, later, by a suitable linear controller containing a sufficient number of iterated integrals of the tracking error or else of the input error [21]. The control scheme proposed for the control of the quadrotor model is illustrated in Figure 2.

Figure 2.

GPI control scheme

From the developments obtained in the previous section the following input-to-highest-derivative of the flat outputs' relations is achieved:

$$\begin{array}{l}\left[\begin{array}{c}\ddot{u}\\ {\tau }_{\theta }\\ {\tau }_{\varphi }\end{array}\right]=N\cdot \left[\begin{array}{c}{x}^{(4)}\\ y^{(4)}\\ z^{(4)}\end{array}\right]+B\\ {\tau }_{\psi }=\ddot{\psi }\end{array}$$

(34)

Now, if we define the virtual input vector,

$$\left[\begin{array}{c}{\Upsilon }_{\ddot{u}}\\ {\Upsilon }_{\theta }\\ {\Upsilon }_{\varphi }\end{array}\right]=N^{-1}\cdot \left[\left[\begin{array}{c}\ddot{u}\\ {\tau }_{\theta }\\ {\tau }_{\varphi }\end{array}\right]-B\right]$$

(35)

and this input transformation is applied to the dynamical system (33), the whole dynamics model is now expressed as

$$x^{(4)}={\Upsilon }_{\ddot{u}}$$

(36)

$$y^{(4)}={\Upsilon }_{\theta }$$

(37)

$$z^{(4)}={\Upsilon }_{\varphi }$$

(38)

$$\ddot{\psi }={\tau }_{\psi }$$

(39)

The GPI-based flat output feedback controller is synthesized as follows: Expression (35) is a fourth-order system in which is regulated the x-position of the quadrotor model towards a given smooth reference trajectory, $x^{*}(t)$ , with $\ddot{u}$ acting as an auxiliary control input. Clearly, if there exists an auxiliary open-loop control input ${ }_{\ddot{u}}^{*}(t)$ that ideally achieves the tracking of $x^{*}(t)$ for suitable initial conditions, it thus satisfies the fourth-order dynamics

$${\left[x^{*}\right]}^{(4)}={\Upsilon }_{\ddot{u}}^{*}$$

(40)

Subtracting (39) from (35) yields the following:

$$e_x^{(4)}=e_{{\Upsilon }_{\ddot{u}}}$$

(41)

where $e_x=x-x^{*}(t)$ and $e_{{ }_{\ddot{u}}}={ }_{\ddot{u}}-{ }_{\ddot{u}}^{*}(t)$ . If we assume that we are able to measure the variables $e_x^{(3)}$ , ${\ddot{e}}_x$ and ${\dot{e}}_x$ , the following control law could be proposed:

$$e_{{\Upsilon }_{\ddot{u}}}=-k_6^x{e}_x^{(3)}-k_5^x{\ddot{e}}_x-k_4^x{\dot{e}}_x-k_3^x{e}_x$$

(42)

In such a case, the closed-loop tracking error for the x-position variable evolves, and is governed by

$$e_x^{(4)}+k_6^x{e}_x^{(3)}+k_5^x{\ddot{e}}_x+k_4^x{\dot{e}}_x+k_3^x{e}_x=0$$

(43)

The design parameters $\{k_6^x,k_5^x,k_4^x,k_3^x\}$ are then chosen so as to render the closed-loop characteristic polynomial into a Hurwitz polynomial with desirable roots. The signals $e_x^{(3)}$ , ${\ddot{e}}_x$ and ${\dot{e}}_x$ have to be either measured or else estimated by means of an observer. In practice, the estimation requires the use of online calculations based on high-frequency samples of the variable's trajectory e_x(t). Using integral reconstructors based on the integration of the system dynamics, such estimations – or time discretizations – are unnecessary. We resort to the following integral reconstructors for such signals,

$${\widehat{e}}_x^{(3)}={\displaystyle {\int }_0^t}{e}_{{\Upsilon }_{\ddot{u}}}(\tau )d\tau$$

(44)

$${\widehat{\ddot{e}}}_x={\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{e}_{{\Upsilon }_{\ddot{u}}}(\lambda )d\lambda d\tau$$

(45)

$${\widehat{\dot{e}}}_x={\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{\displaystyle {\int }_0^{\lambda }}{e}_{{\Upsilon }_{\ddot{u}}}(\sigma )d\sigma d\lambda d\tau$$

(46)

The relation between the structural estimates and the real values of the states of the system are given by

$$e_x^{(3)}={\widehat{e}}_x^{(3)}+e_x^{(3)}(0)$$

(47)

$${\ddot{e}}_x={\widehat{\ddot{e}}}_x+e_x^{(3)}(0)t+{\ddot{e}}_x(0)$$

(48)

$${\dot{e}}_x={\widehat{\dot{e}}}_x+e_x^{(3)}(0)\frac{t^2}{2}+{\ddot{e}}_x(0)t+{\dot{e}}_x(0)$$

(49)

where $e_x^{(3)}(0)$ , ${\ddot{e}}_x(0)$ and ${\dot{e}}_x(0)$ are the unknown initial conditions of the states of the system. It is observed in (46)–(48) that the difference between the structural estimates and the real values exhibits a combination of parabolic error, ramp error and offset error due to the constant initial conditions. This immediately prompts us to consider the possibility of using a modified feedback controller including an integral error term, a double integral error term and a triple integral error term, to reject the possible effects of the possibly non-zero initial conditions $e_x^{(3)}(0)$ , ${\ddot{e}}_x(0)$ and ${\dot{e}}_x(0)$ . We thus proceed to propose the following controller:

$$\begin{array}{l}{e}_{{\Upsilon }_{\ddot{u}}}=-k_6^x{\widehat{e}}_x^{(3)}-k_5^x{\widehat{\ddot{e}}}_x-k_4^x{\widehat{\dot{e}}}_x-k_3^x{e}_x-k_2^x{\displaystyle {\int }_0^t}{e}_x(\tau )d\tau \\ -k_1^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{e}_x(\lambda )d\lambda d\tau -k_0^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{\displaystyle {\int }_0^{\lambda }}{e}_x(\sigma )d\sigma d\lambda d\tau \end{array}$$

(50)

Substituting (43)–(45) in (49), and after some rearrangements, we achieve

$$\begin{array}{l}{e}_{{\Upsilon }_{\ddot{u}}}=-k_6^x{\displaystyle {\int }_0^t}{e}_{{\Upsilon }_{\ddot{u}}}(\tau )d\tau -k_5^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{e}_{{\Upsilon }_{\ddot{u}}}(\lambda )d\lambda d\tau \\ -k_4^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{\displaystyle {\int }_0^{\lambda }}{e}_{{\Upsilon }_{\ddot{u}}}(\sigma )d\sigma d\lambda d\tau -k_3^x{e}_x-k_2^x{\displaystyle {\int }_0^t}{e}_x(\tau )d\tau \\ -k_1^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{e}_x(\lambda )d\lambda d\tau -k_0^x{\displaystyle {\int }_0^t}{\displaystyle {\int }_0^{\tau }}{\displaystyle {\int }_0^{\lambda }}{e}_x(\sigma )d\sigma d\lambda d\tau \end{array}$$

(51)

After using Laplace's transform in (50), we obtain the final form for the auxiliary control input variable $e_{{ }_{\ddot{u}}}$ ,

$$e_{{\Upsilon }_{\ddot{u}}}(s)=-\frac{k_3^x{s}^3+k_2^x{s}^2+k_1^{x}s+k_0^x}{s^3+k_6^x{s}^2+k_5^{x}s+k_4^x}{e}_x(s)$$

(52)

On the other hand, the use of (46)–(48) on the modified controller expression (49) results from the substitution of expression (40) and differentiating on three occasions in the following seventh-order linear tracking-error dynamics:

$$e_x^{(7)}+k_6^x{e}_x^{(6)}+k_5^x{e}_x^{(5)}+k_4^x{e}_x^{(4)}+k_3^x{e}_x^{(3)}+k_2^x{\ddot{e}}_x+k_1^x{\dot{e}}_x+k_0^x{e}_x=0$$

(53)

The design coefficients $\{k_6^x,k_5^x,\dots ,k_1^x,k_0^x\}$ are chosen so as to render the closed-loop characteristic polynomial

$$s^7+k_6^x{s}^6+k_5^x{s}^5+k_4^x{s}^4+k_3^x{s}^3+k_2^x{s}^2+k_1^{x}s+k_0^x=0$$

(54)

into a Hurwitz polynomial with desirable roots. Therefore, the specification of the set of design coefficients $\{k_6^x,k_5^x,\dots ,k_1^x,k_0^x\}$ is chosen so as to locate the desired closed-loop poles in the left half of the complex plane. The control parameters were selected so as to achieve the following desired closed-loop characteristic polynomial,

$$p_x^{des}(s)={\left(s^2+2{\zeta }_x{\omega }_{nx}s+{\omega }_{nx}^2\right)}^3\left(s+p_x\right)$$

(55)

where ζ_x, ω_nx and p_x are positive quantities. Therefore, the design coefficients $\{k_6^x,k_5^x,\dots ,k_1^x,k_0^x\}$ are given by:

$$\begin{array}{l}{k}_6^x=6{\zeta }_x{\omega }_{nx}+p_x\\ k_5^x=12{\zeta }_x^2{\omega }_{nx}^2+3{\omega }_{nx}^2+6{\zeta }_x{\omega }_{nx}{p}_x\\ k_4^x=12{\zeta }_x{\omega }_{nx}^3+8{\zeta }_x^3{\omega }_{nx}^3+12{\zeta }_x^2{\omega }_{nx}^2{p}_x+3{\omega }_{nx}^2{p}_x\\ k_3^x=12{\zeta }_x^2{\omega }_{nx}^2+2{\omega }_{nx}^4+12{\zeta }_x{\omega }_{nx}^3{p}_x+8{\zeta }_x^3{\omega }_{nx}^3{p}_x\\ k_2^x=6{\zeta }_x{\omega }_{nx}^5+12{\zeta }_x^2{\omega }_{nx}^2{p}_x+2{\omega }_{nx}^4{p}_x\\ k_1^x={\omega }_{nx}^6+6{\zeta }_x{\omega }_{nx}^5{p}_x\\ k_0^x={\omega }_{nx}^6{p}_x\end{array}$$

(56)

With a view to avoiding repetition, similar control laws are developed for the position variables y and z (given in expressions (36) and (37)). Substituting the pair $〈x,{ }_{\ddot{u}}〉$ by $〈y,{ }_{\theta }〉$ and $〈z,{ }_{\varphi }〉$ , respectively, the following feedback control laws are obtained for the variables y and z:

$$e_{{\Upsilon }_{\theta }}(s)=-\frac{k_3^y{s}^3+k_2^y{s}^2+k_1^{y}s+k_0^y}{s^3+k_6^y{s}^2+k_5^{y}s+k_4^y}{e}_y(s)$$

(57)

$$e_{{\Upsilon }_{\varphi }}(s)=-\frac{k_3^z{s}^3+k_2^z{s}^2+k_1^{z}s+k_0^z}{s^3+k_6^z{s}^2+k_5^{z}s+k_4^z}{e}_z(s)$$

(58)

On the other hand, the dynamics given in (38) comprise a second-order system in order to control the orientation angle ψ of the quadrotor towards a given smooth reference trajectory, ψ*(t), with τ_ψ acting as the control input. In this case, the open-loop control input ${\tau }_{\psi }^{*}(t)$ that ideally achieves the tracking of ψ*(t) for suitable initial conditions satisfies the following second-order dynamics:

$${\ddot{\psi }}^{*}={\tau }_{\psi }^{*}$$

(59)

Subtracting (58) from (38 yields)

$${\ddot{e}}_{\psi }=e_{{\tau }_{\psi }}$$

(60)

where e_ψ = ψ − ψ* and $e_{{\tau }_{\psi }}={\tau }_{\psi }-{\tau }_{\psi }^{*}(t)$ . A PD feedback controller, specifying the input tracking error, is given by

$$e_{{\tau }_{\psi }}=-k_2^{\psi }{\dot{e}}_{\psi }-k_1^{\psi }{e}_{\psi }$$

(61)

This yield, evidently, a differential equation for the output tracking error, e_ψ, given by

$${\ddot{e}}_{\psi }+k_2^{\psi }{\dot{e}}_{\psi }+k_1^{\psi }{e}_{\psi }=0$$

(62)

The characteristic polynomial, associated with this equation is

$$s^2+k_2^{\psi }s+k_1^{\psi }=0$$

(63)

Thus, the design problem reduces to an appropriate choice of the feedback controller gain so as to make the above polynomial like Hurwitz's. However, the tracking controller (60) requires knowledge of the signal ${\dot{e}}_{\psi }$ . We resort to an integral reconstructor for such a signal, aware of the fact that such a reconstructor differs from the actual signal by an unknown constant quantity fixed by the unchangeable initial condition. We proceed by integrating (38) once and, later on, by disregarding the constant error due to the tracking-error velocity's initial condition ${\dot{e}}_{\psi }(0)$ . The estimated error velocity ${\widehat{\dot{e}}}_{\psi }$ can be computed in the following form:

$${\widehat{\dot{e}}}_{\psi }={\dot{e}}_{\psi }-{\dot{e}}_{\psi }(0)={\displaystyle {\int }_0^t}{e}_{{\tau }_{\psi }}(\tau )d\tau$$

(64)

When the reconstructor is used in the derivative part of the PD controller, the constant error is suitably compensated by the addition of an integral error term to reject the effect of the unknown constant initial conditions ${\dot{e}}_{\psi }(0)$ . We thus proceed to propose the following controller:

$$e_{{\tau }_{\psi }}=-k_2^{\psi }{\widehat{\dot{e}}}_{\psi }-k_1^{\psi }{e}_{\psi }-k_0^{\psi }{\displaystyle {\int }_0^t}{e}_{\psi }(\tau )d\tau$$

(65)

Substituting (63) in (64), and after some rearrangement, we achieve

$$e_{{\tau }_{\psi }}=-k_2^{\psi }{\displaystyle {\int }_0^t}{e}_{{\tau }_{\psi }}(\tau )d\tau -k_1^{\psi }{e}_{\psi }-k_0^{\psi }{\displaystyle {\int }_0^t}{e}_{\psi }(\tau )d\tau$$

(66)

Following a similar procedure to that used for the variables x, y and z, and using Laplace's transform, the following result is yielded for the control input variable $e_{{\tau }_{\psi }}$ :

$$e_{{\tau }_{\psi }}(s)=-\frac{k_1^{\psi }s+k_0^{\psi }}{s+k_2^x}{e}_{\psi }(s)$$

(67)

Now, using (63) in the modified controller (64) results upon substitution of (59), and differentiating once with regard to time for the following third-order linear tracking-error dynamics,

$$e_{\psi }^{(3)}+k_2^{\psi }{\ddot{e}}_{\psi }+k_1^{\psi }{\dot{e}}_{\psi }+k_0^{\psi }{e}_{\psi }=0$$

(68)

the characteristic polynomial, associated with the previous equation, is easily shown to be

$$s^3+k_2^{\psi }{s}^2+k_1^{\psi }s+k_0^{\psi }=0$$

(69)

where the set of design coefficients $\{k_2^{\psi },k_1^{\psi },k_0^{\psi }\}$ are chosen so as to make the above polynomial like Hurwitz's with desirable roots. In order to specify the parameters $\{k_2^{\psi },k_1^{\psi },k_0^{\psi }\}$ , we can choose to locate the desired closed-loop poles in the left half of the complex plane. In particular, they were selected so as to achieve the following desired closed-loop characteristic polynomial,

$$p_{\psi }^{des}(s)=\left(s^2+2{\zeta }_{\psi }{\omega }_{n\psi }s+{\omega }_{n\psi }^2\right)\left(s+p_{\psi }\right)$$

(70)

where ζ_ψ, ω_nψ and p_ψ are positive quantities. Identifying each term of the expression (68) with those of (69), we obtain directly the value of the set of coefficients $\{k_2^{\psi },k_1^{\psi },k_0^{\psi }\}$ , which are given by

$$k_2^{\psi }=2{\zeta }_{\psi }{\omega }_{n\psi }+p_{\psi }$$

(71)

$$k_1^{\psi }=2{\zeta }_{\psi }{\omega }_{n\psi }{p}_{\psi }+{\omega }_{n\psi }^2$$

(72)

$$k_0^{\psi }={\omega }_{n\psi }^2{p}_{\psi }$$

(73)

Next, using (34), the following result is obtained:

$$\left[\begin{array}{c}\ddot{u}\\ {\tau }_{\theta }\\ {\tau }_{\varphi }\end{array}\right]=N\cdot \left[\begin{array}{c}{\Upsilon }_{\ddot{u}}\\ {\Upsilon }_{\theta }\\ {\Upsilon }_{\varphi }\end{array}\right]+B$$

(74)

Finally, the section is concluded by stating the following proposition, as proved in the above developments:

Proposition 2: Given a set of smooth reference trajectories $\left(x^{*},y^{*},z^{*},{\psi }^{*}\right)$ , for the horizontal coordinates x and y, the vertical position z and the orientation variable ψ in the quadrotor dynamics, defined by expressions (1)–(6), then the feedback controller defined by expressions (51), (56), (57), (66) and (71) produces the closed-loop behaviour of the tracking errors, $e_x=x-x^{*}(t)$ , $e_y=y-y^{*}(t)$ , $e_z=z-z^{*}(t)$ and $e_{\psi }=\psi -{\psi }^{*}(t)$ , which is locally governed by the linear dynamics

$$\begin{array}{l}{e}_i^{(7)}+k_6^i{e}_i^{(6)}+k_5^i{e}_i^{(5)}+k_4^i{e}_i^{(4)}+k_3^i{e}_i^{(3)}+k_2^i{\ddot{e}}_i+k_1^i{\dot{e}}_i+k_0^i{e}_i=0\\ e_{\psi }^{(3)}+k_2^{\psi }{\ddot{e}}_{\psi }+k_1^{\psi }{\dot{e}}_{\psi }+k_0^{\psi }{e}_{\psi }=0\end{array}$$

(75)

where the sub-index i = x,y,z, the design coefficients of which can be chosen according to expressions (55) and (70) so as to render the origin of the tracking-error space into an exponential asymptotically equilibrium point.

4. Numerical Simulations

Numerical simulations were carried out in order to verify the efficiency of the proposed approach in terms of the quick convergence of the tracking errors to a small neighbourhood of zero, smooth transient responses and low control effort. In the simulations, in order to evaluate the performance of the proposed controller, it is desirable to track the following sinusoidal trajectories for the variables x and y (defined as x*(t) and y*(t) respectively),

$$x^{*}(t)=R\left(\sin at+\sin bt\right)$$

(76)

$$y^{*}(t)=R\left(\cos at-\cos bt\right)$$

(77)

whereby R = 4[m], a = 0.0625 [rad/s7] and b = 0.1875 [rad/s]. The trajectory defined for the flat output z, defined as z*(t), is a smooth trajectory defined during a finite interval of the form [t_i,t_f], from an initial value $z^{*}(t_i)={\overline{z}}_i$ , to a final desired value $z^{*}(t_f)={\overline{z}}_f$ . We set, for instance,

$$z^{\ast }(t)={\overline{z}}_i+\left({\overline{z}}_f-{\overline{z}}_i\right)\phi (t,t_f,t_i)$$

(78)

with ${\overline{z}}_i=0[m]$ and ${\overline{z}}_f=5[m]$ , $t_i=0[s]$ , and where t_f = 25[s] and ϕ(t,t_f,t_i) are defined as a Bezier polynomial smoothly interpolating between 0 and 1 in the time interval [t_i,t_f]. We choose a 16th-order Bezier polynomial:

$$\begin{array}{l}\phi (t,t_i,t_f)=\\ (\begin{array}{c}\text{0 for }{\text{t<t}}_{\mathrm{i}}\\ {}_{{\left[\frac{\mathrm{t}-{\mathrm{t}}_{\mathrm{i}}}{{\mathrm{t}}_{\mathrm{f}}-{\mathrm{t}}_{\mathrm{i}}}\right]}^8\left[{\mathrm{r}}_1-{\mathrm{r}}_2\left(\frac{\mathrm{t}-{\mathrm{t}}_{\mathrm{i}}}{{\mathrm{t}}_{\mathrm{f}}-{\mathrm{t}}_{\mathrm{i}}}\right)+{\mathrm{r}}_3{\left(\frac{\mathrm{t}-{\mathrm{t}}_{\mathrm{i}}}{{\mathrm{t}}_{\mathrm{f}}-{\mathrm{t}}_{\mathrm{i}}}\right)}^2-\dots +{\mathrm{r}}_9{\left(\frac{\mathrm{t}-{\mathrm{t}}_{\mathrm{i}}}{{\mathrm{t}}_{\mathrm{f}}-{\mathrm{t}}_{\mathrm{i}}}\right)}^8\right]}\\ \text{for }{\mathrm{t}}_{\mathrm{i}}\le \mathrm{t}\le {\mathrm{t}}_{\mathrm{f}}\\ \text{1 for }{\text{t>t}}_{\mathrm{f}}\end{array}\end{array}$$

(79)

where the coefficients r₁,…, r₉ were obtained with polynomial interpolation satisfying the following restrictions:

$$\begin{array}{l}{z}_i=z(t_i);z_f=z(t_f)\\ \dot{z}(t_i)=\ddot{z}(t_i)=z^{(3)}(t_i)=\dots =z^{(7)}(t_i)=0\\ \dot{z}(t_f)=\ddot{z}(t_f)=z^{(3)}(t_f)=\dots =z^{(7)}(t_f)=0\end{array}$$

(80)

We find that:

$$\begin{array}{ccc}{r}_1=12870& r_2=91520& r_3=288288\\ r_4=524160& r_5=600600& r_6=443520\\ r_7=205920& r_8=54912& r_9=6435\end{array}$$

(81)

Finally, the trajectory for the flat output variable ψ, defined as ψ*, is designed as a linear evolution as

$${\psi }^{*}=\frac{\pi }{2}+ct$$

(82)

with c = 0.8[rad / s]

On the other hand, the sets of coefficients of the GPI controller $\{k_6^x,k_5^x,k_4^x,\dots ,k_0^x\}$ , $\{k_6^y,k_5^y,k_4^y,\dots ,k_0^y\}$ , $\{k_6^z,k_5^z,k_4^z,\dots ,k_0^z\}$ and $\{k_2^{\psi },k_1^{\psi },k_0^{\psi }\}$ were designed with the help of the following dominating Hurwitz characteristic polynomials ${(s^2+2{\zeta }_x{\omega }_{nx}s+{\omega }_{nx}^2)}^3(s+p_x)$ , ${(s^2+2{\zeta }_y{\omega }_{ny}s+{\omega }_{ny}^2)}^3(s+p_y)$ , ${(s^2+2{\zeta }_z{\omega }_{nz}s+{\omega }_{nz}^2)}^3(s+p_z)$ and $(s^2+2{\zeta }_{\psi }{\omega }_{n\psi }s+{\omega }_{n\psi }^2)(s+p_{\psi })$ with ${\zeta }_x=1$ , ${\omega }_{nx}=2$ , $p_x=2$ , ${\zeta }_y=1$ , ${\omega }_{ny}=2$ , $p_y=2$ , ${\zeta }_z=1$ , ${\omega }_{nz}=2$ , $p_z=2$ , ${\zeta }_{\psi }=1$ , ${\omega }_{n\psi }=1$ and $p_{\psi }=1$ . The time-sampling used in the simulations is set as 1 · 10⁻³[s].

Two sets of simulations were developed to establish a comparison between the GPI control presented in this paper and a classical PID control. The comparison is carried out on the basis of the following aspects: (1) the stabilization process and trajectory tracking; (2) performance when the signals are noisy; and (3) responses to environmental uncertainties, such as gusty winds. These simulations will be described in detail in what follows.

4.1 Simulation under Ideal Conditions

In this computer simulation, the quadrotor has to track the trajectory defined by expressions (73), (75) and (78) under ideal conditions, which implies that the measured signals are not corrupted by noise and that there are no environmental uncertainties.

The time evolution of the closed-loop centre of mass position-variables using the GPI and the PID control is illustrated in Figure 3, and the controlled evolution of the centre of mass of the quadrotor in 3D is depicted in Figure 4. As can be observed, the tracking of the variables x, y and z of the prescribed trajectory illustrates that the GPI control drives the system towards perfect tracking of the prescribed trajectories, showing an important improvement with respect to the PID control. This fact is demonstrated in full in the tracking trajectory of the quadrotor in 3D. Additionally, the tracking for the ψ variable presents better behaviour using the GPI controller, as is observed in Figure 5. Furthermore, in Figure 5 is shown the evolution of all the closed-loop attitude variables of the quadrotor using the GPI and PID controllers. Finally, the evolution of the control input variables is displayed in Figure 6, using the GPI and the PID control, illustrating the efforts made by the feedback controllers in guiding the errors of the states towards a fairly small neighbourhood close to zero.

Figure 3.

Position variables and reference variable of the center of mass of the quadrotor

Figure 4.

3D centre of mass quadrotor trajectory

4.2 Robustness with Respect to Noisy Signals and Environmental Perturbations

In this computer simulation, the measured controlled variables $\tilde{x}$ , $\tilde{y}$ , $\tilde{z}$ and $\tilde{\psi }$ are affected by an additive, zero mean, high-frequency noise μ_n(t), for n = x,y,z,ψ, such that,

Figure 5.

Attitude variable of the quadrotor

Figure 6.

Applied control inputs

$$\tilde{x}(t)=x(t)+{\mu }_x(t)$$

(83)

$$\tilde{y}(t)=y(t)+{\mu }_y(t)$$

(84)

$$\tilde{z}(t)=z(t)+{\mu }_z(t)$$

(85)

$$\tilde{\psi }(t)=\psi (t)+{\mu }_{\psi }(t)$$

(86)

where the selected standard deviations μ_n(t), for n = x,y,z, were set to be 3.17 · 10⁻⁴, and the standard deviation for μ_ψ(t) was set to be 3.17 · 10⁻⁴. Additionally, we introduced in the simulation the atmospheric disturbances (gusty wind) shown in Figure 7 and defined as in [26]. Figure 8 depicts the controlled evolution of the position variables of the centre of mass of the quadrotor under noisy measurements and atmospheric disturbances, and Figure 9 displays the evolution of the centre of mass of the quadrotor in 3D under these new conditions. Similarly as with the previous simulations, the performance of the quadrotor using the GPI control is improved significantly in comparison to that obtained with the PID control. It is observed that, when the atmospheric disturbances affect the quadrotor, the GPI controller corrects these undesirable effects and again drives the tracking error trajectories to a small neighbourhood on the origin of the tracking-error phase space. In Figure 10 is presented the evolution of all the closed-loop attitude variables of the quadrotor using the GPI and PID controllers as well as the tracking for the variable ψ, showing again improved behaviour with respect to the PID controller. Finally, in Figure 11 is displayed the evolution of the control inputs using the GPI and the PID controllers. In this figure the high robustness of the GPI control is demonstrated in comparison with the PID control when the measured controlled signal are affected by noise.

Figure 7.

Atmospheric disturbances used in the simulations

Figure 8.

Position variable and reference variable of the center of mass of the quadrotor under noisy measurements and with gusts of wind

Figure 9.

3D center of mass quadrotor trajectory under noisy measurements and with gusts of wind

Figure 10.

Attitude variable of the quadrotor under noisy measurements and with gusts of wind

Figure 11.

Applied control inputs noisy measurements and with gusts of wind

5. Conclusions

In this paper, the theoretical applicability of the GPI controller technique for regulation and trajectory tracking problems in a quadrotor has been investigated. The proposed scheme renders state observers and time discretizations completely unnecessary. GPI control sidesteps the need for traditional asymptotic state observers and proceeds directly to use structural state estimates in place of the actual state variables. The effect of such structural estimates in the controller are neglected in the feedback control law by means of suitable integral output tracking-error feedback control actions. Numerical simulations were provided to demonstrate the effectiveness of the proposed approach in comparison with the classical PID control in the following terms: (a) stabilization and trajectory tracking tasks; (b) performance when the measured signals are corrupted by noise; and (c) dynamic response when atmospheric disturbances, such as gusty wind, affect the quadrotor. The results show that the behaviour of the proposed approach improves the behaviour of the system in all aspects in comparison to the PID controller. Future work will be devoted to verifying the effectiveness of the proposed control algorithm through use in experiments with a real platform. This will be the topic of future publications.