Particle Swarm Optimization and Cuckoo Search-Based Approaches for Quadrotor Control and Trajectory Tracking

: This paper explores the full control of a quadrotor Unmanned Aerial Vehicles (UAVs) by exploiting the nature-inspired algorithms of Particle Swarm Optimization (PSO), Cuckoo Search (CS), and the cooperative Particle Swarm Optimization-Cuckoo Search (PSO-CS). The proposed PSO-CS algorithm combines the ability of social thinking in PSO with the local search capability in CS, which helps to overcome the problem of low convergence speed of CS. First, the quadrotor dynamic modeling is deﬁned using Newton-Euler formalism. Second, PID (Proportional, Integral, and Derivative) controllers are optimized by using the intelligent proposed approaches and the classical method of Reference Model (RM) for quadrotor full control. Finally, simulation results prove that PSO and PSO-CS are more e ﬃ cient in tuning of optimal parameters for the quadrotor control. Indeed, the ability of PSO and PSO-CS to track the imposed trajectories is well seen from 3D path tracking simulations and even in presence of wind disturbances.


Introduction
Over the past few years, Unmanned Aerial Vehicles (UAVs) of type quadrotors or quadcopters have seen an increasing interest since their wide range of civilian and military applications. These quadrotors are considered the most complex flying machines due to many physical effects influencing their dynamics including aerodynamic effects, gravity, gyroscopic effects, friction, and inertia. However, they have advantages over conventional helicopters. Given that the left and the right motors rotate clockwise, and the front and the rear motors rotate counterclockwise, gyroscopic effects and aerodynamic torques tend to cancel in trimmed flight.
Many researchers modeled the quadrotor using the formalisms of Newton-Euler [1,2] and Euler-Lagrange [2][3][4]. For both methods, the quadrotor modeling is regarded as a delicate task and the model obtained using these approaches is strongly nonlinear, fully coupled, under-actuated (Six Degrees of Freedom System (6-DOF) with only four actuators) and dynamically unstable with complex behavior. In this work, we adopt the most used formalism of Newton-Euler to define the dynamics of the quadrotor.
Regarding the quadrotor's control, there are several types of research done both in linear and nonlinear control methods [5][6][7]. For practicality, the linear controllers, especially PID (Proportional, Integral, and Derivative)/PD (Proportional and Derivative)/ PI (Proportional and Integral), are easy to design and simple to tune, which perform better in practical implementation [8]. In literature, PIDs are the most used for quadrotor control [2,9,10]. However, tuning the parameters of these controllers is very important to get best performances. In this scope, classical methods can be used. They include the methods of Ziegler-Nichols [11], Graham-Lathrop [12], Naslin [13], Cohen-Coon, Reference Model (RM) [14,15], and many other techniques.
Nature-inspired optimization algorithms can effectively resolve complex problems compared to classical and statistical methods. Some of these algorithms proposed in the literature are Genetic Algorithm (GA) [16,17], Particle Swarm Optimization (PSO) [18], Ant Colony Optimization (ACO) [19], Artificial Bee Colony (ABC) [20], and Cuckoo Search (CS) [21,22]. GA is inspired by the crossover, mutation, and recombination principles of genetics. PSO algorithm is inspired by social behavior patterns of organisms that live and interact within large groups. It incorporates swarming behaviors observed in flocks of birds and schools of fish. ACO is inspired by ants' behavior in seeking food sources. The ABC algorithm is based on the natural behavior of honeybee swarm, and CS algorithm is motivated by the aggressive breeding of a bird called 'cuckoo'. The PSO method has shown superior performances. It has features such as a straightforward algorithm, rapid convergence, easy implementation, and is computationally efficient when applied to a diverse set of optimization problems [15]. In CS, immigration and environmental specifications have the advantage to help cuckoos' groups to converge and reach the best places for breeding and egg laying [21]. Another advantage of CS compared to PSO and GA is that it uses a smaller number of parameters to be tuned, which makes it more adaptable [21]. Due to these advantages, PSO and CS are still successfully applied in recent works for control problems [23][24][25], and PSO is especially used in tuning PID controllers for quadrotor's control [26,27]. Therefore, regarding PSO and CS programs, there are many common points. Both use the same initialization matrices of particles in n-dimensional research space for PSO and nests with D eggs (problem's dimension) for CS. The fitness function evaluates the quality of particles and nests in the same way. Then, the velocities of displacement of particles and nests are calculated in different ways for PSO and CS, which generates differences in the convergence speeds of both algorithms. For these reasons, it would be very interesting to analyze the social thinking and local search abilities of PSO and CS in quadrotor control.
The main contribution of our paper lies on using PSO and CS to improve results concerning the quadrotor full control. However, the CS algorithm suffers from a low convergence speed, since it uses a fixed step size over generations. To overcome this problem, Naik et al. [28] proposed an adaptive cuckoo search algorithm (ACS) with adaptive step size. The ACS converges to a near optimum solution faster than CS algorithm. However, the problem of ACS algorithm is that it does not achieve better solutions compared to CS in most cases. In our work, we propose the cooperative PSO-CS algorithm that combines the ability of social thinking in PSO with the local search capability of CS. Such a cooperative process is likely to offer proper guidance for cuckoos to the global best positions and to ensure a balance between exploitation and exploration of the search space. Consequently, we are interested in examining if it is possible for us to improve results of quadrotor control by effectively integrating these two algorithms with complementary strengths. Our work treats position and attitude control in order to cover all topics in quadrotor control. For this objective, the control design includes six PID controllers: Three for angles (φ, θ, ψ) and three for positions (x, y, z). These controllers are optimally tuned using the heuristics PSO and CS, the cooperative PSO-CS and the classical method of Reference Model (RM). A comparative study is done to highlight the efficiency of the proposed intelligent controllers of PSO and PSO-CS. Indeed, their robustness is demonstrated in presence of wind disturbance.
The rest of this paper is organized as follows. In Section 2, Newton-Euler formalism is used to establish the quadrotor dynamic model. In Section 3, the proposed algorithms of PSO, CS, the cooperative PSO-CS and RM are presented. In Section 4, the quadrotor full control is treated using PID controllers optimally tuned with the proposed methods. In the last section, our conclusions are given.

