Abstract

This paper explores the model parameters estimation of a quadrotor UAV by exploiting the cooperative particle swarm optimization-cuckoo search (PSO-CS). The PSO-CS regulates the convergence velocity benefiting from the capabilities of social thinking and local search in PSO and CS. To evaluate the efficiency of the proposed methods, it is regarded as important to apply these approaches for identifying the autonomous complex and nonlinear dynamics of the quadrotor. After defining the quadrotor dynamic modelling using Newton–Euler formalism, the quadrotor model’s parameters are extracted by using intelligent PSO, CS, PSO-CS, and the statistical least squares (LS) methods. Finally, simulation results prove that PSO and PSO-CS are more efficient in optimal tuning of parameters values for the quadrotor identification.

1. Introduction

Over the past few years, the unmanned aerial vehicles (UAVs) demand has increased dramatically because of the wide range of civilian and military applications. Some of those applications include low cost filming, panoramic picturing, area mapping, surveillance, air pollution monitoring, hostile zones intervention, transmission lines and power distribution inspections, earth science research assistance, etc [1].

In this paper, we are particularly interested in the behavior of one kind of UAVs that has two pairs of rotating rotors attached to the end of a cross named quadrotor or quadcopter. These aircraft are 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 two motors (the left one and the right one) rotate clockwise and the two others rotate counter clockwise, gyroscopic effects and aerodynamic torques tend to cancel in trimmed flight.

The quadrotor modelling is regarded as a delicate task either using Newton–Euler [2] or Euler–Lagrange [3] formalisms. The model obtained using these approaches is strongly nonlinear, fully coupled, underactuated, and dynamically unstable with complex behavior. So, it is of our interest to explore an efficient model parameter estimation technique to realize precise modelling results without using complicated model structures.

By using systems identification, the quadrotor can be represented by an estimated mathematical model, based only on input and output data considering the aircraft as a black-box process. Therefore, the model’s parameter values need to be estimated optimally. In the literature, some classical methods allow the process identification based on the step response such as the Strejc and Broida methods that require aperiodic systems. The integrating process allows identifying systems whose output response corresponds in steady state to a variation in a ramp. However, it is hard to model the unstable divergent behavior of the highly nonlinear quadrotor system using these classical techniques. Other methods raise the problem of identification in statistical terms of parameters estimation. The least squares (LS) method minimizes the squared error between the values predicted by the model and the observed values. This method demonstrated its superiority in parametric identification. In [4], it is shown that the problem of parametric identification of a wiener system could be reduced to a linear parametric estimation problem by a simple input-output data recorded using recursive least squares method (RLS).

The optimization methods of evolutionary algorithms (EAs) and swarm intelligence (SI) techniques can effectively resolve complex optimization problems compared to the described classical and statistical methods. EAs use the survival principle on a set of potential solutions to produce gradual approximations to the optimum, where SI is based on the study of group behavior in decentralized and self-organized societies. Particle swarm optimization (PSO) algorithm [5, 6] and genetic algorithm (GA) [79] seem to be the most successful types of EAs and SI, respectively. PSO deals with problems in which a best solution can be represented as a point or surface in a D-dimensional space. This intelligent method has shown superior performances. First, it can escape from local optimization problems. Second, it has no evolution operators such as mutation. Other advantages of PSO method are its less computational complexity compared to GA and its ease of implementation.

Cuckoo search (CS) algorithm is a novel SI algorithm motivated by the aggressive breeding of a bird called “cuckoo.” An advantage of CS compared to PSO and GA is that it uses less number of parameters to be tuned, which makes it more adaptable [10]. Also, the immigration and environmental specifications of cuckoos’ groups help to converge and reach best places for breeding and egg laying [10]. Many applications of PSO and CS are proposed in recent works of real discrete optimization problems [1114], in identification problems [1518], and for quadrotor’s control [1922]. However, the CS algorithm suffers from its low convergence speed, since it uses a fixed step size over generations. Our proposed cooperative PSO-CS algorithm in [22] combines the social and local search capabilities of PSO and CS. The PSO-CS offers great guidance for cuckoos to the global best positions and ensures a balance between exploitation and exploration of the search space [22].

