# Trajectory Planning of UAVs with Fault Tolerance Based on Monte Carlo Sampling

Authors: Jianwei Wu; Lin Chen; Yang Zhou; Fuyun Liu
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Due to the limitation of system positioning accuracy, it is necessary to correct an unmanned aerial vehicle’s (UAV) error by some correct points arranged in advance, in order to successfully accomplish their tasks. Planning the optimal trajectory of UAVs considering error correction with fault tolerance is a great challenge that almost has not been broken out because corrections may fail by these points. To address the issue, this paper proposes two trajectory planning methods based on “the definite rule” and Monte Carlo sampling. By constructing five calculation models of the two methods in detail, designing the algorithms appropriately, and developing computational programs based on MATLAB, the results of trajectory planning of UAVs from a practical problem are obtained. The results achieve the optimization objectives that the total trajectory length is as short as possible and that the correction points UAV passes are as few as possible, and meanwhile, the optimal trajectory satisfies all the constraints, which illustrate the feasibility and reasonability of trajectory planning based on the definite rule and Monte Carlo sampling. The comparable results show that computation time is 104 and 195 times that of this paper when the Dijkstra algorithm with ant colony algorithm and the greedy algorithm with tabu search algorithm are respectively used, which illustrate the high efficiency of the proposed method. This work provides a feasible solution for UAVs’ trajectory planning that considers the error corrections and failure probability on certain correction points and achieves trajectory planning of UAVs with fault tolerance.


---

## Body

## 1. Introduction

Trajectory planning of UAVs refers to generating an optimal flight trajectory between the starting and ending points, considering factors such as energy consumption, maneuverability, arrival time, terrain, and threat [1–3]. It is one of the key technologies for UAVs and helps to accomplish many missions such as reconnaissance [4, 5], surveillance [6, 7], road patrol [8], aerial photo [9], object detection [10, 11], and cargo transfer [12].It is a very important research topic to realize trajectory planning of a UAV fast in complex environments. Due to the limitation of the UAV’s structure and the influence of the complex environments, its positioning system cannot be precisely located [13]. Once the positioning errors accumulate to a certain extent, the task of a UAV may fail. At present, there is a technology of error correction [14–16], which is used for eliminating UAV’s errors when it flies into the vicinity of some points (called correction points) arranged in advance. The technology can ensure that a UAV flies along the planned trajectory and performs its task successfully. However, affected by complex environments such as weather and terrain, these points could be invalid even after the positioning errors of a UAV are corrected by them. Accordingly, the real challenge is how to ensure that the optimal trajectory for UAVs has the ability of fault tolerance.The past decades have witnessed many studies in trajectory planning for UAVs. To realize the trajectory planning, these studies constructed single or multiple objective optimization with multiconstraints by considering the energy consumption, security (for example, the fewer turning operations), and the shortest path, under constraints such as threat, obstacle, terrain, and arrival time. Researchers have studied its modeling, optimization, and solution in detail, and a large number of algorithms have been presented as follows: Voronoi diagram method [17], A∗ algorithm [18], improved A∗ algorithm [19], artificial potential field method [20], algorithm based on disturbed fluid and trajectory propagation [21], genetic algorithm [22, 23], particle swarm optimization (PSO) algorithm [24], neural network algorithm [25, 26], simulated annealing algorithm [27], and so on. These algorithms having their respective characteristics provide us with a basis for selecting them according to different tasks. For example, graph-based approaches are usually very accurate, but since some missions require several constraints on the resultant path, they face some difficulties, such as a long convergence time. Artificial potential field method has the advantages of simple principle and small amount of computation, yet there may be some oscillation of the trajectory if it enters into the narrow gap between two obstacles; heuristic algorithms, such as genetic algorithm and particle swarm optimization algorithm, are also used for multiobjective optimization but do not guarantee that an optimal solution will be found because there is no an analytical proof [28].In spite of so many studies on trajectory planning for UAVs, they mainly focus on improving performances such as trajectory smoothness, trajectory length, obstacle crossing ability, and computation efficiency. In fact, few studies focus on the trajectory planning for UAVs considering the error correction, let alone on the trajectory planning with fault tolerance. Because there are many correction points arranged in advance and performing the error correction is conditional by them, the trajectory planning of the UAVs considering the error correction is more difficult than without considering it, especially under the three premises that UAVs are with minimal turning radius (i.e., considering the kinematical constraints), that the flying distance is short enough, and that the points used for error correction are as few as possible. Moreover, since correction points have a failure probability, the flight errors could not be entirely eliminated by them. Therefore, a greater difficulty is how to design a trajectory that can satisfy the above premises to arrive at the destination successfully, even if the error correction fails. Maybe “it is for these reasons” that the current studies can still not break through the stage that realizes the trajectory planning of UAVs with fault tolerance.To plan an optimal trajectory considering error correction with fault tolerance, this paper proposes two trajectory planning methods based on the definite rule (i.e., minimizing optimization objectives designed) and Monte Carlo sampling. By constructing five calculation models of the two methods in detail, designing the algorithms appropriately, and developing computational programs based on MATLAB, trajectory planning of UAVs with fault tolerance is achieved. Simulation example illustrates the methods’ feasibility and results’ rationality.The main contributions of this study are as follows. (i) This paper proposes two trajectory planning methods for UAVs with error correction, which overcomes the limitation in traditional methods that cannot deal with the long-distance trajectory planning for UAVs needed to consider error correction. (ii) The proposed method based on Monte Carlo sampling can perform the trajectory planning with error correction existing in failure, which further expands its application and has practical value. (iii) A set of solution programs in MATLAB 2018a programming environment are developed based on the trajectory planning methods, and the optimal trajectory is obtained by a practical problem.The rest of this paper is organized as follows. Section2 provides the necessary background. Section 3 describes in detail the trajectory planning methods with fault tolerance based on definite rule and Monte Carlo sampling and respectively establishes the following five calculation models: the calculation model of the trajectory length considering the turning radius, the optimization model of the maximum spherical surface, the optimal correction point model considering the weighted distance, the model of updating errors, and the model of extracting preferable correction points. The reasons for selecting optimization algorithms and steps of the implementation are given in Section 4. A simulation example from a practical problem is presented in Section 5, and the corresponding simulation results and their discussions are described, including the optimized trajectories under the different weights and the feasible trajectories, the optimal trajectory, and the comparison with different methods. Finally, the conclusions are drawn in Section 6. The framework of this paper is shown in Figure 1.Figure 1
Framework of the paper.


## 2. Background

Figure2 shows the flight space of UAVs with the correction points used for eliminating its positioning errors. However, affected by complex environments such as weather and terrain, some correction points could be invalid even after the positioning errors of a UAV are corrected by them. These points with the failure possibility are expressed in squares, as shown in Figure 2. In the figure, A and B are, respectively, the starting and ending points expressed in red, C2 is the first correction point searched, the vertical and horizontal correction points are shown in green and blue dots, respectively, and the black curve represents a planned trajectory. During planning a trajectory from A to B, the following constraints will be taken into account.(1)
A UAV needs real-time positioning during flight, and its positioning accuracy includes vertical and horizontal errors. When the UAV reaches the ending point, if both the vertical and horizontal errors are less thanθ m (called constraint θ in this paper), it is regarded that the UAV successfully arrives at the destination. For every 1 m flying, it is assumed that both the vertical and horizontal errors of a UAV will be increased by δm, respectively, where δ is a rough estimator based on the previous positioning performance.(2)
A UAV needs to correct the positional error during flight. The correction points include two types, horizontal and vertical corrections. When a UAV reaches a correction point, the vertical or horizontal error can be eliminated based on its type. It is assumed that the positions of correction points in the flight space can be obtained before planning flight trajectory and that these points determined by the terrain of flight space are not a uniform law. If both the vertical and horizontal errors of a UAV are corrected in time, it can fly along the planned trajectory and finally arrive at the intended destination.(3)
At the starting point A, the vertical and horizontal errors of a UAV are both zero.(4)
The vertical error of a UAV will be eliminated after it is successfully corrected by vertical correction points, but its horizontal error will remain unchanged.(5)
The horizontal error of a UAV will be eliminated after it is successfully corrected by horizontal correction points, but its vertical error will remain unchanged.(6)
Vertical correction of a UAV can be performed when its vertical and horizontal errors are no more thanα1and α2, respectively.(7)
Horizontal correction of a UAV can be performed when its vertical and horizontal errors are no more thanβ1and β2, respectively.(8)
A UAV has a minimum turning radiusr. Limited by the structure and control system, it cannot finish an immediate turning (i.e., the direction of a UAV cannot be changed abruptly) and hence has a minimum turning radius.(9)
The probability of a UAV being successfully corrected isσ, when it passes a correction point with the failure probability.(10)
Flight error of a UAV, not more thanζ, will be leftover rather than completely eliminated, when the error correction fails by correction points with the failure probability. Namely, er = min{ef, ζ}, where ef and er are the errors before and after being corrected, respectively.Figure 2
Flight space of UAVs with correction points.From the above descriptions, it can be known that the aim of this study is to find the optimal trajectory for a UAV considering all the above constraints. Even if error corrections fail, the trajectory planning methods with failure tolerance still can ensure that the UAV successfully arrives at the destination. Surely, the optimal trajectory needs to achieve the objectives that trajectory length is as short as possible and that the number of corrections the UAV passes is as small as possible.


## 3. Trajectory Planning Methods with Failure Tolerance

During the flight of a UAV from the starting point to the ending point, for eliminating its positioning errors, it needs to pass some correction points one by one and reaches the destination at last. However, there exists the failure probability after being corrected by certain points. To ensure UAVs have the ability of fault tolerance, there is a conservative handling method assuming that the error corrections all fail by the correction points with the failure probability. Under this assumption, without doubt, if a UAV can fly along a series of trajectories, if one of them subject to all the constraints in Section2 is the shortest, it must be the optimal trajectory along which the UAV can successfully arrive at the destination.Based on the handling method, the key is to search the best correction points by means of a definite rule, in order to obtain the optimal trajectory considering the error correction. To achieve the purpose, “the definite rule” mentioned in this paper is to construct the optimization models and to minimize their objective function under constraints. As a result, the trajectory planning of UAVs with fault tolerance becomes an optimization problem, and its objectives include that the trajectory length is the shortest and that the correction points UAV passes are as small as possible.Therefore, the trajectory planning model with the error correction will be established by using the handling method and by considering the optimization objectives and used for obtaining a series of the best correction points. To achieve trajectory planning considering the error correction with fault tolerance, we establish the calculation model of the trajectory length considering the minimal turning radius, the optimization model of the maximum spherical surface, the optimal correction point model considering the weighted distance, the model of updating errors, and the model of extracting preferable correction points. The trajectory planning methods with failure tolerance for UAVs can be formed by integrating these models, and the optimal trajectory with failure tolerance will be obtained by solving them.


### 3.1. Modeling Based on the Definite Rule

#### 3.1.1. Calculation Model of the Trajectory Length considering Turning Radius

