# Cooperative Search Optimization of an Unknown Dynamic Target Based on the Modified TPM

Authors: Kaixin Cheng; Di Wu; Tao Hu; Jinjin Wei; Zhifu Tian
Journal: International Journal of Aerospace Engineering
 (2022)

---

## Abstract

Aiming to address the problem of unknown dynamic target trajectory prediction and search path optimization in unmanned aerial vehicle (UAV) swarm path planning, this paper proposes a target search algorithm based on a modified target probability map (TPM). First, using the TPM, the proposed algorithm generates a high-probability distribution region of a target with directionality to fit the target trajectory and realizes the trajectory prediction of an unknown dynamic target. Then, the distributed ant colony (ACO) algorithm and the artificial potential field (APF) algorithm are combined to generate and optimize the UAV swarm search result and return path with the goal of maximizing task execution efficiency. Finally, the Monte Carlo simulation method is used to analyze the effectiveness of the proposed algorithm, and the results are evaluated from five aspects, including the number of targets captured. The simulation results show that under the condition of an unknown dynamic target trajectory, the average target captured rate and average unknown region search rate of the MTPM method were higher than that of the traditional TPM method, and the performance was improved by 14.6% and 10.7%, respectively.


---

## Body

## 1. Introduction

In an increasingly complex battlefield environment, reconnaissance and strike of targets in unknown battlefield areas are important means of war. Unmanned aerial vehicles (UAVs) play an increasingly important role in many fields, including reconnaissance, fire strikes, and electromagnetic interference, and can effectively reduce war losses and casualties. However, due to the limited load capacity and flight capacity of a single UAV, when performing large-scale combat missions, multiple UAVs have been commonly used to carry various types of sensors and loads to establish a wireless sensor network [1], which can improve the combat efficiency through cooperative mission planning between multiple UAVs. At present, much research has been conducted on multi-UAV cooperative task allocation and path planning, and good research results have been achieved [2–4].In practical applications, the UAV path needs to be replanned frequently under the conditions of changing terrain data and dynamic target trajectories. However, the traditional exact solution algorithm has the problem of high time complexity, and the traditional heuristic algorithm has the problem of poor adaptability, which cannot quickly generate the optimal path [5]. Therefore, it is crucial to plan an efficient and reasonable search path for UAVs under the uncertain environments [6]. For this reason, it is necessary to solve the problem of dynamic target trajectory prediction and then to generate and optimize the search path. The parameter setting in the traditional TPM algorithm [7, 8] is not comprehensive enough to make full use of the target motion state data, so it cannot predict the target trajectory well. In addition, the traditional APF algorithm [9] is commonly used to deal with the obstacle avoidance problem but is easy to fall into a local optimum under an uncertain environment.In this paper, we study the problem of multi-UAV cooperative mission path optimization for unknown dynamic targets. Aiming at this problem, this paper proposes a target search algorithm based on a modified target probability map (MTPM). First, the prediction model of unknown dynamic targets is improved, and a directional two-dimensional normal probability distribution map is constructed to fit the target trajectory, realizing the trajectory prediction of unknown dynamic targets. Then, the UAV swarm search and return path with the maximized task execution efficiency are generated and optimized by a joint distributed ACO–APF algorithm. Finally, 30 batches of data are generated by the Monte Carlo method to simulate and analyze the efficiency of the proposed algorithm. The simulation results are evaluated from five aspects: the search rate of an unknown region, the search rate of TPM, the number of captured targets, the obstacle avoidance situation of UAV, and the target captured rate of 30 simulations. The experimental results verify that the proposed algorithm can solve the problems of trajectory prediction and search path optimization of unknown dynamic targets. Under the condition of an unknown dynamic target trajectory, the average target captured rate and average unknown region search rate of the MTPM method were higher than that of the traditional TPM method, and the performance was improved by 14.6% and 10.7%, respectively.The main contributions of this paper can be summarized as follows:(1)
The prediction model of an unknown dynamic target is improved, and a directional two-dimensional normal probability distribution diagram is generated to fit the target trajectory(2)
The joint ACO–APF method is introduced to generate the artificial attractive potential field and the repulsive potential field based on the target trajectory prediction(3)
The state transition rules based on the ACO algorithm and MTPM method are designed to generate the optimal path of the UAV swarm search under the constraints of multiple systemsThe rest of this paper is organized as follows. Section2 describes the environment modeling method of UAV path planning. Section 3 introduces the proposed MTPM collaborative search path generation and optimization method. Section 4 presents the experimental simulations and compares the experimental results. Section 5 discusses open issues. Finally, Section 6 concludes the paper.


## 2. Related Work

In the research on the UAV cooperative mission path planning, considering the perspective of uncertainty, environments can be divided into deterministic environments and uncertain environments. In deterministic environments, the size, shape, and position of obstacles are known before planning, so no unknown dynamic changes in obstacle parameters and other environmental factors can occur. In contrast, in uncertain environments, obstacles and environmental factors are completely or partially unknown in a planning system, obstacles can suddenly appear or move, and there are dynamic disturbances in an environment [6].The cooperative path planning method for multiple UAVs in a deterministic environment can be regarded as an NP-hard combinatorial optimization problem with many constraints [10]. In recent years, many scholars have conducted extensive research on how to establish mathematical models and solve the NP-hard combinatorial optimization problem [11–18]. The main mathematical models include the multitrip salesman model (MTSP), mixed linear integer programming model (MILP), and vehicle scheduling and path planning model (VRP). These models can be solved by heuristic algorithms, such as genetic algorithm (GA) [12], evolutionary algorithm (EA) [13, 14], ant colony (ACO) algorithm [15, 16], and particle swarm optimization (PSO) algorithm [17, 18].Considering the problem of trajectory prediction of unknown dynamic targets, in reference [7], the authors established a rasterized environmental uncertainty map, search probability map, and pheromone map model and determined the prediction probability distribution of dynamic time-sensitive targets by calculating the target transition probability. Based on the target probability graph method, a linear model was developed to predict the location of unknown targets and solve the search problem of linear moving targets [8]. However, this method is not suitable for curved moving targets and cannot effectively solve the trajectory prediction problem of unknown dynamic targets.Aiming at the generation and optimization of the search path, currently, the general research goal is to maximize the efficiency of task execution by effectively modeling the obstacle area [19, 20], and then, the artificial potential field (APF) method [9], the improved RRT∗ algorithm [21–24], and heuristic algorithms have been employed to generate an optimal path while solving spatial conflicts. The Circular Arc Trajectory method was used to model the irregular obstacle area, which could simplify the computational complexity of the obstacle avoidance problem [19]. The improved APF algorithm and a path optimization method for multiobstacle avoidance were proposed to solve the problem of path of UAV swarm falling into deadlock [20]. However, the method was not suitable for dynamic environments. The approaches abovementioned can solve the problems of irregular obstacle area modeling and obstacle avoidance path generation and optimization to a certain extent. However, under the condition of an unknown battlefield environment, there have still been a number of problems, such as an unreachable path, easy falling into a local optimum, and slow convergence speed. In conclusion, the abovementioned methods cannot effectively solve the problem of search path generation and optimization in an unknown battlefield environment.


## 3. Problem Description and Mathematical Model

Aiming at the problem of UAV group mission planning, the battlefield area is considered a two-dimensional area with a fixed height. The task content is as follows. Assume there areNT unknown moving targets and NO unknown obstacle areas in a battlefield area. Nuav fixed-wing UAVs with initial flight resources of Euav are sent to enter the unknown battlefield area to perform target reconnaissance and capture tasks. When the target reaches target Ratk, the target is determined to be captured, and the target stops moving. When the UAV swarm captures all targets or the flight resources are insufficient, the UAV swarm returns to the preset return point. After all UAVs return to the return point, the mission ends.This section designs a UAV swarm mission planning system from three aspects: mission environment modeling, target modeling, and constraint modeling, and builds an environment for the search path optimization algorithm in the next section.


### 3.1. System Environment Modeling

