A Cauchy mutant pigeon-inspired optimization–based multi-unmanned aerial vehicle path planning method

To improve the performance of multi-unmanned aerial vehicle path planning in plateau narrow area, a control strategy based on Cauchy mutant pigeon-inspired optimization algorithm is proposed in this article. The Cauchy mutation operator is chosen to improve the pigeon-inspired optimization algorithm by comparing and analyzing the changing trend of fitness function of the local optimum position and the global optimum position when dealing with unmanned aerial vehicle path planning problems. The plateau topography model and plateau wind field model are established. Furthermore, a variety of control constrains of unmanned aerial vehicles are summarized and modeled. By combining with relative positions and total flight duration, a cooperative path planning strategy for unmanned aerial vehicle group is put forward. Finally, the simulation results show that the proposed Cauchy mutant pigeon-inspired optimization method gives better robustness and cooperative path planning strategy which are effective and advanced as compared with traditional pigeon-inspired optimization algorithm.


Introduction
Path planning is the first step in a moving multi-agent mission system.To distinguish from path management, path planning aims to find a set of waypoints from a start position to a target position to instruct agents moving. 1 While path management aims to generate feasible trajectories based on waypoints which are parameterized by time and could be described by twicedifferentiable polynomials.Path management is necessary for most of the intelligent agents, for example, unmanned aerial vehicles (UAVs), which are unable to move along zigzags paths which are simply generated by connecting the waypoints.
Two decades ago, path planning problem focus on the ground moving robot mission with known threat locations.3][4] Unfortunately in some complex environments, for example, when unpredictable obstacles exist, these randomly sampling search algorithms usually cannot construct feasible paths for robots.As a result, this contradiction promotes the motivation to work out more practical path planning strategy.][7][8][9] In recent decade, three-dimensional (3D) UAV cooperative mission path planning problem attracts prime attention since the traditional heuristic-based path planning methods suffer more and more challenges that resulted from realistic mission requirements, physical constrains of devices and complex environmental threats.The traditional heuristic path planning strategies are unable to give an optimal or sub-optimal solution rapidly and accurately in the new difficult mission settings, they have slow search speed and easily fall into local areas.To find paths timely and effectively with minimum fuel consumption and collision free in 3D cluttered complex environment, artificial intelligence (AI)-based path planning method attracts wide publicity.This kind of algorithms are inspired from nature including simulated annealing (SA), genetic strategy (GS) algorithm, ant colony (AC) algorithm, artificial neural network (ANN) algorithm, pigeoninspired optimization (PIO) algorithm, bat strategy algorithm (BA) and artificial bee colony (ABC) algorithm.1][12][13][14][15][16][17][18] To improve the optimality of solutions obtained by AIbased path planning algorithm, many researchers dedicate their efforts.Wu et al.'s 19 metropolis criterion which is based on probability theory is utilized to obtain new individuals for traditional particle swarm optimization (PSO) algorithm.In an advanced path planning strategy, chromosome-encoding strategy is proposed for conventional genetic algorithm to guarantee a feasible path to avoid the searching circulation.While a call back mechanism is designed in for traditional bee colony algorithm to give a shortest and safest path. 20,21n this article, multi-UAV co-operative path planning problem in plateau narrow area is defined.By comparing with other 3D sceneries, mountainous plateau environment provides more challenges for UAV path planning, for example, dangerous peaks and unpredictable meteorological factors like atmospheric turbulence, wind shear and so forth. 22These disadvantages will seriously affect the stability of UAV cooperative flight.When implementing path planning in this extremely complex workspace, an AI path planning algorithm with fast convergence speed is urgently required.A new kind of swarm intelligence optimization algorithm based on pigeon homing behavior algorithm is proposed in Duan and Qiao, 12 called PIO algorithm.The most remarkable advantage of PIO algorithm is its fast convergence speed.
The algorithm consists of two independent computing sections: geomagnetic navigation section and landmark navigation section, which simulates the navigation mechanism of pigeons using different navigation tools at different stages in the process of finding their destination.Specifically, pigeons first adjust their direction through geomagnetic field and solar altitude, and then they will judge the direction and location by the landmark information adjacent to the destination.Inspired by the process by which pigeons find their nest and the similarity of the process with UAV path planning process, that is, both have multi-objective optimization problems with known starting and ending points, it is feasible to apply the pigeon's behavior model for the path planning of UAVs.
However, due to the high complexity of plateau narrow area and its internal attribution, PIO algorithm is easy to fall into local optimum in the early stage of operation, resulting in unsatisfactory results.The best way to avoid local optimum is to increase the diversity of population.In this article, mutation operator is employed for map and compass operator.Normally, mutation operator is used in genetic operation to ensure population diversity, which is similar to the process of biological variation.Mutation operation could create a new individual by changing some information of an original individual.In general, mutation operators 23 include bit element mutation, boundary mutation, uniform mutation, Cauchy mutation and so on.
By comparing and analyzing the changing trend of fitness function of the local optimum position and the global optimum position when solving UAV path planning problems, Cauchy mutation operator is chosen to improve the PIO algorithm, the new algorithm is called Cauchy mutation PIO (CM-PIO) algorithm.This improved PIO algorithm can effectively extend the search area by mutation operation and could ensure that the pigeons can fly to the global optimal solution quickly.It can also reduce the risk that the solution will fall into the local optimal.Cauchy mutation operator aims to change the individual's state value by selecting random numbers obeying Cauchy distribution, and to realize individual variation.In addition, considering the length and location information of each UAV path, a cooperative path planning strategy for multiple UAV is designed.
The main contributions in this article are as follows: 1.The integrated model of multi The organization of this article are as follows: section ''Problem formulation'' defines the complete problem formulation constraints, that is, UAV parameters, terrain constraints, wind field model, the cost of path and multi-objective optimization model which is followed by section ''PIO Algorithm,'' which discusses the PIO algorithm.Section ''CM-PIO-based path planning'' defines the CM-PIO-based path planning of UAVs.In section ''Cooperative path planning strategy,'' the path planning strategy is discussed.Section ''Simulation results and discussion,'' defines the simulation results and its discussion.Finally, section ''Conclusion'' concludes the whole article.