The calculation model considering turning radius is shown in Figure3, where Ci−1, Ci, and Ci+1 denote the previous starting point, the current starting point (i.e., the best correction point searched at the previous iteration), and a new correction point searched, respectively. It can be known that the following flight trajectory between starting point to the correction point can ensure the flight distance is small, when the UAV flies into point Ci+1, its flight direction is just tangent to straight line CiCi+1; namely, during the flying from the point Ci to Ci+1, the UAV flies first with the minimum radius r and then with a large maximum radius R. It should be noted that the new correction point is just the ending point B, when the trajectory errors ex and ez from the current starting point to ending point B are both less than θ.Figure 3
Calculation model of the trajectory length considering turning radius.In Figure3, points Oi and Oi+1 are the centers of the circle corresponding to radii r and R, respectively. Straight line OiQ is perpendicular to straight line CiCi+1. N is the point of intersection between straight lines OiOi+1 and CiCi+1, and P is the tangent point of the two circles. α is the angle between straight lines CiOi and CiCi+1, and β is the angle between straight lines OiOi+1 and OiQ.By the geometrical relationship shown in Figure3, the following expressions can be obtained:(1)γ=arccosCiCi−1⟶⋅CiCi+1⟶CiCi−1⟶CiCi+1⟶−π2,(2)CiCi+1⟶=CiQ⟶+QCi+1⟶=rcosγ+R+rsinψ,(3)R=Oi+1N⟶cosψ=R+r−OiN⟶cosβ=R+rcosψ−rsinγ⇒,R=rcosψ−rsinγ1−cosψ.By substituting equations (3) into (2), the calculation formula of the CiCi+1⟶ can be given as(4)CiCi+1⟶=rcosγ+1−sinγ1−cosψsinψ.According to equation (3), ψ=arccosR+rsinγ/R+r, which manifests that angle ψ belongs to the first quadrant. Since equation (3) is an algebraic equation, it is not hard to get angle ψ by solving it.Next, by using the calculatedψ and relationship between the arc length and central angle, the trajectory length of a UAV from the point Ci to Ci+1 can be written as(5)CiCi+1⌢=rπ2−γ+ψ+rcosψ−sinγ1−cosψψ,where CiCi+1⌢ denotes the sum of the arc lengths CiP⌢ and PCi+1⌢.


#### 3.1.2. Optimization Model of the Maximum Spherical Surface

To find a best correction point at each iteration, a set of correction points satisfying constraints can be determined first. The optimization model of the maximum spherical surface is established by considering constraints (6) and (7). Because the optimization model is preliminarily to obtain a set of correction points, accordingly, it will neglect the influence of the turning radius, but the optimal correction point model (described in Section 3.1.3) considering the weighted distance will not.The optimization model of the maximum spherical surface, including two types of constraint functions caused by the vertical and horizontal correction points, is given as follows:(1)
Define design variable(6)X=x1,x2,x3T.(2)
Determine objection function(7)minfX=−CiX⟶, where Ci is the starting point at iteration i (i.e., the best correction point searched at iteration i − 1), and CiX⟶ denotes the length of straight line from the point Ci to X whose coordinates in the three-dimensional space are x1, x2, and x3.(3)
Define two constraint functions according to the types of correction points.The constraints (constraint (6) mentioned in Section 2) of the maximum errors allowed to be vertically corrected are given by(8)g1X=ezi+CiX⟶δ−α1≤0,g2X=exi+CiX⟶δ−α2≤0,where ezi and exi denote the vertical and horizontal errors (i.e., the errors after being corrected at iteration i − 1) of the starting point at the iteration i, respectively.The constraints (constraint (7) mentioned in Section 2) of the maximum errors allowed to be horizontally corrected are given by(9)g3X=ezi+CiX⟶δ−β1≤0,g4X=exi+CiX⟶δ−β2≤0.The maximum radii areRzi and Rxi by solving the optimization model above with the two constraint functions. Then, we will select a set of correction points belonging to the maximum radius from Rzi and Rxi so that each flight distance of a UAV is farther before its flight errors are corrected as shown in Figures 4 and 5. Obviously, the selection method can ensure that the correction points UAV passes will be as few as possible. In the two figures, the green and blue points represent vertical and horizontal correction points, respectively. After determining the type of correction points and judging whether the distance from the point Ci to a correction point is less than the chosen maximum radius, a set of correction points can be obtained. In the next section, the set of correction points will be used for searching the best correction point based on the optimal correction point model considering the weighted distance.Figure 4
Vertical correction points allowed to be corrected.Figure 5
Horizontal correction points allowed to be corrected.


#### 3.1.3. Optimal Correction Point Model considering the Weighted Distance

The weighted distance model is shown in Figure6 where point Ci+1 is the best point being searched, and point D is the foot of the perpendicular from the point Ci+1 to straight line CiB. In Figure 6, point Xxi (x=1,2,⋯,n) denotes one of the correction points, and n denotes the total number of the correction points selected by the optimization model of the maximum spherical surface at iteration i.Figure 6
The weighted distance model.To realize the optimization objectives that trajectory length is as short as possible and that the correction points a UAV passes are as few as possible, on the one hand, the flight direction should be toward the ending point so that the trajectory length shortens, that is, the distanceXxiD⟶ is as small as possible. That the distance XxiB⟶ shortens, on the other hand, will result in the decrease of total correction points a UAV passes. Therefore, by taking into account the weight of distances XxiD⟶ and XxiB⟶, the optimal correction point model can be established as follows:(1)
Define design variablex(10)x=1,2,⋯,n.(2)
Design objective function(11)minfx=BXxi⟶2BXxi⟶+DXxi⟶+ωDXxi⟶2BXxi⟶+DXxi⟶,whereXxi is one of the correction points selected by the optimization model of the maximum spherical surface, BXxi⟶2/BXxi⟶+DXxi⟶ and DXxi⟶2/BXxi⟶+DXxi⟶ are, respectively, the normalized distances, and ω ≥ 0 is the weight.(3)
Define two types of constraint functions.Vertical correction constraints(12)g5X=ezi+Ci⌢Xxi⋅δ−α1≤0,g6X=exi+Ci⌢Xxi⋅δ−α2≤0.Horizontal correction constraints(13)g7X=ezi+Ci⌢Xxi⋅δ−β1≤0,g8X=exi+Ci⌢Xxi⋅δ−β2≤0.In equations (12) and (13), CiXxi⌢, similar to CiCi+1⌢, denotes the sum of two arc lengths shown in Figure 3 and can be computed by equation (5). In fact, the type of constraint functions is determined after selecting the maximum radius from Rzi and Rxi.


#### 3.1.4. Model of Updating Errors

After obtaining the best correction pointCi+1, the current starting point will be replaced by it as a new starting point, and a UAV’s flight errors will be updated for the next iteration. The handling method assumes that error corrections all fail by those correction points with failure probability. Accordingly, updating errors also needs to consider the failure label of the point Ci+1 besides its type, and the formulas for updating errors can be expressed as follows.Case 1:Ci+1 is a vertical correction point(14)ezi+1=τ⋅minezi+Ci⌢Ci+1⋅δ,ζ,exi+1=exi+CiCi+1⌢⋅δ, where ezi+1 and exi+1 denote the vertical and horizontal errors of the starting point at iteration i + 1, respectively, and τ—either “0 or 1”—represents failure label of point Ci+1.Case 2:Ci+1 is a horizontal correction point(15)ezi+1=ezi+Ci⌢Ci+1⋅δ,exi+1=τ⋅minexi+Ci⌢Ci+1⋅δ,ζ.


### 3.2. Modeling Based on Monte Carlo Sampling

#### 3.2.1. Motivation

Generally, using the definite rule by minimizing the designed objective functions can obtain the optimal trajectory whose length and searched points of the planned trajectory are as small as possible. However, affected by the distribution and number of correction points, any feasible trajectory cannot be found based on the definite rule, and the simulation example in Section5 will give the illustration. Actually, any correction point searched at a certain iteration has a great influence on the next iteration, such as the number, type, failure label of correction points allowed, and trajectory feasibility (i.e., whether a UAV flying along the trajectory can arrive at the destination successfully). Even the total number of iterations is uncertain. So, it is difficult to programing by considering those preferable correction points based on traditional methods—say exhaustive search or Dijkstra algorithm [29, 30]. To address this issue, in this section, we will propose a modeling method based on Monte Carlo Sampling for planning trajectory with failure tolerance. This method can be regarded as an improvement of the definite rule.Monte Carlo method [31–34] is a numerical simulation based on probability theory and mathematical statistics, which can solve many problems difficult to traditional methods. With the development of modern computation technology, the statistical simulation method has become an effective tool to solve or study high-dimensional complex problems. By a large number of samplings, it can get close to definite results, such as an estimation of statistical parameters of response variables [35], theoretical failure probability [36], and posterior distribution [37]. In other words, in this paper, it is nearly equivalent to forming a new definite rule which uses those preferable correction points at each iteration. By the new definite rule, a large number of feasible trajectories can be obtained.


#### 3.2.2. Model of Extracting Preferable Correction Points

To avoid that no correction point can be found based on the definite rule, we use not only the best correction point but also others at each iteration. In this way, the role of weightω mentioned in equation (11) will be weakened. Accordingly, the distance XxiB⟶ can be directly taken as the objective function for extracting preferable correction points by sorting the distances XxiB⟶ (i = 1, 2, 3, …, m). m denotes the number of correction points satisfying the constraints from equations (12) and (13).The model of extracting preferable correction points is shown in Figure7, and it extracts k correction points from m by distances XxiB⟶ (i = 1, 2, 3, …, m) from small to large. It can be regarded as an improvement for the optimal correction point model mentioned in Section 3.1.3, and one of the k points will be randomly selected based on Monte Carlo sampling.Figure 7
Model of extracting preferable correction points.It should be noted: ifm < k, in the next step, one point will be randomly selected from such correction points less than k; if m = 0 and the trajectory errors (from the updated starting point to the ending point) do not satisfy constraint θ, a UAV cannot successfully arrive at the destination when along the trajectory, so it is necessary to automatically return the initial states to restart.


## 3.1. Modeling Based on the Definite Rule

### 3.1.1. Calculation Model of the Trajectory Length considering Turning Radius

The calculation model considering turning radius is shown in Figure3, where Ci−1, Ci, and Ci+1 denote the previous starting point, the current starting point (i.e., the best correction point searched at the previous iteration), and a new correction point searched, respectively. It can be known that the following flight trajectory between starting point to the correction point can ensure the flight distance is small, when the UAV flies into point Ci+1, its flight direction is just tangent to straight line CiCi+1; namely, during the flying from the point Ci to Ci+1, the UAV flies first with the minimum radius r and then with a large maximum radius R. It should be noted that the new correction point is just the ending point B, when the trajectory errors ex and ez from the current starting point to ending point B are both less than θ.Figure 3
Calculation model of the trajectory length considering turning radius.In Figure3, points Oi and Oi+1 are the centers of the circle corresponding to radii r and R, respectively. Straight line OiQ is perpendicular to straight line CiCi+1. N is the point of intersection between straight lines OiOi+1 and CiCi+1, and P is the tangent point of the two circles. α is the angle between straight lines CiOi and CiCi+1, and β is the angle between straight lines OiOi+1 and OiQ.By the geometrical relationship shown in Figure3, the following expressions can be obtained:(1)γ=arccosCiCi−1⟶⋅CiCi+1⟶CiCi−1⟶CiCi+1⟶−π2,(2)CiCi+1⟶=CiQ⟶+QCi+1⟶=rcosγ+R+rsinψ,(3)R=Oi+1N⟶cosψ=R+r−OiN⟶cosβ=R+rcosψ−rsinγ⇒,R=rcosψ−rsinγ1−cosψ.By substituting equations (3) into (2), the calculation formula of the CiCi+1⟶ can be given as(4)CiCi+1⟶=rcosγ+1−sinγ1−cosψsinψ.According to equation (3), ψ=arccosR+rsinγ/R+r, which manifests that angle ψ belongs to the first quadrant. Since equation (3) is an algebraic equation, it is not hard to get angle ψ by solving it.Next, by using the calculatedψ and relationship between the arc length and central angle, the trajectory length of a UAV from the point Ci to Ci+1 can be written as(5)CiCi+1⌢=rπ2−γ+ψ+rcosψ−sinγ1−cosψψ,where CiCi+1⌢ denotes the sum of the arc lengths CiP⌢ and PCi+1⌢.