Using the geometric graph method [25], the battlefield area can be rasterized into a set MAP=m,nm=1,2,⋯,Lx,n=1,2,⋯,Ly, with a size of Lx×Ly corresponding to the cell in row m and column n.The UAV track is connected by adjacent grid nodes. A UAV is defined as a fixed-wing UAV with a fixed speedVuav cruise, a detection radius Ruav, and a target capture radius Ratk, having the maximum steering angle of φmax.A setXuav,Yuav,Euav,Θuav describes the UAV motion state, where Xuav and Yuav represent the abscissa and ordinate of the current UAV position, respectively, Euav is the current flight resource surplus of the UAV, and Θuav denoted the current motion direction of the UAV.To facilitate research, this paper makes the following assumptions:(i)
UAVs can communicate with each other without considering the communication delay and communication distance, namely, a UAV can obtain information on other UAVs immediately(ii)
The flight resource consumption is not considered during the take-off and landing of a UAV(iii)
The UAV airborne processor’s data processing time is not considered


### 3.2. Moving Target Modeling

First, the attributes and target trajectory of a target are modeled. It is assumed that a targetTi with coordinates XTi,YTi moves periodically following the trajectory parameter equations given by
(1)XTi=RTxicoswTit0i+t+CTxi,YTi=RTyisinwTit0i+t+CTyi,where i∈1,⋯,NT; set CTyi,CTxi represents the trajectory center coordinates; RTxi and RTyi represent trajectory equation’s horizontal and vertical axis wheelbases, respectively; ωTi is the angular velocity of a target Ti moving on the trajectory; t0i is the initial time of a target Ti motion; t denotes the current system time.SetXT,YT,VT,ΘT describes the motion state of targets, where XT and YT represent the current abscissa and ordinate values of the target, respectively, VT is the current velocity scalar of the target, and ΘT denotes the current motion direction of the target.


### 3.3. System Constraint Modeling

Reference [18] proposed a circular overfitting method to establish an obstacle area model, which solved the problem of modeling and obstacle avoidance in irregular obstacle areas. According to the obstacle region modeling method based on a circular trajectory proposed in this paper, the obstacle region can be considered a circular region, and COxi,COyi and ROi represent the center coordinates and radius of the Obstaclei region, respectively.The following system constraints are defined:(2)Su:duuij≥Rsafe,So:duoij≥Rsafe,Se:euavik≥λEminduEij,Sp:Sopφsopj−suavi⟶,θuavi⟶≤φmax,where Su and So represent the distance constraint between UAVs and the distance constraint between UAVs and obstacle region, respectively, and they should be larger than the minimum threshold; Se represents the UAV flight resource constraint, which needs to meet the minimum return flight resource threshold during the task execution process; Sp represent the UAV flight performance constraint, and the UAV yaw angle needs to be less than the maximum yaw angle threshold.


## 3.1. System Environment Modeling

Using the geometric graph method [25], the battlefield area can be rasterized into a set MAP=m,nm=1,2,⋯,Lx,n=1,2,⋯,Ly, with a size of Lx×Ly corresponding to the cell in row m and column n.The UAV track is connected by adjacent grid nodes. A UAV is defined as a fixed-wing UAV with a fixed speedVuav cruise, a detection radius Ruav, and a target capture radius Ratk, having the maximum steering angle of φmax.A setXuav,Yuav,Euav,Θuav describes the UAV motion state, where Xuav and Yuav represent the abscissa and ordinate of the current UAV position, respectively, Euav is the current flight resource surplus of the UAV, and Θuav denoted the current motion direction of the UAV.To facilitate research, this paper makes the following assumptions:(i)
UAVs can communicate with each other without considering the communication delay and communication distance, namely, a UAV can obtain information on other UAVs immediately(ii)
The flight resource consumption is not considered during the take-off and landing of a UAV(iii)
The UAV airborne processor’s data processing time is not considered


## 3.2. Moving Target Modeling

First, the attributes and target trajectory of a target are modeled. It is assumed that a targetTi with coordinates XTi,YTi moves periodically following the trajectory parameter equations given by
(1)XTi=RTxicoswTit0i+t+CTxi,YTi=RTyisinwTit0i+t+CTyi,where i∈1,⋯,NT; set CTyi,CTxi represents the trajectory center coordinates; RTxi and RTyi represent trajectory equation’s horizontal and vertical axis wheelbases, respectively; ωTi is the angular velocity of a target Ti moving on the trajectory; t0i is the initial time of a target Ti motion; t denotes the current system time.SetXT,YT,VT,ΘT describes the motion state of targets, where XT and YT represent the current abscissa and ordinate values of the target, respectively, VT is the current velocity scalar of the target, and ΘT denotes the current motion direction of the target.


## 3.3. System Constraint Modeling

Reference [18] proposed a circular overfitting method to establish an obstacle area model, which solved the problem of modeling and obstacle avoidance in irregular obstacle areas. According to the obstacle region modeling method based on a circular trajectory proposed in this paper, the obstacle region can be considered a circular region, and COxi,COyi and ROi represent the center coordinates and radius of the Obstaclei region, respectively.The following system constraints are defined:(2)Su:duuij≥Rsafe,So:duoij≥Rsafe,Se:euavik≥λEminduEij,Sp:Sopφsopj−suavi⟶,θuavi⟶≤φmax,where Su and So represent the distance constraint between UAVs and the distance constraint between UAVs and obstacle region, respectively, and they should be larger than the minimum threshold; Se represents the UAV flight resource constraint, which needs to meet the minimum return flight resource threshold during the task execution process; Sp represent the UAV flight performance constraint, and the UAV yaw angle needs to be less than the maximum yaw angle threshold.


## 4. Search Path Generation and Optimization Algorithm

Aiming at the problem of unknown dynamic target trajectory prediction and search path optimization, based on the system modeling presented in the previous section, this section uses the MTPM method to generate a directional two-dimensional normal probability distribution map to fit the target trajectory and realize the trajectory prediction of unknown dynamic targets. The joint ACO–APF algorithm, where the APF algorithm is used to generate the MTPM attraction field, obstacle area repulsion field, and interaircraft repulsion field, is introduced. Combined with the state transition rules based on the ACO algorithm, the UAV swarm search and return path are generated and optimized under the premise of maximizing task execution efficiency.


### 4.1. Modified Target Probability Map

