Next Article in Journal
On Implementing a Two-Step Interior Point Method for Solving Linear Programs
Next Article in Special Issue
An Improved Iterated Greedy Algorithm for Solving Collaborative Helicopter Rescue Routing Problem with Time Window and Limited Survival Time
Previous Article in Journal
SCMs: Systematic Conglomerated Models for Audio Cough Signal Classification
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Sequential Convex Programming for Nonlinear Optimal Control in UAV Trajectory Planning

College of Intelligent Systems, Science and Engineering, Harbin Engineering University, Harbin 150001, China
*
Authors to whom correspondence should be addressed.
Algorithms 2024, 17(7), 304; https://doi.org/10.3390/a17070304
Submission received: 23 May 2024 / Revised: 25 June 2024 / Accepted: 28 June 2024 / Published: 8 July 2024

Abstract

:
In this paper, an algorithm is proposed to solve the non-convex optimization problem using sequential convex programming. An approximation method was used to solve the collision avoidance constraint. An iterative approach was utilized to estimate the non-convex constraints, replacing them with their linear approximations. Through the simulation, we can see that this method allows for quadcopters to take off from a given initial position and fly to the desired final position within a specified flight time. It is guaranteed that the quadcopters will not collide with each other in different scenarios.

1. Introduction

Quadcopters are the subject of ongoing research due to their agility, mechanical simplicity, and robustness. They are capable of performing a wide range of tasks and have been utilized as platforms for investigating vision-based pose estimation [1], nonlinear control [2], and learning [3], with many other applications. Moreover, they are helpful tools for solving practical issues, such as surveillance [4] and inspection [5]. Trajectory/path planning has obtained more attention as the interest in quadcopters has increased. These machines may interact with each other and the environment. As described in [6,7,8], recent trajectory optimization methods utilize a simple kinematics model and exclusively consider the acceleration vector as the decision variable. Consequently, the resulting trajectories only include information about velocities and positions. These methods assume the presence of a reliable onboard trajectory tracking controller, which justifies using a simple kinematics model. In other words, the trajectory planning algorithm independently controls each motor on the quadrotor to achieve stability. The early works [9,10] introduced traditional control methods for quadrotors, specifically a nonlinear geometric tracking controller on SE (3). These controllers can be classified as PID (proportional–integral–derivative) since they directly control the error dynamics in the attitude and position between the quadrotor and its desired target trajectory. After introducing traditional control methods, trajectory planning techniques have emerged as a distinct optimization issue, generally ignoring the actual dynamics of the drone.
In robotics, the generation of trajectories that adhere to avoidance constraints has been significantly studied. As stated in [11], there are two distinct approaches: planning and reacting. The planned strategy involves generating feasible courses in advance, whereas the reactive approach usually relies on an online collision avoidance system to respond to dangerous situations, as mentioned in [12]. Additionally, the two methods can also be combined. This work focused on a planned methodology. The planning can be performed in a decentralized manner [13] or in a centralized manner [14], as described in this paper. This study was conducted on trajectory generation with avoidance constraints for numerous vehicles, including differential drive robots [15], underwater vehicles [16], and aircraft [17]. Pattern formation research [18,19] also relies on collision-free trajectory generation. Various convex optimization techniques are employed to find collision-free trajectories while considering specific performance criteria. A model predictive control problem is described in [20], which examined highly specialized quadcopter dynamics and subsequently derived the constraints. The generation of trajectories involves constructing the trajectory of the quadcopter based on its jerk and solving a convex optimization problem for each independent axis. A technique for generating trajectories that accounts for simple dynamics constraints was proposed in [6], which is independent of the type of vehicle.
Various trajectory generation strategies have been proposed, depending on the desired objective. Advanced techniques for generating trajectories for quadcopters that consider the agent constraints can be found in [17,21]. This work aims to minimize the total thrust needed to fly the trajectories while using avoidance constraints to synchronize the trajectories of several vehicles. We assume that the agent can follow the generated trajectories (using an appropriate controller) if they satisfy a feasibility check. Separating the lengthy feasibility check from the trajectory planning algorithm results in a simple and flexible approach for rapidly generating collision-free trajectories for many vehicles. A method for planning aircraft trajectories that considers avoidance constraints and utilizes simplified dynamics was proposed in [10]. Two recent reviews [22,23] provide diverse research on motion planning for Unmanned Aerial Vehicles (UAVs). The complexity of trajectory generation is exacerbated by the under-actuated and nonlinear dynamics of the quadcopter and the challenging task of modeling aerodynamics (see [24,25] for studies on aerodynamic effects and [3] for a learning-based compensation technique). In [26], collision-free trajectories are computed to guide a fleet of UAVs from their starting positions to their destination positions, ensuring that the trajectories maintain a minimum separation distance while minimizing the total thrust produced by the quadcopters. The rapid solutions are obtained using sequential convex programming, with times in the range of seconds. In [27], the paper presents a convex programming approach for a fuel-optimal powered landing problem in the presence of aerodynamic drag forces and new types of non-convex control constraints.
Based on the classical methods proposed in the above literature, the mature business software [28] solves the linear problem quickly and effectively. Other methods treat the optimal control problem as a general nonlinear programming problem to better describe the actual model. However, in practice, it is of great difficulty to obtain a solution to nonlinear programming because of its nature, unbounded operation time, and the requirement for an initial guess, which should be provided by humans. Furthermore, the solution has the chance to be a local optimal rather than a global. Due to the above reasons, treating the optimal control problem as a general nonlinear programming approach for solving is not appropriate for onboard applications. To solve nonlinear programming stably and quickly is the primary goal of this paper. Conversely, the problem is solved perfectly if the nonlinear programming can be transformed into convex programming without losing the optimal solution. Here, the local optimal and global optimal are identical, and the solution time is bounded, therefore an initial prediction is unnecessary. Unfortunately, the transformation cannot be identical, and the optimality or the equivalence must be removed. The solution to the approximated problem, which identifies a key point of the original problem, is a feasible alternative. So, we need to find a convex programming formulation or a series of convex programming formulations whose solutions correspond to key points of the original problem. Approximation methodology plays a key role in the process. Regarding the above issue, this paper proposes an algorithm that iteratively uses convex optimization to approximate non-convex optimization.
The trajectory planning problem we considered involves solving the trajectory tracking problem using optimization techniques. There are various existing categories of optimization algorithms, such as convex optimization, non-convex optimization, etc. The primary objective of this work was to design and implement a sequential convex programming algorithm (SCP) for quadrotor trajectory planning. Sequential convex programming can be categorized as a non-convex optimization algorithm. Such methods rely on constructing convex sub-problems via successive linearizations of the original (nonlinear) objective functions and constraints. The resulting sub-problems become convex optimization sub-problems, which are typically solved by existing convex programming solvers. Sequential convex programming (SCP) is a technique employed to plan trajectories for numerous agents in three-dimensional (3D) space [29]. The objective is to transition from an initial set of states, including each vehicle’s position, velocity, and acceleration, to a final set of states. This must be done while maintaining a minimum distance between the agents and satisfying additional trajectory constraints. The methodology is implemented on quadcopters (Figure 1) but can be easily adapted for various platforms by adjusting the constraints. This study was inspired by the research conducted on sequential convex programming (SCP) by Wang et al. [30].
In this paper, we explored the capability of sequential convex programming and applied it to the quadrotor dynamics model, including trajectory dynamics. An iterative approach was used to estimate the non-convex constraints, replacing them with their linear approximations, and showing the simulation results. The simulation demonstrated that the approach can generate collision-free trajectories in different scenarios. The results demonstrate that this method can be used in real-world applications.
This study is outlined in the following order: In Section 2, trajectory dynamics are constructed. In Section 3, the constraints are explained in detail, including formulation of the problem as a convex optimization problem and discussion on using the iterative approach for approximating non-convex constraints. In Section 4, the simulation results are discussed. Section 5 contains the conclusion and an analysis of this paper.