UAV's parameters
Each UAVs have their own parameters, taken from Hebecker et al. 24 and Jian et al., 25 such as maximum range limit, the shortest route limit, maximum turning angle limit, maximum climb/dive angle limit and maximum flight altitude limit.
Maximal flight length.The maximal flight length of UAV is limited by the amount of fuel.Besides, when UAV flies through a wind field, the drag becomes larger and causes an increment in the consumption of fuel.However, when passing through the mountain areas, the UAVs have the action of climbing or descending.In this case, the energy consumption of the UAV is higher.So, the maximal flight length will be limited.Assuming that the maximal flight length is L max and a complete path is divided into n sections.l i , represents the length of section i, the flight length should satisfy Minimal length of leg.When UAV flies in a leg, it needs some time to adjust their altitude to reach the next leg.
To ensure flight safety and quality, the length of leg must be more than the smallest value.The length of every legs should satisfy Maximal turning angle.The maximal angle between the directions in current leg and next leg are called the maximal turning angle.To avoid mountains, sometimes UAVs have to change their flight direction.But, the maximal turning angle of UAV is limited by its maneuverability.The real turning angle must be smaller than their maximal turning angle.Assuming a i and a i + 1 represent the projection of the sections i and i + 1, respectivley.The legs on horizontal plane, u max , represents the maximal turning angle.a i and a i + 1 should satisfy Maximal climb/dive angle.Maximal climb/dive angle refers to the maximal angle to which UAV is able to climb or dive relative to horizontal flight, and the maximal climb/dive angle is expressed as q max .Assuming p i and p i + 1 represent the ith waypoint and (i + 1)th waypoint with the coordinates (x i , y i , z i ) and (x i + 1 , y i + 1 , z i + 1 ), respectively.While, p 0 i and p 0 i + 1 represents the projection of p i and p i + 1 on horizontal plane, the following relationship should be satisfied as Maximal flight altitude.When flying at low temperature environment, which is one of the main features of plateau, performance engine always suffers degradation.
The maximal flight altitude of UAV is lower than that in plain areas.Assuming the maximal altitude is H max , the real-time altitude H i cannot exceed its maximal altitude, that is Terrain constraint A function commonly used in terrain simulation is used to generate datum terrain and peak terrain in Chao. 26he expression of datum terrain function is where a, b, c, d, e, f, g are all constants, x and y are the horizontal coordinates and Z 1 is the altitude at a certain horizontal position.The values of a, b, c, d, e, f, g are relative to the degree of height fluctuation.The expression of peak terrain function is where (x 0i , y 0i ) represents the position of the center of ith peak projecting on the ground plane.h i is the altitude of the peak i. (x si , y si ) represents the decline scale along the x direction and y direction of the peak i.
While, Z 0 represents the base altitude.N represents the number of peaks.
Combining the datum terrain function and peak terrain function, the plateau terrain will be described as which is depicted in Figure 1.