An unknown dynamic target refers to a target with uncertain coordinates, velocity, and direction of motion. To make full use of the limited known information to describe the distribution of unknown dynamic targets in the battlefield area, the target probability map method has usually been used to characterize the existence probability of the target at certain coordinates in the battlefield area. However, due to the limitations on the parameter settings of the two-dimensional normal probability distribution model in the traditional algorithm, the traditional target probability map methods can describe the distribution only of unknown static targets and cannot solve the trajectory prediction problem of unknown dynamic targets.To solve this problem, using the Pauta criterion and the two-dimensional normal probability distribution model, this paper maps the motion state of a target to the variance of the horizontal and vertical axes and generates a directional two-dimensional normal probability distribution map to fit the target motion trajectory and realize the trajectory prediction of unknown dynamic targets.Definepim,n as a distribution probability of a target Ti at coordinates m,n and assume that the motion state of target Ti can be detected as follows. Assume the detection coordinate is m0,n0, the detection interval is τscan, the velocity is vTi, and the direction of motion is vector θTi. Based on the two-dimensional normal distribution model, the probability of the presence pim∗,n∗ of target Ti at coordinates of m∗,n∗~Nm0,n0,σx2,σy2,ρ is given by
(3)pim∗,n∗=12πσx2σy21−ρ2e−1/21−ρ2m∗−m02/2σx2−2ρm∗−m0n∗−n0/σxσy+n∗−n02/2σy2,where a set σx2,σy2,ρ refers to undetermine parameter groups, representing the horizontal axis variance, the vertical axis variance, and horizontal-vertical axes correlation coefficient; τscan represents the time needed for determining the target motion state; θTi represents the angle between vector θTi and the horizontal axis’s positive direction.It is assumed that the motion state of targetTi remains unchanged during the detection interval τscan, namely, target Ti moves from point Am0,n0 to point Bm+cosθTivTiτscan,n+sinθTivTiτscan along the direction of θTi; its motion time is τscan, and its motion distance is δABi=vTiτscan. Substituting the motion state information of target Ti into the variance theorem yields to
(4)σy=tan2θTiσx,θTi∈−π,π.After obtaining the mathematical relationship between variances along the horizontal and vertical axes, the next challenge is to increase the distribution probability of the projection of line segmentAB in the model maximally so that the probability distribution model can fit the trajectory of target Ti and achieve the purpose of trajectory prediction. From the differential point of view, when the line segment AB is short enough, namely, when the detection frequency is high enough, continuous target detection can generate a continuous target trajectory.To make the distribution probability of the projection of the motion trajectory in the model the highest and the trajectory the best fit, the Pauta criterion, namely the three-sigma criterion of normal distribution, is used in this study. Based on this criterion, a mathematical relationship between the variance of three times the motion trajectory and the motion distanceδABi can be obtained as δABi=3σx2+σy2. Then, Equation (6) can be obtained as follows:
(5)σx=vTτscan3tan4θTi+1,σy=tan2θTivTτscan3tan4θTi+1.According to the two-dimensional normal probability distribution theorem and Equation (6), the correlation matrix Φ can be obtained as follows:
(6)Φ=vTτscan3tan4θTi+12ρitanθTivTτscan3tan4θTi+12ρitanθTivTτscan3tan4θTi+12tan2θTivTτscan3tan4θTi+12,ρ∈−1,1,where parameter ρi describes the linear correlation between the horizontal and vertical coordinates, and when θTi∈0,π/2∩−π,−π/2, ρi is positive; otherwise, it is negative.Based on the abovementioned method and the motion state information of multiple targets, multiple directed two-dimensional normal probability distribution models can be generated to fit the transient motion trajectories of multiple targets. Using the normalization method, the target motion state information is brought into the calculation to generate the multitarget probability distribution matrixPMTPM with a size of LX×LY, which can reflect the distribution probability of multitarget in the battlefield area; matrix PMTPM is defined by
(7)PMTPM=p1,1⋯p1,Ly⋮⋱⋮pLx,1⋯pLx,Ly.System initialization and iterative update of matrixPMTPM are performed as follows.Before a task starts, the motion state of an unknown moving target is detectedNscan times, the detection interval is τscan, and the initialization PMTPMk=0 is performed. During the mission execution process, matrix PMTPMk is updated through a multi-UAV cooperative search with a search cycle of τk.


### 4.2. Path Generation and Optimization

The traditional ACO algorithm uses the initialization information for heuristic search through the positive feedback mechanism, which has the advantages of skipping a local optimal solution and approaching the global optimal solution. However, there are also problems of slow convergence and dependence on initialization information, which are specific to carpet search scenarios. The traditional APF algorithm can rapidly generate a direct path through the repulsion attraction mechanism and is suitable for searching for and approaching known targets and rapid escapes from obstacle areas, but there are problems, such as easy falling into a local optimum and unreachable generation path.Combining the two aforementioned algorithms and using the initialization information of matrixPMTPMk=0, a path generation algorithm, which can rapidly generate paths, avoid unreachable paths, and solve the above problems effectively, can be developed. The specific combined method is as follows. When the APF information is strong, the APF method is used to approach targets or escape from obstacles; otherwise, the ACO heuristic transfer method is used to generate an optimal search path.


#### 4.2.1. Pheromone Map and APF Modeling

First, a pheromone matrixAco with a size of LX×LY is constructed to describe the current battlefield area search as follows:
(8)Aco=aco1,1⋯aco1,Ly⋮⋱⋮acoLx,1⋯acoLx,Ly,acom,n∈acomin,acomax,where acom,n represents the pheromone value of a position m,n; the higher the pheromone value is, the higher the search value of the position will be.System initialization and iterative update are performed on the matrixAco as follows. Before the task starts, the pheromone matrix value is initialized to acomin+acomax/2. During the mission, the global pheromone increases with the system iteration number with an increment of Acorise/τk, and UAVs cooperatively detect the pheromone value within the update distance. The update formula is given by
(9)sacojxuavik−a2+yuavik−b2≤Ruavduaij=xuavik−a2+yuavik−b2Acoa,b=1−λAco1−duaij4Ruav4Acoa,b,where λAco represents the pheromone elimination coefficient. Then, by using the APF algorithm, the repulsion potential field of the UAV is constructed with an obstacle area and a UAV as a center, and the collision avoidance paths between the UAV and the obstacle area and adjacent aircraft are generated. The attraction potential field of the UAV is centered on the target and the return point, and the approach search and return path of the UAV are generated.Next, set the farthest distance of the repulsive potential field toRrep, the minimum safe distance around a drone is Rsafe, the modulus of the repulsive potential field force is Fuij,Foij∈Rsafe/Rrep2,1, and it is negatively correlated with the distance; then, the resultant repulsive force is obtained by
(10)Frepxi=∑j=1NOFoxij+∑j=1NTFuxij,Frepyi=∑j=1NOFoyij+∑j=1NTFuyij,where Frepxi and Frepyi represent the transverse and longitudinal components of the resultant repulsive force; Fuxij and Fuyij denote the transverse and longitudinal components of the repulsive force Fuij between UAVi and UAVj; Foxij and Foyij are the transverse and longitudinal components of the repulsive force Foij between UAVi and ObstaclejThe farthest action distance of the target attraction potential field isRuav, indicating that when the drone detects the target, the target can attract the drone, and the attraction mode is given by Fattij∈Ratk/Ruav2,1.Then, set the action distance of the return point attraction potential field to be unlimited; the return attraction is a unit vector, and the potential field action condition is given byeuavik≤λEminduEij, where euavi represents the current flight resource surplus of UAVi, duEij is the distance from UAVi to the return point j, and λE is the return parameter. If the condition is judged to be true, then the UAV flight resources are insufficient and need to return immediately; parameter λE is used to adjust the UAV return time point.When a drone is affected by multiple attraction potential fields at the same time, only the attraction potential field force with the maximum modulus is retained; namely, the drone can be subjected to the attraction potential field force only in one direction at a time, as given by(11)Fatti=FattargmaxjFatti,j,where Fattij represents the attraction vector of the target pair, argmaxFattijj is the attraction potential field force number with the largest modulus, and Fatti denotes the current attraction vector.In summary, the generation formula of the artificial potential forceFi of UAVi at position xuavi,yuavi is defined by
(12)Fxi=Frepxi+Fattxi,Fyi=Frepyi+Fattyi..


#### 4.2.2. Route Generation Algorithm

The optimal decision-making method based on a battlefield environment of a UAV is as follows.When APF information is strong, the APF method is used to approach targets or escape from obstacles; otherwise, the ACO heuristic transfer method is used to determine an optimal search path according to the target probability graph matrix and the pheromone matrix.First, an optional grid setSop can be obtained based on the flight performance constraint set Sp:Sopφsopj−suavi⟶,θuavi⟶≤φmax by Equation (3), where sopj represents the coordinates of the optional grid j; suavi and θuavi denote the current coordinates and motion direction of UAVi; φx,y is the angle between the two vectors. Then, the optimal path grid number jbest is obtained according to the optimal decision method; namely, sopjbest is the next transfer coordinate of UAVi.WhenFi⟶≠0, that is, when the APF information exists, the APF method is used to generate the optimal search path. The grid in Sop with the smallest angle with vector Fi represents an optimal grid, and it is given by
(13)jbest=argminjφsopj−suavi⟶,Fi⟶.WhenFi⟶=0, that is, when the APF information does not exist, the ACO heuristic transfer method is used to generate an optimal search path. The target probability graph matrix PMTPMk and pheromone matrix Acok are substituted into the ACO algorithm transfer probability formula, and the Bayesian probability pijk of sopj is calculated. Then, the roulette algorithm is used to generate an optimal grid as follows:
(14)pijk=PMTPMkn,mijα×Acokn,mijβ∑soj∈SopPMTPMkn,mijα×Acokn,mijβ,where parameter α is a heuristic factor, which represents the relative importance of the target probability distribution matrix PMTPMk, and parameter β is a pheromone factor, which represents the relative importance of pheromone matrix Acok.