2. Dynamic Model

Based on the given set of initial and final positions for N drones, the objective was to ensure that these drones reached their intended positions within a specified time interval, T, while keeping a minimum distance of at least R from each other during their trajectories. An optimization problem was formulated to accomplish this goal, which involved several convex and affine constraints and a non-convex collision avoidance constraint. The solution to this optimization problem guaranteed the achievement of the specified objective.
The goal was to generate trajectories for N agents that do not collide, allowing a transition from an initial set of states to a final set of states (position, velocity, and acceleration) within a specified time T. The trajectories must adhere to many boundaries and avoidance constraints.
The quadcopter was modeled as a rigid body with six degrees of freedom: it can linearly translate along the inertial x 3 axes and rotate around these three inertial axes. The rotation is described by the frame attached to the body relative to the inertial frame, using the proper orthogonal matrix R. The control inputs to the system were taken as the total thrust produced, denoted as f, for simplicity, normalized by the vehicle mass, with units of acceleration. The body rates, expressed as ω = (ω1, ω2, ω3), are considered in the body-fixed frame. x 1 , x 2 , x 3 represents the position of the quadcopter in inertial space. An illustration of this model can be seen in Figure 1.
Euler angles describe the orientation of a rigid body in space. Each angle represents the rotation of the quadrotor in a certain axis. It is used to construct the rotation matrices (Rx1, Rx2, Rx3). These rotation matrices map the body frame onto the inertial frame. Multiplying the three rotation matrices in a x3–x2–x1 order will yield the final orientation of the object. The final rotation matrix can be expressed as follows:
R x 3 x 2 x 1 ϕ , θ , ψ = R x 3 ψ R x 2 θ R x 1 ϕ = c θ c ψ s ϕ s θ c ψ c ϕ s ψ c ϕ s θ c ψ + s ϕ s ψ c θ s ψ s ϕ s θ s ψ + c ϕ c ψ c ϕ s θ s ψ s ϕ c ψ s θ s ϕ c θ c ϕ c θ
where ϕ , θ , ψ are the quadrotor’s roll, pitch, and yaw angles, respectively, s represents sin, and c represents cos.
The rotational dynamics of a quadcopter can be described by Euler’s equations of motion, as follows:
I ω ˙ + ω × I ω = τ
where I is the inertia matrix, ω is the angular velocity vector, and τ represents the torque vector.
The differential equations governing the flight of the quadcopter are now those of a rigid body, as follows:
x ¨ = R e 3 f + g
R ˙ = R ω ×
With e 3 = 0 , 0 , 1 and ω × denoting the skew-symmetric matrix form of the vector cross product, such that:
ω × = 0 ω 3 ω 2 ω 3 0 ω 1 ω 2 ω 1 0
and g = 0 , 0 , g represents the acceleration due to gravity.