Wind field model
When UAVs fly in plateau mountainous area, their speeds and altitudes will be disturbed by wind.As a result, it is difficult to maintain a fixed formation of flight and even brings disadvantage to the flight safety.
In this research, the wind field is taken as a threat area and the spherical model is used to represent the wind field which is written as where (x k0 , y k0 , z k0 ) is the center of k the wind field, r k represents the radius of its spherical model.In the cooperative path planning process, the wind field should be avoided as far as possible.

The cost of path
The cost of path is an important index to evaluate the quality of a planned path.The purpose of path planning is to find a suitable path with the lowest cost of path.Based on the actual flight environment and the UAV's performance, the cost of path is classified into three sectors: distance cost, altitude changing cost and threat cost.
Distance cost.The limited amount of fuel must be taken into account in path planning.So, the flight distance should be shortened as far as possible.Assuming that a complete path is composed of n legs, the distance cost is expressed as where l i represents the length of i, the leg of the path.
Altitude changing cost.When UAVs fly near a peak, some height changing actions have to be implemented to avoid the peaks.But this maneuver will cause fuel consumption.In addition, low temperature makes the engine perform worst, frequent maneuver may bring risk for the flight safety.Therefore, frequent altitude change should be avoided as far as possible.
The altitude changing cost can be expressed as Threat cost.The threat cost is related to the distance between the UAV and the threat area.The closer the UAV is to the threat area, the bigger the threat cost is.
The main threat areas considered in this article are terrain threat and wind field threat.The threat cost is defined as where avgR = where K is a constant whose value is relative to the real flight space, e is penalty coefficient, which is generally a big number, n is the number of legs, r i is the distance between i, the leg of path and the center of the threat area.

Multi-objective optimization model
Considering the parameters of UAV and the cost of path, a multi-objective optimization model for UAV path planning can be obtained using the fitness function that can be expressed as and the constraint function is written as For the multi-objective function, it is difficult to find a solution that satisfies the minimal value of each evaluation function at the same time.It is necessary to synthesize all the cost values by appropriate weights to coordinate every performance.The fitness function of UAV path planning is constructed as The weights in the formula can be adjusted according to the actual mission requirements and flight environment requirements.When the desired mission time is short, the distance cost weight will be increased.While if the desired flight path of UAV is smooth, the altitude changing cost will be increased.

