Unmanned aerial vehicle path planning method based on improved exponential distribution optimization

Abstract

This invention belongs to the field of UAV trajectory planning technology and discloses a UAV trajectory planning method based on improved exponential distribution optimization. The specific steps are as follows: 1. Normalize the three costs of UAV trajectory distance, UAV turning angle, and UAV distance to obstacles, and add weights to each cost to obtain a cost function; 2. The UAV trajectory consists of N+1 sub-trajectories, and determine whether each sub-trajectory intersects with the obstacle; 3. Plan several collision-free trajectories as the initial population for the exponential distribution optimization algorithm; 4. Optimize the UAV trajectory using the exponential distribution optimization algorithm; 5. Smooth the generated trajectory using a third-order B-spline curve to generate a smooth trajectory that satisfies dynamic constraints. The optimization algorithm of this invention has a fast convergence speed and can optimize the UAV trajectory in fewer iterations compared to other intelligent optimization algorithms.

Description

Technical Field

[0001] This invention belongs to the field of UAV trajectory planning technology, specifically relating to a UAV trajectory planning method based on improved exponential distribution optimization. Background Technology

[0002] Unmanned aerial vehicles (UAVs) are playing an increasingly important role in future integrated joint operations due to their characteristics such as low detectability, low cost, lack of risk of casualties, simple takeoff and landing, flexible operation, diverse system configurations, and intelligent automatic control. Current UAV trajectory algorithms can be divided into graph search algorithms, intelligent optimization algorithms, and artificial intelligence algorithms. Among them, intelligent optimization algorithms are widely used due to their wide applicability and short planned trajectory distances. However, current intelligent optimization algorithms generally suffer from slow convergence speeds, failing to optimize trajectories within the shortest number of iterations. Therefore, research on intelligent optimization algorithms for UAV trajectory planning still has significant gaps.

[0003] UAV trajectory planning is an important aspect of UAV research, and swarm intelligence optimization algorithms are one of the key methods for UAV trajectory planning. Their task is to optimize trajectory distance, turning angle, etc., within the fewest iterations to obtain a collision-free shortest trajectory. Currently common trajectory planning algorithms include the Grey Wolf Algorithm, Particle Swarm Optimization Algorithm, Ant Colony Algorithm, and Genetic Algorithm.

[0004] The current algorithm has three problems. First, the current intelligent optimization algorithm generally has a slow convergence speed. Second, the current swarm intelligence optimization algorithm generally considers a relatively simple cost function, which means that the optimization algorithm can only optimize the track distance or turning angle. Third, most of the tracks obtained by the current intelligent optimization algorithm are right-angle turns, which are difficult for UAVs to adapt to. Summary of the Invention

[0005] To address the technical problems existing in the prior art, this invention provides a UAV trajectory planning method based on improved exponential distribution optimization. In two-dimensional and three-dimensional spaces with known irregular obstacles, a novel optimization algorithm is used to optimize the generated trajectory, minimizing the cost of the trajectory in the fewest iterations. The cost function used comprehensively considers the trajectory distance, turning angle, and distance to obstacles, resulting in a shorter optimized trajectory distance, smaller turning angle, and greater distance to obstacles, ensuring safe flight of the UAV. Furthermore, the trajectory generated by optimizing it with a third-order B-spline curve makes the trajectory more suitable for UAV flight.

[0006] To achieve the above objectives, the technical solution adopted by this invention is: a UAV trajectory planning method based on improved exponential distribution optimization, the specific steps of which are as follows:

[0007] Step 1: Construct the cost function: Normalize the three costs of the drone's flight path distance l1, the drone's turning angle θ, and the distance from the drone to the obstacle l2, and add weights to the three costs to obtain the cost function;

[0008] Step 2: Construct a collision detection strategy: The drone's trajectory consists of N+1 sub-trajectories. If the distance from a sub-trajectory to the center of each obstacle is less than or equal to the radius of the obstacle, then the sub-trajectory is determined to intersect with the obstacle, and the drone's trajectory is determined to intersect with the obstacle. If the distance from each sub-trajectory to the center of each obstacle is greater than the radius of the obstacle, then all sub-trajectories do not intersect with the obstacles, and the drone's trajectory is determined to not intersect with the obstacles.

[0009] Step 3: Generate initial collision-free tracks: Plan several collision-free tracks;

[0010] Step 4: Optimize the UAV trajectory using the exponential distribution optimization algorithm;

[0011] Step 5: Use a third-order B-spline curve to smooth the generated track, resulting in a smooth track that satisfies dynamic constraints.

[0012] In two-dimensional space, the map is placed in the first quadrant of the plane coordinate system, and the approximation method is used to adjust the irregularly shaped obstacles into circular obstacles;

[0013] In three-dimensional space, the map is placed in the first quadrant of the spatial coordinate system, and an approximation method is used to adjust irregularly shaped obstacles into spherical obstacles;

[0014] Each UAV's trajectory is defined as consisting of N+1 sub-trajectories, and each UAV trajectory has N+2 endpoints. In two-dimensional space, the N+2 endpoints are {D}. s ,D0,D1,...,D N D g}, the coordinates of each endpoint are {(x s ,y s ),(x0,y0),(x1,y1)...(x N ,y N ),(x g ,y g )}, where (x s ,y s (x) is the starting coordinate, (x) g ,y g Let X be the endpoint coordinates. There are M obstacles in the space. The obstacles are approximated as circular obstacles, and the centers of the obstacles are {(X1,Y1),(X2,Y2),...,(X... M ,Y MThe radii of the obstacles are {R1, R2, ..., R}. M}