#### 4.2.3. System Iterative Update

After obtaining the optimal grid numberjbest, sopjbest represents the next grid coordinates of a UAV, the UAV path is generated, and the system is iteratively updated. The updating process is as follows. Update a UAV’s position suavik+1, flight resource surplus euavik+1, and motion direction θuavik+1 as follows:
(15)suavik+1=sopjbestkeuavik+1=euavik−dsuavik,sopjbestkθuavik+1=θsopjbestk⟶−suavik⟶.The algorithm flow of the path generation and optimization method based on the improved target probability graph is shown in Figure1.Figure 1
The block diagram of the proposed path generation and optimization method based on the MTPM.


## 4.1. Modified Target Probability Map

An unknown dynamic target refers to a target with uncertain coordinates, velocity, and direction of motion. To make full use of the limited known information to describe the distribution of unknown dynamic targets in the battlefield area, the target probability map method has usually been used to characterize the existence probability of the target at certain coordinates in the battlefield area. However, due to the limitations on the parameter settings of the two-dimensional normal probability distribution model in the traditional algorithm, the traditional target probability map methods can describe the distribution only of unknown static targets and cannot solve the trajectory prediction problem of unknown dynamic targets.To solve this problem, using the Pauta criterion and the two-dimensional normal probability distribution model, this paper maps the motion state of a target to the variance of the horizontal and vertical axes and generates a directional two-dimensional normal probability distribution map to fit the target motion trajectory and realize the trajectory prediction of unknown dynamic targets.Definepim,n as a distribution probability of a target Ti at coordinates m,n and assume that the motion state of target Ti can be detected as follows. Assume the detection coordinate is m0,n0, the detection interval is τscan, the velocity is vTi, and the direction of motion is vector θTi. Based on the two-dimensional normal distribution model, the probability of the presence pim∗,n∗ of target Ti at coordinates of m∗,n∗~Nm0,n0,σx2,σy2,ρ is given by
(3)pim∗,n∗=12πσx2σy21−ρ2e−1/21−ρ2m∗−m02/2σx2−2ρm∗−m0n∗−n0/σxσy+n∗−n02/2σy2,where a set σx2,σy2,ρ refers to undetermine parameter groups, representing the horizontal axis variance, the vertical axis variance, and horizontal-vertical axes correlation coefficient; τscan represents the time needed for determining the target motion state; θTi represents the angle between vector θTi and the horizontal axis’s positive direction.It is assumed that the motion state of targetTi remains unchanged during the detection interval τscan, namely, target Ti moves from point Am0,n0 to point Bm+cosθTivTiτscan,n+sinθTivTiτscan along the direction of θTi; its motion time is τscan, and its motion distance is δABi=vTiτscan. Substituting the motion state information of target Ti into the variance theorem yields to
(4)σy=tan2θTiσx,θTi∈−π,π.After obtaining the mathematical relationship between variances along the horizontal and vertical axes, the next challenge is to increase the distribution probability of the projection of line segmentAB in the model maximally so that the probability distribution model can fit the trajectory of target Ti and achieve the purpose of trajectory prediction. From the differential point of view, when the line segment AB is short enough, namely, when the detection frequency is high enough, continuous target detection can generate a continuous target trajectory.To make the distribution probability of the projection of the motion trajectory in the model the highest and the trajectory the best fit, the Pauta criterion, namely the three-sigma criterion of normal distribution, is used in this study. Based on this criterion, a mathematical relationship between the variance of three times the motion trajectory and the motion distanceδABi can be obtained as δABi=3σx2+σy2. Then, Equation (6) can be obtained as follows:
(5)σx=vTτscan3tan4θTi+1,σy=tan2θTivTτscan3tan4θTi+1.According to the two-dimensional normal probability distribution theorem and Equation (6), the correlation matrix Φ can be obtained as follows:
(6)Φ=vTτscan3tan4θTi+12ρitanθTivTτscan3tan4θTi+12ρitanθTivTτscan3tan4θTi+12tan2θTivTτscan3tan4θTi+12,ρ∈−1,1,where parameter ρi describes the linear correlation between the horizontal and vertical coordinates, and when θTi∈0,π/2∩−π,−π/2, ρi is positive; otherwise, it is negative.Based on the abovementioned method and the motion state information of multiple targets, multiple directed two-dimensional normal probability distribution models can be generated to fit the transient motion trajectories of multiple targets. Using the normalization method, the target motion state information is brought into the calculation to generate the multitarget probability distribution matrixPMTPM with a size of LX×LY, which can reflect the distribution probability of multitarget in the battlefield area; matrix PMTPM is defined by
(7)PMTPM=p1,1⋯p1,Ly⋮⋱⋮pLx,1⋯pLx,Ly.System initialization and iterative update of matrixPMTPM are performed as follows.Before a task starts, the motion state of an unknown moving target is detectedNscan times, the detection interval is τscan, and the initialization PMTPMk=0 is performed. During the mission execution process, matrix PMTPMk is updated through a multi-UAV cooperative search with a search cycle of τk.


## 4.2. Path Generation and Optimization

The traditional ACO algorithm uses the initialization information for heuristic search through the positive feedback mechanism, which has the advantages of skipping a local optimal solution and approaching the global optimal solution. However, there are also problems of slow convergence and dependence on initialization information, which are specific to carpet search scenarios. The traditional APF algorithm can rapidly generate a direct path through the repulsion attraction mechanism and is suitable for searching for and approaching known targets and rapid escapes from obstacle areas, but there are problems, such as easy falling into a local optimum and unreachable generation path.Combining the two aforementioned algorithms and using the initialization information of matrixPMTPMk=0, a path generation algorithm, which can rapidly generate paths, avoid unreachable paths, and solve the above problems effectively, can be developed. The specific combined method is as follows. When the APF information is strong, the APF method is used to approach targets or escape from obstacles; otherwise, the ACO heuristic transfer method is used to generate an optimal search path.


### 4.2.1. Pheromone Map and APF Modeling

First, a pheromone matrixAco with a size of LX×LY is constructed to describe the current battlefield area search as follows:
(8)Aco=aco1,1⋯aco1,Ly⋮⋱⋮acoLx,1⋯acoLx,Ly,acom,n∈acomin,acomax,where acom,n represents the pheromone value of a position m,n; the higher the pheromone value is, the higher the search value of the position will be.System initialization and iterative update are performed on the matrixAco as follows. Before the task starts, the pheromone matrix value is initialized to acomin+acomax/2. During the mission, the global pheromone increases with the system iteration number with an increment of Acorise/τk, and UAVs cooperatively detect the pheromone value within the update distance. The update formula is given by
(9)sacojxuavik−a2+yuavik−b2≤Ruavduaij=xuavik−a2+yuavik−b2Acoa,b=1−λAco1−duaij4Ruav4Acoa,b,where λAco represents the pheromone elimination coefficient. Then, by using the APF algorithm, the repulsion potential field of the UAV is constructed with an obstacle area and a UAV as a center, and the collision avoidance paths between the UAV and the obstacle area and adjacent aircraft are generated. The attraction potential field of the UAV is centered on the target and the return point, and the approach search and return path of the UAV are generated.Next, set the farthest distance of the repulsive potential field toRrep, the minimum safe distance around a drone is Rsafe, the modulus of the repulsive potential field force is Fuij,Foij∈Rsafe/Rrep2,1, and it is negatively correlated with the distance; then, the resultant repulsive force is obtained by
(10)Frepxi=∑j=1NOFoxij+∑j=1NTFuxij,Frepyi=∑j=1NOFoyij+∑j=1NTFuyij,where Frepxi and Frepyi represent the transverse and longitudinal components of the resultant repulsive force; Fuxij and Fuyij denote the transverse and longitudinal components of the repulsive force Fuij between UAVi and UAVj; Foxij and Foyij are the transverse and longitudinal components of the repulsive force Foij between UAVi and ObstaclejThe farthest action distance of the target attraction potential field isRuav, indicating that when the drone detects the target, the target can attract the drone, and the attraction mode is given by Fattij∈Ratk/Ruav2,1.Then, set the action distance of the return point attraction potential field to be unlimited; the return attraction is a unit vector, and the potential field action condition is given byeuavik≤λEminduEij, where euavi represents the current flight resource surplus of UAVi, duEij is the distance from UAVi to the return point j, and λE is the return parameter. If the condition is judged to be true, then the UAV flight resources are insufficient and need to return immediately; parameter λE is used to adjust the UAV return time point.When a drone is affected by multiple attraction potential fields at the same time, only the attraction potential field force with the maximum modulus is retained; namely, the drone can be subjected to the attraction potential field force only in one direction at a time, as given by(11)Fatti=FattargmaxjFatti,j,where Fattij represents the attraction vector of the target pair, argmaxFattijj is the attraction potential field force number with the largest modulus, and Fatti denotes the current attraction vector.In summary, the generation formula of the artificial potential forceFi of UAVi at position xuavi,yuavi is defined by
(12)Fxi=Frepxi+Fattxi,Fyi=Frepyi+Fattyi..


### 4.2.2. Route Generation Algorithm

The optimal decision-making method based on a battlefield environment of a UAV is as follows.When APF information is strong, the APF method is used to approach targets or escape from obstacles; otherwise, the ACO heuristic transfer method is used to determine an optimal search path according to the target probability graph matrix and the pheromone matrix.First, an optional grid setSop can be obtained based on the flight performance constraint set Sp:Sopφsopj−suavi⟶,θuavi⟶≤φmax by Equation (3), where sopj represents the coordinates of the optional grid j; suavi and θuavi denote the current coordinates and motion direction of UAVi; φx,y is the angle between the two vectors. Then, the optimal path grid number jbest is obtained according to the optimal decision method; namely, sopjbest is the next transfer coordinate of UAVi.WhenFi⟶≠0, that is, when the APF information exists, the APF method is used to generate the optimal search path. The grid in Sop with the smallest angle with vector Fi represents an optimal grid, and it is given by
(13)jbest=argminjφsopj−suavi⟶,Fi⟶.WhenFi⟶=0, that is, when the APF information does not exist, the ACO heuristic transfer method is used to generate an optimal search path. The target probability graph matrix PMTPMk and pheromone matrix Acok are substituted into the ACO algorithm transfer probability formula, and the Bayesian probability pijk of sopj is calculated. Then, the roulette algorithm is used to generate an optimal grid as follows:
(14)pijk=PMTPMkn,mijα×Acokn,mijβ∑soj∈SopPMTPMkn,mijα×Acokn,mijβ,where parameter α is a heuristic factor, which represents the relative importance of the target probability distribution matrix PMTPMk, and parameter β is a pheromone factor, which represents the relative importance of pheromone matrix Acok.


### 4.2.3. System Iterative Update

After obtaining the optimal grid numberjbest, sopjbest represents the next grid coordinates of a UAV, the UAV path is generated, and the system is iteratively updated. The updating process is as follows. Update a UAV’s position suavik+1, flight resource surplus euavik+1, and motion direction θuavik+1 as follows:
(15)suavik+1=sopjbestkeuavik+1=euavik−dsuavik,sopjbestkθuavik+1=θsopjbestk⟶−suavik⟶.The algorithm flow of the path generation and optimization method based on the improved target probability graph is shown in Figure1.Figure 1
The block diagram of the proposed path generation and optimization method based on the MTPM.


## 4.2.1. Pheromone Map and APF Modeling

First, a pheromone matrixAco with a size of LX×LY is constructed to describe the current battlefield area search as follows:
(8)Aco=aco1,1⋯aco1,Ly⋮⋱⋮acoLx,1⋯acoLx,Ly,acom,n∈acomin,acomax,where acom,n represents the pheromone value of a position m,n; the higher the pheromone value is, the higher the search value of the position will be.System initialization and iterative update are performed on the matrixAco as follows. Before the task starts, the pheromone matrix value is initialized to acomin+acomax/2. During the mission, the global pheromone increases with the system iteration number with an increment of Acorise/τk, and UAVs cooperatively detect the pheromone value within the update distance. The update formula is given by
(9)sacojxuavik−a2+yuavik−b2≤Ruavduaij=xuavik−a2+yuavik−b2Acoa,b=1−λAco1−duaij4Ruav4Acoa,b,where λAco represents the pheromone elimination coefficient. Then, by using the APF algorithm, the repulsion potential field of the UAV is constructed with an obstacle area and a UAV as a center, and the collision avoidance paths between the UAV and the obstacle area and adjacent aircraft are generated. The attraction potential field of the UAV is centered on the target and the return point, and the approach search and return path of the UAV are generated.Next, set the farthest distance of the repulsive potential field toRrep, the minimum safe distance around a drone is Rsafe, the modulus of the repulsive potential field force is Fuij,Foij∈Rsafe/Rrep2,1, and it is negatively correlated with the distance; then, the resultant repulsive force is obtained by
(10)Frepxi=∑j=1NOFoxij+∑j=1NTFuxij,Frepyi=∑j=1NOFoyij+∑j=1NTFuyij,where Frepxi and Frepyi represent the transverse and longitudinal components of the resultant repulsive force; Fuxij and Fuyij denote the transverse and longitudinal components of the repulsive force Fuij between UAVi and UAVj; Foxij and Foyij are the transverse and longitudinal components of the repulsive force Foij between UAVi and ObstaclejThe farthest action distance of the target attraction potential field isRuav, indicating that when the drone detects the target, the target can attract the drone, and the attraction mode is given by Fattij∈Ratk/Ruav2,1.Then, set the action distance of the return point attraction potential field to be unlimited; the return attraction is a unit vector, and the potential field action condition is given byeuavik≤λEminduEij, where euavi represents the current flight resource surplus of UAVi, duEij is the distance from UAVi to the return point j, and λE is the return parameter. If the condition is judged to be true, then the UAV flight resources are insufficient and need to return immediately; parameter λE is used to adjust the UAV return time point.When a drone is affected by multiple attraction potential fields at the same time, only the attraction potential field force with the maximum modulus is retained; namely, the drone can be subjected to the attraction potential field force only in one direction at a time, as given by(11)Fatti=FattargmaxjFatti,j,where Fattij represents the attraction vector of the target pair, argmaxFattijj is the attraction potential field force number with the largest modulus, and Fatti denotes the current attraction vector.In summary, the generation formula of the artificial potential forceFi of UAVi at position xuavi,yuavi is defined by
(12)Fxi=Frepxi+Fattxi,Fyi=Frepyi+Fattyi..


## 4.2.2. Route Generation Algorithm

The optimal decision-making method based on a battlefield environment of a UAV is as follows.When APF information is strong, the APF method is used to approach targets or escape from obstacles; otherwise, the ACO heuristic transfer method is used to determine an optimal search path according to the target probability graph matrix and the pheromone matrix.First, an optional grid setSop can be obtained based on the flight performance constraint set Sp:Sopφsopj−suavi⟶,θuavi⟶≤φmax by Equation (3), where sopj represents the coordinates of the optional grid j; suavi and θuavi denote the current coordinates and motion direction of UAVi; φx,y is the angle between the two vectors. Then, the optimal path grid number jbest is obtained according to the optimal decision method; namely, sopjbest is the next transfer coordinate of UAVi.WhenFi⟶≠0, that is, when the APF information exists, the APF method is used to generate the optimal search path. The grid in Sop with the smallest angle with vector Fi represents an optimal grid, and it is given by
(13)jbest=argminjφsopj−suavi⟶,Fi⟶.WhenFi⟶=0, that is, when the APF information does not exist, the ACO heuristic transfer method is used to generate an optimal search path. The target probability graph matrix PMTPMk and pheromone matrix Acok are substituted into the ACO algorithm transfer probability formula, and the Bayesian probability pijk of sopj is calculated. Then, the roulette algorithm is used to generate an optimal grid as follows:
(14)pijk=PMTPMkn,mijα×Acokn,mijβ∑soj∈SopPMTPMkn,mijα×Acokn,mijβ,where parameter α is a heuristic factor, which represents the relative importance of the target probability distribution matrix PMTPMk, and parameter β is a pheromone factor, which represents the relative importance of pheromone matrix Acok.