PIO algorithm
PIO algorithm is divided into two independent stages.In the first stage, the pigeon updates the position and speed of the individual by map and compass operator principle.In the second stage, the pigeon updates the state by landmark operator principle.The description of PIO algorithm is as follows: Assume that the search area of PIO algorithm is a ''D'' dimensional space, and there are ''M'' pigeons in the group, in which the position of the ith pigeon represents one of the solutions of the optimization problem.Each individual pigeon will constantly adjust its position according to certain rules.With the help of fitness function, the fitness of each individual can be calculated.The individual position will be continuously updated in each iteration, and finally the optimal solution will be found.The position of the ith pigeon in the ''D'' dimension space can be expressed as X i = (X i1 , X i2 , :::, X iD ), velocity can be expressed as V i = (V i1 , V i2 , ::::, V iD ), where i = 1, 2, :::, N.
In the first stage, each individual states are updated according to the information of the best individual in the population by updating position according to formula (16) and updating the speed according to formula (17) where N c represents the current iterations, R represents the map and compass factor.rand is a random value from 0 to 1.After N c À 1 iteration, the optimal position X gbest in the whole population is obtained by comparing the position of all the (pigeons) or UAVs.The maximal iteration in the first stage is set as N c1max .After N c1max iterations obtained by map and compass factor, the position and speed that update the rules of pigeon group will change, and the whole process enters the second stage.
In the second stage, when the pigeons approach target position, their states will be updated according to the landmark information.Some individuals will follow other individuals who are near the landmark, while some far from the landmark that will be abandoned because they have no ability to find the optimal path.In the calculation process, the number of pigeons will be halved in each iteration by abandoning pigeons that have poor adaptability and then look for the central location X center in the remaining population, regarding it as the landmark for flight.The state changes according to the following formulas where F is related to the fitness function.It could be expressed as The maximal iteration in the second stage is N c2 max .After the N c2 max iteration of landmark operator, the optimization process is finished and the optimal solution is obtained as its output.

CM-PIO-based path planning
Cauchy distribution is a continuous probability distribution without variance and mathematical expectation.If the random variable x satisfies the probability density function, it is called Cauchy distribution.
The probability density function of Cauchy distribution is In equation ( 22), x 0 is the peak value of the distribution and g is the width corresponding to half of the maximum value.When g = 1 and x 0 = 0, the random variable satisfies the standard Cauchy distribution, it is recorded as X;C(1, 0).The corresponding cumulative distribution function is an incremental function as

Cauchy mutation map and compass operators
In the traditional PIO, map and compass operators are mainly used to find the global optimal individuals in the search space, and then positions and speeds of the pigeon swarm will be changed by referring to the optimal individuals.While choosing the Cauchy mutation c 1 as the maps and compass operators can not only expand the search area, but also reduce the risk of falling into local optimum.The Cauchy weight coefficient c 1 obeys the Cauchy distribution where rand is a random value from 0 to 1, c 1 could be expressed as In each iteration, the position of the individual pigeon is updated by the following rule where X i0 is the calculated position of i, the pigeon.X gbest is the optimal individual with the optimal fitness after N c À 1 iterations.While X 0 i is the position of i, the pigeon after update.The position of the i pigeon at next iteration will be As for Cauchy mutation map and compass operators c 1 , when it is positive, the updated position X 0 i will be farther from the global optimal position X gbest and vice versa.By adding Cauchy mutation, half of the individuals in the population will diffuse outward to find a better location, while the other half will be not.Comparing the undated positions with the original position, the individual with better fitness can be retained, which can not only ensure the quality of optimization but also improve the diversity of the population.

Cauchy mutation landmark operators
In the tradition PIO algorithm, the scale of population will be reduced by half at each iteration in the landmark operator stage, with pigeon moving toward the central position of the pigeon swarm.Unfortunately, this rapid reduction of population will lead to premature convergence of the algorithm, which will have a negative impact on the optimization of the landmark operator stage.To avoid the premature convergence or missing the optimal solution, the traditional landmark operator is replaced by a Cauchy function, updating the individual position relying on the optimal position in the population.The Cauchy function is described as Cauchy weight coefficient c 2 obeys Cauchy distribution, which is where c 2 is a positive constant.The rule for updating individual position could be expressed as where X N c À1 i0 is the position of ith individual in N c À 1 iteration, X gbest is the position of global optimal individual.
In the Cauchy mutant landmark operator stage, all individuals gradually approach the global optimal solution because of the positive constant c 2 .A suitable Cauchy mutant can make the pigeon swarm effectively move with the appropriate speed and direction and guarantee the stable convergence of the algorithm.