2.1. Trajectory Dynamics

Let P i t , V i t , and a i t represent the location, velocity, and acceleration, respectively, of agent ‘i’ at discrete times t. The trajectory of agent ‘i’ is determined by the following discretized dynamics equations:
P i t + 1 = P i t + Δ t V i t + Δ t 2 2 a i t
V i t + 1 = V i t + Δ t a i t
where Δ t represents the duration of each time step. The variable t represents the current time step, ranging from 0 to T. The variable i represents the agent number from 1 to N. T and N represent the final time step and the total number of agents, respectively. The position and the velocity of each agent can be expressed as affine functions of acceleration, as indicated by Equations (6) and (7). Equations (6) and (7) can be iteratively expanded to be written explicitly as a function of the accelerations of the vehicle at the previous time steps. Therefore, the velocity of quadcopter i at time t is determined as follows:
V i t = V i 0 + Δ t k = 1 t a i k
From Equations (6) to (8), the position of quadcopter i at time t is represented by the following:
P i t = P i 0 + Δ t   V i 0 + Δ t 2 2   k = 1 t a i k 2 t 2 k + 1  
where t ranges from 1 to T, and i ranges from 1 to N. The V i 0 is equal to 0, so Equations (8) and (9) can be simplified as follows:
V i t = Δ t k = 1 t a i k
P i t = P i 0 + Δ t 2 2 k = 1 t a i k 2 t 2 k + 1
Equations (10) and (11) can be utilized to create an optimization problem involving the optimization parameter X . This parameter represents the acceleration of each agent at every time step, i.e., X 3 N T .

2.2. Objective Function

The control inputs for a quadcopter consist of the total thrust produced and the angular velocities. This study examined a normalized thrust force expressed in units of acceleration for the sake of simplicity.
a = R 0 0 f + 0 0 g
where R represents the rotation matrix that accurately describes the rotation between the inertial and body-fixed frames.
The quadcopter is described by six degrees of freedom. The translational position x 1 , x 2 , x 3 is measured in the inertial coordinate system, as shown in Figure 1. The quadcopter attitude b is defined by the rotation matrix R . The rotation matrix is defined such that when multiplying a vector v in the quadcopter coordinate b system with it, the same vector, described in the inertial coordinate system i, is obtained.
v i = R b i v i
The input thrust f can be calculated by applying the Euclidean norm ǁ ǁ to Equation (12), as follows:
f = a g
where g is the acceleration due to gravity and is equal to (0, 0, −g) in its vector form. To determine the jerk, it is necessary to express R ˙ in terms of the known parameters. The expression of R ˙ can be obtained by using the rotation matrix, R , and the skew-symmetric matrix of body rates, as follows:
R ˙ = R 0 ω 3 ω 2 ω 3 0 ω 1 ω 2 ω 1 0
By applying Equations (13) and (15) and taking the derivative of the acceleration, the jerk can be calculated as follows:
j = R ˙ 0 0 f + R 0 0 f ˙
j = R 0 ω 3 ω 2 ω 3 0 ω 1 ω 2 ω 1 0 0 0 f + 0 0 f ˙
Therefore, the cumulative thrust exerted by all agents at every time step can be selected as the cost function, as follows:
f 0 = t = 1 T i = 1 N a i t + g 2 2
Equation (16) can be formulated as a quadratic equation, with G representing the gravity vector of each agent at every time step, i.e., G 3 N T , as follows:
f 0 = X T X + 2 G T X + G T G
Equation (18) is expressed in a quadratic form. Hence, employing this cost function alongside appropriate affine equality/inequality constraints is feasible. Therefore, this issue can quickly be resolved using existing quadratic programming solvers.
Once the acceleration of the agents is determined at each time step, the total thrust force may be calculated using Equation (14). By rearranging Equation (16), we may derive the following equation:
ω 2 ω 1 0 = 1 f 1 0 0 0 1 0 0 0 0 R T j
The absence of the third component of body rates, ω3, can be explained by observing that a rotation about the ω3 axis does not impact the translational acceleration. The system inputs for a trajectory defined in its jerk are given using Equations (14) and (20), with one degree of freedom left in ω3. This could be fully specified if the total attitude of the quadcopter, including rotation around the thrust axis, is also known. For simplicity, we will assume that the value of ω3 is zero and that this rotation is insignificant.
In essence, the output of the optimization problem, which is the acceleration of the agents at each time step, may be used to calculate suitable inputs, body rates, and thrust quickly.

