A bionic intelligent method combining evolutionary game theory with particle swarm optimization for UAV 3D path planning

1 Introduction

Due to its extensive real-world applications, such as rescue, agriculture and military search or attack, unmanned aerial vehicle (UAV) has recently attracted significant attention in robotics community (Xu et al., 2024). As a fundamental and paramount component of the motion planning system, the path planning issue seeks to generate a feasible path for the UAV from a predefined start point to a destination in a complex environment. The quality of the generated path strongly affects the intelligence and automation of the motion planning system of UAV (Lyu and Yang, 2024). Consequently, addressing the UAV path planning issue via different planning methods has become increasingly prevalent in recent years (Yang et al., 2025). Nevertheless, the NP-hard nature of this issue makes it difficult to solve efficently (Tang et al., 2020). To cope with this challenge, numerous researchers have been recently devoted themselves to developing different path planning methods.

To the best knowledge of the authors, the existing path planning methods can be generally classified into two different categories. The first category focuses on traditional planning methods, such as mixed integer linear programming (MILP) (Radmanesh and Kumar, 2016), rapidly-exploring random tree (RRT) (Yang et al., 2014), spares A* (Wang et al., 2014b) and Dijkstra algorithm (Karve and Kapadia, 2020). Despite being mature of the mathematical models and relative easy implementation, the traditional planning algorithms may still struggle with the defect of easily trapping into local optimal when dealing with complicated, nonlinear and non-differential path planning issues (Lyu and Yang, 2024; Zhang et al., 2013).

The aforementioned limitations of classical methods have catalyzed the development of more robust planning approaches based on evolutionary algorithms (EAs). Generally, EAs are bionic intelligent optimization algorithms characterized by their population-based nature and capability for parallel computation (Tang et al., 2017). Owing to these advantages, different EA-based path planning methods, such as ant colony optimization (ACO) algorithm (Hu et al., 2025), genetic algorithm (GA) (Heng and Rahiman, 2025), differential evolution (DE) (Freitas et al., 2024), simulated annealing (SA) (Yu et al., 2024) and black-winged kite algorithm (BKA) (Chen et al., 2025), have been developed and treated as the second but most preferred category in the area of path planning.

As one of the most popular EAs, the particle swarm optimization (PSO) algorithm has been widely employed to address the path planning issue due to its population-based nature. easy of implementation and fast convergence speed (Gupta and Verma, 2023; Gao and Ma, 2024). Capitalizing on these advantages, many different PSO-based path planning algorithms have been recently developed. For instance, an adaptive PSO algorithm by mixing Q-learning with standard PSO has been proposed to handle the UAV path planning problem (Tan et al., 2024). Similarly, a hybrid method combing PSO with GA has been presented to solve the UAV 3D path planning problem (Li et al., 2025). Further examples of PSO applications in UAV path planning can be found in (Liu et al., 2024; Nguyen et al., 2024; Tan et al., 2023). Yet, as a stochastic algorithm, the standard PSO algorithm struggles to dynamically balance exploration and exploitation, which often results in premature convergence of the swarm (Zhang et al., 2025; Raj and Kos, 2023; Tao and Kim, 2024). Therefore, remedying this inherent flaw is essential to further enhancing the performance of PSO.

To efficiently address the UAV 3D path planning problem, this paper focuses on developing a PSO-based path planner that remedies the defects of the standard PSO noted above. To this end, a novel improved self-adaptive PSO algorithm (ISAPSO) is proposed by integrating the standard PSO 2011 (SPSO 2011) (Tang et al., 2020) with evolutionary game theory (EGT) (Leboucher et al., 2016). To achieve a well-balanced exploration-exploitation trade-off, a new self-adaptive parameter updating strategy inspired by the evolutionary stable strategy (ESS) of EGT and the hyperbolic tangent function is developed to tune the three key control parameters of the particle in ISAPSO. As a part of development, this paper completes the establishment of a ISAPSO-based path planning method, in which a modified self-adaptive constraint handling technology based on the feasibility-based rule (Wang et al., 2014a) is presented to easily deal with the constraints of the path planning problem. The performance of the proposed ISAPSO-based path planning method is evaluated under different scenarios by comparing the proposed planner with six well-known heuristic algorithms after a comprehensive benchmark study over 20 benchmark test functions. The non-parameter comparison results over the conducted benchmark study reveal that the ISAPSO algorithm performs significantly better than the six opponents compared at a 90% confidence level. Besides, the developed ISAPSO-based path planning method performs superior to its six peers in terms of the path optimality.

The remainder of this paper is organized as follows. Section 2 presents the problem statement and mathematical model of the path planning issue. Section 3 introduces the proposed ISAPSO algorithm. Section 4 details the ISAPSO-based path planning method, including particle encoding and the modified self-adaptive constraint handling technology. The numerical simulation results and analysis regarding the benchmark study and 3D path planning in different scenarios are presented in Section 5. Section 6 concludes this study and outlines potential directions of future work.

2 Description and modeling of UAV path planning

2.1 Modeling of the terrain environment

As displayed in Figure 1, the primary goal of of the path planning problem is to generate a collision-free path connecting the start and destination points for a UAV in a 3D environment. Modeling of the 3D workspace is crucial since it directly impacts the quality of the generated path. In this work, the mountain model is applied to represent the terrain information of the 3D environment due to its simplicity and universality as following (Li et al., 2025):

z(x,y)=∑i=1Nziexp[-(x-xixsi)2-(y-yiysi)2](1)

where N denotes the total number of mountains in the workspace. x_*i* and y_*i* stand for the coordinates of the ith mountain along the x-axis and y-axis, respectively. z_*i* represents the peak of the ith mountain. z(x, y) indicates the height of the ith mountain. x_*si* and y_*si* are slopes of the ith mountain along the x-axis and y-axis, respectively.

2.2 Mathematical model of path planning