## 4.2.3. System Iterative Update

After obtaining the optimal grid numberjbest, sopjbest represents the next grid coordinates of a UAV, the UAV path is generated, and the system is iteratively updated. The updating process is as follows. Update a UAV’s position suavik+1, flight resource surplus euavik+1, and motion direction θuavik+1 as follows:
(15)suavik+1=sopjbestkeuavik+1=euavik−dsuavik,sopjbestkθuavik+1=θsopjbestk⟶−suavik⟶.The algorithm flow of the path generation and optimization method based on the improved target probability graph is shown in Figure1.Figure 1
The block diagram of the proposed path generation and optimization method based on the MTPM.


## 5. Simulation Experiments

This section presents the numerical simulation results of collaborative search tasks in uncertain dynamic environments. In the simulations, it was assumed that there wereNT unknown moving targets and NO unknown obstacle areas in the battlefield area of a fixed height; Nuav fixed-wing UAVs with the initial flight resources Euav were dispatched to perform reconnaissance and capture tasks. A target was judged to be captured when the distance between the UAV and the target was less than Ratk, and then, the target stopped moving. When the UAV captured all targets or flight resources were insufficient, the UAV returned to the preset return point. After all UAVs returned to the return point, the mission ended.


### 5.1. Parameter Settings

The mission of a UAV group was as follows. It was assumed that 20 patrol vehicles (unknown dynamic targets) patrolled at a constant speed around a certain point in a mission area. Then, four UAVs, which took off from the airports located at the four vertices of the mission area map, were dispatched to enter the mission area and perform a search and capture patrol vehicle mission. To improve the search efficiency of a UAV, the patrol vehicle was detected at intervals ofτscan before the UAV swarm entered the task area, and the detection results were imported into the UAV swarm path planning system as a prior information. There were multiple unknown no-fly zones (unknown obstacle zones) in the mission area, and the minimum safe distance Rsafe between the UAV group and the no-fly zone had to be maintained.The battlefield area size was200km×200km, and it was rasterized into 1,000grid×1,000grid, so 1 grid=0.2km. The UAV path defined the system iteration step size to be τk=5s, and the system clock was tk=τk×k−1. The simulation parameters are given in Table 1.Table 1
Simulation parameters.
ParameterSimulated valueActual valueDetection rangeRuav50 grid10 kmCatching rangeRatk5 grid1 kmUAV speedVuav1,2grid/τk(144,200) km/hFlight resourcesEuav2,400 grid480 kmSafety distanceRsafe4 grid0.8 kmRepulsive force distanceRrep10 grid2 kmSystem step lengthτk1 step5 sSampling periodτscan120 step600 sThe Monte Carlo method was under the condition that the intersection of the target motion trajectory, and the obstacle area was empty. To mobilize the UAV group fully to perform a wide search of the task area, the special situation of the target trajectory concentration should be avoided as much as possible. Therefore, according to the two-dimensional uniform distribution, multiple batches of relevant parameters were generated. Each batch contained 10 moving targets and 20 obstacle areas. The relevant parameters included the trajectory center of a moving target, long-axis wheelbase, short-axis wheelbase, the angular velocity of motion, the initial phase of motion, motion direction (clockwise or counterclockwise), and center and radius of the obstacle area. The ranges of parameters were set as presented in Table2.Table 2
Target and obstacle area parameters.
ParameterSimulation value rangeActual value rangeTarget trajectory radiusRT(20,50) grid(4, 10) kmTarget moving speedVT0.4,0.7grid/τk(57,100) km/hObstacle area radiusRO(20,50) grid(4, 10) kmThe Monte Carlo method was used to generate 30 batches of experimental data, and the experimental results were described and compared. The distribution of unknown dynamic target trajectory and obstacle area is shown in Figure2. As shown in Figure 2, 10 moving targets and 20 obstacle areas were evenly distributed in the task area. The moving target moved in a periodic closed-loop around the center of the trajectory as shown in the blue curve, and the obstacle areas overlapped, forming an irregular obstacle area as shown in the black area. There was no intersection between the moving target trajectory and the obstacle area to avoid the situation of a drone entering the obstacle area to search for a target.Figure 2
Unknown dynamic target trajectory and obstacle area distribution graph.


### 5.2. Experimental Results and Analysis