3. Problem Formulation

3.1. Convex Constraints

In this section, we describe the trajectory planning problem as an optimization problem. Therefore, we define an objective function that we aim to minimize, a variable that may be optimized freely, simple dynamics equations that describe the evolution of the trajectories, and constraints that must be met regarding inequalities and equalities. The proposed method can generate trajectories for vehicles of any arbitrary nature.
Imposing several constraints while minimizing accelerations is necessary to achieve the desired objective. The initial and final positions of all the vehicles are known. Additionally, to accommodate the physical constraints of an agent, velocity, acceleration, and jerk may be restricted. Finally, we need to ensure that the generated trajectories maintain a minimum distance of R between each other at all time steps, thereby confirming collision avoidance.
Next, we will establish and express these constraints in matrix form. To simplify the representation, we analyzed one vehicle at a time and then combined these constraints to form the constraints for all the agents.

3.2. Equality Constraints

To begin the algorithm, it is necessary to select the initial velocity and the position parameters. Furthermore, they must decelerate when they near their individual goals until they reach complete rest. Therefore, their final velocities and accelerations must also be zero. Equations (6) and (7) can be utilized to enforce the constraints on initial velocity and position. Equality constraints can enforce final acceleration, velocity, and position constraints. Considering that the constraints apply to each axis for final acceleration, velocity, and position, there must be nine equality constraints in total for every agent. Therefore, it is necessary to have a total of 9N boundary constraints in the system. These constraints can be expressed as equality constraints in the following manner:
  • Final velocity constraint: v i t = 0 3 × 1
  • Final acceleration constraint: a i t = 0 3 × 1
Final position constraint: the final positions in the form of xf. For a single agent ‘i’, we have the following:
p i t = X f , i
The matrix representation of these constraints can be written as follows:
A e q , f p , i χ i = b e q , f p , i
After establishing the equality constraint matrices for a single agent, our next step was to combine these constraints into a unified matrix. The equality constraints for a single vehicle, ‘i’, can be stacked on top of each other, resulting in a matrix representation of the equality constraints as follows:
A e q , i χ i = b e q , i
We have now created the matrices representing the equality constraints for a single vehicle. We generalized this concept to encompass the matrices that express all the equality constraints for all the agents.
A e q X = b e q
where A e q 9 N × 3 N T + 1 and b e q 9 N .

3.3. Inequality Constraints

Position constraint: The motion range in the workspace is physically limited by maximum position (pmax) and minimum position (pmin). Thus, given a single agent ‘i’ at time step k (>1), we can write the position constraints as follows:
p min p i t p max
Velocity constraint: To ensure safe flight in a confined space and dynamic feasibility, the velocity of the UAVs must be limited by minimum velocity (vmin) and maximum velocity (vmax). Therefore, for a single agent ‘i’ at time step k (>1), we can express velocity constraints as follows:
v min v i t v max
Acceleration constraint: To ensure a smooth and dynamically feasible trajectory, it is necessary to limit the acceleration of all the vehicles at all times within the range defined by minimum acceleration (amin) and maximum acceleration (amax). Thus, for a specific vehicle ‘i’ at time step k (>1), the acceleration constraints can be expressed as follows:
a min a i t a max
Jerk constraint: Limiting a specific vehicle’s jerk (ji) is essential. Jerk can be defined as the first derivative of the acceleration to maintain smooth acceleration in the quadcopter [2,3]. Jerk must be bounded for smooth and feasible trajectories. For simplicity, the jerk can be constrained in each axis individually.
j min j i t j max
The equation below expresses the jerk of any agent at any time jit in discrete time.
j i t = a i t a i t 1 Δ t
By combining all these constraints, we derive the following inequality:
A f e a s i b i l i t y X b f e a s i b i l i t y
where A f e a s i b i l i t y 6 N 2 T + 1 × 3 N T + 1 and b f e a s i b i l i t y 6 N 2 T + 1 .
The inequality constraints established so far constitute a convex feasible set. However, the collision avoidance constraint is formulated as a non-convex constraint. We employed the SCP algorithm to obtain an affine approximation of this constraint. Subsequently, we constructed the total inequality matrices for the true optimization variable X .

3.4. Collision Avoidance Constraint