### 3.1.2. Optimization Model of the Maximum Spherical Surface

To find a best correction point at each iteration, a set of correction points satisfying constraints can be determined first. The optimization model of the maximum spherical surface is established by considering constraints (6) and (7). Because the optimization model is preliminarily to obtain a set of correction points, accordingly, it will neglect the influence of the turning radius, but the optimal correction point model (described in Section 3.1.3) considering the weighted distance will not.The optimization model of the maximum spherical surface, including two types of constraint functions caused by the vertical and horizontal correction points, is given as follows:(1)
Define design variable(6)X=x1,x2,x3T.(2)
Determine objection function(7)minfX=−CiX⟶, where Ci is the starting point at iteration i (i.e., the best correction point searched at iteration i − 1), and CiX⟶ denotes the length of straight line from the point Ci to X whose coordinates in the three-dimensional space are x1, x2, and x3.(3)
Define two constraint functions according to the types of correction points.The constraints (constraint (6) mentioned in Section 2) of the maximum errors allowed to be vertically corrected are given by(8)g1X=ezi+CiX⟶δ−α1≤0,g2X=exi+CiX⟶δ−α2≤0,where ezi and exi denote the vertical and horizontal errors (i.e., the errors after being corrected at iteration i − 1) of the starting point at the iteration i, respectively.The constraints (constraint (7) mentioned in Section 2) of the maximum errors allowed to be horizontally corrected are given by(9)g3X=ezi+CiX⟶δ−β1≤0,g4X=exi+CiX⟶δ−β2≤0.The maximum radii areRzi and Rxi by solving the optimization model above with the two constraint functions. Then, we will select a set of correction points belonging to the maximum radius from Rzi and Rxi so that each flight distance of a UAV is farther before its flight errors are corrected as shown in Figures 4 and 5. Obviously, the selection method can ensure that the correction points UAV passes will be as few as possible. In the two figures, the green and blue points represent vertical and horizontal correction points, respectively. After determining the type of correction points and judging whether the distance from the point Ci to a correction point is less than the chosen maximum radius, a set of correction points can be obtained. In the next section, the set of correction points will be used for searching the best correction point based on the optimal correction point model considering the weighted distance.Figure 4
Vertical correction points allowed to be corrected.Figure 5
Horizontal correction points allowed to be corrected.


### 3.1.3. Optimal Correction Point Model considering the Weighted Distance

The weighted distance model is shown in Figure6 where point Ci+1 is the best point being searched, and point D is the foot of the perpendicular from the point Ci+1 to straight line CiB. In Figure 6, point Xxi (x=1,2,⋯,n) denotes one of the correction points, and n denotes the total number of the correction points selected by the optimization model of the maximum spherical surface at iteration i.Figure 6
The weighted distance model.To realize the optimization objectives that trajectory length is as short as possible and that the correction points a UAV passes are as few as possible, on the one hand, the flight direction should be toward the ending point so that the trajectory length shortens, that is, the distanceXxiD⟶ is as small as possible. That the distance XxiB⟶ shortens, on the other hand, will result in the decrease of total correction points a UAV passes. Therefore, by taking into account the weight of distances XxiD⟶ and XxiB⟶, the optimal correction point model can be established as follows:(1)
Define design variablex(10)x=1,2,⋯,n.(2)
Design objective function(11)minfx=BXxi⟶2BXxi⟶+DXxi⟶+ωDXxi⟶2BXxi⟶+DXxi⟶,whereXxi is one of the correction points selected by the optimization model of the maximum spherical surface, BXxi⟶2/BXxi⟶+DXxi⟶ and DXxi⟶2/BXxi⟶+DXxi⟶ are, respectively, the normalized distances, and ω ≥ 0 is the weight.(3)
Define two types of constraint functions.Vertical correction constraints(12)g5X=ezi+Ci⌢Xxi⋅δ−α1≤0,g6X=exi+Ci⌢Xxi⋅δ−α2≤0.Horizontal correction constraints(13)g7X=ezi+Ci⌢Xxi⋅δ−β1≤0,g8X=exi+Ci⌢Xxi⋅δ−β2≤0.In equations (12) and (13), CiXxi⌢, similar to CiCi+1⌢, denotes the sum of two arc lengths shown in Figure 3 and can be computed by equation (5). In fact, the type of constraint functions is determined after selecting the maximum radius from Rzi and Rxi.


### 3.1.4. Model of Updating Errors

After obtaining the best correction pointCi+1, the current starting point will be replaced by it as a new starting point, and a UAV’s flight errors will be updated for the next iteration. The handling method assumes that error corrections all fail by those correction points with failure probability. Accordingly, updating errors also needs to consider the failure label of the point Ci+1 besides its type, and the formulas for updating errors can be expressed as follows.Case 1:Ci+1 is a vertical correction point(14)ezi+1=τ⋅minezi+Ci⌢Ci+1⋅δ,ζ,exi+1=exi+CiCi+1⌢⋅δ, where ezi+1 and exi+1 denote the vertical and horizontal errors of the starting point at iteration i + 1, respectively, and τ—either “0 or 1”—represents failure label of point Ci+1.Case 2:Ci+1 is a horizontal correction point(15)ezi+1=ezi+Ci⌢Ci+1⋅δ,exi+1=τ⋅minexi+Ci⌢Ci+1⋅δ,ζ.


## 3.1.1. Calculation Model of the Trajectory Length considering Turning Radius

The calculation model considering turning radius is shown in Figure3, where Ci−1, Ci, and Ci+1 denote the previous starting point, the current starting point (i.e., the best correction point searched at the previous iteration), and a new correction point searched, respectively. It can be known that the following flight trajectory between starting point to the correction point can ensure the flight distance is small, when the UAV flies into point Ci+1, its flight direction is just tangent to straight line CiCi+1; namely, during the flying from the point Ci to Ci+1, the UAV flies first with the minimum radius r and then with a large maximum radius R. It should be noted that the new correction point is just the ending point B, when the trajectory errors ex and ez from the current starting point to ending point B are both less than θ.Figure 3
Calculation model of the trajectory length considering turning radius.In Figure3, points Oi and Oi+1 are the centers of the circle corresponding to radii r and R, respectively. Straight line OiQ is perpendicular to straight line CiCi+1. N is the point of intersection between straight lines OiOi+1 and CiCi+1, and P is the tangent point of the two circles. α is the angle between straight lines CiOi and CiCi+1, and β is the angle between straight lines OiOi+1 and OiQ.By the geometrical relationship shown in Figure3, the following expressions can be obtained:(1)γ=arccosCiCi−1⟶⋅CiCi+1⟶CiCi−1⟶CiCi+1⟶−π2,(2)CiCi+1⟶=CiQ⟶+QCi+1⟶=rcosγ+R+rsinψ,(3)R=Oi+1N⟶cosψ=R+r−OiN⟶cosβ=R+rcosψ−rsinγ⇒,R=rcosψ−rsinγ1−cosψ.By substituting equations (3) into (2), the calculation formula of the CiCi+1⟶ can be given as(4)CiCi+1⟶=rcosγ+1−sinγ1−cosψsinψ.According to equation (3), ψ=arccosR+rsinγ/R+r, which manifests that angle ψ belongs to the first quadrant. Since equation (3) is an algebraic equation, it is not hard to get angle ψ by solving it.Next, by using the calculatedψ and relationship between the arc length and central angle, the trajectory length of a UAV from the point Ci to Ci+1 can be written as(5)CiCi+1⌢=rπ2−γ+ψ+rcosψ−sinγ1−cosψψ,where CiCi+1⌢ denotes the sum of the arc lengths CiP⌢ and PCi+1⌢.


## 3.1.2. Optimization Model of the Maximum Spherical Surface

To find a best correction point at each iteration, a set of correction points satisfying constraints can be determined first. The optimization model of the maximum spherical surface is established by considering constraints (6) and (7). Because the optimization model is preliminarily to obtain a set of correction points, accordingly, it will neglect the influence of the turning radius, but the optimal correction point model (described in Section 3.1.3) considering the weighted distance will not.The optimization model of the maximum spherical surface, including two types of constraint functions caused by the vertical and horizontal correction points, is given as follows:(1)
Define design variable(6)X=x1,x2,x3T.(2)
Determine objection function(7)minfX=−CiX⟶, where Ci is the starting point at iteration i (i.e., the best correction point searched at iteration i − 1), and CiX⟶ denotes the length of straight line from the point Ci to X whose coordinates in the three-dimensional space are x1, x2, and x3.(3)
Define two constraint functions according to the types of correction points.The constraints (constraint (6) mentioned in Section 2) of the maximum errors allowed to be vertically corrected are given by(8)g1X=ezi+CiX⟶δ−α1≤0,g2X=exi+CiX⟶δ−α2≤0,where ezi and exi denote the vertical and horizontal errors (i.e., the errors after being corrected at iteration i − 1) of the starting point at the iteration i, respectively.The constraints (constraint (7) mentioned in Section 2) of the maximum errors allowed to be horizontally corrected are given by(9)g3X=ezi+CiX⟶δ−β1≤0,g4X=exi+CiX⟶δ−β2≤0.The maximum radii areRzi and Rxi by solving the optimization model above with the two constraint functions. Then, we will select a set of correction points belonging to the maximum radius from Rzi and Rxi so that each flight distance of a UAV is farther before its flight errors are corrected as shown in Figures 4 and 5. Obviously, the selection method can ensure that the correction points UAV passes will be as few as possible. In the two figures, the green and blue points represent vertical and horizontal correction points, respectively. After determining the type of correction points and judging whether the distance from the point Ci to a correction point is less than the chosen maximum radius, a set of correction points can be obtained. In the next section, the set of correction points will be used for searching the best correction point based on the optimal correction point model considering the weighted distance.Figure 4
Vertical correction points allowed to be corrected.Figure 5
Horizontal correction points allowed to be corrected.


## 3.1.3. Optimal Correction Point Model considering the Weighted Distance

The weighted distance model is shown in Figure6 where point Ci+1 is the best point being searched, and point D is the foot of the perpendicular from the point Ci+1 to straight line CiB. In Figure 6, point Xxi (x=1,2,⋯,n) denotes one of the correction points, and n denotes the total number of the correction points selected by the optimization model of the maximum spherical surface at iteration i.Figure 6
The weighted distance model.To realize the optimization objectives that trajectory length is as short as possible and that the correction points a UAV passes are as few as possible, on the one hand, the flight direction should be toward the ending point so that the trajectory length shortens, that is, the distanceXxiD⟶ is as small as possible. That the distance XxiB⟶ shortens, on the other hand, will result in the decrease of total correction points a UAV passes. Therefore, by taking into account the weight of distances XxiD⟶ and XxiB⟶, the optimal correction point model can be established as follows:(1)
Define design variablex(10)x=1,2,⋯,n.(2)
Design objective function(11)minfx=BXxi⟶2BXxi⟶+DXxi⟶+ωDXxi⟶2BXxi⟶+DXxi⟶,whereXxi is one of the correction points selected by the optimization model of the maximum spherical surface, BXxi⟶2/BXxi⟶+DXxi⟶ and DXxi⟶2/BXxi⟶+DXxi⟶ are, respectively, the normalized distances, and ω ≥ 0 is the weight.(3)
Define two types of constraint functions.Vertical correction constraints(12)g5X=ezi+Ci⌢Xxi⋅δ−α1≤0,g6X=exi+Ci⌢Xxi⋅δ−α2≤0.Horizontal correction constraints(13)g7X=ezi+Ci⌢Xxi⋅δ−β1≤0,g8X=exi+Ci⌢Xxi⋅δ−β2≤0.In equations (12) and (13), CiXxi⌢, similar to CiCi+1⌢, denotes the sum of two arc lengths shown in Figure 3 and can be computed by equation (5). In fact, the type of constraint functions is determined after selecting the maximum radius from Rzi and Rxi.