Cooperative path planning strategy
Assuming three UAVs take off from different locations at a same time and fly to one specific area for a mission.During the flight process, the three UAVs are required to keep safety distance and reach the target area within a similar duration.
Figure 2 depicts a scenario for multiple UAV cooperative path planning.Basically, each UAV must avoid flying in threatening area.The figure shows the shortest safe distance D between each two UAVs.So, the distance between any two planned paths must be less than this safe distance D. In addition, UAV1, UAV2 and UAV3 are required to meet the time constraints 27 by adjusting the lengths of paths.
Therefore, in addition to care about the parameters of UAV, multi-UAV path planning also care about the space and time relationships between all UAVs which may give some constraints to cooperative path planning.
First, a sufficiently safe distance d safe between two adjacent UAVs must be kept during the flight.Assuming that the position of UAV1 is x 1 (t) at the moment t and that of UAV2 is x 2 (t).The flight distance between the two UAVs should satisfy the following relationship When all UAVs generally approach the target area, and their distances get closer and closer, so the safe distance should be set different value in different stages as where T is 80% of the total flight time, D and d are set according to the specific planning area.Moreover, both UAVs are desired to reach their desired target position simultaneously.According to the allowed speed range and length of path, the shortest time and the longest time for each UAV to reach the target point will be calculated.Assuming the allowed speed range of UAV i is V i = ½V i min , V i max and the length of path is L i , the flight duration will satisfy The temporal cooperation rule for two UAVs is In this case, the two UAVs will have the possibility to simultaneously arrival to the target position.
To meet the spatial requirement and the temporal requirement of UAV group, two coordination coefficients are added to each UAV's path evaluation function.They could be expressed as where f space is the spatial coordination coefficient, d is the distance between two UAVs, f time is the temporal coordination coefficient, N is a constant.T 1 and T 2 are the required time for UAV1 and UAV2 to arrive the targeted position.
In terms of the proposed strategy for cooperative path planning problem, each UAV is represented by a sub-population, and each sub-population evolves independently.Only when evaluating the individual fitness, the information sub-populations will be communicated.In the multi-population's co-evolution process, the optimal individuals of each population are selected to communicate with other populations.This co-evolution process is shown in Figure 3.
The path planning process for each UAV will refer to other UAV's paths and calculating the spatial and temporal coordination coefficients between UAVs.Finally, more suitable path will be selected which not only have advanced fitness but also meet the cooperation requirement.
The specific steps of CM-PIO cooperative path planning strategy are presented as following: Step 1. Define the area for path planning of every UAV.
Step 2. Initialize the algorithm parameters including the scale of population M, the dimensions of area D, map and compass factors R, maximal iteration N c1max and N c2 max in the first stage and second stage, respectively, the weight coefficients and so on.
Step 3. Individual coding rules randomly generates M individuals with positions and speeds which are within the allowed range of v min , v max ½ .
Step 4. Find m, the most optimal individuals in each population, by fitness function.And then combining with the information of other populations, find the ith individual with the biggest cooperation coefficient marked as x igbest in each population.
Step 5. Update the individual's positions and speeds in each population by the Cauchy mutation map and compass operator.
Step 6. Compare the current iteration in map and compass operator stage with the maximal iteration N c1max .If the current iteration is bigger than the maximal iteration N c1max , the operation will enter the landmark operator stage.Otherwise, return to step 4.
Step 7. Update the individual's positions by landmark operators.
Step 8. Compare the current iteration in landmark operator stage with the maximal iteration N c2 max .If the current iteration is bigger than the maximum iteration N c2 max , the operation finishes.Otherwise, return to step 7.
The diagram of multi-UAV coordination path planning strategy is shown in Figure 4.

Simulation results and discussion
The simulation results is carried out on MATLAB/ 2014a.Providing a path consisting of 11 waypoints and 10 path sections.Flight space is set as 20 km320 km35 km.There are three fixed-wing UAVs with same parameters involved in the experiment.

