Fractional Order Sliding Mode Control of a Class of Second Order Perturbed Nonlinear Systems : Application to the Trajectory Tracking of a Quadrotor

Arturo Govea-Vargas 1,*, Rafael Castro-Linares 1, Manuel A. Duarte-Mermoud 2,3 , Norelys Aguila-Camacho 4 and Gustavo E. Ceballos-Benavides 5 1 Department of Electrical Engineering, CINVESTAV, Av. IPN 2508, Ciudad de México 07360, México; rcastro@cinvestav.mx 2 Advanced Mining Technology Center, University of Chile, Av. Tupper 2007, Santiago 8370451, Chile; mduartem@ing.uchile.cl 3 Department of Electrical Engineering, University of Chile, Av. Tupper 2007, Santiago 8370451, Chile 4 Departmento de Electricidad, Universidad Tecnológica Metropolitana, Av. José Pedro Alessandri 1242, Santiago 7800002, Chile; norelys.aguila@utem.cl 5 Facultad de Negocios, Ingeniería y Artes Digitales,Universidad Gabriela Mistral, Av. Ricardo Lyon 1177, Providencia, Santiago 7510549, Chile; guceballos@gmail.com * Correspondence: arturogovea@gmail.com; Tel.: +52-(55)-57473791


Introduction
Recently, the design of fractional order controllers has become one of the most exciting topics in control theory leading to interesting applications in the control of physical systems such as suspension systems, permanent magnet synchronous motors, power electronic systems and unmanned aerial vehicles, among others.The idea of fractional order control was first proposed by [1] introducing a robust fractional order control scheme.Also Podlubny [2,3] introduced the well known fractional order proportional-integral-derivative controller.Then, many fractional order controllers have been proposed in the literature, including tilt integral derivative (TID) controllers [4], fractional order lead-lag compensators [5,6], fractional order optimal controllers [7,8], and fractional order adaptive controllers [9,10].
On the other hand sliding mode control (SMC) is a well known control technique that is applied to both linear and nonlinear systems allowing to deal with uncertainties [11].However, the monotonous switching feedback produced by these controllers causes high frequency chattering in the control signals, which leads to undesirable loads on control actuators [12,13].In recent years, fractional order sliding mode control (FOSMC) has been employed to overcome these drawbacks.For example, in [14], a control strategy called fractional order terminal sliding mode control (FOTSMC) is designed for a class of uncertain dynamical systems; a fractional order switching surface is proposed to satisfy a sliding condition and numerical simulations are provided that show the finite time stability of the closed-loop system.Also, in [15] a FOSMC is proposed for the control of a single-link flexible manipulator constructing a switching surface based on fractional derivatives; such a controller achieves better performance with small control chattering and robustness with respect to external load disturbances and parameter variations.A FOSMC for antilock braking systems (ABS) is proposed in [16,17] to regulate the slip to a desired value; these controllers deal with the uncertainties in the ABS and track the desired slip faster than conventional integer order sliding mode control (IOSMC).Aditionally, FOSMC has been used for speed control of permanent magnet synchronous motors [18], for vibration suppression of uncertain structures [19], for control of fractional order chaotic systems [20,21] and so on leading to better control performance.
In [22], the stability of a FOSMC has been studied showing that a stable reaching law, in the fractional order case, corresponds to a stable reaching law, in the integer order case.Also, a sufficient condition is derived to guarantee global attraction to the sliding surface when a fractional order linear system with uncertainties is considered.Motivated by this result, in this paper FOSMC design is proposed to solve a trajectory tracking problem for an integer order nonlinear system where an unknown additive perturbation term is present.A sufficient condition is also given that assures the reaching of a suitable sliding surface; on this surface trajectory tracking is performed.The FOSMC proposed is applied to the trajectory tracking problem of a quadrotor.For vertical take-off, trajectory tracking and landing, an error sliding surface is considered.Numerical simulations show the robustness of the FOSMC with respect to additive uncertainties (the error tracking is bounded to a small region around zero) and the chattering in the control signals is reduced when compared to a integer order sliding mode controller.
The paper is organized as follows.In Section 2, fractional order differentiation is briefly presented together with the FOSMC proposed and the sufficient condition for the attractiveness of the sliding surface.Section 3 presents the FOSMC design for the trajectory control of a quadrotor model.Some simulations results are given in Section 4. Finally, in Section 5 some conclusions are presented.