To prevent the drones from colliding, it is necessary to maintain a minimum distance between them at each time step.
p i t p j t 2 R   i , j   i j , t
where R is the collision threshold distance. Equation (6) is linearized as follows:
p i t , q p j t , q + p i t , q p j t , q T p i t , q p j t , q p i t p j t R
where the superscript pit,q denotes the position vector of agent ‘i’ at time step ‘k’ in the previous outer-loop iteration q.
The Collision Avoidance Constraint can be derived as follows:
A c o l l i s i o n X b c o l l i s i o n
By combining Equations (28) and (31), we obtain the following inequality:
A X b
An optimization problem can be solved in two phases. In the first step, we disregard the collision avoidance constraint and construct our quadratic programming (QP) problem, relying on the position, velocity, acceleration, and jerk constraints. Consequently, a single agent’s total inequality coefficient and residual matrices are obtained by vertically stacking all the coefficients and residual matrices. The result Χ ^ is subsequently employed to commence the second phase. Starting from the second iteration, we can apply the collision avoidance constraints because we now have a previous solution that serves as a reference to simplify the collision avoidance constraint.
(1)
We found an initial value by using the feasibility constraints.
min X f 0
which is subjected to the following:
A e q X = b e q
A f e a s i b i l i t y X b f e a s i b i l i t y
which will determine Χ ^ .
(2)
We utilized the estimated value Χ ^ to create new inequality matrices A and b.
(3)
min X f 0
which is subjected to A e q X = b e q
A X b
During the second step, an iterative approach is used, where the most recent estimation of X ( x ^ ) is applied in each iteration, assuming that f t h r e s h o l d is a specified threshold. Step two is repeated if f0 − (f0)prev > = f t h r e s h o l d .
Because all the constraints are either convex or convexified, quadratic programming tools can be used to solve the optimization problem. The following is the pseudocode of this algorithm. A 1 , f e a s i b i l i t y represents the value in the first iteration of A f e a s i b i l i t y and b 1 , f e a s i b i l i t y represents the value in the first iteration of b f e a s i b i l i t y .

3.5. Sequential Convex Programming

Sequential Convex Programming (SCP) is a local optimization method for non-convex problems that leverage convex optimization.
We considered the non-convex problem as follows:
min f 0 x
s . t . f i x 0 , i = 1 , 2 , , m
h j x = 0 , j = 1 , 2 , , p
with x n , f 0 x and f i x non-convex and h j non-affine.
In the k t h iteration of SCP, we maintained the estimate x k of the solution, and the convex trust region τ k n . We constructed convex approximations f ¯ 0 and f ¯ i of f 0 and f i , respectively, and affine approximations h ¯ j of h j over the trust region τ k . The approximate convex problem obtained is written as follows:
min f ¯ 0 x
f ¯ i x 0 , i = 1 , 2 , , m
h ¯ j x = 0 , j = 1 , 2 , , p
The solution to this approximate problem is x k + 1 .
In summary, the specific algorithm for trajectory optimization using SCP is depicted as Algorithm 1:
Algorithm 1: Trajectory Optimization Using SCP
Require: N,R, T, f t h r e s h o l d , p max , p min , v max , v min , a max , a min , j max , j min
1: X 0 Q P A e q , b e q , A 1 , f e a s i b i l i t y , b 1 , f e a s i b i l i t y
2: f p f X 0
3: n i 0
4: while n i m i do
5: X Q P A e q , b e q , A f e a s i b i l i t y , b f e a s i b i l i t y
6: f c f X
7: if f c f q < f t h r e s h o l d then
8: break
9: end if
10: f p f c
11: n i n i + 1
12: end while

4. Experiment