Aiming at the problem of unknown dynamic target trajectory prediction and search path optimization, the traditional TPM method is compared. The simulation results are shown in Figures3 and 4, where the initial distribution of target probability (k=0) and the cooperative search path are presented, respectively. Then, the effectiveness of the proposed method in an uncertain dynamic environment was evaluated using several evaluation indices, which included the search rate of an unknown region, the search rate of a target probability distribution, the obstacle avoidance of UAV, the number of captured targets, the number of captured targets in 30 times’ simulations, and the unknown region search rate in 30 times’ simulations, shown in Figures 5–8, respectively. The task execution efficiency was comprehensively evaluated in comparison with the traditional TPM method.Figure 3
The initial distribution of target probabilities (k=0).
(a)
MTPM (3D)(b)
MTPM (2D)(c)
TPM (3D)(d)
TPM (2D)Figure 4
The cooperative search results.
(a)
Path of MTPM(b)
Path of TPMFigure 5
Comparison on the search rate with the iteration number between the two methods.
(a)
The changing trend of the unknown region search rate with the iteration number(b)
The changing trend of the TPM search rate with the iteration numberFigure 6
Comparison on the number of captured targets with the iteration number between the two methods.Figure 7
Comparison on the UAV obstacle avoidance result between the two methods.Figure 8
Comparison results in 30 simulations between the two methods.
(a)
Comparison results of targets captured in 30 simulations between the two methods(b)
Comparison results of unknown region search rate in 30 simulations between the two methodsAfter importing the prior information into the system, the initial target probability distribution matrixPTPMk=0 was generated using the proposed MTPM method and the traditional TPM method, as shown in Figure 3.The initial distribution of target probability reflected a target’s movement in the task area at the beginning of the task. The higher the target probability was, the more likely the target was to move around the position. In Figures3(a) and 3(b), the results of the MTPM method are presented. In Figures 3(c) and 3(d), the results of the TPM method are illustrated. Figures 3(a) and 3(c) show the 3D images, and Figures 3(b) and 3(d) display the 2D images. Further, in Figures 3(c) and 3(d), the high-probability area of the target in the distribution map obtained by the traditional TPM method is dotted. Compared with the real motion trajectory of the target presented in Figure 2, the high-probability area of the target is discontinuous and has poor directivity, which cannot predict the target trajectory. Further, as shown in Figures 3(a) and 3(b), the distribution map obtained by the MTPM method showed a banded distribution of the high-probability region of the target. Compared with the actual motion trajectory of the target, which is displayed in Figure 2, the trajectory presented in Figure 3 could accurately reflect the motion trend and process of the target, indicating that the proposed MTPM method is more suitable for the motion trajectory prediction and fitting of the curve trajectory target than the traditional TPM method.After the mission began, the UAV planned the search path according to the initial distribution of the target probability and reduced the target probability value of the explored area as the system iterated, updating the target probability distribution map. Figures4(a) and 4(b) show the drone’s search path. Figure 4(a) shows the results of the MTPM method. Further, Figure 4(b) shows the results of the traditional TPM method.As shown in Figure4(a), under the combined action of the distributed ACO algorithm and the APF algorithm, after the UAV effectively searched the high-probability area of the target, it captured all targets and successfully returned to the predefined return point, and the system iterations ended at k=1639. The results in Figure 4(b) demonstrate that the UAV searched the high-probability area of the target, but due to the lack of effective target trajectory prediction, the drown cannot capture all targets, and the system iteration was ended at k=1686. The comparison of the results showed that the initial distribution of the target probability generated by the MTPM method could accurately describe the motion trend of an unknown target. The joint distributed ACO–APF algorithm could correctly guide the UAV to search the high-probability area of the target effectively and jump out of the path deadlock.The unknown area search rate and the target probability distribution search rate were calculated by Equations (16) and (17).
(16)rcovk=Scovk/SLx×Ly,(17)rTPMk=1−∑b=1Ly∑a=1LxPTPMka,b∑b=1Ly∑a=1LxPTPMk=0a,b.Figure5(a) reflects the comparison of the two methods on the unknown region search rate which reflects the extent of the drone’s search for unknown region, and Figure 5(b) reflects the comparison of the two methods on the TPM search rate which reflects the extent of the drone’s search for the high-probability distribution region of targets. With the iteration of the system, the drone’s search for unknown region gradually deepens based on the TPM matrix. In Figure 5, the end point represents the number of system iterations at the end of the task.Figure6 reflects the comparison of the two methods on the number of target capture which reflects the task efficiency of the drone. With the iteration of the system, the drone continuously discovers and captures the targets when searching the unknown region. In Figure 6, the end point represents the number of system iterations at the end of the task.As shown in Figures5 and 6, the MTPM method completed the iteration process at k=1639, achieving a search rate of an unknown region of 44.7%, a search rate of the target probability distribution of 97.3%, and the number of captured targets of 10, which were better than those of the traditional TPM method. The three mentioned metrics of the traditional TPM method were 41.4%, 88.5%, and 9 at k=1686, indicating that the effective trajectory prediction of the moving target by the MTPM method could accurately guide a UAV to search the path planning and achieve a better task execution effect.Figure7 reflects the comparison of the two methods on the distance to obstacle which reflects the obstacle avoidance situation of the drone. The vertical axis in Figure 5 represents the minimum distance from the drone to the obstacle area. The black solid line that changes with time corresponds to the data of the MTPM method, and the black dotted line corresponds to the data of the TPM method. The constant bold black solid line represents the maximum repulsive force distance Rrep(Rrep=2km) in the APF algorithm. If the distance from the UAV to the obstacle area exceeds Rrep, the repulsive force from the obstacle area would be no longer applied to the UAV. The constant bold black dotted line represents the safety distance Rsafe(Rsafe=0.8km) of the UAV in the APF algorithm. If the distance from the UAV to the obstacle area is less than Rsafe, a collision would be happened between the UAV and the obstacle area.As shown in Figure7, the minimum distance between the UAV and the obstacle area is always longer than the safe distance Rsafe, which shows that the obstacle repulsion potential field can avoid collisions between the UAV and the obstacle area. The reason is that the modulus of the repulsive force is inversely proportional to the square of the distance, so that the UAV can avoid obstacles at the right time and position.Figure8 reflects the comparison results in 30 simulations of the two methods. In Figure 8, the black color corresponds to the data of the MTPM method, and the blue color corresponds to the data of the TPM method.In Figure8(a), the vertical axis represents the number of targets captured by the drone. The dotted lines represent to the number of captured targets in 30 simulations, and the solid lines represent to the average number of captured targets in 30 simulations. After the subsimulation, the average number of captured targets of the MTPM method was 9.63, and the average number of captured targets of the TPM method was 8.4.In Figure8(b), the vertical axis represents the unknown region search rate. The dotted lines represent to the unknown region search rate in 30 simulations, and the solid lines represent to the average unknown region search rate in 30 simulations. After the subsimulation, the average unknown region search rate of the MTPM method was 43.3%, and the average unknown region search rate of the TPM method was 39.1%.According to the above-presented evaluation results, after the subsimulation, the average target captured rate and average unknown region search rate of the MTPM method were higher than that of the traditional TPM method, and the performance was improved by 14.6% and 10.7%, respectively, indicating that the traditional TPM method is deficient in trajectory prediction of dynamic targets. In contrast, the proposed method could effectively solve the search path optimization problem of unknown dynamic targets.


## 5.1. Parameter Settings

The mission of a UAV group was as follows. It was assumed that 20 patrol vehicles (unknown dynamic targets) patrolled at a constant speed around a certain point in a mission area. Then, four UAVs, which took off from the airports located at the four vertices of the mission area map, were dispatched to enter the mission area and perform a search and capture patrol vehicle mission. To improve the search efficiency of a UAV, the patrol vehicle was detected at intervals ofτscan before the UAV swarm entered the task area, and the detection results were imported into the UAV swarm path planning system as a prior information. There were multiple unknown no-fly zones (unknown obstacle zones) in the mission area, and the minimum safe distance Rsafe between the UAV group and the no-fly zone had to be maintained.The battlefield area size was200km×200km, and it was rasterized into 1,000grid×1,000grid, so 1 grid=0.2km. The UAV path defined the system iteration step size to be τk=5s, and the system clock was tk=τk×k−1. The simulation parameters are given in Table 1.Table 1
Simulation parameters.
ParameterSimulated valueActual valueDetection rangeRuav50 grid10 kmCatching rangeRatk5 grid1 kmUAV speedVuav1,2grid/τk(144,200) km/hFlight resourcesEuav2,400 grid480 kmSafety distanceRsafe4 grid0.8 kmRepulsive force distanceRrep10 grid2 kmSystem step lengthτk1 step5 sSampling periodτscan120 step600 sThe Monte Carlo method was under the condition that the intersection of the target motion trajectory, and the obstacle area was empty. To mobilize the UAV group fully to perform a wide search of the task area, the special situation of the target trajectory concentration should be avoided as much as possible. Therefore, according to the two-dimensional uniform distribution, multiple batches of relevant parameters were generated. Each batch contained 10 moving targets and 20 obstacle areas. The relevant parameters included the trajectory center of a moving target, long-axis wheelbase, short-axis wheelbase, the angular velocity of motion, the initial phase of motion, motion direction (clockwise or counterclockwise), and center and radius of the obstacle area. The ranges of parameters were set as presented in Table2.Table 2
Target and obstacle area parameters.
ParameterSimulation value rangeActual value rangeTarget trajectory radiusRT(20,50) grid(4, 10) kmTarget moving speedVT0.4,0.7grid/τk(57,100) km/hObstacle area radiusRO(20,50) grid(4, 10) kmThe Monte Carlo method was used to generate 30 batches of experimental data, and the experimental results were described and compared. The distribution of unknown dynamic target trajectory and obstacle area is shown in Figure2. As shown in Figure 2, 10 moving targets and 20 obstacle areas were evenly distributed in the task area. The moving target moved in a periodic closed-loop around the center of the trajectory as shown in the blue curve, and the obstacle areas overlapped, forming an irregular obstacle area as shown in the black area. There was no intersection between the moving target trajectory and the obstacle area to avoid the situation of a drone entering the obstacle area to search for a target.Figure 2
Unknown dynamic target trajectory and obstacle area distribution graph.


## 5.2. Experimental Results and Analysis