Simulation parameters
The parameters of UAV are presented in Table 1.
Based on the field investigation and topographic analysis of the mountainous area in the Qinghai-Tibet plateau, the parameters of the datum terrain are indicated as a = 10, b = 0:2, c = 0:1, d = 0:6, e = 1, f = 1 and g = 0:1.The parameters of the peak terrain are shown in Table 2.
The parameters of wind field model are shown in Table 3.
The parameters of proposed CM-PIO algorithm are presented in Table 4.
For the fitness function, weights of distance cost, altitude changing cost and threat cost are k 1 = 0:4,k 2 = 0:2 and k 3 = 0:4, respectively.In terms of minimal safe distance in cooperative path planning, setting D = 500 m, d = 100 m.The temporal cooperative coefficient is set as N = 300.
The coordinates of start and target points for each UAV are shown in Table 5.
Results of CM-PIO-based path planning method.To verify the robustness of the proposed path planning algorithm, the paths planned by traditional PIO algorithm are also presented in Figure 5.The cone represents the peak, the sphere areas near the peak represent the wind field, the solid line represents the optimal path planned by CM-PIO algorithm and the dash line represents the optimal path planned by PSO algorithm.Each path includes 11 waypoints which are represented by circles.
It can be seen from Figure 5 that both algorithms can plan feasible paths.Although, the route planned by PIO algorithm can avoid the threat of mountain peaks, the total distances are longer and the altitude near the peak and wind field suffer greatly changing, which will increase the risk of collision.The paths planned by CM-PIO is shorter, smoother, lower cost and safer.Figure 6 illustrates the change in the fitness of the    optimal paths of three UAVs, by PIO and CM-PIO.It can be seen that CM-PIO can find the optimal path faster and more excellent.Besides, the Cauchy mutant could create the most optimal individual when the traditional PIO has fallen in to the local optimum.
Results of cooperative path planning strategy.It can be seen in Figure 7 that without cooperative planning the path of UAV2 and the path of UAV3 are close and have intersection.To ensure the spatial requirement, cooperative path planning is needed.Figure 8 shows the fitness of optimal individuals changing with iteration planned by CM-PIO cooperative planning strategy and Figure 9 shows coordination coefficient f co changing with iteration.Table 6 shows the length of the planned path with and without coordination planning.
Combining the above data it can be seen that proposed CM-PIO cooperative path planning strategy could not only meet the coordination requirements but also finds optimal individuals with advanced fitness.

Conclusion
For the path planning problem of the UAV group in plateau area, first the problem was modeled as a multi-      objective optimization problem and then CM-PIO algorithm is designed to find the optimal path for each UAV.Compared with the traditional PIO algorithm, the result shows that the solution of proposed method has better robustness.Finally, considering the spatial constrain and temporal constrain of multi-UAV cooperation path planning problem, the CM-PIO cooperation path planning strategy was designed.The simulation results show the proposed method could not only find optimal individuals with advanced fitness but also meet the coordination requirements.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship and/or publication of this article.
path not through the threat area e, path through the threat area &ð12Þ

Figure 2 .
Figure 2. A background setting for cooperative path planning.

Figure 3 .
Figure 3. Co-evolution process in coordination path planning.

Figure 4 .
Figure 4. Procedure of CM-PIO cooperative path planning strategy.

Figure 5 .
Figure 5.The paths planned by CM-PIO and PIO.

Figure 6 .
Figure 6.Fitness of optimal individual changing with iteration.

Figure 7 .
Figure 7.The paths with and without cooperative planning.

Figure 8 .
Figure 8. Fitness of CM-PIO cooperative plan changing with iteration.

Figure 9 .
Figure 9. Coordination coefficients of three paths changing with iteration.

Table 2 .
Parameters of the peak terrain.

Table 3 .
Parameters of wind field model.

Table 4 .
Parameters of proposed CM-PIO.

Table 5 .
Coordinates of start point and target point for each UAV.

Table 6 .
Length of planned paths.