## 3.1.4. Model of Updating Errors

After obtaining the best correction pointCi+1, the current starting point will be replaced by it as a new starting point, and a UAV’s flight errors will be updated for the next iteration. The handling method assumes that error corrections all fail by those correction points with failure probability. Accordingly, updating errors also needs to consider the failure label of the point Ci+1 besides its type, and the formulas for updating errors can be expressed as follows.Case 1:Ci+1 is a vertical correction point(14)ezi+1=τ⋅minezi+Ci⌢Ci+1⋅δ,ζ,exi+1=exi+CiCi+1⌢⋅δ, where ezi+1 and exi+1 denote the vertical and horizontal errors of the starting point at iteration i + 1, respectively, and τ—either “0 or 1”—represents failure label of point Ci+1.Case 2:Ci+1 is a horizontal correction point(15)ezi+1=ezi+Ci⌢Ci+1⋅δ,exi+1=τ⋅minexi+Ci⌢Ci+1⋅δ,ζ.


## 3.2. Modeling Based on Monte Carlo Sampling

### 3.2.1. Motivation

Generally, using the definite rule by minimizing the designed objective functions can obtain the optimal trajectory whose length and searched points of the planned trajectory are as small as possible. However, affected by the distribution and number of correction points, any feasible trajectory cannot be found based on the definite rule, and the simulation example in Section5 will give the illustration. Actually, any correction point searched at a certain iteration has a great influence on the next iteration, such as the number, type, failure label of correction points allowed, and trajectory feasibility (i.e., whether a UAV flying along the trajectory can arrive at the destination successfully). Even the total number of iterations is uncertain. So, it is difficult to programing by considering those preferable correction points based on traditional methods—say exhaustive search or Dijkstra algorithm [29, 30]. To address this issue, in this section, we will propose a modeling method based on Monte Carlo Sampling for planning trajectory with failure tolerance. This method can be regarded as an improvement of the definite rule.Monte Carlo method [31–34] is a numerical simulation based on probability theory and mathematical statistics, which can solve many problems difficult to traditional methods. With the development of modern computation technology, the statistical simulation method has become an effective tool to solve or study high-dimensional complex problems. By a large number of samplings, it can get close to definite results, such as an estimation of statistical parameters of response variables [35], theoretical failure probability [36], and posterior distribution [37]. In other words, in this paper, it is nearly equivalent to forming a new definite rule which uses those preferable correction points at each iteration. By the new definite rule, a large number of feasible trajectories can be obtained.


### 3.2.2. Model of Extracting Preferable Correction Points

To avoid that no correction point can be found based on the definite rule, we use not only the best correction point but also others at each iteration. In this way, the role of weightω mentioned in equation (11) will be weakened. Accordingly, the distance XxiB⟶ can be directly taken as the objective function for extracting preferable correction points by sorting the distances XxiB⟶ (i = 1, 2, 3, …, m). m denotes the number of correction points satisfying the constraints from equations (12) and (13).The model of extracting preferable correction points is shown in Figure7, and it extracts k correction points from m by distances XxiB⟶ (i = 1, 2, 3, …, m) from small to large. It can be regarded as an improvement for the optimal correction point model mentioned in Section 3.1.3, and one of the k points will be randomly selected based on Monte Carlo sampling.Figure 7
Model of extracting preferable correction points.It should be noted: ifm < k, in the next step, one point will be randomly selected from such correction points less than k; if m = 0 and the trajectory errors (from the updated starting point to the ending point) do not satisfy constraint θ, a UAV cannot successfully arrive at the destination when along the trajectory, so it is necessary to automatically return the initial states to restart.


## 3.2.1. Motivation

Generally, using the definite rule by minimizing the designed objective functions can obtain the optimal trajectory whose length and searched points of the planned trajectory are as small as possible. However, affected by the distribution and number of correction points, any feasible trajectory cannot be found based on the definite rule, and the simulation example in Section5 will give the illustration. Actually, any correction point searched at a certain iteration has a great influence on the next iteration, such as the number, type, failure label of correction points allowed, and trajectory feasibility (i.e., whether a UAV flying along the trajectory can arrive at the destination successfully). Even the total number of iterations is uncertain. So, it is difficult to programing by considering those preferable correction points based on traditional methods—say exhaustive search or Dijkstra algorithm [29, 30]. To address this issue, in this section, we will propose a modeling method based on Monte Carlo Sampling for planning trajectory with failure tolerance. This method can be regarded as an improvement of the definite rule.Monte Carlo method [31–34] is a numerical simulation based on probability theory and mathematical statistics, which can solve many problems difficult to traditional methods. With the development of modern computation technology, the statistical simulation method has become an effective tool to solve or study high-dimensional complex problems. By a large number of samplings, it can get close to definite results, such as an estimation of statistical parameters of response variables [35], theoretical failure probability [36], and posterior distribution [37]. In other words, in this paper, it is nearly equivalent to forming a new definite rule which uses those preferable correction points at each iteration. By the new definite rule, a large number of feasible trajectories can be obtained.


## 3.2.2. Model of Extracting Preferable Correction Points

To avoid that no correction point can be found based on the definite rule, we use not only the best correction point but also others at each iteration. In this way, the role of weightω mentioned in equation (11) will be weakened. Accordingly, the distance XxiB⟶ can be directly taken as the objective function for extracting preferable correction points by sorting the distances XxiB⟶ (i = 1, 2, 3, …, m). m denotes the number of correction points satisfying the constraints from equations (12) and (13).The model of extracting preferable correction points is shown in Figure7, and it extracts k correction points from m by distances XxiB⟶ (i = 1, 2, 3, …, m) from small to large. It can be regarded as an improvement for the optimal correction point model mentioned in Section 3.1.3, and one of the k points will be randomly selected based on Monte Carlo sampling.Figure 7
Model of extracting preferable correction points.It should be noted: ifm < k, in the next step, one point will be randomly selected from such correction points less than k; if m = 0 and the trajectory errors (from the updated starting point to the ending point) do not satisfy constraint θ, a UAV cannot successfully arrive at the destination when along the trajectory, so it is necessary to automatically return the initial states to restart.


## 4. Algorithms and Implementation

### 4.1. Algorithms

To solve the optimization models, sequential quadratic programming (SQP) and the exhaustive search will be adopted, mainly based on the following reasons.(1)
The optimization model of the maximum spherical surface, whose objective and constraint functions are quadratic, belongs to a quadratic programming problem. On the one hand, the trajectory planning considering error correction involves many optimizations at each iteration. That is to say, the computational cost should be given priority; otherwise, the real-time performance cannot be ensured in practice. On the other hand, our main work is to emphasize the modeling of trajectory planning considering error correction, so solving it can directly call the SQP algorithm provided by MATLAB. As described in the literature [38–41], SQP, being a high-efficiency method, is considered to be one of the most efficient methods to solve the nonlinearly constrained optimization problems, particularly suitable for the quadratic programming problem. Accordingly, “the SQP method will be adopted to solve the optimization model of the maximum spherical surface.”(2)
The optimal correction point model considering the weighted distance will employ the exhaustive search algorithm by programming based on the MATLAB because its feasible points are very finite after a set of correction points obtained by the optimization model of the maximum spherical surface. By using the barbaric algorithm, the best correction point can be found under an assigned weight.


### 4.2. Implementation

The flowchart of trajectory planning methods with fault tolerance is shown in Figure8, and the implementation steps are given as follows.Step 1. Obtain the starting and ending points, as well as all of the correction points in the input data. That means acquiring their information, such as their coordinates in the three-dimensional space, types, and failure labels.Step 2. Judge whether the flight errors (which can be obtained by the calculation model of the trajectory length considering the turning radius) from the starting point to the ending point satisfy the conditionθ and determine whether to continue the next iteration.Step 3. Determine the type of correction points from the two radii by selecting the maximum one and obtain a set of the correction points within the spherical surface of the selected radius, according to the optimization model of the maximum spherical surface with its two constraint functions.Step 4. Perform the following steps separately. Cases 1 and 2 represent the trajectory planning methods based on the definite rule and Monte Carlo sampling, respectively.Case 1: Obtain the best correction point from the set of the correction points, according to the optimal correction point model considering the weighted distance.Case 2: Extractk correction points considering the distances BXxi⟶ at each iteration.Step 5.Case 1: Take the best correction point as a new starting point, update the flight errors after the UAV is corrected, and continue the next iteration.Case 2: Monte Carlo sampling. Randomly select one from thek correction points (preferable correction points) as a new starting point, update the flight error after the UAV is corrected, and continue the next iteration.Step 6. Repeat Steps 1–5, until Step 2 not to search next correction point. That is, the vertical and horizontal errors caused by the flight distance between the new starting point and ending point B are both less thanθ m.Figure 8
Flowchart of the trajectory planning methods with fault tolerance.


## 4.1. Algorithms

To solve the optimization models, sequential quadratic programming (SQP) and the exhaustive search will be adopted, mainly based on the following reasons.(1)
The optimization model of the maximum spherical surface, whose objective and constraint functions are quadratic, belongs to a quadratic programming problem. On the one hand, the trajectory planning considering error correction involves many optimizations at each iteration. That is to say, the computational cost should be given priority; otherwise, the real-time performance cannot be ensured in practice. On the other hand, our main work is to emphasize the modeling of trajectory planning considering error correction, so solving it can directly call the SQP algorithm provided by MATLAB. As described in the literature [38–41], SQP, being a high-efficiency method, is considered to be one of the most efficient methods to solve the nonlinearly constrained optimization problems, particularly suitable for the quadratic programming problem. Accordingly, “the SQP method will be adopted to solve the optimization model of the maximum spherical surface.”(2)
The optimal correction point model considering the weighted distance will employ the exhaustive search algorithm by programming based on the MATLAB because its feasible points are very finite after a set of correction points obtained by the optimization model of the maximum spherical surface. By using the barbaric algorithm, the best correction point can be found under an assigned weight.


## 4.2. Implementation

The flowchart of trajectory planning methods with fault tolerance is shown in Figure8, and the implementation steps are given as follows.Step 1. Obtain the starting and ending points, as well as all of the correction points in the input data. That means acquiring their information, such as their coordinates in the three-dimensional space, types, and failure labels.Step 2. Judge whether the flight errors (which can be obtained by the calculation model of the trajectory length considering the turning radius) from the starting point to the ending point satisfy the conditionθ and determine whether to continue the next iteration.Step 3. Determine the type of correction points from the two radii by selecting the maximum one and obtain a set of the correction points within the spherical surface of the selected radius, according to the optimization model of the maximum spherical surface with its two constraint functions.Step 4. Perform the following steps separately. Cases 1 and 2 represent the trajectory planning methods based on the definite rule and Monte Carlo sampling, respectively.Case 1: Obtain the best correction point from the set of the correction points, according to the optimal correction point model considering the weighted distance.Case 2: Extractk correction points considering the distances BXxi⟶ at each iteration.Step 5.Case 1: Take the best correction point as a new starting point, update the flight errors after the UAV is corrected, and continue the next iteration.Case 2: Monte Carlo sampling. Randomly select one from thek correction points (preferable correction points) as a new starting point, update the flight error after the UAV is corrected, and continue the next iteration.Step 6. Repeat Steps 1–5, until Step 2 not to search next correction point. That is, the vertical and horizontal errors caused by the flight distance between the new starting point and ending point B are both less thanθ m.Figure 8
Flowchart of the trajectory planning methods with fault tolerance.