This section presents four case studies and demonstrates the performance of the suggested strategy through simulation data. The algorithm was executed using MATLAB R2022b on a PC with an Intel i9-13900HX CPU operating at a clock cycle of 2.2 GHz and 16.0 GB RAM. Specifically, we used MATLAB’s quadratic programming function, quadprop. The case study examines two, four, six, and eight quadcopters operating within a three-dimensional space of 5 × 5 × 3, respectively. Initially, the simulation was executed without imposing any distance constraints. Figure 2a, Figure 3a, Figure 4a and Figure 5 illustrate that the algorithm produced the optimal trajectories, i.e., the shortest distances between the initial and final locations. Due to the absence of constraints on collision, the trajectories intersect at some points.
For the first iteration of the optimization loop of SCP, we did not consider the collision avoidance constraint and we formulated our QP problem using only the position, velocity, acceleration, and jerk constraints; therefore, the total inequality coefficient matrices, as well as the residual matrices, for a single drone were obtained by stacking all the coefficient matrices and residual matrices vertically.
In the first iteration, we did not have any previous solution to linearize our collision constraint; therefore, we did not consider the collision constraint in the first iteration. However, from the second iteration onwards, we enforced the collision avoidance constraints, since we had a previous solution around which we could linearize the collision avoidance constraint.
Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10 show the simulation results of two, four, six, and eight quadcopters in the absence of collision avoidance constraints and collision-free.
The parameters used for the case study were chosen as follows:
(1)
Number of agents = 2
Simulation length = 30 s
Discretization step-size, dt = 0.2 s
Cost Threshold, f t h r e s h o l d = 0.05
Run time for the first step of the algorithm = 2.45 s
Total run time = 35 s
Figure 2. (a) Trajectories for two quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for two quadcopters.
Figure 2. (a) Trajectories for two quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for two quadcopters.
Algorithms 17 00304 g002
(2)
Number of agents = 4
Simulation length = 30 s
Discretization step-size, dt = 0.2 s
Cost Threshold, f t h r e s h o l d = 0.05
Run time for the first step of the algorithm = 4.86 s
Total run time = 67 s
Figure 3. (a) Trajectories for four quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for four quadcopters.
Figure 3. (a) Trajectories for four quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for four quadcopters.
Algorithms 17 00304 g003
(3)
Number of agents = 6
Simulation length = 30 s
Discretization step-size, dt = 0.2 s
Cost Threshold, f t h r e s h o l d = 0.05
Run time for the first step of the algorithm = 6.75 s
Total run time = 94 s
Figure 4. (a) Trajectories for six quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for six quadcopters.
Figure 4. (a) Trajectories for six quadcopters in the absence of collision avoidance constraints. (b) Collision-free trajectories for six quadcopters.
Algorithms 17 00304 g004
(4)
Number of agents = 8
Simulation length = 30 s
Discretization step-size, dt = 0.2 s
Number of inequality constraints without distance constraints = 7200
Total number of inequality constraints with distance
constraints = 11,400
Cost Threshold, f t h r e s h o l d = 0.05
Number of iterations for convergence = 2
Run time for the first step of the algorithm = 10.85 s
Total run time = 168 s
Figure 5 and Figure 6 show the variation of position and velocity as time of eight quadcopters in the absence of collision avoidance constraints. Figure 9 and Figure 10 show the variation of position and velocity as time of eight quadcopters with collision-free constraints.
Figure 5. Trajectories for eight quadcopters in the absence of collision avoidance constraints.
Figure 5. Trajectories for eight quadcopters in the absence of collision avoidance constraints.
Algorithms 17 00304 g005
Figure 6. Position vs. time of eight quadcopters in the absence of collision avoidance constraints.
Figure 6. Position vs. time of eight quadcopters in the absence of collision avoidance constraints.
Algorithms 17 00304 g006
Figure 7. Velocity vs. time of eight quadcopters in the absence of collision avoidance constraints.
Figure 7. Velocity vs. time of eight quadcopters in the absence of collision avoidance constraints.
Algorithms 17 00304 g007
During the execution of the second stage of the algorithm, the trajectories were automatically adjusted by the algorithm to prevent collision. Figure 2b, Figure 3b, Figure 4b and Figure 8a,b demonstrate that deviation in the trajectories was minimal while ensuring that the trajectories did not intersect or collide. In these specific case studies, the minimum distance was defined as 1 m (dmin = 1). Depending on the quadcopter’s size and the size of the environment, this parameter was chosen accordingly. The selection of this parameter can be adjusted based on the dimensions of the quadcopter and the surrounding environment. Figure 11 shows the objective value plotted against the number of iterations. It can be seen from Figure 10 that the objective function value rose suddenly after the first iteration and continued to decrease from the second iteration onwards. This behavior was expected since, for the first iteration, the optimization problem was solved without the collision constraint. Following the first iteration, the collision constraint was taken into account, which reduced the objective value and resulted in a collision-free optimal trajectory.
Figure 8. (a) Collision-free trajectories for eight quadcopters. (b) Collision-free trajectories for eight quadcopters.
Figure 8. (a) Collision-free trajectories for eight quadcopters. (b) Collision-free trajectories for eight quadcopters.
Algorithms 17 00304 g008
Figure 9. Position vs. time of eight quadcopters with collision-free constraints.
Figure 9. Position vs. time of eight quadcopters with collision-free constraints.
Algorithms 17 00304 g009
Figure 10. Velocity vs. time of eight quadcopters with collision-free constraints.
Figure 10. Velocity vs. time of eight quadcopters with collision-free constraints.
Algorithms 17 00304 g010
Figure 11. Objective value plotted against the number of iterations.
Figure 11. Objective value plotted against the number of iterations.
Algorithms 17 00304 g011

5. Conclusions

In this study, an algorithm is presented that enables the generation of collision-free trajectories for a multi-agent quadcopter fleet. The non-convex optimization problem was iteratively solved using sequential convex programming that approximates non-convex constraints by using convex constraints. Through the simulation, we can see that this method allows for quadcopters to take off from a given initial position and fly to the desired final position within a specified flight time. It is guaranteed that the quadcopters will not collide with each other. Since the control effort was minimized, these trajectories are the shortest possible trajectories under the collision avoidance constraints. The simulation results show that, in different scenarios, this method can generate collision-free trajectories from the beginning point to the target point. Furthermore, the simulation results demonstrate the computational effort, and the performance indicates that this algorithm can be used in real-life applications with well-tuned hyper-parameters. These collision-free trajectories can be generated online.

Author Contributions

Conceptualization, Q.Z. and Y.L.; methodology, Y.L. and Q.Z.; software, Y.L. and A.E.; validation, Y.L. and A.E.; formal analysis, Y.L.; investigation, Y.L.; resources, Q.Z.; writing-original draft preparation, Y.L.; writing-review and editing, Y.L. and A.E.; visualization, Y.L.; supervision, Q.Z.; project administration, Q.Z.; funding acquisition, Q.Z. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