Quadrotor Dynamic Modeling
The quadrotor is a complex flying system, strongly nonlinear, multivariable, fully coupled, and under-actuated (6-DOF and only four control inputs). Therefore, its modeling is a delicate task. For better understanding of its dynamic modeling, different working hypotheses were assumed: A rigid and a symmetrical quadrotor structure (diagonal matrix of inertia); a rigid propellers (negligible effect of deformation during rotation); the lift and the drag forces are proportional to the square of the rotational speed of the rotors (a very close approximation of the aerodynamic behavior); and the center of mass is exactly the origin of reference related to the structure.
From the configuration shown in Figure 1, the quadrotor has two pairs of rotating rotors attached to the end of a cross and the control electronics is situated in the center of the cross. The front and the rear propellers spin counterclockwise, while the left and the right ones spin clockwise, which effectively neutralizes the unwanted reactive torque and allows the vehicle to fly without overturning.
Appl. Sci. 2019, 9, x FOR PEER REVIEW  3 of 25 and a symmetrical quadrotor structure (diagonal matrix of inertia); a rigid propellers (negligible effect of deformation during rotation); the lift and the drag forces are proportional to the square of the rotational speed of the rotors (a very close approximation of the aerodynamic behavior); and the center of mass is exactly the origin of reference related to the structure. From the configuration shown in Figure 1, the quadrotor has two pairs of rotating rotors attached to the end of a cross and the control electronics is situated in the center of the cross. The front and the rear propellers spin counterclockwise, while the left and the right ones spin clockwise, which effectively neutralizes the unwanted reactive torque and allows the vehicle to fly without overturning. Notice that the absolute position of the mass center is described by the three coordinates (x, y, z) and its attitude by the three Euler's angles (φ, θ, ψ). Using Newton-Euler formalism [2], the dynamic model of the quadrotor can be expressed as: Notice that the absolute position of the mass center is described by the three coordinates (x, y, z) and its attitude by the three Euler's angles (φ, θ, ψ). Using Newton-Euler formalism [2], the dynamic model of the quadrotor can be expressed as: .. .. ..
x = cos ϕ cos ϕ sin θ + sin ϕ sin ψ m Particle Swarm Optimization (PSO) is a class of stochastic algorithms, which was developed by Eberhart and Kennedy in 1995 for hard optimization problems. The basic principle of PSO is inspired by the social behavior of animals moving in swarm as bird flocking. To search for food, each bird called particle, flies in the space of solutions, and determines its speed according to its personal experience and the information gained through interaction with other swarm members [18].
PSO is employed to intelligently select optimal parameters from N particles. The initialization matrix contains N particles dispersed in a D-dimensional search space. Like most optimization techniques, PSO requires a fitness function relevant to the particle's position. Each particle i stores its best position Pb i (t+1) and the best solution in its vicinity Pg(t+1), which is the position of the particle that has the smallest fitness value in the swarm as expressed in Equation (9). The mechanism of displacement of each particle is managed by three rules. First, the particle tends to follow the direction of its current velocity. Second, it wants to move toward its best position. Finally, it tends to move to the best position reached by its neighbors [18,29]. In fact, the new velocity matrix V ij and position matrix X ij of the particle i = {1, 2, . . . , N} in the search space of dimension D, with j = {1, 2, . . . , D} are calculated at iteration (t+1), according to Equations (10) and (11).
where • Pb ij is the best position found by the particle i; • Pg ij is the best position found by the neighborhood; • w, C 1 , and C 2 are weighting coefficients; • R 1 and R 2 are random variables generated from a uniform distribution in [0,1]; • ⊗ means element wise multiplications.
The weighting coefficients govern the process of finding the best solution and reflect the sociability of the particles. In order to move the particle towards the best solutions, selection of these coefficients must ensure a compromise between local and global exploration of the search space. In previous works [30,31], the constants C 1 and C 2 must verify the condition C 1 + C 2 ≤ 4. Also, experimental results in Reference [30] show that [1/2((C 1 + C 2 ) − 1)] < w ≤ 1. PSO search works in G iterations and the coordinates of the obtained global best position Pg are the parameters optimally tuned.