[0015] The specific steps for constructing the cost function are as follows:

[0016] 11) Calculate the track distance l1: Each track consists of N+1 sub-tracks. The track distance l1 is obtained by summing the distances of the sub-tracks. The specific algorithm for solving this problem is shown in equation (1):

[0017]

[0018] 12) Calculate the track offset angle θ: Calculate the offset angle of each subsequent sub-track segment relative to the previous sub-track segment and sum them to obtain the total offset angle θ. 012 It is the offset angle of the two sub-tracks formed by the three points {D0, D1, D2}, where the coordinates of {D0, D1, D2} are {(x0, y0), (x1, y1), (x2, y2)}, respectively. The specific calculation method is shown in Equation (2-6):

[0019]

[0020] 13) Calculate the distance l2 from the sub-track to the obstacle: If the minimum distance l from the sub-track to the obstacle is... min Larger than the radius R of the obstacle i Add safety margin l safe The distance l2 from the sub-track to the obstacle is defined as the length l of the sub-track; if the minimum distance l from the sub-track to the obstacle is... min Less than or equal to the radius R of the obstacle i Add safety margin l safe The distance l2 from the sub-track to the obstacle is determined as the maximum track length l. max Specifically, as shown in equation (12);

[0021] Calculate the distance from each sub-track to the obstacle, i.e., use a vector algorithm to calculate the distance from the center of the obstacle circle to the line segment; first, construct the vector. According to the vector Calculate vectors Specifically, as shown in equation (8); then, the value of r is calculated according to equation (9) for judgment, and the specific judgment method is shown in equation (10); equation (10) calculates the center R of the circle based on the value of r. i to line segment Distance, line segment The representation is shown in equation (7):

[0022]

[0023] The distance from the subtrack to the circular obstacle is shown in equation (10):

[0024]

[0025] l min =min{l o0 ,l o1 ...l oM} (11)

[0026]

[0027] Among them, R i Let D be the center of the obstacle circle, D0 be the beginning of the sub-track, D1 be the end of the sub-track, and D be the center of the circle. i Projection onto a line segment; in equation (7), * represents the product of a number and a variable; in equation (8), · represents the dot product of two vectors, || represents the magnitude of the vector, || 2 The vector modulus is the square of the vector; in equation (10), l oi Let || represent the distance from a certain segment of the flight path to the i-th obstacle, and || represent the absolute value of a number; in equation (12) This represents the minimum distance from the i-th sub-track to all obstacles.

[0028] The normalization of the track distance l1 and the distance to the obstacle l2, and the normalization of the track offset angle θ, are then summed to obtain the final cost function, as shown in equation (14):

[0029]

[0030] In equation (14), N is the number of intermediate endpoints, and α+β+γ=1;

[0031] Where l1 is the track distance, θ is the track offset angle, l2 is the distance from the sub-track to the obstacle, and the maximum track distance is l max .

[0032] In step three, the space is divided into N subspaces, with S as the starting point, G as the ending point, and x... i Let y be the length of the subspace along the x-axis, and y be the length of the subspace along the x-axis. i Let be the length of the subspace along the y-axis. The specific calculation method is shown in the formula:

[0033]

[0034] Based on the specific lengths of the subspaces along the x and y axes, the exact location of the i-th subspace is determined. The coordinates of the four points in the i-th subspace are:

[0035] The specific location of each subspace is determined based on the coordinates of the aforementioned subspaces, thereby planning the initial trajectory.

[0036] In step four, the global optimum is found using a guided solution, which is defined as the average of the three optimal solutions and is calculated as follows:

[0037]

[0038] In the formula, This represents the guiding solution for the time-th iteration. These represent the top three optimal solutions generated in the previous iteration, and time represents the number of iterations.

[0039] The new solution is updated using the mean, exponential rate, and standard deviation. The i-th solution in the initial matrix at the time-th iteration is defined by equation (18):

[0040]

[0041] In equation (18), This represents the i-th solution in the initial matrix during the time-th iteration. This represents the i-th solution in the memoryless matrix during the time-th iteration;

[0042] The specific update strategy during the development phase is shown in equation (19):

[0043]

[0044] a = f 10 b = f 5 f = 2*rand-1 (20)

[0045] It is a random number in the range [0,1], σ 2 The exponential variance is represented by the calculation method shown in equations (21-23). ​​In equation (20), rand is also a random number in the range [0,1].

[0046]

[0047] Parameter a controls the distance between the winner and its variance, while parameter b expands the search space around the guided solution. This represents the ratio of useful information shared from the winners, and b controls the distance between the losers and the exponential variance.

[0048] The specific update strategy for the exploration phase is as follows:

[0049]

[0050] In the formula, It's a strategy of time + 1 times, M time This represents the average of all solutions obtained in the original population;

[0051]

[0052] c = d * f (26)

[0053] In the formula, f is a random number in the range [-1, 1], and d is an adaptive parameter;

[0054]

[0055] In the formula, time is the current iteration number, and Max_time is the maximum iteration number;

[0056] Z1=M-D1+D2 (28)

[0057] Z2=M-D2+D1 (29)

[0058] D1 = M-Xwinners rand1 (30)

[0059] D2 = M-Xwinners rand2 (31)

