<figure><img src="image_1.png" /></figure>

HAL
open science

# Path tracking of a small autonomous airplane in wind gusts

Elie Kahale, Yasmina Bestaoui, Pedro Castillo

► To cite this version:

Elie Kahale, Yasmina Bestaoui, Pedro Castillo. Path tracking of a small autonomous airplane in wind gusts. 25th IEEE/RSJ International Conference on Robotics and Intelligent Systems (IROS 2012), Oct 2012, Vilamoura, Algarve, Portugal. pp.402-407, 10.1109/IROS.2012.6385916 . hal-00868716

HAL Id: hal-00868716

https://hal.science/hal-00868716

Submitted on 7 Sep 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Path tracking of a Small Autonomous Airplane in Wind Gusts

E. Kahale, Y. Bestaoui, P. Castillo

**Abstract**— In this paper, we study the trajectory tracking problem of autonomous small airplanes in presence of unknown wind gusts. The vehicle is represented by its center of gravity, and the mathematical translational model is derived using Newton's second law. The controller is designed using Robust Control Lyapunov Function and Sontag's feedback stabilizing universal. The proposed control laws guarantees, by structure, the robustness of the system with respect to different uncertainties due to model and unknown external parameters. Numerical simulations are performed to demonstrate the performance of proposed control strategy.

## I. INTRODUCTION

A control approach widely used in practice is based on linearization of the nonlinear model around an operating point depending on the current flight regime [1], whereas the linearization may result in the design of sufficiently accurate controllers in the case of stabilization around the operating regime [2]. In the case of tracking a desired trajectory, the problem becomes more difficult because the linearized model is time varying. Using a fixed linear controller in such a case can result in an unacceptable response. Linear Parameter varying control has been presented as a reliable alternative to classical gain scheduling for multi-variable systems [3].

Gain scheduling is a standard method to design controllers for aircrafts over a wide performance envelope. It yields a global controller based on interpolation of a family of locally linearized controllers. The obtained controller comes with no guarantee on its stability or performance other than at the operating points. Linear parameter varying synthesis techniques naturally fit into the gain scheduling framework. In [4], [5], a robust control method for an airplane is introduced. The control strategy consists of an inner $H_\infty$ controller for the dynamics and an outer Single Input Single Output Proportional or Proportional-Integral controller for the remaining states. The classical control algorithms (like PID) have the advantage to be easily implemented and to provide reliable control performance. On the other hand, advanced control methods are mainly developed to improve its control performance in a complex and unstable flight environment.

Nonlinear control methodologies for autonomous aircrafts have been developed to increase performance and development times by dealing with the complete dynamics of the airplane. The most commonly employed control techniques are

feedback linearization, nonlinear predictive control, optimal control, sliding mode, direct and indirect adaptive approaches as well as intelligent approaches such as fuzzy logic, neural networks [6], [7]. Adaptive backstepping type approaches in which the linearized aircraft model parameters are estimated on line have also been explored for flight control. These methods allow a stability analysis but do not provide an easy means of designing the transient response. Model predictive control also referred to as receding horizon control is based on the principle of re-optimizing the controller in real time over a finite horizon [8], [9].

The main motivation for the present work is to produce a control design that is also more robust to model error. The restrictions concerning small autonomous airplanes and their mathematical models have to be considered. A great deal of uncertainty is caused by the constraints of the sensors and their measurements. Low cost components are often used and the measurements obtained are of lower quality when compared with conventional aircraft avionics system. These drawbacks have major influences on the quality of the dynamic model as well as on the estimation of the system states output measurements.

Even though there are a lot of modern control techniques which can be adopted for airplane control, robust controllers remain interesting for real applications of these vehicles when they are exposed to unknown environment (like wind gusts). In this work, we focus on introducing a robust control strategy to make autonomous path tracking of an autonomous airplane in the presence of wind gust. The control algorithm is based on Robust Control Lyapunov Function and Sontag's universal stabilizing feedback.

In order to design the control laws, the airplane equations of motion taking into account the wind disturbances are introduced in section II. For the purpose of guidance and navigation, the equations of translational dynamics are resolved in a reference frame fixed to the Earth, with origin at the center of mass. The control method and the control laws are presented and developed in section III. Some simulations are carried out to evaluate the methodology and the main results (depicted in graphs) are shown in section IV. Finally, conclusions and future work are discussed in section V.