To generate a path similar to that shown in Figure 1, the UAV path is denoted by ordered coordinates in the workspace. To ensure the smoothness of the generated path, this study applies the B-spline curve to model the path due to its simplicity and flexibility in characterizing the path (Zhang et al., 2024; Huang et al., 2024). Let (S, w₁, w₂, w₃, ..., w_*n*, D) represent the waypoints along the path, where S = (x₀, y₀, z₀), D = (x_*n*+1, y_*n*+1, z_*n*+1) and w_*i* = (x_*i*, y_*i*, z_*i*) stand for the start point, destination point and intermediate waypoint, respectively. Then, the path points along the UAV path can be produced through the B-spline curve as follows:

{xi′=∑k=1nBk,d(τi)xkyi′=∑k=1nBk,d(τi)ykzi′=∑k=1nBk,d(τi)zk(2)

where d is regarded as the measure metric of the smoothness of the curve. Since the smoothness and computation complexity of the curve are increased with the growing of d, d is empirically set to be 3 in this study by compromising the path smoothness and computation complexity. B_{*k, d*}(τ_*i*) represents the blending function and can be determined by the following knot vector:

vk={0,if k≤dk-d,if d<k≤n-1n-d,if n≤d(3)

After determining the knot vector through the above model, the blending function can be mathematically obtained as follows:

Bk,1(τi)={1,if vk≤τi<vk+10,otherwise(4)

Bk,d(τi)=τi-vkvk+d-1-vkBk,d-1(τi)+vk+d-τivk+d-vk+1Bk+1,d-1(τi)(5)

where τ_*i* denotes a parameter incrementally increased with a fixed step size within the range of [0, n + d]. n is the total number of points along the path. d is the order of the B-spline curve. Recall that d = 3 in order to simultaneously consider the path smoothness and computation complexity in this paper.

After representing the smooth path using the cubic B-spline curve depicted above, the path planning problem is formulated as a constrained optimization issue. Since the path length correlates with the UAV's overall fuel consumption in real-world applications, the objective function of the concerned path planning issue is then modeled to minimize the path length. Moreover, the constraints of the path planning issue are determined by the physical features of the UAV and mission execution constraint, which strongly influence the mission execution performance of the UAV. Note that the physical constraints of the UAV mainly include the flight distance and turning angle constraints, respectively. The mission conduction constraint requires the flight safety. Thus, the path planning issue can be mathematically represented as follows:

{Generating: ph=[S=w0,w1,...,wn,D=wn+1]Minimize: f=∑i=0n∥wi+1-wi∥(6)

Subject to:

{f≤Lmax(7)ψi≤ψimax, 1≤i≤n−1(8)θi≤θmax, 1≤i≤n−1(9)zi,interp(xi,yi)≤zi≤zmax, 1≤i≤n(10)

where, it has:

f=∑i=0n(xi+1-xi)2+(yi+1-yi)2+(zi+1-zi)2(11)

ψi=arctan(yi+2-yi+1xi+2-xi+1)-arctan(yi+1-yixi+1-xi)(12)

θi=arctan(zi+2-zi+1(xi+2-xi+1)2+(yi+2-yi+1)2) -arctan(zi+1-zi(xi+1-xi)2+(yi+1-yi)2)(13)

where S and D are the start and destination points of the UAV, respectively. f is the objective function denoted by the path length. ||w_*i*+1−w_*i*|| represents the Euclidean distance between waypoint w_*i*+1 and waypoint w_*i*. n is the total number of waypoints of the generated path. L^*max* is the maximum flight path length of the UAV. ψ_*i* and ψ^*max* are the ith turning angle and maximum turning angle of the UAV, respectively. θ_*i* and θ^*max* are the ith climbing slope angle and maximum climbing slope angle of the UAV, respectively. z_{i, interp(_*x**i*, y_*i*)} stands for the terrain height corresponding to (x_*i*, y_*i*). z_*i* represents the actual flight height at the ith waypoint. z^*max* is the flight range constraint along z-axis.

3 The proposed ISAPSO

Inspired by bird flocking and fish schooling, the standard PSO has been developed since 1995 to handle different global optimization issues. Particularly, PSO is well-suited for solving the UAV 3D path planning problem due to its distinct advantages: (1) it operates directly in continuous space, making it naturally compatible with the continuous coordinate representation of 3D paths; (2) its population-based architecture facilitates parallel computation, thereby significantly reducing computational time; (3) as a stochastic evolutionary algorithm, PSO possesses strong global search powers, making it effective in handling the complexity and non-linearity inherent in the NP-hard path planning problem. However, extensive real-world applications have demonstrated that the performance of the standard PSO has been limited since it struggles to balance exploration and exploitation and thus easily fall into iterative stagnation (Gupta and Verma, 2023; Gao and Ma, 2024). To remedy this deficiency, significant efforts have been made to develop different PSO algorithms. Among those enhanced PSO algorithms, SPSO 2011 would be one of the most preferred version (Tang et al., 2020; Liu and Han, 2017; Zambrano-Bigiarini et al., 2013), which is detailed in the herein contents.

3.1 Review of SPSO 2011

SPSO 2011 is a representative bionic intelligent algorithm, where each agent is regarded as a particle characterized by position and velocity vectors. As displayed in Figure 2, each particle dynamically adjusts its trajectory based on its own flight experience and those of the swarm, converging toward the global optimum while navigating around obstacles within the search space. To overcome the flaw of the standard PSO noted above, SPSO 2011 incorporates random disturbances into the velocity vector. This mechanism ensures that particles maintain non-zero velocities, thereby preserving their exploratory capability and effectively preventing premature convergence (Tang et al., 2020; Liu and Han, 2017; Zambrano-Bigiarini et al., 2013). Due to its noted bionic nature and these advantages, SPSO 2011 is applied in this study to handle the UAV path planning issue.

Let Xmt, Pbestmt, and Gbest^*t* denote the current position, the personal best position of the mth particle and the global best position of the swarm at the tth iteration in SPSO 2011, respectively. As illustrated in Figure 3, an iso-barycenter of gravity denoted as Gmt is first defined around the current position Xmt, the personal best position Pbestmt and the global best position Gbest^*t* as follows:

Gmt=Xmt+c1(Pbestmt-Xmt)+c2(Gbestt-Xmt)3(14)

where c₁ and c₂ are the cognitive and social acceleration parameters of the particle, respectively.

Then, a random position X′mt is defined in the hypersphere denoted by H(G_*m*(t), ||G_*m*(t)−X_*m*(t)||) and added to the velocity of each particle as follows:

Vmt+1=ωVmt+X′mt-Xmt(15)

where ω represents the inertia weight parameter of the particle.

Mixing the newly-disturbed velocity defined by Equation 11 with its previous flight experience, each particle updates its position information as follows:

Xmt+1=Vmt+1+Xmt(16)

where Xmt+1 is the position of the mth particle at iteration t+1. In SPSO 2011, the three main control parameters are fixed as constants as follows (Liu and Han, 2017):

{ω=12ln(2)≈0.7213c1=c2=0.5+ln(2)≈1.1931(17)

It is notable from Equation 11 that since a random position is defined in the hypersphere H(G_*m*(t), ||G_*m*(t)−X_*m*(t)||) and added to the velocity of each particle, the particle could be encouraged to keep searching in the solution space with a non-null velocity. This mechanism may prevent SPSO 2011 from prematurely plugging into iterative stagnation, which would promote the optimization performance of this algorithm to some degree. Nevertheless, SPSO 2011 may still suffer from the drawback of unpromising ability in adjusting its exploration and exploitation capabilities since the three main control parameters of the particle are constant, as well as there exists no difference between the cognitive acceleration parameter and the social acceleration parameter (Tang et al., 2020). To overcome the mentioned inherent flaw of SPSO 2011, a newly-developed parameter updating rule is proposed to tune the three key control parameters of the particle in SPSO 2011 by leveraging the EGT. The statements of EGT and novel parameter updating rule are depicted in the herein contents of this paper.

3.2 Statement of EGT

EGT is a bionic computational paradigm inspired by biological competition, where the final outcome is determined by the concept of ESS (Leboucher et al., 2016). In EGT, a strategy is treated as an ESS when all individuals in a swarm follow a specific strategy, no other strategies can invade the swarm under the natural selection. Suppose there are a few individuals or players adopting p₁ and some other individuals using strategy p₂. let f(p1,Sd¯) represent the payoff, that is, the fitness value of a specific player following strategy p₁. For a minimization optimization problem, p₁ is considered as an ESS if f(p2,Sd¯)>f(p1,Sd¯) for p₂ ≠ p₁ and any sufficiently small μ (Leboucher et al., 2016). Here, Sd¯ represents the strategy distribution of the swarm and is computed by:

Sd¯=(1-μ)p1+μp2(18)

where μ is a sufficiently small parameter, representing the frequency of p₂ in the swarm.

As a common way to describe the strategy interactions among players, the matrix game is stated using the following notations: (1) e_*l* stands for the pure strategy l, where l = 1, 2, ..., M and M denotes the total number of pure strategies; (2) K = f(e_*l*, e_*h*) is the M × M payoff matrix; (3) ▽M≡{p=(p1,...,pM)|∑l=1mpl=1,0≤pl≤1} is the set of mixed strategies and represents the probability distribution over the pure strategy e_*l*; (4) f(p, q) = pAq^*T* is the payoff, namely, the fitness value of players following the pth strategy facing other individuals using strategy q. In order to depict the difference between the payoff of strategy l and the average payoff of the swarm, the replicator dynamics equation (RDE) in EGT is defined as follows:

ṗl=-pl(elKpT-pKpT)(19)

where p_*l* denotes the probability distribution associated with the pure strategy e_*l*. K is the payoff matrix.

3.3 Proposed Self-adaptive parameter updating rule in ISAPSO

For enhancing the performance of PSO, this paper proposes a novel ISAPSO algorithm based on SPSO 2011 and EGT. The main motivation of the development concerning to the proposed algorithm is to remedy the inherent defect of SPSO 2011 noted above. To achieve this goal, a novel self-adaptive parameter updating strategy inspired by the ESS of EGT and the hyperbolic tangent function is developed to dynamically tune the three aforementioned main control parameters in ISAPSO. To prevent the swarm from early iterative stagnation, particles in the proposed algorithm first stick to the moving rules defined by Equations 10–12 in SPSO 2011 to update their flight experience. Then, a self-adaptive parameter updating strategy is proposed to adjust the three main control parameters of each particle in the proposed algorithm by leveraging the ESS of EGT and the hyperbolic tangent function.

Prior to introducing the proposed parameter updating strategy, the analogy between EGT and ISAPSO is first described as follows: (1) the players in EGT analogize particles in ISAPSO; (2) each particle in ISAPSO adopts three candidate strategies, these are, moving only according to its inertia weight, just following its personal best memory and merely following the global best memory of the swarm, respectively; and (3) the payoff matrix of EGT is consisted by the average performance obtained by each particle in ISAPSO following a specific strategy.

Let e₁, e₂, and e₃ stand for the three aforementioned strategies, respectively. Then, the payoff matrix applied in this study is computed as follows:

K={K(e1)K(e1)-K(e2)2K(e1)-K(e3)2K(e2)-K(e1)2K(e2)K(e2)-K(e3)2K(e3)-K(e1)2K(e3)-K(e2)2K(e3)}(20)

where K(e_*l*) (l = 1, 2, 3) is the payoff obtained in the case where the particle only adopts the lth strategy.

In this study, the ESS indicates the ratio of each strategy when the swarm reaches a stable point and is mathematically defined as follows:

Ess(t)=[Z1(t),Z2(t),Z3(t)](21)

Z1(t)+Z2(t)+Z3(t)=1(22)

At each iteration, the payoff K(e_*l*) (l = 1, 2, 3) of each particle is gained according to its previous flight experience as follows:

K(el)=∑t1=1t-1Zl(t1)f(x(t1))t(23)

where t is the iteration number at the current stage. f(x(t₁)) is the fitness value of particle x at the previous iteration t₁.

Once the payoff K(e_*l*) (l = 1, 2, 3) of each particle is obtained via Equation 19, the payoff matrix is then constructed by Equation 16. Subsequently, the RDE given in Equation 15 is implemented to compute the ESS, that is, E_*ss*(t) = [Z₁(t), Z₂(t), Z₃(t)]. By integrating the three ratios Z₁(t), Z₂(t) and Z₃(t) contained in E_*ss*(t) with the hyperbolic tangent function, the newly-developed parameter updating strategy in ISAPSO is defined as follows:

ωmt+1=ωs-(ωs-ωf)·tanh(δωtα)(24)

c1mt+1=c1s-(c1s-c1f)·tanh(δc1tα)(25)

c2mt+1=c2s-(c2s-c2f)·tanh(-δc2tα)(26)

where:

δω=ωs-ωftmax(27)

δc1=c1s-c1ftmax(28)

δc2=c2s-c2ftmax(29)

α=Z1(t)+Z2(t)Z3(t)(30)

tanh(x)=ex-e-xex+e-x(31)

where ωmt+1, c1mt+1, and c2mt+1 are the inertial weight, cognitive and social acceleration parameters of mth particle at iteration t+1, respectively. t_*max* is the maximum iteration number. ω_*s* and ω_*f* stand for the initial and final values of the inertial weight parameter, respectively. c_1*s* and c_1*f* are the initial and final values of the cognitive acceleration parameter, respectively. c_2*s* and c_2*f* are the initial and final values of the social acceleration parameter, respectively. tanh(x) indicates the hyperbolic tangent function.

3.4 Analysis of the proposed self-adaptive parameter updating rule

It is apparent from Equations 20–27 that the inertial weight parameter and cognitive acceleration parameter decrease nonlinearly, whereas the the social acceleration parameter of the particle increases nonlinearly with the iteration number increasing. This parameter variation trend indicates that the inertial weight and cognitive acceleration parameter remain relatively large in the early evolution stage, while the social acceleration parameter becomes relatively big in the latter stage of the evolution. Since big values of inertial weight and cognitive acceleration parameter promotes the exploration capability (i.e., the global search power) of the particle and greater social acceleration parameter enhances the exploitation ability (i.e., the local search capability) of the particle (Liu et al., 2024; Nguyen et al., 2024; Tan et al., 2023). Thus, adopting the proposed self-adaptive parameter updating strategy, the global search capability of ISAPSO could be naturally strengthened in the early evolution and its local search capability could be more retained in the late evolution. This natural feature regarding to ISAPSO on one hand could encourage the swarm convergent toward solution space nearby the optimal as soon as possible in the early evolutionary process. On the other hand, the swarm could be encouraged to perform an intensive search in the local solution space containing the optimal solution, thereby enhancing the likelihood of finding the optimal solution.

Apart from the iteration number t, the trade-off between the global and local searche abilities of ISAPSO is also modulated by the parameter α. It can be observed from Equations 20–27 that decreasing rates of the inertial weight parameter and cognitive acceleration parameter diminish, while the increasing tendency of the social acceleration parameter also decreases with α increasing. This may imply that the global search capability of ISAPSO would be more maintained with a larger α. Note that a relatively larger α denotes a bigger value of Z₁(t) + Z₂(t) based on Equation 26. Since a bigger value of Z₁(t) + Z₂(t) in ESS indicates that the search direction of swarm is more stable when most particles adopt the strategies of following their inertia and personal flight experience, it is reasonable to promote the global search capability of ISAPSO for a relatively bigger α. Contrarily, it is natural to promote the local search ability of ISAPSO for a relatively smaller α since a smaller α indicates a bigger value of Z₃(t), which implies the search performance of the proposed algorithm is more stable when the majority of the particles play the strategy based on the global best flight memory of the swarm.

Moreover, since the concerned path planning problem is a nonlinear NP-hard issue, adaptively tuning the three key parameters via the developed parameter updating rule may enhance the adaptability of proposed algorithm over the studied path planning problem. Besides, it is well known that the hyperbolic tangent function is a typical saturated activation function with promising convergence speed due to its exponential function characteristic. Employing the hyperbolic tangent manner to update the three key control parameters of the particle may further promote the convergence speed of ISAPSO. Briefly, on one hand, the iterative stagnation issue regarding to ISAPSO would be remedied in the early evolution since this algorithm adopts the same velocity disturbance mechanism introduced in SPSO 2011. On the other hand, the global and local search capabilities of ISAPSO would be well balanced through the developed self-adaptive parameter updating strategy defined by Equations 20–27. All these improvements may promote the optimization performance and convergence speed of the proposed ISPASO over complicated nonlinear optimization issues. The variation tendencies regarding to the three control parameters through the developed self-adaptive updating rule under different values of α are displayed in Figures 4–6 in the case where ω_*s* = 1, ω_*f* = 0.1, c_1*s* = c_2*f* = 2 and c_1*f* = c_2*s* = 0.2.

3.5 Complexity analysis of ISAPSO

This subsection presents a theoretical analysis of the computational complexity of the proposed ISAPSO algorithm. The computational complexity primarily arises from four sources: (1) initializing the position and velocity of each particle; (2) updating the velocity and position vectors in each iteration; (3) updating the personal best and global best solutions; and (4) adaptively updating the three control parameters via the proposed EGT-based strategy.

Let NP, D, and T_*max* denote the swarm size, problem dimensionality and maximum number of iterations, respectively. The initialization phase requires O(NP·D) operations. In each iteration, updating the velocity and position vectors incurs a complexity of O(NP·D). Additionally, the update of personal best and global best positions contributes O(NP) per iteration. Since the proposed parameter adaptation mechanism employs EGT with M candidate strategies, the complexity of updating control parameters across the entire swarm is O(NP·M²) per iteration.

Consequently, the total computational complexity (tcc) of the proposed algorithm is expressed as:

tcc=O(NP·D)+O(Tmax·NP·D)+O(Tmax·NP)+O(Tmax·NP·M2)(32)

Given that T_*max* > 1 and M is a constant (M = 3 in this study), Equation 28 can be simplified by extracting the dominant terms:

tcc=O(Tmax·NP·D)+O(Tmax·NP·M2)=O(Tmax·NP·(D+M2)(33)

It is evident from Equation 29 that the computational complexity of the proposed ISAPSO is primarily governed by the swarm size, the maximum iteration number and the problem dimensionality. Furthermore, the incorporation of EGT for parameter adaptation introduces additionally computational overhead due to the calculation of the payoff matrix and strategy evolution, which is reflected in the O(NP·M²) term.

4 The ISAPSO-based path planner

This section completes the design of the ISAPSO-based path planning method utilizing the algorithm proposed in Section 3. During this process, particle encoding and constraint handling are two key issues needed to be addressed. These two key technologies are depicted in the following contents of this section.

4.1 Encoding particle

Particle encoding aims to establish a mapping between the particle's position vector and a candidate solution to the optimization problem. Since this study focuses on 3D path planning, the position vector of each particle is encoded as a 3 × n vector, where n denotes the total number of waypoints. The column number of the position vector corresponds to the number of the waypoint. The three position coordinates of each waypoint are orderly encoded based on each column of the corresponding particle's position vector. For guaranteeing each particle to search within the given 3D workspace of UAV, the following saturation method is used to modify each element contained in the position vector of the particle at each iteration as:

xi={xmax, if xi>xmaxxmin, if xi<xminxi, otherwise(34)

yi={ymax, if yi>ymaxymin, if yi<yminyi, otherwise(35)

zi={zmax, if zi>zmaxzmin, if zi<zminzi, otherwise(36)

where x_*i*, y_*i*, and z_*i* (i = 1, 2, .., n) are the three position coordinates of the ith waypoint along the genrated path, that is, the three elements in the ith column of the position vector of the particle. x^*max*, y^*max*, and z^*max* are the maximum values of three axes along the 3D workspace system of the UAV, respectively. x^*min*, y^*min*, and z^*min* are the minimum values of three axes along the 3D workspace system of the UAV, respectively.

4.2 Handling constraints

As defined by Equations 7–10, since the path planning problem is subject to physical and flight safety constraints of the UAV, handling these constraints remains the second paramount issue in the ISAPSO-based path planner to guarantee the feasibility of the generated path. To address this issue, a newly-developed constraint handling technology based on the feasibility-based rule is proposed in this paper. In this technology, the violation degree of each particle at each iteration is first computed as follows:

CVmt=Clmt+Cψmt+Cθmt+Csmt(37)

where CVmt is the total violation degree of the mth particle at iteration t. Clmt,Cψmt, Cθmt, and Csmt are the violation degrees of the flight distance, turning angle, climbing slope angle and the path safety constraints, respectively.

Given the maximum physical constraints of the UAV and the terrain information of the 3D workspace, Clmt,Cψmt, Cθmt, and Csmt can be calculated as follows:

Clmt={0if Lmt≤LmaxLmt-LmaxLmaxotherwise(38)

where Lmt is the path length of particle m at iteration t. L^*max* is a predefined constant, denoting the maximum flight distance of the UAV.

Cψmt=1n∑i=1nAimAim={0if |ψimt|≤ψmax1otherwise(39)

where n is the total number of waypoints along the generated path. ψimt represents the ith turning angle of the path searched by the mth particle at iteration t, which can be computed via Equation 8. ψ^*max* is the maximum turning angle of the UAV.

Cθmt=1n∑i=1nBimBim={0if |θimt|≤θmax1otherwise(40)

where θimt indicates the ith climbing slope angle of the path found by the mth particle at iteration t, which can be obtained by Equation 9. θ^*max* is the given maximum climbing slope angle of the UAV.

Csmt=1n∑i=1nCimCim={0if zimt>zi,interp(xi,yi)1otherwise(41)

where zimt is the height of the ith waypoint along the path yielded by the mth particle at iteration t. z_{i, interp(_*x**i*, y_*i*)} is the terrain height corresponding to point (x_*i*, y_*i*) and can be gained by the interpolation calculation according to the terrain model defined by Equation 1. (x_*i*, y_*i*) is the position coordinates of the ith waypoint of the path searched by the mth particle at iteration t.

After calculating the total violation degree of the particle using Equations 33–37, an adaptive relaxation variation method is introduced within the proposed constraint-handling scheme of the ISAPSO-based path planner, so as to preserve population diversity. To implement this strategy, a relaxation variable is first defined and computed as follows:

CVrt=CV0·(1-NftNP)·exp(-t)(42)

where CVrt is the value of the relaxation variable at iteration t. Nft is the total number of feasible particles at iteration t. NP is the total number of particles in the ISAPSO-based planner. CV₀ is the average violation constraint degree regarding to the particle swarm at the initial iteration. After initializing the swarm and evaluating each particle's total constraint violation degree via Equations 33–37, CV₀ can be computed as follows:

CV0=1NP∑m=1NPCVm0(43)

Figure 7 displays the variation curves of the relaxation variable in the case where NftNP=0.3,0.5,0.7 and CV₀ = 8, respectively. After iteratively computing the relaxation variable following Equations 38, 39, each particle at each iteration is considered and defined as a feasible particle (or solution) in the case where CVmt≤CVrt. Followafter redefining the feasibility of the particle at each iteration, the following feasibility-based rule is applied in the developed constraint handling technology to select the elite solution between any two candidate solutions in the ISAPSO-based planner as: (1) for any two feasible solutions, the solution owing the smaller fitness value dominates the one having greater fitness value; (2) for any feasible solution and any infeasible solution, the feasible solution dominates the infeasible solution; (3) for any two infeasible solutions, the solution with smaller violation degree dominates the solution with greater violation degree.

It is worth noting that the developed self-adaptive constraint handling technology evaluates the fitness values and constraint violation degrees of any two candidate solutions independently, without introducing additional control parameters. Thus, the computational burden of the ISAPSO-based planner in managing constraints could be effectively reduced. Besides, since some real infeasible solutions may still contain valuable information about the solution space, permitting them to advance to the next iteration may help to diversify the swarm, thereby increasing the likelihood of locating the global optimum.

To achieve the aforementioned objective, a self-adaptive variation strategy defined by Equations 38, 39 is introduced within the proposed constraint handling technology to dynamically redefine the feasibility of each particle. Adopting this strategy, some real infeasible solutions could be reclassified as feasible and permitted into next generations. This would naturally contribute to maintaining population diversity and enhance the capability of generating high-quality paths. Moreover, as evidenced by Equations 38, 39 and Figure 7, the relaxation variable decreases nonlinearly with the iteration number increasing and converges to zero in the latter stage of the evolution. This may indicate that the feasibility criteria become progressively stricter. Therefore, infeasible solutions are more likely to be retained during the early evolutionary phase to promote diversity, whereas most are discarded in the later stages to ensure that the final path strictly satisfies the UAV's physical and safety constraints. To summarize, the proposed self-adaptive constraint handling method may not only address constraints without imposing additional computation burden on the optimizer, but also diversify the swarm, therefore enhancing the probability of obtaining high-quality path for the UAV.

4.3 The framework of the ISAPSO-based path planner

This subsection outlines the overall framework of the ISAPSO-based path planning method for solving the path planning problem defined by Equations 6–10. Let NP and t_*max* denote the swarm size and maximum iteration number of ISAPSO, respectively. The algorithmic steps of ISAPSO-based path planner are summarized and shown in Table 1. For a minimization problem, the ways of updating the personal best and global best solutions at the following algorithmic steps can be mathematically given as follows:

Pbestmt={xmt,if f(xmt)≼f(Pbestm(t-1))Pbestm(t-1), otherwise(44)

Gbestt={Pbestmt,if f(Pbestmt)≼f(Gbest(t-1))Gbest(t-1), otherwise(45)

where xmt is the position vector of particle m at iteration t. Pbestmt denotes the personal best solution of particle m at iteration t. Gbest^*t* is the global best solution of the swarm at iteration t. f(·) is the fitness value corresponding to different position vectors. f(xmt)≼f(Pbestm(t-1)) indicates that solution xmt dominates solution Pbestm(t-1) according to the developed constraint handling technology.

Table 1. The algorithmic steps of the ISAPSO-based path planning method.

5 Numerical simulation and result analysis

To verify the proposed method, three distinct sets of numerical simulations are conducted in this paper. The first numerical simulation employs ablation experiments to validate the efficacy of the velocity perturbation mechanism in SPSO 2011 and the self-adaptive parameter tuning strategy in the proposed ISAPSO. The second set focuses on evaluating the performance of the ISAPSO algorithm across 20 benchmark test functions. Following the benchmark tests, the ISAPSO-based path planning method is further assessed under two different 3D path planning scenarios. All the conducted simulation studies are described in detail in the following sections.

5.1 Ablation experiment

The purpose of this ablation study is to simply verify the effectiveness of the random velocity perturbation mechanism in SPSO 2011 and the adaptive parameter adjustment mechanism in the proposed ISAPSO. It is worth noting that the primary distinction between SPSO 2011 and the standard PSO lies in the introduction of the random velocity perturbation mechanism. Furthermore, the key difference between ISAPSO and SPSO 2011 is that ISAPSO incorporates an adaptive parameter tuning law into the SPSO 2011 framework. Therefore, in the conducted ablation experiments, the performance of the standard PSO, SPSO 2011 and ISAPSO is comparatively evaluated through numerical simulations on two benchmark functions. The characteristics of these benchmark functions are summarized in Table 2.

Table 2. Information of the benchmarks used in the ablation study (“D” stands for the dimension).

In the conducted ablation experiments, each algorithm is independently run 10 times on each benchmark function. For each run, the population size and the maximum iteration number of each algorithm are set to 50 and 500, respectively. After completing the runs for each benchmark function, the average solution accuracy, denoted by E_*mean*, is adopted as the metric to evaluate the performance of the considered algorithms. To ensure the convergence rate of the proposed method, the initial and final values of the control parameters in the proposed self-adaptive parameter updating strategy are empirically set as ω_*s* = 0.9, ω_*f* = 0.3, c_1*s* = c_2*f* = 2 and c_1*f* = c_2*s* = 0.2.

The statistical results of E_mean obtained by each algorithm on each benchmark function in the conducted ablation study are given in Table 3. The average fitness value of each algorithm on each benchmark function is shown in Figure 8. As shown in Table Table 3, the results of the ablation experiments shown in Table 3 indicate that the performance ranking of the three algorithms is ISAPSO, SPSO 2011 and Standard PSO in terms of average solution quality. It is worth noting that the primary distinction between SPSO 2011 and Standard PSO is the introduction of a velocity perturbation mechanism designed to prevent premature convergence of the swarm. Furthermore, the key difference between ISAPSO and SPSO 2011 lies in the incorporation of the proposed self-adaptive parameter updating strategy. Therefore, the data presented in Table 3 can reflect the efficacy of the velocity perturbation mechanism in SPSO 2011 and the parameter updating strategy proposed in ISAPSO to some extent.

Table 3. Statistical results of E_mean obtained by each algorithm on each benchmark function in the ablation experiment.

The bold values provided in these tables denote the best value obtained by each method of the corresponding item.

5.2 Benchmark study

Following the ablation study to simply verify the effectiveness of the proposed parameter updating strategy in ISAPSO, the performance of this algorithm is further evaluated using 20 benchmark test functions, as listed in Table 4. For a rigorous comparison, ISAPSO is benchmarked against six well-known evolutionary algorithms: IPSO (Raj and Kos, 2023), SPSO 2011 (Liu and Han, 2017), EGPSO (Leboucher et al., 2016), GWO (Yang et al., 2025), FDE (Gupta et al., 2025) and quantum technology-based genetic algorithm (QGA) (Wu, 2024). For comprehensive comparison, a Monte-carlo experiment test with 25 runs is conducted for each algorithm on every benchmark function. In each run, the swarm size and maximum iteration number of each method are empirically set to be 50 and 1E+05, respectively. The evolutionary process terminates when either the maximum iteration number is reached or the solution error meets the predefined accuracy threshold specified in Table 4. For each benchmark test function, the solution accuracy metric denoted by E_*mean* is applied as a metric to verify the performance of each considered algorithm. Furthermore, the initial and final values of the control parameters in ISAPSO are set consistent with those shown in SubSection 5.1.

Table 4. Information of selected benchmark test functions (“D” stands for the dimension and “ξ” indicates the accuracy level).

Following the Monte Carlo simulations described above, the stochastic results regarding solution accuracy obtained by different algorithms for each benchmark function are summarized in Tables 5, 6, where the best average result of E_*mean* yielded by each method for each test function is highlighted in boldface. The average fitness curve of each method for each test function is illustrated in Appendix.

Table 5. The statistical results regarding to E_*mean* gained by each algorithm for each benchmark.

The bold values provided in these tables denote the best value obtained by each method of the corresponding item.

Table 6. The statistical results regarding to E_*mean* gained by each algorithm for each benchmark.

The bold values provided in these tables denote the best value obtained by each method of the corresponding item.

It is clear from Table 5, 6 that the proposed algorithm demonstrates superior performance over the majority of the 20 benchmark functions (specifically excluding F₄ and F₈) compared to the other six algorithms in terms of E_*mean*. Note that the proposed algorithm is ranked the second among the seven algorithms when addressing F₄ and F₈. Thus, this overall analysis can be intuitively claimed that the proposed algorithm is highly promising over the most of the 20 test functions in terms of the solution optimality.

It is important to note that although ISAPSO performs well overall, one cannot strictly claim that it universally dominates its competitors across all 20 test functions. This is primarily due to the significant variability in the E_*mean* metric among different algorithms and benchmark functions. To rigorously quantify the significant improvements offered by the proposed ISAPSO, a statistical comparison needs to be conducted. In the conducted statistical comparison, a rank-based analysis is first executed to check the mean rank value of E_*mean* quality of each algorithm over the 20 test function. Based on the stochastic results displayed in Table 5, the rank value of E_*mean* obtained by each method for each test function and the average rank value regarding to E_*mean* of each method over the 20 test functions are summarized and displayed in Table 7.

Table 7. Rank values with respect to E_*mean* obtained by different methods for the 20 benchmarks.

As shown in Table 7, the proposed ISAPSO achieves the highest average rank based on E_*mean* across the 20 benchmarks, followed by IPSO, GWO, FDE, EGPSO, QGA, and SPSO 2011. However, it is notable that this observation merely indicates the relative efficiency of ISAPSO,rather than concluding that ISAPSO performs significantly differently from or better than its counterparts. To address this, a non-parametric Friedman test (Demšar, 2006) is conducted at a 90% confidence level to evaluate the significance of the differences among the algorithms.Given that K = 7 algorithms are tested on N = 20 functions, the critical value F_*stastic* is calculated as 1.8265 using the MATLAB command finv(α, K−1, (K−1)(N−1)), where α = 0.1 represents the significance level. Based on the data in Table 7, the calculated Friedman statistic F_*s* is 80.8124, which is greater than the value of F_*stastic*. Therefore, we can claim that the seven compared algorithms perform significantly different across the 20 benchmarks at the given confidence level in terms of E_*mean* metric (Demšar, 2006).

Although the non-parametric Friedman test confirms significant differences among the seven algorithms, it does not specify which algorithms differ from each other. To identify specific pairwise improvements, a post-hoc Bonferroni–Dunn test (Demšar, 2006) is conducted at a 90% confidence level. Given that K = 7 algorithms are evaluated on N = 20 test functions, the critical difference (CD) of pairwise post hoc Bonferroni–Dunn test equals to 1.6354 based on qα(K(k+1)/(6N)), where q_α is constantly set to be 2.394 (Demšar, 2006). Based on the summarized results shown in Table 7, the absolute differences in average ranks between ISAPSO and IPSO, GWO, FDE, EGPSO, QGA, and SPSO2011 are 1.65, 1.70, 3.00, 3.75, 4.45, and 5.75, respectively. Since all these values are greater than that of the critical difference of pairwise post hoc Bonferroni–Dunn test (i.e., CD = 1.6354), it can be sufficiently conclusive that the proposed ISAPSO significantly dominates its six contenders over the 20 test functions at the confidence level of 90% in terms of the average solution optimality. Thus, the proposed ISAPSO algorithm can be efficient in solving different global optimization issues.

5.3 Verification of the ISAPSO-based path planner

After verifying the proposed ISAPSO via the benchmark study conducted in the above subsection, this subsection evaluates the performance of the ISAPSO-based path planning method on the UAV 3D path planning issue. For a rigorous comparison, the performance of the ISAPSO-based path planning approach is compared with those of the same six evolutionary algorithms depicted in SubSection 5.2. In each simulation scenario, the swarm size and maximum iteration number of each approach are empirically set to be 40 and 100, respectively. The rest simulation parameters for each approach are referred to SubSection 5.2. For a comprehensive comparison and reducing impacts of randomness, a Monte-Carlo experiment with 30 runs is conducted in each simulation scenario for each considered method. Moreover, the fitness value (FV) and the computation time (CT) of each method are reported in each numerical simulation in order to quantitatively assess the performance of each approach.

In each numerical simulation, the number of waypoints is set to be 5 to compromise the solution optimality and computation time. The physical constraints of the UAV in the two numerical simulations are displayed in Table 8. The first numerical simulation is conducted in a 100km × 100km × 100km workspace with 10 mountains in the case where the UAV is navigated from the start point (1, 1, 1) to the destination point (100, 100, 40). The second scenario simulates to navigate the UAV from the start point (10, 10, 0) to the destination point (200, 200, 100) in a 200km × 200km × 200km workspace owing 12 mountains. The terrain information of each mountain in the two conducted simulations is illustrated in Tables 9, 10.

Table 8. The physical constraints of the UAV in the two numerical simulations.

Table 9. The terrain information of each mountain in the first numerical simulation (where mn denotes mountain number).

Table 10. The terrain information of each mountain in the second numerical simulation (where mn denotes mountain number).

After conducting the described Monte-Carlo experiment for each simulation scenario, the statistic results regarding FV and CT gained by each approach are summarized in Tables 11, 12. The best generated path yielded by each method among the Monte-Carlo study for the two numerical simulations is demonstrated in Figures 9,10. The fitness curve corresponding to the best generated path for each studied case is illustrated in Figures 11, 12.

Table 11. The statistical results of FV and CT gained by approach in the first simulation.

The bold values provided in these tables denote the best value obtained by each method of the corresponding item.

Table 12. The statistical results of FV and CT gained by approach in the second simulation.

The bold values provided in these tables denote the best value obtained by each method of the corresponding item.

It is apparent from Figures 9, 10 that each considered approach can generate a collision-avoidance flyable path for the UAV in each numerical simulation. This can reflect the feasibility and efficiency of each method over the UAV 3D path planning issue to some extent. However, it is obvious from Tables Tables 11-12 that the developed ISAPSO-based path planning approach dominates the six compared approaches in terms of the average fitness value in the two numerical simulations, which can confirm that the developed approach performs generally superior to its six peers in handling UAV 3D path planning problem. Moreover, it can be also observed from Tables 11, 12 that the developed path planner can gain the the best values regarding to the best and worst fitness values among the seven considered approaches in the Monte-Carlo experiment conducted for the two numerical simulations. This could further reflect the high efficiency and performance of the developed path planning approach.

Moreover, it can be easily inferred from Tables 11, 12 that SPSO 2011 and the ISAPSO-based planning method are respectively ranked the first and second among the seven approaches in the Monte-Carlo experiment conducted for the two numerical simulations in terms of the computation time.This indicates that the computation complexity of SPSO 2011 and ISAPSO-based approach are cheaper than those of the other compared approaches. It is important to note that although SPSO 2011 is computationally cheaper than our developed path planning approach probably due to the fact that a novel self-adaptive parameter updating strategy is adopted in ISAPSO, the distinction of the computation time between SPSO 2011 and ISAPSO-based method is not outstanding, which could be neglected since the path planning issue considered in this paper is a global path planning issue, which requires to generate path off-line.

It is notable from Figures 11, 12 that the fitness curves of some approaches keep raising in the early of the evolution process and maintain dropping or stable in the latter evolution stage. This changing tendency of the fitness curve could be likely interpreted by the fact that the self-adaptive constraint handling technology is implemented in these path planning methods. Since both unfeasible and feasible solutions are allowed to participate into the next iterations in the early evolution stage via our developed constraint handling technology, the feasible solution with greater fitness value would dominate the unfeasible solution having smaller fitness value, which thus leads the fitness curves of some approaches to raising in the early evolution. In the latter evolution stage, since the redefinition of feasible solution becomes stricter though the developed constraint handling technology, only the real feasible solutions are permitted to join in the next iteration, which results in the feasible solution having smaller fitness value vanquishing the one with bigger fitness value. Probably due to this mechanism, as displayed in Figures 11, 12, the fitness curves of some considered approaches keep dropping or stable in the latter of evolution process.

To summarize the analysis described above, it can be conclusive that the developed ISAPSO-based path panning approach is highly powerful and promising in solving UAV 3D path planning problem in terms of the solution optimality. Moreover, the computation complexity of this developed approach is highly competitive and comparable with those of its peer. Thus, the developed ISAPSO-based path planning method could be regarded as a paramount alternative in the field of UAV path planning.

6 Conclusion and future work

To efficiently handle the UAV 3D path planning problem, this paper proposes a novel ISAPSO algorithm based on SPSO 2011 and EGT. For well balancing the global and local search abilities of particles, a novel self-adaptive parameter updating principle inspired by the ESS of EGT hyperbolic tangent function is proposed in ISAPSO to tune the three key control parameters of the particle. Moreover, leveraging the proposed PSO algorithm, a new ISAPSO-based path planing approach is completed in this paper. Aiming at easily handling constraints of the path planning problem and diversifying the particle swarm, a modified self-adaptive constraint handling technology according to the feasibility-based principle is adopted in the ISAPSO-based planner. Finally, the performance of the proposed ISAPSO algorithm is assessed through a benchmark study over 20 benchmark test functions. Followed the benchmark study, the ISAPSO-based path planning method is evaluated under different simulation scenarios. The simulation results regarding to the benchmark study reveal that the proposed ISAPSO algorithm performs significantly better than its peers at the confidence level of 90% over the 20 selected test functions. Moreover, the ISAPSO-based path planning method generally dominates its competitors in terms of the path optimality.

The results and method shown in this paper raise some works deserved future study. Firstly, as a stochastic algorithm, the convergence issue of the proposed PSO remains unaddressed. To fully guarantee convergence, an analytical investigation of the algorithm is needed in the near future. Secondly, the path planning issue considered in this paper belongs to the community of global path planing. The feasibility and efficiency of the developed path planning approach over online path planning problem with uncertain or moving obstacles could be considered in the near future. Next, this paper mainly devotes its interest on single UAV path planning. As an extended study of this paper, we are considering the possibility of applying the developed approach to handle multiple-UAVs cooperative path planning problem. Last but not least, since the proposed method is a heuristic algorithm with swarm-based nature, the proposed method could be considered to solve the multiple-objective path planning issue by adopting the non-dominated definition of Pareto front.

Generative AI statement

The author(s) declared that Generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnbot.2026.1857152/full#supplementary-material

Supplementary material Supplementary material at the source

Footnotes

Edited by: Zhengtai Xie, Lanzhou University, China

Reviewed by: János Botzheim, Eótvós Loránd University, Hungary

Baichun Gong, Nanjing University of Aeronautics and Astronautics, China

Keming Jiao, Tongji University, China