[0060] In the formula, Xwinners rand1 Xwinners rand2 D1 and D2 represent the distances between the average solution and the two random winners, respectively.

[0061] Formula (24) is the specific method for updating the scheme during the exploration phase of the exponential distribution optimization algorithm. The final strategy is confirmed by two random winners, and the global optimal scheme is explored.

[0062] The exponential distribution optimization algorithm uses the following steps to simply explore the global optimum search space, resulting in an improved exponential distribution optimization:

[0063] The specific improvement strategy is as follows: three optimal strategies are given adaptive weights to obtain the final guided solution. The specific calculation method is shown in equation (32):

[0064]

[0065] In the formula, A, B, and C represent the costs of the first three optimal strategies, respectively. These are three optimal strategies;

[0066] The coefficients in the improved equation (19) are specifically improved as shown in equation (33):

[0067]

[0068] In the formula, σ 2 Let a represent the standard deviation, where a is a random number in (0, 1), and a + b = 1;

[0069] The improved exponential distribution optimization is used for trajectory planning. The planned UAV trajectories are all composed of N+1 sub-trajectories, that is, each trajectory has N+2 endpoints.

[0070] In two-dimensional space, the improved exponential distribution optimization is used to update the x and y coordinates of each endpoint one by one. The process of updating all endpoints of a track is one iteration. After each iteration, the cost of the new track is calculated. All tracks are sorted in ascending order of cost. The tracks with the lowest cost are retained as the initial population of the optimization algorithm. Multiple iterations are performed until the maximum number of iterations of the optimization algorithm is reached. The cost of the UAV is optimized to the minimum and the algorithm ends.

[0071] In three-dimensional space, this method requires updating the x, y, z coordinates of each endpoint one by one during the endpoint update process, while the rest of the process is the same as in two-dimensional space.

[0072] In step four, there are n+1 control points D. i (i = 0, 1, ..., n) and node vector T = {t0, t1, ..., t} m}, then the expression for the k+1 order (k) B-spline is shown in equation (34), where the number of control points, degree, and node vector must satisfy equation (35):

[0073]

[0074] m = n + k + 1 (35)

[0075] Among them, F i,k (t) is a k-th degree B-spline basis function, where k represents the degree of the basis function, and the range of k is: 2 <= k <= n+1;

[0076] The basis function F is obtained using the deBoor-Cox recursion. i,k The expression for (t) is shown in equation (36):

[0077]

[0078] Among them, t i It is a continuously changing sequence of non-decreasing values ​​called the node vector, and this sequence is t. i , t i+1 ,...,t i+k , interval [t i , t i+k ] is called F i,k The support range of (t);

[0079] Let the overlap between the starting point and the target point be k, that is, the node vector satisfies equation (37):

[0080]

[0081] Generally, a higher degree results in a higher degree of derivative of the curve, leading to more zeros. More zeros in the derivative result in more extreme values ​​on the original curve, causing more peaks and valleys. Conversely, a lower degree results in a better approximation of the control points by the spline curve. Third-order B-spline curves achieve continuity in the second derivative, making them the most suitable choice for trajectory planning.

[0082] Compared with the prior art, the specific beneficial effects of this invention are reflected in:

[0083] 1. The optimization algorithm of this invention has a fast convergence speed and can optimize the UAV's trajectory in fewer iterations compared to other intelligent optimization algorithms.

[0084] 2. The cost function of the present invention is more comprehensive. The algorithm cost function mentioned in the present invention takes into account the track distance, turning angle and distance to obstacles, so that the optimized track distance is short, the turning angle is large and it is far away from obstacles.

[0085] 3. This invention optimizes the flight path of the UAV, resulting in a smoother flight path that better reflects the actual flight conditions of the UAV.

[0086] 4. This invention applies a new optimization algorithm to UAV trajectory planning and makes relevant improvements to the optimization algorithm, which is beneficial to the research of UAV trajectory planning.

[0087] 5. The algorithm proposed in this invention can be used in a variety of scenarios. For example, in two-dimensional space, this invention can plan the flight path for delivery drones to improve the efficiency of drone delivery. In three-dimensional space, when a drone malfunctions during a drone formation performance, this invention can quickly plan an optimal flight path for a backup drone, allowing the backup drone to arrive quickly and ensuring the performance effect. Attached Figure Description

[0088] Figure 1 This is the overall flowchart of the present invention.

[0089] Figure 2 This is a schematic diagram of the offset angle between two sub-tracks.

[0090] Figure 3 A schematic diagram of the sub-track to the center of the obstacle circle. Figure 3 (a) is a schematic diagram of the projection of the center of the circle onto the line segment. Figure 3(b) is a schematic diagram showing the projection of the center of the circle onto the extension line at the end of the line segment. Figure 3 (c) is a schematic diagram of the projection of the center of the circle onto the extension line at the beginning of the line segment.

[0091] Figure 4 A flowchart for planning a collision-free trajectory.

[0092] Figure 5 This is a schematic diagram of subspace partitioning.

[0093] Figure 6 This is a flowchart of the exponential distribution optimization algorithm.

[0094] Figure 7 The planned flight path map is optimized for two algorithms. Figure 7 (a) The track map planned to improve the exponential distribution optimization. Figure 7 (b) is the trajectory map planned by the Grey Wolf Optimization Algorithm.

[0095] Figure 8 This is a comparison of the convergence curves of the two algorithms, namely the one proposed in this invention and the one optimized by Grey Wolf, in two-dimensional space.

[0096] Figure 9 To improve the exponential distribution, the track map obtained by optimization in three-dimensional space is used. Figure 9 (a) Front view of the flight path. Figure 9 (b) is a side view of the flight path. Figure 9 (c) is a top view of the flight path.

[0097] Figure 10 Optimize the trajectory map of the gray wolf in three-dimensional space. Figure 10 (a) Front view of the flight path. Figure 10 (b) is a side view of the flight path. Figure 10 (c) is a top view of the flight path.

[0098] Figure 11 This is a comparison of the convergence curves of the two algorithms, namely the one proposed in this invention and the one optimized by Grey Wolf, in three-dimensional space. Detailed Implementation

[0099] To make the technical problems, solutions, and beneficial effects of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0100] Based on an improved exponential distribution optimization method for UAV trajectory planning, this invention uses an improved exponential distribution optimization algorithm to optimize UAV trajectories. First, multiple collision-free trajectories with fixed trajectory segments are planned and sorted in ascending order. Then, the improved exponential distribution optimization algorithm is used to optimize the trajectories, ultimately generating collision-free trajectories with short distances, small turning angles, and high distances from obstacles. The specific implementation process is as follows: Figure 1 As shown.

[0101] The algorithm of this invention can be used in two-dimensional and three-dimensional spaces where obstacle conditions are known. In two-dimensional space, the map is placed in the first quadrant of the planar coordinate system. Obstacles can be adjusted from irregular shapes to circular obstacles using an approximation method. The positions of the irregular obstacles are known, so the center and radius of the circular obstacles are also known. In three-dimensional space, the map is placed in the first quadrant of the spatial coordinate system. Obstacles can be adjusted from irregular shapes to spherical obstacles using an approximation method. The positions of the irregular obstacles are known, so the center and radius of the spherical obstacles are also known. Finally, to facilitate track optimization, each track is defined as consisting of N+1 sub-tracks, i.e., each track has N+2 endpoints. In two-dimensional space, these endpoints are described as {D}. s ,D0,D1,...,D N D g}, specifically {(x s ,y s ),(x0,y0),(x1,y1)...(x N ,y N ),(x g ,y g )}, where (x s ,y s (x) is the starting coordinate, (x) g ,y g The coordinates of the endpoint are given. There are M obstacles in the space, which have been approximated as circular obstacles with centers {(X1,Y1),(X2,Y2),...,(X... M ,Y M )}, with radii {R1,R2,...,R}. M}

[0102] 1. Construct the cost function:

[0103] The cost function C proposed in this invention comprehensively considers the flight path distance l1, the turning angle θ, and the distance to the obstacle l2. These three costs are normalized, and then a certain weight is given to each of the three costs and they are added together to obtain the final cost function.

[0104] (1) Calculation of track distance l1: Each track consists of N+1 sub-tracks. This method sums the distances of the sub-tracks to obtain the track distance l1. The specific algorithm is shown in equation (1):

[0105]

[0106] (2) Calculation of track offset angle θ: Each track consists of N+1 sub-tracks. This method calculates the offset angle of each subsequent sub-track relative to the previous sub-track and sums them to obtain the total offset angle. For example... Figure 2 As shown, θ 012 It is the offset angle between the two sub-tracks formed by the three points {D0, D1, D2}, and the specific calculation method is shown in equation (2-6):

[0107]

[0108] (3) Calculation of distance l2 to obstacles: There are M obstacles in the space that are adjusted to be circular, and the center and radius of the circle are known. Each planned path consists of N+1 sub-paths. It is necessary to calculate the distance from each sub-path to each obstacle and determine the minimum distance l from the sub-path to the obstacle. min Is it greater than the radius R of the obstacle? i Add safety margin l safe If true, the distance l2 to the obstacle is determined as the length l of the sub-track; if false, the distance l2 to the obstacle is determined as the maximum track length l. max Specifically, as shown in equation (12).

[0109] Calculating the distance from each sub-track to the obstacle is equivalent to calculating the distance from the center of the circle to the line segment. This method uses a vector algorithm to calculate this distance. The distance from the center of the circle to the line segment can be divided into three cases, as detailed below. Figure 3 As shown, R i Let D be the center of the obstacle circle, D0 be the beginning of the sub-track, D1 be the end of the sub-track, and D be the center of the circle. i Projecting onto the line segment, 'a' shows the projection of the center of the circle onto the line segment, with the distance from the center to the line segment being point R. i to line segment Distance, line segment The representation is shown in equation (7); the projection of the center of the circle in b is on the extension of the end of the line segment, and the distance from the center of the circle to the line segment is the center R. i The distance to the endpoint D1; the projection of the center of circle c onto the extension of the beginning of the line segment, and the distance from the center of the circle to the line segment is the center R. i Distance to the first endpoint D0. Construct a vector. Based on vectors Calculate vectors Specifically, as shown in equation (7), where · represents the dot product of two vectors. Based on the magnitude of r in equation (8), we can determine which of a, b, or c it belongs to. The specific judgment method is shown in equation (9), thereby calculating the distance from the center of the circle to the line segment.

[0110]

[0111] l min =min{l o0 ,l o1 ...l oM} (11)

[0112]

[0113] The above process has yielded the track distance l1, track offset angle θ, and distance to the obstacle l2. Now, these three costs need to be normalized and then weighted and summed to obtain the final cost function. Normalizing the track distance l1 and the distance to the obstacle l2 requires dividing by the maximum track distance l. max To normalize the track offset angle θ, it is necessary to divide by the maximum offset angle θ. max Normalization and summation yield the final cost function, as shown in equation (14):

[0114]

[0115] In equation (14), N is the number of intermediate endpoints, and α+β+γ=1.

[0116] 2. Collision detection strategy:

[0117] In two-dimensional space, the obstacle is a circular obstacle, and the trajectory consists of N+1 sub-trajectories. Collision detection is performed by determining whether the distance from each sub-trajectory segment to the center of each obstacle circle is greater than the radius of the obstacle. If a sub-trajectory segment intersects with the obstacle circle, then the trajectory intersects with the obstacle; otherwise, they do not intersect. The distance from the sub-trajectory to the circular obstacle is calculated as shown in equation (10).

[0118] 3. Generation of the original collision-free trajectory:

[0119] This invention uses an exponential distribution optimization algorithm to optimize UAV trajectories. Therefore, it is necessary to plan several collision-free trajectories first. The specific process for planning collision-free trajectories is as follows: Figure 4 As shown. The above process plans a collision-free path. By repeating this process multiple times, multiple collision-free paths can be planned.

[0120] The specific process of subspace partitioning is as follows: Figure 5 As shown in the diagram, S is the starting point, G is the ending point, and x... i ,yi Let x and y be the lengths of the subspace along the x-axis and y-axis, respectively. The specific calculation methods are shown in equations (15) and (16).

[0121]

[0122] Based on the specific lengths of the subspaces along the x and y axes, the exact location of the i-th subspace can be determined. The coordinates of the four points in the i-th subspace are as follows:

[0123]

[0124] Based on the coordinates of the subspaces mentioned above, the specific location of each subspace can be determined, thereby planning the initial trajectory.

[0125] 4. Improved exponential distribution optimization algorithm:

[0126] (1) Basic principles of exponential distribution optimization

[0127] The Exponential Distribution Optimization (EDO) algorithm is a novel population-based metaheuristic algorithm inspired by the exponential probability distribution in mathematics. In the initialization phase, EDO generates a population of random solutions representing multiple exponential distribution models. The position of each solution is considered as a random variable following the exponential distribution, and the set of these random variables constitutes the potential solutions in the algorithm's search process. The EDO algorithm operates based on two phases: development and exploration. The development phase utilizes the memoryless property, exponential rate, standard deviation, and mean of the exponential distribution model, and employs a bootstrapping solution to guide the search process to the global optimum. The exploration phase uses the two winners from the original population to establish the EDO optimization model.

[0128] The development phase revolves around finding the global optimum using the guiding solution (Xguide). The guiding solution, containing useful information about the global optimum, is defined as the average of the three optimal solutions and is calculated as follows:

[0129]

[0130] In the formula, time represents the number of iterations. However, relying solely on the guided solution can lead the algorithm into local optima. This method uses the mean, exponential rate, and standard deviation to update the new solution, and also introduces the concepts of winners and losers. If the new solution achieves better accuracy than the corresponding solution in the overall population, it is defined as a winner, which can also be defined as expression (18); solutions that do not achieve better accuracy are defined as losers, and losers also play a role in updating the solution.

[0131]

[0132] In the formula, This represents the i-th solution in the initial matrix during the time-th iteration. Let represent the i-th solution in the memoryless matrix at the time-th iteration. It can be concluded that if the solution is updated in both the x-winner and memoryless matrix, then it is a winning solution. Therefore, the specific update strategy during the development phase is shown in equation (19):

[0133]

[0134] In equation (19), It is a random number in the range [0,1], σ 2 The exponential variance is represented by the formula (21-23), and rand in formula (20) is also a random number in the range [0,1].

[0135]

[0136] The first branch of equation (19) revolves around searching for a guiding solution to obtain a solution superior to the current winner. Parameter a controls the distance between the winner and its variance, and b expands the search space around the guiding solution. The second branch shows the loser moving in the direction of the winner. represents the ratio of useful information that can be shared from the winners, and 'b' controls the distance between the losers and the exponential variance. The development phase of the exponential distribution optimization algorithm seeks the optimal solution through the aforementioned update strategy.

[0137] The exploration phase of the algorithm uses the two winners from the original population to find the globally optimal strategy. The specific update strategy for the exploration phase is as follows:

[0138]

[0139] In the formula, It's a strategy of time + 1 times, M time This represents the average of all solutions obtained in the original population.

[0140]

[0141] c = d * f (26)

[0142] In the formula, f is a random number in the range [-1, 1], and d is an adaptive parameter.

[0143]

[0144] In the formula, time is the current iteration number, and Max_time is the maximum iteration number.

[0145] Z1=M-D1+D2 (28)

[0146] Z2=M-D2+D1 (29)

[0147] D1 = M-Xwinners rand1 (30)

[0148] D2 = M-Xwinners rand2 (31)

[0149] In the formula, Xwinners rand1 Xwinners rand2 D1 and D2 represent the distances between the average solution and the two random winners.

[0150] Formula (24) is the specific method for updating the scheme in the exploration phase of the exponential distribution optimization algorithm. The final strategy is confirmed by two random winners, and the global optimal scheme is explored.

[0151] The exponential distribution optimization algorithm uses the following steps to simply explore the global optimum search space, the specific implementation process of which is as follows: Figure 6 As shown.

[0152] (2) Improved exponential distribution optimization:

[0153] The original exponential distribution optimization does not achieve good results when applied to trajectory planning. Therefore, relevant improvements have been made to the exponential distribution optimization algorithm, which is more suitable for UAV trajectory planning. The specific improvement strategy is as follows: First, the calculation method of the guiding solution in the development stage has been improved. Instead of using the average of the three optimal strategies as the guiding solution, the three optimal strategies are given adaptive weights to obtain the final guiding solution. The specific calculation method is shown in Equation (32):

[0154]

[0155] In the formula, A, B, and C represent the costs of the first three optimal strategies, respectively. These are three optimal strategies.

[0156] Secondly, during the development phase of the exponential distribution optimization, the updated points may exceed the boundary. To address this issue, this method makes relevant improvements to the coefficients in equation (19), as shown in equation (33):

[0157]

[0158] In the formula, σ 2Let represent the mean squared error, and 'a' be a random number in the range (0, 1), where a + b = 1. This method makes relevant improvements to the exponential distribution optimization, and the improved exponential distribution optimization algorithm is more suitable for UAV trajectory planning.

[0159] (3) Improved exponential distribution optimization for trajectory planning:

[0160] The UAV trajectories planned by this method consist of N+1 sub-trajectories, meaning each trajectory has N+2 endpoints. In two-dimensional space, an improved exponential distribution optimization is used to update the x, y coordinates of each endpoint one by one. Updating all endpoints of a trajectory constitutes one iteration. After each iteration, the cost of the new trajectory is calculated, and all trajectories are sorted in ascending order of cost. The trajectories with the lowest costs are retained as the initial population for the optimization algorithm. Multiple iterations are performed until the maximum number of iterations is reached, at which point the UAV cost is minimized, and the algorithm terminates. In three-dimensional space, this method requires updating the x, y, z coordinates of each endpoint one by one during the endpoint update process; the other processes are the same as in two-dimensional space.

[0161] 5. Optimization of third-order B-spline curves:

[0162] The above optimization algorithm can obtain a collision-free optimal trajectory. However, the turning point of this trajectory is an angle formed by two line segments, which does not match the actual flight trajectory of the UAV. Therefore, the generated trajectory needs to be smoothed. This method uses a third-order B-spline curve to smooth the generated trajectory, producing a smooth trajectory that satisfies dynamic constraints. B-spline curves offer advantages such as smoothness and local controllability. By locally influencing the shape of the curve through control points, changing a control point only affects a local area of ​​the curve, not the entire curve. Furthermore, B-spline curves are piecewise defined polynomial curves, with smooth transitions at connection points, resulting in a very smooth overall curve. These advantages make B-spline curves widely used in trajectory planning. Therefore, this invention uses a uniform B-spline curve to smooth the generated trajectory, ensuring the stability and safety of the UAV's flight.

[0163] There are n+1 control points D i (i = 0, 1, ..., n) and node vector T = {t0, t1, ..., t} m}, then the expression for the (k)th order B-spline is shown in equation (34), where the number of control points, degree, and node vector must satisfy equation (35).

[0164]

[0165] m = n + k + 1 (35)

[0166] Among them, F i,k(t) is a k-th degree B-spline basis function, where k represents the degree of the basis function, and the range of k is 2 <= k <= n+1. This method uses the de Boor-Cox recursion to obtain the basis function F. i,k The expression for (t) is shown in equation (36):

[0167]

[0168] In the formula, t i It is a continuously changing sequence of non-decreasing values ​​called a node vector, whose first and last values ​​are generally defined as 0 and 1, and this sequence is t. i , t i+1 ,...,t i+k , interval [t i , t i+k ] is called F i,k The support range of (t).

[0169] In order for the curve equation to satisfy the constraints of the starting point and the target point, let the repetition degree of the starting point and the target point be k, that is, their nodal vectors satisfy equation (37):

[0170]

[0171] Generally, a higher degree results in a higher degree of derivative of the curve, leading to more zeros. More zeros in the derivative result in more extreme values ​​on the original curve, causing more peaks and valleys. Conversely, a lower degree results in a better approximation of the control points by the spline curve. Third-order B-spline curves achieve continuity in the second derivative, making them the most suitable choice for trajectory planning.

[0172] To verify the algorithm mentioned in this invention, a simulation experiment was conducted using Python 3.10. A 100×100 map was set in a two-dimensional environment, with the unit set to meters (m). In the scene, black circles represent obstacles, "—" represents the planned trajectory, and "---" represents the trajectory optimized by exponential distribution. The starting coordinates of the UAV were set to [2,3], the target point to [90,94], and the obstacles were four irregular shapes of different sizes, all of which had been expanded into circular obstacles. This method was iterated 500 times to observe the optimization of the trajectory, and the algorithm of this invention was compared with the Grey Wolf optimization algorithm. Figure 7 (a) is the trajectory planned by the improved exponential distribution optimization. Figure 7 (b) is the trajectory planned by the Grey Wolf Optimization Algorithm.

[0173] After multiple iterations, this method plots the convergence curves of the two algorithms, as shown in the following figure. Figure 8As shown in the figure, compared to the Grey Wolf optimization algorithm, the algorithm proposed in this invention is more suitable for UAV trajectory planning in two-dimensional space, and can minimize trajectory cost in fewer iterations.

[0174] A 100×100×100 map was set up in a 3D environment, with the unit set to meters (m). Black spheres in the scene represent obstacles, "—" represents the planned flight path, and "---" represents the exponentially optimized flight path. The starting coordinates of the UAV were set to (2, 3, 6), and the target point was set to (90, 94, 96). The starting point was represented by a triangle in the map, and the ending point by a pentagram. The obstacles consisted of five irregular shapes of varying sizes, all of which had been expanded into spherical obstacles. This method iterated 500 times to observe the optimization of the flight path, and the algorithm of this invention was compared with the Grey Wolf optimization algorithm.

[0175] Figure 9 It is a trajectory obtained by optimizing the improved exponential distribution in three-dimensional space, and the obtained trajectory is displayed from three perspectives: front view, side view, and top view. Figure 9 (a), 9(b), and 9(c) show the flight paths in the three cases, respectively.

[0176] Figure 10 This is the trajectory optimized by Grey Wolf in three-dimensional space, displayed from three perspectives: front view, side view, and top view. Figure 10 (a), 10(b), and 10(c) show the tracks in the three cases, respectively.

[0177] After multiple iterations, this method plots the convergence curves of the two algorithms, as shown in the following figure. Figure 11 As shown in the figure, compared to the Grey Wolf optimization algorithm, the algorithm proposed in this invention is also more suitable for UAV trajectory planning in three-dimensional space, and can minimize trajectory cost in fewer iterations.

[0178] The improved exponential distribution optimization algorithm was applied to UAV trajectory planning in two-dimensional and three-dimensional spaces. Simulation results show that, in both two-dimensional and three-dimensional spaces, the algorithm can minimize the cost of the UAV trajectory in a few iterations. Furthermore, the optimized trajectory is far from obstacles, ensuring the safe flight of the UAV. The optimized trajectory has a short distance, reducing the flight energy consumption of the UAV. Optimizing the optimized trajectory with a third-order B-spline curve results in a smoother trajectory that is more suitable for actual UAV flight conditions.

[0179] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included within the scope of the present invention.

Claims

1. A UAV trajectory planning method based on improved exponential distribution optimization, characterized in that, The specific steps are as follows: Step 1: Construct the cost function: for the UAV's flight path distance l1 and the UAV's turning angle... The three costs of the drone, the distance l2 from the obstacle, and the distance l2 are normalized, and weights are added to each cost to obtain the cost function. In step one, the specific steps for constructing the cost function are as follows: 11) Calculate the track distance l1: Each track consists of N+1 sub-tracks. The track distance l1 is obtained by summing the distances of the sub-tracks. The specific algorithm for solving this problem is shown in equation (1): （1） 12) Calculate the track offset angle : Calculate the offset angle of each subsequent sub-track segment relative to the previous sub-track segment and sum them to obtain the total offset angle. yes The offset angles of the two sub-tracks formed by these three points. The coordinates are respectively The specific calculation method is shown in equation (2-6): （2） （3） （4） （5） （6） 13) Calculate the distance l2 from the sub-track to the obstacle: If the minimum distance l from the sub-track to the obstacle is... min Larger than the radius R of the obstacle i Add safety margin l safe The distance l2 from the sub-track to the obstacle is defined as the length l of the sub-track; if the minimum distance l from the sub-track to the obstacle is... min Less than or equal to the radius R of the obstacle i Add safety margin l safe The distance l2 from the sub-track to the obstacle is determined as the maximum track length l. max Specifically, as shown in equation (12); Calculate the distance from each sub-track to the obstacle, i.e., use a vector algorithm to calculate the distance from the center of the obstacle circle to the line segment; first, construct the vector. , , According to vector , Calculate vectors Specifically, as shown in equation (8); then, the value of r is calculated according to equation (9) for judgment, and the specific judgment method is shown in equation (10); equation (10) calculates the center R of the circle based on the value of r. i to line segment Distance, line segment The representation is shown in equation (7): （7） （8） （9） The distance from the subtrack to the circular obstacle is shown in equation (10): （10） （11） （12） （13） Among them, R i Let D be the center of the obstacle circle, D0 be the beginning of the sub-track, D1 be the end of the sub-track, and D be the center of the circle. i The projection onto the line segment, in equation (7) * represents the product of a number and a variable; in equation (8) · represents the dot product of two vectors, || represents the magnitude of the vector, || 2 In equation (10), the square of the vector magnitude is expressed. Let | represent the distance from a certain segment of the flight path to the i-th obstacle, and | represent the absolute value of a certain number; in equation (12) This represents the minimum distance from the i-th sub-track to all obstacles; For the normalization of the track distance l1 and the distance to the obstacle l2, for the track offset angle The normalization and summation yield the final cost function, as shown in equation (14): （14） In equation (14), N is the number of intermediate endpoints. ; Where l1 is the track distance, l is the track offset angle, l2 is the distance from the sub-track to the obstacle, and the maximum track distance is l. max ; Step 2: Construct a collision detection strategy: The drone's trajectory consists of N+1 sub-trajectories. If the distance from a sub-trajectory to the center of each obstacle is less than or equal to the radius of the obstacle, then the sub-trajectory is determined to intersect with the obstacle, and the drone's trajectory is determined to intersect with the obstacle. If the distance from each sub-trajectory to the center of each obstacle is greater than the radius of the obstacle, then all sub-trajectories do not intersect with the obstacles, and the drone's trajectory is determined to not intersect with the obstacles. Step 3: Generate initial collision-free tracks: Plan several collision-free tracks; Step 4: Optimize the UAV trajectory using the exponential distribution optimization algorithm; Step 5: Use a third-order B-spline curve to smooth the generated track, resulting in a smooth track that satisfies dynamic constraints.

2. The UAV trajectory planning method based on improved exponential distribution optimization according to claim 1, characterized in that, In two-dimensional space, the map is placed in the first quadrant of the plane coordinate system, and the approximation method is used to adjust the irregularly shaped obstacles into circular obstacles; In three-dimensional space, the map is placed in the first quadrant of the spatial coordinate system, and an approximation method is used to adjust irregularly shaped obstacles into spherical obstacles; Each UAV's trajectory is defined as consisting of N+1 sub-trajectories, and each UAV trajectory has N+2 endpoints. In two-dimensional space, the N+2 endpoints are... The coordinates of each endpoint are { },in, As the starting coordinates, Given the endpoint coordinates, there are M obstacles in the space. These obstacles are approximated as circular obstacles, with their centers being {...} The radii of the obstacles are { } 3. The UAV trajectory planning method based on improved exponential distribution optimization according to claim 2, characterized in that, In step three, the space is divided into N subspaces, with S as the starting point, G as the ending point, and x... i Let y be the length of the subspace along the x-axis, and y be the length of the subspace along the x-axis. i Let be the length of the subspace along the y-axis. The specific calculation method is shown in the formula: （15） （16） Based on the specific lengths of the subspaces along the x and y axes, the exact location of the i-th subspace is determined. The coordinates of the four points in the i-th subspace are: ; }； The specific location of each subspace is determined based on the coordinates of the aforementioned subspaces, thereby planning the initial trajectory.

4. The UAV trajectory planning method based on improved exponential distribution optimization according to claim 3, characterized in that, In step four, the global optimum is found using a guided solution, which is defined as the average of the three optimal solutions and is calculated as follows: （17） In the formula, This represents the guiding solution for the time-th iteration. , , These represent the top three optimal solutions generated in the previous iteration, and time represents the number of iterations. The new solution is updated using the mean, exponential rate, and standard deviation. The i-th solution in the initial matrix at the time-th iteration is defined by equation (18): （18） In equation (18), This represents the i-th solution in the initial matrix during the time-th iteration. This represents the i-th solution in the memoryless matrix during the time-th iteration; The specific update strategy during the development phase is shown in equation (19): （19） （20） It is a random number in the range [0,1]. The exponential variance is represented by the calculation method shown in equations (21-23). ​​In equation (20), rand is also a random number in the range [0,1]. （21） （22） （23） Parameter a controls the distance between the winner and its variance, while parameter b expands the search space around the guided solution. This represents the ratio of useful information shared from the winners, and b controls the distance between the losers and the exponential variance. The specific update strategy for the exploration phase is as follows: （24） In the formula, It's a strategy of time + 1 times. This represents the average of all solutions obtained in the original population; （25） （26） In the formula, f is a random number in the range [-1, 1], and d is an adaptive parameter; （27） In the formula, time is the current iteration number, and Max_time is the maximum iteration number; （28） （29） （30） （31） In the formula, , D1 and D2 represent the distances between the average solution and the two random winners, respectively. Formula (24) is the specific method for updating the scheme during the exploration phase of the exponential distribution optimization algorithm. The final strategy is confirmed by two random winners, and the global optimal scheme is explored. The exponential distribution optimization algorithm uses the following steps to simply explore the global optimum search space, resulting in an improved exponential distribution optimization: The specific improvement strategy is as follows: three optimal strategies are given adaptive weights to obtain the final guided solution. The specific calculation method is shown in equation (32): （32） In the formula, A, B, and C represent the costs of the first three optimal strategies, respectively. , , These are three optimal strategies; The coefficients in the improved equation (19) are specifically improved as shown in equation (33): （33） In the formula, Let a represent the mean squared error, where a is a random number in the range (0, 1), and a + b = 1. The improved exponential distribution optimization is used for trajectory planning. The planned UAV trajectories are all composed of N+1 sub-trajectories, that is, each trajectory has N+2 endpoints. In two-dimensional space, the coordinates of each endpoint are updated one by one using an improved exponential distribution optimization. The process of updating all endpoints of a track is one iteration. After each iteration, the cost of the new track is calculated. All tracks are sorted in ascending order of cost, and the tracks with the lowest cost are retained as the initial population of the optimization algorithm. Multiple iterations are performed until the maximum number of iterations of the optimization algorithm is reached, the cost of the UAV is optimized to the minimum, and the algorithm ends. In three-dimensional space, the coordinates of each endpoint are updated one by one during the process of updating the endpoints, and the rest of the process is the same as in two-dimensional space.

5. The UAV trajectory planning method based on improved exponential distribution optimization according to claim 4, characterized in that, In step four, there are n+1 control points. and node vectors Then the expression for the (k)th order B-spline is shown in equation (34), where the number of control points, degree, and node vector must satisfy equation (35): （34） （35） in, It is a k-th degree B-spline basis function, where k represents the degree of the basis function, and the range of k is: 2<=k<=n+1; The basis functions are obtained using the deBoor-Cox recursion. The expression is shown in equation (36): （36） Among them, t i It is a continuously changing sequence of non-decreasing values ​​called the node vector, and this sequence is t. i , t i + 1 ,...,t i + k , interval [t i , t i + k ] is called The support range; Let the overlap between the starting point and the target point be k, that is, the node vector satisfies equation (37): （37）。

Images