In this paper, we propose the intelligent PSO and CS methods, the cooperative PSO-CS, and the statistical LS to optimally identify the quadrotor’s dynamics under determined operating conditions. The model used to estimate the variations of roll (ϕ), pitch (θ), yaw (ψ) angles, and altitude z during the flying process is composed of four subsystems of second order with the same structure and whose coefficients are adjusted with PSO, CS, PSO-CS, and LS to represent as best as possible the quadrotor divergent unstable behaviors. A comparative study is done to highlight the efficiency of the proposed intelligent methods in quadrotor identification.

This paper is organized as follows. In Section 2, Newton–Euler formalism is used to derive the motion equations of the quadrotor. In Section 3, the identification strategy is explained. In Sections 4 and 5, the quadrotor identification is detailed using the proposed intelligent PSO and CS techniques, the cooperative PSO-CS, and the statistical LS. In section 6, the simulation results are given. In the last section, our conclusions are provided.

2. Quadrotor Dynamic Modelling

The quadrotor is a complex flying machine, which is strongly nonlinear, fully coupled, and underactuated (6-DOF and only four control inputs). Therefore, its aerodynamic is affected by many physical effects, including gravity, gyroscopic effects, friction, and moment of inertia. As shown in Figure 1, the quadrotor consists of 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 two pairs of propellers must spin in opposite directions to prevent the vehicle’s overturning. In Figure 1, 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 (ϕ, θ, ψ).

For simplifying the delicate dynamic modelling of the quadrotor, various working hypotheses have been assumed [1]: the quadrotor structure is assumed rigid and symmetrical (diagonal matrix of inertia); the propellers are assumed rigid (negligible effect of deformation during rotation); a very close approximation of the aerodynamic behavior is assumed (the lift and the drag forces are proportional to the square of the rotational speed of the rotors); and the origin of reference related to the structure is fixed on the center of mass of the quadrotor.

By using Newton–Euler formalism, the equations can be written as follows:where is the linear velocity and m is the total mass of the quadrotor.

R is the rotation matrix:

Ω is the angular velocity expressed in the fixed reference:

S(Ω) denotes the oblique symmetric matrix:

J is the inertia of the system:

Ff is the total thrust force generated by the four rotors:where Fi is the force generated by the rotor i and means the vector product.

Ft is the drag force along the axes x, y, and z:where Kftx, Kfty, and Kftz are the translational drag coefficients.

is the gravity vector:

Mf is the moment caused by the thrust and the drag forces following the three rotations:

Ma is the moment resulting from aerodynamic friction:where Kfax, Kfay, and Kfaz are coefficients of aerodynamic friction.

Mgh is the gyroscopic moment of the propellers:where Jr is the inertia of the rotors and Ωr=ω1 − ω2+ω3 − ω4;

Mgm is the gyroscopic moment due to the movements of the quadrotor:

Then, the dynamic model of the quadrotor can be expressed by the following equations, as given in [1]:with

Table 1 shows the definitions of the quadrotor’s parameters.

3. Identification Strategy

Simulations of the quadrotor dynamical model, expressed by equations (13)–(18), with a unit step reference signal for U1, U2, U3, and U4, give unstable responses. The instability observed for the four flight variables is of divergent type, which can be simply modelled by second-order systems. In addition, the flight parameters of the quadrotor can be separated, thus giving four subsystems of order 2 having the same form Gmi (p), with i = {ϕ, θ, ψ, z}.

When the model of system is fixed, the identification task can be treated as an optimization problem. The basic idea of parameter estimation is to compare the time-dependent responses of the system and the model based only on inputs and outputs data [5]. Considering Figure 2, Uk are the excitation inputs and Ek are the errors of identification that characterize the difference in behavior between the quadrotor system and model, with k = {1, 2, 3, 4}. Then, PSO, CS, PSO-CS, and LS are used as adequate methods to optimally estimate the quadrotor’s model parameters (Kmi, ami, and bmi).

4. Quadrotor Identification Using Intelligent Methods

4.1. Intelligent PSO

For hard optimization problems, particle swarm optimization (PSO) was developed by Eberhart and Kennedy in 1995. The basic principle of PSO was inspired by the social behavior of animals moving in a swarm as bird flocking. To search for food, each bird flies in the space of solutions and determines its speed according to its personal experience and the information obtained through interaction with other swarm members [4].

The initialization matrix contains N particles dispersed in a search space along dimension j for j = {1, 2, …, D}. Each particle Pi (i= 1, 2, …, N) stores its best position Pbi (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 (22). The mechanism of displacement of each particle is managed by three rules. Firstly, the particle tends to follow the direction of its current velocity. Secondly, it tends to move towards its best position. Finally, it tends to move to the best position reached by its neighbors [5, 23]. In fact, the new velocity matrix Vij and position matrix Xij of particles are calculated at iteration (t+ 1), according to equations (23) and (24):where Pbij is the best position found by the particle i; Pgij is the best position found by the neighborhood; , C1, and C2 are weighting coefficients; and R1 and R2 are random variables generated from a uniform distribution in [0, 1].

4.2. Intelligent CS

Cuckoo search (CS) algorithm, proposed by Yang and Deb in 2009, is based on the life of a “cuckoo” bird. The basic principle of CS 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 are not discovered and removed by host birds. Reproduction and breeding are favoured by cuckoo groups’ immigration, converging and reaching the best places [10].

The primary population of CS contains N nests, and each nest is composed of D eggs. The best nests with a high quality of eggs (solutions) carry over to the next generations, where the quality evaluation is based on the fitness function F of the habitat (array of 1 × D). The host can discover an alien egg with probability Pa from [0, 1], which is approximated by a fraction Pa of the N nests being replaced by new nests, having new random solutions [24].

To explore the search space when replacing solutions in the nests with new solutions, Lévy flight mechanism is used. The step length S from Mantegna algorithm (based on Gaussian normal distribution denoted by Norm) can be written as represented by equation (25), and σu2 is the variance of the distributions given by equation (26). Therefore, a new solution Xi (t+ 1) for cuckoo i is given by equation (27), and the fraction Pa of worse solutions is generated as given by equation (28):where Γ is the gamma function; α > 0 is the step size; R and r are random variable generated from a uniform distribution in the interval [0, 1]; Xj (t) and Xk (t) are two random solutions chosen by random permutation; and H is a Heaviside function.

4.3. Cooperative PSO-CS

The initialization matrix of the cooperative PSO-CS is of dimension D × N, and the solutions’ quality is evaluated as in PSO and CS. The global best particle (or the best nest) is the particle (or the nest) that has the smallest fitness value among all potential solutions [22].

To overcome the fast convergence speed of PSO and the low convergence speed of CS, the PSO-CS combines the capacities of social thinking in PSO and local search in CS. Thus, the displacement equation is modified, by combining the Lévy flight random walks of cuckoos and the velocity of particles toward the global best solution Pg. New solutions Xi (t+ 1) are given by [22]

Both algorithms’ capabilities are combined to increase the particles’ diversification. PSO-CS guides the cuckoos toward the global best positions (global intelligence of the swarm). In fact, the search ability increases during iterations, and the exploration of the local and global places is achieved by the Lévy flight displacement of cuckoos and the velocity of particles toward the global best solution (Pg), as given by equation (29). The flow chart of PSO-CS is considered the same as that of CS. It is given in Figure 3. The only modification is in the expression of the velocity of displacement, which helps to search at local and global scales in order to move all cuckoos toward best environment and to quickly converge at later stage.

4.4. Parameter Setting for PSO, CS, and PSO-CS

Intelligent PSO, CS, and PSO-CS are applied to optimally select the model parameters from N = 200 solutions. So, the search space interval is chosen sufficient to contain all possible solutions [0, 200], and its dimension D is set to 12. The compromise between local and global exploration in PSO is achieved for  = 0.8. C1 takes random values in [0, 0.8] to avoid the problem of fast convergence, while C2 takes random values in [0, 1.2] to give more importance to the global best solution Pg [22]. While in CS, the parameters used in experiments are as follows: abandon probability, Pa = 0.25, and the Lévy flight settings, α = 0.1 and λ = 1.5 [22]. In PSO-CS, the same settings of PSO and CS are kept [22]: Pa = 0.25 and the Lévy flight settings α = 0.1, λ = 1.5, and C2 with random values in [0, 1.2].

These algorithms are evaluated using a profit defined in a similar way in order to minimize the differences between the output responses of the estimated model and the quadrotor system. This fitness function is defined in equation (30) as the sum of the quadratic errors Ek previously mentioned, with l = 4 and k = {1, 2, 3, 4}:

An appropriate set of PSO, CS, and PSO-CS parameters can yield model responses close to those of the quadrotor. The maximum number of generations for the three programs (PSO, CS, and PSO-CS) is fixed as the stop criterion and set to 20.

5. Quadrotor Identification Using the Statistical LS Method

The method of least squares (LS) provides the parameters of a model so that the sum of squared errors (between the predicted and observed values) is minimal [25].

For the quadrotor’s identification, we consider the same identification scheme represented in Figure 2. The vectors of measurements are extracted from the temporal responses of the quadrotor’s attitude (ϕ, θ, ψ) and altitude z.

The recurring equation in equation (31) is obtained by discretizing at Te the temporal response of each transfer function Gmi (p) (in equation (20)), where uk and yk are, respectively, discrete inputs and outputs of Gmi (p):

The aj and bk coefficients or θ (equation (32)) to be estimated, in the sense of the LS criterion given in equation (33), is computed following equation (34). The matrix F and the vector Y are constituted from the values taken by u and y at the different sampling times:

The transition between continuous (km, am, bm) and discrete parameters (a0, a1, b0b1) is grouped in the following equations:

The application of this approach to the identification of the quadrotor’s movements is discussed in the following section.

6. Simulation Results

The quadrotor system is identified using a model composed of four subsystems with the same structure, expressed by Gmi (p) in equation (20). The quadrotor model is generated in a simulation time of 500 seconds. Table 2 summarizes the parameters of these subsystems (Kmi, ami, bmi) estimated by using PSO, CS, PSO-CS, and LS methods, with i= {ϕ, θ, ψ, z}. Figures 47 show the responses of the quadrotor’s system and models with parameters obtained using PSO, CS, PSO-CS, and LS, in an observation time of 500 seconds. Table 3 shows the results of statistical analysis based on the integral of a positive term involving the error: the integral absolute error (IAE), the integral square error (ISE), the integral time absolute error (ITAE), and the integral time square error (ITSE) [6].

PSO and PSO-CS perform the search with high precision to avoid all solutions that are far from the desired responses of the quadrotor system (angles of roll, pitch, yaw, and altitude z). From Figures 47, the model responses with parameters adjusted using PSO and PSO-CS are very close to those of the quadrotor’s responses, while CS and LS give the estimated model parameters that allow fitting a divergent behavior near to that of the quadrotor system.

These observations spread to 500 seconds prove the validity of the chosen model to identify the divergent responses of the quadrotor, as well as the efficiency of PSO and PSO-CS for seeking for the best solution, even if identifying unstable responses is a very difficult task. During the simulation time, PSO and PSO-CS effectively avoid the divergence that can appear between the model (20) with the calculated parameters and the quadrotor system defined by equations (13)–(18). PSO and PSO-CS reduce the fitness function and result in an optimal search for the parameters of the models. To confirm the effectiveness of PSO and PSO-CS in reducing the errors between the two responses for the four outputs, we present in the following table a statistical analysis based on the integral of the errors IAE, ISE, ITAE, and ITSE.

Statistical results in Table 3 show that the errors between the model outputs and the system outputs, such as IAE, ISE, ITAE, and ITSE are all smaller (higher performance) when the model is established by PSO-CS than PSO, CS, and LS.

In [26], the authors proposed a modified CS algorithm named oriented cuckoo search algorithm (OCS), where they have tested its performance for different probability distributions (DDICS1, DDICS2, MCS1, MCS2, DACS1, DACS2, and OCS-LG). The proposed oriented cuckoo search algorithm with Lévy distribution and standard Gaussian distribution (OCS-LG) has shown better performance in means of mean error. The lower value means the better performance; therefore, OCS-LG is superior to other six algorithms. From the results of Wilcoxon test, there are significant differences (p value below 0.05) among OCS-LG, DACS1, MCS1, and MCS2. In other words, OCS-LG is significantly better than DACS1, MCS1, and MCS2.

The proposed hybrid many-objective cuckoo search (HMaOCS) for many-objective optimization problems (MaOPs) in [27] show that HMaOCS is promising in dealing with most many-objective optimization problems. It was noted that in CS, there is no guidance information (the global best individual) and each individual in the population is not affected by any other individuals. Then, the proposed changes were focused in the population updating method, where the new population was generated from the combined parent and offspring populations to ensure better convergence and diversity. However, the computation time necessary to generate new solutions in this way was not considered.

Comparisons to OCS-LG and HMaOPs show that our proposed PSO-CS algorithm used Lévy flight mechanism with Gaussian normal distribution when replacing solutions in the nests with new solutions, which is the best solution attested in [26]. In addition, PSO-CS benefits from the global intelligence already calculated in each iteration in order to manage the displacements of solutions.

Our proposed programs are performed on a computer with a processor Intel® Core™ i7-3770 CPU of 3.40 GHz and 8.0 GB of RAM using Matlab R2016a. The implemented PSO-CS program consumes an additional CPU time (6 min 15 s) than PSO and CS programs that spend (5 min 07 s) and (5 min 50 s), respectively. This additional simulation time of PSO-CS is produced due to the sum of both matrices generating Lévy flight random walks of cuckoos and the movement of particles toward the global best solution . This updating method is faster than generating solutions from the combined parent using selection and crossover operators.

7. Conclusions

In this paper, recent heuristics optimization of particle swarm optimization (PSO) and cuckoo search (CS), the proposed cooperative particle swarm optimization-cuckoo search (PSO-CS), and the statistical least squares (LS) has been applied to identify the variations of the rotational movements (ϕ, θ, and ψ angles) and the translational movement along z-axis of a quadrotor UAV. Simulation results proved the efficiency of the PSO approach and the cooperative PSO-CS in identifying optimally the quadrotor’s outputs (ϕ, θ, ψ, z) compared to CS and LS methods. The proposed PSO-CS seeks for the best solution by exploiting the local search capacity of CS and benefiting from the global intelligence offered by PSO.

These achievements are the results of the good choice of the model structure to identify the quadrotor system responses with unstable behaviors and the good adjustment of PSO, CS, and PSO-CS parameters (the fitness function, the weighting coefficients , C1, and C2, the probability Pa, and the Lévy flight settings for α and λ). However, PSO-CS program consumes an additional CPU time to calculate the velocity update equation, and a wrong setting of the acceleration constant C2 can cause premature convergence.

Data Availability

No data were used to support this study.

Conflicts of Interest

The authors declare that they have no conflicts of interest.