FOSMC of an Integer Second Order Perturbed Nonlinear System
Fractional calculus is mainly the calculus of derivatives and integrals of fractional order that can be real or complex [23] and it has many applications in science and engineering.In this article a Sliding Mode Control of Fractional Order is proposed for the trajectory tracking of a quadrotor.Among the definitions of fractional derivative the most cited are the Riemman-Liouville and Caputo definitions.Here, the Caputo definition is used which is expressed as where Γ is the so called gamma function defined as Γ(n) = ∞ 0 t z−1 e −t dt and n − 1 ≤ α ≤ n with n being an integer.From this definition, one can notice that when α < 0 one has integration and when α > 0 one has differentiation; when α = 1.0 one has the first order derivative of a function.The gamma function generalizes the notion of a number raised to an arbitrary real number (or complex number).Thus the gamma function somehow involves the computation of a fractional derivative.The definition of a fractional derivative, as stated above, can not be used in practice, thus numeric methods such as the one based on the Grünwald-Letnikov approach is commonly used.
Let us consider a perturbed nonlinear system described by a second order differential equation of the form where χ(•) ∈ R, f χ (•) and g χ (•) are bounded functions of their arguments and u(•) ∈ R is the system input.∆ χ (•, •) ∈ R is a disturbance term which includes uncertainties associated to the modelling of a process and external perturbations.The following assumption is made on for all χ(•) which are solutions of ( 2) and all t ≥ 0 with κ ∆ being a real positive constant different from 0.
It is desired that the signal χ tracks a given reference signal χ d , despite the presence of the disturbance term ∆ χ .Thus a tracking error e is defined as follows: In accordance to the sliding mode control methodology [24] a so called switching function s χ is chosen as where λ χ is a real constant parameter.The switching function (5) defines the sliding surface Thus, λ χ is selected in such a way that the first order differential Equation ( 6) has a solution that exponentially converges to zero; as a consequence χ also converges to χ d exponentially.
Let us also consider the nominal nonlinear system associated to the perturbed system (2), this is with ∆ χ (χ, t) = 0, for all χ ∈ R and all t > 0, given by In order to attract the dynamics of system (2) to the sliding surface (6) and based on the strategy proposed in [25], the (1 + β) fractional order derivative of s χ is set to be where σ χ and µ χ are positive real constants and sgn(•) is the signum function.The derivative of (8) to the order (−β) is now taken (this is equivalent to integrate (8) to the order β), leading to On the other hand, from (5) and the nominal nonlinear system (7), one has that Thus by combining Equations ( 9) and (10), the following FOSMC is obtained: where In fact, the FOSMC Equations ( 11) and ( 12) can achieve the attraction of the perturbed dynamics (2) to the sliding surface (6) under a sufficient condition on the size of the perturbation term ∆ χ , the parameters σ χ ,µ χ and the discontinuous term D −β sgn(s χ ), as stated in the following result.
Theorem 1.Consider the perturbed nonlinear system (2) together with Assumption A1.If the parameters σ χ , µ χ and the bound κ ∆ satisfy then the FOSMC Equations ( 11) and ( 12) assures the attractiveness of the perturbed system to the switching surface s χ = 0 defined by (6).
Proof.Let us consider the following Lyapunov function candidate: which is positive definite.The time derivative of V χ has the form substituting the dynamics of the perturbed system (2) into (15) together with the FOSMC Equations ( 11) and ( 12) leads to Since s χ = |s χ |sgn(s χ ), Vχ can also be written as Also 0 < β < 1, thus sgn(D −β sgn(s χ )) = sgn(s χ ) [25], then (17) takes the form Majoring this last expression and using the bound κ ∆ in (3), one gets Thus if Condition (3) is satisfied, Vχ < 0 and the convergence to the surface s χ = 0 is accomplished.
Remark 1.The FOSMC ( 11)-( 12) compensate the effect of additive perturbations in a second order nonlinear system described by (2).In fact, the main idea to define a sliding surface (6) together with its (1 + β) fractional order derivative is that as soon as the disturbance appears, s χ becomes nonzero and, at that time, the discontinuous control term becomes active and rejects the disturbance.Also, and similar to the control scheme proposed in [26], expression (19) shows that a stronger attraction to s χ = 0 is achieved and that the attraction is higher for any s χ with µ χ = 0 than with µ χ = 0.In particular, the parameter µ χ is not used in [26] where a fractional order linear system is controlled.

Quadrotor's Dynamic Model
Let us consider an inertial reference frame and a body fixed reference frame to specify the position, velocity and acceleration of the quadrotor.The world frame W, is defined by axes X, Y and Z, with Z pointing upward.The body frame B, xB, yB, zB, is attached to the center of mass of the quadrotor.A picture of a quadrotor with the corresponding reference frames is shown in Figure 1.The Newton-Euler approach was used to obtain the dynamical behaviour of the quadrotor [27], which considers a quadrotor as a rigid body.Let us consider the generalized coordinates of the quadrotor denoted as q = (ξ, η), where ξ = (x, y, z) is the position vector, which goes from the center of mass of the quadrotor relative to the fixed reference frame (X, Y, Z), and η = (φ, θ, ψ) are the Euler angles of the quadrotor.These angles have the following names and they represent a rotation around their corresponding axis: φ is the pitch angle around xB axis, θ is the roll angle around yB axis, and ψ is the yaw angle around zB axis.The angular velocity vector Ω = (p, q, r) is related to η by the equivalent form Ω = W η η, where is a standard kinematic matrix.Also R(φ, θ, ψ) SO(3), is the so called rotation matrix that represents the quadrotor's orientation that is relative to the fixed reference frame.Using the compact notation c(α) = cos(α) and s(α) = sin(α), this matrix is given by .
Additionally, the inertia matrix is given by where I xx , I yy and I zz are the moments of inertia of the quadrotor about its corresponding axis, m is the mass of the quadrotor and g = 9.81 [m/s 2 ] is the gravitational constant.Finally, the dynamics of a rigid body under external forces take the form where Ω is the anti-symmetric matrix of Ω, also called the wedge operator defined as and × denotes the cross product.T f is the total thrust applied to the quadrotor, this is i=1 ω 2 i and represents the first control signal, furthermore τ = (τ φ , τ θ , τ ψ ) is the vector of control signals for roll, pitch and yaw angles, respectively.
It is possible to stablish a mapping between the control signals u 1 , u 2 , u 3 and u 4 and the angular velocity of each quadrotor's electric motor ω 1 , ω 2 , ω 3 and ω 4 .In fact, each motor produces a vertical force F i and a moment M i that is related to the motor angular velocity ω i in the form [29] where k F and k M are positive constants that can be experimentally determined.The relation between the motor angular velocities and the control signals is then given by [29] with K MF being an invertible matrix of the form where L is the distance from the rotation axis of the motors to the center of the quadrotor.