## 5. Simulation Example

In this section, we apply the trajectory planning methods proposed in Section3 to solve a practical problem. This problem is very typical because the number of correction points is small and the requirements are very strict, which may lead to the failure of the definite rule method.In MATLAB 2018a programming environment, a set of solution programs is developed based on the trajectory planning methods, and the optimal trajectory is obtained. It is noted that the developed program is not limited to this problem to be solved, if these data of a flight space change, and the new results of trajectory planning for UAVs can be obtained after importing the new data.


### 5.1. Problem to Be Solved

Table1 shows the partial data of a flight space, which includes the coordinates of the starting point, ending point, and correction points. The serial numbers 0 and 326 are the starting point and ending point, respectively, and the rest are correction points. In the column containing the type of correction points, “0” and “1” represent the horizontal and vertical error correction points, respectively. In the last column, “0” and “1” represent the correction points with and without the failure probability, respectively; the probability of a successful correction is only σ when the flight error is corrected by those having the failure probability.Table 1
Partial data of a flight space [42].
Serial numbersCoordinatesType of correction pointsFailure labelX (m)Y (m)Z (m)Starting point A00.0050000.005000.00——Correction points176009.859788.119121.901022448.2071599.881877.1300327800.9315218.078345.4100450056.7391668.091465.0900525311.1551225.791330.3501………………32478846.7073662.653906.981032563645.5671435.226263.1301Ending point B326100000.0074860.555499.61——For this problem, the parameters mentioned in Section2 are given as follows: α1 = 20 m, α2 = 10 m, β1 = 15 m, β2 = 20 m, θ = 20 m, δ = 0.001 m, σ = 80%, ζ = 5 m, and r = 200 m.


### 5.2. Results and Discussion Based on the Definite Rule

After considering that certain correction points have the failure probability and assuming that all corrections fail by these points at each iteration, the results can be obtained by running the programs developed in this paper. However, the simulation results show that, no matter which weight we select, the serial numbers of the correction points searched are all [163, termination]. This is a special phenomenon, which presents that after serial number 163 is searched, solving the optimal correction point model cannot find any correction point at the next iteration. This causes the cease of the programs unexpectedly. In other words, if there exists the failure probability on certain correction points, no feasible trajectory can be found based on the definite method for the practical problem whose distributions on correction points are relatively rigorous. “It is for this reason that we further propose a trajectory planning method based on Monte Carlo Sampling,” and the simulation results are as shown in Section5.3.To illustrate the simulation process and results based on the definite rule, the restriction of the correction points for the problem will be loosened, and it is assumed that all correction points will not fail because most problems fall into this situation. In other words, in equations (14) and (15), τ will remain 0, and the corresponding results are as follows.


#### 5.2.1. Optimized Trajectories under the Different Weights

Because correction points are finite and not uniformly distributed in the whole space, weightω from 0 to infinite can be discretely picked to obtain optimized trajectories as shown in Table 2.Table 2
Optimized trajectory under the different weights.
WeightωTotal number of the best correction pointsTrajectory length (km)Serial numbers of the points searched0∼0.2——[163, 114, 8, 309, 255, 148, 90, 156, 320, termination]0.3∼0.514137.7666[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 187, 61, 292, 326]0.6∼0.714136.0649[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 323, 61, 292, 326]0.8∼314135.0177[163, 114, 8, 309, 255, 148, 90, 156, 33, 289, 50, 323, 61, 292, 326]4∼514127.5342[163, 114, 8, 309, 121, 123, 49, 160, 138, 19, 50, 323, 61, 292, 326]6∼2512111.7787[163, 114, 8, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]27∼2915123.6882[163, 114, 234, 188, 222, 230, 225, 170, 123, 49, 160, 92, 93, 61, 292, 326]30∼30.515122.4127[163, 114, 234, 188, 238, 227, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]31∼4816120.0604[163, 114, 234, 188, 238, 227, 309, 305, 123, 45, 160, 92, 93, 61, 292, 326]>49No solution[163, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, …]If weightω is very small, the best correction points searched mainly depend on the distance BXxi⟶. To minimize the objective function (i.e., equation (11)), a correction point farthest from the starting point will be found. In other words, a UAV will fly a longer distance at each iteration, which may cause that the optimal correction point model cannot find any correction point because flight errors are so large that they cannot satisfy the constraints from equations (12) and (13) at next iteration. When weight ω belongs to [0, 0.2], the series number of the correction points searched just illustrates the point as shown in Table 2. It should be noted that the “termination” shown in Table 2 presents that after serial number 320 is searched, solving the optimal correction point model cannot find any correction point at the next iteration, which causes the cease of the programs unexpectedly.With the increase of weightω, the importance of BXxi⟶ shown in Figure 6 is gradually weakened. Under a larger weight ω, it is easy to understand that if there are three or more correction points (satisfying constraints (6) and (7)) nearby the straight line between the starting and ending points, the optimal correction point model may always search the same points so that running the programs will fall into an endless loop. When ω > 49, the series number of the correction points searched just illustrates the point as shown in Table 2.It can be seen from Table2 that, the optimal trajectory is obtained when weight ω belongs to [6, 25]; in this case, trajectory length is the shortest, and the correction points UAV passes are very few. If there is no failure of all the correction points, therefore, the total length of the optimal trajectory is 111.7787 km, and the total number of the best correction points searched is 12.


#### 5.2.2. Optimal Trajectory

Table3 shows the flight data of the UAV along the optimal trajectory, which includes the correction pints, vertical and horizontal errors before being corrected, type of correction points, and flight distance from the starting point updated continually to the best correction point searched at each iteration, respectively.Table 3
Flight data of the UAV along the optimal trajectory.
Serial numbersVertical error before being correctedHorizontal error before being correctedType of correction pointsFlight distance of each segment (km)Starting point A000——Best correction points searched16313.2878913.2879013.287911418.62955.341615.346813.957319.2989013.957330919.52155.564215.564212111.252316.8166011.252312316.66335.411015.41104911.79817.2090011.798016018.32196.523916.5239925.787812.311705.78789315.2739.485219.4852619.834719.3198809.834729216.3896.554316.5543Ending point B3266.980513.5348—6.9805Total trajectory length111.7787Figure9 shows the optimal trajectory in the space. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 9 and 10. The two figures clearly show the locations of the best correction points searched at each iteration. It is noted that in reality, the planned trajectory is a curve; however, Figures 9–11, the purpose of which is to describe the locations of the best correction points and to illustrate the process of searching correction points at each iteration, hence use a series of straight lines to express the optimal trajectory. Actually, owing to the small turning radius, it is not apparent to express the optimal trajectory with a curve.Figure 9
Optimal trajectory without the failure probability in the space.Figure 10
Optimal trajectory without the failure probability in the XOY plane.Figure 11
Optimal trajectory without the failure probability in the XOZ plane.It can be seen from Table3 that the two types (“0 and 1”) of the best correction points searched appear alternately, and the flight distance of each segment appears with the long and short one. That is coincident with the principle that the best correction points are as few as possible. In other words, if a best correction point searched is the vertical one, the next must be the horizontal one in order to maximize the flight distance before being corrected.Figure9 clearly shows that the flight trajectory is near the straight line AB when its three axes have the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory whose total length is only 111.7787 km illustrates that it is quite short. As shown in Table 3, all the flight errors are limited to the allowed range of the vertical and horizontal correction before being corrected. So, few are the best correction points searched, additionally, that it is sure that the optimal trajectory of the UAV is reasonable.Furthermore, some correction points, which satisfy all constraints, can be randomly selected to generate many feasible trajectories; it can be verified that, regardless of the trajectory length or number of correction points, they both are greater than the optimal trajectory’s counterparts.From what has been analyzed and discussed in the three paragraphs above, it is worth believing that the results obtained by the trajectory planning method based on the definite rule are feasible and reasonable, and this work provides an effective solution for the trajectory planning of UAVs considering the error correction when there is no failure on the correction points.


### 5.3. Results and Discussions Based on Monte Carlo Sampling

In this section, the trajectory planning method based on Monte Carlo sampling will be adopted, considering that there is a failure on certain correction points, and the error corrections all fail by these points. By repeatedly running programs whose flowchart is shown in Figure8, a series of feasible trajectories can be obtained based on Monte Carlo sampling. The optimal trajectory can be selected from them by considering the trajectory length and number of correction points a UAV passes.


#### 5.3.1. Feasible Trajectories

Some feasible trajectories obtained by Monte Carlo sampling are shown in Table4; in fact, many others are not listed due to the large trajectory lengths.Table 4
Feasible trajectories obtained by Monte Carlo sampling.
Number of correction points searchedTrajectory length (km)Serial numbers of points searched (including the ending point 326)122184.0919[252, 266, 270, 89, 236, 133, 3, 55, 144, 69, 298, 117, 85, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]223191.0631[252, 266, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 9, 253, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]325212.61575[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 167, 249, 274, 211, 321, 16, 282, 141, 291, 161, 326]426221.16501[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 73, 207, 70, 211, 321, 16, 282, 141, 291, 84, 124, 326]527226.8188[169, 266, 270, 89, 236, 185, 3, 281, 197, 112, 103, 250, 243, 73, 249, 274, 16, 186, 233, 120, 202, 71, 141, 282, 84, 287, 99, 326]628223.2609[169, 266, 100, 146, 270, 89, 236, 133, 3, 281, 144, 55, 197, 112, 143, 250, 9, 253, 249, 274, 12, 219, 16, 101, 141, 291, 84, 124, 326]729229.7933[169, 266, 270, 89, 236, 133, 3, 281, 144, 55, 69, 298, 117, 85, 76, 6, 211, 44, 191, 160, 324, 19, 50, 177, 87, 246, 313, 210, 212, 326]830233.0981[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 243, 73, 39, 274, 51, 211, 321, 16, 71, 233, 101, 141, 291, 84, 59, 326]931248.6503[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 55, 144, 261, 298, 143, 250, 86, 119, 243, 167, 249, 274, 12, 219, 16, 71, 141, 291, 84, 287, 61, 326]1032259.0214[169, 157, 252, 266, 270, 89, 236, 185, 78, 281, 3, 14, 214, 179, 111, 144, 55, 69, 298, 143, 250, 243, 167, 249, 274, 12, 321, 16, 282, 84, 291, 161, 326]It can be seen from Table4 that the more the correction points searched are, the larger the trajectory length becomes. Apparently, “trajectory 1” is the optimal trajectory; in this case, trajectory length is the shortest, and the correction points the UAV passes are very few. Considering that there exists the failure probability on certain correction points, therefore, the optimal trajectory length is 184.0919 km and its total number of correction points searched is 22. When a UAV flies along the optimal trajectory, even if error corrections fail, it can still arrive at the destination successfully. “It means that the trajectory planning method based on Monte Carlo sampling has the ability of fault tolerance.”


#### 5.3.2. Optimal Trajectory

Table5 shows the flight data of a UAV along the optimal trajectory with the failure probability on certain correction points. The flight data include the correction points, vertical and horizontal errors before and after being corrected, type of correction points, failure label, and flight distance from the starting point updated continually to the point searched at each iteration.Table 5
Flight data of a UAV along the optimal trajectory with the failure probability on certain correction points.
PointsSerial numbersVertical error before being correctedHorizontal error before being correctedVertical error after being correctedHorizontal error after being correctedType of correction pointsFailure labelFlight distance of each segment (km)A00000———Correction points searched25211.107611.107611.10760.00000011.107626619.10717.99950.00007.9995107.99952706.437014.43656.43700.0000006.43708914.08987.65280.00007.6528107.652823610.137917.790710.13790.00000010.137913312.77042.63250.00002.6325102.632538.838311.47088.83830.0000008.83835518.02079.18240.00009.1824109.18241442.716711.89922.71670.0000002.71676914.589611.872914.58960.00000011.872929819.67395.08430.00005.0843105.08431179.905514.98989.90550.0000009.90558515.54455.63900.00005.6390105.63907614.323519.962514.32350.00000014.3235619.54905.22550.00005.2255105.22551210.702715.928210.70270.00000010.702721918.58867.88590.00007.8859107.8859167.916815.80277.91680.0000007.91687116.08428.16730.00008.1673108.16731418.413516.58088.41350.0000008.413529113.99835.58480.00005.5848105.58481616.517812.10276.51785.0000016.5178B32616.665415.147616.66540.00000010.1476Total trajectory length184.0919Figure12 shows the optimal trajectory in the space with the failure probability on certain correction points. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 12 and 13. The two figures clearly show the locations of correction points searched at each iteration. Like Figures 9–11, Figures 12–14 also use a series of straight lines to express the optimal trajectory.Figure 12
Optimal trajectory with the failure probability in the space.Figure 13
Optimal trajectory with the failure probability in the XOY plane.Figure 14
Optimal trajectory with the failure probability in the XOZ plane.It can be seen from Table5 that the two types (“0 and 1”) of correction points searched alternately appear like Table 3, which is in conformity with the principle that correction points are as few as possible. From the vertical and horizontal errors before and after being corrected, it is easy to verify that the optimal trajectory satisfies all the constraints. That is, before flight errors are corrected by vertical or horizontal correction points, they satisfy the constraint ((12) or (13)), respectively. Among correction points searched, only one point, whose serial number is 161, has the failure probability. That is reasonable because those points with the failure probability, not completely eliminating flight errors, should be maximally avoided.It can be shown from Figure12 that the optimal trajectory is not better near straight line AB when its three axes keep the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory (with the failure probability on certain correction points) whose length is 184.0919 km has a large increase in length, further revealing that for the practical problem, the correction points with the failure probability is so severe so that the shorter trajectories cannot avoid them.It is true that both the length and number of correction points belonging to the optimal trajectory are relatively large. In spite of this, the optimal trajectory satisfies all the constraints, and its length and number of correction points are smaller than other trajectories obtained by Monte Carlo sampling. Therefore, it is worth believing that the trajectory planning method based on Monte Carlo sampling is feasible and reasonable and provides a solution for trajectory planning with fault tolerance for UAVs.


### 5.4. Comparison from Different Methods

The computation time is a crucial factor that decides whether the method proposed in this paper can be used in practice. Besides, the computation time based on Monte Carlo sampling may sometimes be unacceptable, so it is necessary to evaluate the time cost to obtain an optimized trajectory. To study the computation time, we give the computer configurations for simulation as follows.Operating system: Windows 10 with 64 bitCPU: AMD Ryzen 7 4800H with Radeon GraphicsRAM: 16 GB with 3200 MHzMATLAB version: MATLAB R2018aFrom the results with the failure probability on certain correction points, it can be seen that using the definite rule cannot obtain a feasible trajectory. Hence, the computation time will be analyzed only for the results without the failure probability.Table6 shows under the different weights the computation time, which is the average time by repeatedly performing 20 times. By observing the correction points searched under different weights shown in Table 2, it can be seen from Table 6 that the time computation mainly depends on their number. That is, the more correction points searched are, the longer the computation time is. Table 6 also reveals that the computation time during each iteration is about 0.031 s (i.e., 0.404/13 where “13” is the number of points searched). However, regardless of which weight is selected, computation time is short enough, which indicates that the method based on the definite rule has high efficiency and will promise to realize trajectory planning in real time under the hardware.Table 6
Computation time under the different weights.
Weightω0.3∼0.50.6∼0.70.8∼34∼56∼2527∼2930∼30.531∼48Computation time (s)0.44630.46900.46930.49010.40420.55890.58860.5869The research in this paper comes from the mathematical modeling competition of Chinese graduate students (including doctoral students) [43]. There are two excellent documents selected in this competition [44]. Their serial numbers are F19102510081 and F19103360018, respectively. As described in the two documents, the adaptive and improved Dijkstra algorithm and ant colony algorithm are used in the first document, and the greedy algorithm and tabu search algorithm are used in the second document.It can be seen after referring to the two documents that the numbers (i.e., 12) of correction points in the optimal trajectory are the same. The computation time is 42.23 s in the first document, 78.87 s in the second document, and however only 0.404 s in this paper. Their computation time is 104 and 195 times that of this paper, respectively. It is noted that the total lengths of the optimal trajectory are very close, and they are 109.464, 109.324, and 111.779 km, respectively, but our computational efficiency is much higher than the two documents.Although our computers may differ, the configuration is generally little different. Moreover, we just use the laptop whose configurations have been given in this paper. Therefore, this work provides an efficient trajectory planning method, which is expected to realize real-time trajectory planning considering the error correction. This efficient method also lays a foundation for the realization of trajectory planning with error correction with failure tolerance based on Monte Carlo sampling.


## 5.1. Problem to Be Solved

Table1 shows the partial data of a flight space, which includes the coordinates of the starting point, ending point, and correction points. The serial numbers 0 and 326 are the starting point and ending point, respectively, and the rest are correction points. In the column containing the type of correction points, “0” and “1” represent the horizontal and vertical error correction points, respectively. In the last column, “0” and “1” represent the correction points with and without the failure probability, respectively; the probability of a successful correction is only σ when the flight error is corrected by those having the failure probability.Table 1
Partial data of a flight space [42].
Serial numbersCoordinatesType of correction pointsFailure labelX (m)Y (m)Z (m)Starting point A00.0050000.005000.00——Correction points176009.859788.119121.901022448.2071599.881877.1300327800.9315218.078345.4100450056.7391668.091465.0900525311.1551225.791330.3501………………32478846.7073662.653906.981032563645.5671435.226263.1301Ending point B326100000.0074860.555499.61——For this problem, the parameters mentioned in Section2 are given as follows: α1 = 20 m, α2 = 10 m, β1 = 15 m, β2 = 20 m, θ = 20 m, δ = 0.001 m, σ = 80%, ζ = 5 m, and r = 200 m.


## 5.2. Results and Discussion Based on the Definite Rule

After considering that certain correction points have the failure probability and assuming that all corrections fail by these points at each iteration, the results can be obtained by running the programs developed in this paper. However, the simulation results show that, no matter which weight we select, the serial numbers of the correction points searched are all [163, termination]. This is a special phenomenon, which presents that after serial number 163 is searched, solving the optimal correction point model cannot find any correction point at the next iteration. This causes the cease of the programs unexpectedly. In other words, if there exists the failure probability on certain correction points, no feasible trajectory can be found based on the definite method for the practical problem whose distributions on correction points are relatively rigorous. “It is for this reason that we further propose a trajectory planning method based on Monte Carlo Sampling,” and the simulation results are as shown in Section5.3.To illustrate the simulation process and results based on the definite rule, the restriction of the correction points for the problem will be loosened, and it is assumed that all correction points will not fail because most problems fall into this situation. In other words, in equations (14) and (15), τ will remain 0, and the corresponding results are as follows.


### 5.2.1. Optimized Trajectories under the Different Weights

Because correction points are finite and not uniformly distributed in the whole space, weightω from 0 to infinite can be discretely picked to obtain optimized trajectories as shown in Table 2.Table 2
Optimized trajectory under the different weights.
WeightωTotal number of the best correction pointsTrajectory length (km)Serial numbers of the points searched0∼0.2——[163, 114, 8, 309, 255, 148, 90, 156, 320, termination]0.3∼0.514137.7666[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 187, 61, 292, 326]0.6∼0.714136.0649[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 323, 61, 292, 326]0.8∼314135.0177[163, 114, 8, 309, 255, 148, 90, 156, 33, 289, 50, 323, 61, 292, 326]4∼514127.5342[163, 114, 8, 309, 121, 123, 49, 160, 138, 19, 50, 323, 61, 292, 326]6∼2512111.7787[163, 114, 8, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]27∼2915123.6882[163, 114, 234, 188, 222, 230, 225, 170, 123, 49, 160, 92, 93, 61, 292, 326]30∼30.515122.4127[163, 114, 234, 188, 238, 227, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]31∼4816120.0604[163, 114, 234, 188, 238, 227, 309, 305, 123, 45, 160, 92, 93, 61, 292, 326]>49No solution[163, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, …]If weightω is very small, the best correction points searched mainly depend on the distance BXxi⟶. To minimize the objective function (i.e., equation (11)), a correction point farthest from the starting point will be found. In other words, a UAV will fly a longer distance at each iteration, which may cause that the optimal correction point model cannot find any correction point because flight errors are so large that they cannot satisfy the constraints from equations (12) and (13) at next iteration. When weight ω belongs to [0, 0.2], the series number of the correction points searched just illustrates the point as shown in Table 2. It should be noted that the “termination” shown in Table 2 presents that after serial number 320 is searched, solving the optimal correction point model cannot find any correction point at the next iteration, which causes the cease of the programs unexpectedly.With the increase of weightω, the importance of BXxi⟶ shown in Figure 6 is gradually weakened. Under a larger weight ω, it is easy to understand that if there are three or more correction points (satisfying constraints (6) and (7)) nearby the straight line between the starting and ending points, the optimal correction point model may always search the same points so that running the programs will fall into an endless loop. When ω > 49, the series number of the correction points searched just illustrates the point as shown in Table 2.It can be seen from Table2 that, the optimal trajectory is obtained when weight ω belongs to [6, 25]; in this case, trajectory length is the shortest, and the correction points UAV passes are very few. If there is no failure of all the correction points, therefore, the total length of the optimal trajectory is 111.7787 km, and the total number of the best correction points searched is 12.


### 5.2.2. Optimal Trajectory