CS Algorithm
Cuckoo Search (CS) algorithm, proposed by Yang and Deb in 2009, is based on the life of the "cuckoo" bird. The basic principle of this optimization algorithm is the specific breeding and egg laying of this bird. In the habitat of other host birds, adult cuckoos lay some eggs that grow and become mature cuckoos if they are not discovered and removed by host birds. The immigration of groups of cuckoos and environmental specifications hopefully lead them to converge and reach the best places for reproduction and breeding [21].
The primary population of CS is composed of N nests representing a set of solutions, where each nest has multiple eggs (D eggs). As in nature, each cuckoo dedicates from five to 20 eggs. These values are used as the limits of the dimension D. Each cuckoo lays one egg at a time and dumps it in a randomly chosen nest. The best nests with high-quality eggs (solutions) will carry over to the next generations, where the profit of a solution is obtained by evaluation of a fitness function F of the habitat, so the profit is an array of 1 × D. The number of available host nests is fixed, and a host can discover an alien egg with probability P a from [0, 1]. In this case, the host bird can either throw the egg away or abandon the nest to build a completely new nest in a new location. For simplicity, this assumption can be approximated by a fraction P a of the N nests being replaced by new nests, having new random solutions [21].
To replace solutions in the nests with new solutions, Lévy flights mechanism is used. A new solution X i (t + 1) for cuckoo I is given by Equation (12), where ⊗ is entry wise product, similar to that used in PSO, α > 0 is the step size, and R is a random variable generated from a uniform distribution in the interval [0, 1] in order to provide a stochastic weight.
The random walk via Lévy flight is more efficient in exploring the search space, as its step length is drawn from a Lévy distribution and it is much longer in the long run. The Mantegna algorithm is used as Lévy stable distribution in most problems to decide the step length. Thus, the step length S from Mantegna algorithm can be written as represented by Equation (13), where u and v are samples drawn from Gaussian normal distribution given by Equation (14), and σ u 2 is the variance of the distributions given by Equation (15), where Γ is the gamma function [22].
The fraction P a of worse solutions are generated as given by Equation (16), where X j (t) and X k (t) are two random solutions chosen by random permutation, H is a Heaviside function, and r is a random number drawn from a normal distribution.

Cooperative PSO-CS Algorithm
The initialization matrix in PSO contains N particles dispersed in a D-dimensional search space, when in CS, this matrix is composed of N nests and each nest has D eggs (problem's dimension). In both PSO and CS, this primary population can be regarded of dimension D × N and the quality of these solutions is evaluated in the same way. The global best particle in PSO (or the best nest in CS) is the particle (or the nest) that has the smallest fitness value among all potential solutions. Then, the velocities of displacement of particles and nests are calculated in different ways in PSO and CS. However, due to the fast convergence of PSO algorithm and the low convergence speed of the CS algorithm, the cooperation between PSO and CS would be very interesting to combine the ability of social thinking in PSO with the local search capability of CS. This cooperative PSO-CS aims to modify the speed equation of displacement by combining the Lévy flights random walks of cuckoos and the movement of particles toward the global best solution in their vicinity Pg. New solutions X i (t + 1) are given by: The use of cooperative PSO-CS helps to improve searching at global and local scales, which allows a balance between exploration and exploitation of the search space. Both algorithms are combined and involved to produce results. This combination of PSO and CS capabilities contribute to increase particle diversification in the entire search space and ensure good coverage during iterations. Thanks to this cooperation, CS algorithm can be proposed to control the quadrotor. Indeed, to evaluate the efficiency of the proposed PSO-CS method, comparisons with the intelligent approaches of PSO and CS are presented in the next sections for quadrotor's control and trajectory tracking.

RM Method
Reference Model (RM) technique is a linear method that approaches the behavior of any system of "n" order to the desired behavior of a reference system of first or second order. If the reference behavior is that of a first order, the desired settling time T s is fixed and the value of τ is deduced, such that T s = 3τ. The dominant pole is then to place at −1/τ and the (n−1) other poles on left of −1/τ in the complex plane [14,15]. If the reference behavior is that of a second order, with imposed values for damping ζ and natural pulsation ω n , the two dominant poles are to place at −ζω n ± jω n √ (1 − ζ 2 ) and the (n−2) other poles on left of −ζω n in the complex plane [14,15].
The closed loop transfer function H(p), of the system including the controller and its denominator D H (p) are computed. By taking into account some desired specifications (T s , τ or ζ, ω n ), the characteristic polynomial R s (p) of the reference system is evaluated.
i f a f irst order behavior desired The parameter "α" in (18) is made much greater than 1 to ensure the dominance of the poles −1/τ or −ζω n ± jω n √ (1 − ζ 2 ). By equating the coefficients of D H (p) and D s (p), the controller parameters can be extracted [14,15].

Virtual Control
To manage the system dynamics, a simplified model is needed. By considering small rotations (cos x = 1 and sin x = x), the equations become linear and the PID controllers' inputs can be clearly seen. Consequently, the control inputs for φ, θ, and ψ are U 2 , U 3 , and U 4 , respectively. To deal with the under-actuated part, φ, θ, and U 1 are used to create three virtual control inputs so that every output will be controlled separately. Physically, these virtual controls mean that the translation motion along x, y, and z are controlled indirectly by the three common inputs (φ, θ, and U 1 ). These virtual controls are given by Equation (19) and the desired trajectories of roll (φ d ), pitch (θ d ) and U 1d are obtained by inverting Equation (19).

Control Law Design
For an appropriate control of the UAV quadrotor system described by Equations (1)-(6), PID controllers are used. These controllers are simple to design and robust. The complete control scheme for the quadrotor is represented in Figure 2. It can be clearly seen from this figure that three PID controllers are used for positions (x, y, z) and three PIDs for attitudes (φ, θ, ψ). The inputs are divided to two parts, the desired and the sensor signals. Desired signals for positions and yaw angle (x d , y d , z d , ψ d ) are fixed by the pilot or the autopilot program. Sensor signals for desired roll and pitch angles (φ d , θ d , and U 1d ) are computed from the three virtual inputs U x , U y , and U z .
Appl. Sci. 2019, 9, x FOR PEER REVIEW 7 of 25 the under-actuated part, ϕ, θ, and U1 are used to create three virtual control inputs so that every output will be controlled separately. Physically, these virtual controls mean that the translation motion along x, y, and z are controlled indirectly by the three common inputs (ϕ, θ, and U1). These virtual controls are given by Equation (19) and the desired trajectories of roll (ϕd), pitch (θd) and U1d are obtained by inverting Equation (19).

Control Law Design
For an appropriate control of the UAV quadrotor system described by Equations (1-6), PID controllers are used. These controllers are simple to design and robust. The complete control scheme for the quadrotor is represented in Figure 2. It can be clearly seen from this figure that three PID controllers are used for positions (x, y, z) and three PIDs for attitudes (ϕ, θ, ѱ). The inputs are divided to two parts, the desired and the sensor signals. Desired signals for positions and yaw angle (xd, yd, zd, ψd) are fixed by the pilot or the autopilot program. Sensor signals for desired roll and pitch angles (ϕd, θd, and U1d) are computed from the three virtual inputs Ux, Uy, and Uz. The transfer functions Cs(p) of PID controllers are given in Equation (21), where s = {ϕ, θ, ψ, x, y, z}. Their gains (Kps, Kis, Kds) should be optimally tuned for the quadrotor in order to ensure stable, fast, precise, and robust responses. For this objective, PSO, CS, PSO-CS, and RM methods are applied to realize an optimal setting of these parameters. The transfer functions C s (p) of PID controllers are given in Equation (21), where s = {φ, θ, ψ, x, y, z}. Their gains (K ps , K is , K ds ) should be optimally tuned for the quadrotor in order to ensure stable, fast, precise, and robust responses. For this objective, PSO, CS, PSO-CS, and RM methods are applied to realize an optimal setting of these parameters.

Quadrotor Intelligent Control Using PSO, CS, and PSO-CS
To control the nonlinear and fully coupled quadrotor defined by Equations (1)- (6), it is suitable to use a method directly applicable to nonlinear systems as PSO, CS, and the cooperative PSO-CS. Instead of linear PID controllers, the intelligent and adaptive controllers can look for optimal PIDs parameters for the six quadrotors' outputs [φ, θ, ψ, x, y, z].
The initial populations (N = 200) of PSO, CS, and PSO-CS are randomly dispersed in a search space that is chosen sufficient to contain all possible solutions [0, 200] and its dimension D is set to 18. To define the exploration capacity of each particle in order to improve the convergence of the process and to not exceed the velocity limits, we tested numerous values of w, and the compromise between local and global exploration in PSO was achieved for w = 0.8. The confidence coefficients C 1 and C 2 define the exploration capacity. They determine the degree of influence of the past positions of the particle itself and those of the other particles of the swarm in order to regulate the velocities relative to the best local and global positions, thus making it possible to determine the movements relative to these best positions. C 1 R 1 and C 2 R 2 are represented by matrices with dimension of 200 × swarm size. We considered that C 1 R 1 takes random values in [0, 0.8] to avoid the problem of fast convergence, when C 2 R 2 takes random values in [0, 1.2] to give more importance to the global best solution Pg (the interval of C 2 R 2 is large to that of C 1 R 1 ). In CS, the parameters used in experiments are as follows: Abandon probability P a = 0.25, and the Lévy flights settings, α = 0.1 and λ = 1.5. In PSO-CS, the same settings of PSO and CS are conserved: P a = 0.25, and the Lévy flights settings, α = 0.1, λ = 1.5 and C 2 with random values in [0, 1.2].
To evaluate the optimization performance of the quadrotor's responses using the proposed approaches, the fitness function (profit) is defined in a similar way in order to minimize the differences between the desired and the controlled outputs responses of the quadrotor. This fitness function is defined in Equation (22) as the sum of the errors E k that characterize the difference in behavior between the inputs and the outputs of the quadrotor based on the Integral Square Error (ISE), with l = 6 and k = {1, 2, 3, 4, 5, 6}.
The performances found with ISE are compared to those established using three other performance indices: Integral Absolute Error (IAE), Integral Time Absolute Error (ITAE), and Integral Time Square Error (ITSE).
The maximum number of generations for the three programs (PSO, CS and PSO-CS) is fixed as the stop criterion and set to 20.

Quadrotor Classical Control Using RM
To use the linear RM technique for quadrotor control, the complete nonlinear Equations (1)- (6) should be simplified. The most used technique to simplify these nonlinear equations is to consider weak angles variations (reduced gyroscopic effects) and retain the linear portion, the Equations (1)-(6) become: The inputs U 2 , U 3 , U 4 , U x , U y , and U z affect respectively angles of roll, pitch, and yaw, and positions x, y, and z, the reason for their usage as the PID controllers' inputs. By applying the Laplace transform, we obtain the following transfer functions for angles (φ, θ, ψ) and positions (x, y, z).
For quadrotor's attitudes and positions control, the behavior of the linear equations defined in (28) is approached to a desired reference one of a first order system. The desired behaviors for angles and positions correspond to the characteristic polynomials R s (p), where s = {φ, θ, ψ, x, y, z} and α, α 1 , α 2 are taken equal to 30, 2 and 30, respectively.
For all cases, we calculate the closed loop transfer functions Hs(p) of the six systems including the PID controllers, for which we desire search their parameters. The computed denominators D Hs (p) are given in Equation (30).

Results and Discussions
Simulations of the proposed intelligent and classical PID controllers are carried out on the quadrotor's attitudes (φ, θ, ψ) and positions (x, y, z). These PIDs are tested firstly in the case of hovering to a determined position coordinates in order study the performances of the proposed intelligent PSO, CS and PSO-CS and the classical RM controllers. Then, for two different desired trajectories, the controllers' gains are optimized through the intelligent proposed approaches. The main objective is to confirm the ability of these intelligent nature inspired methods to control the quadrotor and to efficiently realize different imposed trajectories even under the presence of disturbances. Table 1 gives the parameters of the quadrotor and Table 2 presents the suggested trajectories' dynamics and their conditions.

Comparisons of Intelligent and Classical Control Performances
To get the quadrotor to hover at a determined position coordinates, desired signals for positions [x d , y d , z d ] are considered step signals of 1 m and desired yaw angle (ψ d ) is fixed to 1 degree. Table 3 summarizes the PID controllers gains computed online using PSO, CS, and PSO-CS, and offline with the conventional RM for the quadrotor's outputs (φ, θ, ψ, x, y, z). Figures 3-6 show the step responses for the corresponding outputs [x, y, z, ψ] obtained with these methods. Figures 7 and 8 show the various obtained curves for the desired sensor signals as well as the controlled responses. Table 4 gives the obtained quadrotor's responses performances: Settling time (T s ), maximum overshoot (Mp), and the errors previously mentioned (ISE, IAE, ITAE, and ITSE).
For control systems, the RM method often respects the fixed set of specifications and gives acceptable results [15]. However, this method is linear, and it cannot be applied directly to quadrotor control (the nonlinear quadrotor model). It can be observed from the obtained responses in Figures 3  and 4 Table 4 that the search for the parameters of PIDs by RM is not optimal: Responses with overshoots, settling times, and some error values exceed the amplitude of the desired signal. These unsatisfactory performances of the linear RM method in search of PID gains can be explained by the fact that the linear transfer functions defining the system (26) already have integrations.           CS results present a considerable overshoot in yaw angle and the errors values of ISE and IAE are around 50% of the desired altitude z, which proves that ISE value weight more significantly the quality of solutions and more adequate to evaluate the fitness function. However, the PSO approach is suitable and results in better performances. Comparisons of PSO results to those of Genetic Algorithms presented in Reference [16] show that PSO reduces the settling time Ts from 2.5 s to 0.172 s for yaw angle. For altitude z, the value of Ts is increased to 1.07 s instead of 0.4 s for GA. These results illustrate that PSO and GA techniques are highly competitive to explore the optimal solution   CS results present a considerable overshoot in yaw angle and the errors values of ISE and IAE are around 50% of the desired altitude z, which proves that ISE value weight more significantly the quality of solutions and more adequate to evaluate the fitness function. However, the PSO approach is suitable and results in better performances. Comparisons of PSO results to those of Genetic Algorithms presented in Reference [16] show that PSO reduces the settling time Ts from 2.5 s to 0.172 s for yaw angle. For altitude z, the value of Ts is increased to 1.07 s instead of 0.4 s for GA. These results illustrate that PSO and GA techniques are highly competitive to explore the optimal solution  CS results present a considerable overshoot in yaw angle and the errors values of ISE and IAE are around 50% of the desired altitude z, which proves that ISE value weight more significantly the quality of solutions and more adequate to evaluate the fitness function. However, the PSO approach is suitable and results in better performances. Comparisons of PSO results to those of Genetic Algorithms presented in Reference [16] show that PSO reduces the settling time T s from 2.5 s to 0.172 s for yaw angle. For altitude z, the value of T s is increased to 1.07 s instead of 0.4 s for GA. These results illustrate that PSO and GA techniques are highly competitive to explore the optimal solution in quadrotor control. In the same project of quadrotor control, the work presented in Reference [32] suggests three intelligent controllers, multilayer neural networks controller based on reference model (RMNN), Neural Networks with PD's controller behavior (CPDNN), and Fuzzy Logic to adjust PID controllers' gains (Fuzzy-PID). The proposed Fuzzy-PID and CPDNN achieves better performances compared to RMNN. Thus, comparing PSO results to those of Fuzzy-PID and CPDNN confirm that PSO is highly competitive to soft computing methods.      Table 4, the controllers tuned using PSO-CS have good performances for positions and attitudes. Thus, PSO-CS results are seen concurrent to those obtained using PSO, where responses lead to fast settling times with no overshoot and reduced errors. The desired sensor signal curves in Figures 7 and 8 (desired roll and pitch angles (φ d and θ d )) depend on the outputs of the controllers PID-x, PID-y, and PID-z. It is observed from this figure that these intelligent approaches, and especially PSO and PSO-CS, allow to follow perfectly the sensor signals for roll and pitch (φ d and θ d ), and give controlled responses close to the desired curves.

Considering the Figures 3-8 and
These programs are implemented using Matlab R2016a and performed on a computer with a processor Intel ® Core™ i7-3770 CPU of 3.40 GHz and 8.0 GB of RAM. From the conducted experiments, we can get conclusions that for the same goal, the PSO algorithm spends less CPU time (10 min 45 s) than the PSO-CS algorithm (14 min 12 s) and the CS algorithm (12 min 15 s). During the calculation of the velocity update equation, PSO-CS consumes an additional CPU time as the algorithm combines the two matrices generating Lévy flights random walks of cuckoos and that generating the movement of particles toward the global best solution Pg. However, the PSO-CS and PSO algorithms get higher performances (reduced errors) than the CS algorithm and RM method. These results are achieved thanks to the good adjustment of PSO, CS and PSO-CS parameters: Position and velocity initializations, the appropriate choice of the fitness function, the good tuning of the weighting coefficients w, C 1 , and C 2 , the chosen probability P a , and the Lévy flights settings for α and λ.

Rectangular Trajectory Tracking
This experimental simulation on the quadrotor engages to test the PIDs' efficiency when the vehicle tracks a square trajectory without overturning. In this mode, inputs for desired positions (x d , y d , z d ) are considered pulse signals and desired yaw attitude is set to 0 • . Figure 9 shows the tracking of each position coordinate (x, y, z) and the tracking trajectories in 3D plan obtained using the proposed methods. Responses for desired and controlled roll (φ) and pitch (θ) attitudes are illustrated in Figure 6. Table 5 resumes the errors ISE used to compute fitness function.
It can be noticeably observed from Figure 9 and Table 5 that the adaptive autotuned controllers using PSO, CS, and PSO-CS have good performances (reduced errors) for positions (x, y, z) control. PSO-CS and PSO responses always lead to fast settling times with no overshoot for x position. RM presents overshoots for positions and CS shows few overshoots for y position and altitude z. The intelligent PSO and the cooperative PSO-CS drive the system toward perfect tracking of the prescribed square trajectory. Also, all these methods allowed maintaining the yaw angle in zero degrees: no overturning of the aircraft. Figures 10 and 11 confirm that these intelligent approaches of PSO-CS nd PSO give controlled responses closed to the desired roll and pitch attitudes. and λ.

Rectangular Trajectory Tracking
This experimental simulation on the quadrotor engages to test the PIDs' efficiency when the vehicle tracks a square trajectory without overturning. In this mode, inputs for desired positions (xd, yd, zd) are considered pulse signals and desired yaw attitude is set to 0°. Figure 9 shows the tracking of each position coordinate (x, y, z) and the tracking trajectories in 3D plan obtained using the proposed methods. Responses for desired and controlled roll (ϕ) and pitch (θ) attitudes are illustrated in Figure 6. Table 5 resumes the errors ISE used to compute fitness function.
It can be noticeably observed from Figure 9 and Table 5 that the adaptive autotuned controllers using PSO, CS, and PSO-CS have good performances (reduced errors) for positions (x, y, z) control. PSO-CS and PSO responses always lead to fast settling times with no overshoot for x position. RM presents overshoots for positions and CS shows few overshoots for y position and altitude z. The intelligent PSO and the cooperative PSO-CS drive the system toward perfect tracking of the prescribed square trajectory. Also, all these methods allowed maintaining the yaw angle in zero degrees: no overturning of the aircraft. Figures 10 and 11 confirm that these intelligent approaches of PSO-CS nd PSO give controlled responses closed to the desired roll and pitch attitudes.       To test the robustness of these PID controllers, external disturbance was introduced. The most likely disturbance acting on the quadrotor is the wind disturbance along z axis. The external disturbance is modeled in this simulation by a random signal of variable amplitude inserted in the closed loop scheme of altitude z. It is noticed from Figures 12-14 and the errors ISE presented in Table 6 that when random disturbance affects the quadrotor, PSO, CS, and PSO-CS controllers reject these undesirable effects and provide outputs responses close to the desired ones. Therefore, we can affirm that PSO, CS, and especially the cooperative PSO-CS, are able to tune optimally the six PID controllers and drive the quadrotor to the desired altitude z, positions x and y, and yaw angle.  To test the robustness of these PID controllers, external disturbance was introduced. The most likely disturbance acting on the quadrotor is the wind disturbance along z axis. The external disturbance is modeled in this simulation by a random signal of variable amplitude inserted in the closed loop scheme of altitude z. It is noticed from Figures 12-14 and the errors ISE presented in Table  6 that when random disturbance affects the quadrotor, PSO, CS, and PSO-CS controllers reject these undesirable effects and provide outputs responses close to the desired ones. Therefore, we can affirm that PSO, CS, and especially the cooperative PSO-CS, are able to tune optimally the six PID controllers and drive the quadrotor to the desired altitude z, positions x and y, and yaw angle.

Helical Trajectory Tracking
Through the last simulation, we will investigate more the efficiency of the proposed intelligent control methods (PSO, CS, and PSO-CS) to search for optimal PID autotuned controllers for the quadrotor, when the desired trajectory is a helix. Simulation results of each position's tracking (x, y, z), and the tracking in 3D plan are presented in Figure 15. The tracking of roll and pitch attitudes (φ and θ) are presented in Figures 16 and 17. Table 7 summarizes the obtained ISE errors. Figure 15 clearly shows the possibility of the quadrotor system to follow a linear movement in altitude z in the first mode from 0 to 18 s and to keep perfectly a stable attitude z of 20 m in the second mode from 18 to 30 s. For translation motion, the quadrotor makes a complete circle of 2 m of radius in 30 s, where we notice that the errors in x, y, and z axis are very small (from Table 7). The same figure illustrates in 3D plan that the quadrotor tracks the helix trajectory and ends by a circle form. Also, the quadrotor avoids with high precision the overturning and maintains exactly 0 • in yaw attitude.

Helical Trajectory Tracking
Through the last simulation, we will investigate more the efficiency of the proposed intelligent control methods (PSO, CS, and PSO-CS) to search for optimal PID autotuned controllers for the quadrotor, when the desired trajectory is a helix. Simulation results of each position's tracking (x, y, z), and the tracking in 3D plan are presented in Figure 15. The tracking of roll and pitch attitudes (ϕ and θ) are presented in Figures 16, 17. Table 7 summarizes the obtained ISE errors. Figure 15 clearly shows the possibility of the quadrotor system to follow a linear movement in altitude z in the first mode from 0 to 18 s and to keep perfectly a stable attitude z of 20 m in the second mode from 18 to 30 s. For translation motion, the quadrotor makes a complete circle of 2 m of radius in 30 s, where we notice that the errors in x, y, and z axis are very small (from Table 7). The same figure illustrates in 3D plan that the quadrotor tracks the helix trajectory and ends by a circle form. Also, the quadrotor avoids with high precision the overturning and maintains exactly 0° in yaw attitude.        Table 7 illustrate that PSO, CS, and PSO-CS help to follow perfectly the sensor signals (φ d and θ d ). These results prove the possibility of our system to follow the imposed combined trajectory form with precision, for all variables (positions (x and y), altitude z, and attitudes (φ, θ, and ψ).
These controllers are tested in order to evaluate their performances under the same disturbance (a random input of variable amplitude) inserted in the altitude z closed loop scheme of the quadrotor. It can be observed from     Table 7 illustrate that PSO, CS, and PSO-CS help to follow perfectly the sensor signals (ϕd and θd). These results prove the possibility of our system to follow the imposed combined trajectory form with precision, for all variables (positions (x and y), altitude z, and attitudes (ϕ, θ, and ψ).
These controllers are tested in order to evaluate their performances under the same disturbance (a random input of variable amplitude) inserted in the altitude z closed loop scheme of the quadrotor. It can be observed from   Table 8 that PSO and PSO-CS controllers reject the disturbance effects and drive all variables to the desired ones. Therefore, we can affirm that despite the usage of a combined trajectory form and the presence of forcible disturbance, adaptive PSO-CS and PSO are proven as remarkably robust and efficient for quadrotor control.

Conclusions
In this paper, Particle Swarm Optimization (PSO), Cuckoo Search (CS), and the proposed cooperative PSO-CS were applied to control the quadrotor system. Results proved the efficiency of PSO approach and the cooperative PSO-CS to control optimally the quadrotor's outputs (φ, θ, ψ, x, y, z) compared to CS and the classical Reference Model (RM) methods. Experimental simulations also confirm the effectiveness of the aircraft using the proposed intelligent PID controllers to track the imposed rectangular and helical trajectories. Indeed, the robustness and performances of the intelligent controllers are validated in presence of disturbance. The performances of PSO-CS showed not only the stability robustness against the external disturbance, but also the fast response and the excellent tracking capacity for the different trajectories. These achievements are the results of the good settings of PSO, CS, and PSO-CS parameters. The advantages of the proposed PSO-CS in seeking for the best solutions in control problem are the results of exploiting the local search capacity of CS and benefiting from the global intelligence offered by PSO.
Concerning the limitations of our PSO-CS, one more parameter is added to CS program, which is the acceleration constant C 2 . A wrong setting of this parameter can cause premature convergence.

Conflicts of Interest:
The authors declare no conflict of interest.