Trajectory Tracking of a Quadrotor
Let us first consider the z position dynamics in the quadrotor model ( 24)-( 29), more precisely Equation (26).These dynamics have the form (2) with x = z, u = u 1 and The tracking error e also takes the form e = e z = z − z d with χ d = z d being the reference signal.In accordance with the methodology described in the previous section, the FOSMC that assures to have convergence to z towards z d despite the presence of the perturbation term ∆ z is given by where with s z = s χ = ė + λ χ e = ėz + λ z e z as switching function (λ χ = λ z ) and σ z = σ χ , µ z = µ χ in (35).∆ z = ∆ χ is assumed to satisfy (13), with a bound κ ∆ = κ ∆z .
In order to design the control signals u 2 and u 3 , a desired reference signal is first defined for the angles φ and θ, this is φ d and θ d , respectively.For doing this, we adopt the small angle approximation widely used in the literature of unmanned aerial vehicles (see for example [28,30]).Under this assumption, the dynamics ( 24), ( 25), ( 27)-( 29) take the form ẍ ≈ tan(θ)(g ÿ ≈ −tan(φ)(g where ∆ x and ∆ y are unknown perturbation terms bounded by the real positive constants κ ∆ x and κ ∆ y respectively affecting the dynamics of x and y. It is now considered that tan(φ) and tan(θ) are "virtual" input signals, this is u x = tan(θ) and u y = tan(φ) in (36) and (37).Similar to the previous reasoning for the design of the control signal u 1 , the tracking errors e x and e y are defined as where x d and y d are the reference signals for x and y, respectively.The following switching functions are then chosen: where, as before, λ x and λ y are real constant parameters selected in such a way that the first order linear differential equations defined by the sliding surfaces have solutions that exponentially converge to zero.In this same way, the (1 + β) fractional order derivative of s x and s y are set as x , s where σ x , σ y , µ x and µ y are positive real constants.The (−β) order derivative of Equations ( 42) have the form (43), this is Also from the dynamics (36) and (37), with ∆ x = ∆ y = 0, one has that Equation ( 10) is given by ṡy = ëy + λ y ėy = −tan(φ)(g The desired signals θ d and φ d can then be computed from ( 44) and (45) when the signals tan(θ) and tan(φ) are considered to be virtual input signals, as mentioned before.This consideration allows to have the following expressions for θ d and φ d : where Notice now that, when θ and φ are managed to approach θ d and φ d as soon as possible, the FOSMC strategy discussed above allows to reach the sliding surface s x = 0 and s y = 0 where the tracking errors e x and e y tend to zero exponentially.In this paper, a simple proportional and derivation action is implemented to allow that θ → θ d , φ → φ d and ψ → ψ d , this is where e φ = φ − φ d , e θ = θ − θ d and e ψ = ψ − ψ d are the tracking errors in the angles φ, θ and ψ respectively.Notice in particular, that φ d is set to 0 in (51).

Simulation Results
Some simulations were carried out to evaluate the performance of the control scheme proposed.The simulations were made in MATLAB by integrating the dynamic equations of the quadrotor ( 24)-( 29) using a script that implements a Runge-Kutta method of fourth order.In order to calculate the fractional derivative used in Equations ( 35), ( 47) and (48) of the control scheme, the Grünwald-Letnikov numeric method is used, implementing the adaptive time step memory algorithm described in [31].Different values of the fractional order β were tested, in order to show the performance of the control scheme.Also the signum function was implemented using the expression where the small quantity > 0 is a real positive constant.
A piecewise trajectory was used as a tracking reference in the space as follows: the first 5 s are used for the take-off of the quadrotor until it reaches hover at 1 m, then a 2 m long lemniscata trajectory is followed once, with a period of 30 s per return, finally landing from 1 m to the floor is performed within the last 5 s.
A perturbation term was introduced at the right hand side of the quadrotor dynamic model ( 24)-( 29) at t = 13 s in order to simulate the effect of weather conditions, such as wind variations.This perturbation term used was given by . The model parameters used in the simulation are shown in Table 1 while the control parameters are given in Table 2.The model parameters correspond to a real quadrotor that will be used in real experimentation.The quadrotor has 23 cm of arm longitud; this is the length from the center of the quadrotor to the place where each motor is mounted.The motors used are the Emax MT2213-935KV; each motor can thrust up to 670 grams.An IMU um7 from RedShiftLabs is used and an Optitrack System will be used to locate the quadrotor in the space.The trajectory of the quadrotor in the space is shown in Figure 2. Nine fractional orders are being plotted (from β = 0.1 to 0.9); the colors used in the figures for each value of β are listed in Table 3.Notice that the evolution of the variables with some values of β are indistinguishable, and the light green color, corresponding to the fractional order β = 0.9, is predominant.A projection of the trajectory in the x-y plane is depicted in Figure 3.The desired trajectory is not perfectly tracked, due to the presence of the disturbance, however the deviation error is small (see also Figures 4-9), having a better performance with β = 0.5.Besides, by making a zoom into the projection, for a time interval from t = 0 s to t = 40 s, as it is shown in Figure 10, one can notice that the distinct trajectories for the distinct values of fractional order are very close one from another.The position tracking error in the x, y and z directions, this is e x = x − x d , e y = y − y d and e z = z − z d are shown in Figures 4-6.In fact, the error raises to almost 2 cm after the introduction of the disturbance at second 13, in the coordinates x and y, and 5 mm in the z direction when β = 0.5.Notice that a significant overshoot and a much greater tracking error is obtained with β = 0.9.The Euler angle errors are depicted in Figures 7-9.One can notice that they are within a 0.6 radians absolute error bound for values of β between 0.1 and 0.8; the corresponding Euler angles evolve within a bound of 0.3 radians for β = 0.5, as it is shown in Figure 11.Thus, the small angle approximation is assured through a great part the trajectory, except for take off and landing.The four control signals are depicted in Figures 12-15.It is important to mention that these control signals were limited to the values shown in Table 4 [32].This values are similar to the ones that correspond to the quadrotor that will be used in the real experiments.Since in a real application the absolute value of those signals are limited due to the actuators (motors).
Roll Torque (u2) Roll Torque (u2) with =0.5  Besides, the nine control signals, corresponding to the nine fractional orders, are indistinguishable one from the other, but when β = 0.1 there is not chattering at all; higher fractional orders lead to a chattering to be present.This is even more noticeable when the fractional order tends to 1 (β → 1).A comparison between the vertical thrust for a fractional order of β = 0.5 and for a fractional order of β = 0.9 is depicted in Figure 16.

Simulations Varying the Parameters σ χ
Some simulations were carried out varying the parameters σ z , σ y and σ x in the range [0.1, 0.01], and fixing the fractional order to β = 0.5, in order to analyze the behaviour of the quadrotor system when the control parameters are given different values; the parameter σ χ was selected since the closed-loop behaviour is more sensitive to changes in σ χ .From the results obtained (Figures 17-19) one can observe that with lower values of σ z , σ y and σ x , the chattering in the position error signals is reduced.
With respect to Euler Angles, the simulations show some chattering with higher values of σ χ while the disturbance does not significantly affect the Euler Angles; the disturbance appears at second 13 and it is kept for the rest of the simulation, the errors in Euler Angles are shown in Figures 20-22.It can also be noticed that the control acts in the vicinity of the second 13 where the correction needs to be made.

Comparison with the IOSMC
In order to compare the performance of the FOSMC design presented here, an IOSMC was implemented for the quadrotor's simplified model.Using the small angle approximation discussed above this controller is given by Equations ( 26), ( 27)-( 29) with where, as before, λ z , λ y , λ x , σ z , σ y , σ x , µ z , µ y and µ x , are real coefficients chosen in such a way that the corresponding solutions of the first order linear differential equations defined by the sliding surfaces asymptotically converge to zero.The simulation carried out with this controller corresponds to the one using the FOSMC with the control parameters given in Table 5.The tracking position errors and the Euler angles tracking errors are shown in Figures 27 and 28, respectively, while the 3dPlot is depicted in Figure 29.One can notice that the tracking errors have a good performance when the disturbance appears.However, undesirable chattering appears in the control signals in particular with the vertical thrust u 1 , as it can be seen in Figures 30-33.The errors in Euler Angles seems to have a good performance, again when the disturbance is presented, some deviations appears, see Figure 28.
The MATLAB code generated for the simulations shown in this section are available at the link: https://drive.google.com/drive/folders/1pszDBTfKHeGENOETDpS5Kutq1zGPX1qc?usp= sharing.It should be pointed out that the initial tests which have been carried out with a real quadrotor (the model parameters used in the simulations correspond to this quadrotor) give indications of a behaviour similar to the one shown here.It is expected to have complete experiments in the near future.

Conclusions
A fractional order sliding mode control (FOSMC) was proposed for the trajectory tracking control for an integer second order nonlinear system with an unknown perturbation term.A sufficient condition was given to assure the attraction of the sliding surface where trajectory tracking is achieved.The sufficient condition found allows to adjust two parameters (σ χ and µ χ ) of the FOSMC in order to assure the attraction to the sliding surface proposed.Besides, an application of the methodology to the trajectory tracking of a quadrotor dynamical model is given.Simulation results were obtained for a given trajectory, including take-off and landing.These results were also compared with an IOSMC showing that chattering can be diminished in some of the control signals, when using FOSMC, and in the presence of additive perturbations.The results obtained are promising and further research is being carried out together with experimentation with a real quadrotor.

Figure 1 .
Figure 1.Coordinate Systems and moment/forces acting on the quadrotor.

Figure 2 .
Figure 2. 3d Plot of the position of the trajectory being tracked.

Figure 3 .
Figure 3. Projection in the x-y plane of the trajectory.Orange-Desired trajectory, Else-Trajectory being tracked.

Figure 10 .Figure 11 .
Figure 10.Zoom of the projection in the x-y plane of the trajectory.

Table 4 .
Limitation of the control signals.