Table3 shows the flight data of the UAV along the optimal trajectory, which includes the correction pints, vertical and horizontal errors before being corrected, type of correction points, and flight distance from the starting point updated continually to the best correction point searched at each iteration, respectively.Table 3
Flight data of the UAV along the optimal trajectory.
Serial numbersVertical error before being correctedHorizontal error before being correctedType of correction pointsFlight distance of each segment (km)Starting point A000——Best correction points searched16313.2878913.2879013.287911418.62955.341615.346813.957319.2989013.957330919.52155.564215.564212111.252316.8166011.252312316.66335.411015.41104911.79817.2090011.798016018.32196.523916.5239925.787812.311705.78789315.2739.485219.4852619.834719.3198809.834729216.3896.554316.5543Ending point B3266.980513.5348—6.9805Total trajectory length111.7787Figure9 shows the optimal trajectory in the space. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 9 and 10. The two figures clearly show the locations of the best correction points searched at each iteration. It is noted that in reality, the planned trajectory is a curve; however, Figures 9–11, the purpose of which is to describe the locations of the best correction points and to illustrate the process of searching correction points at each iteration, hence use a series of straight lines to express the optimal trajectory. Actually, owing to the small turning radius, it is not apparent to express the optimal trajectory with a curve.Figure 9
Optimal trajectory without the failure probability in the space.Figure 10
Optimal trajectory without the failure probability in the XOY plane.Figure 11
Optimal trajectory without the failure probability in the XOZ plane.It can be seen from Table3 that the two types (“0 and 1”) of the best correction points searched appear alternately, and the flight distance of each segment appears with the long and short one. That is coincident with the principle that the best correction points are as few as possible. In other words, if a best correction point searched is the vertical one, the next must be the horizontal one in order to maximize the flight distance before being corrected.Figure9 clearly shows that the flight trajectory is near the straight line AB when its three axes have the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory whose total length is only 111.7787 km illustrates that it is quite short. As shown in Table 3, all the flight errors are limited to the allowed range of the vertical and horizontal correction before being corrected. So, few are the best correction points searched, additionally, that it is sure that the optimal trajectory of the UAV is reasonable.Furthermore, some correction points, which satisfy all constraints, can be randomly selected to generate many feasible trajectories; it can be verified that, regardless of the trajectory length or number of correction points, they both are greater than the optimal trajectory’s counterparts.From what has been analyzed and discussed in the three paragraphs above, it is worth believing that the results obtained by the trajectory planning method based on the definite rule are feasible and reasonable, and this work provides an effective solution for the trajectory planning of UAVs considering the error correction when there is no failure on the correction points.


## 5.2.1. Optimized Trajectories under the Different Weights

Because correction points are finite and not uniformly distributed in the whole space, weightω from 0 to infinite can be discretely picked to obtain optimized trajectories as shown in Table 2.Table 2
Optimized trajectory under the different weights.
WeightωTotal number of the best correction pointsTrajectory length (km)Serial numbers of the points searched0∼0.2——[163, 114, 8, 309, 255, 148, 90, 156, 320, termination]0.3∼0.514137.7666[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 187, 61, 292, 326]0.6∼0.714136.0649[163, 114, 8, 309, 255, 148, 90, 156, 13, 164, 50, 323, 61, 292, 326]0.8∼314135.0177[163, 114, 8, 309, 255, 148, 90, 156, 33, 289, 50, 323, 61, 292, 326]4∼514127.5342[163, 114, 8, 309, 121, 123, 49, 160, 138, 19, 50, 323, 61, 292, 326]6∼2512111.7787[163, 114, 8, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]27∼2915123.6882[163, 114, 234, 188, 222, 230, 225, 170, 123, 49, 160, 92, 93, 61, 292, 326]30∼30.515122.4127[163, 114, 234, 188, 238, 227, 309, 121, 123, 49, 160, 92, 93, 61, 292, 326]31∼4816120.0604[163, 114, 234, 188, 238, 227, 309, 305, 123, 45, 160, 92, 93, 61, 292, 326]>49No solution[163, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, 114, 234, 188, …]If weightω is very small, the best correction points searched mainly depend on the distance BXxi⟶. To minimize the objective function (i.e., equation (11)), a correction point farthest from the starting point will be found. In other words, a UAV will fly a longer distance at each iteration, which may cause that the optimal correction point model cannot find any correction point because flight errors are so large that they cannot satisfy the constraints from equations (12) and (13) at next iteration. When weight ω belongs to [0, 0.2], the series number of the correction points searched just illustrates the point as shown in Table 2. It should be noted that the “termination” shown in Table 2 presents that after serial number 320 is searched, solving the optimal correction point model cannot find any correction point at the next iteration, which causes the cease of the programs unexpectedly.With the increase of weightω, the importance of BXxi⟶ shown in Figure 6 is gradually weakened. Under a larger weight ω, it is easy to understand that if there are three or more correction points (satisfying constraints (6) and (7)) nearby the straight line between the starting and ending points, the optimal correction point model may always search the same points so that running the programs will fall into an endless loop. When ω > 49, the series number of the correction points searched just illustrates the point as shown in Table 2.It can be seen from Table2 that, the optimal trajectory is obtained when weight ω belongs to [6, 25]; in this case, trajectory length is the shortest, and the correction points UAV passes are very few. If there is no failure of all the correction points, therefore, the total length of the optimal trajectory is 111.7787 km, and the total number of the best correction points searched is 12.


## 5.2.2. Optimal Trajectory

Table3 shows the flight data of the UAV along the optimal trajectory, which includes the correction pints, vertical and horizontal errors before being corrected, type of correction points, and flight distance from the starting point updated continually to the best correction point searched at each iteration, respectively.Table 3
Flight data of the UAV along the optimal trajectory.
Serial numbersVertical error before being correctedHorizontal error before being correctedType of correction pointsFlight distance of each segment (km)Starting point A000——Best correction points searched16313.2878913.2879013.287911418.62955.341615.346813.957319.2989013.957330919.52155.564215.564212111.252316.8166011.252312316.66335.411015.41104911.79817.2090011.798016018.32196.523916.5239925.787812.311705.78789315.2739.485219.4852619.834719.3198809.834729216.3896.554316.5543Ending point B3266.980513.5348—6.9805Total trajectory length111.7787Figure9 shows the optimal trajectory in the space. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 9 and 10. The two figures clearly show the locations of the best correction points searched at each iteration. It is noted that in reality, the planned trajectory is a curve; however, Figures 9–11, the purpose of which is to describe the locations of the best correction points and to illustrate the process of searching correction points at each iteration, hence use a series of straight lines to express the optimal trajectory. Actually, owing to the small turning radius, it is not apparent to express the optimal trajectory with a curve.Figure 9
Optimal trajectory without the failure probability in the space.Figure 10
Optimal trajectory without the failure probability in the XOY plane.Figure 11
Optimal trajectory without the failure probability in the XOZ plane.It can be seen from Table3 that the two types (“0 and 1”) of the best correction points searched appear alternately, and the flight distance of each segment appears with the long and short one. That is coincident with the principle that the best correction points are as few as possible. In other words, if a best correction point searched is the vertical one, the next must be the horizontal one in order to maximize the flight distance before being corrected.Figure9 clearly shows that the flight trajectory is near the straight line AB when its three axes have the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory whose total length is only 111.7787 km illustrates that it is quite short. As shown in Table 3, all the flight errors are limited to the allowed range of the vertical and horizontal correction before being corrected. So, few are the best correction points searched, additionally, that it is sure that the optimal trajectory of the UAV is reasonable.Furthermore, some correction points, which satisfy all constraints, can be randomly selected to generate many feasible trajectories; it can be verified that, regardless of the trajectory length or number of correction points, they both are greater than the optimal trajectory’s counterparts.From what has been analyzed and discussed in the three paragraphs above, it is worth believing that the results obtained by the trajectory planning method based on the definite rule are feasible and reasonable, and this work provides an effective solution for the trajectory planning of UAVs considering the error correction when there is no failure on the correction points.


## 5.3. Results and Discussions Based on Monte Carlo Sampling

In this section, the trajectory planning method based on Monte Carlo sampling will be adopted, considering that there is a failure on certain correction points, and the error corrections all fail by these points. By repeatedly running programs whose flowchart is shown in Figure8, a series of feasible trajectories can be obtained based on Monte Carlo sampling. The optimal trajectory can be selected from them by considering the trajectory length and number of correction points a UAV passes.


### 5.3.1. Feasible Trajectories

Some feasible trajectories obtained by Monte Carlo sampling are shown in Table4; in fact, many others are not listed due to the large trajectory lengths.Table 4
Feasible trajectories obtained by Monte Carlo sampling.
Number of correction points searchedTrajectory length (km)Serial numbers of points searched (including the ending point 326)122184.0919[252, 266, 270, 89, 236, 133, 3, 55, 144, 69, 298, 117, 85, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]223191.0631[252, 266, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 9, 253, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]325212.61575[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 167, 249, 274, 211, 321, 16, 282, 141, 291, 161, 326]426221.16501[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 73, 207, 70, 211, 321, 16, 282, 141, 291, 84, 124, 326]527226.8188[169, 266, 270, 89, 236, 185, 3, 281, 197, 112, 103, 250, 243, 73, 249, 274, 16, 186, 233, 120, 202, 71, 141, 282, 84, 287, 99, 326]628223.2609[169, 266, 100, 146, 270, 89, 236, 133, 3, 281, 144, 55, 197, 112, 143, 250, 9, 253, 249, 274, 12, 219, 16, 101, 141, 291, 84, 124, 326]729229.7933[169, 266, 270, 89, 236, 133, 3, 281, 144, 55, 69, 298, 117, 85, 76, 6, 211, 44, 191, 160, 324, 19, 50, 177, 87, 246, 313, 210, 212, 326]830233.0981[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 243, 73, 39, 274, 51, 211, 321, 16, 71, 233, 101, 141, 291, 84, 59, 326]931248.6503[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 55, 144, 261, 298, 143, 250, 86, 119, 243, 167, 249, 274, 12, 219, 16, 71, 141, 291, 84, 287, 61, 326]1032259.0214[169, 157, 252, 266, 270, 89, 236, 185, 78, 281, 3, 14, 214, 179, 111, 144, 55, 69, 298, 143, 250, 243, 167, 249, 274, 12, 321, 16, 282, 84, 291, 161, 326]It can be seen from Table4 that the more the correction points searched are, the larger the trajectory length becomes. Apparently, “trajectory 1” is the optimal trajectory; in this case, trajectory length is the shortest, and the correction points the UAV passes are very few. Considering that there exists the failure probability on certain correction points, therefore, the optimal trajectory length is 184.0919 km and its total number of correction points searched is 22. When a UAV flies along the optimal trajectory, even if error corrections fail, it can still arrive at the destination successfully. “It means that the trajectory planning method based on Monte Carlo sampling has the ability of fault tolerance.”


### 5.3.2. Optimal Trajectory

Table5 shows the flight data of a UAV along the optimal trajectory with the failure probability on certain correction points. The flight data include the correction points, vertical and horizontal errors before and after being corrected, type of correction points, failure label, and flight distance from the starting point updated continually to the point searched at each iteration.Table 5
Flight data of a UAV along the optimal trajectory with the failure probability on certain correction points.
PointsSerial numbersVertical error before being correctedHorizontal error before being correctedVertical error after being correctedHorizontal error after being correctedType of correction pointsFailure labelFlight distance of each segment (km)A00000———Correction points searched25211.107611.107611.10760.00000011.107626619.10717.99950.00007.9995107.99952706.437014.43656.43700.0000006.43708914.08987.65280.00007.6528107.652823610.137917.790710.13790.00000010.137913312.77042.63250.00002.6325102.632538.838311.47088.83830.0000008.83835518.02079.18240.00009.1824109.18241442.716711.89922.71670.0000002.71676914.589611.872914.58960.00000011.872929819.67395.08430.00005.0843105.08431179.905514.98989.90550.0000009.90558515.54455.63900.00005.6390105.63907614.323519.962514.32350.00000014.3235619.54905.22550.00005.2255105.22551210.702715.928210.70270.00000010.702721918.58867.88590.00007.8859107.8859167.916815.80277.91680.0000007.91687116.08428.16730.00008.1673108.16731418.413516.58088.41350.0000008.413529113.99835.58480.00005.5848105.58481616.517812.10276.51785.0000016.5178B32616.665415.147616.66540.00000010.1476Total trajectory length184.0919Figure12 shows the optimal trajectory in the space with the failure probability on certain correction points. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 12 and 13. The two figures clearly show the locations of correction points searched at each iteration. Like Figures 9–11, Figures 12–14 also use a series of straight lines to express the optimal trajectory.Figure 12
Optimal trajectory with the failure probability in the space.Figure 13
Optimal trajectory with the failure probability in the XOY plane.Figure 14
Optimal trajectory with the failure probability in the XOZ plane.It can be seen from Table5 that the two types (“0 and 1”) of correction points searched alternately appear like Table 3, which is in conformity with the principle that correction points are as few as possible. From the vertical and horizontal errors before and after being corrected, it is easy to verify that the optimal trajectory satisfies all the constraints. That is, before flight errors are corrected by vertical or horizontal correction points, they satisfy the constraint ((12) or (13)), respectively. Among correction points searched, only one point, whose serial number is 161, has the failure probability. That is reasonable because those points with the failure probability, not completely eliminating flight errors, should be maximally avoided.It can be shown from Figure12 that the optimal trajectory is not better near straight line AB when its three axes keep the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory (with the failure probability on certain correction points) whose length is 184.0919 km has a large increase in length, further revealing that for the practical problem, the correction points with the failure probability is so severe so that the shorter trajectories cannot avoid them.It is true that both the length and number of correction points belonging to the optimal trajectory are relatively large. In spite of this, the optimal trajectory satisfies all the constraints, and its length and number of correction points are smaller than other trajectories obtained by Monte Carlo sampling. Therefore, it is worth believing that the trajectory planning method based on Monte Carlo sampling is feasible and reasonable and provides a solution for trajectory planning with fault tolerance for UAVs.