Aiming at the problem of unknown dynamic target trajectory prediction and search path optimization, the traditional TPM method is compared. The simulation results are shown in Figures3 and 4, where the initial distribution of target probability (k=0) and the cooperative search path are presented, respectively. Then, the effectiveness of the proposed method in an uncertain dynamic environment was evaluated using several evaluation indices, which included the search rate of an unknown region, the search rate of a target probability distribution, the obstacle avoidance of UAV, the number of captured targets, the number of captured targets in 30 times’ simulations, and the unknown region search rate in 30 times’ simulations, shown in Figures 5–8, respectively. The task execution efficiency was comprehensively evaluated in comparison with the traditional TPM method.Figure 3
The initial distribution of target probabilities (k=0).
(a)
MTPM (3D)(b)
MTPM (2D)(c)
TPM (3D)(d)
TPM (2D)Figure 4
The cooperative search results.
(a)
Path of MTPM(b)
Path of TPMFigure 5
Comparison on the search rate with the iteration number between the two methods.
(a)
The changing trend of the unknown region search rate with the iteration number(b)
The changing trend of the TPM search rate with the iteration numberFigure 6
Comparison on the number of captured targets with the iteration number between the two methods.Figure 7
Comparison on the UAV obstacle avoidance result between the two methods.Figure 8
Comparison results in 30 simulations between the two methods.
(a)
Comparison results of targets captured in 30 simulations between the two methods(b)
Comparison results of unknown region search rate in 30 simulations between the two methodsAfter importing the prior information into the system, the initial target probability distribution matrixPTPMk=0 was generated using the proposed MTPM method and the traditional TPM method, as shown in Figure 3.The initial distribution of target probability reflected a target’s movement in the task area at the beginning of the task. The higher the target probability was, the more likely the target was to move around the position. In Figures3(a) and 3(b), the results of the MTPM method are presented. In Figures 3(c) and 3(d), the results of the TPM method are illustrated. Figures 3(a) and 3(c) show the 3D images, and Figures 3(b) and 3(d) display the 2D images. Further, in Figures 3(c) and 3(d), the high-probability area of the target in the distribution map obtained by the traditional TPM method is dotted. Compared with the real motion trajectory of the target presented in Figure 2, the high-probability area of the target is discontinuous and has poor directivity, which cannot predict the target trajectory. Further, as shown in Figures 3(a) and 3(b), the distribution map obtained by the MTPM method showed a banded distribution of the high-probability region of the target. Compared with the actual motion trajectory of the target, which is displayed in Figure 2, the trajectory presented in Figure 3 could accurately reflect the motion trend and process of the target, indicating that the proposed MTPM method is more suitable for the motion trajectory prediction and fitting of the curve trajectory target than the traditional TPM method.After the mission began, the UAV planned the search path according to the initial distribution of the target probability and reduced the target probability value of the explored area as the system iterated, updating the target probability distribution map. Figures4(a) and 4(b) show the drone’s search path. Figure 4(a) shows the results of the MTPM method. Further, Figure 4(b) shows the results of the traditional TPM method.As shown in Figure4(a), under the combined action of the distributed ACO algorithm and the APF algorithm, after the UAV effectively searched the high-probability area of the target, it captured all targets and successfully returned to the predefined return point, and the system iterations ended at k=1639. The results in Figure 4(b) demonstrate that the UAV searched the high-probability area of the target, but due to the lack of effective target trajectory prediction, the drown cannot capture all targets, and the system iteration was ended at k=1686. The comparison of the results showed that the initial distribution of the target probability generated by the MTPM method could accurately describe the motion trend of an unknown target. The joint distributed ACO–APF algorithm could correctly guide the UAV to search the high-probability area of the target effectively and jump out of the path deadlock.The unknown area search rate and the target probability distribution search rate were calculated by Equations (16) and (17).
(16)rcovk=Scovk/SLx×Ly,(17)rTPMk=1−∑b=1Ly∑a=1LxPTPMka,b∑b=1Ly∑a=1LxPTPMk=0a,b.Figure5(a) reflects the comparison of the two methods on the unknown region search rate which reflects the extent of the drone’s search for unknown region, and Figure 5(b) reflects the comparison of the two methods on the TPM search rate which reflects the extent of the drone’s search for the high-probability distribution region of targets. With the iteration of the system, the drone’s search for unknown region gradually deepens based on the TPM matrix. In Figure 5, the end point represents the number of system iterations at the end of the task.Figure6 reflects the comparison of the two methods on the number of target capture which reflects the task efficiency of the drone. With the iteration of the system, the drone continuously discovers and captures the targets when searching the unknown region. In Figure 6, the end point represents the number of system iterations at the end of the task.As shown in Figures5 and 6, the MTPM method completed the iteration process at k=1639, achieving a search rate of an unknown region of 44.7%, a search rate of the target probability distribution of 97.3%, and the number of captured targets of 10, which were better than those of the traditional TPM method. The three mentioned metrics of the traditional TPM method were 41.4%, 88.5%, and 9 at k=1686, indicating that the effective trajectory prediction of the moving target by the MTPM method could accurately guide a UAV to search the path planning and achieve a better task execution effect.Figure7 reflects the comparison of the two methods on the distance to obstacle which reflects the obstacle avoidance situation of the drone. The vertical axis in Figure 5 represents the minimum distance from the drone to the obstacle area. The black solid line that changes with time corresponds to the data of the MTPM method, and the black dotted line corresponds to the data of the TPM method. The constant bold black solid line represents the maximum repulsive force distance Rrep(Rrep=2km) in the APF algorithm. If the distance from the UAV to the obstacle area exceeds Rrep, the repulsive force from the obstacle area would be no longer applied to the UAV. The constant bold black dotted line represents the safety distance Rsafe(Rsafe=0.8km) of the UAV in the APF algorithm. If the distance from the UAV to the obstacle area is less than Rsafe, a collision would be happened between the UAV and the obstacle area.As shown in Figure7, the minimum distance between the UAV and the obstacle area is always longer than the safe distance Rsafe, which shows that the obstacle repulsion potential field can avoid collisions between the UAV and the obstacle area. The reason is that the modulus of the repulsive force is inversely proportional to the square of the distance, so that the UAV can avoid obstacles at the right time and position.Figure8 reflects the comparison results in 30 simulations of the two methods. In Figure 8, the black color corresponds to the data of the MTPM method, and the blue color corresponds to the data of the TPM method.In Figure8(a), the vertical axis represents the number of targets captured by the drone. The dotted lines represent to the number of captured targets in 30 simulations, and the solid lines represent to the average number of captured targets in 30 simulations. After the subsimulation, the average number of captured targets of the MTPM method was 9.63, and the average number of captured targets of the TPM method was 8.4.In Figure8(b), the vertical axis represents the unknown region search rate. The dotted lines represent to the unknown region search rate in 30 simulations, and the solid lines represent to the average unknown region search rate in 30 simulations. After the subsimulation, the average unknown region search rate of the MTPM method was 43.3%, and the average unknown region search rate of the TPM method was 39.1%.According to the above-presented evaluation results, after the subsimulation, the average target captured rate and average unknown region search rate of the MTPM method were higher than that of the traditional TPM method, and the performance was improved by 14.6% and 10.7%, respectively, indicating that the traditional TPM method is deficient in trajectory prediction of dynamic targets. In contrast, the proposed method could effectively solve the search path optimization problem of unknown dynamic targets.


## 6. Conclusion

Under the condition of an uncertain dynamic environment, this paper studies the cooperative search problem of UAV swarm. In this study, an uncertain dynamic environment is defined as an unknown region, including the multiobstacle region and dynamic target with the periodic motion of curve trajectory. In addition, a cooperative search algorithm for UAVs based on the MTPM is proposed. From the aspect of prior information processing of a battlefield environment, aiming at improving the TPM method and making full use of limited target motion state information, a directional target probability distribution map is constructed and used to predict the target motion trajectory. To generate a UAV swarm’s path, this paper combines the distributed ACO algorithm and the APF algorithm to develop the path generation method according to environmental information under the multisystem constraint. When the target information is sufficient, the APF method is used to obtain the target trajectory attraction potential field and the capture path. However, when the target information is insufficient, the ACO algorithm is combined with the MTPM method into a heuristic ACO-based probability transfer method, which is used to determine a search path. After the mission ends, the UAV can return to the predefined return point. Finally, the Monte Carlo method is used to perform simulation experiments, and five indicators are used for evaluation, namely, the unknown region search rate, the TPM search rate variations, the number of captured targets, the UAV obstacle avoidance result, and the comparison results of targets captured in 30 simulations using different methods. The experimental results verify that the proposed search path planning method based on the improved target probability graph can achieve high efficiency under the condition of limited target motion information. The average target captured rate and average unknown region search rate of the MTPM method were higher than that of the traditional TPM method, and the performance was improved by 14.6% and 10.7%, respectively, indicating that the proposed algorithm can effectively solve the problem of UAV cooperative search task planning in uncertain dynamic environments. At present, there are some shortcomings in the design method of the UAV flight model. We are reading relevant literature on the mobile robots and quadrotors control algorithm [26–30], hoping to improve our work by combining these technologies in the next step. Future research will focus on how to improve the UAV flight model and smooth its path. In addition, a new UAV swam path planning algorithm will be designed for the improved UAV flight model to generate a better multi-UAV cooperative search path that can better meet actual needs.



---

*Source: 8561245-2022-12-12.xml*