## II. AIRPLANE EQUATIONS OF MOTION

From the point of view of guidance and control systems, an airplane can be modeled using a point mass model [10], [11], [12]. A geometric and mass symmetry about the body axis is assumed. The equations of motion are the result of Newton's laws of motion assuming that the rotational velocity of the

E. Kahale and Y. Bestaoui are with the IBISC EA 4526, University of Evry, France (e-mail: (kahale,bestaoui)@iup.univ-evry.fr)

P. Castillo is with Heudiasyc laboratory, UMR 7253, Université de Technologie de Compiegne, BP 20529-60205 Compiegne cedex, France (e-mail: castillo@hds.utc.fr).
earth is neglected (also so-called flat earth assumption). The symmetric flight case is considered, hence ignoring the sideslip angle : the angle between the forward velocity and the nose of the airplane.

The vehicle's motion is defined by the following coordinate systems :

- **Earth Fixed Inertial frame** $R_I(O, x, y, h)$ related to the geographic coordinate system **NED** which has the unit vectors $\vec{n}$ on $x$ axis, $\vec{e}$ on $y$ axis and $\vec{d}$ on $h$ axis pointing respectively to the North, East and Down. The origin of this frame $O$, chosen arbitrarily, may be for example the initial point of the vehicle.

- **Local Horizon frame** $R_h(CG, x_h, y_h, h_h)$. This frame moves with the airplane, but its axes remain parallel to the earth fixed inertial frame axes and its origin CG is the center of gravity of the vehicle.

- **Body Fixed frame** $R_b(CG, x_b, y_b, h_b)$ related to the geometry of the airplane : $x_b$ axis is directed along the axis of symmetry, $y_b$ axis to the right and $h_b$ axis to the down. The origin of this frame CG is the center of gravity of the vehicle, considered as an invariant point in the vehicle.

- **Wind Relative frame** $R_w(CG, x_w, y_w, h_w)$. The origin of this frame is also the center of gravity of the vehicle (CG). The $x_w$ axis is aligned with the relative flight speed vector $\vec{V}$. $h_w$ axis belongs to the symmetry plan ($x_b - h_b$). The angle of attack ($\alpha$) and the sideslip angle ($\beta$) are used to defined the position of relative speed vector with respect to $R_b$. $\alpha$ is the angle that the straight line drawn from the leading edge to the trailing edge of the airfoil makes with the freestream velocity vector.

The translational velocity vector is given by :

$$ \vec{r} = \vec{V}_I \qquad (1) $$

and from Newton's law

$$ m \vec{V}_I = \sum \vec{F} \qquad (2) $$

where $\vec{r}$ defines the position vector, $m$ represents the mass of the airplane, $\vec{F}$ describes the external forces acting on the vehicle, and $\vec{V}_I$ denotes the inertial velocity vector which has the form

$$ \vec{V}_I = \vec{V} + \vec{V}_W \qquad (3) $$

where $\vec{V}$ represents the relative velocity vector and $\vec{V}_W$ describes the wind velocity vector. In addition, $\vec{V} \in R_h$ is defined by its magnitude $V$, the heading angle $\chi$ (measured from the north to the projection of $\vec{V}$ in the local horizon plane), and the flight path angle $\gamma$ (vertically up to $\vec{V}$). $\vec{V}_W \in R_I$ is composed by $[W_x \quad W_y \quad -W_z]^T$. Therefore, we can write

$$ \vec{V}_I = (V \cos \chi \cos \gamma + W_x) \vec{n} + (V \sin \chi \cos \gamma + W_y) \vec{e} + (-V \sin \gamma - W_z) \vec{d} \qquad (4) $$

Then, the three-dimensional translational kinematic equations of an aerial vehicle taking into account the wind effect

can be expressed as following

$$ \begin{aligned} \dot{x} &= V \cos \chi \cos \gamma + W_x \\ \dot{y} &= V \sin \chi \cos \gamma + W_y \\ \dot{z} &= -\dot{h} = V \sin \gamma + W_z \end{aligned} \qquad (5) $$

where $x$ and $y$ define the horizontal position of the vehicle and $z$ denotes its altitude. The dynamic of the system is introduced using Newton-Euler equation (2). This equation is described in the wind relative frame ($R_w$). A problem is that coordinate system $R_w$ is a rotating coordinate system in which Newton's Laws do not apply. In order to face this problem we have to apply the following transformation relationship

$$ \left. \frac{d\vec{A}}{dt} \right|_{\text{fixed}} = \left. \frac{\partial \vec{A}}{\partial t} \right|_{\text{rotating}} + \vec{w} \times \vec{A} \qquad (6) $$

where $\vec{w}$ is the angular rotation of $R_w$ relative to $R_h$. This vector is also referred to as the angular velocity of the airplane relative to the earth. Then, the right hand of equation (2) can be expressed by

$$ \frac{dV_I}{dt} = \frac{dV}{dt} + \frac{dV_W}{dt} \qquad (7) $$

where

$$ \begin{aligned} \frac{dV}{dt} &= \begin{bmatrix} \dot{V} \\ 0 \\ 0 \end{bmatrix} + \begin{bmatrix} \dot{\chi} \sin \gamma \\ \dot{\gamma} \\ \chi \cos \gamma \end{bmatrix} \times \begin{bmatrix} V \\ 0 \\ 0 \end{bmatrix} \\ &= \begin{bmatrix} \dot{V} \\ \dot{\chi}V \cos \gamma \\ -\dot{\gamma}V \end{bmatrix} \end{aligned} \qquad (8) $$

and

$$ \left( \frac{dV_W}{dt} \right)_W = C_h^w \left( \frac{dV_W}{dt} \right)_I = \begin{pmatrix} W_x \cos \chi \cos \gamma + \dot{W}_x \sin \chi \cos \gamma + \dot{W}_z \sin \gamma \\ -\dot{W}_x \sin \chi + \dot{W}_y \cos \chi \\ W_x \cos \chi \sin \gamma + W_y \sin \chi \sin \gamma - W_z \cos \gamma \end{pmatrix} \qquad (9) $$

with $C_h^w$ presents the transformation matrix from the local horizon frame $R_h$ to the wind relative frame $R_w$, and is defined as

$$ C_h^w = \begin{bmatrix} \cos \chi \cos \gamma & \sin \chi \cos \gamma & -\sin \gamma \\ -\sin \chi & \cos \chi & 0 \\ \cos \chi \sin \gamma & \sin \chi \sin \gamma & \cos \gamma \end{bmatrix} \qquad (10) $$

The external forces acting on the aircraft are : the gravity force ($f_g$), the aerodynamical force ($f_a$) and the thrust force ($f_T$). These forces are described in wind relative frame $R_w$. $f_g$ is presented by

$$ f_g = C_h^w \begin{pmatrix} 0 \\ 0 \\ mg \end{pmatrix} = \begin{pmatrix} -mg \sin \gamma \\ 0 \\ mg \cos \gamma \end{pmatrix} \qquad (11) $$

whilst $f_a$ and $f_T$ are described as

$$ \begin{aligned} f_{aT} &= f_a + f_T = C_{w'}^w \begin{bmatrix} T \\ 0 \\ 0 \end{bmatrix} + \begin{bmatrix} -D \\ 0 \\ -L \end{bmatrix} \\ &= \begin{bmatrix} T \cos \alpha - D \\ (T \sin \alpha + L) \sin \sigma \\ -(T \sin \alpha + L) \cos \sigma \end{bmatrix} \end{aligned} \qquad (12) $$
Notice that, the aerodynamic forces are produced by the relative aircraft motion with respect to the air flow. In the above equation we have that, $C_b^w$ represents the transformation matrix from $R_b$ to $R_w$ and $C_{w'}^w$ defines the transformation matrix from the frame $R_{w'}$, which is the result of the rotation of lift vector $\vec{L}$ about the velocity vector $\vec{V}$ through an angle $\sigma$ known as bank angle, to wind relative frame $R_w$. These matrices are stated as

$$C_b^w = \begin{bmatrix} \cos \alpha & 0 & -\sin \alpha \\ 0 & 1 & 0 \\ \sin \alpha & 0 & \cos \alpha \end{bmatrix} \quad (13)$$

$$C_{w'}^w = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \sigma & -\sin \sigma \\ 0 & \sin \sigma & \cos \sigma \end{bmatrix} \quad (14)$$

from equations (7), (8), (9), (11) and (12) we obtain

$$\dot{\gamma} = \frac{(T \sin \alpha + L) \cos \sigma}{mV} - \frac{g \cos \gamma}{V} + \frac{\dot{W}_x \cos \chi \sin \gamma}{V} + \frac{\dot{W}_y \sin \chi \sin \gamma + \dot{W}_z \cos \gamma}{V} \quad (15a)$$

$$\dot{\chi} = \frac{(T \sin \alpha + L) \sin \sigma}{mV \cos \gamma} + \frac{\dot{W}_x \sin \chi - \dot{W}_y \cos \chi}{V \cos \gamma} \quad (15b)$$

$$\dot{V} = \frac{T \cos \alpha - D}{m} - g \sin \gamma - \dot{W}_x \cos \gamma \cos \chi - \dot{W}_y \sin \chi \cos \gamma + \dot{W}_z \sin \gamma \quad (15c)$$

The lift $L$ and drag $D$ forces are expressed by the following equations

$$L = \frac{1}{2} C_L(M, \alpha) V^2 A_r \rho \quad (16)$$

$$D = \frac{1}{2} C_D(M, \alpha) V^2 A_r \rho$$

with the parameter $A_r$ is the reference area or characteristic area, $C_L(M, \alpha)$ and $C_D(M, \alpha)$ denote respectively the lift and drag parameters depending on the Mach number ($M$) and the angle of attack ($\alpha$). These parameters are estimated via wind tunnel experiments, more details, see [10].

### III. CONTROL STRATEGY USING RCLF

The control goal is to realize path tracking control of an airplane in presence of wind gusts. The main idea is to use the RCLF (Robust Control Lyapunov Functions) properties to derive the robust control strategy.

#### A. Inverse Optimal Robust Stabilization Problem

The task of constructing control laws using RCLF for uncertain nonlinear control systems is introduced by Freemann and Kokotovic in [13]. The approach can be summarized in the following paragraphs.

Assuming the following affine system with disturbances

$$\dot{x} = f(x, u, w) = f_0(x) + f_1(x)u + f_2(x)w \quad (17)$$

where $x \in \mathcal{X}$ represents the state variables, $u \in \mathcal{U}$ describes the control inputs, $w \in \mathcal{W}$ introduces the disturbances in a convex space $\mathcal{W}$ and $f_0(x)$, $f_1(x)$, $f_2(x)$ defines continuous functions.

Given a Robust Control Lyapunov Function $L$ for the system (17), we define $D: \mathcal{X} \times \mathcal{U} \to \mathbb{R}$ and $K: \mathcal{X} \rightsquigarrow \mathcal{U}$ by

$$D(x, u) := \max_{w \in \mathcal{W}(x)} [L_f L(x, u, w) + \varphi_v(x)]$$

$$K(x) := \{u \in \mathcal{U}(x) : D(x, u) < 0\}$$

Where, $L_f L(x, u, w)$ denote the Lyapunov derivative which is given by the equation:

$$L_f L(x, u, w) := \frac{\partial L(x)}{\partial x} \cdot f(x, u, w)$$

and $\varphi_v(x)$ is a class $K_\infty$ function. Then, the above implies that

$$\begin{aligned} D(x, u) &= \psi_0(x) + \psi_1^T(x)u \\ K(x) &= \{u \in \mathcal{U} : \psi_0(x) + \psi_1^T(x)u < 0\} \end{aligned}$$

where

$$\psi_0(x) = \nabla L(x) \cdot f_0(x) + \|\nabla L(x) \cdot f_2(x)\| + \varphi_v(x)$$

$$\psi_1(x) = [\nabla L(x) \cdot f_1(x)]^T$$

and $\varphi_v(x) > 0$.

Therefore, the control law that stabilize the system (17) is

$$u(x) := \underset{\|u\|: u \in K(x)}{\operatorname{argmin}} \{\|u\| : u \in K(x)\} \quad (18)$$

or

$$u = \begin{cases} \frac{-\psi_0(x)\psi_1(x)}{\psi_1^T(x)\psi_1(x)} & ; \quad \psi_0(x) > 0 \\ 0 & ; \quad \psi_0(x) \le 0 \end{cases} \quad (19)$$

Observe that (19) depends on $\varphi_v$ through the $\psi_0$ function. Note also that, there is never division by zero because the set $K(x)$ is nonempty [15]. Remark that the function $\varphi_v$ represents the desired negativity of the Lyapunov derivative, and it can be adjusted to achieve a tradeoff between the control effort and the rate of convergence of the state to zero.

#### B. Airplane Control

Aircrafts are characterized by nonlinear models that are non affine in the control input [6]. One nonlinear approach to this problem is that based on directly inverting the nonlinear function of the control on a domain. Although, the existence of an inverse function can be guaranteed by the implicit function theorem, it is generally difficult to prescribe a technique to actually obtain such an inverse.

The nonlinear state equation non affine in the control is transformed into an augmented state space model in which the new control appears in a linear fashion [7], [14]. For this purpose, the control variables are regarded as a state variables while its derivatives are used as virtual control inputs. So, the new state variables are $[x \ U]$. Whilst, the virtual control inputs are $\dot{U}$.

Consider the nonlinear system affine in control given by

$$\dot{X} = f_0 + f_1 U + f_2 W \quad (20)$$

where

$$X_1 = (x \ y \ z)^T, \quad X_2 = (\gamma \ \chi \ V)^T$$

$$X_3 = (T \ \sigma \ \alpha)^T$$
and $U = \dot{X}_3$, $W = [W_1 \ W_2]^T$; $W_1 = [W_x \ W_y \ W_z]$, $W_2 = \dot{W}_1$. with

$$f_0(X) = \begin{bmatrix} B(X_2) \\ \Phi(X_2, X_3) \\ 0_{3 \times 1} \end{bmatrix}$$

$$f_1(X) = \begin{bmatrix} 0_{3 \times 3} \\ 0_{3 \times 3} \\ I_{3 \times 3} \end{bmatrix}$$

$$f_2(X) = \begin{bmatrix} I_{3 \times 3} & 0_{3 \times 3} \\ 0_{3 \times 3} & B_w \\ 0_{3 \times 3} & 0_{3 \times 3} \end{bmatrix}$$

Likewise

$$\Phi(X_2, X_3) = \begin{pmatrix} \frac{L \cos \sigma - m g \cos \gamma}{m V} + \frac{T \sin \alpha \cos \sigma}{m V \sin \sigma} \\ \frac{T \sin \alpha \cos \sigma}{m V \sin \sigma} \\ \frac{-D - m g \sin \gamma}{m} + \frac{T \cos \alpha}{m} \end{pmatrix}$$

Indeed,

$$B(X_2) = \begin{pmatrix} V \cos \gamma \cos \chi & V \cos \gamma \sin \chi & V \sin \gamma \end{pmatrix}^T$$

$$B_w = \begin{bmatrix} \frac{\cos \chi \sin \gamma}{V} & \frac{\sin \chi \sin \gamma}{V} & \frac{\cos \gamma}{V} \\ \frac{\sin \chi}{V \cos \gamma} & \frac{-\cos \chi}{V \cos \gamma} & 0 \\ -\cos \chi \cos \gamma & -\sin \chi \cos \gamma & \sin \gamma \end{bmatrix} \quad (21)$$

Let us consider the following Lyapunov function

$$S(X) = \frac{X^T P X}{2}$$

where $P_{9 \times 9}$ is a symmetric and positive definite matrix. Thus, the derivative of the Lyapunov function is given by

$$\begin{aligned} \nabla S(X) &= \frac{1}{2} \frac{\partial [X^T P X]}{\partial X} = P X \\ &= \begin{bmatrix} \nabla_{S_1} & \nabla_{S_2} & \nabla_{S_3} \end{bmatrix}^T \end{aligned}$$

where $\dim(\nabla_{S_i}) = 3$, $i = 1, 2, 3$

Then $\psi_0(X)$ can be rewritten as

$$\psi_0(X) = \nabla S(x) \cdot f_0(X) + \| \nabla S(X) \cdot f_2(X) \| + \varphi_v(X) \quad (22)$$

with

$$\begin{aligned} \nabla S(x) \cdot f_0(X) &= \nabla_{S_1} B_1 X_2 + \nabla_{S_2} \Phi_0(X_2) \\ &+ \nabla_{S_2} \Phi_1(X_3) \end{aligned} \quad (23)$$

$$\| \nabla S(X) \cdot f_2(X) \| = \| [\nabla_{S_1} \ \nabla_{S_2} B_w] \|$$

where $\varphi_v(X) = X^T M X$ and $M_{9 \times 9}$ defines a diagonal positive matrix. On the other hand, $\psi_1(X)$ is computed by

$$\begin{aligned} \psi_1(X) &= [\nabla S(X) \cdot f_1(X)]^T \\ &= \nabla_{S_3} \end{aligned}$$

Once these tasks completed, the robust control laws are given by:

For $\psi_0(X) > 0$

$$u_1 = \frac{-\psi_0(X) \cdot \psi_1(1)}{\psi_1(1)^2 + \psi_1(2)^2 + \psi_1(3)^2}$$

$$u_2 = \frac{-\psi_0(X) \cdot \psi_1(2)}{\psi_1(1)^2 + \psi_1(2)^2 + \psi_1(3)^2}$$

$$u_3 = \frac{-\psi_0(X) \cdot \psi_1(3)}{\psi_1(1)^2 + \psi_1(2)^2 + \psi_1(3)^2}$$

whilst for $\psi_0(X) \le 0$, $u_1 = u_2 = u_3 = 0$.

## IV. SIMULATION RESULTS

In order to validate the performance of the closed-loop system, we carried out some simulations, using first the kinematic model. The wind model used in simulations is a combination of white noise and one minus cosine formulation, see figure 1. The gust velocity is defined as,

$$V_{wind} = \frac{1}{2} V_m \left[ 1 - \cos \left( \frac{\pi x}{H} \right) \right] \quad (24)$$

where $V_m$ presents gust amplitude, $x$ is the traveled distance and $H$ is the distance between the start point of the gust and its maximum value. $P$ is obtained by linearizing the nominal system (without perturbations) and solving the algebraic Riccati equation.

<figure><img src="image_3.png" /><figcaption>Fig. 1. Wind gust applied to the airplane</figcaption></figure>

On the other hand, $M$ can be regarded as a weight matrix which specify the relative cost of each state variable. A large value of $M$ makes $\varphi_v(x)$ larger. That guarantee a large negativity of Lyapunov derivative, which affect the response characteristics. In this paper, $M$ is chosen in a such way to have a suitable results concerning the overshoot which has to respect the vehicle limitations: $|\gamma_{max}| = 45^\circ$ and $|V_{max}| = 50$ m/s. The selected value is $M = 10^8 I_{9 \times 9}$.

### A. Trajectory with a line form

Our objective in this case is to stabilize the airplane with a fixed orientation and velocity. The initial conditions for the airplane are $\gamma_0 = \chi_0 = 0^\circ$ and $V_0 = 10$ m/s, whilst the final conditions for the vehicle are $\gamma_f = 10^\circ$, $\chi_f = 90^\circ$ and $V_f = 20$ m/s. The initial coordinates for the desired trajectory are $x_0 = y_0 = z_0 = 0$ m while the final coordinates are $x_f = 5000$ m, $y_f = 13340$ m and $z_f = 4250$ m.

The performance of the system when applying the control strategies is illustrated in Figures 2 - 4. The Positional errors along $x$, $y$ and $z$ axes are shown in Figure 2. On the other hand, the responses of the airplane's orientation and its velocity are presented in Figure 3. In these figures, the solid line represents the system response and the dashed line describes the desired value or trajectory. In Figure 4, the control inputs obtained are depicted.
<figure><img src="image_4.png" /><figcaption>Fig. 2. Positional error along x, y and z axis (first case).</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 3. The responses of the airplane's orientation and its velocity, $γ$, $χ$ and $V$, in closed-loop (first case).</figcaption></figure>

## B. Trajectory with a curve form

In this case, the orientation and the velocity are changing. For simulation purposes, we assume that all $γ$, $χ$ and $V$

<figure><img src="image_6.png" /><figcaption>Fig. 4. Control inputs $u_1$, $u_2$ and $u_3$ (first case).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 5. Positional error along x, y and z axis (second case).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 6. Flight path angle responses in the closed-loop (second case)</figcaption></figure>

are a piecewise constant functions which allow us to reach $γ_f = 0°$, $χ_f = 90°$ and $V_f = 25$ m/s. The initial conditions for the vehicle are $γ_0 = χ_0 = 0°$ and $V_0 = 15$ m/s. The desired trajectory is proposed to relate the initial position : $x_0 = y_0 = z_0 = 0$ m with the final one : $x_f = 1000$ m, $y_f = 14000$ m and $z_f = 2600$ m.

Figures 5-9 illustrate the performance of the closed-loop system. The positional errors along x, y and z axes are shown in Figure 5. In Figures 6, 7 and 8, we introduce the responses of the airplane's orientation and its velocity. Finally, the control inputs responses are depicted in Figure 9.

From the previous graphs obtained in simulation (either case one or case two) it can be observed that the proposed control strategy has a good performance and the airplane remains stable even in the presence of unknown wind gust. An example of this performance can be clearly illustrated in Figures 5 and 6. In these figures, notice that in the time period $50 \rightarrow 150$ [sec] the wind velocity (see also Figure 1) pushes the airplane strongly and the displacement of the vehicle on z axis differs significantly from the desired one.
<figure><img src="image_9.png" /><figcaption>Fig. 7. Heading angle responses in the closed-loop (second case)</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8. Velocity responses in the closed-loop (second case)</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 9. Control inputs $u_1, u_2$ and $u_3$ (second case)</figcaption></figure>

In this case, the system reaction was to change the flight path angle which allowed the vehicle to rejoin the desired trajectory.

**Remark :** We mentioned previously that the matrix *P* is obtained by linearization around the operating point. So, in the trajectory tracking problem, the selected value for the matrix *P* varies as $\gamma, \chi$ and *V* vary.

## V. CONCLUSIONS

A nonlinear control strategy using a Robust Control Lyapunov Function to stabilize a small autonomous aircraft in the presence of unknown wind gusts was proposed in this paper. The nonlinear model of the airplane was written using a point mass model. Some numerical simulations were carried out and some graphs were presented to illustrate the good performance of the closed-loop system even in presence of unknown disturbances. The drawback is sometimes it can be difficult to select the design parameters to specify the transient response.

## REFERENCES

[1] Stevens B. L., and Lewis J. D., Aircraft Control and Simulation, Wiley Press, 2003.

[2] Tewari A., Automatic Control of Atmospheric and Space Flight, Birkhauser, 2011.

[3] Marcos A., and Balas G. J., Development of linear parameter varying model, In AIAA journal of guidance, control and dynamics, vol. 27, pp. 218-223, April 2004.

[4] Wimmer D. A., Bildstein M., Well K. H., Schlenker M., Kungl P., and Kroplin B. H., Development and operation controllers for autonomous flight phases, InWorkshop on Aerial Robotics, IEEE International Conference on Intelligent Robots and Systems, pages 55-68, Lausanne, Switzerland, October 2002.

[5] Colgren R. D., Applications of robust control to nonlinear systems, AIAA Press 2004.

[6] Boskovic J. V., Chen L., Mehra R. K., Adaptive control design for non affine models arising in flight control, In AIAA Journal of guidance, control and dynamics, vol. 27, pp. 209-217, 2004.

[7] Singh S. N., Steinberg M.L., Page A. B., Nonlinear adaptive and sliding mode flight path control of FA 18 model, IEEE Trans. On aerospace and electronic systems, vol. 39, pp. 1250-1262, 2003.

[8] Scholte E., Campbell M. E., Robust nonlinear model predictive control with partial state information, In IEEE Transactions on Control Systems technology, volume 6 pp : 636-651, July 2008.

[9] Bestaoui, Y., Kahale, E., Time optimal trajectories for an autonomous airship' IEEE Workshop on Robot Motion Control (ROMOCO 2011), Bukowy Dworek, Poland, June 2011.

[10] Vinh, Nguyen X., Flight Mechanics of High-Performance Aircraft, Cambridge University Press, 1993.

[11] Hull, David G., Fundamentals of Airplane Flight Mechanics, Springer, 2007.

[12] Imado F., Heike Y., and Kinoshita T., Research on a New Aircraft Point-Mass Model, Journal of Aircraft Vol.48, No. 4, July-August 2011

[13] Freemann R.A., and Kokotovic P.V., Robust Nonlinear Control Design : State-Space and Lyapunov Techniques, Birkhauser, 1996.

[14] Levine W.S. The Control Handbook, Second Edition : Control System, Advanced Methods, CRC Press, Taylor & Francis Group, 2011.

[15] Munoz L.E., Santos O., and Castillo P., Robust nonlinear real-time control strategy to stabilize a PVTOL aircraft in crosswind, IEEE/RSJ International Conference on Intelligent Robots and Systems, 2010.

[16] Sontag, E. D., A universal construction of Artein's theorem on nonlinear stabilization, System and control letters, vol. 13, no. 2, pp. 117-123, 1989.