## 5.3.1. Feasible Trajectories

Some feasible trajectories obtained by Monte Carlo sampling are shown in Table4; in fact, many others are not listed due to the large trajectory lengths.Table 4
Feasible trajectories obtained by Monte Carlo sampling.
Number of correction points searchedTrajectory length (km)Serial numbers of points searched (including the ending point 326)122184.0919[252, 266, 270, 89, 236, 133, 3, 55, 144, 69, 298, 117, 85, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]223191.0631[252, 266, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 9, 253, 76, 6, 12, 219, 16, 71, 141, 291, 161, 326]325212.61575[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 167, 249, 274, 211, 321, 16, 282, 141, 291, 161, 326]426221.16501[169, 266, 270, 89, 236, 133, 147, 278, 144, 55, 197, 112, 143, 250, 86, 73, 207, 70, 211, 321, 16, 282, 141, 291, 84, 124, 326]527226.8188[169, 266, 270, 89, 236, 185, 3, 281, 197, 112, 103, 250, 243, 73, 249, 274, 16, 186, 233, 120, 202, 71, 141, 282, 84, 287, 99, 326]628223.2609[169, 266, 100, 146, 270, 89, 236, 133, 3, 281, 144, 55, 197, 112, 143, 250, 9, 253, 249, 274, 12, 219, 16, 101, 141, 291, 84, 124, 326]729229.7933[169, 266, 270, 89, 236, 133, 3, 281, 144, 55, 69, 298, 117, 85, 76, 6, 211, 44, 191, 160, 324, 19, 50, 177, 87, 246, 313, 210, 212, 326]830233.0981[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 281, 197, 112, 143, 250, 243, 73, 39, 274, 51, 211, 321, 16, 71, 233, 101, 141, 291, 84, 59, 326]931248.6503[252, 157, 169, 264, 322, 270, 89, 236, 133, 3, 55, 144, 261, 298, 143, 250, 86, 119, 243, 167, 249, 274, 12, 219, 16, 71, 141, 291, 84, 287, 61, 326]1032259.0214[169, 157, 252, 266, 270, 89, 236, 185, 78, 281, 3, 14, 214, 179, 111, 144, 55, 69, 298, 143, 250, 243, 167, 249, 274, 12, 321, 16, 282, 84, 291, 161, 326]It can be seen from Table4 that the more the correction points searched are, the larger the trajectory length becomes. Apparently, “trajectory 1” is the optimal trajectory; in this case, trajectory length is the shortest, and the correction points the UAV passes are very few. Considering that there exists the failure probability on certain correction points, therefore, the optimal trajectory length is 184.0919 km and its total number of correction points searched is 22. When a UAV flies along the optimal trajectory, even if error corrections fail, it can still arrive at the destination successfully. “It means that the trajectory planning method based on Monte Carlo sampling has the ability of fault tolerance.”


## 5.3.2. Optimal Trajectory

Table5 shows the flight data of a UAV along the optimal trajectory with the failure probability on certain correction points. The flight data include the correction points, vertical and horizontal errors before and after being corrected, type of correction points, failure label, and flight distance from the starting point updated continually to the point searched at each iteration.Table 5
Flight data of a UAV along the optimal trajectory with the failure probability on certain correction points.
PointsSerial numbersVertical error before being correctedHorizontal error before being correctedVertical error after being correctedHorizontal error after being correctedType of correction pointsFailure labelFlight distance of each segment (km)A00000———Correction points searched25211.107611.107611.10760.00000011.107626619.10717.99950.00007.9995107.99952706.437014.43656.43700.0000006.43708914.08987.65280.00007.6528107.652823610.137917.790710.13790.00000010.137913312.77042.63250.00002.6325102.632538.838311.47088.83830.0000008.83835518.02079.18240.00009.1824109.18241442.716711.89922.71670.0000002.71676914.589611.872914.58960.00000011.872929819.67395.08430.00005.0843105.08431179.905514.98989.90550.0000009.90558515.54455.63900.00005.6390105.63907614.323519.962514.32350.00000014.3235619.54905.22550.00005.2255105.22551210.702715.928210.70270.00000010.702721918.58867.88590.00007.8859107.8859167.916815.80277.91680.0000007.91687116.08428.16730.00008.1673108.16731418.413516.58088.41350.0000008.413529113.99835.58480.00005.5848105.58481616.517812.10276.51785.0000016.5178B32616.665415.147616.66540.00000010.1476Total trajectory length184.0919Figure12 shows the optimal trajectory in the space with the failure probability on certain correction points. The views in the two coordinate planes projected by the optimal trajectory are shown in Figures 12 and 13. The two figures clearly show the locations of correction points searched at each iteration. Like Figures 9–11, Figures 12–14 also use a series of straight lines to express the optimal trajectory.Figure 12
Optimal trajectory with the failure probability in the space.Figure 13
Optimal trajectory with the failure probability in the XOY plane.Figure 14
Optimal trajectory with the failure probability in the XOZ plane.It can be seen from Table5 that the two types (“0 and 1”) of correction points searched alternately appear like Table 3, which is in conformity with the principle that correction points are as few as possible. From the vertical and horizontal errors before and after being corrected, it is easy to verify that the optimal trajectory satisfies all the constraints. That is, before flight errors are corrected by vertical or horizontal correction points, they satisfy the constraint ((12) or (13)), respectively. Among correction points searched, only one point, whose serial number is 161, has the failure probability. That is reasonable because those points with the failure probability, not completely eliminating flight errors, should be maximally avoided.It can be shown from Figure12 that the optimal trajectory is not better near straight line AB when its three axes keep the same range [0, 100 km]. Compared to the length of straight line AB being 103.0451 km, the optimal trajectory (with the failure probability on certain correction points) whose length is 184.0919 km has a large increase in length, further revealing that for the practical problem, the correction points with the failure probability is so severe so that the shorter trajectories cannot avoid them.It is true that both the length and number of correction points belonging to the optimal trajectory are relatively large. In spite of this, the optimal trajectory satisfies all the constraints, and its length and number of correction points are smaller than other trajectories obtained by Monte Carlo sampling. Therefore, it is worth believing that the trajectory planning method based on Monte Carlo sampling is feasible and reasonable and provides a solution for trajectory planning with fault tolerance for UAVs.


## 5.4. Comparison from Different Methods

The computation time is a crucial factor that decides whether the method proposed in this paper can be used in practice. Besides, the computation time based on Monte Carlo sampling may sometimes be unacceptable, so it is necessary to evaluate the time cost to obtain an optimized trajectory. To study the computation time, we give the computer configurations for simulation as follows.Operating system: Windows 10 with 64 bitCPU: AMD Ryzen 7 4800H with Radeon GraphicsRAM: 16 GB with 3200 MHzMATLAB version: MATLAB R2018aFrom the results with the failure probability on certain correction points, it can be seen that using the definite rule cannot obtain a feasible trajectory. Hence, the computation time will be analyzed only for the results without the failure probability.Table6 shows under the different weights the computation time, which is the average time by repeatedly performing 20 times. By observing the correction points searched under different weights shown in Table 2, it can be seen from Table 6 that the time computation mainly depends on their number. That is, the more correction points searched are, the longer the computation time is. Table 6 also reveals that the computation time during each iteration is about 0.031 s (i.e., 0.404/13 where “13” is the number of points searched). However, regardless of which weight is selected, computation time is short enough, which indicates that the method based on the definite rule has high efficiency and will promise to realize trajectory planning in real time under the hardware.Table 6
Computation time under the different weights.
Weightω0.3∼0.50.6∼0.70.8∼34∼56∼2527∼2930∼30.531∼48Computation time (s)0.44630.46900.46930.49010.40420.55890.58860.5869The research in this paper comes from the mathematical modeling competition of Chinese graduate students (including doctoral students) [43]. There are two excellent documents selected in this competition [44]. Their serial numbers are F19102510081 and F19103360018, respectively. As described in the two documents, the adaptive and improved Dijkstra algorithm and ant colony algorithm are used in the first document, and the greedy algorithm and tabu search algorithm are used in the second document.It can be seen after referring to the two documents that the numbers (i.e., 12) of correction points in the optimal trajectory are the same. The computation time is 42.23 s in the first document, 78.87 s in the second document, and however only 0.404 s in this paper. Their computation time is 104 and 195 times that of this paper, respectively. It is noted that the total lengths of the optimal trajectory are very close, and they are 109.464, 109.324, and 111.779 km, respectively, but our computational efficiency is much higher than the two documents.Although our computers may differ, the configuration is generally little different. Moreover, we just use the laptop whose configurations have been given in this paper. Therefore, this work provides an efficient trajectory planning method, which is expected to realize real-time trajectory planning considering the error correction. This efficient method also lays a foundation for the realization of trajectory planning with error correction with failure tolerance based on Monte Carlo sampling.


## 6. Conclusions

To achieve trajectory planning of UAVs with fault tolerance, this paper respectively proposes two trajectory planning methods based on the definite rule and Monte Carlo sampling. By establishing five models, designing corresponding algorithms, and developing computational programs based on MATLAB, the results of trajectory planning of UAVs are obtained for a practical problem. According to the analyses and discussions on the results, the conclusions are drawn as follows.(1)
The results obtained by trajectory planning methods based on the definite rule and Monte Carlo sampling are feasible and reasonable, which achieve the planning objectives that the total trajectory length is as short as possible and that the correction points UAV passes are as few as possible.(2)
The trajectory planning method based on the definite rule proposed in this paper has high efficiency and will promise to realize trajectory planning in real time, suitable for those problems without the failure probability of correction points. However, the trajectory planning method of UAVs with fault tolerance (allowing error corrections fails) based on Monte Carlo sampling is more suitable for off-line trajectory planning due to its expensive time cost.(3)
This work provides not only an effective solution for the trajectory planning of UAVs considering error correction without the failure probability but also a feasible solution for the trajectory planning of UAVs which considers the error correction, failure probability on certain correction point after being corrected, and arriving at the destination successfully.



---

*Source: 6741227-2022-07-11.xml*