Data are contained within the article.

Conflicts of Interest

The authors declare no conflicts of interest.

References

  1. Weiss, S.; Achtelik, M.; Chli, M.; Siegwart, R. Versatile distributed pose estimation and sensor self-calibration for an autonomous MAV. In Proceedings of the IEEE International Conference on Robotics and Automation, Saint Paul, MN, USA, 14–18 May 2012; pp. 31–38. [Google Scholar]
  2. Cowling, I.D.; Yakimenko, O.A.; Whidborne, J.F.; Cooke, A.K. A prototype of an autonomous controller for a quadrotor UAV. In Proceedings of the European Control Conference, Kos, Greece, 2–5 July 2007; pp. 1–8. [Google Scholar]
  3. Lupashin, S.; Schollig, A.; Sherback, M.; D’Andrea, R. A simple learning strategy for high-speed quadrocopter multi-flips. In Proceedings of the IEEE International Conference on Robotics and Automation, Anchorage, AK, USA, 3–7 May 2010; pp. 1642–1648. [Google Scholar]
  4. Alexis, K.; Nikolakopoulos, G.; Tzes, A.; Dritsas, L. Coordination of helicopter UAVs for aerial forest-fire surveillance. In Applications of Intelligent Control to Engineering Systems; Intelligent Systems, Control, and Automation: Science and Engineering; Springer: Dordrecht, The Netherlands, 2009; Volume 39. [Google Scholar] [CrossRef]
  5. Serrano, N.E. Autonomous Quadrotor Unmanned Aerial Vehicle for Culvert Inspection. Ph.D. Dissertation, Massachusetts Institute of Technology, Cambridge, MA, USA, 2011. [Google Scholar]
  6. Augugliaro, F.; Schoellig, A.P.; D’Andrea, R. Generation of collision-free trajectories for a quadrocopter fleet: A sequential convex programming approach. In Proceedings of the 2012 IEEE/RSJ International Conference on Intelligent Robots and Systems, Vilamoura-Algarve, Portugal, 7–12 October 2012; IEEE: Piscataway, NJ, USA, 2012; pp. 1917–1922. [Google Scholar]
  7. Luis, C.E.; Schoellig, A.P. Trajectory generation for multiagent point-to-point transitions via distributed model predictive control. IEEE Robot. Autom. Lett. 2019, 4, 375–382. [Google Scholar] [CrossRef]
  8. Chen, Y.; Cutler, M.; How, J.P. Decoupled multiagent path planning via incremental sequential convex programming. In Proceedings of the 2015 IEEE International Conference on Robotics and Automation (ICRA), Seattle, WA, USA, 26–30 May 2015; pp. 5954–5961. [Google Scholar]
  9. Lee, T.; Leok, M.; McClamroch, N.H. Geometric tracking control of a quadrotor UAV on SE(3). In Proceedings of the 49th IEEE Conference on Decision and Control (CDC), Atlanta, GA, USA, 15–17 December 2010; pp. 5420–5425. [Google Scholar]
  10. Mellinger, D.; Kumar, V. Minimum snap trajectory generation and control for quadrotors. In Proceedings of the 2011 IEEE International Conference on Robotics and Automation, Shanghai, China, 9–13 May 2011; pp. 2520–2525. [Google Scholar]
  11. Siegwart, R.; Nourbakhsh, I.R.; Scaramuzza, D. Introduction to Autonomous Mobile Robots; MIT Press: Cambridge, MA, USA, 2011. [Google Scholar]
  12. Gillula, J.H.; Hoffmann, G.M.; Vitus, M.P.; Tomlin, C.J. Applications of hybrid reachability analysis to robotic aerial vehicles. Int. J. Robot. Res. 2011, 30, 335–354. [Google Scholar] [CrossRef]
  13. Kuwata, Y.; Richards, A.; Schouwenaars, T.; How, J.P. Distributed Robust Receding Horizon Control for Multivehicle Guidance. IEEE Trans. Control Syst. Technol. 2007, 15, 627–641. [Google Scholar] [CrossRef]
  14. Raghunathan, A.U.; Gopal, V.; Subramanian, D.; Biegler, L.T.; Samad, T. Dynamic optimization strategies for three-dimensional conflict resolution of multiple aircraft. J. Guid. Control. Dyn. 2004, 27, 586–594. [Google Scholar] [CrossRef]
  15. Snape, J.; van den Berg, J.; Guy, S.J.; Manocha, D. Smooth and collision-free navigation for multiple robots under differential-drive constraints. In Proceedings of the EEE/RSJ International Conference on Intelligent Robots and Systems, Taipei, Taiwan, 18–22 October 2010; pp. 4584–4589. [Google Scholar]
  16. Yuan, H.; Qu, Z. Optimal real-time collision-free motion planning for autonomous underwater vehicles in a 3D underwater space. IET Control Theory Appl. 2009, 3, 712. [Google Scholar] [CrossRef]
  17. Hehn, M.; D’Andrea, R. Quadrocopter Trajectory Generation and Control. In Proceedings of the IFAC World Congress, Milan, Italy, 28 August–2 September 2011; pp. 1485–1491. [Google Scholar]
  18. Schouwenaars, T.; How, J.; Feron, E. Decentralized cooperative trajectory planning of multiple aircraft with hard safety guarantees. In Proceedings of the AIAA Guidance, Navigation and Control Conference, Providence, RI, USA, 16–19 August 2004. [Google Scholar]
  19. Varghese, B.; McKee, G.; Beji, L.; Otmane, S.; Abichou, A. Towards a Unifying Framework for Pattern Transformation in Swarm Systems. In AIP Conference Proceedings; AIP: College Park, MD, USA, 2009; Volume 1107, pp. 65–70. [Google Scholar]
  20. Alonso-Mora, J.; Breitenmoser, A.; Rufli, M.; Siegwart, R.; Beardsley, P. Multi-robot system for artistic pattern formation. In Proceedings of the IEEE International Conference on Robotics and Automation, Shanghai, China, 9–13 May 2011; pp. 4512–4517. [Google Scholar]
  21. Mueller, M.W.; D’Andrea, R. A model predictive controller for quadrocopter state interception. In Proceedings of the 2013 European Control Conference (ECC), Zurich, Switzerland, 17–19 July 2013; IEEE: Piscataway, NJ, USA, 2013; pp. 1383–1389. [Google Scholar]
  22. Richards, A.; How, J. Aircraft trajectory planning with collision avoidance using mixed integer linear programming. In Proceedings of the American Control Conference, Anchorage, AK, USA, 8–10 May 2002; Volume 3, pp. 1936–1941. [Google Scholar]
  23. Goerzen, C.; Kong, Z.; Mettler, B. A Survey of Motion Planning Algorithms from the Perspective of Autonomous UAV Guidance. J. Intell. Robot. Syst. 2009, 57, 65–100. [Google Scholar] [CrossRef]
  24. Kendoul, F. Survey of advances in guidance, navigation, and control of unmanned rotorcraft systems. J. Field Robot. 2012, 29, 315–378. [Google Scholar] [CrossRef]
  25. Hoffmann, G.M.; Huang, H.; Waslander, S.L.; Tomlin, C.J. Quadrotor helicopter flight dynamics and control: Theory and experiment. In Proceedings of the AIAA Guidance, Navigation, and Control Conference, Hilton Head, SC, USA, 20–23 August 2007; pp. 1–20. [Google Scholar]
  26. Martin, P.; Salaun, E. The true role of accelerometer feedback in quadrotor control. In Proceedings of the IEEE International Conference on Robotics and Automation, Anchorage, AK, USA, 3–7 May 2010; pp. 1623–1629. [Google Scholar]
  27. Szmuk, M.; Acikmese, B.; Berning, A.W. Successive convexification for fuel-optimal powered landing with aerodynamic drag and non-convex constraints. In Proceedings of the AIAA Guidance, Navigation, & Control Conference, Kissimmee, FL, USA, 5–9 January 2015; p. 0378. [Google Scholar]
  28. Jennings, L.; Teo, K.; Goh, C. Miser 3.2 Optimal Control Software: Theory and User Manual; Department of Mathematics, University of Western Australia: Perth, Australia, 1997. [Google Scholar]
  29. Boyd, S. Sequential Convex Programming; Lecture Notes; Stanford University: Stanford, CA, USA, 2008. [Google Scholar]
  30. Wang, Y.; Akuiyibo, E.; Boyd, S. Applications of Convex Optimization in Control; Stanford University: Stanford, CA, USA, 2011. [Google Scholar]
Figure 1. Dynamic model of a quadcopter.
Figure 1. Dynamic model of a quadcopter.
Algorithms 17 00304 g001
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Li, Y.; Zhu, Q.; Elahi, A. Sequential Convex Programming for Nonlinear Optimal Control in UAV Trajectory Planning. Algorithms 2024, 17, 304. https://doi.org/10.3390/a17070304

AMA Style

Li Y, Zhu Q, Elahi A. Sequential Convex Programming for Nonlinear Optimal Control in UAV Trajectory Planning. Algorithms. 2024; 17(7):304. https://doi.org/10.3390/a17070304

Chicago/Turabian Style

Li, Yong, Qidan Zhu, and Ahsan Elahi. 2024. "Sequential Convex Programming for Nonlinear Optimal Control in UAV Trajectory Planning" Algorithms 17, no. 7: 304. https://doi.org/10.3390/a17070304

APA Style

Li, Y., Zhu, Q., & Elahi, A. (2024). Sequential Convex Programming for Nonlinear Optimal Control in UAV Trajectory Planning. Algorithms, 17(7), 304. https://doi.org/10.3390/a17070304

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop