Currently, Unmanned Aerial Vehicles (UAVs) are used for many tasks: search, rescue, environmental monitoring, surveillance, and inspection. With the help of UAVs, it is possible to easily enter environments that are difficult for people or manned vehicles to enter. This is especially true for extreme situations: firefighting, reconnaissance etc^[1]. UAVs are becoming more popular in the form of multi-rotor devices, more commonly a quadcopter, which is a platform with four rotors. Compared with helicopter-type UAVs with main and tail rotors, quadcopters offer many advantages—reliability and simplicity of design, greater stability, compactness and maneuverability, light take-off weight, and greater payload quality. A quadcopter has four engines with propellers that generate thrust. The axis of the propeller and the angle of the blades are fixed, and only the rotational speed is adjusted, which significantly simplifies the design. Vertical movement is carried out by a synchronized change in the rotational speed of all propellers; for horizontal movement, it is necessary to tilt the quadcopter, which is achieved by a corresponding change in the rotational speed of the different propellers, which produces the moment required for tilting. The opposite direction of rotation of the screw pair ensures the compensation of the drag torque. Currently, such devices are quite widespread and versatile, but this use is primarily limited to the “manual” remote control mode from the operator console.

Researchers have been working on the problem of trajectory tracking control in UAVs along known trajectories. Various methods are discussed and proposed. Kumar et al^[2] use a linear quadratic controller to successfully solve the problem of stability during hover, but not suitable for significant nonlinearity and cross-coupling. In literature [3], an algorithm for tracking the trajectory of a quadcopter based on a PD controller and backstepping with a correction filter is proposed. In literature [4], a neural network-based nonlinear controller is proposed, and in literature [5], a “real-time” algorithm for trajectory planning and collision avoidance is proposed. In literature [6], an adaptive position tracking algorithm for vertical take-off and landing under limited external disturbance is considered. The authors in literature [7] propose a control system that realizes the movement of the quadcopter along a straight line connecting the start points and end points of the trajectory at a given speed and a given altitude above the Earth's surface. In literature [8], an algorithm for moving a quadcopter along any trajectory in space is considered, where the quadcopter is able to follow a moving object with known coordinates and velocities.

In addition to the above methods, a commonly used method for quadrotor UAV trajectory tracking is the sliding mode control. Sliding Mode Control (SMC) is a powerful nonlinear strategy that ensures robustness despite parameter uncertainty in the model^[9-10]. A number of papers have been published applying SMC methods to the position and attitude tracking problems of quadcopters and ensuring the robustness to external disturbances. The author in literature [11] applies this method to the quadrocopter with the motivation to theoretically exploit the robustness of the method with regard to model uncertainties and disturbances. The problem of chattering, which occurs with SMC by switching along with the sliding mode, is prevented by a smooth approximation of the signum function in literature [12]. In literature [13], the authors compare the SMC with the backstepping control of a miniature quadrotor, and the simulation results show a robust and practical behavior of SMC. In literature [14], an integral sliding mode is induced by a Newton-Euler under-driven dynamic model of a quadrotor aircraft, which in turn leads to a robust SMC-based controller. A robust terminal sliding mode control algorithm for the full-motion subsystem of the controller is presented in literature [15]. SMC is also combined with other control algorithms^[16-17]. For example, in literature [18], the authors use adaptive fuzzy hierarchical sliding mode control to solve the trajectory tracking problem of a quadrotor UAV. In addition, the sliding mode technique is also designed as an estimator of external disturbances. A control method for implementing the entire observer-estimator using SMC is presented in literature [19]. For some other related works, see literatures [20-23].

However, in most of these aforementioned SMC-based algorithms, the authors do not consider the first and second-order derivatives of rotor dynamics, 3D trajectories, UAV positions, and angles, which makes these algorithms ineffective in solving the 3D trajectory tracking problem of quadrotor UAVs in noisy environments (note that the noisy environment refers to any external disturbances. In the later experiment part of the paper, it means the perturbations applied to the horizontal and vertical axes of the UAV). The novelty of this work is that we propose a robust sliding mode controller for the quadrotor UAV in the presence of perturbations and parameter uncertainties (specifically, with mass, inertia, and stiffness uncertainties), and verify our algorithm in both simulation (MATLAB Simulink and Universal Mechanism) and real experiment (Pelican quadrotor platform). The manuscript is structured as follows. In Section 2, the considered UAV is described and the flight dynamics model is derived. The design of the sliding mode controller is presented in Section 3. The results of the numerical simulations are shown in Section 4. And the real experimental validation of the proposed controller is presented in Section 5. The review and discussion in Section 6 conclude the paper.

Notation. For an uncertain model parameter $ p $ , the notations $ \overline{p} $ , $ \underline{p} $ and $ \hat p $ stand the maximum, minimum and nominal value of $ p $ , respectively. Moreover we assume that they satisfy $ \hat p = (\overline{p}+\underline{p})/{2} $ .

1 System Dynamics

When modeling the flight of a quadrocopter, the earth's motion can be ignored due to the very limited range and duration of flight, that is, the earth's coordinate system can be regarded as inertia. At the initial stages of development, one can also consider the non-load-bearing part of the apparatus as a solid body, and consider the wind only as an external disturbance.

Then, the $ 6 $ -DOF nonlinear dynamic model of UAV is given by (see Fig. 1):

$ \tag{1a} \ddot x = \frac{F}{m}(\cos\psi \sin\theta \cos\phi+\sin\psi \sin\phi)-\frac{k_x}{m}\dot x $

$ \tag{1b}\ddot y = \frac{F}{m}(\sin\psi \sin\theta \cos\phi-\sin\phi \cos\psi)-\frac{k_y}{m}\dot y $

$ \tag{1c}\ddot z = \frac{F}{m}(\cos\theta \cos\phi)-\frac{k_z}{m}\dot z-g $

$ \tag{1d}\ddot \theta = \frac{1}{I_y}(\tau_1+(I_z-I_x)\dot\phi \dot\psi) $

$\tag{1e} \ddot \phi = \frac{1}{I_x}(\tau_2+(I_y-I_z)\dot\theta \dot\psi) $

$ \tag{1f}\ddot \psi = \frac{1}{I_z}(\tau_3+(I_x-I_y)\dot\phi \dot\theta) $

where $ x $ , $ y $ , $ z $ are positions of the UAV, $ \theta $ is the pitch angle, $ \phi $ is the roll angle, and $ \psi $ is the yaw angle. Parameter $ m $ is the mass, $ I_i $ is the rotational inertia of the $ i $ -axis, $ k_i $ and ${\tau _i} $ are the coefficient of air resistance and the torque of the i-axis, respectively ( $ i\in\{x,y,z\} $ ). Constant $ g $ is the gravitational acceleration. If we neglect the inertia of the screws when changing the angular velocities of their rotation, then these moments (i.e., the control inputs) can be expressed as follows:

$ \begin{bmatrix} F\\ \tau_1\\ \tau_2\\ \tau_3 \end{bmatrix} = \begin{bmatrix} 1&1&1&1\\ -l &l &0& 0\\ 1& 0 &-l& l\\ -c_\tau& -c_\tau& c_\tau & c_\tau \end{bmatrix} {\begin{bmatrix} f_1\\ f_2\\ f_3\\ f_4 \end{bmatrix}} $

(2)

where $ l $ is the length between the center of the UAV and the rotors and $ c_\tau $ is the torque coefficient and $ f_{1}{\text{-}}f_{4} $ are the lifting forces of the rotors. Taking into account the symmetry of the spacecraft and assuming that the center of mass is located at the origin of the coupled system, we can define the following auxiliary inputs:

$\tag{3a} F_x = F(\cos\psi \sin\theta \cos\phi+\sin\psi \sin\phi) $

$ \tag{3b}F_y = F(\sin\psi \sin\theta \cos\phi-\sin\phi \cos\psi) $

$ \tag{3c}F_z = F(\cos\theta \cos\phi) $

It implies that the changes in the Euler angles are determined through the projections of the angular velocity by the kinematic Euler equations:

$ \tag{4a}F = \frac{F_z}{\cos\theta \cos\phi} $

$ \tag{4b}\theta = \arctan\left(\frac{F_x \cos\psi+F_y \sin\psi}{F_z}\right) $

$ \tag{4c}\phi = \arctan\left(\cos\theta\frac{F_x \sin\psi-F_y \sin\psi}{F_z}\right) $

Then we consider $ F_x $ , $ F_y $ and $ F_z $ as control inputs and they result in $ \theta $ and $ \phi $ . In fact, we do not need to feed desired $ \theta $ and desired $ \phi $ and they can be obtained from equation (4). Using equation (1) and equation (3) we have:

$ \tag{5a}\ddot x = \frac{F_x}{m}-\frac{k_x}{m}\dot x $

$ \tag{5b}\ddot y = \frac{F_y}{m}-\frac{k_y}{m}\dot y $

$ \tag{5c}\ddot z = \frac{F_z}{m}-\frac{k_z}{m}\dot z-g $

$ \tag{5d}\ddot\theta = \frac{1}{I_y}(\tau_1+(I_z-I_x)\dot\phi \dot\psi) $

$ \tag{5e}\ddot\phi = \frac{1}{I_x}(\tau_2+(I_y-I_z)\dot\theta \dot\psi) $

$ \tag{5f}\ddot\psi = \frac{1}{I_z}(\tau_3+(I_x-I_y)\dot\phi \dot\theta) $

2 SMC Design for the Quadcopter UAV

In this section, we design the sliding mode controller for system (5). First we define the following tracking errors $ e_{\{1,\cdots,12\}} $ as:

$ \begin{matrix} e_1 = x-x_d & e_2 = \dot x-\dot x_d\\ e_3 = y-y_d & e_4 = \dot y-\dot y_d\\ e_5 = z-z_d & e_6 = \dot z-\dot z_d\\ e_7 = \theta-\theta_d & e_8 = \dot \theta-\dot \theta_d\\ e_9 = \phi-\phi_d & e_{10} = \dot \phi-\dot \psi_d\\ e_{11} = \psi-\psi_d & e_{12} = \dot \psi-\dot \psi_d \end{matrix} $

(6)

Then their corresponding sliding surfaces $ s_j $ , $ j\in \{x,y,z,\theta,\phi,\psi\} $ can be designed as follows^[25]:

$\tag{7a} s_x = e_2+\lambda_1 e_1 $

$\tag{7b} s_y = e_4+\lambda_2 e_3 $

$ \tag{7c}s_z = e_6+\lambda_3 e_5 $

$ \tag{7d}s_\theta = e_8+\lambda_4 e_7 $

$\tag{7e} s_\phi = e_{10}+\lambda_5 e_{9} $

$ \tag{7f}s_\psi = e_{12}+\lambda_6 e_{11} $

where $ \lambda_{\{1,\cdots,6\}} $ are positive constant. By taking derivation of sliding surfaces with respect the time:

$ \tag{8a}\dot s_x = \frac{{{F_x}}}{m}-\frac{{{k_x}}}{m}\dot{x} + {\lambda _1}{e_2}- {{\ddot x}_d} $

$\tag{8b} \dot s_y = \frac{{{F_y}}}{m}-\frac{{{k_y}}}{m}\dot{y} + {\lambda _2}{e_4}- {{\ddot y}_d} $

$ \tag{8c}\dot s_z = \frac{{{F_z}}}{m}-\frac{{{k_z}}}{m}\dot{z} + {\lambda _3}{e_6}-g- {{\ddot z}_d} $

$ \tag{8d}\dot s_\theta = \frac{1}{I_y}(\tau_1+(I_z-I_x)\dot\phi \dot\psi)+ {\lambda _4}{e_8}- {{\ddot \theta}_d} $

$ \tag{8e}\dot s_\phi = \frac{1}{I_x}(\tau_2+(I_y-I_z)\dot\theta \dot\psi)+ {\lambda _5}{e_{10}}- {{\ddot \phi}_d} $

$ \tag{8f}\dot s_\psi = \frac{1}{I_z}(\tau_3+(I_x-I_y)\dot\phi \dot\theta)+ {\lambda _6}{e_{12}}- {{\ddot \psi}_d} $

For a nominal controller denoted by $ \hat{F}_j $ , we have $ \dot s_j = 0 $ . i.e:

$ \tag{9a}\hat{F}_x = \hat m(\ddot x_d-\lambda_1e_2)+{{{\hat k_x}}}\dot{x} $

$ \tag{9b}\hat{F}_y = \hat m(\ddot y_d-\lambda_2e_4)+{{{\hat k_y}}}\dot{y} $

$ \tag{9c}\hat{F}_z = \hat m(g+\ddot z_d-\lambda_3e_6)+{{{\hat k_z}}}\dot{z} $

$ \tag{10a}\hat{\tau}_1 = \hat I_y(\ddot \theta_d-\lambda_4e_8)+(\hat I_x-\hat I_z)\dot\phi \dot\psi $

$ \tag{10b}\hat{\tau}_2 = \hat I_x(\ddot \phi_d-\lambda_5e_{10})+(\hat I_z-\hat I_y)\dot\theta \dot\psi $

$ \tag{10c}\hat{\tau}_3 = \hat I_z(\ddot \psi_d-\lambda_6e_{12})+(\hat I_y-\hat I_x)\dot\phi \dot\theta $

Then we take the controllers as follows:

$ \tag{11a}F_x = \hat{F}_x-\beta_x\mathrm{sign}(s_x) $

$\tag{11b} F_y = \hat{F}_y-\beta_y\mathrm{sign}(s_y) $

$ \tag{11c}F_z = \hat{F}_z-\beta_z\mathrm{sign}(s_z) $

$ \tag{11d}\tau_1 = \hat{\tau}_1-\beta_\theta\mathrm{sign}(s_\theta) $

$ \tag{11e}\tau_2 = \hat{\tau}_2-\beta_\phi\mathrm{sign}(s_\phi) $

$ \tag{11f}\tau_3 = \hat{\tau}_3-\beta_\psi\mathrm{sign}(s_\psi) $

For instance in order to obtain $ \beta_x $ we use the Lyapunov function in the form of $ V(s_x) = \dfrac{1}{2}s_x^2 $ . Using equation (11) we have:

$ \begin{split} \dot V(s_x) =& s_x \dot s_x = s_x\left(\frac{{{F_x}}}{m}-\frac{{{k_x}}}{m}\dot{x} + {\lambda _1}{e_2}- {{\ddot x}_d}\right) =\\ & s_x\left((\frac{\hat m}{m}-1)(\ddot x_d-\lambda_1e_2)+\frac{\hat k_x-k_x}{m}\dot{x}-\frac{\beta_x}{m}\mathrm{sign}(s_x) \right)=\\ & \frac{1}{m}(\hat m-m)(\ddot x_d-\lambda_1e_2)s_x+\frac{\hat k_x-k_x}{m}\dot{x}s_x-\frac{1}{m}\beta_x|s_x| \leqslant\\ & \frac{1}{m} (\overline{ m}-\hat m)(\lambda_1|e_2|+\bar a_x)|s_x|+\frac{\overline{ k_x}-k_x}{m}\bar{v}_x |s_x|-\frac{1}{m}\beta_x|s_x| \leqslant\\ & \frac{\beta^0_x}{m} |s_x| \end{split} $

(12)

if $ \beta_x $ satisfies the following inequality:

$ \beta^0_x+(\overline{ m}-\hat m)(\lambda_1|e_2|+\bar a_x)+(\overline{ k_x}-\hat k_x)\bar v_x\leqslant \beta_x $

(13)

for some $ 0<\beta^0_x $ , where $ \bar a_x $ and $ \bar v_x $ are the absolute value of the maximum of desired acceleration $ \ddot x_d $ and velocity $ \dot x_d $ , respectively. I.e. $ \bar a_x: = \max |\ddot x_d| $ and $ \bar v_x: = \max |\dot x_d| $ . The same procedure is applied for $ F_y,F_z, \tau_1,\tau_2,\tau_3 $ .

3 Simulation Results

In this section, we provide simulation results for trajectory tracking of the UAV. The software environments Matlab Simulink and Universal Mechanism are used for mathematical modeling. In Matlab Simulink, it is convenient to model the dynamics by integrating the equations of motion Universal Mechanism allows. The joint use of these packages not only simplifies the implementation of the model but also allows one to obtain parallel results, the comparison of which can confirm the correctness of the proposed algorithm.

The model in Matlab Simulink (see Fig. 2) consists of four blocks: input voltages, voltage limiters, the model of motors with propellers, the model of the quadrocopter. The input signals of the model are the control voltages to the motors corresponding to control inputs $ (F, \tau_1, \tau_2, \tau_3) $ , the output signals are flight coordinates $ (x, y, z) $ and angles $ (\theta, \phi, \psi) $ . The control voltages for the motors are supplied via the signal distributor and voltage limiters. For the selected motor (mode: $ {\rm{X}}2212 {\rm{KV}}980 $ ) with a linear dependence of the speed on the control voltage in the operating speed range, the input voltages are limited to values from $ 0 \;{\rm{V}} $ to $ 11.1\; {\rm{V}} $ . In the Universal Mechanism package, the quadcopter is presented as a rigid cruciform body with four motors and four propellers with $ 6 $ DOF. And the inertial characteristics are the same as those specified in the Matlab package.

The structure of the control system is shown in Fig. 3, where the numbers indicate: 1—given trajectory; 2—corrective devices (CU) of the trajectory control subsystem; 3—coordinate converter; 4—regulator of the orientation and stabilization subsystem; 5—signal distributor; 6—voltage limiter; 7—model of the propeller group; 8—quadrocopter model.

For the reference trajectory, we use the following desired path:

$ \tag{14a}x_d = 0.05\cos(0.2t) $

$ \tag{14b}y_d = 0.05\sin(0.2t) $

$ \tag{14c}z_d = 0.2+0.001t $

$ \tag{14d}\psi_d = {\text{π}}/3 $

In order to show robustness, we use $ 20\% $ uncertainty in the mass and $ 10\% $ uncertainty in the stiffness constants (the stiffness constant is defined as the rotational inertial of the $ i $ axis $ I_i $ times the torque coefficients $ c_\tau $ ). The parameters values of the UAV and the SMC controller constants are given in Table 1.

Fig. 4 shows the desired trajectories and UAV trajectories. As it can be seen, the UAV trajectories track the desired trajectories. Fig. 5 illustrates the desired path and UAV path and Fig. 6 shows the control inputs. Fig. 7 shows sliding surfaces $ s_x $ , $ s_y $ , $ s_z $ and $ s_{\psi} $ . As it can be seen, all surfaces are stable and they converge to zero. Finally, Fig. 8 illustrates the tracking errors.

4 Real Implementation

To further validate the reliability of the proposed algorithm, we use an actual quadrotor UAV for indoor flight experiments.

The detailed system configuration is described in Fig. 9. In this flight test, four VICON cameras with a capture range of $ 3\;{\rm{m}}\times3\;{\rm{m}}\times1.5\;{\rm{m}} $ are used. A quadrotor-type Pelican platform from AscTec (Ascending Technology) is used as a multi-rotor flying robot platform. Pelican consists of HLP (High-Level Processor) and LLP (Low-Level Processor). In the HLP, guidance and SMC control algorithms are executed based on position data such as GPS, user design code can be compiled using the SDK, and motor control signals are transferred from the compiled implementation code to the LLP. The LLP transfers IMU data and motor control signal data from the HLP to the motor drivers. The Pelican flight platform is equipped with sensors such as a gyroscope, accelerometer, geomagnetometer, GPS, and barometer. In this study, the gyroscope and accelerometer inertial data (IMU Data) are used to calculate the velocity and angular velocity information using the proposed SMC algorithm, and the multi-rotor position and altitude control information is obtained by processing VICON data. The information obtained from VICON is transmitted from the GCS (Ground Control System) to the quadrotor UAV via ZigBee wireless communication and used as input values for flight control. Each IMU, GPS, PC, UART, and motor data are communicated in real-time via GCS and ZigBee. The IMU data includes measured angular velocity and acceleration data, and the GPS includes latitude, longitude, altitude, and speed data. The PC receives signals from the remote control, and UART receives command values from the GCS input. Position and yaw information from VICON are transferred from MATLAB to the Simulink-based GCS, which classifies the data and transmits them to the vehicle. The parameter values for the real-time test flight are: total weight is $ 1.63 \;{\rm{kg}} $ , center distance from the rotor is $ 0.225 \;{\rm{m}} $ , thrust efficiency is $ 3.6 \times 10^{-2} $ , and drag efficiency equals to $ 1.8 \times 10^{-2} $ .

For the real experimental testing, we apply the proposed SMC algorithm to simulate a complete flight mission, including the takeoff and landing of the UAV. Fig. 10 shows the results of the results when tracking the trajectory consisting of sections: AB—takeoff; BC—uniform rectilinear motion; CD—uniform accelerated rectilinear motion; DE—uniform movement around the point $ (5,10,15) $ in a horizontal circle with a radius of $ 5 \;{\rm{m}} $ ; EF—uniform rectilinear motion; FG—uniform movement around the point $ (-5,10,10) $ in a horizontal circle with a radius of $ 10 \;{\rm{m}} $ ; GH—landing. As can be seen, the results are generally satisfactory as well.

As previously described, in the case of multi-rotor UAVs, which are susceptible to external environmental uncertainties, stable flight performance is required even if there are unforeseen external disturbances. Therefore, we also need to validate the robustness of the designed SMC controller. An experiment is conducted using a digital dynamometer to generate disturbances on the quadcopter from the outside. Perturbations are added in the form of impulse shocks through the dynamometer to produce the effect of sudden gusts of wind, and the disturbance forces are measured as digital data values. The perturbation experiments are performed in the $ X $ and $ Z $ directions (horizontal and vertical directions) of the vehicle, respectively.

The perturbation to the $ X $ -direction is input about $ 15\;{\rm{s}} $ after the flight, and the maximum peak of the force is set to $ 4.7\; {\rm{N}} $ . The applied perturbation results in a position error of up to $ 0.45\;{\rm{m}} $ in the $ X $ -direction but is confirmed to return to the reference command value of $ 0 \;{\rm{m}} $ within $ 1 \;{\rm{s}} $ by the SMC controller. The results are shown in Fig. 11. The perturbation to the $ Z $ -direction is input at around $ 41 \;{\rm{s}} $ . The maximum peak of the force is set to $ 7.6\; {\rm{N}} $ . The perturbation to the $ Z $ -direction results in a height error of up to $ 0.3 \;{\rm{m}} $ , but is confirmed to return to the existing standard command value of $ 0.7 \;{\rm{m}} $ within $ 1 \;{\rm{s}} $ . The results are shown in Fig. 12. In the case of the horizontal disturbance, an instantaneous error of $ 0.45\;{\rm{m}} $ occurs when a disturbance force of $ 4.7\;{\rm{N}} $ is applied, and in the case of the vertical direction, an error of about $ 0.3\;{\rm{m}} $ is confirmed to have occurred at a disturbance force of $ 7.6\;{\rm{N}} $ . As shown in the experiment, the SMC controller can stabilize the quadrotor UAV again within $ 1 $ second in the presence of perturbations, showing the excellent robustness in the $ X $ and $ Z $ directions.

5 Conclusion

In this work, we develop a simple sliding mode control algorithm for robust 3D trajectory tracking of a quadcopter UAV, where the controller is able to track desired paths in the presence of uncertainties. Simulation experiments with Matlab Simulink and real flight experiments are conducted to verify the performance of the proposed SMC algorithm. In particular, to verify the robustness against disturbances such as gusts of wind in an uncertain environment of quadcopter operation, disturbances are applied to the horizontal and vertical axes of the UAV and the performance is verified in terms of attitude stability. The designed SMC controller shows immediate re-stabilization of the UAV with excellent resilience in both the horizontal and vertical directions. Both simulation and real testing results show the efficiency of the SMC algorithm and the possibility of its implementation.